14 import modularAnalysis
as ma
15 from variables
import variables
18 from ROOT
import Belle2
20 _TrainingMode = Belle2.ChargedPidMVAWeights.ChargedPidMVATrainingMode
24 _chargednames = [
'pi',
'K',
'p',
'e',
'mu']
25 _pidnames = [
'pionID',
'kaonID',
'protonID',
'electronID',
'muonID']
26 _stdnames = [
'all',
'loose',
'loosepid',
'good',
'higheff']
27 _effnames = [
'95eff',
'90eff',
'85eff']
30 _mostLikelyList =
'mostlikely'
33 def _stdChargedEffCuts(particletype, listtype):
35 Provides the PID cut corresponding to a given efficiency percentile
37 @param particletype type of charged particle (pi, K, p, e, mu)
38 @param listtype efficiency percentile for the list (95eff, 90eff, 85eff)
41 particleindex = _chargednames.index(particletype)
42 effindex = _effnames.index(listtype)
45 effcuts = [[0.001, 0.019, 0.098],
47 [0.000, 0.043, 0.251],
48 [0.093, 0.301, 0.709],
49 [0.187, 0.418, 0.909]]
51 return effcuts[particleindex][effindex]
54 def stdCharged(particletype, listtype, path, writeOut=True):
56 Function to prepare one of several standardized types of charged particle lists:
57 - 'all' with no cuts on track
58 - 'good' high purity lists for data studies
59 - 'loosepid' loose selections for skimming, PID cut only
60 - 'loose' loose selections for skimming
61 - 'higheff' high efficiency list with loose global ID cut for data studies
62 - 'mostlikely' list with the highest PID likelihood
63 Also the following lists, which may or may not be available depending on the release
64 - '99eff' with 99% selection efficiency (calculated for 1<p<4 GeV) and good track (MC only)
65 - '95eff' with 95% selection efficiency (calculated for 1<p<4 GeV) and good track (MC only)
66 - '90eff' with 90% selection efficiency (calculated for 1<p<4 GeV) and good track (MC only)
67 - '85eff' with 85% selection efficiency (calculated for 1<p<4 GeV) and good track (MC only)
69 @param particletype type of charged particle to make a list of
70 @param listtype name of standard list
71 @param path modules are added to this path
72 @param writeOut whether RootOutput module should save the created ParticleList
76 trackQuality =
'thetaInCDCAcceptance and nCDCHits>20'
77 ipCut =
'dr < 0.5 and abs(dz) < 2'
78 goodTrack = trackQuality +
' and ' + ipCut
80 if particletype
not in _chargednames:
81 b2.B2ERROR(
"The requested list is not a standard charged particle. Use one of pi, K, e, mu, p.")
84 ma.fillParticleList(particletype +
'+:all',
'', writeOut=writeOut, path=path)
85 elif listtype ==
'good':
87 particletype +
'+:good',
88 _pidnames[_chargednames.index(particletype)] +
' > 0.5 and ' + goodTrack,
91 elif listtype ==
'loose':
93 particletype +
'+:loose',
94 _pidnames[_chargednames.index(particletype)] +
' > 0.1 and ' + goodTrack,
97 elif listtype ==
'loosepid':
99 particletype +
'+:loosepid',
100 _pidnames[_chargednames.index(particletype)] +
' > 0.1',
103 elif listtype ==
'higheff':
105 particletype +
'+:higheff',
106 _pidnames[_chargednames.index(particletype)] +
' > 0.002 and ' + goodTrack,
109 elif listtype
not in _effnames:
110 b2.B2ERROR(
"The requested list is not defined. Please refer to the stdCharged documentation.")
112 pidcut = _stdChargedEffCuts(particletype, listtype)
113 if 0.0 < pidcut < 1.0:
118 _pidnames[_chargednames.index(particletype)] +
126 b2.B2ERROR(
'The requested standard particle list ' + particletype +
127 '+:' + listtype +
' is not available in this release.')
130 def stdPi(listtype=_defaultlist, path=None, writeOut=True):
132 Function to prepare standard pion lists, refer to `stdCharged` for details
134 @param listtype name of standard list
135 @param path modules are added to this path
136 @param writeOut whether RootOutput module should save the created ParticleList
141 def stdK(listtype=_defaultlist, path=None, writeOut=True):
143 Function to prepare standard kaon lists, refer to `stdCharged` for details
145 @param listtype name of standard list
146 @param path modules are added to this path
147 @param writeOut whether RootOutput module should save the created ParticleList
152 def stdPr(listtype=_defaultlist, path=None, writeOut=True):
154 Function to prepare standard proton lists, refer to `stdCharged` for details
156 @param listtype name of standard list
157 @param path modules are added to this path
158 @param writeOut whether RootOutput module should save the created ParticleList
169 channel_eff="combination",
170 channel_misid_pi="combination",
171 channel_misid_K="combination",
173 outputListLabel=None,
174 trainingModeMulticlass=_TrainingMode.c_Multiclass,
175 trainingModeBinary=_TrainingMode.c_Classification,
178 Function to prepare one of several standardized types of lepton (:math:`e,\\mu`) lists, with the following working points:
180 * 'FixedThresh05', PID cut of > 0.5 for each particle in the list.
181 * 'FixedThresh09', PID cut of > 0.9 for each particle in the list.
182 * 'FixedThresh095', PID cut of > 0.95 for each particle in the list.
183 * 'FixedThresh099', PID cut of > 0.99 for each particle in the list.
184 * 'UniformEff60' 60% lepton efficiency list, uniform in a given multi-dimensional parametrisation.
185 * 'UniformEff70' 70% lepton efficiency list, uniform in a given multi-dimensional parametrisation.
186 * 'UniformEff80' 80% lepton efficiency list, uniform in a given multi-dimensional parametrisation.
187 * 'UniformEff90' 90% lepton efficiency list, uniform in a given multi-dimensional parametrisation.
188 * 'UniformEff95' 95% lepton efficiency list, uniform in a given multi-dimensional parametrisation.
189 * 'UniformPiFR5EM1' 50% pion to lepton fake rate, uniform in a given multi-dimensional parametrisation.
190 * 'UniformPiFR1EM1' 10% pion to lepton fake rate, uniform in a given multi-dimensional parametrisation.
191 * 'UniformPiFR5EM2' 5% pion to lepton fake rate, uniform in a given multi-dimensional parametrisation.
192 * 'UniformPiFR1EM2' 1% pion to lepton fake rate, uniform in a given multi-dimensional parametrisation.
193 * 'UniformPiFR5EM3' 0.5% pion to lepton fake rate, uniform in a given multi-dimensional parametrisation.
194 * 'UniformPiFR1EM3' 0.1% pion to lepton fake rate, uniform in a given multi-dimensional parametrisation.
196 The function creates a ``ParticleList``, selecting particles according to the chosen ``working_point``,
197 and decorates each candidate in the list with the nominal Data/MC :math:`\\ell` ID efficiency and
198 :math:`\\pi,K` fake rate correction factors and their stat, syst uncertainty, reading the info
199 from the Conditions Database (CDB) according to the chosen input global tag (GT).
202 Particles will **not** be selected if they are outside the Data/MC *efficiency* corrections' phase space coverage
203 for the given working point.
204 In fact, the threshold value for the PID cut in such cases is set to NaN.
207 At the moment, the only supported *binary* lepton identification is against the **pion** hypothesis.
208 This implies that, if binary classification is chosen, only :math:`\\pi` fake rate corrections will be added to the
209 resulting particle list.
212 pdgId (int): the lepton pdg code.
213 working_point (str): name of the chosen working point that defines the content of the list. Choose among the above values.
214 method (str): the PID method: 'likelihood' or 'bdt'.
215 classification (str): the type of classifier: 'binary' (one-vs-pion) or 'global' (one-vs-all).
216 lid_weights_gt (str): the name identifier of the global tag with the recommended Data/MC correction weights.
220 `Lepton ID Confluence page <https://confluence.desy.de/display/BI/Lepton+ID+Performance>`_
221 for info about the recommended global tags.
223 release (Optional[int]): the major release number of the data and MC campaigns considered.
224 If specified, it ensures the correct :math:`\\ell` ID variables are used.
228 `Lepton ID Confluence page <https://confluence.desy.de/display/BI/Lepton+ID+Performance>`_
229 for info about lepton identification variables and campaigns.
231 channel_eff (Optional[str]): the channel used to derive the :math:`\\ell` ID efficiency corrections.
232 By default, 'combination' is set, meaning they are obtained by combining results
233 of several hadronic and low multiplicity channels, wherever they overlap.
237 `Lepton ID Confluence page <https://confluence.desy.de/display/BI/Lepton+ID+Performance>`_
238 for other possible choices (if any).
240 channel_misid_pi (Optional[str]): the channel used to derive the :math:`\\pi` fake rate corrections.
241 channel_misid_K (Optional[str]): the channel used to derive the :math:`K` fake rate corrections.
242 inputListName (Optional[str]): the name of a pre-existing ``ParticleList`` object (defined as a full ``decayString``,
243 e.g. 'e-:my_input_electrons') of which the standard lepton list will be a subset.
244 For instance, users might want to apply a Bremsstrahlung correction to electrons first,
245 which modifies their 4-momentum, and only later define the subset passing the PID selection,
246 including the appropriate PID weights and uncertainties (which are :math:`p`-dependent).
247 By default, the standard lepton list is created from all ``Track`` objects in the event.
250 Do **not** apply any PID selection on the input list, otherwise results could be biased.
252 outputListLabel (Optional[str]): the name of the output lepton list label, i.e.,
253 the string that follows the particle identifier ('e-:', 'mu-:').
254 By default, it is assigned as:
255 ``'{method}_{classification}_{working_point}'``.
257 trainingModeMulticlass (Optional[``Belle2.ChargedPidMVAWeights.ChargedPidMVATrainingMode``]): enum identifier
258 of the multi-class (global PID) training mode.
259 See `modularAnalysis.applyChargedPidMVA` docs for available options.
260 trainingModeBinary (Optional[``Belle2.ChargedPidMVAWeights.ChargedPidMVATrainingMode``]): enum identifier
261 of the classification (binary PID) training mode.
262 See `modularAnalysis.applyChargedPidMVA` docs for available options.
263 path (basf2.Path): modules are added to this path.
266 tuple(str, list(str)): the alias for the lepton ID variable, and the list of aliases for the weights.
287 available_methods = (
"likelihood",
"bdt")
288 available_classificators = (
"global",
"binary")
290 if working_point
not in working_points:
291 b2.B2ERROR(f
"The requested lepton list working point: {working_point} is not defined. \
292 Please refer to the stdLep and stdCharged documentation.")
295 if method
not in available_methods:
296 b2.B2ERROR(f
"method: {method}. Must be any of: {available_methods}.")
299 if classification
not in available_classificators:
300 b2.B2ERROR(f
"classification: {classification}. Must be any of: {available_classificators}.")
303 b2.B2INFO(f
"Prepending GT with LID corrections: {lid_weights_gt}")
304 b2.conditions.prepend_globaltag(lid_weights_gt)
310 electron = Const.electron.getPDGCode()
311 muon = Const.muon.getPDGCode()
312 pion = Const.pion.getPDGCode()
314 if lepton
not in (electron, muon):
315 b2.B2ERROR(f
"{pdgId} is not that of a light charged lepton.")
323 electron:
"electronID",
327 electron:
"electronID",
333 electron: f
"binaryPID({electron}, {pion})",
334 muon: f
"binaryPID({muon}, {pion})",
337 electron:
"binaryPID_e_pi",
338 muon:
"binaryPID_mu_pi",
345 lepton: f
"pidChargedBDTScore({lepton}, ALL)",
348 lepton: re.sub(
r"\W+",
"", f
"pidChargedBDTScore_{lepton_name}"),
353 lepton: f
"pidPairChargedBDTScore({lepton}, {pion}, ALL)",
356 lepton: re.sub(
r"\W+",
"", f
"pidPairChargedBDTScore_{lepton_name}_pi"),
364 if int(release)
in [5, 6]:
365 if lepton == electron:
366 b2.B2INFO(f
"The likelihood-based electron ID in release {release} samples is defined w/o the SVD and the TOP")
367 pid_variables[
"likelihood"][
"global"][
"var"][electron] =
"electronID_noSVD_noTOP"
368 pid_variables[
"likelihood"][
"global"][
"alias"][electron] =
"electronID_noSVD_noTOP"
369 pid_variables[
"likelihood"][
"binary"][
"var"][electron] = f
"binaryElectronID_noSVD_noTOP({pion})"
370 pid_variables[
"likelihood"][
"binary"][
"alias"][electron] =
"binaryElectronID_noSVD_noTOP_pi"
372 b2.B2INFO(f
"The likelihood-based muon ID in release {release} samples is defined w/o the SVD")
373 pid_variables[
"likelihood"][
"global"][
"var"][muon] =
"muonID_noSVD"
374 pid_variables[
"likelihood"][
"global"][
"alias"][muon] =
"muonID_noSVD"
375 pid_variables[
"likelihood"][
"binary"][
"var"][muon] = f
"binaryPID_noSVD({muon}, {pion})"
376 pid_variables[
"likelihood"][
"binary"][
"alias"][muon] =
"binaryMuonID_noSVD_pi"
379 pid_var = pid_variables[method][classification][
"var"][lepton]
380 pid_alias = pid_variables[method][classification][
"alias"][lepton]
381 if pid_alias != pid_var:
382 variables.addAlias(pid_alias, pid_var)
385 outputListName = f
"{lepton_name}:{method}_{classification}_{working_point}"
386 if outputListLabel
is not None:
387 outputListName = f
"{lepton_name}:{outputListLabel}"
389 if inputListName
is None:
390 ma.fillParticleList(outputListName,
"", path=path)
393 f
"The standard lepton list: '{outputListName}' will be created as a subset \
394 of the following ParticleList: '{inputListName}'")
395 ma.copyList(outputListName, inputListName, path=path)
399 if classification ==
"global":
400 ma.applyChargedPidMVA(particleLists=[outputListName],
402 trainingMode=trainingModeMulticlass)
403 elif classification ==
"binary":
404 ma.applyChargedPidMVA(particleLists=[outputListName],
406 binaryHypoPDGCodes=(lepton, pion),
407 trainingMode=trainingModeBinary)
410 payload_eff = f
"ParticleReweighting:{pid_alias}_eff_{channel_eff}_{working_point}"
411 payload_misid_pi = f
"ParticleReweighting:{pid_alias}_misid_pi_{channel_misid_pi}_{working_point}"
414 path.add_module(
"ParticleWeighting",
415 particleList=outputListName,
416 tableName=payload_eff).set_name(f
"ParticleWeighting_eff_{outputListName}")
417 path.add_module(
"ParticleWeighting",
418 particleList=outputListName,
419 tableName=payload_misid_pi,
420 allowToSkip=
True).set_name(f
"ParticleWeighting_misid_pi_{outputListName}")
422 if classification ==
"global":
423 payload_misid_K = f
"ParticleReweighting:{pid_alias}_misid_K_{channel_misid_K}_{working_point}"
424 path.add_module(
"ParticleWeighting",
425 particleList=outputListName,
426 tableName=payload_misid_K,
427 allowToSkip=
True).set_name(f
"ParticleWeighting_misid_K_{outputListName}")
431 cut = f
"[{pid_alias} >= extraInfo({payload_eff}_threshold)]"
432 ma.applyCuts(outputListName, cut, path=path)
435 weight_aliases_to_var = {
436 f
"weight_{pid_alias}_eff_{working_point}": f
"extraInfo({payload_eff}_data_MC_ratio)",
437 f
"weight_{pid_alias}_misid_pi_{working_point}": f
"extraInfo({payload_misid_pi}_data_MC_ratio)",
439 f
"weight_{pid_alias}_eff_{working_point}_stat_up": f
"extraInfo({payload_eff}_data_MC_uncertainty_stat_up)",
440 f
"weight_{pid_alias}_eff_{working_point}_stat_dn": f
"extraInfo({payload_eff}_data_MC_uncertainty_stat_dn)",
441 f
"weight_{pid_alias}_eff_{working_point}_sys_up": f
"extraInfo({payload_eff}_data_MC_uncertainty_sys_up)",
442 f
"weight_{pid_alias}_eff_{working_point}_sys_dn": f
"extraInfo({payload_eff}_data_MC_uncertainty_sys_dn)",
443 f
"weight_{pid_alias}_misid_pi_{working_point}_stat_up": f
"extraInfo({payload_misid_pi}_data_MC_uncertainty_stat_up)",
444 f
"weight_{pid_alias}_misid_pi_{working_point}_stat_dn": f
"extraInfo({payload_misid_pi}_data_MC_uncertainty_stat_dn)",
445 f
"weight_{pid_alias}_misid_pi_{working_point}_sys_up": f
"extraInfo({payload_misid_pi}_data_MC_uncertainty_sys_up)",
446 f
"weight_{pid_alias}_misid_pi_{working_point}_sys_dn": f
"extraInfo({payload_misid_pi}_data_MC_uncertainty_sys_dn)",
455 f
"weight_{pid_alias}_eff_{working_point}_rel_stat_up":
456 f
"formula(1+[extraInfo({payload_eff}_data_MC_uncertainty_stat_up)/extraInfo({payload_eff}_data_MC_ratio)])",
457 f
"weight_{pid_alias}_eff_{working_point}_rel_stat_dn":
458 f
"formula(1-[extraInfo({payload_eff}_data_MC_uncertainty_stat_dn)/extraInfo({payload_eff}_data_MC_ratio)])",
459 f
"weight_{pid_alias}_eff_{working_point}_rel_sys_up":
460 f
"formula(1+[extraInfo({payload_eff}_data_MC_uncertainty_sys_up)/extraInfo({payload_eff}_data_MC_ratio)])",
461 f
"weight_{pid_alias}_eff_{working_point}_rel_sys_dn":
462 f
"formula(1-[extraInfo({payload_eff}_data_MC_uncertainty_sys_dn)/extraInfo({payload_eff}_data_MC_ratio)])",
463 f
"weight_{pid_alias}_misid_pi_{working_point}_rel_stat_up":
464 f
"formula(1+[extraInfo({payload_misid_pi}_data_MC_uncertainty_stat_up)/extraInfo({payload_misid_pi}_data_MC_ratio)])",
465 f
"weight_{pid_alias}_misid_pi_{working_point}_rel_stat_dn":
466 f
"formula(1-[extraInfo({payload_misid_pi}_data_MC_uncertainty_stat_dn)/extraInfo({payload_misid_pi}_data_MC_ratio)])",
467 f
"weight_{pid_alias}_misid_pi_{working_point}_rel_sys_up":
468 f
"formula(1+[extraInfo({payload_misid_pi}_data_MC_uncertainty_sys_up)/extraInfo({payload_misid_pi}_data_MC_ratio)])",
469 f
"weight_{pid_alias}_misid_pi_{working_point}_rel_sys_dn":
470 f
"formula(1-[extraInfo({payload_misid_pi}_data_MC_uncertainty_sys_dn)/extraInfo({payload_misid_pi}_data_MC_ratio)])",
472 if classification ==
"global":
473 weight_aliases_to_var.update({
474 f
"weight_{pid_alias}_misid_K_{working_point}": f
"extraInfo({payload_misid_K}_data_MC_ratio)",
476 f
"weight_{pid_alias}_misid_K_{working_point}_stat_up": f
"extraInfo({payload_misid_K}_data_MC_uncertainty_stat_up)",
477 f
"weight_{pid_alias}_misid_K_{working_point}_stat_dn": f
"extraInfo({payload_misid_K}_data_MC_uncertainty_stat_dn)",
478 f
"weight_{pid_alias}_misid_K_{working_point}_sys_up": f
"extraInfo({payload_misid_K}_data_MC_uncertainty_sys_up)",
479 f
"weight_{pid_alias}_misid_K_{working_point}_sys_dn": f
"extraInfo({payload_misid_K}_data_MC_uncertainty_sys_dn)",
481 f
"weight_{pid_alias}_misid_K_{working_point}_rel_stat_up":
482 f
"formula(1+[extraInfo({payload_misid_K}_data_MC_uncertainty_stat_up)/extraInfo({payload_misid_K}_data_MC_ratio)])",
483 f
"weight_{pid_alias}_misid_K_{working_point}_rel_stat_dn":
484 f
"formula(1-[extraInfo({payload_misid_K}_data_MC_uncertainty_stat_dn)/extraInfo({payload_misid_K}_data_MC_ratio)])",
485 f
"weight_{pid_alias}_misid_K_{working_point}_rel_sys_up":
486 f
"formula(1+[extraInfo({payload_misid_K}_data_MC_uncertainty_sys_up)/extraInfo({payload_misid_K}_data_MC_ratio)])",
487 f
"weight_{pid_alias}_misid_K_{working_point}_rel_sys_dn":
488 f
"formula(1-[extraInfo({payload_misid_K}_data_MC_uncertainty_sys_dn)/extraInfo({payload_misid_K}_data_MC_ratio)])",
491 for alias, var
in weight_aliases_to_var.items():
492 variables.addAlias(alias, var)
494 return pid_alias, list(weight_aliases_to_var.keys())
497 def stdE(listtype=_defaultlist,
502 channel_eff="combination",
503 channel_misid_pi="combination",
504 channel_misid_K="combination",
506 outputListLabel=None,
507 trainingModeMulticlass=_TrainingMode.c_Multiclass,
508 trainingModeBinary=_TrainingMode.c_Classification,
510 """ Function to prepare one of several standardized types of electron lists.
511 See the documentation of `stdLep` for details.
513 It also accepts any of the legacy definitions
514 for the ``listtype`` parameter (aka ``working_point`` in `stdLep`) to fall back to the `stdCharged` behaviour:
526 tuple(str, list(str)): the alias for the electron ID variable, and the list of aliases for the weights.
529 if listtype
in _stdnames + _effnames:
531 return _pidnames[_chargednames.index(
"e")],
None
533 return stdLep(Const.electron.getPDGCode(), listtype, method, classification, lid_weights_gt,
535 channel_eff=channel_eff,
536 channel_misid_pi=channel_misid_pi,
537 channel_misid_K=channel_misid_K,
538 inputListName=inputListName,
539 outputListLabel=outputListLabel,
540 trainingModeMulticlass=trainingModeMulticlass,
541 trainingModeBinary=trainingModeBinary,
545 def stdMu(listtype=_defaultlist,
550 channel_eff="combination",
551 channel_misid_pi="combination",
552 channel_misid_K="combination",
554 outputListLabel=None,
555 trainingModeMulticlass=_TrainingMode.c_Multiclass,
556 trainingModeBinary=_TrainingMode.c_Classification,
558 """ Function to prepare one of several standardized types of muon lists.
559 See the documentation of `stdLep` for details.
561 It also accepts any of the legacy definitions
562 for the ``listtype`` parameter (aka ``working_point`` in `stdLep`) to fall back to the `stdCharged` behaviour:
574 tuple(str, list(str)): the alias for the muon ID variable, and the list of aliases for the weights.
577 if listtype
in _stdnames + _effnames:
579 return _pidnames[_chargednames.index(
"mu")],
None
581 return stdLep(Const.muon.getPDGCode(), listtype, method, classification, lid_weights_gt,
583 channel_eff=channel_eff,
584 channel_misid_pi=channel_misid_pi,
585 channel_misid_K=channel_misid_K,
586 inputListName=inputListName,
587 outputListLabel=outputListLabel,
588 trainingModeMulticlass=trainingModeMulticlass,
589 trainingModeBinary=trainingModeBinary,
593 def stdMostLikely(pidPriors=None, suffix='', custom_cuts='', path=None, writeOut=True):
595 Function to prepare most likely particle lists according to PID likelihood, refer to stdCharged for details
597 @param pidPriors list of 6 float numbers used to reweight PID likelihoods, for e, mu, pi, K, p and d
598 @param suffix string added to the end of particle list names
599 @param custom_cuts custom selection cut string, if empty, standard track quality cuts will be applied
600 @param path modules are added to this path
601 @param writeOut whether RootOutput module should save the created ParticleList
607 if pidPriors
is not None:
608 args = str(pidPriors)[1:-1]
609 trackQuality =
'thetaInCDCAcceptance and nCDCHits>20'
610 if custom_cuts !=
'':
611 trackQuality = custom_cuts
612 for name
in _chargednames:
613 ma.fillParticleList(f
'{name}+:{_mostLikelyList}{suffix}',
614 f
'pidIsMostLikely({args}) > 0 and {trackQuality}', writeOut=writeOut, path=path)
This class provides a set of constants for the framework.