Skip to content

Introducing Particle Combiners

Learning Objectives

  • Understanding how particle combiners work in LHCb
  • Learning how to filter through and combine tracks to form decay candidates

We have already seen how we can use DaVinci to process the data we have collected (or simulated) into a ROOT file we can then use for analysis. However, the DaVinci software's capabilities go much further.

In this lesson, we will see how the particles are combined to form decay descriptors, and how we can apply different filters on the tracks or combination of, to build our physics cases. This is the basis for writing a Sprucing or HLT2 line (as will be explained in their dedicated lesson), for tupling NoBias data (as will be explored in the next lesson), and for checking the effect on different selections in simulated samples before any cuts are applied, among other capabilities.

In Run 3, DaVinci is built on top of Moore, the application that configures HLT2. An advantage of this is that it allows us to use Moore's functions to build our own decays and then write the output into a ROOT file for inspection.

Take the example from the First Analysis Steps lesson, which took simulated candidates that had passed an HLT2 line selecting \(B_s^0\to (D_s^- \to K^- K^+ \pi^-) \pi^+\). The line essentially selects three hadron tracks -- two kaons and a pion -- and fits them to a common vertex, which we associate to a \(D_s^-\). Then, it combines the particle associated to that vertex to another pion, to fit both of them into another vertex, associated to a \(B_s^0\).

In addition, we can set some requirements for each of these tracks. For instance, we may be interested in kaons and pions with a high IPCHI2, such that they are less likely to originate from the PV; transverse momentum can also be used for this purpose. We might also be interested in some particle identification on the kaons, like a cut on PID_K.

Then we can also set cuts on the combination and the vertex. Since we want our hadrons to come from a \(D_s^-\), we can ask for this by setting the invariant mass of the three tracks to be in the range of the parent particle's mass. We can also ask for the fit to the vertex to be good by setting a cut on CHI2DOF, and we increase our chances that the three particles come from the same point by cutting on the distance of closest approach between pairs.

We now have an object ds associated to a good vertex where two kaons and a pion originated. We can now combine this object with another pion in the event and form the \(B_s^0\) (or \(B^0\)) vertex.

This is how we combine the three hadrons to form a \(D_s^-\), and then the \(D_s^-\) with another pion to form the \(B\) meson:

def setup_selections() -> algorithms_thor.ParticleCombiner:
    '''
    Builds the D_s from two kaons and a pion with some basic selection
    Then builds the B_s from the D_s and a pion, also with basic selection
    Finally returns the B_s combination
    '''
    import Functors as F
    from Functors.math import in_range
    from RecoConf import standard_particles, algorithms_thor
    from RecoConf.reconstruction_objects import make_pvs

    from GaudiKernel.SystemOfUnits import GeV, MeV, mm, picosecond

    make_pions = standard_particles.make_long_pions()
    make_kaons = standard_particles.make_long_kaons()

    ## We need to import the PVs to calculate associated variables
    pvs = make_pvs()

    # Apply some cuts to the pions and kaons
    code_pions = F.require_all(F.PT > 250 * MeV, F.MINIPCHI2(pvs) > 4,)
    code_kaons = F.require_all(F.PT > 250 * MeV, F.MINIPCHI2(pvs) > 4, F.PID_K > 5) # Kaons can have an extra PID cut

    pions = algorithms_thor.ParticleFilter(make_pions, F.FILTER(code_pions))
    kaons = algorithms_thor.ParticleFilter(make_kaons, F.FILTER(code_kaons))

    # We now also apply cuts on the parent D_s
    ds_mass_bounds = in_range(1920 * MeV, F.MASS, 2010 * MeV)

    ds_combination = F.require_all(ds_mass_bounds, F.SDOCA(1, 3) < 0.2 * mm, 
                                F.SDOCA(2, 3) < 0.2 * mm)
    ds_vertex = F.require_all(ds_mass_bounds, F.CHI2DOF < 9)

    # Finally we combine the child particles into the parent
    ds = algorithms_thor.ParticleCombiner(
        [kaons, kaons, pions], ## The order is important here, must match that of the DecayDescriptor
        name="DsToKpKmPip",
        DecayDescriptor="[D_s- -> K+ K- pi-]cc",
        CombinationCut=ds_combination,
        CompositeCut=ds_vertex)

    # Cuts for the B_s (which we make from a pion and the D_s)
    b_mass_bounds = in_range(5050 * MeV, F.MASS, 5650 * MeV)

    b_combination = F.require_all(b_mass_bounds)
    b_vertex = F.require_all(b_mass_bounds, F.BPVIPCHI2(pvs) < 16,
                            F.BPVLTIME(pvs) > 0.2 * picosecond)

    # Combine the pion and the D_s
    bmeson = algorithms_thor.ParticleCombiner(
        [ds, pions],
        name="BToDsmPip",
        DecayDescriptor="[B0 -> D_s- pi+]cc",
        CombinationCut=b_combination,
        CompositeCut=b_vertex
    )

    return bmeson

We can then use this to build our own decay on simulated data:

'''
lb-run -c x86_64_v3-el9-gcc13+detdesc-opt+g DaVinci/v65r0 lbexec DV_advanced1:main DVoptions.yaml
'''
import Functors as F
import FunTuple.functorcollections as FC
from DaVinci import Options, make_config
from DaVinci.algorithms import create_lines_filter
from FunTuple import FunctorCollection
from FunTuple import FunTuple_Particles as Funtuple
from PyConf.reading import get_particles, get_pvs
from RecoConf.event_filters import require_pvs


def main(options: Options):
    bmeson = setup_selections() # Taken from the above snippet
    pvs = get_pvs()

    fields = {
        "Bs": "[B0 -> (D_s- -> K- K+ pi-) pi+]CC", #Note here, we call the parent a B0 even though it's a Bs, because it needs to match the reconstruction in the trigger line
        "Ds": "[B0 -> ^(D_s- -> K- K+ pi-) pi+]CC",
        "Bs_pi": "[B0 -> (D_s- -> K- K+ pi-) ^pi+]CC",
        "Ds_pi": "[B0 -> (D_s- -> K- K+ ^pi-) pi+]CC",
        "Km": "[B0 -> (D_s- -> ^K- K+ pi-) pi+]CC",
        "Kp": "[B0 -> (D_s- -> K- ^K+ pi-) pi+]CC",
    }

    all_vars = FunctorCollection({
        "M": F.MASS,
        "P": F.P,
        "PT": F.PT
    })

    variables = {"ALL": all_vars}

    funtuple = Funtuple(
        name="MyTuple",
        tuple_name="DecayTree",
        fields=fields,
        variables=variables,
        inputs=bmeson,
    )

    return make_config(options, [require_pvs(pvs), funtuple])

The options file is the same as in the previous DaVinci lesson:

input_files:
- root://eoslhcb.cern.ch//eos/lhcb/wg/dpa/wp7/Run3SK/exampleDST/00257716_00000011_1.hlt2.dst
input_type: ROOT
output_type: ROOT
input_raw_format: 0.5
simulation: true
input_process: Hlt2
input_stream: b2oc
conddb_tag: sim10-2024.W37.39-v00.00-md100
dddb_tag: dddb-20240427
lumi: false
data_type: '2024'
evt_max: 1000
print_freq: 100
ntuple_file: tuple.root