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]
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
75 trackQuality =
'thetaInCDCAcceptance and nCDCHits>20'
76 ipCut =
'dr < 0.5 and abs(dz) < 2'
77 goodTrack = trackQuality +
' and ' + ipCut
79 if particletype
not in _chargednames:
80 b2.B2ERROR(
"The requested list is not a standard charged particle. Use one of pi, K, e, mu, p.")
83 ma.fillParticleList(particletype +
'+:all',
'',
True, path=path)
84 elif listtype ==
'good':
86 particletype +
'+:good',
87 _pidnames[_chargednames.index(particletype)] +
' > 0.5 and ' + goodTrack,
90 elif listtype ==
'loose':
92 particletype +
'+:loose',
93 _pidnames[_chargednames.index(particletype)] +
' > 0.1 and ' + goodTrack,
96 elif listtype ==
'loosepid':
98 particletype +
'+:loosepid',
99 _pidnames[_chargednames.index(particletype)] +
' > 0.1',
102 elif listtype ==
'higheff':
104 particletype +
'+:higheff',
105 _pidnames[_chargednames.index(particletype)] +
' > 0.002 and ' + goodTrack,
108 elif listtype
not in _effnames:
109 b2.B2ERROR(
"The requested list is not defined. Please refer to the stdCharged documentation.")
111 pidcut = _stdChargedEffCuts(particletype, listtype)
112 if 0.0 < pidcut < 1.0:
117 _pidnames[_chargednames.index(particletype)] +
125 b2.B2ERROR(
'The requested standard particle list ' + particletype +
126 '+:' + listtype +
' is not available in this release.')
129 def stdPi(listtype=_defaultlist, path=None):
131 Function to prepare standard pion lists, refer to `stdCharged` for details
133 @param listtype name of standard list
134 @param path modules are added to this path
139 def stdK(listtype=_defaultlist, path=None):
141 Function to prepare standard kaon lists, refer to `stdCharged` for details
143 @param listtype name of standard list
144 @param path modules are added to this path
149 def stdPr(listtype=_defaultlist, path=None):
151 Function to prepare standard proton lists, refer to `stdCharged` for details
153 @param listtype name of standard list
154 @param path modules are added to this path
165 channel_eff="combination",
166 channel_misid_pi="combination",
167 channel_misid_K="combination",
169 outputListLabel=None,
170 trainingModeMulticlass=_TrainingMode.c_Multiclass,
171 trainingModeBinary=_TrainingMode.c_Classification,
174 Function to prepare one of several standardized types of lepton (:math:`e,\\mu`) lists, with the following working points:
176 * 'FixedThresh05', PID cut of > 0.5 for each particle in the list.
177 * 'FixedThresh09', PID cut of > 0.9 for each particle in the list.
178 * 'FixedThresh095', PID cut of > 0.95 for each particle in the list.
179 * 'UniformEff60' 60% lepton efficiency list, uniform in a given multi-dimensional parametrisation.
180 * 'UniformEff70' 70% lepton efficiency list, uniform in a given multi-dimensional parametrisation.
181 * 'UniformEff80' 80% lepton efficiency list, uniform in a given multi-dimensional parametrisation.
182 * 'UniformEff90' 90% lepton efficiency list, uniform in a given multi-dimensional parametrisation.
183 * 'UniformEff95' 95% lepton efficiency list, uniform in a given multi-dimensional parametrisation.
185 The function creates a ``ParticleList``, selecting particles according to the chosen ``working_point``,
186 and decorates each candidate in the list with the nominal Data/MC :math:`\\ell` ID efficiency and
187 :math:`\\pi,K` fake rate correction factors and their stat, syst uncertainty, reading the info
188 from the Conditions Database (CDB) according to the chosen input global tag (GT).
191 Particles will **not** be selected if they are outside the Data/MC *efficiency* corrections' phase space coverage
192 for the given working point.
193 In fact, the threshold value for the PID cut in such cases is set to NaN.
196 At the moment, the only supported *binary* lepton identification is against the **pion** hypothesis.
197 This implies that, if binary classification is chosen, only :math:`\\pi` fake rate corrections will be added to the
198 resulting particle list.
201 pdgId (int): the lepton pdg code.
202 working_point (str): name of the chosen working point that defines the content of the list. Choose among the above values.
203 method (str): the PID method: 'likelihood' or 'bdt'.
204 classification (str): the type of classifier: 'binary' (one-vs-pion) or 'global' (one-vs-all).
205 lid_weights_gt (str): the name identifier of the global tag with the recommended Data/MC correction weights.
209 `Lepton ID Confluence page <https://confluence.desy.de/display/BI/Lepton+ID+Performance>`_
210 for info about the recommended global tags.
212 release (Optional[int]): the major release number of the data and MC campaigns considered.
213 If specified, it ensures the correct :math:`\\ell` ID variables are used.
217 `Lepton ID Confluence page <https://confluence.desy.de/display/BI/Lepton+ID+Performance>`_
218 for info about lepton identification variables and campaigns.
220 channel_eff (Optional[str]): the channel used to derive the :math:`\\ell` ID efficiency corrections.
221 By default, 'combination' is set, meaning they are obtained by combining results
222 of several hadronic and low multiplicity channels, wherever they overlap.
226 `Lepton ID Confluence page <https://confluence.desy.de/display/BI/Lepton+ID+Performance>`_
227 for other possible choices (if any).
229 channel_misid_pi (Optional[str]): the channel used to derive the :math:`\\pi` fake rate corrections.
230 channel_misid_K (Optional[str]): the channel used to derive the :math:`K` fake rate corrections.
231 inputListName (Optional[str]): the name of a pre-existing ``ParticleList`` object (defined as a full ``decayString``,
232 e.g. 'e-:my_input_electrons') of which the standard lepton list will be a subset.
233 For instance, users might want to apply a Bremsstrahlung correction to electrons first,
234 which modifies their 4-momentum, and only later define the subset passing the PID selection,
235 including the appropriate PID weights and uncertainties (which are :math:`p`-dependent).
236 By default, the standard lepton list is created from all ``Track`` objects in the event.
239 Do **not** apply any PID selection on the input list, otherwise results could be biased.
241 outputListLabel (Optional[str]): the name of the output lepton list label, i.e.,
242 the string that follows the particle identifier ('e-:', 'mu-:').
243 By default, it is assigned as:
244 ``'{method}_{classification}_{working_point}'``.
246 trainingModeMulticlass (Optional[``Belle2.ChargedPidMVAWeights.ChargedPidMVATrainingMode``]): enum identifier
247 of the multi-class (global PID) training mode.
248 See `modularAnalysis.applyChargedPidMVA` docs for available options.
249 trainingModeBinary (Optional[``Belle2.ChargedPidMVAWeights.ChargedPidMVATrainingMode``]): enum identifier
250 of the classification (binary PID) training mode.
251 See `modularAnalysis.applyChargedPidMVA` docs for available options.
252 path (basf2.Path): modules are added to this path.
255 tuple(str, list(str)): the alias for the lepton ID variable, and the list of aliases for the weights.
269 available_methods = (
"likelihood",
"bdt")
270 available_classificators = (
"global",
"binary")
272 if working_point
not in working_points:
273 b2.B2ERROR(
"The requested lepton list working point is not defined. \
274 Please refer to the stdLep and stdCharged documentation.")
277 if method
not in available_methods:
278 b2.B2ERROR(f
"method: {method}. Must be any of: {available_methods}.")
281 if classification
not in available_classificators:
282 b2.B2ERROR(f
"classification: {classification}. Must be any of: {available_classificators}.")
285 b2.B2INFO(f
"Prepending GT with LID corrections: {lid_weights_gt}")
286 b2.conditions.prepend_globaltag(lid_weights_gt)
292 electron = Const.electron.getPDGCode()
293 muon = Const.muon.getPDGCode()
294 pion = Const.pion.getPDGCode()
296 if lepton
not in (electron, muon):
297 b2.B2ERROR(f
"{pdgId} is not that of a light charged lepton.")
305 electron:
"electronID",
309 electron:
"electronID",
315 electron: f
"binaryPID({electron}, {pion})",
316 muon: f
"binaryPID({muon}, {pion})",
319 electron:
"binaryPID_e_pi",
320 muon:
"binaryPID_mu_pi",
327 lepton: f
"pidChargedBDTScore({lepton}, ALL)",
330 lepton: re.sub(
r"\W+",
"", f
"pidChargedBDTScore_{lepton_name}"),
335 lepton: f
"pidPairChargedBDTScore({lepton}, {pion}, ALL)",
338 lepton: re.sub(
r"\W+",
"", f
"pidPairChargedBDTScore_{lepton_name}_pi"),
346 if int(release) == 5:
347 if lepton == electron:
348 b2.B2INFO(
"The likelihood-based electron ID in release 5 samples is defined w/o the SVD and the TOP")
349 pid_variables[
"likelihood"][
"global"][
"var"][electron] =
"electronID_noSVD_noTOP"
350 pid_variables[
"likelihood"][
"global"][
"alias"][electron] =
"electronID_noSVD_noTOP"
351 pid_variables[
"likelihood"][
"binary"][
"var"][electron] = f
"binaryElectronID_noSVD_noTOP({pion})"
352 pid_variables[
"likelihood"][
"binary"][
"alias"][electron] =
"binaryElectronID_noSVD_noTOP_pi"
354 b2.B2INFO(
"The likelihood-based muon ID in release 5 samples is defined w/o the SVD")
355 pid_variables[
"likelihood"][
"global"][
"var"][muon] =
"muonID_noSVD"
356 pid_variables[
"likelihood"][
"global"][
"alias"][muon] =
"muonID_noSVD"
357 pid_variables[
"likelihood"][
"binary"][
"var"][muon] = f
"binaryPID_noSVD({muon}, {pion})"
358 pid_variables[
"likelihood"][
"binary"][
"alias"][muon] =
"binaryMuonID_noSVD_pi"
359 if int(release) == 6:
360 if lepton == electron:
361 b2.B2INFO(
"The likelihood-based electron ID in release 6 samples is defined w/o the TOP")
362 pid_variables[
"likelihood"][
"global"][
"var"][electron] =
"electronID_noTOP"
363 pid_variables[
"likelihood"][
"global"][
"alias"][electron] =
"electronID_noTOP"
364 pid_variables[
"likelihood"][
"binary"][
"var"][electron] = f
"binaryElectronID_noTOP({pion})"
365 pid_variables[
"likelihood"][
"binary"][
"alias"][electron] =
"binaryElectronID_noTOP_pi"
368 pid_var = pid_variables[method][classification][
"var"][lepton]
369 pid_alias = pid_variables[method][classification][
"alias"][lepton]
370 if pid_alias != pid_var:
371 variables.addAlias(pid_alias, pid_var)
374 outputListName = f
"{lepton_name}:{method}_{classification}_{working_point}"
375 if outputListLabel
is not None:
376 outputListName = f
"{lepton_name}:{outputListLabel}"
378 if inputListName
is None:
379 ma.fillParticleList(outputListName,
"", path=path)
382 f
"The standard lepton list: '{outputListName}' will be created as a subset \
383 of the following ParticleList: '{inputListName}'")
384 ma.copyList(outputListName, inputListName, path=path)
388 if classification ==
"global":
389 ma.applyChargedPidMVA(particleLists=[outputListName],
391 trainingMode=trainingModeMulticlass)
392 elif classification ==
"binary":
393 ma.applyChargedPidMVA(particleLists=[outputListName],
395 binaryHypoPDGCodes=(lepton, pion),
396 trainingMode=trainingModeBinary)
399 payload_eff = f
"ParticleReweighting:{pid_alias}_eff_{channel_eff}_{working_point}"
400 payload_misid_pi = f
"ParticleReweighting:{pid_alias}_misid_pi_{channel_misid_pi}_{working_point}"
401 payload_misid_K = f
"ParticleReweighting:{pid_alias}_misid_K_{channel_misid_K}_{working_point}"
404 path.add_module(
"ParticleWeighting",
405 particleList=outputListName,
406 tableName=payload_eff).set_name(f
"ParticleWeighting_eff_{outputListName}")
407 path.add_module(
"ParticleWeighting",
408 particleList=outputListName,
409 tableName=payload_misid_pi).set_name(f
"ParticleWeighting_misid_pi_{outputListName}")
410 if classification ==
"global":
411 path.add_module(
"ParticleWeighting",
412 particleList=outputListName,
413 tableName=payload_misid_K).set_name(f
"ParticleWeighting_misid_K_{outputListName}")
417 cut = f
"[{pid_alias} >= extraInfo({payload_eff}_threshold)]"
418 ma.applyCuts(outputListName, cut, path=path)
421 weight_aliases_to_var = {
422 f
"weight_{pid_alias}_eff_{working_point}": f
"extraInfo({payload_eff}_data_MC_ratio)",
423 f
"weight_{pid_alias}_misid_pi_{working_point}": f
"extraInfo({payload_misid_pi}_data_MC_ratio)",
425 f
"weight_{pid_alias}_eff_{working_point}_stat_up": f
"extraInfo({payload_eff}_data_MC_uncertainty_stat_up)",
426 f
"weight_{pid_alias}_eff_{working_point}_stat_dn": f
"extraInfo({payload_eff}_data_MC_uncertainty_stat_dn)",
427 f
"weight_{pid_alias}_eff_{working_point}_sys_up": f
"extraInfo({payload_eff}_data_MC_uncertainty_sys_up)",
428 f
"weight_{pid_alias}_eff_{working_point}_sys_dn": f
"extraInfo({payload_eff}_data_MC_uncertainty_sys_dn)",
429 f
"weight_{pid_alias}_misid_pi_{working_point}_stat_up": f
"extraInfo({payload_misid_pi}_data_MC_uncertainty_stat_up)",
430 f
"weight_{pid_alias}_misid_pi_{working_point}_stat_dn": f
"extraInfo({payload_misid_pi}_data_MC_uncertainty_stat_dn)",
431 f
"weight_{pid_alias}_misid_pi_{working_point}_sys_up": f
"extraInfo({payload_misid_pi}_data_MC_uncertainty_sys_up)",
432 f
"weight_{pid_alias}_misid_pi_{working_point}_sys_dn": f
"extraInfo({payload_misid_pi}_data_MC_uncertainty_sys_dn)",
441 f
"weight_{pid_alias}_eff_{working_point}_rel_stat_up":
442 f
"formula(1+[extraInfo({payload_eff}_data_MC_uncertainty_stat_up)/extraInfo({payload_eff}_data_MC_ratio)])",
443 f
"weight_{pid_alias}_eff_{working_point}_rel_stat_dn":
444 f
"formula(1-[extraInfo({payload_eff}_data_MC_uncertainty_stat_dn)/extraInfo({payload_eff}_data_MC_ratio)])",
445 f
"weight_{pid_alias}_eff_{working_point}_rel_sys_up":
446 f
"formula(1+[extraInfo({payload_eff}_data_MC_uncertainty_sys_up)/extraInfo({payload_eff}_data_MC_ratio)])",
447 f
"weight_{pid_alias}_eff_{working_point}_rel_sys_dn":
448 f
"formula(1-[extraInfo({payload_eff}_data_MC_uncertainty_sys_dn)/extraInfo({payload_eff}_data_MC_ratio)])",
449 f
"weight_{pid_alias}_misid_pi_{working_point}_rel_stat_up":
450 f
"formula(1+[extraInfo({payload_misid_pi}_data_MC_uncertainty_stat_up)/extraInfo({payload_misid_pi}_data_MC_ratio)])",
451 f
"weight_{pid_alias}_misid_pi_{working_point}_rel_stat_dn":
452 f
"formula(1-[extraInfo({payload_misid_pi}_data_MC_uncertainty_stat_dn)/extraInfo({payload_misid_pi}_data_MC_ratio)])",
453 f
"weight_{pid_alias}_misid_pi_{working_point}_rel_sys_up":
454 f
"formula(1+[extraInfo({payload_misid_pi}_data_MC_uncertainty_sys_up)/extraInfo({payload_misid_pi}_data_MC_ratio)])",
455 f
"weight_{pid_alias}_misid_pi_{working_point}_rel_sys_dn":
456 f
"formula(1-[extraInfo({payload_misid_pi}_data_MC_uncertainty_sys_dn)/extraInfo({payload_misid_pi}_data_MC_ratio)])",
458 if classification ==
"global":
459 weight_aliases_to_var.update({
460 f
"weight_{pid_alias}_misid_K_{working_point}": f
"extraInfo({payload_misid_K}_data_MC_ratio)",
462 f
"weight_{pid_alias}_misid_K_{working_point}_stat_up": f
"extraInfo({payload_misid_K}_data_MC_uncertainty_stat_up)",
463 f
"weight_{pid_alias}_misid_K_{working_point}_stat_dn": f
"extraInfo({payload_misid_K}_data_MC_uncertainty_stat_dn)",
464 f
"weight_{pid_alias}_misid_K_{working_point}_sys_up": f
"extraInfo({payload_misid_K}_data_MC_uncertainty_sys_up)",
465 f
"weight_{pid_alias}_misid_K_{working_point}_sys_dn": f
"extraInfo({payload_misid_K}_data_MC_uncertainty_sys_dn)",
467 f
"weight_{pid_alias}_misid_K_{working_point}_rel_stat_up":
468 f
"formula(1+[extraInfo({payload_misid_K}_data_MC_uncertainty_stat_up)/extraInfo({payload_misid_K}_data_MC_ratio)])",
469 f
"weight_{pid_alias}_misid_K_{working_point}_rel_stat_dn":
470 f
"formula(1-[extraInfo({payload_misid_K}_data_MC_uncertainty_stat_dn)/extraInfo({payload_misid_K}_data_MC_ratio)])",
471 f
"weight_{pid_alias}_misid_K_{working_point}_rel_sys_up":
472 f
"formula(1+[extraInfo({payload_misid_K}_data_MC_uncertainty_sys_up)/extraInfo({payload_misid_K}_data_MC_ratio)])",
473 f
"weight_{pid_alias}_misid_K_{working_point}_rel_sys_dn":
474 f
"formula(1-[extraInfo({payload_misid_K}_data_MC_uncertainty_sys_dn)/extraInfo({payload_misid_K}_data_MC_ratio)])",
477 for alias, var
in weight_aliases_to_var.items():
478 variables.addAlias(alias, var)
480 return pid_alias, list(weight_aliases_to_var.keys())
483 def stdE(listtype=_defaultlist,
488 channel_eff="combination",
489 channel_misid_pi="combination",
490 channel_misid_K="combination",
492 outputListLabel=None,
493 trainingModeMulticlass=_TrainingMode.c_Multiclass,
494 trainingModeBinary=_TrainingMode.c_Classification,
496 """ Function to prepare one of several standardized types of electron lists.
497 See the documentation of `stdLep` for details.
499 It also accepts any of the legacy definitions
500 for the ``listtype`` parameter (aka ``working_point`` in `stdLep`) to fall back to the `stdCharged` behaviour:
512 tuple(str, list(str)): the alias for the electron ID variable, and the list of aliases for the weights.
515 if listtype
in _stdnames + _effnames:
517 return _pidnames[_chargednames.index(
"e")],
None
519 return stdLep(Const.electron.getPDGCode(), listtype, method, classification, lid_weights_gt,
521 channel_eff=channel_eff,
522 channel_misid_pi=channel_misid_pi,
523 channel_misid_K=channel_misid_K,
524 inputListName=inputListName,
525 outputListLabel=outputListLabel,
526 trainingModeMulticlass=trainingModeMulticlass,
527 trainingModeBinary=trainingModeBinary,
531 def stdMu(listtype=_defaultlist,
536 channel_eff="combination",
537 channel_misid_pi="combination",
538 channel_misid_K="combination",
540 outputListLabel=None,
541 trainingModeMulticlass=_TrainingMode.c_Multiclass,
542 trainingModeBinary=_TrainingMode.c_Classification,
544 """ Function to prepare one of several standardized types of muon lists.
545 See the documentation of `stdLep` for details.
547 It also accepts any of the legacy definitions
548 for the ``listtype`` parameter (aka ``working_point`` in `stdLep`) to fall back to the `stdCharged` behaviour:
560 tuple(str, list(str)): the alias for the muon ID variable, and the list of aliases for the weights.
563 if listtype
in _stdnames + _effnames:
565 return _pidnames[_chargednames.index(
"mu")],
None
567 return stdLep(Const.muon.getPDGCode(), listtype, method, classification, lid_weights_gt,
569 channel_eff=channel_eff,
570 channel_misid_pi=channel_misid_pi,
571 channel_misid_K=channel_misid_K,
572 inputListName=inputListName,
573 outputListLabel=outputListLabel,
574 trainingModeMulticlass=trainingModeMulticlass,
575 trainingModeBinary=trainingModeBinary,
579 def stdMostLikely(pidPriors=None, suffix='', custom_cuts='', path=None):
581 Function to prepare most likely particle lists according to PID likelihood, refer to stdCharged for details
583 @param pidPriors list of 6 float numbers used to reweight PID likelihoods, for e, mu, pi, K, p and d
584 @param suffix string added to the end of particle list names
585 @param custom_cuts custom selection cut string, if empty, standard track quality cuts will be applied
586 @param path modules are added to this path
592 if pidPriors
is not None:
593 args = str(pidPriors)[1:-1]
594 trackQuality =
'thetaInCDCAcceptance and nCDCHits>20'
595 if custom_cuts !=
'':
596 trackQuality = custom_cuts
597 for name
in _chargednames:
598 ma.fillParticleList(f
'{name}+:{_mostLikelyList}{suffix}',
599 f
'pidIsMostLikely({args}) > 0 and {trackQuality}',
True, path=path)
This class provides a set of constants for the framework.