B2A403-KFit-VertexFit.py

#!/usr/bin/env python3
 
 
 
 
 
import basf2 as b2
from modularAnalysis import inputMdst
from modularAnalysis import reconstructDecay
from modularAnalysis import matchMCTruth
from vertex import kFit
from stdCharged import stdPi, stdK
from modularAnalysis import variablesToNtuple
import variables.collections as vc
import variables.utils as vu
 
# create path
my_path = b2.create_path()
 
# load input ROOT file
inputMdst(filename=b2.find_file('B02pi0D0_D2kpi_B2Dstarpi_Dstar2Dpi_D2kpi.root', 'examples', False),
          path=my_path)
 
 
# use standard final state particle lists
#
# creates "pi+:all" ParticleList (and c.c.)
stdPi('all', path=my_path)
# creates "pi+:loose" ParticleList (and c.c.)
stdPi('loose', path=my_path)
# creates "K+:loose" ParticleList (and c.c.)
stdK('loose', path=my_path)
 
# reconstruct D0 -> K- pi+ decay
# keep only candidates with 1.8 < M(Kpi) < 1.9 GeV
reconstructDecay('D0:kpi -> K-:loose pi+:loose', '1.8 < M < 1.9', path=my_path)
 
# Perform D0 vertex fit.
# Reject the candidates with failed fit.
kFit('D0:kpi', 0.0, path=my_path)
 
# reconstruct D*+ -> D0 pi+ decay
# keep only candidates with Q = M(D0pi) - M(D0) - M(pi) < 20 MeV
# and D* CMS momentum > 2.5 GeV
reconstructDecay('D*+ -> D0:kpi pi+:all', '0.0 <= Q < 0.02 and 2.5 < useCMSFrame(p) < 5.5', path=my_path)
 
# perform MC matching (MC truth association)
matchMCTruth('D*+', path=my_path)
 
# perform D*+ vertex fit
# keep candidates only passing C.L. value of the fit > 0.0 (no cut)
kFit('D*+', 0.0, path=my_path)
 
# Select variables that we want to store to ntuple
 
dstar_vars = vc.inv_mass + vc.mc_truth
 
fs_hadron_vars = vu.create_aliases_for_selected(
    vc.pid + vc.track + vc.mc_truth,
    'D*+ -> [D0 -> ^K- ^pi+] ^pi+')
 
d0_vars = vu.create_aliases_for_selected(
    vc.inv_mass + vc.mc_truth,
    'D*+ -> ^D0 pi+', 'D0')
 
 
# Saving variables to ntuple
output_file = 'B2A403-KFit-VertexFit.root'
variablesToNtuple('D*+', dstar_vars + d0_vars + fs_hadron_vars,
                  filename=output_file, treename='dsttree', path=my_path)
 
# Process the events
b2.process(my_path)
 
# print out the summary
print(b2.statistics)