Skip to content
Snippets Groups Projects

Add SL WG production for Run 3 charged isolation studies

Merged Mark Smith requested to merge masmith/SL_run3_isolation into master
2 files
+ 454
0
Compare changes
  • Side-by-side
  • Inline
Files
2
+ 291
0
import sys
from DaVinci import Options
from PyConf.reading import get_pvs, get_rec_summary, get_particles
import Functors as F
from FunTuple import FunctorCollection
from FunTuple import FunTuple_Particles as Funtuple
from DaVinci.algorithms import create_lines_filter
from DaVinci import make_config
import FunTuple.functorcollections as FC
from DaVinciMCTools import MCTruthAndBkgCat
from Hlt2Conf.algorithms_thor import ParticleCombiner, ParticleFilter
from Hlt2Conf.algorithms_thor import ParticleContainersMerger
from Functors.grammar import Functor
_BPVCORRM = Functor('_BPVCORRM', 'Composite::CorrectedMass','Compute the corrected mass')
_allpv_FDVEC = (F.ENDVERTEX_POS @ F.FORWARDARG1 - F.TOLINALG @ F.POSITION @ F.FORWARDARG0)
_allpv_NORMEDDOT = F.NORMEDDOT.bind(F.THREEMOMENTUM @ F.FORWARDARG1, _allpv_FDVEC)
ALLPV_IPCHI2 = lambda Vertices: F.MAP(F.IPCHI2).bind(F.TES(Vertices), F.FORWARDARGS)
ALLPV_FDCHI2 = lambda Vertices: F.MAP(F.VTX_FDCHI2).bind(F.TES(Vertices), F.FORWARDARGS)
ALLPV_DIRA = lambda Vertices: F.MAP(_allpv_NORMEDDOT).bind(F.TES(Vertices), F.FORWARDARGS)
ALLPV_CORRM = lambda Vertices: F.MAP(_BPVCORRM).bind(F.TES(Vertices), F.FORWARDARGS)
ALLPV_LTIME = lambda Vertices: F.MAP(F.VTX_LTIME).bind(F.TES(Vertices), F.FORWARDARGS)
ALLPV_DLS = lambda Vertices: F.MAP(F.VTX_DLS).bind(F.TES(Vertices), F.FORWARDARGS)
ALLPV_FD_COORDINATE = lambda coordinate, Vertices: F.MAP(coordinate @ _allpv_FDVEC).bind(F.TES(Vertices), F.FORWARDARGS)
TRUE_ID_IS = lambda id, mctruth: (F.VALUE_OR(0) @ mctruth(F.ABS @ F.PARTICLE_ID) == id)
def make_Bst(B, pions, hadron_candidate_id):
#filter B keeping only true D0 and true leptons
MCTRUTH_B = MCTruthAndBkgCat(B, name='MCTRUTH_Bu')
Dz_child = F.CHILD(1, F.FORWARDARG0) #get the first child of B- which is D0
lep_child = F.CHILD(2, F.FORWARDARG0) #get the second child of B- which is lepton
Dz_truth_cut = TRUE_ID_IS(hadron_candidate_id, MCTRUTH_B) @ Dz_child #require that D0 is true
lep_truth_cut = F.require_any(TRUE_ID_IS(11, MCTRUTH_B) @ lep_child , TRUE_ID_IS(13, MCTRUTH_B) @ lep_child) #require that lepton is true i.e. electron or muon
cut_b = F.require_all(Dz_truth_cut, lep_truth_cut) #make a cut
signal_B = ParticleFilter(Input=B, Cut=F.FILTER(cut_b), name='signal_Dzmu') #filter the B candidates
#combine to make [B*0 -> B- pi+]cc
Bst_1 = ParticleCombiner(
Inputs=[signal_B, pions],
name='Bst_1',
DecayDescriptor='[B*0 -> B- pi+]cc',
CombinationCut=F.ALL,
CompositeCut=F.ALL,
ParticleCombiner="ParticleVertexFitter", #for neutrals need to use different combiner e.g. ParticleAdder
#OutputLevel=1
)
#combine to make [B*0 -> B- pi-]cc
Bst_2 = ParticleCombiner(
Inputs=[signal_B, pions],
name='Bst_2',
DecayDescriptor='[B*0 -> B- pi-]cc',
CombinationCut=F.ALL,
CompositeCut=F.ALL,
ParticleCombiner="ParticleVertexFitter", #for neutrals need to use different combiner e.g. ParticleAdder
#OutputLevel=1
)
#merge the two Bst candidates
Bst = ParticleContainersMerger([Bst_1, Bst_2], name = 'Bst_combiner')
return Bst
def get_functors(MCTRUTH, v2_pvs, rec_summary):
#define common variables for all fields (composite and basic)
vars_common = FunctorCollection()
#all pvs: ip, ipchi2. The PV positions are stored in event variables
vars_common['ALLPV_IP[pv_indx]'] = F.ALLPV_IP(v2_pvs)
vars_common['ALLPV_IPCHI2[pv_indx]'] = ALLPV_IPCHI2(v2_pvs)
#best pv: x, y, z, ip, ipchi2
vars_common['BPV_X'] = F.BPVX(v2_pvs)
vars_common['BPV_Y'] = F.BPVY(v2_pvs)
vars_common['BPV_Z'] = F.BPVZ(v2_pvs)
vars_common['BPV_IP'] = F.BPVIP(v2_pvs)
vars_common['BPV_IPCHI2'] = F.BPVIPCHI2(v2_pvs)
vars_common['BPV_CHI2'] = F.CHI2 @ F.BPV(v2_pvs)
vars_common['BPV_CHI2DOF'] = F.CHI2DOF @ F.BPV(v2_pvs)
#particle id, key, truekey. The trueid is stored in MCHierarchy
vars_common['ID'] = F.PARTICLE_ID
vars_common['KEY'] = F.OBJECT_KEY
vars_common['TRUEKEY'] = F.VALUE_OR(-1) @ MCTRUTH(F.OBJECT_KEY)
#get charge, min ip and min ipchi2
vars_common['CHARGE'] = F.CHARGE
vars_common['MINIP'] = F.MINIP(v2_pvs)
vars_common['MINIPCHI2'] = F.MINIPCHI2(v2_pvs)
#reco kinematics
vars_common += FC.Kinematics()
vars_common['ETA'] = F.ETA
vars_common['PHI'] = F.PHI
#mc truth kinematics
vars_common += FC.MCKinematics(mctruth_alg = MCTRUTH)
vars_common['TRUEETA'] = MCTRUTH(F.ETA)
vars_common['TRUEPHI'] = MCTRUTH(F.PHI)
#type of the origin vertex
vars_common['MC_VTX_TYPE'] = MCTRUTH(F.VALUE_OR(-1) @ F.MC_VTX_TYPE @ F.MC_ORIGINVERTEX)
##mc truth hierarchy (done seperately)
#vars_common += FC.MCHierarchy(mctruth_alg = MCTRUTH)
#make some helper functions
MCMOTHER_ID = lambda n: F.VALUE_OR(0) @ MCTRUTH(F.MC_MOTHER(n, F.PARTICLE_ID))
MCMOTHER_KEY = lambda n: F.VALUE_OR(-1) @ MCTRUTH(F.MC_MOTHER(n, F.OBJECT_KEY ))
vars_common["TRUEID"] = F.VALUE_OR(0) @ MCTRUTH(F.PARTICLE_ID)
vars_common["MCKEY"] = F.VALUE_OR(0) @ MCTRUTH(F.OBJECT_KEY)
vars_common["MC_MOTHER_ID"] = MCMOTHER_ID(1)
vars_common["MC_MOTHER_KEY"] = MCMOTHER_KEY(1)
for i in range(2,14):
prefix = "MC_GD_MOTHER_{}".format(i)
vars_common[f"{prefix}_ID"] = MCMOTHER_ID(i)
vars_common[f"{prefix}_KEY"] = MCMOTHER_KEY(i)
#variables for composite particles
vars_composite = FunctorCollection()
#end vertex position and end vertex chi2
vars_composite['ENDVERTEX_POS_'] = F.ENDVERTEX_POS
vars_composite['ENDVERTEX_CHI2'] = F.CHI2 @ F.ENDVERTEX
vars_composite['ENDVERTEX_CHI2DOF'] = F.CHI2DOF @ F.ENDVERTEX
#all pvs: dira, fd, fdchi2
vars_composite['ALLPV_DIRA[pv_indx]'] = ALLPV_DIRA(v2_pvs)
vars_composite['ALLPV_FD_X[pv_indx]'] = ALLPV_FD_COORDINATE(F.X_COORDINATE, v2_pvs)
vars_composite['ALLPV_FD_Y[pv_indx]'] = ALLPV_FD_COORDINATE(F.Y_COORDINATE, v2_pvs)
vars_composite['ALLPV_FD_Z[pv_indx]'] = ALLPV_FD_COORDINATE(F.Z_COORDINATE, v2_pvs)
vars_composite['ALLPV_FDCHI2[pv_indx]'] = ALLPV_FDCHI2(v2_pvs)
vars_composite['ALLPV_CORRM[pv_indx]'] = ALLPV_CORRM(v2_pvs)
vars_composite['ALLPV_LTIME[pv_indx]'] = ALLPV_LTIME(v2_pvs)
vars_composite['ALLPV_DLS[pv_indx]'] = ALLPV_DLS(v2_pvs)
#best pv: dira, fd, fdchi2, corrm, ltime, dls
vars_composite['BPV_DIRA'] = F.BPVDIRA(v2_pvs)
vars_composite['BPV_FD_'] = F.BPVFDVEC(v2_pvs)
vars_composite['BPV_FDCHI2'] = F.BPVFDCHI2(v2_pvs)
vars_composite['BPV_CORRM'] = F.BPVCORRM(v2_pvs)
vars_composite['BPV_LTIME'] = F.BPVLTIME(v2_pvs)
vars_composite['BPV_DLS'] = F.BPVDLS(v2_pvs)
#mc composite vertex information
vars_composite += FC.MCVertexInfo(mctruth_alg = MCTRUTH)
#Compute maximum DOCA and maximum DOCACHI2. Since there are
# only two daughters of Sb particles, the maximum DOCA/DOCACHI2 is
# the DOCA/DOCACHI2 see below.
vars_composite['MAX_DOCA'] = F.MAXDOCA
vars_composite['MAX_DOCACHI2'] = F.MAXDOCACHI2
vars_composite['MAX_SDOCA'] = F.MAXSDOCA
vars_composite['MAX_SDOCACHI2'] = F.MAXSDOCACHI2
#variables for Lb field
var_B = FunctorCollection()
#variables for basics
vars_basic = FunctorCollection()
vars_basic['TRACKTYPE'] = F.VALUE_OR(-1) @ F.CAST_TO_INT @ F.TRACKTYPE @ F.TRACK
vars_basic['TRACKHISTORY'] = F.VALUE_OR(-1) @ F.CAST_TO_INT @ F.TRACKHISTORY @ F.TRACK
vars_basic['TRACKFLAG'] = F.VALUE_OR(-1) @ F.CAST_TO_INT @ F.TRACKFLAG @ F.TRACK
vars_basic['TRACKNDOF'] = F.VALUE_OR(-1) @ F.NDOF @ F.TRACK
vars_basic['TRACKCHI2'] = F.CHI2 @ F.TRACK
vars_basic['TRACKCHI2DOF'] = F.CHI2DOF @ F.TRACK
vars_basic['GHOSTPROB'] = F.GHOSTPROB
vars_basic['PIDpi'] = F.PID_PI
vars_basic['PIDk'] = F.PID_K
vars_basic['PIDp'] = F.PID_P
vars_basic['PIDe'] = F.PID_E
vars_basic['PIDmu'] = F.PID_MU
vars_basic['TRUEORIGINVERTEX_X'] = MCTRUTH(F.ORIGIN_VX)
vars_basic['TRUEORIGINVERTEX_Y'] = MCTRUTH(F.ORIGIN_VY)
vars_basic['TRUEORIGINVERTEX_Z'] = MCTRUTH(F.ORIGIN_VZ)
#variables for event
evt_vars = FunctorCollection()
evt_vars['nTracks'] = F.VALUE_OR(-1) @ F.RECSUMMARY_INFO(rec_summary, "nTracks")
evt_vars['nLongTracks'] = F.VALUE_OR(-1) @ F.RECSUMMARY_INFO(rec_summary, "nLongTracks")
evt_vars['nPVs'] = F.VALUE_OR(-1) @ F.RECSUMMARY_INFO(rec_summary, "nPVs")
#return all functors
functors = (vars_common, vars_composite, var_B, vars_basic, evt_vars)
return functors
def tuple_Bst(Bst, MCTRUTH, v2_pvs, rec_summary, hadron_name, get_electron = False):
#get functors
functors = get_functors(MCTRUTH, v2_pvs, rec_summary)
vars_common = functors[0]
vars_composite = functors[1]
var_B = functors[2]
vars_basic = functors[3]
evt_vars = functors[4]
#variables for Sb field
vars_Bst = FunctorCollection()
B_child = F.CHILD(1, F.FORWARDARG0)
Epi_child = F.CHILD(2, F.FORWARDARG0)
bpv_pi = F.BPV(v2_pvs).bind(Epi_child)
#diff in vtx chi2 with and without extra particle
vars_Bst['DELTA_ENDVERTEX_CHI2'] = (F.CHI2 @ F.ENDVERTEX) - (F.CHI2 @ F.ENDVERTEX.bind(B_child))
#IP and IPChi2 of extra particle wrt to Sb end vertex
vars_Bst['Epi_IP_WRT_SbENDVERTEX'] = F.IP.bind(F.TOLINALG @ F.POSITION @ F.ENDVERTEX @ F.FORWARDARG0 , Epi_child)
vars_Bst['Epi_IPCHI2_WRT_SbENDVERTEX'] = F.IPCHI2.bind(F.ENDVERTEX @ F.FORWARDARG0 , Epi_child)
#IP and IPChi2 of extra particle wrt to Lb end vertex
vars_Bst['Epi_IP_WRT_LbENDVERTEX'] = F.IP.bind(F.TOLINALG @ F.POSITION @ F.ENDVERTEX @ B_child , Epi_child)
vars_Bst['Epi_IPCHI2_WRT_LbENDVERTEX'] = F.IPCHI2.bind(F.ENDVERTEX @ B_child , Epi_child)
#angle between extra particle and Lb
vars_Bst['Epi_COSANGLE_Lb'] = F.COSANGLE.bind(F.THREEMOMENTUM @ B_child, F.THREEMOMENTUM @ Epi_child)
#diff in fd chi2 with and without extra particle
vars_Bst['Delta_BPV_of_Epi_FDCHI2'] = F.VTX_FDCHI2.bind(bpv_pi, F.FORWARDARGS) - F.VTX_FDCHI2.bind(bpv_pi, B_child)
vars_Bst['BPV_of_Epi_FDCHI2'] = F.VTX_FDCHI2.bind(bpv_pi, F.FORWARDARGS)
#IP chi2 of extra particle wrt PV and SV of Sb
vars_Bst['Epi_IP_WRT_SbBPV'] = F.IP.bind(F.TOLINALG @ F.POSITION @ F.BPV(v2_pvs) @ F.FORWARDARG0, Epi_child)
vars_Bst['Epi_IPCHI2_WRT_SbBPV'] = F.IPCHI2.bind(F.BPV(v2_pvs) @ F.FORWARDARG0, Epi_child)
#IP chi2 of extra particle wrt PV and SV of Bb
vars_Bst['Epi_IP_WRT_LbBPV'] = F.IP.bind(F.TOLINALG @ F.POSITION @ F.BPV(v2_pvs) @ B_child, Epi_child)
vars_Bst['Epi_IPCHI2_WRT_LbBPV'] = F.IPCHI2.bind(F.BPV(v2_pvs) @ B_child, Epi_child)
#DOCA and DOCACHI2 b/w Lb and extra particle
vars_Bst['DOCA12'] = F.DOCA(1, 2)
vars_Bst['DOCA12_CHI2_12'] = F.DOCACHI2(1, 2)
#DOCACHI2 of extra particle wrt to mother i.e. Sb
vars_Bst['MTDOCACHI2_1'] = F.MTDOCACHI2(1, v2_pvs)
vars_Bst['MTDOCACHI2_2'] = F.MTDOCACHI2(2, v2_pvs)
#define fields
fields_mu = {}
fields_mu['Sb'] = f'[B*0 -> (B- -> {hadron_name} mu-) [pi+]CC]CC'
fields_mu['Lb'] = f'[B*0 -> ^(B- -> {hadron_name} mu-) [pi+]CC]CC'
fields_mu['Lc'] = f'[B*0 -> (B- -> ^{hadron_name} mu-) [pi+]CC]CC'
fields_mu['Lep'] = f'[B*0 -> (B- -> {hadron_name} ^mu-) [pi+]CC]CC'
fields_mu['Epi'] = f'[B*0 -> (B- -> {hadron_name} mu-) ^[pi+]CC]CC'
#add variables
variables = {}
variables["ALL"] = vars_common
variables["Sb"] = vars_composite + vars_Bst
variables["Lb"] = vars_composite + var_B
variables["Lc"] = vars_composite
variables["Lep"] = vars_basic
variables["Epi"] = vars_basic
#define tuple_mu
my_tuple_mu = Funtuple(
name="Tuple_mu",
tuple_name="DecayTree",
fields=fields_mu,
variables=variables,
event_variables=evt_vars,
inputs=Bst)
if get_electron:
#field for electron mode
fields_e = {}
fields_e['Sb'] = fields_mu['Sb'].replace('mu', 'e')
fields_e['Lb'] = fields_mu['Lb'].replace('mu', 'e')
fields_e['Lc'] = fields_mu['Lc'].replace('mu', 'e')
fields_e['Lep'] = fields_mu['Lep'].replace('mu', 'e')
fields_e['Epi'] = fields_mu['Epi'].replace('mu', 'e')
#define tuple_e
my_tuple_e = Funtuple(
name="Tuple_e",
tuple_name="DecayTree",
fields=fields_e,
variables=variables,
event_variables=evt_vars,
inputs=Bst)
return my_tuple_mu, my_tuple_e
return my_tuple_mu
def main(options : Options):
#line name
SUFFIX_NAME = '_BToDzl'
line_name = f'Hlt2SLB{SUFFIX_NAME}'
hadron_candidate_id = 421 #B -> Hadron mu, Here Hadron is D0 so absid is 421
hadron_name = 'D0'
store_electron = True
#define filer
my_filter = create_lines_filter(name="Filter", lines=[line_name])
#get data and extra particles
lb = get_particles(f"/Event/HLT2/{line_name}/Particles")
extra_particles = get_particles(f"/Event/HLT2/{line_name}/ExtraPions/Particles")
#get v2_pvs and rec_summary
v2_pvs = get_pvs()
rec_summary = get_rec_summary()
#make sb candidate
Bst = make_Bst(lb, extra_particles, hadron_candidate_id)
MCTRUTH_Bst = MCTruthAndBkgCat(Bst, name='MCTRUTH_Bst')
tuple_file = tuple_Bst(Bst, MCTRUTH_Bst, v2_pvs, rec_summary, hadron_name, get_electron=store_electron)
#define algorithms
user_algorithms = {}
user_algorithms['Alg_mu'] = [my_filter, tuple_file[0]]
if store_electron:
user_algorithms['Alg_e'] = [my_filter, tuple_file[1]]
return make_config(options, user_algorithms)
Loading