diff --git a/Liv_Lb0NeutronLines/hlt1.py b/Liv_Lb0NeutronLines/hlt1.py
new file mode 100644
index 0000000000000000000000000000000000000000..8765323b88da3b8075162c50e6d36dfba2c8defe
--- /dev/null
+++ b/Liv_Lb0NeutronLines/hlt1.py
@@ -0,0 +1,14 @@
+"""
+Stolen verbatim from Hlt/Hlt1Conf/options/hlt1_pp_expected_24_without_UT.py from https://gitlab.cern.ch/lhcb/Moore/-/merge_requests/3620/diffs
+ - Added @cburr's AgeLimit to reduce memory consumption
+"""
+from Moore import Options
+from Moore.production import hlt1
+
+
+def alg_config(options: Options):
+    config = hlt1(options, "--sequence=hlt1_pp_matching_no_ut_1000KHz",
+                  "--flagging")
+    # Will be the default in releases including https://gitlab.cern.ch/lhcb/LHCb/-/merge_requests/4641
+    config["Gaudi::IODataManager/IODataManager"].AgeLimit = 0
+    return config
diff --git a/Liv_Lb0NeutronLines/hlt2.py b/Liv_Lb0NeutronLines/hlt2.py
new file mode 100644
index 0000000000000000000000000000000000000000..dd304ca2d9ac84ade0e83d6cd3f677887e6f78b4
--- /dev/null
+++ b/Liv_Lb0NeutronLines/hlt2.py
@@ -0,0 +1,36 @@
+"""
+Stolen and adapted from Hlt/Hlt2Conf/options/hlt2_pp_expected_24_without_UT.py from v55r11.
+ - Added @cburr's AgeLimit to reduce memory consumption
+"""
+from Moore import Options
+from RecoConf.reconstruction_objects import reconstruction
+
+from Hlt2Conf.lines.mc.mc_lines import MC_HLT2_PASSTHROUGH_LINE_NAME
+from Hlt2Conf.lines.b_to_open_charm import all_lines
+
+
+def _make_regex():
+    """Use some QEE lines for the sake of example."""
+    # Here we 'brute-force' this regex via a regex of each individual line name we want to consider.
+    # A neater format is possible via knowing how the lines are named but this feels almost less robust.
+    lnames = [
+        linename for line_dict in [all_lines] for linename in line_dict.keys()
+    ]
+    lnames.append(MC_HLT2_PASSTHROUGH_LINE_NAME)
+    return f"({'|'.join(lnames)})"
+
+
+def alg_config(options: Options):
+    with reconstruction.bind(from_file=False):
+        from Moore.production import hlt2
+        # Rather than processing via Moore:run_moore, using Moore.production:hlt2 emulates MC centralised processing.
+        hlt2_extra_args = [
+            '--flagging',  # prescales+postscales removed + passthrough_line added
+            '--settings=hlt2_pp_2024',  # Decides which settings to create the streams from
+            f'--lines-regex={_make_regex()}',  # Filters the settings to only run the lines you wish
+        ]
+        # instantiates public_tools and run_moore
+        config = hlt2(options, *hlt2_extra_args)
+    # Will be the default in releases including https://gitlab.cern.ch/lhcb/LHCb/-/merge_requests/4641
+    config["Gaudi::IODataManager/IODataManager"].AgeLimit = 0
+    return config
diff --git a/Liv_Lb0NeutronLines/info.yaml b/Liv_Lb0NeutronLines/info.yaml
new file mode 100644
index 0000000000000000000000000000000000000000..4435551a2c1fbd740ef4847f3a8ffbddd04a61f5
--- /dev/null
+++ b/Liv_Lb0NeutronLines/info.yaml
@@ -0,0 +1,220 @@
+defaults:
+  inform:
+    - juan.baptista.leite@cern.ch
+    - shuqi.sheng@cern.ch
+  wg: B2OC
+
+###################
+## MC HLT2:
+###################
+
+# 2024
+{%- set polarities = [
+    'Down', 'Up'
+]%}
+
+{%- set mc_datasets = [
+ ('Lb0ToLcpLcmN0', 'Sim10c' , '7.6', 'DIGI', '15196810'),
+ ('BpToLcpLcmKp', 'Sim10c' , '7.6', 'DIGI', '12197003'),
+ ('Lb0ToPbarPN0','Sim10c','7.6','DIGI','15102831'),
+]%}
+
+
+# {%- for p in polarities %}
+# {%- for sigtype, simul, nu, inptype, evttype in mc_datasets %}
+
+# {{sigtype}}_Mag{{p}}_MC_HLT1:
+#   application: "Moore/v55r11p2@x86_64_v2-el9-gcc13+detdesc-opt"
+#   input: 
+#     bk_query: /MC/Dev/Beam6800GeV-expected-2024-Mag{{p}}-Nu{{nu}}-25ns-Pythia8/{{simul}}/{{evttype}}/{{inptype}}
+#     dq_flags:
+#       - OK
+#   output: MC_HLT1.DST
+#   options: 
+#     entrypoint: Liv_Lb0NeutronLines.hlt1:alg_config
+#     extra_options:
+#       input_raw_format: 0.5
+#       conddb_tag: sim-20231017-vc-md100
+#       dddb_tag: dddb-20231017
+#       input_type: ROOT
+#       output_type: ROOT
+#       simulation: True
+#       data_type: "Upgrade"
+#       scheduler_legacy_mode: False
+#       compression:
+#         algorithm: ZSTD
+#         level: 1
+#         max_buffer_size: 1048576
+
+# {{sigtype}}_Mag{{p}}_MC_HLT2:
+#   application: "Moore/v55r11p2@x86_64_v2-el9-gcc13+detdesc-opt"
+#   input:
+#     job_name: {{sigtype}}_Mag{{p}}_MC_HLT1
+#   output: MC_HLT2_full.DST
+#   options:
+#     entrypoint: Liv_Lb0NeutronLines.hlt2:alg_config
+#     extra_options:
+#       conddb_tag: sim-20231017-vc-md100
+#       dddb_tag: dddb-20231017
+#       input_type: "ROOT"
+#       output_type: "ROOT"
+#       simulation: True
+#       input_process: "Hlt2"
+#       persistreco_version: 0.0
+#       data_type: "Upgrade"
+#       evt_max: -1
+#       output_manifest_file: "HLT2.tck.json"
+#       scheduler_legacy_mode: False
+#       compression:
+#         algorithm: ZSTD
+#         level: 1
+#         max_buffer_size: 1048576
+        
+# {{sigtype}}_Mag{{p}}_MC_Tuple:
+#   application: "DaVinci/v64r8@x86_64_v2-el9-clang16+detdesc-opt"
+#   input:
+#     job_name: {{sigtype}}_Mag{{p}}_MC_HLT2
+#   output: MC_{{sigtype}}.ROOT
+#   options:
+#     entrypoint: Liv_Lb0NeutronLines.main_Liv_{{sigtype}}:main
+#     extra_options:
+#       input_raw_format: 0.5
+#       input_type: ROOT
+#       simulation: True
+#       data_type: Upgrade
+#       conddb_tag: sim-20231017-vc-md100
+#       dddb_tag: dddb-20231017
+#       input_process: "Hlt2"
+#       conditions_version: master
+#       lumi: false
+#       output_type: ROOT
+#       persistreco_version: 0.0
+#       write_decoding_keys_to_git: True
+#       input_manifest_file: "HLT2.tck.json"
+#       input_stream: default
+
+
+# {%- endfor %}
+# {%- endfor %}
+
+
+###################
+## Data HLT2:
+###################
+
+{%- for p in polarities %}
+{%- for sigtype, simul, nu, inptype, evttype in mc_datasets %}
+{{sigtype}}_Mag{{p}}_Data_Hlt2_Spc2_Tuple:
+  application: "DaVinci/v64r8@x86_64_v2-el9-clang16-opt"
+  input:
+    bk_query: /LHCb/Collision24/Beam6800GeV-VeloClosed-Mag{{p}}/Real Data/Sprucing24c2/94000000/B2OC.DST
+    dq_flags:
+      - UNCHECKED
+      - OK
+    keep_running: True
+  output: RealData_{{sigtype}}.ROOT
+  options:
+    entrypoint: Liv_Lb0NeutronLines.main_Liv_{{sigtype}}:main
+    extra_options:
+      input_raw_format: 0.5
+      input_type: ROOT
+      simulation: False
+      data_type: Upgrade
+      geometry_version: run3/trunk
+      conditions_version: master
+      input_process: "TurboPass"
+      input_stream: "b2oc"
+{%- endfor %}
+{%- endfor %}
+
+
+# I don't have HLT2 24c1 for Lb0 channels here
+{%- set mc_datasets = [
+ ('BpToLcpLcmKp', '' , '', '', '')
+]%}
+
+{%- for p in polarities %}
+{%- for sigtype, simul, nu, inptype, evttype in mc_datasets %}
+{{sigtype}}_Mag{{p}}_Data_Hlt2_ExclUT_Spc1_Tuple:
+  application: "DaVinci/v64r8@x86_64_v2-el9-clang16-opt"
+  input:
+    bk_query: /LHCb/Collision24/Beam6800GeV-VeloClosed-Mag{{p}}-Excl-UT/Real Data/Sprucing24c1/94000000/B2OC.DST
+    dq_flags:
+      - UNCHECKED
+      - OK
+    keep_running: True
+  output: RealData_{{sigtype}}.ROOT
+  options:
+    entrypoint: Liv_Lb0NeutronLines.main_Liv_{{sigtype}}:main
+    extra_options:
+      input_raw_format: 0.5
+      input_type: ROOT
+      simulation: False
+      data_type: Upgrade
+      geometry_version: run3/trunk
+      conditions_version: master
+      input_process: "TurboPass"
+      input_stream: "b2oc"
+{%- endfor %}
+{%- endfor %}
+
+###################
+## Data Spruce:
+###################
+
+{%- set mc_datasets = [
+ ('Lb0ToLcpLcmN0', '' , '', '', ''),
+ ('BpToLcpLcmKp', '' , '', '', ''),
+ ('Lb0ToPbarPN0','' , '', '', ''),
+]%}
+
+
+{%- for p in polarities %} 
+{%- for sigtype, simul, nu, inptype, evttype in mc_datasets %}
+{{sigtype}}_Mag{{p}}_Data_Spruce_Spc2_Tuple:
+  application: "DaVinci/v64r8@x86_64_v2-el9-clang16-opt"
+  input:
+    bk_query: /LHCb/Collision24/Beam6800GeV-VeloClosed-Mag{{p}}/Real Data/Sprucing24c2/90000000/B2OC.DST
+    dq_flags:
+      - UNCHECKED
+      - OK
+    keep_running: True
+  output: RealData_{{sigtype}}.ROOT
+  options:
+    entrypoint: Liv_Lb0NeutronLines.main_Liv_{{sigtype}}_Spruce:main
+    extra_options:
+      input_raw_format: 0.5
+      input_type: ROOT
+      simulation: False
+      data_type: Upgrade
+      geometry_version: run3/trunk
+      conditions_version: master
+      input_process: Spruce
+      input_stream: "b2oc"
+{%- endfor %}
+{%- endfor %}
+
+{%- for p in polarities %}
+{%- for sigtype, simul, nu, inptype, evttype in mc_datasets %}
+{{sigtype}}_Mag{{p}}_Data_Spruce_ExclUT_Spc1_Tuple:
+  application: "DaVinci/v64r8@x86_64_v2-el9-clang16-opt"
+  input:
+    bk_query: /LHCb/Collision24/Beam6800GeV-VeloClosed-Mag{{p}}-Excl-UT/Real Data/Sprucing24c1/90000000/B2OC.DST
+    dq_flags:
+      - UNCHECKED
+      - OK
+    keep_running: True
+  output: RealData_{{sigtype}}.ROOT
+  options:
+    entrypoint: Liv_Lb0NeutronLines.main_Liv_{{sigtype}}_Spruce:main
+    extra_options:
+      input_raw_format: 0.5
+      input_type: ROOT
+      simulation: False
+      data_type: Upgrade
+      geometry_version: run3/trunk
+      conditions_version: master
+      input_process: Spruce
+      input_stream: b2oc
+{%- endfor %}
+{%- endfor %}
diff --git a/Liv_Lb0NeutronLines/main_Liv_BpToLcpLcmKp.py b/Liv_Lb0NeutronLines/main_Liv_BpToLcpLcmKp.py
new file mode 100644
index 0000000000000000000000000000000000000000..e0c4c999ad6cb47ecc99a2ee65b0fa1ea558b00b
--- /dev/null
+++ b/Liv_Lb0NeutronLines/main_Liv_BpToLcpLcmKp.py
@@ -0,0 +1,32 @@
+###############################################################################
+# (c) Copyright 2021-2022 CERN for the benefit of the LHCb Collaboration      #
+#                                                                             #
+# This software is distributed under the terms of the GNU General Public      #
+# Licence version 3 (GPL Version 3), copied verbatim in the file "COPYING".   #
+#                                                                             #
+# In applying this licence, CERN does not waive the privileges and immunities #
+# granted to it by virtue of its status as an Intergovernmental Organization  #
+# or submit itself to any jurisdiction.                                       #
+###############################################################################
+
+from .options.Lb_NeutronLines_HLT2 import *
+
+from DaVinci import Options, make_config
+from PyConf.reading import get_pvs, get_rec_summary
+
+def main(options: Options):
+        
+    # get ODIN and DecReports location
+    pvs = get_pvs()
+    rec_summary = get_rec_summary()
+
+    """ CORREGGERE TUPLES """
+
+    tuples = {    
+        "BpToLcpLcmKp" : maketuple_BpToLcpLcmKp(options, pvs, rec_summary)
+    }
+        
+    config = make_config(options, tuples)
+
+    return config
+
diff --git a/Liv_Lb0NeutronLines/main_Liv_BpToLcpLcmKp_Spruce.py b/Liv_Lb0NeutronLines/main_Liv_BpToLcpLcmKp_Spruce.py
new file mode 100644
index 0000000000000000000000000000000000000000..e58d39a826dc9278b1f8886a8d33c843a7d84e58
--- /dev/null
+++ b/Liv_Lb0NeutronLines/main_Liv_BpToLcpLcmKp_Spruce.py
@@ -0,0 +1,32 @@
+###############################################################################
+# (c) Copyright 2021-2022 CERN for the benefit of the LHCb Collaboration      #
+#                                                                             #
+# This software is distributed under the terms of the GNU General Public      #
+# Licence version 3 (GPL Version 3), copied verbatim in the file "COPYING".   #
+#                                                                             #
+# In applying this licence, CERN does not waive the privileges and immunities #
+# granted to it by virtue of its status as an Intergovernmental Organization  #
+# or submit itself to any jurisdiction.                                       #
+###############################################################################
+
+from .options.Lb_NeutronLines_Spruce import *
+
+from DaVinci import Options, make_config
+from PyConf.reading import get_pvs, get_rec_summary
+
+def main(options: Options):
+        
+    # get ODIN and DecReports location
+    pvs = get_pvs()
+    rec_summary = get_rec_summary()
+
+    """ CORREGGERE TUPLES """
+
+    tuples = {    
+        "BpToLcpLcmKp" : maketuple_BpToLcpLcmKp(options, pvs, rec_summary)
+    }
+        
+    config = make_config(options, tuples)
+
+    return config
+
diff --git a/Liv_Lb0NeutronLines/main_Liv_Lb0ToLcpLcmN0.py b/Liv_Lb0NeutronLines/main_Liv_Lb0ToLcpLcmN0.py
new file mode 100644
index 0000000000000000000000000000000000000000..9cf6e161d9bf7345be4d572b882e0fc7fe336072
--- /dev/null
+++ b/Liv_Lb0NeutronLines/main_Liv_Lb0ToLcpLcmN0.py
@@ -0,0 +1,32 @@
+###############################################################################
+# (c) Copyright 2021-2022 CERN for the benefit of the LHCb Collaboration      #
+#                                                                             #
+# This software is distributed under the terms of the GNU General Public      #
+# Licence version 3 (GPL Version 3), copied verbatim in the file "COPYING".   #
+#                                                                             #
+# In applying this licence, CERN does not waive the privileges and immunities #
+# granted to it by virtue of its status as an Intergovernmental Organization  #
+# or submit itself to any jurisdiction.                                       #
+###############################################################################
+
+from .options.Lb_NeutronLines_HLT2 import *
+
+from DaVinci import Options, make_config
+from PyConf.reading import get_pvs, get_rec_summary
+
+def main(options: Options):
+        
+    # get ODIN and DecReports location
+    pvs = get_pvs()
+    rec_summary = get_rec_summary()
+
+    """ CORREGGERE TUPLES """
+
+    tuples = {    
+        "Lb0ToLcpLcmN0" : maketuple_Lb0ToLcpLcnN0(options, pvs, rec_summary)
+    }
+        
+    config = make_config(options, tuples)
+
+    return config
+
diff --git a/Liv_Lb0NeutronLines/main_Liv_Lb0ToLcpLcmN0_Spruce.py b/Liv_Lb0NeutronLines/main_Liv_Lb0ToLcpLcmN0_Spruce.py
new file mode 100644
index 0000000000000000000000000000000000000000..8fa00a3470eb384299c6a5e42c1b1b4c5fbded7c
--- /dev/null
+++ b/Liv_Lb0NeutronLines/main_Liv_Lb0ToLcpLcmN0_Spruce.py
@@ -0,0 +1,32 @@
+###############################################################################
+# (c) Copyright 2021-2022 CERN for the benefit of the LHCb Collaboration      #
+#                                                                             #
+# This software is distributed under the terms of the GNU General Public      #
+# Licence version 3 (GPL Version 3), copied verbatim in the file "COPYING".   #
+#                                                                             #
+# In applying this licence, CERN does not waive the privileges and immunities #
+# granted to it by virtue of its status as an Intergovernmental Organization  #
+# or submit itself to any jurisdiction.                                       #
+###############################################################################
+
+from .options.Lb_NeutronLines_Spruce import *
+
+from DaVinci import Options, make_config
+from PyConf.reading import get_pvs, get_rec_summary
+
+def main(options: Options):
+        
+    # get ODIN and DecReports location
+    pvs = get_pvs()
+    rec_summary = get_rec_summary()
+
+    """ CORREGGERE TUPLES """
+
+    tuples = {    
+        "Lb0ToLcpLcmN0" : maketuple_Lb0ToLcpLcnN0(options, pvs, rec_summary)
+    }
+        
+    config = make_config(options, tuples)
+
+    return config
+
diff --git a/Liv_Lb0NeutronLines/main_Liv_Lb0ToPbarPN0.py b/Liv_Lb0NeutronLines/main_Liv_Lb0ToPbarPN0.py
new file mode 100644
index 0000000000000000000000000000000000000000..65a044e77dcbdf928d495928ba7effb83584b56c
--- /dev/null
+++ b/Liv_Lb0NeutronLines/main_Liv_Lb0ToPbarPN0.py
@@ -0,0 +1,32 @@
+###############################################################################
+# (c) Copyright 2021-2022 CERN for the benefit of the LHCb Collaboration      #
+#                                                                             #
+# This software is distributed under the terms of the GNU General Public      #
+# Licence version 3 (GPL Version 3), copied verbatim in the file "COPYING".   #
+#                                                                             #
+# In applying this licence, CERN does not waive the privileges and immunities #
+# granted to it by virtue of its status as an Intergovernmental Organization  #
+# or submit itself to any jurisdiction.                                       #
+###############################################################################
+
+from .options.Lb_NeutronLines_HLT2 import *
+
+from DaVinci import Options, make_config
+from PyConf.reading import get_pvs, get_rec_summary
+
+def main(options: Options):
+        
+    # get ODIN and DecReports location
+    pvs = get_pvs()
+    rec_summary = get_rec_summary()
+
+    """ CORREGGERE TUPLES """
+
+    tuples = {    
+        "Lb0ToPpPmN0" : maketuple_Lb0ToPpPmN0(options, pvs, rec_summary)
+    }
+        
+    config = make_config(options, tuples)
+
+    return config
+
diff --git a/Liv_Lb0NeutronLines/main_Liv_Lb0ToPbarPN0_Spruce.py b/Liv_Lb0NeutronLines/main_Liv_Lb0ToPbarPN0_Spruce.py
new file mode 100644
index 0000000000000000000000000000000000000000..c7acf2a6d8d065d9c61036949806c94d61b9e5db
--- /dev/null
+++ b/Liv_Lb0NeutronLines/main_Liv_Lb0ToPbarPN0_Spruce.py
@@ -0,0 +1,32 @@
+###############################################################################
+# (c) Copyright 2021-2022 CERN for the benefit of the LHCb Collaboration      #
+#                                                                             #
+# This software is distributed under the terms of the GNU General Public      #
+# Licence version 3 (GPL Version 3), copied verbatim in the file "COPYING".   #
+#                                                                             #
+# In applying this licence, CERN does not waive the privileges and immunities #
+# granted to it by virtue of its status as an Intergovernmental Organization  #
+# or submit itself to any jurisdiction.                                       #
+###############################################################################
+
+from .options.Lb_NeutronLines_Spruce import *
+
+from DaVinci import Options, make_config
+from PyConf.reading import get_pvs, get_rec_summary
+
+def main(options: Options):
+        
+    # get ODIN and DecReports location
+    pvs = get_pvs()
+    rec_summary = get_rec_summary()
+
+    """ CORREGGERE TUPLES """
+
+    tuples = {    
+        "Lb0ToPpPmN0" : maketuple_Lb0ToPpPmN0(options, pvs, rec_summary)
+    }
+        
+    config = make_config(options, tuples)
+
+    return config
+
diff --git a/Liv_Lb0NeutronLines/options/Lb_NeutronLines_HLT2.py b/Liv_Lb0NeutronLines/options/Lb_NeutronLines_HLT2.py
new file mode 100644
index 0000000000000000000000000000000000000000..94b99244cf852c99dee34731039f8470eff8e41b
--- /dev/null
+++ b/Liv_Lb0NeutronLines/options/Lb_NeutronLines_HLT2.py
@@ -0,0 +1,216 @@
+from .tupling import *
+
+import Functors as F
+from Functors.math import log
+from DaVinci import Options, make_config
+from DaVinci.algorithms import create_lines_filter
+from PyConf.reading import get_particles
+from FunTuple import FunctorCollection
+from PyConf.reading import get_particles, get_pvs, get_extended_pvs
+from PyConf.Algorithms import ParticleUnbiasedPVAdder
+import FunTuple.functorcollections as FC
+from FunTuple import FunTuple_Particles as Funtuple
+from DecayTreeFitter import DecayTreeFitter
+
+
+def maketuple_Lb0ToLcpLcnN0(options, pvs, rec_summary):
+
+    name = "Lb0ToLcpLcmN0"
+    turbo_line = "Hlt2B2OC_Lb0ToLcpLcmN0"
+    input_data = get_particles(f"/Event/HLT2/{turbo_line}/Particles")
+
+    myfilter = create_lines_filter(
+        f"LineFilter_{turbo_line}_{{hash}}",
+        [turbo_line],
+    )
+
+    branches = {
+        "Lb0":
+        "[Lambda_b0 -> (Lambda_c+ -> p+ K- pi+) (Lambda_c~- -> p~- K+ pi-)]CC",
+        "Lcp":
+        "[Lambda_b0 -> ^(Lambda_c+ -> p+ K- pi+) (Lambda_c~- -> p~- K+ pi-)]CC",
+        "Lcm":
+        "[Lambda_b0 -> (Lambda_c+ -> p+ K- pi+) ^(Lambda_c~- -> p~- K+ pi-)]CC",
+        "pp":
+        "[Lambda_b0 -> (Lambda_c+ -> ^p+ K- pi+) (Lambda_c~- -> p~- K+ pi-)]CC",
+        "pm":
+        "[Lambda_b0 -> (Lambda_c+ -> p+ K- pi+) (Lambda_c~- -> ^p~- K+ pi-)]CC",
+        "pip":
+        "[Lambda_b0 -> (Lambda_c+ -> p+ K- ^pi+) (Lambda_c~- -> p~- K+ pi-)]CC",
+        "pim":
+        "[Lambda_b0 -> (Lambda_c+ -> p+ K- pi+) (Lambda_c~- -> p~- K+ ^pi-)]CC",
+        "Kp":
+        "[Lambda_b0 -> (Lambda_c+ -> p+ K- pi+) (Lambda_c~- -> p~- ^K+ pi-)]CC",
+        "Km":
+        "[Lambda_b0 -> (Lambda_c+ -> p+ ^K- pi+) (Lambda_c~- -> p~- K+ pi-)]CC",
+    }
+
+    DTF_MASS = DecayTreeFitter(
+        name='DTF_{hash}',
+        input_particles=input_data,
+        mass_constraints=["Lambda_c+"],
+    )
+
+    DTF_MASS_OWNPV = DecayTreeFitter(
+        name='DTF_{hash}',
+        input_particles=input_data,
+        mass_constraints=["Lambda_c+"],
+        constrain_to_ownpv=True,
+    )
+
+    basic_variables = make_basic_variables(
+        options, pvs, input_data, add_tistos=False)
+
+    composite_variables = make_composite_variables(
+        options, pvs, input_data, add_tistos=True)
+
+    composite_variables_3body = make_composite_variables_3body(
+        options, pvs, input_data, add_tistos=False)
+
+    DTF_TwoBody = make_composite_dtf_variables(options, pvs, input_data, add_truth=True, DTF=DTF_MASS, pv_constraint=False, mass_constraint=True, particle_name="_Lc", three_body=False) +\
+               make_composite_dtf_variables(options, pvs, input_data, add_truth=True, DTF=DTF_MASS_OWNPV, pv_constraint=True, three_body=False, mass_constraint=True,particle_name="_Lc") #same as BEST
+
+    variables = {
+        "Lb0": composite_variables + DTF_TwoBody,
+        "Lcp": composite_variables_3body,
+        "Lcm": composite_variables_3body,
+        "pp": basic_variables,
+        "pm": basic_variables,
+        "pip": basic_variables,
+        "pim": basic_variables,
+        "Kp": basic_variables,
+        "Km": basic_variables,
+    }
+
+    mytuple = Funtuple(
+        name=name,
+        tuple_name="DecayTree",
+        fields=branches,
+        variables=variables,
+        event_variables=make_hlt2_event_variables(options, pvs, rec_summary),
+        inputs=input_data,
+        store_multiple_cand_info=True)
+
+    return [myfilter, mytuple]
+
+
+def maketuple_BpToLcpLcmKp(options, pvs, rec_summary):
+
+    name = "BuToLcpLcmK"
+    turbo_line = "Hlt2B2OC_BuToLcpLcmK_LcpToPKPi"
+    input_data = get_particles(f"/Event/HLT2/{turbo_line}/Particles")
+
+    myfilter = create_lines_filter(
+        f"LineFilter_{turbo_line}_{{hash}}",
+        [turbo_line],
+    )
+
+    branches = {
+        "Bp":
+        "[B+ -> (Lambda_c+ -> p+ K- pi+) (Lambda_c~- -> p~- K+ pi-) K+]CC",
+        "Lcp":
+        "[B+ -> ^(Lambda_c+ -> p+ K- pi+) (Lambda_c~- -> p~- K+ pi-) K+]CC",
+        "Lcm":
+        "[B+ -> (Lambda_c+ -> p+ K- pi+) ^(Lambda_c~- -> p~- K+ pi-) K+]CC",
+        "pp":
+        "[B+ -> (Lambda_c+ -> ^p+ K- pi+) (Lambda_c~- -> p~- K+ pi-) K+]CC",
+        "pm":
+        "[B+ -> (Lambda_c+ -> p+ K- pi+) (Lambda_c~- -> ^p~- K+ pi-) K+]CC",
+        "pip":
+        "[B+ -> (Lambda_c+ -> p+ K- ^pi+) (Lambda_c~- -> p~- K+ pi-) K+]CC",
+        "pim":
+        "[B+ -> (Lambda_c+ -> p+ K- pi+) (Lambda_c~- -> p~- K+ ^pi-) K+]CC",
+        "Kp":
+        "[B+ -> (Lambda_c+ -> p+ K- pi+) (Lambda_c~- -> p~- ^K+ pi-) K+]CC",
+        "Km":
+        "[B+ -> (Lambda_c+ -> p+ ^K- pi+) (Lambda_c~- -> p~- K+ pi-) K+]CC",
+        "Kp_ext":
+        "[B+ -> (Lambda_c+ -> p+ K- pi+) (Lambda_c~- -> p~- K+ pi-) ^K+]CC",
+    }
+
+    DTF_MASS = DecayTreeFitter(
+        name='DTF_{hash}',
+        input_particles=input_data,
+        mass_constraints=["Lambda_c+"],
+    )
+
+    DTF_MASS_OWNPV = DecayTreeFitter(
+        name='DTF_{hash}',
+        input_particles=input_data,
+        mass_constraints=["Lambda_c+"],
+        constrain_to_ownpv=True,
+    )
+
+    basic_variables = make_basic_variables(options, pvs, input_data)
+
+    composite_variables_3body = make_composite_variables_3body(
+        options, pvs, input_data)
+
+    DTFVars = make_composite_dtf_variables(options, pvs, input_data, add_truth=True, DTF=DTF_MASS, pv_constraint=False, three_body=True, mass_constraint=True,particle_name="_Lc") +\
+        make_composite_dtf_variables(options, pvs, input_data, add_truth=True, DTF=DTF_MASS_OWNPV, pv_constraint=True, three_body=True, mass_constraint=True,particle_name="_Lc") #same as BEST
+
+    MVA_variables = make_MVA_variables()
+
+    variables = {
+        "Bp": composite_variables_3body + DTFVars + MVA_variables,
+        "Lcp": composite_variables_3body,
+        "Lcm": composite_variables_3body,
+        "pp": basic_variables,
+        "pm": basic_variables,
+        "pip": basic_variables,
+        "pim": basic_variables,
+        "Kp": basic_variables,
+        "Km": basic_variables,
+        "Kp_ext": basic_variables,
+    }
+
+    mytuple = Funtuple(
+        name=name,
+        tuple_name="DecayTree",
+        fields=branches,
+        variables=variables,
+        event_variables=make_hlt2_event_variables(options, pvs, rec_summary),
+        inputs=input_data,
+        store_multiple_cand_info=True)
+
+    return [myfilter, mytuple]
+
+
+def maketuple_Lb0ToPpPmN0(options, pvs, rec_summary):
+
+    name = "Lb0ToPbarPN0"
+    turbo_line = "Hlt2B2OC_Lb0ToPbarPN0"
+    input_data = get_particles(f"/Event/HLT2/{turbo_line}/Particles")
+
+    myfilter = create_lines_filter(
+        f"LineFilter_{turbo_line}_{{hash}}",
+        [turbo_line],
+    )
+
+    branches = {
+        "Lb0": "[Lambda_b0 -> p+ p~-]CC",
+        "pp": "[Lambda_b0 -> ^p+ p~-]CC",
+        "pm": "[Lambda_b0 -> p+ ^p~-]CC",
+    }
+
+    basic_variables = make_basic_variables(
+        options, pvs, input_data, add_tistos=False)
+    composite_variables = make_composite_variables(
+        options, pvs, input_data, add_tistos=True)
+
+    variables = {
+        "Lb0": composite_variables,
+        "pp": basic_variables,
+        "pm": basic_variables,
+    }
+
+    mytuple = Funtuple(
+        name=name,
+        tuple_name="DecayTree",
+        fields=branches,
+        variables=variables,
+        event_variables=make_hlt2_event_variables(options, pvs, rec_summary),
+        inputs=input_data,
+        store_multiple_cand_info=True)
+
+    return [myfilter, mytuple]
diff --git a/Liv_Lb0NeutronLines/options/Lb_NeutronLines_Spruce.py b/Liv_Lb0NeutronLines/options/Lb_NeutronLines_Spruce.py
new file mode 100644
index 0000000000000000000000000000000000000000..a0504dfbe23b474bcd4bb1bcffeb59dea1153586
--- /dev/null
+++ b/Liv_Lb0NeutronLines/options/Lb_NeutronLines_Spruce.py
@@ -0,0 +1,210 @@
+from .tupling import *
+
+import Functors as F
+from Functors.math import log
+from DaVinci import Options, make_config
+from DaVinci.algorithms import create_lines_filter
+from PyConf.reading import get_particles
+from FunTuple import FunctorCollection
+from PyConf.reading import get_particles, get_pvs, get_extended_pvs
+from PyConf.Algorithms import ParticleUnbiasedPVAdder
+import FunTuple.functorcollections as FC
+from FunTuple import FunTuple_Particles as Funtuple
+from DecayTreeFitter import DecayTreeFitter
+
+
+def maketuple_Lb0ToLcpLcnN0(options, pvs, rec_summary):
+
+    name = "Lb0ToLcpLcmN0"
+    turbo_line = "SpruceB2OC_Lb0ToLcpLcmN0"
+    input_data = get_particles(f"/Event/Spruce/{turbo_line}/Particles")
+
+    myfilter = create_lines_filter(
+        f"LineFilter_{turbo_line}_{{hash}}",
+        [turbo_line],
+    )
+
+    branches = {
+        "Lb0":
+        "[Lambda_b0 -> (Lambda_c+ -> p+ K- pi+) (Lambda_c~- -> p~- K+ pi-)]CC",
+        "Lcp":
+        "[Lambda_b0 -> ^(Lambda_c+ -> p+ K- pi+) (Lambda_c~- -> p~- K+ pi-)]CC",
+        "Lcm":
+        "[Lambda_b0 -> (Lambda_c+ -> p+ K- pi+) ^(Lambda_c~- -> p~- K+ pi-)]CC",
+        "pp":
+        "[Lambda_b0 -> (Lambda_c+ -> ^p+ K- pi+) (Lambda_c~- -> p~- K+ pi-)]CC",
+        "pm":
+        "[Lambda_b0 -> (Lambda_c+ -> p+ K- pi+) (Lambda_c~- -> ^p~- K+ pi-)]CC",
+        "pip":
+        "[Lambda_b0 -> (Lambda_c+ -> p+ K- ^pi+) (Lambda_c~- -> p~- K+ pi-)]CC",
+        "pim":
+        "[Lambda_b0 -> (Lambda_c+ -> p+ K- pi+) (Lambda_c~- -> p~- K+ ^pi-)]CC",
+        "Kp":
+        "[Lambda_b0 -> (Lambda_c+ -> p+ K- pi+) (Lambda_c~- -> p~- ^K+ pi-)]CC",
+        "Km":
+        "[Lambda_b0 -> (Lambda_c+ -> p+ ^K- pi+) (Lambda_c~- -> p~- K+ pi-)]CC",
+    }
+
+    DTF_MASS = DecayTreeFitter(
+        name='DTF_{hash}',
+        input_particles=input_data,
+        mass_constraints=["Lambda_c+"],
+    )
+
+    DTF_MASS_OWNPV = DecayTreeFitter(
+        name='DTF_{hash}',
+        input_particles=input_data,
+        mass_constraints=["Lambda_c+"],
+        constrain_to_ownpv=True,
+    )
+
+    basic_variables = make_basic_variables(options, pvs, input_data)
+    composite_variables = make_composite_variables(options, pvs, input_data)
+    composite_variables_3body = make_composite_variables_3body(
+        options, pvs, input_data)
+
+    DTF_TwoBody = make_composite_dtf_variables(options, pvs, input_data, add_truth=True, DTF=DTF_MASS, pv_constraint=False, mass_constraint=True, particle_name="_Lc", three_body=False) +\
+               make_composite_dtf_variables(options, pvs, input_data, add_truth=True, DTF=DTF_MASS_OWNPV, pv_constraint=True, three_body=False, mass_constraint=True,particle_name="_Lc") #same as BEST
+
+    variables = {
+        "Lb0": composite_variables + DTF_TwoBody,
+        "Lcp": composite_variables_3body,
+        "Lcm": composite_variables_3body,
+        "pp": basic_variables,
+        "pm": basic_variables,
+        "pip": basic_variables,
+        "pim": basic_variables,
+        "Kp": basic_variables,
+        "Km": basic_variables,
+    }
+
+    mytuple = Funtuple(
+        name=name,
+        tuple_name="DecayTree",
+        fields=branches,
+        variables=variables,
+        event_variables=make_hlt2_event_variables(options, pvs, rec_summary),
+        inputs=input_data,
+        store_multiple_cand_info=True)
+
+    return [myfilter, mytuple]
+
+
+def maketuple_BpToLcpLcmKp(options, pvs, rec_summary):
+
+    name = "BuToLcpLcmK"
+    turbo_line = "SpruceB2OC_BuToLcpLcmK_LcpToPKPi"
+    input_data = get_particles(f"/Event/Spruce/{turbo_line}/Particles")
+
+    myfilter = create_lines_filter(
+        f"LineFilter_{turbo_line}_{{hash}}",
+        [turbo_line],
+    )
+
+    branches = {
+        "Bp":
+        "[B+ -> (Lambda_c+ -> p+ K- pi+) (Lambda_c~- -> p~- K+ pi-) K+]CC",
+        "Lcp":
+        "[B+ -> ^(Lambda_c+ -> p+ K- pi+) (Lambda_c~- -> p~- K+ pi-) K+]CC",
+        "Lcm":
+        "[B+ -> (Lambda_c+ -> p+ K- pi+) ^(Lambda_c~- -> p~- K+ pi-) K+]CC",
+        "pp":
+        "[B+ -> (Lambda_c+ -> ^p+ K- pi+) (Lambda_c~- -> p~- K+ pi-) K+]CC",
+        "pm":
+        "[B+ -> (Lambda_c+ -> p+ K- pi+) (Lambda_c~- -> ^p~- K+ pi-) K+]CC",
+        "pip":
+        "[B+ -> (Lambda_c+ -> p+ K- ^pi+) (Lambda_c~- -> p~- K+ pi-) K+]CC",
+        "pim":
+        "[B+ -> (Lambda_c+ -> p+ K- pi+) (Lambda_c~- -> p~- K+ ^pi-) K+]CC",
+        "Kp":
+        "[B+ -> (Lambda_c+ -> p+ K- pi+) (Lambda_c~- -> p~- ^K+ pi-) K+]CC",
+        "Km":
+        "[B+ -> (Lambda_c+ -> p+ ^K- pi+) (Lambda_c~- -> p~- K+ pi-) K+]CC",
+        "Kp_ext":
+        "[B+ -> (Lambda_c+ -> p+ K- pi+) (Lambda_c~- -> p~- K+ pi-) ^K+]CC",
+    }
+
+    DTF_MASS = DecayTreeFitter(
+        name='DTF_{hash}',
+        input_particles=input_data,
+        mass_constraints=["Lambda_c+"],
+    )
+
+    DTF_MASS_OWNPV = DecayTreeFitter(
+        name='DTF_{hash}',
+        input_particles=input_data,
+        mass_constraints=["Lambda_c+"],
+        constrain_to_ownpv=True,
+    )
+
+    basic_variables = make_basic_variables(options, pvs, input_data)
+
+    composite_variables_3body = make_composite_variables_3body(options, pvs, input_data) +\
+        make_composite_dtf_variables(options, pvs, input_data, add_truth=True, DTF=DTF_MASS, pv_constraint=False, three_body=True, mass_constraint=True,particle_name="_Lc") +\
+        make_composite_dtf_variables(options, pvs, input_data, add_truth=True, DTF=DTF_MASS_OWNPV, pv_constraint=True, three_body=True, mass_constraint=True,particle_name="_Lc") #same as BEST
+
+    MVA_variables = make_MVA_variables()
+
+    variables = {
+        "Bp": composite_variables_3body + MVA_variables,
+        "Lcp": composite_variables_3body,
+        "Lcm": composite_variables_3body,
+        "pp": basic_variables,
+        "pm": basic_variables,
+        "pip": basic_variables,
+        "pim": basic_variables,
+        "Kp": basic_variables,
+        "Km": basic_variables,
+        "Kp_ext": basic_variables,
+    }
+
+    mytuple = Funtuple(
+        name=name,
+        tuple_name="DecayTree",
+        fields=branches,
+        variables=variables,
+        event_variables=make_hlt2_event_variables(options, pvs, rec_summary),
+        inputs=input_data,
+        store_multiple_cand_info=True)
+
+    return [myfilter, mytuple]
+
+
+def maketuple_Lb0ToPpPmN0(options, pvs, rec_summary):
+
+    name = "Lb0ToPbarPN0"
+    turbo_line = "SpruceB2OC_Lb0ToPbarPN0"
+    input_data = get_particles(f"/Event/Spruce/{turbo_line}/Particles")
+
+    myfilter = create_lines_filter(
+        f"LineFilter_{turbo_line}_{{hash}}",
+        [turbo_line],
+    )
+
+    branches = {
+        "Lb0": "[Lambda_b0 -> p+ p~-]CC",
+        "pp": "[Lambda_b0 -> ^p+ p~-]CC",
+        "pm": "[Lambda_b0 -> p+ ^p~-]CC",
+    }
+
+    basic_variables = make_basic_variables(
+        options, pvs, input_data, add_tistos=False)
+    composite_variables = make_composite_variables(
+        options, pvs, input_data, add_tistos=True)
+
+    variables = {
+        "Lb0": composite_variables,
+        "pp": basic_variables,
+        "pm": basic_variables,
+    }
+
+    mytuple = Funtuple(
+        name=name,
+        tuple_name="DecayTree",
+        fields=branches,
+        variables=variables,
+        event_variables=make_hlt2_event_variables(options, pvs, rec_summary),
+        inputs=input_data,
+        store_multiple_cand_info=True)
+
+    return [myfilter, mytuple]
diff --git a/Liv_Lb0NeutronLines/options/tupling.py b/Liv_Lb0NeutronLines/options/tupling.py
new file mode 100644
index 0000000000000000000000000000000000000000..c0ec073625f65453ec6e0e23dddf0cc852589014
--- /dev/null
+++ b/Liv_Lb0NeutronLines/options/tupling.py
@@ -0,0 +1,831 @@
+###############################################################################
+# (c) Copyright 2022 CERN for the benefit of the LHCb Collaboration           #
+#                                                                             #
+# This software is distributed under the terms of the GNU General Public      #
+# Licence version 3 (GPL Version 3), copied verbatim in the file "COPYING".   #
+#                                                                             #
+# In applying this licence, CERN does not waive the privileges and immunities #
+# granted to it by virtue of its status as an Intergovernmental Organization  #
+# or submit itself to any jurisdiction.                                       #
+###############################################################################
+"""Common configuration functions
+
+"""
+
+import Functors as F
+from Functors.math import log
+import Functors.math as fmath
+from PyConf.reading import get_pvs
+
+from FunTuple import FunctorCollection
+from FunTuple.functorcollections import (
+    MCHierarchy,
+    MCPromptDecay,
+    Kinematics,
+    SelectionInfo,
+    HltTisTos,
+    MCVertexInfo,
+    MCKinematics,
+    ParticleID,  #wrong variables PID_PI = 0, PROBNN_D = nan
+    EventInfo,
+    ParticleIsolation,
+    MCPrimaries,
+)
+
+from DaVinciMCTools import MCTruthAndBkgCat
+from PyConf.Algorithms import ParticleToSubcombinationsAlg
+from DecayTreeFitter import DecayTreeFitter
+
+Hlt1_global_lines = [
+    "Hlt1GECPassthroughDecision",
+    "Hlt1BeamGasDecision",
+    "Hlt1PassthroughDecision",
+    "Hlt1NoBeamDecision",
+    "Hlt1BeamOneDecision",
+    "Hlt1BeamTwoDecision",
+    "Hlt1BothBeamsDecision",
+    "Hlt1ODINLumiDecision",
+    "Hlt1ODINVeloOpenDecision",
+    "Hlt1ODINNoBiasDecision",
+    "Hlt1VeloMicroBiasDecision",
+    "Hlt1RICH1AlignmentDecision",
+    "Hlt1RICH2AlignmentDecision",
+    "Hlt1BeamGasDecision",
+    "Hlt1L02PPiDecision",
+    "Hlt1LowMassNoipDielectron_massSlice1_promptDecision",
+    "Hlt1LowMassNoipDielectron_SS_massSlice1_promptDecision",
+    "Hlt1LowMassNoipDielectron_massSlice2_promptDecision",
+    "Hlt1LowMassNoipDielectron_SS_massSlice2_promptDecision",
+    "Hlt1LowMassNoipDielectron_massSlice3_promptDecision",
+    "Hlt1LowMassNoipDielectron_SS_massSlice3_promptDecision",
+    "Hlt1LowMassNoipDielectron_massSlice4_promptDecision",
+    "Hlt1LowMassNoipDielectron_SS_massSlice4_promptDecision",
+    "Hlt1LowMassNoipDielectron_massSlice1_displacedDecision",
+    "Hlt1LowMassNoipDielectron_SS_massSlice1_displacedDecision",
+    "Hlt1LowMassNoipDielectron_massSlice2_displacedDecision",
+    "Hlt1LowMassNoipDielectron_SS_massSlice2_displacedDecision",
+    "Hlt1LowMassNoipDielectron_massSlice3_displacedDecision",
+    "Hlt1LowMassNoipDielectron_SS_massSlice3_displacedDecision",
+    "Hlt1LowMassNoipDielectron_massSlice4_displacedDecision",
+    "Hlt1LowMassNoipDielectron_SS_massSlice4_displacedDecision",
+]
+Hlt1_1track_lines = [
+    "Hlt1TrackMVADecision",
+    # "Hlt1LowPtMuonDecision",
+    # "Hlt1SingleHighPtMuonDecision",
+    # "Hlt1SingleHighPtMuonNoMuIDDecision",
+    # "Hlt1TrackMuonMVADecision",
+    # "Hlt1OneMuonTrackLineDecision",
+    # "Hlt1TrackElectronMVADecision",
+    # "Hlt1SingleHighPtElectronDecision",
+    # "Hlt1SingleHighEtDecision",
+]
+Hlt1_lines = Hlt1_1track_lines + [
+    # "Hlt1TwoTrackMVACharmXSecDecision",
+    "Hlt1TwoTrackMVADecision",
+    # "Hlt1TwoTrackKsDecision",
+    # "Hlt1D2KPiDecision",
+    # "Hlt1D2KKDecision",
+    # "Hlt1D2PiPiDecision",
+    # "Hlt1KsToPiPiDecision",
+    # "Hlt1DiMuonNoIPDecision",
+    # "Hlt1DiMuonNoIP_ssDecision",
+    # "Hlt1DiMuonHighMassDecision",
+    # "Hlt1DiMuonSoftDecision",
+    # "Hlt1DiMuonDisplacedDecision",
+    # "Hlt1TwoKsDecision",
+    # "Hlt1D2KPiAlignmentDecision",
+    # "Hlt1DiMuonHighMassAlignmentDecision",
+    # "Hlt1DisplacedDiMuonAlignmentDecision",
+    # "Hlt1DisplacedDielectronDecision",
+]
+
+Hlt2_lines = [
+    "Hlt2Topo2BodyDecision", "Hlt2Topo3BodyDecision",
+    "Hlt2B2OC_Lb0ToLcpLcmN0Decision", "Hlt2B2OC_BuToLcpLcmK_LcpToPKPiDecision",
+    "Hlt2B2OC_Lb0ToPbarPN0Decision", "SpruceB2OC_Lb0ToLcpLcmN0Decision",
+    "SpruceB2OC_BuToLcpLcmK_LcpToPKPiDecision", "SpruceB2OC_Lb0ToPbarPN0Decision",
+]
+
+
+def make_MVA_variables():
+
+    pvs = get_pvs()
+
+    variables = (FunctorCollection({
+        "MVA":
+        F.MVA(
+            MVAType="SigmaNet",
+            Config={
+                "File":
+                "paramfile://data/Hlt2B2OC_B_SigmaNet_Run3-v2.json",
+                "Name":
+                "B2OC_SigmaNet_Generic",
+                "Lambda":
+                "2.0",
+                "NLayers":
+                "3",
+                "InputSize":
+                "6",
+                "Monotone_Constraints":
+                "[1,-1,-1,-1,-1,-1]",
+                "Variables":
+                "log_B_PT,B_ETA,log_B_DIRA,log_B_ENDVERTEX_CHI2,log_B_IPCHI2_OWNPV,log_B_IP_OWNPV",
+            },
+            Inputs={
+                "log_B_PT": fmath.log(F.PT),
+                "B_ETA": F.ETA,
+                "log_B_DIRA": fmath.log(1. + 1.e-6 - F.BPVDIRA(pvs)),
+                "log_B_ENDVERTEX_CHI2": fmath.log(F.CHI2DOF),
+                "log_B_IPCHI2_OWNPV": fmath.log(F.BPVIPCHI2(pvs)),
+                "log_B_IP_OWNPV": fmath.log(F.BPVIP(pvs)),
+            }),
+        "MVAold":
+        F.MVA(
+            MVAType="SigmaNet",
+            Config={
+                "File":
+                "paramfile://data/Hlt2B2OC_B_SigmaNet_Run3-v1.json",
+                "Name":
+                "B2OC_SigmaNet_Generic",
+                "Lambda":
+                "2.0",
+                "NLayers":
+                "3",
+                "InputSize":
+                "6",
+                "Monotone_Constraints":
+                "[1,-1,-1,-1,-1,-1]",
+                "Variables":
+                "log_B_PT,B_ETA,log_B_DIRA,log_B_ENDVERTEX_CHI2,log_B_IPCHI2_OWNPV,log_B_IP_OWNPV",
+            },
+            Inputs={
+                "log_B_PT": fmath.log(F.PT),
+                "B_ETA": F.ETA,
+                "log_B_DIRA": fmath.log(1. - F.BPVDIRA(pvs)),
+                "log_B_ENDVERTEX_CHI2": fmath.log(F.CHI2DOF),
+                "log_B_IPCHI2_OWNPV": fmath.log(F.BPVIPCHI2(pvs)),
+                "log_B_IP_OWNPV": fmath.log(F.BPVIP(pvs)),
+            }),
+    }))
+
+    return variables
+
+
+def make_composite_variables(options,
+                             pvs,
+                             data,
+                             add_truth=True,
+                             add_tistos=True):
+    if not options.simulation:
+        add_truth = False
+    variables = (
+        FunctorCollection({
+            "ID": F.PARTICLE_ID,
+            "CHARGE": F.CHARGE,
+            "VCHI2DOF": F.CHI2DOF,  #CHI2VXNDOF,
+            "MCORR":
+            F.BPVCORRM(pvs),  #"MCorr": F.BPVCORRM(939.56542052), for lb0
+            "MCORRERR": F.BPVCORRMERR(pvs),
+            "M12": F.SUBCOMB(Functor=F.MASS, Indices=(1, 2)),
+            "MAXPT": F.MAX(F.PT),
+            "MINPT": F.MIN(F.PT),
+            "SUMPT": F.SUM(F.PT),
+            "MAXP": F.MAX(F.P),
+            "MINP": F.MIN(F.P),
+            "SUMP": F.SUM(F.P),
+            "BPVDIRA": F.BPVDIRA(pvs),
+            "BPVFDCHI2": F.BPVFDCHI2(pvs),
+            "BPVFD": F.BPVFD(pvs),
+            "BPVVDRHO": F.BPVVDRHO(pvs),
+            "BPVX": F.BPVX(pvs),
+            "BPVY": F.BPVY(pvs),
+            "BPVZ": F.BPVZ(pvs),
+            "BPVIPCHI2": F.BPVIPCHI2(pvs),
+            "LOGBPVIPCHI2": log(F.BPVIPCHI2(pvs)),
+            "BPVIP": F.BPVIP(pvs),
+            "MINIPCHI2": F.MINIPCHI2(pvs),
+            "MINIP": F.MINIP(pvs),
+            "BPVLTIME": F.BPVLTIME(pvs),
+            "MAXBPVIPCHI2": F.MAX(F.BPVIPCHI2(pvs)),  #MAX_
+            "MINBPVIPCHI2": F.MIN(F.BPVIPCHI2(pvs)),
+            "MAXBPVIP": F.MAX(F.BPVIP(pvs)),
+            "MINBPVIP": F.MIN(F.BPVIP(pvs)),
+            "MAXDOCACHI2": F.MAXDOCACHI2,
+            "MAXDOCA": F.MAXDOCA,
+            "DOCA_12": F.DOCA(Child1=1, Child2=2),
+            "DOCACHI2_12": F.DOCACHI2(Child1=1, Child2=2),
+            "ETA": F.ETA,
+            "PHI": F.PHI,
+            "END_VX": F.END_VX,  #END_
+            "END_VY": F.END_VY,
+            "END_VZ": F.END_VZ,
+            "END_VCHI2DOF": F.CHI2DOF @ F.ENDVERTEX,
+        }) + Kinematics())
+
+    if add_tistos:
+        variables += HltTisTos(
+            selection_type="Hlt1", trigger_lines=Hlt1_lines, data=data)
+        # variables += HltTisTos(
+        #     selection_type="Hlt2", trigger_lines=[x for x in Hlt2_lines if x.startswith("Hlt2")], data=data
+        # )
+
+    if add_truth:
+        variables = add_truth_matching_functors(
+            options,
+            variables,
+            data,
+            hierarchy=True,
+            kinematics=True,
+            vertex_info=True,
+            bkgcat=True,
+        )
+
+    return variables
+
+
+def make_tistoscombinations(options, pvs, data):
+    algo_output = ParticleToSubcombinationsAlg(Input=data)
+    relations = algo_output.OutputRelations
+
+    tistos_variables_extra = HltTisTos(
+        selection_type="Hlt1",
+        trigger_lines=["Hlt1TwoTrackMVADecision"],
+        data=algo_output.OutputParticles)
+    tistos_combinations = {}
+    for k, v in tistos_variables_extra.get_thor_functors().items():
+        tistos_combinations[k + "_SUBCOMB"] = F.MAP_INPUT_ARRAY(v, relations)
+    tistos_combinations = FunctorCollection(tistos_combinations)
+
+    return tistos_combinations
+
+
+def make_composite_variables_3body(options,
+                                   pvs,
+                                   data,
+                                   add_truth=True,
+                                   add_tistos=True):
+
+    variables = (
+        FunctorCollection({
+            #13
+            "M13": F.SUBCOMB(Functor=F.MASS, Indices=(1, 3)),
+            "DOCA13": F.DOCA(Child1=1, Child2=3),
+            "DOCACHI213": F.DOCACHI2(Child1=1, Child2=3),
+            "COS13": F.ALV(1, 3),
+            #23
+            "M23": F.SUBCOMB(Functor=F.MASS, Indices=(2, 3)),
+            "DOCA23": F.DOCA(Child1=2, Child2=3),
+            "DOCACHI223": F.DOCACHI2(Child1=2, Child2=3),
+            "COS23": F.ALV(2, 3),
+        }))
+
+    variables += make_composite_variables(options, pvs, data, add_truth,
+                                          add_tistos)
+
+    if add_tistos:
+        variables += make_tistoscombinations(options, pvs, data)
+
+    return variables
+
+
+def make_composite_variables_4body(options, pvs, data, add_truth=True):
+
+    variables = (
+        FunctorCollection({
+            #14
+            "M14": F.SUBCOMB(Functor=F.MASS, Indices=(1, 4)),
+            "DOCA14": F.DOCA(Child1=1, Child2=4),
+            "DOCACHI214": F.DOCACHI2(Child1=1, Child2=4),
+            "COS14": F.ALV(1, 4),
+            #24
+            "M24": F.SUBCOMB(Functor=F.MASS, Indices=(2, 4)),
+            "DOCA24": F.DOCA(Child1=2, Child2=4),
+            "DOCACHI224": F.DOCACHI2(Child1=2, Child2=4),
+            "COS24": F.ALV(2, 4),
+            #34
+            "M34": F.SUBCOMB(Functor=F.MASS, Indices=(3, 4)),
+            "DOCA34": F.DOCA(Child1=3, Child2=4),
+            "DOCACHI234": F.DOCACHI2(Child1=3, Child2=4),
+            "COS34": F.ALV(3, 4),
+        }))
+    return make_composite_variables_3body(options, pvs, data,
+                                          add_truth) + variables
+
+
+def make_DeltaM_variable(options):
+    return FunctorCollection({"DM": F.MASS - F.CHILD(1, F.MASS)})
+
+
+def make_basic_variables(options, pvs, data, add_truth=True, add_tistos=True):
+    if not options.simulation:
+        add_truth = False
+    variables = (FunctorCollection(
+        {
+            "ID": F.PARTICLE_ID,
+            "CHARGE": F.CHARGE,
+            "TRCHI2DOF": F.CHI2DOF @ F.TRACK,
+            "VCHI2DOF": F.CHI2DOF,
+            "TRGHOSTPROB": F.GHOSTPROB,
+            "TX": F.TX,
+            "TY": F.TY,
+            "QOVERP": F.QOVERP @ F.TRACK,
+            "ETA": F.ETA,
+            "PHI": F.PHI,
+            "BPVIPCHI2": F.BPVIPCHI2(pvs),
+            "BPVIP": F.BPVIP(pvs),
+            "MINIPCHI2": F.MINIPCHI2(pvs),
+            "MINIP": F.MINIP(pvs),
+            "ISMUON": F.ISMUON,
+            "NFTHITS": F.VALUE_OR(-1) @ F.NFTHITS @ F.TRACK,
+            "NHITS": F.VALUE_OR(-1) @ F.NHITS @ F.TRACK,
+            "NUTHITS": F.VALUE_OR(-1) @ F.NUTHITS @ F.TRACK,
+            "NVPHITS": F.VALUE_OR(-1) @ F.NVPHITS @ F.TRACK,
+            "TRACKHASVELO": F.VALUE_OR(-1) @ F.TRACKHASVELO @ F.TRACK,
+        }) + Kinematics() + ParticleID(extra_info=True))
+    if add_tistos:
+        variables += HltTisTos(
+            selection_type="Hlt1", trigger_lines=Hlt1_1track_lines, data=data)
+        # variables += HltTisTos(
+        #     selection_type="Hlt2", trigger_lines=[ x for x in Hlt2_lines if x.startswith("Hlt2") ], data=data
+        # )
+
+    if add_truth:
+        variables = add_truth_matching_functors(
+            options,
+            variables,
+            data,
+            hierarchy=True,
+            kinematics=True,
+            vertex_info=True,
+            bkgcat=False,
+        )
+
+    return variables
+
+
+def make_hlt2_event_variables(options, pvs, rec_summary):
+    # define event level variables
+    evt_variables = EventInfo()
+    evt_variables += SelectionInfo(
+        selection_type="Hlt2", trigger_lines=Hlt2_lines)
+    evt_variables += SelectionInfo(
+        selection_type="Hlt1", trigger_lines=Hlt1_lines)
+    evt_variables += FunctorCollection({
+        "ALLPVX":
+        F.ALLPVX(pvs),
+        "ALLPVY":
+        F.ALLPVY(pvs),
+        "ALLPVZ":
+        F.ALLPVZ(pvs),
+        "nPVs":
+        F.SIZE(pvs),
+        "nTracks":
+        F.VALUE_OR(-1) @ F.RECSUMMARY_INFO(rec_summary,
+                                           "nTracks"),  # but always 0 ...
+        "nLongTracks":
+        F.VALUE_OR(-1) @ F.RECSUMMARY_INFO(rec_summary, "nLongTracks"),
+        "nMuonTracks":
+        F.VALUE_OR(-1) @ F.RECSUMMARY_INFO(rec_summary, "nMuonTracks"),
+        "nFTClusters":
+        F.VALUE_OR(-1) @ F.RECSUMMARY_INFO(rec_summary, "nFTClusters"),
+        "nVPClusters":
+        F.VALUE_OR(-1) @ F.RECSUMMARY_INFO(rec_summary, "nVPClusters"),
+        "nUTClusters":
+        F.VALUE_OR(-1) @ F.RECSUMMARY_INFO(rec_summary, "nUTClusters"),
+        "nSPDhits":
+        F.VALUE_OR(-1) @ F.RECSUMMARY_INFO(rec_summary, "nSPDhits"),
+        "nEcalClusters":
+        F.VALUE_OR(-1) @ F.RECSUMMARY_INFO(rec_summary, "nEcalClusters"),
+        "nEcalClusters":
+        F.VALUE_OR(-1) @ F.RECSUMMARY_INFO(rec_summary, "nEcalClusters"),
+        "eCalTot":
+        F.VALUE_OR(-1) @ F.RECSUMMARY_INFO(rec_summary, "eCalTot"),
+        "hCalTot":
+        F.VALUE_OR(-1) @ F.RECSUMMARY_INFO(rec_summary, "hCalTot"),
+        "nRich1Hits":
+        F.VALUE_OR(-1) @ F.RECSUMMARY_INFO(rec_summary, "nRich1Hits"),
+        "nRich2Hits":
+        F.VALUE_OR(-1) @ F.RECSUMMARY_INFO(rec_summary, "nRich2Hits"),
+    })
+    return evt_variables
+
+
+def make_unpack_only_mc_variables():
+    MC_MOTHER_ID = lambda gen: F.VALUE_OR(0) @ F.MC_MOTHER(gen, F.PARTICLE_ID)
+    MC_MOTHER_KEY = lambda gen: F.VALUE_OR(-1) @ F.MC_MOTHER(gen, F.OBJECT_KEY)
+    return (FunctorCollection({
+        "TRUEID": F.PARTICLE_ID,
+        "MC_MOTHER_ID": MC_MOTHER_ID(1),
+        "MC_MOTHER_KEY": MC_MOTHER_KEY(1),
+        "MC_GD_MOTHER_ID": MC_MOTHER_ID(2),
+        "MC_GD_MOTHER_KEY": MC_MOTHER_KEY(2),
+        "MC_GD_GD_MOTHER_ID": MC_MOTHER_ID(3),
+        "MC_GD_GD_MOTHER_KEY": MC_MOTHER_KEY(3),
+        "TRUEORIGIN_VX": F.ORIGIN_VX,
+        "TRUEORIGIN_VY": F.ORIGIN_VY,
+        "TRUEORIGIN_VZ": F.ORIGIN_VZ,
+        "TRUEEND_VX": F.END_VX,
+        "TRUEEND_VY": F.END_VY,
+        "TRUEEND_VZ": F.END_VZ,
+    }) + Kinematics() + MCPromptDecay())
+
+
+####
+
+
+def add_truth_matching_functors(
+        options,
+        variables,
+        data,
+        hierarchy=True,
+        kinematics=True,
+        vertex_info=True,
+        bkgcat=False,
+        prompt_info=True,
+):
+    """Add MC truth matching functors to an existing FunctorCollection.
+
+    Args:
+        options (DaVinci.Options): DaVinci options object.
+        variables (FunctorCollection): FunctorCollection to modify
+        data: data handle
+        hierarchy (bool): Add MCHierarchy info (default True)
+        kinematics (bool): Add MCKinematics info (default True)
+        vertex_info (bool): Add MCVertexInfo (default True)
+        bkgcat (bool): Add TRUEKEY and BKGCAT functors - intended for composite particles (default False)
+        prompt_info (bool): Add MCPromptDecay info (default True)
+
+    Returns:
+        FunctorCollection: modified FunctorCollection with truth matched variables.
+    """
+    assert (
+        options.simulation
+    ), "options.simulation is set to False - it doesn't make sense to add truth matching."
+
+    MCTRUTH = MCTruthAndBkgCat(data)
+
+    if bkgcat:
+        variables += FunctorCollection({
+            # Important note: specify an invalid value for integer functors if there exists no truth info.
+            #                 The invalid value for floating point functors is set to nan.
+            "TRUEKEY":
+            F.VALUE_OR(-1) @ MCTRUTH(F.OBJECT_KEY),
+            "BKGCAT":
+            MCTRUTH.BkgCat,
+        })
+    if hierarchy:
+        variables += MCHierarchy(
+            mctruth_alg=MCTRUTH)  # likely to change again!!
+    if vertex_info:
+        variables += MCVertexInfo(
+            mctruth_alg=MCTRUTH)  # likely to change again!!
+    if kinematics:
+        variables += MCKinematics(
+            mctruth_alg=MCTRUTH)  # likely to change again!!
+    if prompt_info:
+        variables += MCPromptDecay(mctruth_alg=MCTRUTH)
+
+    return variables
+
+
+### DTF variables ###
+
+
+def make_basic_dtf_variables(options,
+                             pvs,
+                             data,
+                             add_truth=True,
+                             DTF=None,
+                             pv_constraint=False,
+                             mass_constraint=False,
+                             particle_name=""):
+    if not options.simulation:
+        add_truth = False
+    variables = (
+        FunctorCollection({
+            "ETA": F.ETA,
+            "PHI": F.PHI,
+            "BPVIPCHI2": F.BPVIPCHI2(pvs),
+            "BPVIP": F.BPVIP(pvs),
+            "TX": F.TX,
+            "TY": F.TY,
+            "MINIPCHI2": F.MINIPCHI2(pvs),
+            "MINIP": F.MINIP(pvs),
+        }) + Kinematics()
+        #+ ParticleID(extra_info=True)
+    )
+
+    dtf_extras = FunctorCollection({
+        "DTF_NITER": DTF.NITER,
+        "DTF_CHI2": DTF.CHI2,
+        "DTF_NDOF": DTF.NDOF,
+        "DTF_CHI2DOF": DTF.CHI2DOF
+    })
+
+    if (mass_constraint):
+        if (pv_constraint):  # MASS + PV
+            dtf_variables_mass_pv = FunctorCollection({
+                'DTF_PV_M' + particle_name + '_' + k: DTF(v)
+                for k, v in variables.get_thor_functors().items()
+            })
+            return dtf_variables_mass_pv
+        else:  # MASS
+            dtf_variables_mass = FunctorCollection({
+                'DTF_M' + particle_name + '_' + k: DTF(v)
+                for k, v in variables.get_thor_functors().items()
+            })
+        return dtf_variables_mass
+    elif (pv_constraint):  # PV
+        dtf_variables_pv = FunctorCollection({
+            'DTF_PV_' + k: DTF(v)
+            for k, v in variables.get_thor_functors().items()
+        })
+        return dtf_variables_pv
+
+    else:  # NO MASS/PV
+        dtf_variables = FunctorCollection({
+            'DTF_' + k: DTF(v)
+            for k, v in variables.get_thor_functors().items()
+        })
+        return dtf_variables
+
+
+def make_composite_dtf_variables(options,
+                                 pvs,
+                                 data,
+                                 add_truth=True,
+                                 DTF=None,
+                                 pv_constraint=False,
+                                 mass_constraint=False,
+                                 particle_name="",
+                                 three_body=False):
+    if not options.simulation:
+        add_truth = False
+    variables = (
+        FunctorCollection({
+            # "MAXPT": F.MAX(F.PT),
+            # "MINPT": F.MIN(F.PT),
+            # "SUMPT": F.SUM(F.PT),
+            # "MAXP": F.MAX(F.P),
+            # "MINP": F.MIN(F.P),
+            # "BPVDIRA": F.BPVDIRA(pvs),
+            # "BPVFDCHI2": F.BPVFDCHI2(pvs),
+            #"BPVFD": F.BPVFD(pvs),
+            # "BPVVDRHO": F.BPVVDRHO(pvs),
+            # "BPVVDZ": F.BPVVDZ(pvs),
+            # "BPVIPCHI2": F.BPVIPCHI2(pvs),
+            # "BPVIP": F.BPVIP(pvs),
+            # "LOGBPVIPCHI2": log(F.BPVIPCHI2(pvs)),
+            #"BPVLTIME": F.BPVLTIME(pvs),
+            # "MAXBPVIPCHI2": F.MAX(F.BPVIPCHI2(pvs)), #MAX_
+            # "MINBPVIPCHI2": F.MIN(F.BPVIPCHI2(pvs)),
+            # "MAXBPVIP": F.MAX(F.BPVIP(pvs)),
+            # "MINBPVIP": F.MIN(F.BPVIP(pvs)),
+            "MCORR":
+            F.BPVCORRM(pvs),  #"MCorr": F.BPVCORRM(939.56542052), for lb0
+            "MCORRERR": F.BPVCORRMERR(pvs),
+            "M12": F.SUBCOMB(Functor=F.MASS, Indices=(1, 2)),
+            # "ETA": F.ETA,
+            # "PHI": F.PHI
+        })
+        # + Kinematics()
+
+        #+ ParticleID(extra_info=True) #only for daughters
+    )
+
+    if three_body:
+        variables += FunctorCollection({
+            "M13":
+            F.SUBCOMB(Functor=F.MASS, Indices=(1, 3)),
+            "M23":
+            F.SUBCOMB(Functor=F.MASS, Indices=(2, 3)),
+        })
+
+
+    if (mass_constraint):
+        if (pv_constraint):  # MASS + PV
+            dtf_extras = FunctorCollection({
+                "DTF_PV_M" + particle_name + "_NITER":
+                DTF.NITER,
+                "DTF_PV_M" + particle_name + "_CHI2DOF" :
+                DTF.CHI2DOF,
+                "DTF_PV_M" + particle_name + "_M":
+                DTF.MASS,
+                "DTF_PV_M" + particle_name + "_MASSERR":
+                DTF.MASSERR,
+            })
+            dtf_variables_mass_pv = FunctorCollection({
+                'DTF_PV_M' + particle_name + '_' + k: DTF(v)
+                for k, v in variables.get_thor_functors().items()
+            })
+            return dtf_variables_mass_pv + dtf_extras
+        else:  # MASS
+            dtf_extras = FunctorCollection({
+                "DTF_M" + particle_name + "_NITER":
+                DTF.NITER,
+                "DTF_M" + particle_name + "_CHI2DOF" :
+                DTF.CHI2DOF,
+                "DTF_M" + particle_name + "_M":
+                DTF.MASS,
+                "DTF_M" + particle_name + "_MASSERR":
+                DTF.MASSERR,
+            })
+            dtf_variables_mass = FunctorCollection({
+                'DTF_M' + particle_name + '_' + k: DTF(v)
+                for k, v in variables.get_thor_functors().items()
+            })
+        return dtf_variables_mass + dtf_extras
+
+    elif (pv_constraint):  # PV
+        dtf_extras = FunctorCollection({
+                "DTF_PV" + particle_name + "_NITER":
+                DTF.NITER,
+                "DTF_PV" + particle_name + "_CHI2DOF" :
+                DTF.CHI2DOF,
+                "DTF_PV" + particle_name + "_M":
+                DTF.MASS,
+                "DTF_PV" + particle_name + "_MASSERR":
+                DTF.MASSERR,
+            })
+        dtf_variables_pv = FunctorCollection({
+            'DTF_PV_' + k: DTF(v)
+            for k, v in variables.get_thor_functors().items()
+        })
+        return dtf_variables_pv + dtf_extras
+
+    else:  # NO MASS/PV
+        dtf_extras = FunctorCollection({
+                "DTF" + particle_name + "_NITER":
+                DTF.NITER,
+                "DTF" + particle_name + "_CHI2DOF" :
+                DTF.CHI2DOF,
+                "DTF" + particle_name + "_M":
+                DTF.MASS,
+                "DTF" + particle_name + "_MASSERR":
+                DTF.MASSERR,
+            })
+        dtf_variables = FunctorCollection({
+            'DTF_' + k: DTF(v)
+            for k, v in variables.get_thor_functors().items()
+        })
+        return dtf_variables + dtf_extras
+
+
+def make_MC_basic_variables():
+    variables = (FunctorCollection({
+        "OBJECT_KEY": F.OBJECT_KEY,
+        "ETA": F.ETA,
+        "PHI": F.PHI,
+        "ORIGIN_VX": F.ORIGIN_VX,
+        "ORIGIN_VY": F.ORIGIN_VY,
+        "ORIGIN_VZ": F.ORIGIN_VZ,
+        "FOURMOMENTUM": F.FOURMOMENTUM,
+    }) + Kinematics())
+
+    return variables
+
+
+def make_MC_composite_variables():
+    variables = (
+        FunctorCollection({
+            "END_VX": F.END_VX,  #END_
+            "END_VY": F.END_VY,
+            "END_VZ": F.END_VZ,
+            "LTIME": F.MC_LIFETIME,
+            "OBJECT_KEY": F.OBJECT_KEY,
+            "ETA": F.ETA,
+            "PHI": F.PHI,
+            "ORIGIN_VX": F.ORIGIN_VX,
+            "ORIGIN_VY": F.ORIGIN_VY,
+            "ORIGIN_VZ": F.ORIGIN_VZ,
+            "FOURMOMENTUM": F.FOURMOMENTUM,
+        }) + Kinematics())
+
+    return variables
+
+
+def make_MC_event_variables(mc_header):
+    # define event level variables
+    evt_variables = EventInfo()
+    evt_variables += MCPrimaries(mc_header=mc_header)
+
+    return evt_variables
+
+
+def make_top_isolation_variables(hlt2_line,
+                                 input_data,
+                                 locations=["LongTrackIso", "NeutralIso"]):
+    from PyConf.reading import get_particles
+    from IsolationTools import VertexAndConeIsolation
+
+    possible_charm_locations = [
+        "LongTrackIso", "TTrackIso", "DownstreamTrackIso", "UpstreamTrackIso",
+        "NeutralIso", "PizIso"
+    ]
+    coneangles = [0.25, 0.5, 1., 1.5, 2.]
+
+    count = 0
+    for location in locations:
+        extra_particles = get_particles(
+            f"/Event/HLT2/{hlt2_line}/{location}/Particles")
+
+        for coneangle in coneangles:
+            top_RTAlg = VertexAndConeIsolation(
+                name=location + "_" + str(coneangle),
+                reference_particles=input_data,
+                related_particles=extra_particles,
+                cut=F.DR2 < coneangle)
+
+            if count == 0:
+                top_iso_variables = ParticleIsolation(
+                    isolation_alg=top_RTAlg, array_indx_name='indx')
+            else:
+                top_iso_variables += ParticleIsolation(
+                    isolation_alg=top_RTAlg, array_indx_name='indx')
+            count += 1
+
+    return top_iso_variables
+
+
+def make_basic_isolation_variables(hlt2_line,
+                                   input_data,
+                                   locations=["LongTrackIso", "NeutralIso"]):
+    from PyConf.reading import get_particles
+    from IsolationTools import VertexAndConeIsolation
+    from PyConf.Algorithms import ThOrParticleSelection
+
+    basic_code = (F.FILTER(F.ALL) @ F.GET_ALL_BASICS())
+    basic_particles = ThOrParticleSelection(
+        InputParticles=input_data, Functor=basic_code).OutputSelection
+
+    possible_charm_locations = [
+        "LongTrackIso", "TTrackIso", "DownstreamTrackIso", "UpstreamTrackIso",
+        "NeutralIso", "PizIso"
+    ]
+    coneangles = [0.25, 0.5, 1., 1.5, 2.]
+
+    count = 0
+    for location in locations:
+        extra_particles = get_particles(
+            f"/Event/HLT2/{hlt2_line}/{location}/Particles")
+
+        for coneangle in coneangles:
+            basic_RTAlg = VertexAndConeIsolation(
+                name=location + "_" + str(coneangle),
+                reference_particles=basic_particles,
+                related_particles=extra_particles,
+                cut=F.DR2 < coneangle)
+
+            if count == 0:
+                basic_iso_variables = ParticleIsolation(
+                    isolation_alg=basic_RTAlg, array_indx_name='indx')
+            else:
+                basic_iso_variables += ParticleIsolation(
+                    isolation_alg=basic_RTAlg, array_indx_name='indx')
+            count += 1
+
+    return basic_iso_variables
+
+
+def make_intermediate_isolation_variables(
+        hlt2_line,
+        input_data,
+        locations=["LongTrackIso", "NeutralIso"],
+        composite_ID="J/psi(1S)"):
+    from PyConf.reading import get_particles
+    from IsolationTools import VertexAndConeIsolation
+    from PyConf.Algorithms import ThOrParticleSelection
+
+    intermediate_code = F.FILTER(
+        F.IS_ABS_ID(composite_ID)) @ F.GET_ALL_DESCENDANTS()
+    intermediate_particles = ThOrParticleSelection(
+        InputParticles=input_data, Functor=intermediate_code).OutputSelection
+
+    possible_charm_locations = [
+        "LongTrackIso", "TTrackIso", "DownstreamTrackIso", "UpstreamTrackIso",
+        "NeutralIso", "PizIso"
+    ]
+    coneangles = [0.25, 0.5, 1., 1.5, 2.]
+
+    count = 0
+    for location in locations:
+        extra_particles = get_particles(
+            f"/Event/HLT2/{hlt2_line}/{location}/Particles")
+
+        for coneangle in coneangles:
+            intermediate_RTAlg = VertexAndConeIsolation(
+                name=location + "_" + str(coneangle),
+                reference_particles=intermediate_particles,
+                related_particles=extra_particles,
+                cut=F.DR2 < coneangle)
+
+            if count == 0:
+                intermediate_iso_variables = ParticleIsolation(
+                    isolation_alg=intermediate_RTAlg, array_indx_name='indx')
+            else:
+                intermediate_iso_variables += ParticleIsolation(
+                    isolation_alg=intermediate_RTAlg, array_indx_name='indx')
+            count += 1
+
+    return intermediate_iso_variables