Belle II Software development
modularAnalysis.py
1#!/usr/bin/env python3
2
3
10
11"""
12This module defines wrapper functions around the analysis modules.
13"""
14
15import b2bii
16from basf2 import register_module, create_path
17from basf2 import B2INFO, B2WARNING, B2ERROR, B2FATAL
18import basf2
19import subprocess
20
21
22def setAnalysisConfigParams(configParametersAndValues, path):
23 """
24 Sets analysis configuration parameters.
25
26 These are:
27
28 - 'tupleStyle': 'Default' (default) or 'Laconic'
29
30 - defines the style of the branch name in the ntuple
31
32 - 'mcMatchingVersion': Specifies what version of mc matching algorithm is going to be used:
33
34 - 'Belle' - analysis of Belle MC
35 - 'BelleII' (default) - all other cases
36
37 @param configParametersAndValues dictionary of parameters and their values of the form {param1: value, param2: value, ...)
38 @param modules are added to this path
39 """
40
41 conf = register_module('AnalysisConfiguration')
42
43 allParameters = ['tupleStyle', 'mcMatchingVersion']
44
45 keys = configParametersAndValues.keys()
46 for key in keys:
47 if key not in allParameters:
48 allParametersString = ', '.join(allParameters)
49 B2ERROR('Invalid analysis configuration parameter: ' + key + '.\n'
50 'Please use one of the following: ' + allParametersString)
51
52 for param in allParameters:
53 if param in configParametersAndValues:
54 conf.param(param, configParametersAndValues.get(param))
55
56 path.add_module(conf)
57
58
59def inputMdst(filename, path, environmentType='default', skipNEvents=0, entrySequence=None, *, parentLevel=0, **kwargs):
60 """
61 Loads the specified :ref:`mDST <mdst>` (or :ref:`uDST <analysis_udstoutput>`) file with the RootInput module.
62
63 The correct environment (e.g. magnetic field settings) is determined from
64 ``environmentType``. Options are either: 'default' (for Belle II MC and
65 data: falls back to database), 'Belle': for analysis of converted Belle 1
66 data and MC.
67
68 Parameters:
69 filename (str): the name of the file to be loaded
70 path (basf2.Path): modules are added to this path
71 environmentType (str): type of the environment to be loaded (either 'default' or 'Belle')
72 skipNEvents (int): N events of the input file are skipped
73 entrySequence (str): The number sequences (e.g. 23:42,101) defining the entries which are processed.
74 parentLevel (int): Number of generations of parent files (files used as input when creating a file) to be read
75 """
76
77 # FIXME remove this check of "filename" at release-07
78 if filename == 'default':
79 B2FATAL("""
80We have simplified the arguments to inputMdst! If you are running on Belle II
81data or MC, you don't have to use "default" any more.
82Please replace:
83 inputMdst("default", "/your/input/file.root", path=mypath)
84With:
85 inputMdst("/your/input/file.root", path=mypath)
86 """)
87 elif filename == "Belle":
88 B2FATAL("""
89We have reordered the arguments to inputMdst! If you are running on Belle 1
90data or MC, you need to specify the 'environmentType'.
91Please replace:
92 inputMdst("Belle", "/your/input/file.root", path=mypath)
93With:
94 inputMdst("/your/input/file.root", path=mypath, environmentType='Belle')
95 """)
96 elif filename in [f"MC{i}" for i in range(5, 10)]:
97 B2FATAL(f"We no longer support the MC version {filename}. Sorry.")
98
99 if entrySequence is not None:
100 entrySequence = [entrySequence]
101
102 inputMdstList([filename], path, environmentType, skipNEvents, entrySequence, parentLevel=parentLevel, **kwargs)
103
104
105def inputMdstList(
106 filelist,
107 path,
108 environmentType='default',
109 skipNEvents=0,
110 entrySequences=None,
111 *,
112 parentLevel=0,
113 useB2BIIDBCache=True):
114 """
115 Loads the specified list of :ref:`mDST <mdst>` (or :ref:`uDST <analysis_udstoutput>`) files with the RootInput module.
116
117 The correct environment (e.g. magnetic field settings) is determined from
118 ``environmentType``. Options are either: 'default' (for Belle II MC and
119 data: falls back to database), 'Belle': for analysis of converted Belle 1
120 data and MC.
121
122 Parameters:
123 filelist (list(str)): the filename list of files to be loaded
124 path (basf2.Path): modules are added to this path
125 environmentType (str): type of the environment to be loaded (either 'default' or 'Belle')
126 skipNEvents (int): N events of the input files are skipped
127 entrySequences (list(str)): The number sequences (e.g. 23:42,101) defining
128 the entries which are processed for each inputFileName.
129 parentLevel (int): Number of generations of parent files (files used as input when creating a file) to be read
130 useB2BIIDBCache (bool): Loading of local KEKCC database (only to be deactivated in very special cases)
131 """
132
133 # FIXME remove this check of "filename" at release-07
134 if filelist == 'default':
135 B2FATAL("""
136We have simplified the arguments to inputMdstList! If you are running on
137Belle II data or MC, you don't have to use "default" any more.
138Please replace:
139 inputMdstList("default", list_of_your_files, path=mypath)
140With:
141 inputMdstList(list_of_your_files, path=mypath)
142 """)
143 elif filelist == "Belle":
144 B2FATAL("""
145We have reordered the arguments to inputMdstList! If you are running on
146Belle 1 data or MC, you need to specify the 'environmentType'.
147Please replace:
148 inputMdstList("Belle", list_of_your_files, path=mypath)
149With:
150 inputMdstList(list_of_your_files, path=mypath, environmentType='Belle')
151 """)
152 elif filelist in [f"MC{i}" for i in range(5, 10)]:
153 B2FATAL(f"We no longer support the MC version {filelist}. Sorry.")
154
155 roinput = register_module('RootInput')
156 roinput.param('inputFileNames', filelist)
157 roinput.param('skipNEvents', skipNEvents)
158 if entrySequences is not None:
159 roinput.param('entrySequences', entrySequences)
160 roinput.param('parentLevel', parentLevel)
161
162 path.add_module(roinput)
163 path.add_module('ProgressBar')
164
165 if environmentType == 'Belle':
166 # Belle 1 constant magnetic field
167 # -------------------------------
168 # n.b. slightly unfortunate syntax: the MagneticField is a member of the
169 # Belle2 namespace but will be set to the Belle 1 values
170 from ROOT import Belle2 # reduced scope of potentially-misbehaving import
171 from ROOT.Math import XYZVector
172 belle1_field = Belle2.MagneticField()
173 belle1_field.addComponent(Belle2.MagneticFieldComponentConstant(XYZVector(0, 0, 1.5 * Belle2.Unit.T)))
174 Belle2.DBStore.Instance().addConstantOverride("MagneticField", belle1_field, False)
175 # also set the MC matching for Belle 1
176 setAnalysisConfigParams({'mcMatchingVersion': 'Belle'}, path)
178 if useB2BIIDBCache:
179 basf2.conditions.metadata_providers = ["/sw/belle/b2bii/database/conditions/b2bii.sqlite"]
180 basf2.conditions.payload_locations = ["/sw/belle/b2bii/database/conditions/"]
181
182
183def outputMdst(filename, path):
184 """
185 Saves mDST (mini-Data Summary Tables) to the output root file.
186
187 .. warning::
188
189 This function is kept for backward-compatibility.
190 Better to use `mdst.add_mdst_output` directly.
191
192 """
193
194 import mdst
195 mdst.add_mdst_output(path, mc=True, filename=filename)
196
197
198def outputUdst(filename, particleLists=None, includeArrays=None, path=None, dataDescription=None):
199 """
200 Save uDST (user-defined Data Summary Tables) = MDST + Particles + ParticleLists
201 The charge-conjugate lists of those given in particleLists are also stored.
202 Additional Store Arrays and Relations to be stored can be specified via includeArrays
203 list argument.
204
205 Note:
206 This does not reduce the amount of Particle objects saved,
207 see `udst.add_skimmed_udst_output` for a function that does.
208
209 """
210
211 import udst
213 path=path, filename=filename, particleLists=particleLists,
214 additionalBranches=includeArrays, dataDescription=dataDescription)
215
216
217def outputIndex(filename, path, includeArrays=None, keepParents=False, mc=True):
218 """
219 Write out all particle lists as an index file to be reprocessed using parentLevel flag.
220 Additional branches necessary for file to be read are automatically included.
221 Additional Store Arrays and Relations to be stored can be specified via includeArrays
222 list argument.
223
224 @param str filename the name of the output index file
225 @param str path modules are added to this path
226 @param list(str) includeArrays: datastore arrays/objects to write to the output
227 file in addition to particle lists and related information
228 @param bool keepParents whether the parents of the input event will be saved as the parents of the same event
229 in the output index file. Useful if you are only adding more information to another index file
230 @param bool mc whether the input data is MC or not
231 """
232
233 if includeArrays is None:
234 includeArrays = []
235
236 # Module to mark all branches to not be saved except particle lists
237 onlyPLists = register_module('OnlyWriteOutParticleLists')
238 path.add_module(onlyPLists)
239
240 # Set up list of all other branches we need to make index file complete
241 partBranches = [
242 'Particles',
243 'ParticlesToMCParticles',
244 'ParticlesToPIDLikelihoods',
245 'ParticleExtraInfoMap',
246 'EventExtraInfo'
247 ]
248 branches = ['EventMetaData']
249 persistentBranches = ['FileMetaData']
250 if mc:
251 branches += []
252 # persistentBranches += ['BackgroundInfos']
253 branches += partBranches
254 branches += includeArrays
255
256 r1 = register_module('RootOutput')
257 r1.param('outputFileName', filename)
258 r1.param('additionalBranchNames', branches)
259 r1.param('branchNamesPersistent', persistentBranches)
260 r1.param('keepParents', keepParents)
261 path.add_module(r1)
262
263
264def setupEventInfo(noEvents, path):
265 """
266 Prepare to generate events. This function sets up the EventInfoSetter.
267 You should call this before adding a generator from generators.
268 The experiment and run numbers are set to 0 (run independent generic MC in phase 3).
269 https://xwiki.desy.de/xwiki/rest/p/59192
270
271 Parameters:
272 noEvents (int): number of events to be generated
273 path (basf2.Path): modules are added to this path
274 """
275
276 evtnumbers = register_module('EventInfoSetter')
277 evtnumbers.param('evtNumList', [noEvents])
278 evtnumbers.param('runList', [0])
279 evtnumbers.param('expList', [0])
280 path.add_module(evtnumbers)
281
282
283def loadGearbox(path, silence_warning=False):
284 """
285 Loads Gearbox module to the path.
286
287 Warning:
288 Should be used in a job with *cosmic event generation only*
289
290 Needed for scripts which only generate cosmic events in order to
291 load the geometry.
292
293 @param path modules are added to this path
294 @param silence_warning stops a verbose warning message if you know you want to use this function
295 """
296
297 if not silence_warning:
298 B2WARNING("""You are overwriting the geometry from the database with Gearbox.
299 This is fine if you're generating cosmic events. But in most other cases you probably don't want this.
300
301 If you're really sure you know what you're doing you can suppress this message with:
302
303 >>> loadGearbox(silence_warning=True)
304
305 """)
306
307 paramloader = register_module('Gearbox')
308 path.add_module(paramloader)
309
310
311def printPrimaryMCParticles(path, **kwargs):
312 """
313 Prints all primary MCParticles, that is particles from
314 the physics generator and not particles created by the simulation
315
316 This is equivalent to `printMCParticles(onlyPrimaries=True, path=path) <printMCParticles>` and additional
317 keyword arguments are just forwarded to that function
318 """
319
320 return printMCParticles(onlyPrimaries=True, path=path, **kwargs)
321
322
323def printMCParticles(onlyPrimaries=False, maxLevel=-1, path=None, *,
324 showProperties=False, showMomenta=False, showVertices=False, showStatus=False, suppressPrint=False):
325 """
326 Prints all MCParticles or just primary MCParticles up to specified level. -1 means no limit.
327
328 By default this will print a tree of just the particle names and their pdg
329 codes in the event, for example ::
330
331 [INFO] Content of MCParticle list
332 ├── e- (11)
333 ├── e+ (-11)
334 ╰── Upsilon(4S) (300553)
335 ├── B+ (521)
336 │ ├── anti-D_0*0 (-10421)
337 │ │ ├── D- (-411)
338 │ │ │ ├── K*- (-323)
339 │ │ │ │ ├── anti-K0 (-311)
340 │ │ │ │ │ ╰── K_S0 (310)
341 │ │ │ │ │ ├── pi+ (211)
342 │ │ │ │ │ │ ╰╶╶ p+ (2212)
343 │ │ │ │ │ ╰── pi- (-211)
344 │ │ │ │ │ ├╶╶ e- (11)
345 │ │ │ │ │ ├╶╶ n0 (2112)
346 │ │ │ │ │ ├╶╶ n0 (2112)
347 │ │ │ │ │ ╰╶╶ n0 (2112)
348 │ │ │ │ ╰── pi- (-211)
349 │ │ │ │ ├╶╶ anti-nu_mu (-14)
350 │ │ │ │ ╰╶╶ mu- (13)
351 │ │ │ │ ├╶╶ nu_mu (14)
352 │ │ │ │ ├╶╶ anti-nu_e (-12)
353 │ │ │ │ ╰╶╶ e- (11)
354 │ │ │ ╰── K_S0 (310)
355 │ │ │ ├── pi0 (111)
356 │ │ │ │ ├── gamma (22)
357 │ │ │ │ ╰── gamma (22)
358 │ │ │ ╰── pi0 (111)
359 │ │ │ ├── gamma (22)
360 │ │ │ ╰── gamma (22)
361 │ │ ╰── pi+ (211)
362 │ ├── mu+ (-13)
363 │ │ ├╶╶ anti-nu_mu (-14)
364 │ │ ├╶╶ nu_e (12)
365 │ │ ╰╶╶ e+ (-11)
366 │ ├── nu_mu (14)
367 │ ╰── gamma (22)
368 ...
369
370
371 There's a distinction between primary and secondary particles. Primary
372 particles are the ones created by the physics generator while secondary
373 particles are ones generated by the simulation of the detector interaction.
374
375 Secondaries are indicated with a dashed line leading to the particle name
376 and if the output is to the terminal they will be printed in red. If
377 ``onlyPrimaries`` is True they will not be included in the tree.
378
379 On demand, extra information on all the particles can be displayed by
380 enabling any of the ``showProperties``, ``showMomenta``, ``showVertices``
381 and ``showStatus`` flags. Enabling all of them will look like
382 this::
383
384 ...
385 ╰── pi- (-211)
386 │ mass=0.14 energy=0.445 charge=-1 lifetime=6.36
387 │ p=(0.257, -0.335, 0.0238) |p|=0.423
388 │ production vertex=(0.113, -0.0531, 0.0156), time=0.00589
389 │ status flags=PrimaryParticle, StableInGenerator, StoppedInDetector
390 │ list index=48
391
392 ╰╶╶ n0 (2112)
393 mass=0.94 energy=0.94 charge=0 lifetime=5.28e+03
394 p=(-0.000238, -0.0127, 0.0116) |p|=0.0172
395 production vertex=(144, 21.9, -1.29), time=39
396 status flags=StoppedInDetector
397 creation process=HadronInelastic
398 list index=66
399
400 The first line of extra information is enabled by ``showProperties``, the
401 second line by ``showMomenta``, the third line by ``showVertices`` and the
402 last two lines by ``showStatus``. Note that all values are given in Belle II
403 standard units, that is GeV, centimeter and nanoseconds.
404
405 The depth of the tree can be limited with the ``maxLevel`` argument: If it's
406 bigger than zero it will limit the tree to the given number of generations.
407 A visual indicator will be added after each particle which would have
408 additional daughters that are skipped due to this limit. An example event
409 with ``maxLevel=3`` is given below. In this case only the tau neutrino and
410 the pion don't have additional daughters. ::
411
412 [INFO] Content of MCParticle list
413 ├── e- (11)
414 ├── e+ (-11)
415 ╰── Upsilon(4S) (300553)
416 ├── B+ (521)
417 │ ├── anti-D*0 (-423) → …
418 │ ├── tau+ (-15) → …
419 │ ╰── nu_tau (16)
420 ╰── B- (-521)
421 ├── D*0 (423) → …
422 ├── K*- (-323) → …
423 ├── K*+ (323) → …
424 ╰── pi- (-211)
425
426 The same information will be stored in the branch ``__MCDecayString__`` of
427 TTree created by `VariablesToNtuple` or `VariablesToEventBasedTree` module.
428 This branch is automatically created when `PrintMCParticles` modules is called.
429 Printing the information on the log message can be suppressed if ``suppressPrint``
430 is True, while the branch ``__MCDecayString__``. This option helps to reduce the
431 size of the log message.
432
433 Parameters:
434 onlyPrimaries (bool): If True show only primary particles, that is particles coming from
435 the generator and not created by the simulation.
436 maxLevel (int): If 0 or less print the whole tree, otherwise stop after n generations
437 showProperties (bool): If True show mass, energy and charge of the particles
438 showMomenta (bool): if True show the momenta of the particles
439 showVertices (bool): if True show production vertex and production time of all particles
440 showStatus (bool): if True show some status information on the particles.
441 For secondary particles this includes creation process.
442 suppressPrint (bool): if True printing the information on the log message is suppressed.
443 Even if True, the branch ``__MCDecayString__`` is created.
444 """
445
446 return path.add_module(
447 "PrintMCParticles",
448 onlyPrimaries=onlyPrimaries,
449 maxLevel=maxLevel,
450 showProperties=showProperties,
451 showMomenta=showMomenta,
452 showVertices=showVertices,
453 showStatus=showStatus,
454 suppressPrint=suppressPrint,
455 )
456
457
458def correctBrems(outputList,
459 inputList,
460 gammaList,
461 maximumAcceptance=3.0,
462 multiplePhotons=False,
463 usePhotonOnlyOnce=True,
464 writeOut=False,
465 path=None):
466 """
467 For each particle in the given ``inputList``, copies it to the ``outputList`` and adds the
468 4-vector of the photon(s) in the ``gammaList`` which has(have) a weighted named relation to
469 the particle's track, set by the ``ECLTrackBremFinder`` module during reconstruction.
470
471 Warning:
472 So far, there haven't been any comprehensive comparisons of the performance of the `BremsFinder` module, which
473 is called in this function, with the `BelleBremRecovery` module, which is called via the `correctBremsBelle`
474 function. If your analysis is very sensitive to the Bremsstrahlung corrections, it is currently advised to use
475 `correctBremsBelle`.
476
477 The reason is that studies by the tau WG revealed that in the past the cuts applied by the
478 ``ECLTrackBremFinder`` module were too tight. They were only loosened for proc16 and MC16. New performance
479 studies are needed to verify that now this module outperforms the Belle-like approach.
480
481 Information:
482 A detailed description of how the weights are set can be found directly at the documentation of the
483 `BremsFinder` module.
484
485 Please note that a new particle is always generated, with the old particle and -if found- one or more
486 photons as daughters.
487
488 The ``inputList`` should contain particles with associated tracks. Otherwise, the module will exit with an error.
489
490 The ``gammaList`` should contain photons. Otherwise, the module will exit with an error.
491
492 @param outputList The output particle list name containing the corrected particles
493 @param inputList The initial particle list name containing the particles to correct. *It should already exist.*
494 @param gammaList The photon list containing possibly bremsstrahlung photons; *It should already exist.*
495 @param maximumAcceptance Maximum value of the relation weight. Should be a number between [0,3)
496 @param multiplePhotons Whether to use only one photon (the one with the smallest acceptance) or as many as possible
497 @param usePhotonOnlyOnce If true, each brems candidate is used to correct only the track with the smallest relation weight
498 @param writeOut Whether `RootOutput` module should save the created ``outputList``
499 @param path The module is added to this path
500 """
501
502 import b2bii
503 if b2bii.isB2BII():
504 B2ERROR("The BremsFinder can only be run over Belle II data.")
505
506 bremscorrector = register_module('BremsFinder')
507 bremscorrector.set_name('bremsCorrector_' + outputList)
508 bremscorrector.param('inputList', inputList)
509 bremscorrector.param('outputList', outputList)
510 bremscorrector.param('gammaList', gammaList)
511 bremscorrector.param('maximumAcceptance', maximumAcceptance)
512 bremscorrector.param('multiplePhotons', multiplePhotons)
513 bremscorrector.param('usePhotonOnlyOnce', usePhotonOnlyOnce)
514 bremscorrector.param('writeOut', writeOut)
515 path.add_module(bremscorrector)
516
517
518def copyList(outputListName, inputListName, writeOut=False, path=None):
519 """
520 Copy all Particle indices from input ParticleList to the output ParticleList.
521 Note that the Particles themselves are not copied. The original and copied
522 ParticleLists will point to the same Particles.
523
524 @param ouputListName copied ParticleList
525 @param inputListName original ParticleList to be copied
526 @param writeOut whether RootOutput module should save the created ParticleList
527 @param path modules are added to this path
528 """
529
530 copyLists(outputListName, [inputListName], writeOut, path)
531
532
533def correctBremsBelle(outputListName,
534 inputListName,
535 gammaListName,
536 multiplePhotons=True,
537 angleThreshold=0.05,
538 usePhotonOnlyOnce=False,
539 writeOut=False,
540 path=None):
541 """
542 Run the Belle - like brems finding on the ``inputListName`` of charged particles.
543 Adds all photons in ``gammaListName`` to a copy of the charged particle that are within
544 ``angleThreshold``.
545
546 Tip:
547 Studies by the tau WG show that using a rather wide opening angle (up to
548 0.2 rad) and rather low energetic photons results in good correction.
549 However, this should only serve as a starting point for your own studies
550 because the optimal criteria are likely mode-dependent
551
552 Parameters:
553 outputListName (str): The output charged particle list containing the corrected charged particles
554 inputListName (str): The initial charged particle list containing the charged particles to correct.
555 gammaListName (str): The gammas list containing possibly radiative gammas, should already exist.
556 multiplePhotons (bool): How many photons should be added to the charged particle? nearest one -> False,
557 add all the photons within the cone -> True
558 angleThreshold (float): The maximum angle in radians between the charged particle and the (radiative)
559 gamma to be accepted.
560 writeOut (bool): whether RootOutput module should save the created ParticleList
561 usePhotonOnlyOnce (bool): If true, a photon is used for correction of the closest charged particle in the inputList.
562 If false, a photon is allowed to be used for correction multiple times (Default).
563
564 Warning:
565 One cannot use a photon twice to reconstruct a composite particle. Thus, for example, if ``e+`` and ``e-`` are corrected
566 with a ``gamma``, the pair of ``e+`` and ``e-`` cannot form a ``J/psi -> e+ e-`` candidate.
567
568 path (basf2.Path): modules are added to this path
569 """
570
571 fsrcorrector = register_module('BelleBremRecovery')
572 fsrcorrector.set_name('BelleFSRCorrection_' + outputListName)
573 fsrcorrector.param('inputListName', inputListName)
574 fsrcorrector.param('outputListName', outputListName)
575 fsrcorrector.param('gammaListName', gammaListName)
576 fsrcorrector.param('multiplePhotons', multiplePhotons)
577 fsrcorrector.param('angleThreshold', angleThreshold)
578 fsrcorrector.param('usePhotonOnlyOnce', usePhotonOnlyOnce)
579 fsrcorrector.param('writeOut', writeOut)
580 path.add_module(fsrcorrector)
581
582
583def copyLists(outputListName, inputListNames, writeOut=False, path=None):
584 """
585 Copy all Particle indices from all input ParticleLists to the
586 single output ParticleList.
587 Note that the Particles themselves are not copied.
588 The original and copied ParticleLists will point to the same Particles.
589
590 Duplicates are removed based on the first-come, first-served principle.
591 Therefore, the order of the input ParticleLists matters.
592
593 .. seealso::
594 If you want to select the best duplicate based on another criterion, have
595 a look at the function `mergeListsWithBestDuplicate`.
596
597 .. note::
598 Two particles that differ only by the order of their daughters are
599 considered duplicates and one of them will be removed.
600
601 @param ouputListName copied ParticleList
602 @param inputListName vector of original ParticleLists to be copied
603 @param writeOut whether RootOutput module should save the created ParticleList
604 @param path modules are added to this path
605 """
606
607 pmanipulate = register_module('ParticleListManipulator')
608 pmanipulate.set_name('PListCopy_' + outputListName)
609 pmanipulate.param('outputListName', outputListName)
610 pmanipulate.param('inputListNames', inputListNames)
611 pmanipulate.param('writeOut', writeOut)
612 path.add_module(pmanipulate)
613
614
615def copyParticles(outputListName, inputListName, writeOut=False, path=None):
616 """
617 Create copies of Particles given in the input ParticleList and add them to the output ParticleList.
618
619 The existing relations of the original Particle (or it's (grand-)^n-daughters)
620 are copied as well. Note that only the relation is copied and that the related
621 object is not. Copied particles are therefore related to the *same* object as
622 the original ones.
623
624 @param ouputListName new ParticleList filled with copied Particles
625 @param inputListName input ParticleList with original Particles
626 @param writeOut whether RootOutput module should save the created ParticleList
627 @param path modules are added to this path
628 """
629
630 # first copy original particles to the new ParticleList
631 pmanipulate = register_module('ParticleListManipulator')
632 pmanipulate.set_name('PListCopy_' + outputListName)
633 pmanipulate.param('outputListName', outputListName)
634 pmanipulate.param('inputListNames', [inputListName])
635 pmanipulate.param('writeOut', writeOut)
636 path.add_module(pmanipulate)
637
638 # now replace original particles with their copies
639 pcopy = register_module('ParticleCopier')
640 pcopy.param('inputListNames', [outputListName])
641 path.add_module(pcopy)
642
643
644def cutAndCopyLists(outputListName, inputListNames, cut, writeOut=False, path=None):
645 """
646 Copy candidates from all lists in ``inputListNames`` to
647 ``outputListName`` if they pass ``cut`` (given selection criteria).
648
649 Note:
650 Note that the Particles themselves are not copied.
651 The original and copied ParticleLists will point to the same Particles.
652
653 Example:
654 Require energetic pions safely inside the cdc
655
656 .. code-block:: python
657
658 cutAndCopyLists("pi+:energeticPions", ["pi+:good", "pi+:loose"], "[E > 2] and thetaInCDCAcceptance", path=mypath)
659
660 Warning:
661 You must use square braces ``[`` and ``]`` for conditional statements.
662
663 Parameters:
664 outputListName (str): the new ParticleList name
665 inputListName (list(str)): list of input ParticleList names
666 cut (str): Candidates that do not pass these selection criteria are removed from the ParticleList
667 writeOut (bool): whether RootOutput module should save the created ParticleList
668 path (basf2.Path): modules are added to this path
669 """
670
671 pmanipulate = register_module('ParticleListManipulator')
672 pmanipulate.set_name('PListCutAndCopy_' + outputListName)
673 pmanipulate.param('outputListName', outputListName)
674 pmanipulate.param('inputListNames', inputListNames)
675 pmanipulate.param('cut', cut)
676 pmanipulate.param('writeOut', writeOut)
677 path.add_module(pmanipulate)
678
679
680def cutAndCopyList(outputListName, inputListName, cut, writeOut=False, path=None):
681 """
682 Copy candidates from ``inputListName`` to ``outputListName`` if they pass
683 ``cut`` (given selection criteria).
684
685 Note:
686 Note the Particles themselves are not copied.
687 The original and copied ParticleLists will point to the same Particles.
688
689 Example:
690 require energetic pions safely inside the cdc
691
692 .. code-block:: python
693
694 cutAndCopyList("pi+:energeticPions", "pi+:loose", "[E > 2] and thetaInCDCAcceptance", path=mypath)
695
696 Warning:
697 You must use square braces ``[`` and ``]`` for conditional statements.
698
699 Parameters:
700 outputListName (str): the new ParticleList name
701 inputListName (str): input ParticleList name
702 cut (str): Candidates that do not pass these selection criteria are removed from the ParticleList
703 writeOut (bool): whether RootOutput module should save the created ParticleList
704 path (basf2.Path): modules are added to this path
705 """
706
707 cutAndCopyLists(outputListName, [inputListName], cut, writeOut, path)
708
709
710def removeTracksForTrackingEfficiencyCalculation(inputListNames, fraction, path=None):
711 """
712 Randomly remove tracks from the provided particle lists to estimate the tracking efficiency.
713 Takes care of the duplicates, if any.
714
715 Parameters:
716 inputListNames (list(str)): input particle list names
717 fraction (float): fraction of particles to be removed randomly
718 path (basf2.Path): module is added to this path
719 """
720
721 trackingefficiency = register_module('TrackingEfficiency')
722 trackingefficiency.param('particleLists', inputListNames)
723 trackingefficiency.param('frac', fraction)
724 path.add_module(trackingefficiency)
725
726
727def scaleTrackMomenta(inputListNames, scale=float('nan'), payloadName="", scalingFactorName="SF", path=None):
728 """
729 Scale momenta of the particles according to a scaling factor scale.
730 This scaling factor can either be given as constant number or as the name of the payload which contains
731 the variable scale factors.
732 If the particle list contains composite particles, the momenta of the track-based daughters are scaled.
733 Subsequently, the momentum of the mother particle is updated as well.
734
735 Parameters:
736 inputListNames (list(str)): input particle list names
737 scale (float): scaling factor (1.0 -- no scaling)
738 payloadName (str): name of the payload which contains the phase-space dependent scaling factors
739 scalingFactorName (str): name of scaling factor variable in the payload.
740 path (basf2.Path): module is added to this path
741 """
742
743 import b2bii
744 if b2bii.isB2BII():
745 B2ERROR("The tracking momentum scaler can only be run over Belle II data.")
746
747 TrackingMomentumScaleFactors = register_module('TrackingMomentumScaleFactors')
748 TrackingMomentumScaleFactors.param('particleLists', inputListNames)
749 TrackingMomentumScaleFactors.param('scale', scale)
750 TrackingMomentumScaleFactors.param('payloadName', payloadName)
751 TrackingMomentumScaleFactors.param('scalingFactorName', scalingFactorName)
752
753 path.add_module(TrackingMomentumScaleFactors)
754
755
756def correctTrackEnergy(inputListNames, correction=float('nan'), payloadName="", correctionName="SF", path=None):
757 """
758 Correct the energy loss of tracks according to a 'correction' value.
759 This correction can either be given as constant number or as the name of the payload which contains
760 the variable corrections.
761 If the particle list contains composite particles, the momenta of the track-based daughters are corrected.
762 Subsequently, the momentum of the mother particle is updated as well.
763
764 Parameters:
765 inputListNames (list(str)): input particle list names
766 correction (float): correction value to be subtracted to the particle energy (0.0 -- no correction)
767 payloadName (str): name of the payload which contains the phase-space dependent scaling factors
768 correctionName (str): name of correction variable in the payload.
769 path (basf2.Path): module is added to this path
770 """
771
772 import b2bii
773 if b2bii.isB2BII():
774 B2ERROR("The tracking energy correction can only be run over Belle II data.")
775
776 TrackingEnergyLossCorrection = register_module('TrackingEnergyLossCorrection')
777 TrackingEnergyLossCorrection.param('particleLists', inputListNames)
778 TrackingEnergyLossCorrection.param('correction', correction)
779 TrackingEnergyLossCorrection.param('payloadName', payloadName)
780 TrackingEnergyLossCorrection.param('correctionName', correctionName)
781
782 path.add_module(TrackingEnergyLossCorrection)
783
784
785def smearTrackMomenta(inputListNames, payloadName="", smearingFactorName="smear", path=None):
786 """
787 Smear the momenta of the particles according the values read from the given payload.
788 If the particle list contains composite particles, the momenta of the track-based daughters are smeared.
789 Subsequently, the momentum of the mother particle is updated as well.
790
791 Parameters:
792 inputListNames (list(str)): input particle list names
793 payloadName (str): name of the payload which contains the smearing values
794 smearingFactorName (str): name of smearing factor variable in the payload.
795 path (basf2.Path): module is added to this path
796 """
797
798 TrackingMomentumScaleFactors = register_module('TrackingMomentumScaleFactors')
799 TrackingMomentumScaleFactors.param('particleLists', inputListNames)
800 TrackingMomentumScaleFactors.param('payloadName', payloadName)
801 TrackingMomentumScaleFactors.param('smearingFactorName', smearingFactorName)
802
803 path.add_module(TrackingMomentumScaleFactors)
804
805
806def mergeListsWithBestDuplicate(outputListName,
807 inputListNames,
808 variable,
809 preferLowest=True,
810 writeOut=False,
811 ignoreMotherFlavor=False,
812 path=None):
813 """
814 Merge input ParticleLists into one output ParticleList. Only the best
815 among duplicates is kept. The lowest or highest value (configurable via
816 preferLowest) of the provided variable determines which duplicate is the
817 best.
818
819 @param ouputListName name of merged ParticleList
820 @param inputListName vector of original ParticleLists to be merged
821 @param variable variable to determine best duplicate
822 @param preferLowest whether lowest or highest value of variable should be preferred
823 @param writeOut whether RootOutput module should save the created ParticleList
824 @param ignoreMotherFlavor whether the flavor of the mother particle is ignored when trying to find duplicates
825 @param path modules are added to this path
826 """
827
828 pmanipulate = register_module('ParticleListManipulator')
829 pmanipulate.set_name('PListMerger_' + outputListName)
830 pmanipulate.param('outputListName', outputListName)
831 pmanipulate.param('inputListNames', inputListNames)
832 pmanipulate.param('variable', variable)
833 pmanipulate.param('preferLowest', preferLowest)
834 pmanipulate.param('writeOut', writeOut)
835 pmanipulate.param('ignoreMotherFlavor', ignoreMotherFlavor)
836 path.add_module(pmanipulate)
837
838
839def fillSignalSideParticleList(outputListName, decayString, path):
840 """
841 This function should only be used in the ROE path, that is a path
842 that is executed for each ROE object in the DataStore.
843
844 Example: fillSignalSideParticleList('gamma:sig','B0 -> K*0 ^gamma', roe_path)
845
846 Function will create a ParticleList with name 'gamma:sig' which will be filled
847 with the existing photon Particle, being the second daughter of the B0 candidate
848 to which the ROE object has to be related.
849
850 @param ouputListName name of the created ParticleList
851 @param decayString specify Particle to be added to the ParticleList
852 """
853
854 pload = register_module('SignalSideParticleListCreator')
855 pload.set_name('SSParticleList_' + outputListName)
856 pload.param('particleListName', outputListName)
857 pload.param('decayString', decayString)
858 path.add_module(pload)
859
860
861def fillParticleLists(decayStringsWithCuts, writeOut=False, path=None, enforceFitHypothesis=False,
862 loadPhotonsFromKLM=False):
863 """
864 Creates Particles of the desired types from the corresponding ``mdst`` dataobjects,
865 loads them to the ``StoreArray<Particle>`` and fills the ParticleLists.
866
867 The multiple ParticleLists with their own selection criteria are specified
868 via list tuples (decayString, cut), for example
869
870 .. code-block:: python
871
872 kaons = ('K+:mykaons', 'kaonID>0.1')
873 pions = ('pi+:mypions','pionID>0.1')
874 fillParticleLists([kaons, pions], path=mypath)
875
876 If you are unsure what selection you want, you might like to see the
877 :doc:`StandardParticles` functions.
878
879 The type of the particles to be loaded is specified via the decayString module parameter.
880 The type of the ``mdst`` dataobject that is used as an input is determined from the type of
881 the particle. The following types of the particles can be loaded:
882
883 * charged final state particles (input ``mdst`` type = Tracks)
884 - e+, mu+, pi+, K+, p, deuteron (and charge conjugated particles)
885
886 * neutral final state particles
887 - "gamma" (input ``mdst`` type = ECLCluster)
888 - "K_S0", "Lambda0" (input ``mdst`` type = V0)
889 - "K_L0", "n0" (input ``mdst`` type = KLMCluster or ECLCluster)
890
891 Note:
892 For "K_S0" and "Lambda0" you must specify the daughter ordering.
893
894 For example, to load V0s as :math:`\\Lambda^0\\to p^+\\pi^-` decays from V0s:
895
896 .. code-block:: python
897
898 v0lambdas = ('Lambda0 -> p+ pi-', '0.9 < M < 1.3')
899 fillParticleLists([kaons, pions, v0lambdas], path=mypath)
900
901 Tip:
902 Gammas can also be loaded from KLMClusters by explicitly setting the
903 parameter ``loadPhotonsFromKLM`` to True. However, this should only be
904 done in selected use-cases and the effect should be studied carefully.
905
906 Tip:
907 For "K_L0" it is now possible to load from ECLClusters, to revert to
908 the old (Belle) behavior, you can require ``'isFromKLM > 0'``.
909
910 .. code-block:: python
911
912 klongs = ('K_L0', 'isFromKLM > 0')
913 fillParticleLists([kaons, pions, klongs], path=mypath)
914
915 * Charged kinks final state particles (input ``mdst`` type = Kink)
916
917 Note:
918 To reconstruct charged particle kink you must specify the daughter.
919
920 For example, to load Kinks as :math:`K^- \\to \\pi^-\\pi^0` decays from Kinks:
921
922 .. code-block:: python
923
924 kinkKaons = ('K- -> pi-', yourCut)
925 fillParticleLists([kaons, pions, v0lambdas, kinkKaons], path=mypath)
926
927
928 Parameters:
929 decayStringsWithCuts (list): A list of python ntuples of (decayString, cut).
930 The decay string determines the type of Particle
931 and the name of the ParticleList.
932 If the input MDST type is V0 the whole
933 decay chain needs to be specified, so that
934 the user decides and controls the daughters
935 ' order (e.g. ``K_S0 -> pi+ pi-``).
936 If the input MDST type is Kink the decay chain needs to be specified
937 with only one daughter (e.g. ``K- -> pi-``).
938 The cut is the selection criteria
939 to be added to the ParticleList. It can be an empty string.
940 writeOut (bool): whether RootOutput module should save the created ParticleList
941 path (basf2.Path): modules are added to this path
942 enforceFitHypothesis (bool): If true, Particles will be created only for the tracks which have been fitted
943 using a mass hypothesis of the exact type passed to fillParticleLists().
944 If enforceFitHypothesis is False (the default) the next closest fit hypothesis
945 in terms of mass difference will be used if the fit using exact particle
946 type is not available.
947 loadPhotonsFromKLM (bool): If true, photon candidates will be created from KLMClusters as well.
948 """
949
950 pload = register_module('ParticleLoader')
951 pload.set_name('ParticleLoader_' + 'PLists')
952 pload.param('decayStrings', [decayString for decayString, cut in decayStringsWithCuts])
953 pload.param('writeOut', writeOut)
954 pload.param("enforceFitHypothesis", enforceFitHypothesis)
955 path.add_module(pload)
956
957 from ROOT import Belle2
958 decayDescriptor = Belle2.DecayDescriptor()
959 for decayString, cut in decayStringsWithCuts:
960 if not decayDescriptor.init(decayString):
961 raise ValueError("Invalid decay string")
962 # need to check some logic to unpack possible scenarios
963 if decayDescriptor.getNDaughters() > 0:
964 # ... then we have an actual decay in the decay string which must be a V0 (if more than 1 daughter)
965 # or a kink (if 1 daughter)
966 # the particle loader automatically calls this "V0" or "kink", respectively, so we have to copy over
967 # the list to name/format that user wants
968 if (decayDescriptor.getNDaughters() == 1) and (decayDescriptor.getMother().getLabel() != 'kink'):
969 copyList(decayDescriptor.getMother().getFullName(), decayDescriptor.getMother().getName() + ':kink',
970 writeOut, path)
971 if (decayDescriptor.getNDaughters() > 1) and (decayDescriptor.getMother().getLabel() != 'V0'):
972 copyList(decayDescriptor.getMother().getFullName(), decayDescriptor.getMother().getName() + ':V0', writeOut, path)
973 elif (decayDescriptor.getMother().getLabel() != 'all' and
974 abs(decayDescriptor.getMother().getPDGCode()) != Belle2.Const.neutron.getPDGCode()):
975 # then we have a non-V0/kink particle which the particle loader automatically calls "all"
976 # as with the special V0 and kink cases we have to copy over the list to the name/format requested
977 copyList(decayString, decayDescriptor.getMother().getName() + ':all', writeOut, path)
978
979 # optionally apply a cut
980 if cut != "":
981 applyCuts(decayDescriptor.getMother().getFullName(), cut, path)
982
983 if decayString.startswith("gamma"):
984 # keep KLM-source photons as a experts-only for now: they are loaded by the particle loader,
985 # but the user has to explicitly request them.
986 if not loadPhotonsFromKLM:
987 applyCuts(decayString, 'isFromECL', path)
988
989
990def fillParticleList(decayString, cut, writeOut=False, path=None, enforceFitHypothesis=False,
991 loadPhotonsFromKLM=False):
992 """
993 Creates Particles of the desired type from the corresponding ``mdst`` dataobjects,
994 loads them to the StoreArray<Particle> and fills the ParticleList.
995
996 See also:
997 the :doc:`StandardParticles` functions.
998
999 The type of the particles to be loaded is specified via the decayString module parameter.
1000 The type of the ``mdst`` dataobject that is used as an input is determined from the type of
1001 the particle. The following types of the particles can be loaded:
1002
1003 * charged final state particles (input ``mdst`` type = Tracks)
1004 - e+, mu+, pi+, K+, p, deuteron (and charge conjugated particles)
1005
1006 * neutral final state particles
1007 - "gamma" (input ``mdst`` type = ECLCluster)
1008 - "K_S0", "Lambda0" (input ``mdst`` type = V0)
1009 - "K_L0", "n0" (input ``mdst`` type = KLMCluster or ECLCluster)
1010
1011 Note:
1012 For "K_S0" and "Lambda0" you must specify the daughter ordering.
1013
1014 For example, to load V0s as :math:`\\Lambda^0\\to p^+\\pi^-` decays from V0s:
1015
1016 .. code-block:: python
1017
1018 fillParticleList('Lambda0 -> p+ pi-', '0.9 < M < 1.3', path=mypath)
1019
1020 Tip:
1021 Gammas can also be loaded from KLMClusters by explicitly setting the
1022 parameter ``loadPhotonsFromKLM`` to True. However, this should only be
1023 done in selected use-cases and the effect should be studied carefully.
1024
1025 Tip:
1026 For "K_L0" it is now possible to load from ECLClusters, to revert to
1027 the old (Belle) behavior, you can require ``'isFromKLM > 0'``.
1028
1029 .. code-block:: python
1030
1031 fillParticleList('K_L0', 'isFromKLM > 0', path=mypath)
1032
1033 * Charged kinks final state particles (input ``mdst`` type = Kink)
1034
1035 .. note::
1036 To reconstruct charged particle kink you must specify the daughter.
1037
1038 For example, to load Kinks as :math:`K^- \\to \\pi^-\\pi^0` decays from Kinks:
1039
1040 .. code-block:: python
1041
1042 fillParticleList('K- -> pi-', yourCut, path=mypath)
1043
1044
1045 Parameters:
1046 decayString (str): Type of Particle and determines the name of the ParticleList.
1047 If the input MDST type is V0 the whole decay chain needs to be specified, so that
1048 the user decides and controls the daughters' order (e.g. ``K_S0 -> pi+ pi-``).
1049 If the input MDST type is Kink the decay chain needs to be specified
1050 with only one daughter (e.g. ``K- -> pi-``).
1051 cut (str): Particles need to pass these selection criteria to be added to the ParticleList
1052 writeOut (bool): whether RootOutput module should save the created ParticleList
1053 path (basf2.Path): modules are added to this path
1054 enforceFitHypothesis (bool): If true, Particles will be created only for the tracks which have been fitted
1055 using a mass hypothesis of the exact type passed to fillParticleLists().
1056 If enforceFitHypothesis is False (the default) the next closest fit hypothesis
1057 in terms of mass difference will be used if the fit using exact particle
1058 type is not available.
1059 loadPhotonsFromKLM (bool): If true, photon candidates will be created from KLMClusters as well.
1060 """
1061
1062 pload = register_module('ParticleLoader')
1063 pload.set_name('ParticleLoader_' + decayString)
1064 pload.param('decayStrings', [decayString])
1065 pload.param('writeOut', writeOut)
1066 pload.param("enforceFitHypothesis", enforceFitHypothesis)
1067 path.add_module(pload)
1068
1069 # need to check some logic to unpack possible scenarios
1070 from ROOT import Belle2
1071 decayDescriptor = Belle2.DecayDescriptor()
1072 if not decayDescriptor.init(decayString):
1073 raise ValueError("Invalid decay string")
1074 if decayDescriptor.getNDaughters() > 0:
1075 # ... then we have an actual decay in the decay string which must be a V0 (if more than 1 daughter)
1076 # or a kink (if 1 daughter)
1077 # the particle loader automatically calls this "V0" or "kink", respectively, so we have to copy over
1078 # the list to name/format that user wants
1079 if (decayDescriptor.getNDaughters() == 1) and (decayDescriptor.getMother().getLabel() != 'kink'):
1080 copyList(decayDescriptor.getMother().getFullName(), decayDescriptor.getMother().getName() + ':kink',
1081 writeOut, path)
1082 if (decayDescriptor.getNDaughters() > 1) and (decayDescriptor.getMother().getLabel() != 'V0'):
1083 copyList(decayDescriptor.getMother().getFullName(), decayDescriptor.getMother().getName() + ':V0', writeOut,
1084 path)
1085 elif (decayDescriptor.getMother().getLabel() != 'all' and
1086 abs(decayDescriptor.getMother().getPDGCode()) != Belle2.Const.neutron.getPDGCode()):
1087 # then we have a non-V0/kink particle which the particle loader automatically calls "all"
1088 # as with the special V0 and kink cases we have to copy over the list to the name/format requested
1089 copyList(decayString, decayDescriptor.getMother().getName() + ':all', writeOut, path)
1090
1091 # optionally apply a cut
1092 if cut != "":
1093 applyCuts(decayDescriptor.getMother().getFullName(), cut, path)
1094
1095 if decayString.startswith("gamma"):
1096 # keep KLM-source photons as a experts-only for now: they are loaded by the particle loader,
1097 # but the user has to explicitly request them.
1098 if not loadPhotonsFromKLM:
1099 applyCuts(decayString, 'isFromECL', path)
1100
1101
1102def fillParticleListWithTrackHypothesis(decayString,
1103 cut,
1104 hypothesis,
1105 writeOut=False,
1106 enforceFitHypothesis=False,
1107 path=None):
1108 """
1109 As fillParticleList, but if used for a charged FSP, loads the particle with the requested hypothesis if available
1110
1111 @param decayString specifies type of Particles and determines the name of the ParticleList
1112 @param cut Particles need to pass these selection criteria to be added to the ParticleList
1113 @param hypothesis the PDG code of the desired track hypothesis
1114 @param writeOut whether RootOutput module should save the created ParticleList
1115 @param enforceFitHypothesis If true, Particles will be created only for the tracks which have been fitted
1116 using a mass hypothesis of the exact type passed to fillParticleLists().
1117 If enforceFitHypothesis is False (the default) the next closest fit hypothesis
1118 in terms of mass difference will be used if the fit using exact particle
1119 type is not available.
1120 @param path modules are added to this path
1121 """
1122
1123 pload = register_module('ParticleLoader')
1124 pload.set_name('ParticleLoader_' + decayString)
1125 pload.param('decayStrings', [decayString])
1126 pload.param('trackHypothesis', hypothesis)
1127 pload.param('writeOut', writeOut)
1128 pload.param("enforceFitHypothesis", enforceFitHypothesis)
1129 path.add_module(pload)
1130
1131 from ROOT import Belle2
1132 decayDescriptor = Belle2.DecayDescriptor()
1133 if not decayDescriptor.init(decayString):
1134 raise ValueError("Invalid decay string")
1135 if decayDescriptor.getMother().getLabel() != 'all':
1136 # the particle loader automatically calls particle lists of charged FSPs "all"
1137 # so we have to copy over the list to the name/format requested
1138 copyList(decayString, decayDescriptor.getMother().getName() + ':all', writeOut, path)
1139
1140 # apply a cut if a non-empty cut string is provided
1141 if cut != "":
1142 applyCuts(decayString, cut, path)
1143
1144
1145def fillConvertedPhotonsList(decayString, cut, writeOut=False, path=None):
1146 """
1147 Creates photon Particle object for each e+e- combination in the V0 StoreArray.
1148
1149 Note:
1150 You must specify the daughter ordering.
1151
1152 .. code-block:: python
1153
1154 fillConvertedPhotonsList('gamma:converted -> e+ e-', '', path=mypath)
1155
1156 Parameters:
1157 decayString (str): Must be gamma to an e+e- pair. You must specify the daughter ordering.
1158 Will also determine the name of the particleList.
1159 cut (str): Particles need to pass these selection criteria to be added to the ParticleList
1160 writeOut (bool): whether RootOutput module should save the created ParticleList
1161 path (basf2.Path): modules are added to this path
1162
1163 """
1164
1165 import b2bii
1166 if b2bii.isB2BII():
1167 B2ERROR('For Belle converted photons are available in the pre-defined list "gamma:v0mdst".')
1168
1169 pload = register_module('ParticleLoader')
1170 pload.set_name('ParticleLoader_' + decayString)
1171 pload.param('decayStrings', [decayString])
1172 pload.param('addDaughters', True)
1173 pload.param('writeOut', writeOut)
1174 path.add_module(pload)
1175
1176 from ROOT import Belle2
1177 decayDescriptor = Belle2.DecayDescriptor()
1178 if not decayDescriptor.init(decayString):
1179 raise ValueError("Invalid decay string")
1180 if decayDescriptor.getMother().getLabel() != 'V0':
1181 # the particle loader automatically calls converted photons "V0" so we have to copy over
1182 # the list to name/format that user wants
1183 copyList(decayDescriptor.getMother().getFullName(), decayDescriptor.getMother().getName() + ':V0', writeOut, path)
1184
1185 # apply a cut if a non-empty cut string is provided
1186 if cut != "":
1187 applyCuts(decayDescriptor.getMother().getFullName(), cut, path)
1188
1189
1190def fillParticleListFromROE(decayString,
1191 cut,
1192 maskName='all',
1193 sourceParticleListName='',
1194 useMissing=False,
1195 writeOut=False,
1196 path=None):
1197 """
1198 Creates Particle object for each ROE of the desired type found in the
1199 StoreArray<RestOfEvent>, loads them to the StoreArray<Particle>
1200 and fills the ParticleList. If useMissing is True, then the missing
1201 momentum is used instead of ROE.
1202
1203 The type of the particles to be loaded is specified via the decayString module parameter.
1204
1205 @param decayString specifies type of Particles and determines the name of the ParticleList.
1206 Source ROEs can be taken as a daughter list, for example:
1207 'B0:tagFromROE -> B0:signal'
1208 @param cut Particles need to pass these selection criteria to be added to the ParticleList
1209 @param maskName Name of the ROE mask to use
1210 @param sourceParticleListName Use related ROEs to this particle list as a source
1211 @param useMissing Use missing momentum instead of ROE momentum
1212 @param writeOut whether RootOutput module should save the created ParticleList
1213 @param path modules are added to this path
1214 """
1215
1216 pload = register_module('ParticleLoader')
1217 pload.set_name('ParticleLoader_' + decayString)
1218 pload.param('decayStrings', [decayString])
1219 pload.param('writeOut', writeOut)
1220 pload.param('roeMaskName', maskName)
1221 pload.param('useMissing', useMissing)
1222 pload.param('sourceParticleListName', sourceParticleListName)
1223 pload.param('useROEs', True)
1224 path.add_module(pload)
1225
1226 from ROOT import Belle2
1227 decayDescriptor = Belle2.DecayDescriptor()
1228 if not decayDescriptor.init(decayString):
1229 raise ValueError("Invalid decay string")
1230
1231 # apply a cut if a non-empty cut string is provided
1232 if cut != "":
1233 applyCuts(decayDescriptor.getMother().getFullName(), cut, path)
1234
1235
1236def fillParticleListFromDummy(decayString,
1237 mdstIndex=0,
1238 covMatrix=10000.,
1239 treatAsInvisible=True,
1240 writeOut=False,
1241 path=None):
1242 """
1243 Creates a ParticleList and fills it with dummy Particles. For self-conjugated Particles one dummy
1244 Particle is created, for Particles that are not self-conjugated one Particle and one anti-Particle is
1245 created. The four-momentum is set to zero.
1246
1247 The type of the particles to be loaded is specified via the decayString module parameter.
1248
1249 @param decayString specifies type of Particles and determines the name of the ParticleList
1250 @param mdstIndex sets the mdst index of Particles
1251 @param covMatrix sets the value of the diagonal covariance matrix of Particles
1252 @param treatAsInvisible whether treeFitter should treat the Particles as invisible
1253 @param writeOut whether RootOutput module should save the created ParticleList
1254 @param path modules are added to this path
1255 """
1256
1257 pload = register_module('ParticleLoader')
1258 pload.set_name('ParticleLoader_' + decayString)
1259 pload.param('decayStrings', [decayString])
1260 pload.param('useDummy', True)
1261 pload.param('dummyMDSTIndex', mdstIndex)
1262 pload.param('dummyCovMatrix', covMatrix)
1263 pload.param('dummyTreatAsInvisible', treatAsInvisible)
1264 pload.param('writeOut', writeOut)
1265 path.add_module(pload)
1266
1267
1268def fillParticleListFromMC(decayString,
1269 cut,
1270 addDaughters=False,
1271 skipNonPrimaryDaughters=False,
1272 writeOut=False,
1273 path=None,
1274 skipNonPrimary=False,
1275 skipInitial=True):
1276 """
1277 Creates Particle object for each MCParticle of the desired type found in the StoreArray<MCParticle>,
1278 loads them to the StoreArray<Particle> and fills the ParticleList.
1279
1280 The type of the particles to be loaded is specified via the decayString module parameter.
1281
1282 @param decayString specifies type of Particles and determines the name of the ParticleList
1283 @param cut Particles need to pass these selection criteria to be added to the ParticleList
1284 @param addDaughters adds the bottom part of the decay chain of the particle to the datastore and
1285 sets mother-daughter relations
1286 @param skipNonPrimaryDaughters if true, skip non primary daughters, useful to study final state daughter particles
1287 @param writeOut whether RootOutput module should save the created ParticleList
1288 @param path modules are added to this path
1289 @param skipNonPrimary if true, skip non primary particle
1290 @param skipInitial if true, skip initial particles
1291 """
1292
1293 pload = register_module('ParticleLoader')
1294 pload.set_name('ParticleLoader_' + decayString)
1295 pload.param('decayStrings', [decayString])
1296 pload.param('addDaughters', addDaughters)
1297 pload.param('skipNonPrimaryDaughters', skipNonPrimaryDaughters)
1298 pload.param('writeOut', writeOut)
1299 pload.param('useMCParticles', True)
1300 pload.param('skipNonPrimary', skipNonPrimary)
1301 pload.param('skipInitial', skipInitial)
1302 path.add_module(pload)
1303
1304 from ROOT import Belle2
1305 decayDescriptor = Belle2.DecayDescriptor()
1306 if not decayDescriptor.init(decayString):
1307 raise ValueError("Invalid decay string")
1308
1309 # apply a cut if a non-empty cut string is provided
1310 if cut != "":
1311 applyCuts(decayString, cut, path)
1312
1313
1314def fillParticleListsFromMC(decayStringsWithCuts,
1315 addDaughters=False,
1316 skipNonPrimaryDaughters=False,
1317 writeOut=False,
1318 path=None,
1319 skipNonPrimary=False,
1320 skipInitial=True):
1321 """
1322 Creates Particle object for each MCParticle of the desired type found in the StoreArray<MCParticle>,
1323 loads them to the StoreArray<Particle> and fills the ParticleLists.
1324
1325 The types of the particles to be loaded are specified via the (decayString, cut) tuples given in a list.
1326 For example:
1327
1328 .. code-block:: python
1329
1330 kaons = ('K+:gen', '')
1331 pions = ('pi+:gen', 'pionID>0.1')
1332 fillParticleListsFromMC([kaons, pions], path=mypath)
1333
1334 .. tip::
1335 Daughters of ``Lambda0`` are not primary, but ``Lambda0`` is not final state particle.
1336 Thus, when one reconstructs a particle from ``Lambda0``, that is created with
1337 ``addDaughters=True`` and ``skipNonPrimaryDaughters=True``, the particle always has ``isSignal==0``.
1338 Please set options for ``Lambda0`` to use MC-matching variables properly as follows,
1339 ``addDaughters=True`` and ``skipNonPrimaryDaughters=False``.
1340
1341 @param decayString specifies type of Particles and determines the name of the ParticleList
1342 @param cut Particles need to pass these selection criteria to be added to the ParticleList
1343 @param addDaughters adds the bottom part of the decay chain of the particle to the datastore and
1344 sets mother-daughter relations
1345 @param skipNonPrimaryDaughters if true, skip non primary daughters, useful to study final state daughter particles
1346 @param writeOut whether RootOutput module should save the created ParticleList
1347 @param path modules are added to this path
1348 @param skipNonPrimary if true, skip non primary particle
1349 @param skipInitial if true, skip initial particles
1350 """
1351
1352 pload = register_module('ParticleLoader')
1353 pload.set_name('ParticleLoader_' + 'PLists')
1354 pload.param('decayStrings', [decayString for decayString, cut in decayStringsWithCuts])
1355 pload.param('addDaughters', addDaughters)
1356 pload.param('skipNonPrimaryDaughters', skipNonPrimaryDaughters)
1357 pload.param('writeOut', writeOut)
1358 pload.param('useMCParticles', True)
1359 pload.param('skipNonPrimary', skipNonPrimary)
1360 pload.param('skipInitial', skipInitial)
1361 path.add_module(pload)
1362
1363 from ROOT import Belle2
1364 decayDescriptor = Belle2.DecayDescriptor()
1365 for decayString, cut in decayStringsWithCuts:
1366 if not decayDescriptor.init(decayString):
1367 raise ValueError("Invalid decay string")
1368
1369 # apply a cut if a non-empty cut string is provided
1370 if cut != "":
1371 applyCuts(decayString, cut, path)
1372
1373
1374def fillParticleListFromChargedCluster(outputParticleList,
1375 inputParticleList,
1376 cut,
1377 useOnlyMostEnergeticECLCluster=True,
1378 writeOut=False,
1379 path=None):
1380 """
1381 Creates the Particle object from ECLCluster and KLMCluster that are being matched with the Track of inputParticleList.
1382
1383 @param outputParticleList The output ParticleList. Only neutral final state particles are supported.
1384 @param inputParticleList The input ParticleList that is required to have the relation to the Track object.
1385 @param cut Particles need to pass these selection criteria to be added to the ParticleList
1386 @param useOnlyMostEnergeticECLCluster If True, only the most energetic ECLCluster among ones that are matched with the Track is
1387 used. If False, all matched ECLClusters are loaded. The default is True. Regardless of
1388 this option, the KLMCluster is loaded.
1389 @param writeOut whether RootOutput module should save the created ParticleList
1390 @param path modules are added to this path
1391 """
1392
1393 pload = register_module('ParticleLoader')
1394 pload.set_name('ParticleLoader_' + outputParticleList)
1395
1396 pload.param('decayStrings', [outputParticleList])
1397 pload.param('sourceParticleListName', inputParticleList)
1398 pload.param('writeOut', writeOut)
1399 pload.param('loadChargedCluster', True)
1400 pload.param('useOnlyMostEnergeticECLCluster', useOnlyMostEnergeticECLCluster)
1401 path.add_module(pload)
1402
1403 # apply a cut if a non-empty cut string is provided
1404 if cut != "":
1405 applyCuts(outputParticleList, cut, path)
1406
1407
1408def extractParticlesFromROE(particleLists,
1409 signalSideParticleList=None,
1410 maskName='all',
1411 writeOut=False,
1412 path=None):
1413 """
1414 Extract Particle objects that belong to the Rest-Of-Events and fill them into the ParticleLists.
1415 The types of the particles other than those specified by ``particleLists`` are not stored.
1416 If one creates a ROE with ``fillWithMostLikely=True`` via `buildRestOfEvent`, for example,
1417 one should create particleLists for not only ``pi+``, ``gamma``, ``K_L0`` but also other charged final state particles.
1418
1419 When one calls the function in the main path, one has to set the argument ``signalSideParticleList`` and the signal side
1420 ParticleList must have only one candidate.
1421
1422 .. code-block:: python
1423
1424 buildRestOfEvent('B0:sig', fillWithMostLikely=True, path=mypath)
1425
1426 roe_path = create_path()
1427 deadEndPath = create_path()
1428 signalSideParticleFilter('B0:sig', '', roe_path, deadEndPath)
1429
1430 plists = ['%s:in_roe' % ptype for ptype in ['pi+', 'gamma', 'K_L0', 'K+', 'p+', 'e+', 'mu+']]
1431 extractParticlesFromROE(plists, maskName='all', path=roe_path)
1432
1433 # one can analyze these ParticleLists in the roe_path
1434
1435 mypath.for_each('RestOfEvent', 'RestOfEvents', roe_path)
1436
1437 rankByLowest('B0:sig', 'deltaE', numBest=1, path=mypath)
1438 extractParticlesFromROE(plists, signalSideParticleList='B0:sig', maskName='all', path=mypath)
1439
1440 # one can analyze these ParticleLists in the main path
1441
1442
1443 @param particleLists (str or list(str)) Name of output ParticleLists
1444 @param signalSideParticleList (str) Name of signal side ParticleList
1445 @param maskName (str) Name of the ROE mask to be applied on Particles
1446 @param writeOut (bool) whether RootOutput module should save the created ParticleList
1447 @param path (basf2.Path) modules are added to this path
1448 """
1449
1450 if isinstance(particleLists, str):
1451 particleLists = [particleLists]
1452
1453 pext = register_module('ParticleExtractorFromROE')
1454 pext.set_name('ParticleExtractorFromROE_' + '_'.join(particleLists))
1455 pext.param('outputListNames', particleLists)
1456 if signalSideParticleList is not None:
1457 pext.param('signalSideParticleListName', signalSideParticleList)
1458 pext.param('maskName', maskName)
1459 pext.param('writeOut', writeOut)
1460 path.add_module(pext)
1461
1462
1463def applyCuts(list_name, cut, path):
1464 """
1465 Removes particle candidates from ``list_name`` that do not pass ``cut``
1466 (given selection criteria).
1467
1468 Example:
1469 require energetic pions safely inside the cdc
1470
1471 .. code-block:: python
1472
1473 applyCuts("pi+:mypions", "[E > 2] and thetaInCDCAcceptance", path=mypath)
1474
1475 Warning:
1476 You must use square braces ``[`` and ``]`` for conditional statements.
1477
1478 Parameters:
1479 list_name (str): input ParticleList name
1480 cut (str): Candidates that do not pass these selection criteria are removed from the ParticleList
1481 path (basf2.Path): modules are added to this path
1482 """
1483
1484 pselect = register_module('ParticleSelector')
1485 pselect.set_name('ParticleSelector_applyCuts_' + list_name)
1486 pselect.param('decayString', list_name)
1487 pselect.param('cut', cut)
1488 path.add_module(pselect)
1489
1490
1491def applyEventCuts(cut, path, metavariables=None):
1492 """
1493 Removes events that do not pass the ``cut`` (given selection criteria).
1494
1495 Example:
1496 continuum events (in mc only) with more than 5 tracks
1497
1498 .. code-block:: python
1499
1500 applyEventCuts("[nTracks > 5] and [isContinuumEvent], path=mypath)
1501
1502 .. warning::
1503 Only event-based variables are allowed in this function
1504 and only square brackets ``[`` and ``]`` for conditional statements.
1505
1506 Parameters:
1507 cut (str): Events that do not pass these selection criteria are skipped
1508 path (basf2.Path): modules are added to this path
1509 metavariables (list(str)): List of meta variables to be considered in decomposition of cut
1510 """
1511
1512 import b2parser
1513 from variables import variables
1514
1515 def find_vars(t: tuple, var_list: list, meta_list: list) -> None:
1516 """ Recursive helper function to find variable names """
1517 if not isinstance(t, tuple):
1518 return
1519 if t[0] == b2parser.B2ExpressionParser.node_types['IdentifierNode']:
1520 var_list += [t[1]]
1521 return
1522 if t[0] == b2parser.B2ExpressionParser.node_types['FunctionNode']:
1523 meta_list.append(list(t[1:]))
1524 return
1525 for i in t:
1526 if isinstance(i, tuple):
1527 find_vars(i, var_list, meta_list)
1528
1529 def check_variable(var_list: list, metavar_ids: list) -> None:
1530 for var_string in var_list:
1531 # Check if the var_string is alias
1532 orig_name = variables.resolveAlias(var_string)
1533 if orig_name != var_string:
1534 var_list_temp = []
1535 meta_list_temp = []
1536 find_vars(b2parser.parse(orig_name), var_list_temp, meta_list_temp)
1537
1538 check_variable(var_list_temp, metavar_ids)
1539 check_meta(meta_list_temp, metavar_ids)
1540 else:
1541 # Get the variable
1542 var = variables.getVariable(var_string)
1543 if event_var_id not in var.description:
1544 B2ERROR(f'Variable {var_string} is not an event-based variable! "\
1545 "Please check your inputs to the applyEventCuts method!')
1546
1547 def check_meta(meta_list: list, metavar_ids: list) -> None:
1548 for meta_string_list in meta_list:
1549 var_list_temp = []
1550 while meta_string_list[0] in metavar_ids:
1551 # remove special meta variable
1552 meta_string_list.pop(0)
1553 for meta_string in meta_string_list[0].split(","):
1554 find_vars(b2parser.parse(meta_string), var_list_temp, meta_string_list)
1555 if len(meta_string_list) > 0:
1556 meta_string_list.pop(0)
1557 if len(meta_string_list) == 0:
1558 break
1559 if len(meta_string_list) > 1:
1560 meta_list += meta_string_list[1:]
1561 if isinstance(meta_string_list[0], list):
1562 meta_string_list = [element for element in meta_string_list[0]]
1563
1564 check_variable(var_list_temp, metavar_ids)
1565
1566 if len(meta_string_list) == 0:
1567 continue
1568 elif len(meta_string_list) == 1:
1569 var = variables.getVariable(meta_string_list[0])
1570 else:
1571 var = variables.getVariable(meta_string_list[0], meta_string_list[1].split(","))
1572 # Check if the variable's description contains event-based marker
1573 if event_var_id in var.description:
1574 continue
1575 # Throw an error message if non event-based variable is used
1576 B2ERROR(f'Variable {var.name} is not an event-based variable! Please check your inputs to the applyEventCuts method!')
1577
1578 event_var_id = '[Eventbased]'
1579 metavar_ids = ['formula', 'abs',
1580 'cos', 'acos',
1581 'tan', 'atan',
1582 'sin', 'asin',
1583 'exp', 'log', 'log10',
1584 'min', 'max',
1585 'isNAN', 'ifNANgiveX']
1586 if metavariables:
1587 metavar_ids += metavariables
1588
1589 var_list = []
1590 meta_list = []
1591 find_vars(b2parser.parse(cut), var_list=var_list, meta_list=meta_list)
1592
1593 if len(var_list) == 0 and len(meta_list) == 0:
1594 B2WARNING(f'Cut string "{cut}" has no variables for applyEventCuts helper function!')
1595
1596 check_variable(var_list, metavar_ids)
1597 check_meta(meta_list, metavar_ids)
1598
1599 eselect = register_module('VariableToReturnValue')
1600 eselect.param('variable', 'passesEventCut(' + cut + ')')
1601 path.add_module(eselect)
1602 empty_path = create_path()
1603 eselect.if_value('<1', empty_path)
1604
1605
1606def reconstructDecay(decayString,
1607 cut,
1608 dmID=0,
1609 writeOut=False,
1610 path=None,
1611 candidate_limit=None,
1612 ignoreIfTooManyCandidates=True,
1613 chargeConjugation=True,
1614 allowChargeViolation=False):
1615 r"""
1616 Creates new Particles by making combinations of existing Particles - it reconstructs unstable particles via their specified
1617 decay mode, e.g. in form of a :ref:`DecayString`: :code:`D0 -> K- pi+` or :code:`B+ -> anti-D0 pi+`, ... All possible
1618 combinations are created (particles are used only once per candidate) and combinations that pass the specified selection
1619 criteria are saved to a newly created (mother) ParticleList. By default the charge conjugated decay is reconstructed as well
1620 (meaning that the charge conjugated mother list is created as well) but this can be deactivated.
1621
1622 One can use an ``@``-sign to mark a particle as unspecified for inclusive analyses,
1623 e.g. in a DecayString: :code:`'@Xsd -> K+ pi-'`.
1624
1625 .. seealso:: :ref:`Marker_of_unspecified_particle`
1626
1627 .. warning::
1628 The input ParticleLists are typically ordered according to the upstream reconstruction algorithm.
1629 Therefore, if you combine two or more identical particles in the decay chain you should not expect to see the same
1630 distribution for the daughter kinematics as they may be sorted by geometry, momentum etc.
1631
1632 For example, in the decay :code:`D0 -> pi0 pi0` the momentum distributions of the two ``pi0`` s are not identical.
1633 This can be solved by manually randomising the lists before combining.
1634
1635 See Also:
1636
1637 * `Particle combiner how does it work? <https://questions.belle2.org/question/4318/particle-combiner-how-does-it-work/>`_
1638 * `Identical particles in decay chain <https://questions.belle2.org/question/5724/identical-particles-in-decay-chain/>`_
1639
1640 @param decayString :ref:`DecayString` specifying what kind of the decay should be reconstructed
1641 (from the DecayString the mother and daughter ParticleLists are determined)
1642 @param cut created (mother) Particles are added to the mother ParticleList if they
1643 pass give cuts (in VariableManager style) and rejected otherwise
1644 @param dmID user specified decay mode identifier
1645 @param writeOut whether RootOutput module should save the created ParticleList
1646 @param path modules are added to this path
1647 @param candidate_limit Maximum amount of candidates to be reconstructed. If
1648 the number of candidates is exceeded a Warning will be
1649 printed.
1650 By default, all these candidates will be removed and event will be ignored.
1651 This behaviour can be changed by \'ignoreIfTooManyCandidates\' flag.
1652 If no value is given the amount is limited to a sensible
1653 default. A value <=0 will disable this limit and can
1654 cause huge memory amounts so be careful.
1655 @param ignoreIfTooManyCandidates whether event should be ignored or not if number of reconstructed
1656 candidates reaches limit. If event is ignored, no candidates are reconstructed,
1657 otherwise, number of candidates in candidate_limit is reconstructed.
1658 @param chargeConjugation boolean to decide whether charge conjugated mode should be reconstructed as well (on by default)
1659 @param allowChargeViolation whether the decay string needs to conserve the electric charge
1660 """
1661
1662 pmake = register_module('ParticleCombiner')
1663 pmake.set_name('ParticleCombiner_' + decayString)
1664 pmake.param('decayString', decayString)
1665 pmake.param('cut', cut)
1666 pmake.param('decayMode', dmID)
1667 pmake.param('writeOut', writeOut)
1668 if candidate_limit is not None:
1669 pmake.param("maximumNumberOfCandidates", candidate_limit)
1670 pmake.param("ignoreIfTooManyCandidates", ignoreIfTooManyCandidates)
1671 pmake.param('chargeConjugation', chargeConjugation)
1672 pmake.param("allowChargeViolation", allowChargeViolation)
1673 path.add_module(pmake)
1674
1675
1676def combineAllParticles(inputParticleLists, outputList, cut='', writeOut=False, path=None):
1677 """
1678 Creates a new Particle as the combination of all Particles from all
1679 provided inputParticleLists. However, each particle is used only once
1680 (even if duplicates are provided) and the combination has to pass the
1681 specified selection criteria to be saved in the newly created (mother)
1682 ParticleList.
1683
1684 @param inputParticleLists List of input particle lists which are combined to the new Particle
1685 @param outputList Name of the particle combination created with this module
1686 @param cut created (mother) Particle is added to the mother ParticleList if it passes
1687 these given cuts (in VariableManager style) and is rejected otherwise
1688 @param writeOut whether RootOutput module should save the created ParticleList
1689 @param path module is added to this path
1690 """
1691
1692 pmake = register_module('AllParticleCombiner')
1693 pmake.set_name('AllParticleCombiner_' + outputList)
1694 pmake.param('inputListNames', inputParticleLists)
1695 pmake.param('outputListName', outputList)
1696 pmake.param('cut', cut)
1697 pmake.param('writeOut', writeOut)
1698 path.add_module(pmake)
1699
1700
1701def reconstructMissingKlongDecayExpert(decayString,
1702 cut,
1703 dmID=0,
1704 writeOut=False,
1705 path=None,
1706 recoList="_reco"):
1707 """
1708 Creates a list of K_L0's and of B -> K_L0 + X, with X being a fully-reconstructed state.
1709 The K_L0 momentum is determined from kinematic constraints of the two-body B decay into K_L0 and X
1710
1711 @param decayString DecayString specifying what kind of the decay should be reconstructed
1712 (from the DecayString the mother and daughter ParticleLists are determined)
1713 @param cut Particles are added to the K_L0 and B ParticleList if the B candidates
1714 pass the given cuts (in VariableManager style) and rejected otherwise
1715 @param dmID user specified decay mode identifier
1716 @param writeOut whether RootOutput module should save the created ParticleList
1717 @param path modules are added to this path
1718 @param recoList suffix appended to original K_L0 and B ParticleList that identify the newly created K_L0 and B lists
1719 """
1720
1721 pcalc = register_module('KlongMomentumCalculatorExpert')
1722 pcalc.set_name('KlongMomentumCalculatorExpert_' + decayString)
1723 pcalc.param('decayString', decayString)
1724 pcalc.param('writeOut', writeOut)
1725 pcalc.param('recoList', recoList)
1726 path.add_module(pcalc)
1727
1728 rmake = register_module('KlongDecayReconstructorExpert')
1729 rmake.set_name('KlongDecayReconstructorExpert_' + decayString)
1730 rmake.param('decayString', decayString)
1731 rmake.param('cut', cut)
1732 rmake.param('decayMode', dmID)
1733 rmake.param('writeOut', writeOut)
1734 rmake.param('recoList', recoList)
1735 path.add_module(rmake)
1736
1737
1738def setBeamConstrainedMomentum(particleList, decayStringTarget, decayStringDaughters, path=None):
1739 """
1740 Replace the four-momentum of the target Particle by p(beam) - p(selected daughters).
1741 The momentum of the mother Particle will not be changed.
1742
1743 @param particleList mother Particlelist
1744 @param decayStringTarget DecayString specifying the target particle whose momentum
1745 will be updated
1746 @param decayStringDaughters DecayString specifying the daughter particles used to replace
1747 the momentum of the target particle by p(beam)-p(daughters)
1748 """
1749
1750 mod = register_module('ParticleMomentumUpdater')
1751 mod.set_name('ParticleMomentumUpdater' + particleList)
1752 mod.param('particleList', particleList)
1753 mod.param('decayStringTarget', decayStringTarget)
1754 mod.param('decayStringDaughters', decayStringDaughters)
1755 path.add_module(mod)
1756
1757
1758def updateKlongKinematicsExpert(particleList,
1759 writeOut=False,
1760 path=None):
1761 """
1762 Calculates and updates the kinematics of B->K_L0 + something else with same method as
1763 `reconstructMissingKlongDecayExpert`. This helps to revert the kinematics after the vertex fitting.
1764
1765 @param particleList input ParticleList of B meson that decays to K_L0 + X
1766 @param writeOut whether RootOutput module should save the ParticleList
1767 @param path modules are added to this path
1768 """
1769
1770 mod = register_module('KlongMomentumUpdaterExpert')
1771 mod.set_name('KlongMomentumUpdaterExpert_' + particleList)
1772 mod.param('listName', particleList)
1773 mod.param('writeOut', writeOut)
1774 path.add_module(mod)
1775
1776
1777def replaceMass(replacerName, particleLists=None, pdgCode=22, path=None):
1778 """
1779 replaces the mass of the particles inside the given particleLists
1780 with the invariant mass of the particle corresponding to the given pdgCode.
1781
1782 @param particleLists new ParticleList filled with copied Particles
1783 @param pdgCode PDG code for mass reference
1784 @param path modules are added to this path
1785 """
1786
1787 if particleLists is None:
1788 particleLists = []
1789
1790 # first copy original particles to the new ParticleList
1791 pmassupdater = register_module('ParticleMassUpdater')
1792 pmassupdater.set_name('ParticleMassUpdater_' + replacerName)
1793 pmassupdater.param('particleLists', particleLists)
1794 pmassupdater.param('pdgCode', pdgCode)
1795 path.add_module(pmassupdater)
1796
1797
1798def reconstructRecoil(decayString,
1799 cut,
1800 dmID=0,
1801 writeOut=False,
1802 path=None,
1803 candidate_limit=None,
1804 allowChargeViolation=False):
1805 """
1806 Creates new Particles that recoil against the input particles.
1807
1808 For example the decay string M -> D1 D2 D3 will:
1809 - create mother Particle M for each unique combination of D1, D2, D3 Particles
1810 - Particles D1, D2, D3 will be appended as daughters to M
1811 - the 4-momentum of the mother Particle M is given by
1812 p(M) = p(HER) + p(LER) - Sum_i p(Di)
1813
1814 @param decayString DecayString specifying what kind of the decay should be reconstructed
1815 (from the DecayString the mother and daughter ParticleLists are determined)
1816 @param cut created (mother) Particles are added to the mother ParticleList if they
1817 pass give cuts (in VariableManager style) and rejected otherwise
1818 @param dmID user specified decay mode identifier
1819 @param writeOut whether RootOutput module should save the created ParticleList
1820 @param path modules are added to this path
1821 @param candidate_limit Maximum amount of candidates to be reconstructed. If
1822 the number of candidates is exceeded no candidate will be
1823 reconstructed for that event and a Warning will be
1824 printed.
1825 If no value is given the amount is limited to a sensible
1826 default. A value <=0 will disable this limit and can
1827 cause huge memory amounts so be careful.
1828 @param allowChargeViolation whether the decay string needs to conserve the electric charge
1829 """
1830
1831 pmake = register_module('ParticleCombiner')
1832 pmake.set_name('ParticleCombiner_' + decayString)
1833 pmake.param('decayString', decayString)
1834 pmake.param('cut', cut)
1835 pmake.param('decayMode', dmID)
1836 pmake.param('writeOut', writeOut)
1837 pmake.param('recoilParticleType', 1)
1838 if candidate_limit is not None:
1839 pmake.param("maximumNumberOfCandidates", candidate_limit)
1840 pmake.param('allowChargeViolation', allowChargeViolation)
1841 path.add_module(pmake)
1842
1843
1844def reconstructRecoilDaughter(decayString,
1845 cut,
1846 dmID=0,
1847 writeOut=False,
1848 path=None,
1849 candidate_limit=None,
1850 allowChargeViolation=False):
1851 """
1852 Creates new Particles that are daughters of the particle reconstructed in the recoil (always assumed to be the first daughter).
1853
1854 For example the decay string M -> D1 D2 D3 will:
1855 - create mother Particle M for each unique combination of D1, D2, D3 Particles
1856 - Particles D1, D2, D3 will be appended as daughters to M
1857 - the 4-momentum of the mother Particle M is given by
1858 p(M) = p(D1) - Sum_i p(Di), where i>1
1859
1860 @param decayString DecayString specifying what kind of the decay should be reconstructed
1861 (from the DecayString the mother and daughter ParticleLists are determined)
1862 @param cut created (mother) Particles are added to the mother ParticleList if they
1863 pass give cuts (in VariableManager style) and rejected otherwise
1864 @param dmID user specified decay mode identifier
1865 @param writeOut whether RootOutput module should save the created ParticleList
1866 @param path modules are added to this path
1867 @param candidate_limit Maximum amount of candidates to be reconstructed. If
1868 the number of candidates is exceeded no candidate will be
1869 reconstructed for that event and a Warning will be
1870 printed.
1871 If no value is given the amount is limited to a sensible
1872 default. A value <=0 will disable this limit and can
1873 cause huge memory amounts so be careful.
1874 @param allowChargeViolation whether the decay string needs to conserve the electric charge taking into account that the first
1875 daughter is actually the mother
1876 """
1877
1878 pmake = register_module('ParticleCombiner')
1879 pmake.set_name('ParticleCombiner_' + decayString)
1880 pmake.param('decayString', decayString)
1881 pmake.param('cut', cut)
1882 pmake.param('decayMode', dmID)
1883 pmake.param('writeOut', writeOut)
1884 pmake.param('recoilParticleType', 2)
1885 if candidate_limit is not None:
1886 pmake.param("maximumNumberOfCandidates", candidate_limit)
1887 pmake.param('allowChargeViolation', allowChargeViolation)
1888 path.add_module(pmake)
1889
1890
1891def rankByHighest(particleList,
1892 variable,
1893 numBest=0,
1894 outputVariable='',
1895 allowMultiRank=False,
1896 cut='',
1897 overwriteRank=False,
1898 path=None):
1899 """
1900 Ranks particles in the input list by the given variable (highest to lowest), and stores an integer rank for each Particle
1901 in an :b2:var:`extraInfo` field ``${variable}_rank`` starting at 1 (best).
1902 The list is also sorted from best to worst candidate
1903 (each charge, e.g. B+/B-, separately).
1904 This can be used to perform a best candidate selection by cutting on the corresponding rank value, or by specifying
1905 a non-zero value for 'numBest'.
1906
1907 .. tip::
1908 Extra-info fields can be accessed by the :b2:var:`extraInfo` metavariable.
1909 These variable names can become clunky, so it's probably a good idea to set an alias.
1910 For example if you rank your B candidates by momentum,
1911
1912 .. code:: python
1913
1914 rankByHighest("B0:myCandidates", "p", path=mypath)
1915 vm.addAlias("momentumRank", "extraInfo(p_rank)")
1916
1917
1918 @param particleList The input ParticleList
1919 @param variable Variable to order Particles by.
1920 @param numBest If not zero, only the $numBest Particles in particleList with rank <= numBest are kept.
1921 @param outputVariable Name for the variable that will be created which contains the rank, Default is '${variable}_rank'.
1922 @param allowMultiRank If true, candidates with the same value will get the same rank.
1923 @param cut Only candidates passing the cut will be ranked. The others will have rank -1
1924 @param overwriteRank If true, the extraInfo of rank is overwritten when the particle has already the extraInfo.
1925 @param path modules are added to this path
1926 """
1927
1928 bcs = register_module('BestCandidateSelection')
1929 bcs.set_name('BestCandidateSelection_' + particleList + '_' + variable)
1930 bcs.param('particleList', particleList)
1931 bcs.param('variable', variable)
1932 bcs.param('numBest', numBest)
1933 bcs.param('outputVariable', outputVariable)
1934 bcs.param('allowMultiRank', allowMultiRank)
1935 bcs.param('cut', cut)
1936 bcs.param('overwriteRank', overwriteRank)
1937 path.add_module(bcs)
1938
1939
1940def rankByLowest(particleList,
1941 variable,
1942 numBest=0,
1943 outputVariable='',
1944 allowMultiRank=False,
1945 cut='',
1946 overwriteRank=False,
1947 path=None):
1948 """
1949 Ranks particles in the input list by the given variable (lowest to highest), and stores an integer rank for each Particle
1950 in an :b2:var:`extraInfo` field ``${variable}_rank`` starting at 1 (best).
1951 The list is also sorted from best to worst candidate
1952 (each charge, e.g. B+/B-, separately).
1953 This can be used to perform a best candidate selection by cutting on the corresponding rank value, or by specifying
1954 a non-zero value for 'numBest'.
1955
1956 .. tip::
1957 Extra-info fields can be accessed by the :b2:var:`extraInfo` metavariable.
1958 These variable names can become clunky, so it's probably a good idea to set an alias.
1959 For example if you rank your B candidates by :b2:var:`dM`,
1960
1961 .. code:: python
1962
1963 rankByLowest("B0:myCandidates", "dM", path=mypath)
1964 vm.addAlias("massDifferenceRank", "extraInfo(dM_rank)")
1965
1966
1967 @param particleList The input ParticleList
1968 @param variable Variable to order Particles by.
1969 @param numBest If not zero, only the $numBest Particles in particleList with rank <= numBest are kept.
1970 @param outputVariable Name for the variable that will be created which contains the rank, Default is '${variable}_rank'.
1971 @param allowMultiRank If true, candidates with the same value will get the same rank.
1972 @param cut Only candidates passing the cut will be ranked. The others will have rank -1
1973 @param overwriteRank If true, the extraInfo of rank is overwritten when the particle has already the extraInfo.
1974 @param path modules are added to this path
1975 """
1976
1977 bcs = register_module('BestCandidateSelection')
1978 bcs.set_name('BestCandidateSelection_' + particleList + '_' + variable)
1979 bcs.param('particleList', particleList)
1980 bcs.param('variable', variable)
1981 bcs.param('numBest', numBest)
1982 bcs.param('selectLowest', True)
1983 bcs.param('allowMultiRank', allowMultiRank)
1984 bcs.param('outputVariable', outputVariable)
1985 bcs.param('cut', cut)
1986 bcs.param('overwriteRank', overwriteRank)
1987 path.add_module(bcs)
1988
1989
1990def applyRandomCandidateSelection(particleList, path=None):
1991 """
1992 If there are multiple candidates in the provided particleList, all but one of them are removed randomly.
1993 This is done on a event-by-event basis.
1994
1995 @param particleList ParticleList for which the random candidate selection should be applied
1996 @param path module is added to this path
1997 """
1998
1999 rcs = register_module('BestCandidateSelection')
2000 rcs.set_name('RandomCandidateSelection_' + particleList)
2001 rcs.param('particleList', particleList)
2002 rcs.param('variable', 'random')
2003 rcs.param('selectLowest', False)
2004 rcs.param('allowMultiRank', False)
2005 rcs.param('numBest', 1)
2006 rcs.param('cut', '')
2007 rcs.param('outputVariable', '')
2008 path.add_module(rcs)
2009
2010
2011def printDataStore(eventNumber=-1, path=None):
2012 """
2013 Prints the contents of DataStore in the first event (or a specific event number or all events).
2014 Will list all objects and arrays (including size).
2015
2016 See also:
2017 The command line tool: ``b2file-size``.
2018
2019 Parameters:
2020 eventNumber (int): Print the datastore only for this event. The default
2021 (-1) prints only the first event, 0 means print for all events (can produce large output)
2022 path (basf2.Path): the PrintCollections module is added to this path
2023
2024 Warning:
2025 This will print a lot of output if you print it for all events and process many events.
2026
2027 """
2028
2029 printDS = register_module('PrintCollections')
2030 printDS.param('printForEvent', eventNumber)
2031 path.add_module(printDS)
2032
2033
2034def printVariableValues(list_name, var_names, path):
2035 """
2036 Prints out values of specified variables of all Particles included in given ParticleList. For debugging purposes.
2037
2038 @param list_name input ParticleList name
2039 @param var_names vector of variable names to be printed
2040 @param path modules are added to this path
2041 """
2042
2043 prlist = register_module('ParticlePrinter')
2044 prlist.set_name('ParticlePrinter_' + list_name)
2045 prlist.param('listName', list_name)
2046 prlist.param('fullPrint', False)
2047 prlist.param('variables', var_names)
2048 path.add_module(prlist)
2049
2050
2051def printList(list_name, full, path):
2052 """
2053 Prints the size and executes Particle->print() (if full=True)
2054 method for all Particles in given ParticleList. For debugging purposes.
2055
2056 @param list_name input ParticleList name
2057 @param full execute Particle->print() method for all Particles
2058 @param path modules are added to this path
2059 """
2060
2061 prlist = register_module('ParticlePrinter')
2062 prlist.set_name('ParticlePrinter_' + list_name)
2063 prlist.param('listName', list_name)
2064 prlist.param('fullPrint', full)
2065 path.add_module(prlist)
2066
2067
2068def variablesToNtuple(decayString, variables, treename='variables', filename='ntuple.root', path=None, basketsize=1600,
2069 signalSideParticleList="", filenameSuffix="", useFloat=False, storeEventType=True,
2070 ignoreCommandLineOverride=False):
2071 """
2072 Creates and fills a flat ntuple with the specified variables from the VariableManager.
2073 If a decayString is provided, then there will be one entry per candidate (for particle in list of candidates).
2074 If an empty decayString is provided, there will be one entry per event (useful for trigger studies, etc).
2075
2076 Parameters:
2077 decayString (str): specifies type of Particles and determines the name of the ParticleList
2078 variables (list(str)): the list of variables (which must be registered in the VariableManager)
2079 treename (str): name of the ntuple tree
2080 filename (str): which is used to store the variables
2081 path (basf2.Path): the basf2 path where the analysis is processed
2082 basketsize (int): size of baskets in the output ntuple in bytes
2083 signalSideParticleList (str): The name of the signal-side ParticleList.
2084 Only valid if the module is called in a for_each loop over the RestOfEvent.
2085 filenameSuffix (str): suffix to be appended to the filename before ``.root``.
2086 useFloat (bool): Use single precision (float) instead of double precision (double)
2087 for floating-point numbers.
2088 storeEventType (bool) : if true, the branch __eventType__ is added for the MC event type information.
2089 The information is available from MC16 on.
2090 ignoreCommandLineOverride (bool) : if true, ignore override of file name via command line argument ``-o``.
2091
2092 .. tip:: The output filename can be overridden using the ``-o`` argument of basf2.
2093 """
2094
2095 output = register_module('VariablesToNtuple')
2096 output.set_name('VariablesToNtuple_' + decayString)
2097 output.param('particleList', decayString)
2098 output.param('variables', variables)
2099 output.param('fileName', filename)
2100 output.param('treeName', treename)
2101 output.param('basketSize', basketsize)
2102 output.param('signalSideParticleList', signalSideParticleList)
2103 output.param('fileNameSuffix', filenameSuffix)
2104 output.param('useFloat', useFloat)
2105 output.param('storeEventType', storeEventType)
2106 output.param('ignoreCommandLineOverride', ignoreCommandLineOverride)
2107 path.add_module(output)
2108
2109
2110def variablesToHistogram(decayString,
2111 variables,
2112 variables_2d=None,
2113 filename='ntuple.root',
2114 path=None, *,
2115 directory=None,
2116 prefixDecayString=False,
2117 filenameSuffix="",
2118 ignoreCommandLineOverride=False):
2119 """
2120 Creates and fills a flat ntuple with the specified variables from the VariableManager
2121
2122 Parameters:
2123 decayString (str): specifies type of Particles and determines the name of the ParticleList
2124 variables (list(tuple))): variables + binning which must be registered in the VariableManager
2125 variables_2d (list(tuple)): pair of variables + binning for each which must be registered in the VariableManager
2126 filename (str): which is used to store the variables
2127 path (basf2.Path): the basf2 path where the analysis is processed
2128 directory (str): directory inside the output file where the histograms should be saved.
2129 Useful if you want to have different histograms in the same file to separate them.
2130 prefixDecayString (bool): If True the decayString will be prepended to the directory name to allow for more
2131 programmatic naming of the structure in the file.
2132 filenameSuffix (str): suffix to be appended to the filename before ``.root``.
2133 ignoreCommandLineOverride (bool) : if true, ignore override of file name via command line argument ``-o``.
2134
2135 .. tip:: The output filename can be overridden using the ``-o`` argument of basf2.
2136 """
2137
2138 if variables_2d is None:
2139 variables_2d = []
2140 output = register_module('VariablesToHistogram')
2141 output.set_name('VariablesToHistogram_' + decayString)
2142 output.param('particleList', decayString)
2143 output.param('variables', variables)
2144 output.param('variables_2d', variables_2d)
2145 output.param('fileName', filename)
2146 output.param('fileNameSuffix', filenameSuffix)
2147 output.param('ignoreCommandLineOverride', ignoreCommandLineOverride)
2148 if directory is not None or prefixDecayString:
2149 if directory is None:
2150 directory = ""
2151 if prefixDecayString:
2152 directory = decayString + "_" + directory
2153 output.param("directory", directory)
2154 path.add_module(output)
2155
2156
2157def variablesToExtraInfo(particleList, variables, option=0, path=None):
2158 """
2159 For each particle in the input list the selected variables are saved in an extra-info field with the given name.
2160 Can be used when wanting to save variables before modifying them, e.g. when performing vertex fits.
2161
2162 Parameters:
2163 particleList (str): The input ParticleList
2164 variables (dict[str,str]): Dictionary of Variables (key) and extraInfo names (value).
2165 option (int): Option to overwrite an existing extraInfo. Choose among -1, 0, 1, 2.
2166 An existing extra info with the same name will be overwritten if the new
2167 value is lower / will never be overwritten / will be overwritten if the
2168 new value is higher / will always be overwritten (option = -1/0/1/2).
2169 path (basf2.Path): modules are added to this path
2170 """
2171
2172 mod = register_module('VariablesToExtraInfo')
2173 mod.set_name('VariablesToExtraInfo_' + particleList)
2174 mod.param('particleList', particleList)
2175 mod.param('variables', variables)
2176 mod.param('overwrite', option)
2177 path.add_module(mod)
2178
2179
2180def variablesToDaughterExtraInfo(particleList, decayString, variables, option=0, path=None):
2181 """
2182 For each daughter particle specified via decay string the selected variables (estimated for the mother particle)
2183 are saved in an extra-info field with the given name. In other words, the property of mother is saved as extra-info
2184 to specified daughter particle.
2185
2186 Parameters:
2187 particleList (str): The input ParticleList
2188 decayString (str): Decay string that specifies to which daughter the extra info should be appended
2189 variables (dict[str,str]): Dictionary of Variables (key) and extraInfo names (value).
2190 option (int): Option to overwrite an existing extraInfo. Choose among -1, 0, 1, 2.
2191 An existing extra info with the same name will be overwritten if the new
2192 value is lower / will never be overwritten / will be overwritten if the
2193 new value is higher / will always be overwritten (option = -1/0/1/2).
2194 path (basf2.Path): modules are added to this path
2195 """
2196
2197 mod = register_module('VariablesToExtraInfo')
2198 mod.set_name('VariablesToDaughterExtraInfo_' + particleList)
2199 mod.param('particleList', particleList)
2200 mod.param('decayString', decayString)
2201 mod.param('variables', variables)
2202 mod.param('overwrite', option)
2203 path.add_module(mod)
2204
2205
2206def variablesToEventExtraInfo(particleList, variables, option=0, path=None):
2207 """
2208 For each particle in the input list the selected variables are saved in an event-extra-info field with the given name,
2209 Can be used to save MC truth information, for example, in a ntuple of reconstructed particles.
2210
2211 .. tip::
2212 When the function is called first time not in the main path but in a sub-path e.g. ``roe_path``,
2213 the eventExtraInfo cannot be accessed from the main path because of the shorter lifetime of the event-extra-info field.
2214 If one wants to call the function in a sub-path, one has to call the function in the main path beforehand.
2215
2216 Parameters:
2217 particleList (str): The input ParticleList
2218 variables (dict[str,str]): Dictionary of Variables (key) and extraInfo names (value).
2219 option (int): Option to overwrite an existing extraInfo. Choose among -1, 0, 1, 2.
2220 An existing extra info with the same name will be overwritten if the new
2221 value is lower / will never be overwritten / will be overwritten if the
2222 new value is higher / will always be overwritten (option = -1/0/1/2).
2223 path (basf2.Path): modules are added to this path
2224 """
2225
2226 mod = register_module('VariablesToEventExtraInfo')
2227 mod.set_name('VariablesToEventExtraInfo_' + particleList)
2228 mod.param('particleList', particleList)
2229 mod.param('variables', variables)
2230 mod.param('overwrite', option)
2231 path.add_module(mod)
2232
2233
2234def variableToSignalSideExtraInfo(particleList, varToExtraInfo, path):
2235 """
2236 Write the value of specified variables estimated for the single particle in the input list (has to contain exactly 1
2237 particle) as an extra info to the particle related to current ROE.
2238 Should be used only in the for_each roe path.
2239
2240 Parameters:
2241 particleList (str): The input ParticleList
2242 varToExtraInfo (dict[str,str]): Dictionary of Variables (key) and extraInfo names (value).
2243 path (basf2.Path): modules are added to this path
2244 """
2245
2246 mod = register_module('SignalSideVariablesToExtraInfo')
2247 mod.set_name('SigSideVarToExtraInfo_' + particleList)
2248 mod.param('particleListName', particleList)
2249 mod.param('variableToExtraInfo', varToExtraInfo)
2250 path.add_module(mod)
2251
2252
2253def signalRegion(particleList, cut, path=None, name="isSignalRegion", blind_data=True):
2254 """
2255 Define and blind a signal region.
2256 Per default, the defined signal region is cut out if ran on data.
2257 This function will provide a new variable 'isSignalRegion' as default, which is either 0 or 1 depending on the cut
2258 provided.
2259
2260 Example:
2261 .. code-block:: python
2262
2263 ma.reconstructDecay("B+:sig -> D+ pi0", "Mbc>5.2", path=path)
2264 ma.signalRegion("B+:sig",
2265 "Mbc>5.27 and abs(deltaE)<0.2",
2266 blind_data=True,
2267 path=path)
2268 ma.variablesToNtuples("B+:sig", ["isSignalRegion"], path=path)
2269
2270 Parameters:
2271 particleList (str): The input ParticleList
2272 cut (str): Cut string describing the signal region
2273 path (basf2.Path):: Modules are added to this path
2274 name (str): Name of the Signal region in the variable manager
2275 blind_data (bool): Automatically exclude signal region from data
2276
2277 """
2278
2279 from variables import variables
2280 mod = register_module('VariablesToExtraInfo')
2281 mod.set_name(f'{name}_' + particleList)
2282 mod.param('particleList', particleList)
2283 mod.param('variables', {f"passesCut({cut})": name})
2284 variables.addAlias(name, f"extraInfo({name})")
2285 path.add_module(mod)
2286
2287 # Check if we run on Data
2288 if blind_data:
2289 applyCuts(particleList, f"{name}==0 or isMC==1", path=path)
2290
2291
2292def removeExtraInfo(particleLists=None, removeEventExtraInfo=False, path=None):
2293 """
2294 Removes the ExtraInfo of the given particleLists. If specified (removeEventExtraInfo = True) also the EventExtraInfo is removed.
2295 """
2296
2297 if particleLists is None:
2298 particleLists = []
2299 mod = register_module('ExtraInfoRemover')
2300 mod.param('particleLists', particleLists)
2301 mod.param('removeEventExtraInfo', removeEventExtraInfo)
2302 path.add_module(mod)
2303
2304
2305def signalSideParticleFilter(particleList, selection, roe_path, deadEndPath):
2306 """
2307 Checks if the current ROE object in the for_each roe path (argument roe_path) is related
2308 to the particle from the input ParticleList. Additional selection criteria can be applied.
2309 If ROE is not related to any of the Particles from ParticleList or the Particle doesn't
2310 meet the selection criteria the execution of deadEndPath is started. This path, as the name
2311 suggests should be empty and its purpose is to end the execution of for_each roe path for
2312 the current ROE object.
2313
2314 @param particleList The input ParticleList
2315 @param selection Selection criteria that Particle needs meet in order for for_each ROE path to continue
2316 @param for_each roe path in which this filter is executed
2317 @param deadEndPath empty path that ends execution of or_each roe path for the current ROE object.
2318 """
2319
2320 mod = register_module('SignalSideParticleFilter')
2321 mod.set_name('SigSideParticleFilter_' + particleList)
2322 mod.param('particleLists', [particleList])
2323 mod.param('selection', selection)
2324 roe_path.add_module(mod)
2325 mod.if_false(deadEndPath)
2326
2327
2328def signalSideParticleListsFilter(particleLists, selection, roe_path, deadEndPath):
2329 """
2330 Checks if the current ROE object in the for_each roe path (argument roe_path) is related
2331 to the particle from the input ParticleList. Additional selection criteria can be applied.
2332 If ROE is not related to any of the Particles from ParticleList or the Particle doesn't
2333 meet the selection criteria the execution of deadEndPath is started. This path, as the name
2334 suggests should be empty and its purpose is to end the execution of for_each roe path for
2335 the current ROE object.
2336
2337 @param particleLists The input ParticleLists
2338 @param selection Selection criteria that Particle needs meet in order for for_each ROE path to continue
2339 @param for_each roe path in which this filter is executed
2340 @param deadEndPath empty path that ends execution of or_each roe path for the current ROE object.
2341 """
2342
2343 mod = register_module('SignalSideParticleFilter')
2344 mod.set_name('SigSideParticleFilter_' + particleLists[0])
2345 mod.param('particleLists', particleLists)
2346 mod.param('selection', selection)
2347 roe_path.add_module(mod)
2348 mod.if_false(deadEndPath)
2349
2350
2352 decayString,
2353 cut,
2354 dmID=0,
2355 writeOut=False,
2356 path=None,
2357 chargeConjugation=True,
2358):
2359 r"""
2360 Finds and creates a ``ParticleList`` from given decay string.
2361 ``ParticleList`` of daughters with sub-decay is created.
2362
2363 Only the particles made from MCParticle, which can be loaded by `fillParticleListFromMC`, are accepted as daughters.
2364
2365 Only signal particle, which means :b2:var:`isSignal` is equal to 1, is stored. One can use the decay string grammar
2366 to change the behavior of :b2:var:`isSignal`. One can find detailed information in :ref:`DecayString`.
2367
2368 .. tip::
2369 If one uses same sub-decay twice, same particles are registered to a ``ParticleList``. For example,
2370 ``K_S0:pi0pi0 =direct=> [pi0:gg =direct=> gamma:MC gamma:MC] [pi0:gg =direct=> gamma:MC gamma:MC]``.
2371 One can skip the second sub-decay, ``K_S0:pi0pi0 =direct=> [pi0:gg =direct=> gamma:MC gamma:MC] pi0:gg``.
2372
2373 .. tip::
2374 It is recommended to use only primary particles as daughter particles unless you want to explicitly study the secondary
2375 particles. The behavior of MC-matching for secondary particles from a stable particle decay is not guaranteed.
2376 Please consider to use `fillParticleListFromMC` with ``skipNonPrimary=True`` to load daughter particles.
2377 Moreover, it is recommended to load ``K_S0`` and ``Lambda0`` directly from MCParticle by `fillParticleListFromMC` rather
2378 than reconstructing from two pions or a proton-pion pair, because their direct daughters can be the secondary particle.
2379
2380
2381 @param decayString :ref:`DecayString` specifying what kind of the decay should be reconstructed
2382 (from the DecayString the mother and daughter ParticleLists are determined)
2383 @param cut created (mother) Particles are added to the mother ParticleList if they
2384 pass given cuts (in VariableManager style) and rejected otherwise
2385 isSignal==1 is always required by default.
2386 @param dmID user specified decay mode identifier
2387 @param writeOut whether RootOutput module should save the created ParticleList
2388 @param path modules are added to this path
2389 @param chargeConjugation boolean to decide whether charge conjugated mode should be reconstructed as well (on by default)
2390 """
2391
2392 pmake = register_module('ParticleCombinerFromMC')
2393 pmake.set_name('ParticleCombinerFromMC_' + decayString)
2394 pmake.param('decayString', decayString)
2395 pmake.param('cut', cut)
2396 pmake.param('decayMode', dmID)
2397 pmake.param('writeOut', writeOut)
2398 pmake.param('chargeConjugation', chargeConjugation)
2399 path.add_module(pmake)
2400
2401
2402def findMCDecay(
2403 list_name,
2404 decay,
2405 writeOut=False,
2406 appendAllDaughters=False,
2407 skipNonPrimaryDaughters=True,
2408 path=None,
2409):
2410 """
2411 Finds and creates a ``ParticleList`` for all ``MCParticle`` decays matching a given :ref:`DecayString`.
2412 The decay string is required to describe correctly what you want.
2413 In the case of inclusive decays, you can use :ref:`Grammar_for_custom_MCMatching`
2414
2415 The output particles has only the daughter particles written in the given decay string, if
2416 ``appendAllDaughters=False`` (default). If ``appendAllDaughters=True``, all daughters of the matched MCParticle are
2417 appended in the order defined at the MCParticle level. For example,
2418
2419 .. code-block:: python
2420
2421 findMCDecay('B0:Xee', 'B0 -> e+ e- ... ?gamma', appendAllDaughters=False, path=mypath)
2422
2423 The output ParticleList ``B0:Xee`` will match the inclusive ``B0 -> e+ e-`` decays (but neutrinos are not included),
2424 in both cases of ``appendAllDaughters`` is false and true.
2425 If the ``appendAllDaughters=False`` as above example, the ``B0:Xee`` has only two electrons as daughters.
2426 While, if ``appendAllDaughters=True``, all daughters of the matched MCParticles are appended. When the truth decay mode of
2427 the MCParticle is ``B0 -> [K*0 -> K+ pi-] [J/psi -> e+ e-]``, the first daughter of ``B0:Xee`` is ``K*0`` and ``e+``
2428 will be the first daughter of second daughter of ``B0:Xee``.
2429
2430 The option ``skipNonPrimaryDaughters`` only has an effect if ``appendAllDaughters=True``. If ``skipNonPrimaryDaughters=True``,
2431 all primary daughters are appended but the secondary particles are not.
2432
2433 .. tip::
2434 Daughters of ``Lambda0`` are not primary, but ``Lambda0`` is not a final state particle.
2435 In order for the MCMatching to work properly, the daughters of ``Lambda0`` are appended to
2436 ``Lambda0`` regardless of the value of the option ``skipNonPrimaryDaughters``.
2437
2438
2439 @param list_name The output particle list name
2440 @param decay The decay string which you want
2441 @param writeOut Whether `RootOutput` module should save the created ``outputList``
2442 @param skipNonPrimaryDaughters if true, skip non primary daughters, useful to study final state daughter particles
2443 @param appendAllDaughters if true, not only the daughters described in the decay string but all daughters are appended
2444 @param path modules are added to this path
2445 """
2446
2447 decayfinder = register_module('MCDecayFinder')
2448 decayfinder.set_name('MCDecayFinder_' + list_name)
2449 decayfinder.param('listName', list_name)
2450 decayfinder.param('decayString', decay)
2451 decayfinder.param('appendAllDaughters', appendAllDaughters)
2452 decayfinder.param('skipNonPrimaryDaughters', skipNonPrimaryDaughters)
2453 decayfinder.param('writeOut', writeOut)
2454 path.add_module(decayfinder)
2455
2456
2457def summaryOfLists(particleLists, outputFile=None, path=None):
2458 """
2459 Prints out Particle statistics at the end of the job: number of events with at
2460 least one candidate, average number of candidates per event, etc.
2461 If an output file name is provided the statistics is also dumped into a json file with that name.
2462
2463 @param particleLists list of input ParticleLists
2464 @param outputFile output file name (not created by default)
2465 """
2466
2467 particleStats = register_module('ParticleStats')
2468 particleStats.param('particleLists', particleLists)
2469 if outputFile is not None:
2470 particleStats.param('outputFile', outputFile)
2471 path.add_module(particleStats)
2472
2473
2474def matchMCTruth(list_name, path):
2475 """
2476 Performs MC matching (sets relation Particle->MCParticle) for
2477 all particles (and its (grand)^N-daughter particles) in the specified
2478 ParticleList.
2479
2480 @param list_name name of the input ParticleList
2481 @param path modules are added to this path
2482 """
2483
2484 mcMatch = register_module('MCMatcherParticles')
2485 mcMatch.set_name('MCMatch_' + list_name)
2486 mcMatch.param('listName', list_name)
2487 path.add_module(mcMatch)
2488
2489
2490def looseMCTruth(list_name, path):
2491 """
2492 Performs loose MC matching for all particles in the specified
2493 ParticleList.
2494 The difference between loose and normal mc matching algorithm is that
2495 the loose algorithm will find the common mother of the majority of daughter
2496 particles while the normal algorithm finds the common mother of all daughters.
2497 The results of loose mc matching algorithm are stored to the following extraInfo
2498 items:
2499
2500 - looseMCMotherPDG: PDG code of most common mother
2501 - looseMCMotherIndex: 1-based StoreArray<MCParticle> index of most common mother
2502 - looseMCWrongDaughterN: number of daughters that don't originate from the most common mother
2503 - looseMCWrongDaughterPDG: PDG code of the daughter that doesn't originate from the most common mother (only if
2504 looseMCWrongDaughterN = 1)
2505 - looseMCWrongDaughterBiB: 1 if the wrong daughter is Beam Induced Background Particle
2506
2507 @param list_name name of the input ParticleList
2508 @param path modules are added to this path
2509 """
2510
2511 mcMatch = register_module('MCMatcherParticles')
2512 mcMatch.set_name('LooseMCMatch_' + list_name)
2513 mcMatch.param('listName', list_name)
2514 mcMatch.param('looseMCMatching', True)
2515 path.add_module(mcMatch)
2516
2517
2518def buildRestOfEvent(target_list_name, inputParticlelists=None,
2519 fillWithMostLikely=True,
2520 chargedPIDPriors=None, path=None):
2521 """
2522 Creates for each Particle in the given ParticleList a RestOfEvent
2523 dataobject and makes basf2 relation between them. User can provide additional
2524 particle lists with a different particle hypothesis like ['K+:good, e+:good'], etc.
2525
2526 @param target_list_name name of the input ParticleList
2527 @param inputParticlelists list of user-defined input particle list names, which serve
2528 as source of particles to build the ROE, the FSP particles from
2529 target_list_name are automatically excluded from the ROE object
2530 @param fillWithMostLikely By default the module uses the most likely particle mass hypothesis for charged particles
2531 based on the PID likelihood. Turn this behavior off if you want to configure your own
2532 input particle lists.
2533 @param chargedPIDPriors The prior PID fractions, that are used to regulate the
2534 amount of certain charged particle species, should be a list of
2535 six floats if not None. The order of particle types is
2536 the following: [e-, mu-, pi-, K-, p+, d+]
2537 @param path modules are added to this path
2538 """
2539
2540 if inputParticlelists is None:
2541 inputParticlelists = []
2542 fillParticleList('pi+:all', '', path=path)
2543 if fillWithMostLikely:
2544 from stdCharged import stdMostLikely
2545 stdMostLikely(chargedPIDPriors, '_roe', path=path)
2546 inputParticlelists = [f'{ptype}:mostlikely_roe' for ptype in ['K+', 'p+', 'e+', 'mu+']]
2547 import b2bii
2548 if not b2bii.isB2BII():
2549 fillParticleList('gamma:all', '', path=path)
2550 fillParticleList('K_L0:roe_default', 'isFromKLM > 0', path=path)
2551 inputParticlelists += ['pi+:all', 'gamma:all', 'K_L0:roe_default']
2552 else:
2553 inputParticlelists += ['pi+:all', 'gamma:mdst']
2554 roeBuilder = register_module('RestOfEventBuilder')
2555 roeBuilder.set_name('ROEBuilder_' + target_list_name)
2556 roeBuilder.param('particleList', target_list_name)
2557 roeBuilder.param('particleListsInput', inputParticlelists)
2558 roeBuilder.param('mostLikely', fillWithMostLikely)
2559 path.add_module(roeBuilder)
2560
2561
2562def buildNestedRestOfEvent(target_list_name, maskName='all', path=None):
2563 """
2564 Creates for each Particle in the given ParticleList a RestOfEvent
2565 @param target_list_name name of the input ParticleList
2566 @param mask_name name of the ROEMask to be used
2567 @param path modules are added to this path
2568 """
2569
2570 roeBuilder = register_module('RestOfEventBuilder')
2571 roeBuilder.set_name('NestedROEBuilder_' + target_list_name)
2572 roeBuilder.param('particleList', target_list_name)
2573 roeBuilder.param('nestedROEMask', maskName)
2574 roeBuilder.param('createNestedROE', True)
2575 path.add_module(roeBuilder)
2576
2577
2578def buildRestOfEventFromMC(target_list_name, inputParticlelists=None, path=None):
2579 """
2580 Creates for each Particle in the given ParticleList a RestOfEvent
2581 @param target_list_name name of the input ParticleList
2582 @param inputParticlelists list of input particle list names, which serve
2583 as a source of particles to build ROE, the FSP particles from
2584 target_list_name are excluded from ROE object
2585 @param path modules are added to this path
2586 """
2587
2588 if inputParticlelists is None:
2589 inputParticlelists = []
2590 if (len(inputParticlelists) == 0):
2591 # Type of particles to use for ROEBuilder
2592 # K_S0 and Lambda0 are added here because some of them have interacted
2593 # with the detector material
2594 types = ['gamma', 'e+', 'mu+', 'pi+', 'K+', 'p+', 'K_L0',
2595 'n0', 'nu_e', 'nu_mu', 'nu_tau',
2596 'K_S0', 'Lambda0']
2597 for t in types:
2598 fillParticleListFromMC(f"{t}:roe_default_gen", 'mcPrimary > 0 and nDaughters == 0',
2599 True, True, path=path)
2600 inputParticlelists += [f"{t}:roe_default_gen"]
2601 roeBuilder = register_module('RestOfEventBuilder')
2602 roeBuilder.set_name('MCROEBuilder_' + target_list_name)
2603 roeBuilder.param('particleList', target_list_name)
2604 roeBuilder.param('particleListsInput', inputParticlelists)
2605 roeBuilder.param('fromMC', True)
2606 path.add_module(roeBuilder)
2607
2608
2609def appendROEMask(list_name,
2610 mask_name,
2611 trackSelection,
2612 eclClusterSelection,
2613 klmClusterSelection='',
2614 path=None):
2615 """
2616 Loads the ROE object of a particle and creates a ROE mask with a specific name. It applies
2617 selection criteria for tracks and eclClusters which will be used by variables in ROEVariables.cc.
2618
2619 - append a ROE mask with all tracks in ROE coming from the IP region
2620
2621 .. code-block:: python
2622
2623 appendROEMask('B+:sig', 'IPtracks', '[dr < 2] and [abs(dz) < 5]', path=mypath)
2624
2625 - append a ROE mask with only ECL-based particles that pass as good photon candidates
2626
2627 .. code-block:: python
2628
2629 goodPhotons = 'inCDCAcceptance and clusterErrorTiming < 1e6 and [clusterE1E9 > 0.4 or E > 0.075]'
2630 appendROEMask('B+:sig', 'goodROEGamma', '', goodPhotons, path=mypath)
2631
2632
2633 @param list_name name of the input ParticleList
2634 @param mask_name name of the appended ROEMask
2635 @param trackSelection decay string for the track-based particles in ROE
2636 @param eclClusterSelection decay string for the ECL-based particles in ROE
2637 @param klmClusterSelection decay string for the KLM-based particles in ROE
2638 @param path modules are added to this path
2639 """
2640
2641 roeMask = register_module('RestOfEventInterpreter')
2642 roeMask.set_name('RestOfEventInterpreter_' + list_name + '_' + mask_name)
2643 roeMask.param('particleList', list_name)
2644 roeMask.param('ROEMasks', [(mask_name, trackSelection, eclClusterSelection, klmClusterSelection)])
2645 path.add_module(roeMask)
2646
2647
2648def appendROEMasks(list_name, mask_tuples, path=None):
2649 """
2650 Loads the ROE object of a particle and creates a ROE mask with a specific name. It applies
2651 selection criteria for track-, ECL- and KLM-based particles which will be used by ROE variables.
2652
2653 The multiple ROE masks with their own selection criteria are specified
2654 via list of tuples (mask_name, trackParticleSelection, eclParticleSelection, klmParticleSelection) or
2655 (mask_name, trackSelection, eclClusterSelection) in case with fractions.
2656
2657 - Example for two tuples, one with and one without fractions
2658
2659 .. code-block:: python
2660
2661 ipTracks = ('IPtracks', '[dr < 2] and [abs(dz) < 5]', '', '')
2662 goodPhotons = 'inCDCAcceptance and [clusterErrorTiming < 1e6] and [clusterE1E9 > 0.4 or E > 0.075]'
2663 goodROEGamma = ('ROESel', '[dr < 2] and [abs(dz) < 5]', goodPhotons, '')
2664 goodROEKLM = ('IPtracks', '[dr < 2] and [abs(dz) < 5]', '', 'nKLMClusterTrackMatches == 0')
2665 appendROEMasks('B+:sig', [ipTracks, goodROEGamma, goodROEKLM], path=mypath)
2666
2667 @param list_name name of the input ParticleList
2668 @param mask_tuples array of ROEMask list tuples to be appended
2669 @param path modules are added to this path
2670 """
2671
2672 compatible_masks = []
2673 for mask in mask_tuples:
2674 # add empty KLM-based selection if it's absent:
2675 if len(mask) == 3:
2676 compatible_masks += [(*mask, '')]
2677 else:
2678 compatible_masks += [mask]
2679 roeMask = register_module('RestOfEventInterpreter')
2680 roeMask.set_name('RestOfEventInterpreter_' + list_name + '_' + 'MaskList')
2681 roeMask.param('particleList', list_name)
2682 roeMask.param('ROEMasks', compatible_masks)
2683 path.add_module(roeMask)
2684
2685
2686def updateROEMask(list_name,
2687 mask_name,
2688 trackSelection,
2689 eclClusterSelection='',
2690 klmClusterSelection='',
2691 path=None):
2692 """
2693 Update an existing ROE mask by applying additional selection cuts for
2694 tracks and/or clusters.
2695
2696 See function `appendROEMask`!
2697
2698 @param list_name name of the input ParticleList
2699 @param mask_name name of the ROEMask to update
2700 @param trackSelection decay string for the track-based particles in ROE
2701 @param eclClusterSelection decay string for the ECL-based particles in ROE
2702 @param klmClusterSelection decay string for the KLM-based particles in ROE
2703 @param path modules are added to this path
2704 """
2705
2706 roeMask = register_module('RestOfEventInterpreter')
2707 roeMask.set_name('RestOfEventInterpreter_' + list_name + '_' + mask_name)
2708 roeMask.param('particleList', list_name)
2709 roeMask.param('ROEMasks', [(mask_name, trackSelection, eclClusterSelection, klmClusterSelection)])
2710 roeMask.param('update', True)
2711 path.add_module(roeMask)
2712
2713
2714def updateROEMasks(list_name, mask_tuples, path):
2715 """
2716 Update existing ROE masks by applying additional selection cuts for tracks
2717 and/or clusters.
2718
2719 The multiple ROE masks with their own selection criteria are specified
2720 via list tuples (mask_name, trackSelection, eclClusterSelection, klmClusterSelection)
2721
2722 See function `appendROEMasks`!
2723
2724 @param list_name name of the input ParticleList
2725 @param mask_tuples array of ROEMask list tuples to be appended
2726 @param path modules are added to this path
2727 """
2728
2729 compatible_masks = []
2730 for mask in mask_tuples:
2731 # add empty KLM-based selection if it's absent:
2732 if len(mask) == 3:
2733 compatible_masks += [(*mask, '')]
2734 else:
2735 compatible_masks += [mask]
2736
2737 roeMask = register_module('RestOfEventInterpreter')
2738 roeMask.set_name('RestOfEventInterpreter_' + list_name + '_' + 'MaskList')
2739 roeMask.param('particleList', list_name)
2740 roeMask.param('ROEMasks', compatible_masks)
2741 roeMask.param('update', True)
2742 path.add_module(roeMask)
2743
2744
2745def keepInROEMasks(list_name, mask_names, cut_string, path=None):
2746 """
2747 This function is used to apply particle list specific cuts on one or more ROE masks (track or eclCluster).
2748 With this function one can KEEP the tracks/eclclusters used in particles from provided particle list.
2749 This function should be executed only in the for_each roe path for the current ROE object.
2750
2751 To avoid unnecessary computation, the input particle list should only contain particles from ROE
2752 (use cut 'isInRestOfEvent == 1'). To update the ECLCluster masks, the input particle list should be a photon
2753 particle list (e.g. 'gamma:someLabel'). To update the Track masks, the input particle list should be a charged
2754 pion particle list (e.g. 'pi+:someLabel').
2755
2756 Updating a non-existing mask will create a new one.
2757
2758 - keep only those tracks that were used in provided particle list
2759
2760 .. code-block:: python
2761
2762 keepInROEMasks('pi+:goodTracks', 'mask', '', path=mypath)
2763
2764 - keep only those clusters that were used in provided particle list and pass a cut, apply to several masks
2765
2766 .. code-block:: python
2767
2768 keepInROEMasks('gamma:goodClusters', ['mask1', 'mask2'], 'E > 0.1', path=mypath)
2769
2770
2771 @param list_name name of the input ParticleList
2772 @param mask_names array of ROEMasks to be updated
2773 @param cut_string decay string with which the mask will be updated
2774 @param path modules are added to this path
2775 """
2776
2777 updateMask = register_module('RestOfEventUpdater')
2778 updateMask.set_name('RestOfEventUpdater_' + list_name + '_masks')
2779 updateMask.param('particleList', list_name)
2780 updateMask.param('updateMasks', mask_names)
2781 updateMask.param('cutString', cut_string)
2782 updateMask.param('discard', False)
2783 path.add_module(updateMask)
2784
2785
2786def discardFromROEMasks(list_name, mask_names, cut_string, path=None):
2787 """
2788 This function is used to apply particle list specific cuts on one or more ROE masks (track or eclCluster).
2789 With this function one can DISCARD the tracks/eclclusters used in particles from provided particle list.
2790 This function should be executed only in the for_each roe path for the current ROE object.
2791
2792 To avoid unnecessary computation, the input particle list should only contain particles from ROE
2793 (use cut 'isInRestOfEvent == 1'). To update the ECLCluster masks, the input particle list should be a photon
2794 particle list (e.g. 'gamma:someLabel'). To update the Track masks, the input particle list should be a charged
2795 pion particle list (e.g. 'pi+:someLabel').
2796
2797 Updating a non-existing mask will create a new one.
2798
2799 - discard tracks that were used in provided particle list
2800
2801 .. code-block:: python
2802
2803 discardFromROEMasks('pi+:badTracks', 'mask', '', path=mypath)
2804
2805 - discard clusters that were used in provided particle list and pass a cut, apply to several masks
2806
2807 .. code-block:: python
2808
2809 discardFromROEMasks('gamma:badClusters', ['mask1', 'mask2'], 'E < 0.1', path=mypath)
2810
2811
2812 @param list_name name of the input ParticleList
2813 @param mask_names array of ROEMasks to be updated
2814 @param cut_string decay string with which the mask will be updated
2815 @param path modules are added to this path
2816 """
2817
2818 updateMask = register_module('RestOfEventUpdater')
2819 updateMask.set_name('RestOfEventUpdater_' + list_name + '_masks')
2820 updateMask.param('particleList', list_name)
2821 updateMask.param('updateMasks', mask_names)
2822 updateMask.param('cutString', cut_string)
2823 updateMask.param('discard', True)
2824 path.add_module(updateMask)
2825
2826
2827def optimizeROEWithV0(list_name, mask_names, cut_string, path=None):
2828 """
2829 This function is used to apply particle list specific cuts on one or more ROE masks for Tracks.
2830 It is possible to optimize the ROE selection by treating tracks from V0's separately, meaning,
2831 taking V0's 4-momentum into account instead of 4-momenta of tracks. A cut for only specific V0's
2832 passing it can be applied.
2833
2834 The input particle list should be a V0 particle list: K_S0 ('K_S0:someLabel', ''),
2835 Lambda ('Lambda:someLabel', '') or converted photons ('gamma:someLabel').
2836
2837 Updating a non-existing mask will create a new one.
2838
2839 - treat tracks from K_S0 inside mass window separately, replace track momenta with K_S0 momentum
2840
2841 .. code-block:: python
2842
2843 optimizeROEWithV0('K_S0:opt', 'mask', '0.450 < M < 0.550', path=mypath)
2844
2845 @param list_name name of the input ParticleList
2846 @param mask_names array of ROEMasks to be updated
2847 @param cut_string decay string with which the mask will be updated
2848 @param path modules are added to this path
2849 """
2850
2851 updateMask = register_module('RestOfEventUpdater')
2852 updateMask.set_name('RestOfEventUpdater_' + list_name + '_masks')
2853 updateMask.param('particleList', list_name)
2854 updateMask.param('updateMasks', mask_names)
2855 updateMask.param('cutString', cut_string)
2856 path.add_module(updateMask)
2857
2858
2859def updateROEUsingV0Lists(target_particle_list, mask_names, default_cleanup=True, selection_cuts=None,
2860 apply_mass_fit=False, fitter='treefit', path=None):
2861 """
2862 This function creates V0 particle lists (photons, :math:`K^0_S` and :math:`\\Lambda^0`)
2863 and it uses V0 candidates to update the Rest Of Event, which is associated to the target particle list.
2864 It is possible to apply a standard or customized selection and mass fit to the V0 candidates.
2865
2866
2867 @param target_particle_list name of the input ParticleList
2868 @param mask_names array of ROE masks to be applied
2869 @param default_cleanup if True, predefined cuts will be applied on the V0 lists
2870 @param selection_cuts a single string of selection cuts or tuple of three strings (photon_cuts, K_S0_cuts, Lambda0_cuts),
2871 which will be applied to the V0 lists. These cuts will have a priority over the default ones.
2872 @param apply_mass_fit if True, a mass fit will be applied to the V0 particles
2873 @param fitter string, that represent a fitter choice: "treefit" for TreeFitter and "kfit" for KFit
2874 @param path modules are added to this path
2875 """
2876
2877 roe_path = create_path()
2878 deadEndPath = create_path()
2879 signalSideParticleFilter(target_particle_list, '', roe_path, deadEndPath)
2880
2881 if (default_cleanup and selection_cuts is None):
2882 B2INFO("Using default cleanup in updateROEUsingV0Lists.")
2883 selection_cuts = 'abs(dM) < 0.1 '
2884 selection_cuts += 'and daughter(0,particleID) > 0.2 and daughter(1,particleID) > 0.2 '
2885 selection_cuts += 'and daughter(0,thetaInCDCAcceptance) and daughter(1,thetaInCDCAcceptance)'
2886 if (selection_cuts is None or selection_cuts == ''):
2887 B2INFO("No cleanup in updateROEUsingV0Lists.")
2888 selection_cuts = ('True', 'True', 'True')
2889 if (isinstance(selection_cuts, str)):
2890 selection_cuts = (selection_cuts, selection_cuts, selection_cuts)
2891 # The isInRestOfEvent variable will be applied on FSPs of composite particles automatically:
2892 roe_cuts = 'isInRestOfEvent > 0'
2893 fillConvertedPhotonsList('gamma:v0_roe -> e+ e-', f'{selection_cuts[0]} and {roe_cuts}',
2894 path=roe_path)
2895 fillParticleList('K_S0:v0_roe -> pi+ pi-', f'{selection_cuts[1]} and {roe_cuts}',
2896 path=roe_path)
2897 fillParticleList('Lambda0:v0_roe -> p+ pi-', f'{selection_cuts[2]} and {roe_cuts}',
2898 path=roe_path)
2899 fitter = fitter.lower()
2900 if (fitter != 'treefit' and fitter != 'kfit'):
2901 B2WARNING('Argument "fitter" in updateROEUsingV0Lists has only "treefit" and "kfit" options, '
2902 f'but "{fitter}" was provided! TreeFitter will be used instead.')
2903 fitter = 'treefit'
2904 from vertex import kFit, treeFit
2905 for v0 in ['gamma:v0_roe', 'K_S0:v0_roe', 'Lambda0:v0_roe']:
2906 if (apply_mass_fit and fitter == 'kfit'):
2907 kFit(v0, conf_level=0.0, fit_type='massvertex', path=roe_path)
2908 if (apply_mass_fit and fitter == 'treefit'):
2909 treeFit(v0, conf_level=0.0, massConstraint=[v0.split(':')[0]], path=roe_path)
2910 optimizeROEWithV0(v0, mask_names, '', path=roe_path)
2911 path.for_each('RestOfEvent', 'RestOfEvents', roe_path)
2912
2913
2914def printROEInfo(mask_names=None, full_print=False,
2915 unpackComposites=True, path=None):
2916 """
2917 This function prints out the information for the current ROE, so it should only be used in the for_each path.
2918 It prints out basic ROE object info.
2919
2920 If mask names are provided, specific information for those masks will be printed out.
2921
2922 It is also possible to print out all particles in a given mask if the
2923 'full_print' is set to True.
2924
2925 @param mask_names array of ROEMask names for printing out info
2926 @param unpackComposites if true, replace composite particles by their daughters
2927 @param full_print print out particles in mask
2928 @param path modules are added to this path
2929 """
2930
2931 if mask_names is None:
2932 mask_names = []
2933 printMask = register_module('RestOfEventPrinter')
2934 printMask.set_name('RestOfEventPrinter')
2935 printMask.param('maskNames', mask_names)
2936 printMask.param('fullPrint', full_print)
2937 printMask.param('unpackComposites', unpackComposites)
2938 path.add_module(printMask)
2939
2940
2941def buildContinuumSuppression(list_name, roe_mask, ipprofile_fit=False, path=None):
2942 """
2943 Creates for each Particle in the given ParticleList a ContinuumSuppression
2944 dataobject and makes basf2 relation between them.
2945
2946 :param list_name: name of the input ParticleList
2947 :param roe_mask: name of the ROE mask
2948 :param ipprofile_fit: turn on vertex fit of input tracks with IP profile constraint
2949 :param path: modules are added to this path
2950 """
2951
2952 qqBuilder = register_module('ContinuumSuppressionBuilder')
2953 qqBuilder.set_name('QQBuilder_' + list_name)
2954 qqBuilder.param('particleList', list_name)
2955 qqBuilder.param('ROEMask', roe_mask)
2956 qqBuilder.param('performIPProfileFit', ipprofile_fit)
2957 path.add_module(qqBuilder)
2958
2959
2960def removeParticlesNotInLists(lists_to_keep, path):
2961 """
2962 Removes all Particles that are not in a given list of ParticleLists (or daughters of those).
2963 All relations from/to Particles, daughter indices, and other ParticleLists are fixed.
2964
2965 @param lists_to_keep Keep the Particles and their daughters in these ParticleLists.
2966 @param path modules are added to this path
2967 """
2968
2969 mod = register_module('RemoveParticlesNotInLists')
2970 mod.param('particleLists', lists_to_keep)
2971 path.add_module(mod)
2972
2973
2974def inclusiveBtagReconstruction(upsilon_list_name, bsig_list_name, btag_list_name, input_lists_names, path):
2975 """
2976 Reconstructs Btag from particles in given ParticleLists which do not share any final state particles (mdstSource) with Bsig.
2977
2978 @param upsilon_list_name Name of the ParticleList to be filled with 'Upsilon(4S) -> B:sig anti-B:tag'
2979 @param bsig_list_name Name of the Bsig ParticleList
2980 @param btag_list_name Name of the Bsig ParticleList
2981 @param input_lists_names List of names of the ParticleLists which are used to reconstruct Btag from
2982 """
2983
2984 btag = register_module('InclusiveBtagReconstruction')
2985 btag.set_name('InclusiveBtagReconstruction_' + bsig_list_name)
2986 btag.param('upsilonListName', upsilon_list_name)
2987 btag.param('bsigListName', bsig_list_name)
2988 btag.param('btagListName', btag_list_name)
2989 btag.param('inputListsNames', input_lists_names)
2990 path.add_module(btag)
2991
2992
2993def selectDaughters(particle_list_name, decay_string, path):
2994 """
2995 Redefine the Daughters of a particle: select from decayString
2996
2997 @param particle_list_name input particle list
2998 @param decay_string for selecting the Daughters to be preserved
2999 """
3000
3001 seld = register_module('SelectDaughters')
3002 seld.set_name('SelectDaughters_' + particle_list_name)
3003 seld.param('listName', particle_list_name)
3004 seld.param('decayString', decay_string)
3005 path.add_module(seld)
3006
3007
3008def markDuplicate(particleList, prioritiseV0, path):
3009 """
3010 Call DuplicateVertexMarker to find duplicate particles in a list and
3011 flag the ones that should be kept
3012
3013 @param particleList input particle list
3014 @param prioritiseV0 if true, give V0s a higher priority
3015 """
3016
3017 markdup = register_module('DuplicateVertexMarker')
3018 markdup.param('particleList', particleList)
3019 markdup.param('prioritiseV0', prioritiseV0)
3020 path.add_module(markdup)
3021
3022
3023PI0ETAVETO_COUNTER = 0
3024
3025
3026def oldwritePi0EtaVeto(
3027 particleList,
3028 decayString,
3029 workingDirectory='.',
3030 pi0vetoname='Pi0_Prob',
3031 etavetoname='Eta_Prob',
3032 downloadFlag=True,
3033 selection='',
3034 path=None
3035):
3036 """
3037 Give pi0/eta probability for hard photon.
3038
3039 In the default weight files a value of 1.4 GeV is set as the lower limit for the hard photon energy in the CMS frame.
3040
3041 The current default weight files are optimised using MC9.
3042 The input variables are as below. Aliases are set to some variables during training.
3043
3044 * M: pi0/eta candidates Invariant mass
3045 * lowE: soft photon energy in lab frame
3046 * cTheta: soft photon ECL cluster's polar angle
3047 * Zmva: soft photon output of MVA using Zernike moments of the cluster
3048 * minC2Hdist: soft photon distance from eclCluster to nearest point on nearest Helix at the ECL cylindrical radius
3049
3050 If you don't have weight files in your workingDirectory,
3051 these files are downloaded from database to your workingDirectory automatically.
3052 Please refer to analysis/examples/tutorials/B2A306-B02RhoGamma-withPi0EtaVeto.py
3053 about how to use this function.
3054
3055 NOTE:
3056 Please don't use following ParticleList names elsewhere:
3057
3058 ``gamma:HARDPHOTON``, ``pi0:PI0VETO``, ``eta:ETAVETO``,
3059 ``gamma:PI0SOFT + str(PI0ETAVETO_COUNTER)``, ``gamma:ETASOFT + str(PI0ETAVETO_COUNTER)``
3060
3061 Please don't use ``lowE``, ``cTheta``, ``Zmva``, ``minC2Hdist`` as alias elsewhere.
3062
3063 @param particleList The input ParticleList
3064 @param decayString specify Particle to be added to the ParticleList
3065 @param workingDirectory The weight file directory
3066 @param downloadFlag whether download default weight files or not
3067 @param pi0vetoname extraInfo name of pi0 probability
3068 @param etavetoname extraInfo name of eta probability
3069 @param selection Selection criteria that Particle needs meet in order for for_each ROE path to continue
3070 @param path modules are added to this path
3071 """
3072
3073 import b2bii
3074 if b2bii.isB2BII():
3075 B2ERROR("The old pi0 / eta veto is not suitable for Belle analyses.")
3076
3077 import os
3078 import basf2_mva
3079
3080 global PI0ETAVETO_COUNTER
3081
3082 if PI0ETAVETO_COUNTER == 0:
3083 from variables import variables
3084 variables.addAlias('lowE', 'daughter(1,E)')
3085 variables.addAlias('cTheta', 'daughter(1,clusterTheta)')
3086 variables.addAlias('Zmva', 'daughter(1,clusterZernikeMVA)')
3087 variables.addAlias('minC2Tdist', 'daughter(1,minC2TDist)')
3088 variables.addAlias('cluNHits', 'daughter(1,clusterNHits)')
3089 variables.addAlias('E9E21', 'daughter(1,clusterE9E21)')
3090
3091 PI0ETAVETO_COUNTER = PI0ETAVETO_COUNTER + 1
3092
3093 roe_path = create_path()
3094
3095 deadEndPath = create_path()
3096
3097 signalSideParticleFilter(particleList, selection, roe_path, deadEndPath)
3098
3099 fillSignalSideParticleList('gamma:HARDPHOTON', decayString, path=roe_path)
3100
3101 pi0softname = 'gamma:PI0SOFT'
3102 etasoftname = 'gamma:ETASOFT'
3103 softphoton1 = pi0softname + str(PI0ETAVETO_COUNTER)
3104 softphoton2 = etasoftname + str(PI0ETAVETO_COUNTER)
3105
3106 fillParticleList(
3107 softphoton1,
3108 '[clusterReg==1 and E>0.025] or [clusterReg==2 and E>0.02] or [clusterReg==3 and E>0.02]',
3109 path=roe_path)
3110 applyCuts(softphoton1, 'abs(clusterTiming)<120', path=roe_path)
3111 fillParticleList(
3112 softphoton2,
3113 '[clusterReg==1 and E>0.035] or [clusterReg==2 and E>0.03] or [clusterReg==3 and E>0.03]',
3114 path=roe_path)
3115 applyCuts(softphoton2, 'abs(clusterTiming)<120', path=roe_path)
3116
3117 reconstructDecay('pi0:PI0VETO -> gamma:HARDPHOTON ' + softphoton1, '', path=roe_path)
3118 reconstructDecay('eta:ETAVETO -> gamma:HARDPHOTON ' + softphoton2, '', path=roe_path)
3119
3120 if not os.path.isdir(workingDirectory):
3121 os.mkdir(workingDirectory)
3122 B2INFO('oldwritePi0EtaVeto: ' + workingDirectory + ' has been created as workingDirectory.')
3123
3124 if not os.path.isfile(workingDirectory + '/pi0veto.root'):
3125 if downloadFlag:
3126 basf2_mva.download('Pi0VetoIdentifier', workingDirectory + '/pi0veto.root')
3127 B2INFO('oldwritePi0EtaVeto: pi0veto.root has been downloaded from database to workingDirectory.')
3128
3129 if not os.path.isfile(workingDirectory + '/etaveto.root'):
3130 if downloadFlag:
3131 basf2_mva.download('EtaVetoIdentifier', workingDirectory + '/etaveto.root')
3132 B2INFO('oldwritePi0EtaVeto: etaveto.root has been downloaded from database to workingDirectory.')
3133
3134 roe_path.add_module('MVAExpert', listNames=['pi0:PI0VETO'], extraInfoName='Pi0Veto',
3135 identifier=workingDirectory + '/pi0veto.root')
3136 roe_path.add_module('MVAExpert', listNames=['eta:ETAVETO'], extraInfoName='EtaVeto',
3137 identifier=workingDirectory + '/etaveto.root')
3138
3139 rankByHighest('pi0:PI0VETO', 'extraInfo(Pi0Veto)', numBest=1, path=roe_path)
3140 rankByHighest('eta:ETAVETO', 'extraInfo(EtaVeto)', numBest=1, path=roe_path)
3141
3142 variableToSignalSideExtraInfo('pi0:PI0VETO', {'extraInfo(Pi0Veto)': pi0vetoname}, path=roe_path)
3143 variableToSignalSideExtraInfo('eta:ETAVETO', {'extraInfo(EtaVeto)': etavetoname}, path=roe_path)
3144
3145 path.for_each('RestOfEvent', 'RestOfEvents', roe_path)
3146
3147
3148def writePi0EtaVeto(
3149 particleList,
3150 decayString,
3151 mode='standardMC16rd',
3152 selection='',
3153 path=None,
3154 suffix='',
3155 hardParticle='gamma',
3156 pi0PayloadNameOverride=None,
3157 pi0SoftPhotonCutOverride=None,
3158 etaPayloadNameOverride=None,
3159 etaSoftPhotonCutOverride=None,
3160 requireSoftPhotonIsInROE=False,
3161 pi0Selection='',
3162 etaSelection=''
3163):
3164 """
3165 Give pi0/eta probability for hard photon.
3166
3167 In the default weight files a value of 1.4 GeV is set as the lower limit for the hard photon energy in the CMS frame.
3168 For MC15rd/MC16rd weight files, the BtoXGamma skim is applied during the MVA training.
3169
3170 The current default weight files are for MC16rd. The weight files for MC15rd/MC12 are still available.
3171
3172 The input variables of the mva training for pi0 veto using MC16rd are:
3173
3174 * M: Invariant mass of pi0 candidates
3175 * cosHelicityAngleMomentum: Cosine of angle between momentum difference of the photons in the pi0 rest frame
3176 and momentum of pi0 in lab frame
3177 * daughter(1,E): soft photon energy in lab frame
3178 * daughter(1,clusterTheta): soft photon ECL cluster's polar angle
3179 * daughter(1,clusterLAT): soft photon lateral energy distribution
3180 * daughter(1,beamBackgroundSuppression): soft photon beam background suppression MVA output
3181 * daughter(1,fakePhotonSuppression): soft photon fake photon suppression MVA output
3182
3183 The input variables of the mva training for eta veto using MC16rd are:
3184
3185 * M: Invariant mass of eta candidates
3186 * cosHelicityAngleMomentum: Cosine of angle between momentum difference of the photons in the eta rest frame
3187 and momentum of eta in lab frame
3188 * daughter(1,E): soft photon energy in lab frame
3189 * daughter(1,clusterTheta): soft photon ECL cluster's polar angle
3190 * daughter(1,clusterLAT): soft photon lateral energy distribution
3191 * daughter(1,clusterNHits): soft photon total crystal weights sum(w_i) with w_i<=1
3192 * daughter(1,clusterE1E9): soft photon ratio between energies of central crystal and inner 3x3 crystals
3193 * daughter(1,clusterE9E21): soft photon ratio of energies in inner 3x3 crystals and 5x5 crystals without corners
3194 * daughter(1,clusterSecondMoment): soft photon second moment
3195 * daughter(1,clusterAbsZernikeMoment40): soft photon Zernike moment 40
3196 * daughter(1,clusterAbsZernikeMoment51): soft photon Zernike moment 51
3197 * daughter(1,beamBackgroundSuppression): soft photon beam background suppression MVA output
3198 * daughter(1,fakePhotonSuppression): soft photon fake photon suppression MVA output
3199
3200
3201 The input variables of the mva training for pi0 veto using MC15rd are:
3202
3203 * M: Invariant mass of pi0 candidates
3204 * cosHelicityAngleMomentum: Cosine of angle between momentum difference of the photons in the pi0 rest frame
3205 and momentum of pi0 in lab frame
3206 * daughter(1,E): soft photon energy in lab frame
3207 * daughter(1,clusterTheta): soft photon ECL cluster's polar angle
3208 * daughter(1,clusterLAT): soft photon lateral energy distribution
3209
3210 The input variables of the mva training for eta veto using MC15rd are:
3211
3212 * M: Invariant mass of eta candidates
3213 * cosHelicityAngleMomentum: Cosine of angle between momentum difference of the photons in the eta rest frame
3214 and momentum of eta in lab frame
3215 * daughter(1,E): soft photon energy in lab frame
3216 * daughter(1,clusterTheta): soft photon ECL cluster's polar angle
3217 * daughter(1,clusterLAT): soft photon lateral energy distribution
3218 * daughter(1,clusterNHits): soft photon total crystal weights sum(w_i) with w_i<=1
3219 * daughter(1,clusterE1E9): soft photon ratio between energies of central crystal and inner 3x3 crystals
3220 * daughter(1,clusterE9E21): soft photon ratio of energies in inner 3x3 crystals and 5x5 crystals without corners
3221 * daughter(1,clusterSecondMoment): soft photon second moment
3222 * daughter(1,clusterAbsZernikeMoment40): soft photon Zernike moment 40
3223 * daughter(1,clusterAbsZernikeMoment51): soft photon Zernike moment 51
3224
3225 The input variables of the mva training using MC12 are:
3226
3227 * M: Invariant mass of pi0/eta candidates
3228 * daughter(1,E): soft photon energy in lab frame
3229 * daughter(1,clusterTheta): soft photon ECL cluster's polar angle
3230 * daughter(1,minC2TDist): soft photon distance from eclCluster to nearest point on nearest Helix at the ECL cylindrical radius
3231 * daughter(1,clusterZernikeMVA): soft photon output of MVA using Zernike moments of the cluster
3232 * daughter(1,clusterNHits): soft photon total crystal weights sum(w_i) with w_i<=1
3233 * daughter(1,clusterE9E21): soft photon ratio of energies in inner 3x3 crystals and 5x5 crystals without corners
3234 * cosHelicityAngleMomentum: Cosine of angle between momentum difference of the photons in the pi0/eta rest frame
3235 and momentum of pi0/eta in lab frame
3236
3237 The following strings are available for mode:
3238
3239 * standard: loose energy cut and no clusterNHits cut are applied to soft photon
3240 * tight: tight energy cut and no clusterNHits cut are applied to soft photon
3241 * cluster: loose energy cut and clusterNHits cut are applied to soft photon
3242 * both: tight energy cut and clusterNHits cut are applied to soft photon
3243 * standardMC15rd: loose energy cut is applied to soft photon and the weight files are trained using MC15rd
3244 * tightMC15rd: tight energy cut is applied to soft photon and the weight files are trained using MC15rd
3245 * standardMC16rd: loose energy cut is applied to soft photon and the weight files are trained using MC16rd
3246 * tightMC16rd: tight energy cut is applied to soft photon and the weight files are trained using MC16rd
3247
3248 The final probability of the pi0/eta veto is stored as an extraInfo. If no suffix is set it can be obtained from the variables
3249 `pi0Prob`/`etaProb`. Otherwise, it is available as '{Pi0, Eta}ProbOrigin', '{Pi0, Eta}ProbTightEnergyThreshold', '{Pi0,
3250 Eta}ProbLargeClusterSize', '{Pi0, Eta}ProbTightEnergyThresholdAndLargeClusterSize', '{Pi0, Eta}ProbOriginMC15rd', or
3251 '{Pi0, Eta}ProbTightEnergyThresholdMC15rd' for the six modes described above, with the chosen suffix appended. If one would
3252 like to call this veto twice in one script, add suffix in the second time!
3253 The second highest probability of the pi0/eta veto also is stored as an extraInfo, with a prefix of 'second' to the previous
3254 ones, e.g. secondPi0ProbOrigin{suffix}. This can be used to do validation/systematics study.
3255
3256 NOTE:
3257 Please don't use following ParticleList names elsewhere:
3258
3259 ``gamma:HardPhoton``,
3260 ``gamma:Pi0Soft + ListName + '_' + particleList.replace(':', '_')``,
3261 ``gamma:EtaSoft + ListName + '_' + particleList.replace(':', '_')``,
3262 ``pi0:EtaVeto + ListName``,
3263 ``eta:EtaVeto + ListName``
3264
3265 @param particleList the input ParticleList
3266 @param decayString specify Particle to be added to the ParticleList
3267 @param mode choose one mode out of 'standardMC16rd', 'tightMC16rd', 'standardMC15rd', 'tightMC15rd',
3268 'standard', 'tight', 'cluster' and 'both'
3269 @param selection selection criteria that Particle needs meet in order for for_each ROE path to continue
3270 @param path modules are added to this path
3271 @param suffix optional suffix to be appended to the usual extraInfo name
3272 @param hardParticle particle name which is used to calculate the pi0/eta probability (default is gamma)
3273 @param pi0PayloadNameOverride specify the payload name of pi0 veto only if one wants to use non-default one. (default is None)
3274 @param pi0SoftPhotonCutOverride specify the soft photon selection criteria of pi0 veto only if one wants to use non-default one.
3275 (default is None)
3276 @param etaPayloadNameOverride specify the payload name of eta veto only if one wants to use non-default one. (default is None)
3277 @param etaSoftPhotonCutOverride specify the soft photon selection criteria of eta veto only if one wants to use non-default one.
3278 (default is None)
3279 @param requireSoftPhotonIsInROE specify if the soft photons used to build pi0 and eta candidates have to be in the current ROE
3280 or not. Default is False, i.e. all soft photons in the event are used.
3281 @param pi0Selection Selection for the pi0 reconstruction. Default is "".
3282 @param etaSelection Selection for the eta reconstruction. Default is "".
3283 """
3284
3285 import b2bii
3286 if b2bii.isB2BII():
3287 B2ERROR("The pi0 / eta veto is not suitable for Belle analyses.")
3288
3289 if (requireSoftPhotonIsInROE):
3290 B2WARNING("Requiring the soft photon to being in the ROE was not done for the MVA training. "
3291 "Please check the results carefully.")
3292 showWarning = False
3293
3294 if (mode == 'standardMC15rd' or mode == 'tightMC15rd'):
3295 if (pi0Selection != '[0.03 < M < 0.23]' or etaSelection != '[0.25 < M < 0.75]'):
3296 showWarning = True
3297 else:
3298 if (pi0Selection != '' or etaSelection != ''):
3299 showWarning = True
3300 if showWarning:
3301 B2WARNING(
3302 "Selection criteria for the pi0 or the eta during reconstructDecay differ from those used during the MVA training. "
3303 "You may get NAN value. Please check the results carefully.")
3304
3305 renameSuffix = False
3306
3307 for module in path.modules():
3308 if module.type() == "SubEvent" and not renameSuffix:
3309 for subpath in [p.values for p in module.available_params() if p.name == "path"]:
3310 if renameSuffix:
3311 break
3312 for submodule in subpath.modules():
3313 if f'{hardParticle}:HardPhoton{suffix}' in submodule.name():
3314 suffix += '_0'
3315 B2WARNING("Same extension already used in writePi0EtaVeto, append '_0'")
3316 renameSuffix = True
3317 break
3318
3319 roe_path = create_path()
3320 deadEndPath = create_path()
3321 signalSideParticleFilter(particleList, selection, roe_path, deadEndPath)
3322 fillSignalSideParticleList(f'{hardParticle}:HardPhoton{suffix}', decayString, path=roe_path)
3323
3324 dictListName = {'standard': 'Origin',
3325 'tight': 'TightEnergyThreshold',
3326 'cluster': 'LargeClusterSize',
3327 'both': 'TightEnrgyThresholdAndLargeClusterSize',
3328 'standardMC15rd': 'OriginMC15rd',
3329 'tightMC15rd': 'TightEnergyThresholdMC15rd',
3330 'standardMC16rd': 'OriginMC16rd',
3331 'tightMC16rd': 'TightEnergyThresholdMC16rd'}
3332
3333 dictPi0EnergyCut = {
3334 'standard': '[[clusterReg==1 and E>0.025] or [clusterReg==2 and E>0.02] or [clusterReg==3 and E>0.02]]',
3335 'tight': '[[clusterReg==1 and E>0.03] or [clusterReg==2 and E>0.03] or [clusterReg==3 and E>0.04]]',
3336 'cluster': '[[clusterReg==1 and E>0.025] or [clusterReg==2 and E>0.02] or [clusterReg==3 and E>0.02]]',
3337 'both': '[[clusterReg==1 and E>0.03] or [clusterReg==2 and E>0.03] or [clusterReg==3 and E>0.04]]',
3338 'standardMC15rd': '[[clusterReg==1 and E>0.0225] or [clusterReg==2 and E>0.02] or [clusterReg==3 and E>0.02]]',
3339 'tightMC15rd': '[[clusterReg==1 and E>0.03] or [clusterReg==2 and E>0.03] or [clusterReg==3 and E>0.04]]',
3340 'standardMC16rd': '[[clusterReg==1 and E>0.0225] or [clusterReg==2 and E>0.02] or [clusterReg==3 and E>0.02]]',
3341 'tightMC16rd': '[[clusterReg==1 and E>0.03] or [clusterReg==2 and E>0.03] or [clusterReg==3 and E>0.04]]'}
3342
3343 dictEtaEnergyCut = {
3344 'standard': '[[clusterReg==1 and E>0.035] or [clusterReg==2 and E>0.03] or [clusterReg==3 and E>0.03]]',
3345 'tight': '[[clusterReg==1 and E>0.06] or [clusterReg==2 and E>0.06] or [clusterReg==3 and E>0.06]]',
3346 'cluster': '[[clusterReg==1 and E>0.035] or [clusterReg==2 and E>0.03] or [clusterReg==3 and E>0.03]]',
3347 'both': '[[clusterReg==1 and E>0.06] or [clusterReg==2 and E>0.06] or [clusterReg==3 and E>0.06]]',
3348 'standardMC15rd': '[[clusterReg==1 and E>0.0225] or [clusterReg==2 and E>0.02] or [clusterReg==3 and E>0.02]]',
3349 'tightMC15rd': '[[clusterReg==1 and E>0.03] or [clusterReg==2 and E>0.03] or [clusterReg==3 and E>0.04]]',
3350 'standardMC16rd': '[[clusterReg==1 and E>0.0225] or [clusterReg==2 and E>0.02] or [clusterReg==3 and E>0.02]]',
3351 'tightMC16rd': '[[clusterReg==1 and E>0.03] or [clusterReg==2 and E>0.03] or [clusterReg==3 and E>0.04]]'}
3352
3353 dictNHitsTimingCut = {'standard': 'clusterNHits >= 0 and abs(clusterTiming)<clusterErrorTiming',
3354 'tight': 'clusterNHits >= 0 and abs(clusterTiming)<clusterErrorTiming',
3355 'cluster': 'clusterNHits >= 2 and abs(clusterTiming)<clusterErrorTiming',
3356 'both': 'clusterNHits >= 2 and abs(clusterTiming)<clusterErrorTiming',
3357 'standardMC15rd': 'clusterNHits > 1.5 and abs(clusterTiming) < 200',
3358 'tightMC15rd': 'clusterNHits > 1.5 and abs(clusterTiming) < 200',
3359 'standardMC16rd': 'clusterNHits > 1.5 and abs(clusterTiming) < 200',
3360 'tightMC16rd': 'clusterNHits > 1.5 and abs(clusterTiming) < 200'}
3361
3362 dictPi0PayloadName = {'standard': 'Pi0VetoIdentifierStandard',
3363 'tight': 'Pi0VetoIdentifierWithHigherEnergyThreshold',
3364 'cluster': 'Pi0VetoIdentifierWithLargerClusterSize',
3365 'both': 'Pi0VetoIdentifierWithHigherEnergyThresholdAndLargerClusterSize',
3366 'standardMC15rd': 'Pi0VetoIdentifierStandardMC15rd',
3367 'tightMC15rd': 'Pi0VetoIdentifierWithHigherEnergyThresholdMC15rd',
3368 'standardMC16rd': 'Pi0VetoIdentifierStandardMC16rd',
3369 'tightMC16rd': 'Pi0VetoIdentifierWithHigherEnergyThresholdMC16rd'}
3370
3371 dictEtaPayloadName = {'standard': 'EtaVetoIdentifierStandard',
3372 'tight': 'EtaVetoIdentifierWithHigherEnergyThreshold',
3373 'cluster': 'EtaVetoIdentifierWithLargerClusterSize',
3374 'both': 'EtaVetoIdentifierWithHigherEnergyThresholdAndLargerClusterSize',
3375 'standardMC15rd': 'EtaVetoIdentifierStandardMC15rd',
3376 'tightMC15rd': 'EtaVetoIdentifierWithHigherEnergyThresholdMC15rd',
3377 'standardMC16rd': 'EtaVetoIdentifierStandardMC16rd',
3378 'tightMC16rd': 'EtaVetoIdentifierWithHigherEnergyThresholdMC16rd'}
3379
3380 dictPi0ExtraInfoName = {'standard': 'Pi0ProbOrigin',
3381 'tight': 'Pi0ProbTightEnergyThreshold',
3382 'cluster': 'Pi0ProbLargeClusterSize',
3383 'both': 'Pi0ProbTightEnergyThresholdAndLargeClusterSize',
3384 'standardMC15rd': 'Pi0ProbOriginMC15rd',
3385 'tightMC15rd': 'Pi0ProbTightEnergyThresholdMC15rd',
3386 'standardMC16rd': 'Pi0ProbOriginMC16rd',
3387 'tightMC16rd': 'Pi0ProbTightEnergyThresholdMC16rd'}
3388
3389 dictEtaExtraInfoName = {'standard': 'EtaProbOrigin',
3390 'tight': 'EtaProbTightEnergyThreshold',
3391 'cluster': 'EtaProbLargeClusterSize',
3392 'both': 'EtaProbTightEnergyThresholdAndLargeClusterSize',
3393 'standardMC15rd': 'EtaProbOriginMC15rd',
3394 'tightMC15rd': 'EtaProbTightEnergyThresholdMC15rd',
3395 'standardMC16rd': 'EtaProbOriginMC16rd',
3396 'tightMC16rd': 'EtaProbTightEnergyThresholdMC16rd'}
3397
3398 ListName = dictListName[mode]
3399 Pi0EnergyCut = dictPi0EnergyCut[mode]
3400 EtaEnergyCut = dictEtaEnergyCut[mode]
3401 NHitsTimingCut = dictNHitsTimingCut[mode]
3402 Pi0PayloadName = dictPi0PayloadName[mode]
3403 EtaPayloadName = dictEtaPayloadName[mode]
3404 Pi0ExtraInfoName = dictPi0ExtraInfoName[mode]
3405 EtaExtraInfoName = dictEtaExtraInfoName[mode]
3406
3407 # pi0 veto
3408 if pi0PayloadNameOverride is not None:
3409 Pi0PayloadName = pi0PayloadNameOverride
3410 B2WARNING("You're using personal weight files, be careful. ")
3411 if pi0SoftPhotonCutOverride is None:
3412 Pi0SoftPhotonCut = Pi0EnergyCut + ' and ' + NHitsTimingCut
3413 else:
3414 Pi0SoftPhotonCut = pi0SoftPhotonCutOverride
3415 B2WARNING("You're applying personal cuts on the soft photon candidates, be careful. ")
3416
3417 if requireSoftPhotonIsInROE:
3418 Pi0SoftPhotonCut += ' and isInRestOfEvent==1'
3419
3420 # define the particleList name for soft photon
3421 pi0soft = f'gamma:Pi0Soft{suffix}' + ListName + '_' + particleList.replace(':', '_')
3422 # fill the particleList for soft photon with energy, timing and clusterNHits cuts
3423 fillParticleList(pi0soft, Pi0SoftPhotonCut, path=roe_path)
3424 # reconstruct pi0
3425 reconstructDecay('pi0:Pi0Veto' + ListName + suffix + f' -> {hardParticle}:HardPhoton{suffix} ' + pi0soft, pi0Selection,
3426 allowChargeViolation=True, path=roe_path)
3427 # MVA training is conducted.
3428 roe_path.add_module('MVAExpert', listNames=['pi0:Pi0Veto' + ListName + suffix],
3429 extraInfoName=Pi0ExtraInfoName, identifier=Pi0PayloadName)
3430 # Pick up the pi0/eta candidate with the highest pi0/eta probability.
3431 rankByHighest(
3432 'pi0:Pi0Veto' + ListName + suffix,
3433 'extraInfo(' + Pi0ExtraInfoName + ')',
3434 numBest=2,
3435 outputVariable="Pi0VetoRank",
3436 path=roe_path)
3437 cutAndCopyList(outputListName='pi0:Pi0VetoFirst' + ListName + suffix,
3438 inputListName='pi0:Pi0Veto' + ListName + suffix,
3439 cut='extraInfo(Pi0VetoRank)==1',
3440 path=roe_path)
3441 variableToSignalSideExtraInfo('pi0:Pi0VetoFirst' + ListName + suffix,
3442 {'extraInfo(' + Pi0ExtraInfoName + ')': Pi0ExtraInfoName + suffix}, path=roe_path)
3443 # Pick up the pi0/eta candidate with the second highest pi0/eta probability.
3444 cutAndCopyList(outputListName='pi0:Pi0VetoSecond' + ListName + suffix,
3445 inputListName='pi0:Pi0Veto' + ListName + suffix,
3446 cut='extraInfo(Pi0VetoRank)==2',
3447 path=roe_path)
3448 variableToSignalSideExtraInfo('pi0:Pi0VetoSecond' + ListName + suffix,
3449 {'extraInfo(' + Pi0ExtraInfoName + ')': 'second' + Pi0ExtraInfoName + suffix}, path=roe_path)
3450
3451 # eta veto
3452 if etaPayloadNameOverride is not None:
3453 EtaPayloadName = etaPayloadNameOverride
3454 B2WARNING("You're using personal weight files, be careful. ")
3455 if etaSoftPhotonCutOverride is None:
3456 EtaSoftPhotonCut = EtaEnergyCut + ' and ' + NHitsTimingCut
3457 else:
3458 EtaSoftPhotonCut = etaSoftPhotonCutOverride
3459 B2WARNING("You're applying personal cuts on the soft photon candidates, be careful. ")
3460
3461 if requireSoftPhotonIsInROE:
3462 EtaSoftPhotonCut += ' and isInRestOfEvent==1'
3463
3464 etasoft = f'gamma:EtaSoft{suffix}' + ListName + '_' + particleList.replace(':', '_')
3465 fillParticleList(etasoft, EtaSoftPhotonCut, path=roe_path)
3466 reconstructDecay('eta:EtaVeto' + ListName + suffix + f' -> {hardParticle}:HardPhoton{suffix} ' + etasoft, etaSelection,
3467 allowChargeViolation=True, path=roe_path)
3468 roe_path.add_module('MVAExpert', listNames=['eta:EtaVeto' + ListName + suffix],
3469 extraInfoName=EtaExtraInfoName, identifier=EtaPayloadName)
3470 rankByHighest(
3471 'eta:EtaVeto' + ListName + suffix,
3472 'extraInfo(' + EtaExtraInfoName + ')',
3473 numBest=2,
3474 outputVariable="EtaVetoRank",
3475 path=roe_path)
3476 cutAndCopyList(outputListName='eta:EtaVetoFirst' + ListName + suffix,
3477 inputListName='eta:EtaVeto' + ListName + suffix,
3478 cut='extraInfo(EtaVetoRank)==1',
3479 path=roe_path)
3480 variableToSignalSideExtraInfo('eta:EtaVetoFirst' + ListName + suffix,
3481 {'extraInfo(' + EtaExtraInfoName + ')': EtaExtraInfoName + suffix}, path=roe_path)
3482 cutAndCopyList(outputListName='eta:EtaVetoSecond' + ListName + suffix,
3483 inputListName='eta:EtaVeto' + ListName + suffix,
3484 cut='extraInfo(EtaVetoRank)==2',
3485 path=roe_path)
3486 variableToSignalSideExtraInfo('eta:EtaVetoSecond' + ListName + suffix,
3487 {'extraInfo(' + EtaExtraInfoName + ')': 'second' + EtaExtraInfoName + suffix}, path=roe_path)
3488
3489 path.for_each('RestOfEvent', 'RestOfEvents', roe_path)
3490
3491
3492def lowEnergyPi0Identification(pi0List, gammaList, payloadNameSuffix,
3493 path=None):
3494 """
3495 Calculate low-energy pi0 identification.
3496 The result is stored as ExtraInfo ``lowEnergyPi0Identification`` for
3497 the list pi0List.
3498
3499 Parameters:
3500 pi0List (str): Pi0 list.
3501
3502 gammaList (str): Gamma list. First, an energy cut E > 0.2 is applied to the photons from this list.
3503 Then, all possible combinations with a pi0 daughter photon are formed except the one
3504 corresponding to the reconstructed pi0.
3505 The maximum low-energy pi0 veto value is calculated for such photon pairs
3506 and used as one of the input variables for the identification classifier.
3507
3508 payloadNameSuffix (str): Payload name suffix. The weight payloads are stored in the analysis global
3509 tag and have the following names:\n
3510 * ``'LowEnergyPi0Veto' + payloadNameSuffix``
3511 * ``'LowEnergyPi0Identification' + payloadNameSuffix``\n
3512 The possible suffixes are:\n
3513 * ``'Belle1'`` for Belle data.
3514 * ``'Belle2Release5'`` for Belle II release 5 data (MC14, proc12, buckets 16 - 25).
3515 * ``'Belle2Release6'`` for Belle II release 6 data (MC15, proc13, buckets 26 - 36).
3516
3517 path (basf2.Path): Module path.
3518 """
3519
3520 # Select photons with higher energy for formation of veto combinations.
3521 gammaListVeto = f'{gammaList}_pi0veto'
3522 cutAndCopyList(gammaListVeto, gammaList, 'E > 0.2', path=path)
3523 import b2bii
3524 payload_name = 'LowEnergyPi0Veto' + payloadNameSuffix
3525 path.add_module('LowEnergyPi0VetoExpert', identifier=payload_name,
3526 VetoPi0Daughters=True, GammaListName=gammaListVeto,
3527 Pi0ListName=pi0List, Belle1=b2bii.isB2BII())
3528 payload_name = 'LowEnergyPi0Identification' + payloadNameSuffix
3529 path.add_module('LowEnergyPi0IdentificationExpert',
3530 identifier=payload_name, Pi0ListName=pi0List,
3531 Belle1=b2bii.isB2BII())
3532
3533
3534def getNeutralHadronGeomMatches(
3535 particleLists,
3536 addKL=True,
3537 addNeutrons=False,
3538 efficiencyCorrectionKl=0.83,
3539 efficiencyCorrectionNeutrons=1.0,
3540 path=None):
3541 """
3542 For an ECL-based list, assign the mcdistanceKL and mcdistanceNeutron variables that correspond
3543 to the distance to the closest MC KL and neutron, respectively.
3544 @param particleLists the input ParticleLists, must be ECL-based lists (e.g. photons)
3545 @param addKL (default True) add distance to MC KL
3546 @param addNeutrons (default False) add distance to MC neutrons
3547 @param efficiencyCorrectionKl (default 0.83) apply overall efficiency correction
3548 @param efficiencyCorrectionNeutrons (default 1.0) apply overall efficiency correction
3549 @param path modules are added to this path
3550 """
3551 from ROOT import Belle2
3552 Const = Belle2.Const
3553
3554 if addKL:
3555 path.add_module(
3556 "NeutralHadronMatcher",
3557 particleLists=particleLists,
3558 mcPDGcode=Const.Klong.getPDGCode(),
3559 efficiencyCorrection=efficiencyCorrectionKl)
3560 if addNeutrons:
3561 path.add_module(
3562 "NeutralHadronMatcher",
3563 particleLists=particleLists,
3564 mcPDGcode=Const.neutron.getPDGCode(),
3565 efficiencyCorrection=efficiencyCorrectionNeutrons)
3566
3567
3568def getBeamBackgroundProbability(particleList, weight, path=None):
3569 """
3570 Assign a probability to each ECL cluster as being signal like (1) compared to beam background like (0)
3571 @param particleList the input ParticleList, must be a photon list
3572 @param weight type of weight file to use
3573 @param path modules are added to this path
3574 """
3575
3576 import b2bii
3577 if b2bii.isB2BII() and weight != "Belle":
3578 B2WARNING("weight type must be 'Belle' for b2bii.")
3579
3580 path.add_module('MVAExpert',
3581 listNames=particleList,
3582 extraInfoName='beamBackgroundSuppression',
3583 identifier=f'BeamBackgroundMVA_{weight}')
3584
3585
3586def getFakePhotonProbability(particleList, weight, path=None):
3587 """
3588 Assign a probability to each ECL cluster as being signal like (1) compared to fake photon like (0)
3589 @param particleList the input ParticleList, must be a photon list
3590 @param weight type of weight file to use
3591 @param path modules are added to this path
3592 """
3593
3594 import b2bii
3595 if b2bii.isB2BII() and weight != "Belle":
3596 B2WARNING("weight type must be 'Belle' for b2bii.")
3597
3598 path.add_module('MVAExpert',
3599 listNames=particleList,
3600 extraInfoName='fakePhotonSuppression',
3601 identifier=f'FakePhotonMVA_{weight}')
3602
3603
3604def buildEventKinematics(inputListNames=None, default_cleanup=True, custom_cuts=None,
3605 chargedPIDPriors=None, fillWithMostLikely=False, path=None):
3606 """
3607 Calculates the global kinematics of the event (visible energy, missing momentum, missing mass...)
3608 using ParticleLists provided. If no ParticleList is provided, default ParticleLists are used
3609 (all track and all hits in ECL without associated track).
3610
3611 The visible energy missing values are
3612 stored in a EventKinematics dataobject.
3613
3614 @param inputListNames list of ParticleLists used to calculate the global event kinematics.
3615 If the list is empty, default ParticleLists pi+:evtkin and gamma:evtkin are filled.
3616 @param fillWithMostLikely if True, the module uses the most likely particle mass hypothesis for charged particles
3617 according to the PID likelihood and the option inputListNames will be ignored.
3618 @param chargedPIDPriors The prior PID fractions, that are used to regulate
3619 amount of certain charged particle species, should be a list of
3620 six floats if not None. The order of particle types is
3621 the following: [e-, mu-, pi-, K-, p+, d+]
3622 @param default_cleanup if True and either inputListNames empty or fillWithMostLikely True, default clean up cuts are applied
3623 @param custom_cuts tuple of selection cut strings of form (trackCuts, photonCuts), default is None,
3624 which would result in a standard predefined selection cuts
3625 @param path modules are added to this path
3626 """
3627
3628 if inputListNames is None:
3629 inputListNames = []
3630 trackCuts = 'pt > 0.1'
3631 trackCuts += ' and thetaInCDCAcceptance'
3632 trackCuts += ' and abs(dz) < 3'
3633 trackCuts += ' and dr < 0.5'
3634
3635 gammaCuts = 'E > 0.05'
3636 gammaCuts += ' and thetaInCDCAcceptance'
3637 if not b2bii.isB2BII():
3638 gammaCuts += ' and abs(clusterTiming) < 200'
3639 if (custom_cuts is not None):
3640 trackCuts, gammaCuts = custom_cuts
3641
3642 if fillWithMostLikely:
3643 from stdCharged import stdMostLikely
3644 stdMostLikely(chargedPIDPriors, '_evtkin', path=path)
3645 inputListNames = [f'{ptype}:mostlikely_evtkin' for ptype in ['K+', 'p+', 'e+', 'mu+', 'pi+']]
3646 if b2bii.isB2BII():
3647 copyList('gamma:evtkin', 'gamma:mdst', path=path)
3648 else:
3649 fillParticleList('gamma:evtkin', '', path=path)
3650 inputListNames += ['gamma:evtkin']
3651 if default_cleanup:
3652 B2INFO("Using default cleanup in EventKinematics module.")
3653 for ptype in ['K+', 'p+', 'e+', 'mu+', 'pi+']:
3654 applyCuts(f'{ptype}:mostlikely_evtkin', trackCuts, path=path)
3655 applyCuts('gamma:evtkin', gammaCuts, path=path)
3656 else:
3657 B2INFO("No cleanup in EventKinematics module.")
3658 if not inputListNames:
3659 B2INFO("Creating particle lists pi+:evtkin and gamma:evtkin to get the global kinematics of the event.")
3660 fillParticleList('pi+:evtkin', '', path=path)
3661 if b2bii.isB2BII():
3662 copyList('gamma:evtkin', 'gamma:mdst', path=path)
3663 else:
3664 fillParticleList('gamma:evtkin', '', path=path)
3665 particleLists = ['pi+:evtkin', 'gamma:evtkin']
3666 if default_cleanup:
3667 if (custom_cuts is not None):
3668 B2INFO("Using default cleanup in EventKinematics module.")
3669 applyCuts('pi+:evtkin', trackCuts, path=path)
3670 applyCuts('gamma:evtkin', gammaCuts, path=path)
3671 else:
3672 B2INFO("No cleanup in EventKinematics module.")
3673 else:
3674 particleLists = inputListNames
3675
3676 eventKinematicsModule = register_module('EventKinematics')
3677 eventKinematicsModule.set_name('EventKinematics_reco')
3678 eventKinematicsModule.param('particleLists', particleLists)
3679 path.add_module(eventKinematicsModule)
3680
3681
3682def buildEventKinematicsFromMC(inputListNames=None, selectionCut='', path=None):
3683 """
3684 Calculates the global kinematics of the event (visible energy, missing momentum, missing mass...)
3685 using generated particles. If no ParticleList is provided, default generated ParticleLists are used.
3686
3687 @param inputListNames list of ParticleLists used to calculate the global event kinematics.
3688 If the list is empty, default ParticleLists are filled.
3689 @param selectionCut optional selection cuts
3690 @param path Path to append the eventKinematics module to.
3691 """
3692
3693 if inputListNames is None:
3694 inputListNames = []
3695 if (len(inputListNames) == 0):
3696 # Type of particles to use for EventKinematics
3697 # K_S0 and Lambda0 are added here because some of them have interacted
3698 # with the detector material
3699 types = ['gamma', 'e+', 'mu+', 'pi+', 'K+', 'p+',
3700 'K_S0', 'Lambda0']
3701 for t in types:
3702 fillParticleListFromMC(f"{t}:evtkin_default_gen", 'mcPrimary > 0 and nDaughters == 0',
3703 True, True, path=path)
3704 if (selectionCut != ''):
3705 applyCuts(f"{t}:evtkin_default_gen", selectionCut, path=path)
3706 inputListNames += [f"{t}:evtkin_default_gen"]
3707
3708 eventKinematicsModule = register_module('EventKinematics')
3709 eventKinematicsModule.set_name('EventKinematics_gen')
3710 eventKinematicsModule.param('particleLists', inputListNames)
3711 eventKinematicsModule.param('usingMC', True)
3712 path.add_module(eventKinematicsModule)
3713
3714
3715def buildEventShape(inputListNames=None,
3716 default_cleanup=True,
3717 custom_cuts=None,
3718 allMoments=False,
3719 cleoCones=True,
3720 collisionAxis=True,
3721 foxWolfram=True,
3722 harmonicMoments=True,
3723 jets=True,
3724 sphericity=True,
3725 thrust=True,
3726 checkForDuplicates=False,
3727 path=None):
3728 """
3729 Calculates the event-level shape quantities (thrust, sphericity, Fox-Wolfram moments...)
3730 using the particles in the lists provided by the user. If no particle list is provided,
3731 the function will internally create a list of good tracks and a list of good photons
3732 with (optionally) minimal quality cuts.
3733
3734
3735 The results of the calculation are then stored into the EventShapeContainer dataobject,
3736 and are accessible using the variables of the EventShape group.
3737
3738 The user can switch the calculation of certain quantities on or off to save computing
3739 time. By default the calculation of the high-order moments (5-8) is turned off.
3740 Switching off an option will make the corresponding variables not available.
3741
3742 Info:
3743 The user can provide as many particle lists as needed, using also composite particles.
3744 In these cases, it is recommended to activate the checkForDuplicates flag since it
3745 will eliminate duplicates, e.g., if the same track is provided multiple times
3746 (either with different mass hypothesis or once as an independent particle and once
3747 as daughter of a composite particle). The first occurrence will be used in the
3748 calculations so the order in which the particle lists are given as well as within
3749 the particle lists matters.
3750
3751 @param inputListNames List of ParticleLists used to calculate the
3752 event shape variables. If the list is empty the default
3753 particleLists pi+:evtshape and gamma:evtshape are filled.
3754 @param default_cleanup If True, applies standard cuts on pt and cosTheta when
3755 defining the internal lists. This option is ignored if the
3756 particleLists are provided by the user.
3757 @param custom_cuts tuple of selection cut strings of form (trackCuts, photonCuts), default is None,
3758 which would result in a standard predefined selection cuts
3759 @param path Path to append the eventShape modules to.
3760 @param thrust Enables the calculation of thrust-related quantities (CLEO
3761 cones, Harmonic moments, jets).
3762 @param collisionAxis Enables the calculation of the quantities related to the
3763 collision axis .
3764 @param foxWolfram Enables the calculation of the Fox-Wolfram moments.
3765 @param harmonicMoments Enables the calculation of the Harmonic moments with respect
3766 to both the thrust axis and, if collisionAxis = True, the collision axis.
3767 @param allMoments If True, calculates also the FW and harmonic moments from order
3768 5 to 8 instead of the low-order ones only.
3769 @param cleoCones Enables the calculation of the CLEO cones with respect to both the thrust
3770 axis and, if collisionAxis = True, the collision axis.
3771 @param jets Enables the calculation of the hemisphere momenta and masses.
3772 Requires thrust = True.
3773 @param sphericity Enables the calculation of the sphericity-related quantities.
3774 @param checkForDuplicates Perform a check for duplicate particles before adding them. Regardless of the value of this option,
3775 it is recommended to consider sanitizing the lists you are passing to the function since this will
3776 speed up the processing.
3777
3778 """
3779
3780 if inputListNames is None:
3781 inputListNames = []
3782 trackCuts = 'pt > 0.1'
3783 trackCuts += ' and thetaInCDCAcceptance'
3784 trackCuts += ' and abs(dz) < 3.0'
3785 trackCuts += ' and dr < 0.5'
3786
3787 gammaCuts = 'E > 0.05'
3788 gammaCuts += ' and thetaInCDCAcceptance'
3789 if not b2bii.isB2BII():
3790 gammaCuts += ' and abs(clusterTiming) < 200'
3791 if (custom_cuts is not None):
3792 trackCuts, gammaCuts = custom_cuts
3793
3794 if not inputListNames:
3795 B2INFO("Creating particle lists pi+:evtshape and gamma:evtshape to get the event shape variables.")
3796 fillParticleList('pi+:evtshape', '', path=path)
3797 if b2bii.isB2BII():
3798 copyList('gamma:evtshape', 'gamma:mdst', path=path)
3799 else:
3800 fillParticleList(
3801 'gamma:evtshape',
3802 '',
3803 path=path)
3804 particleLists = ['pi+:evtshape', 'gamma:evtshape']
3805
3806 if default_cleanup:
3807 if (custom_cuts is not None):
3808 B2INFO("Applying standard cuts")
3809 applyCuts('pi+:evtshape', trackCuts, path=path)
3810
3811 applyCuts('gamma:evtshape', gammaCuts, path=path)
3812 else:
3813 B2WARNING("Creating the default lists with no cleanup.")
3814 else:
3815 particleLists = inputListNames
3816
3817 eventShapeModule = register_module('EventShapeCalculator')
3818 eventShapeModule.set_name('EventShape')
3819 eventShapeModule.param('particleListNames', particleLists)
3820 eventShapeModule.param('enableAllMoments', allMoments)
3821 eventShapeModule.param('enableCleoCones', cleoCones)
3822 eventShapeModule.param('enableCollisionAxis', collisionAxis)
3823 eventShapeModule.param('enableFoxWolfram', foxWolfram)
3824 eventShapeModule.param('enableJets', jets)
3825 eventShapeModule.param('enableHarmonicMoments', harmonicMoments)
3826 eventShapeModule.param('enableSphericity', sphericity)
3827 eventShapeModule.param('enableThrust', thrust)
3828 eventShapeModule.param('checkForDuplicates', checkForDuplicates)
3829
3830 path.add_module(eventShapeModule)
3831
3832
3833def labelTauPairMC(printDecayInfo=False, path=None, TauolaBelle=False, mapping_minus=None, mapping_plus=None):
3834 """
3835 Search tau leptons into the MC information of the event. If confirms it's a generated tau pair decay,
3836 labels the decay generated of the positive and negative leptons using the ID of KKMC tau decay table.
3837
3838 @param printDecayInfo: If true, prints ID and prong of each tau lepton in the event.
3839 @param path: module is added to this path
3840 @param TauolaBelle: if False, TauDecayMode is set. If True, TauDecayMarker is set.
3841 @param mapping_minus: if None, the map is the default one, else the path for the map is given by the user for tau-
3842 @param mapping_plus: if None, the map is the default one, else the path for the map is given by the user for tau+
3843 """
3844
3845 from basf2 import find_file
3846 if not TauolaBelle:
3847
3848 if printDecayInfo:
3849 m_printmode = 'all'
3850 else:
3851 m_printmode = 'default'
3852
3853 if mapping_minus is None:
3854 mp_file_minus = find_file('data/analysis/modules/TauDecayMode/map_tauminus.txt')
3855 else:
3856 mp_file_minus = mapping_minus
3857
3858 if mapping_plus is None:
3859 mp_file_plus = find_file('data/analysis/modules/TauDecayMode/map_tauplus.txt')
3860 else:
3861 mp_file_plus = mapping_plus
3862
3863 path.add_module('TauDecayMode', printmode=m_printmode, file_minus=mp_file_minus, file_plus=mp_file_plus)
3864
3865 else:
3866 tauDecayMarker = register_module('TauDecayMarker')
3867 tauDecayMarker.set_name('TauDecayMarker_')
3868
3869 path.add_module(tauDecayMarker, printDecayInfo=printDecayInfo)
3870
3871
3872def tagCurlTracks(particleLists,
3873 mcTruth=False,
3874 responseCut=-1.0,
3875 selectorType='cut',
3876 ptCut=0.5,
3877 expert_train=False,
3878 expert_filename="",
3879 path=None):
3880 """
3881 Warning:
3882 The cut selector is not calibrated with Belle II data and should not be used without extensive study.
3883
3884 Identifies curl tracks and tags them with extraInfo(isCurl=1) for later removal.
3885 For Belle data with a `b2bii` analysis the available cut based selection is described in `BN1079`_.
3886
3887 .. _BN1079: https://belle.kek.jp/secured/belle_note/gn1079/bn1079.pdf
3888
3889
3890 The module loops over all particles in a given list with a transverse momentum below the pre-selection **ptCut**
3891 and assigns them to bundles based on the response of the chosen **selector** and the required minimum response set by the
3892 **responseCut**. Once all particles are assigned they are ranked by 25dr^2+dz^2. All but the lowest are tagged
3893 with extraInfo(isCurl=1) to allow for later removal by cutting the list or removing these from ROE as
3894 applicable.
3895
3896
3897 @param particleLists: list of particle lists to check for curls.
3898 @param mcTruth: bool flag to additionally assign particles with extraInfo(isTruthCurl) and
3899 extraInfo(truthBundleSize). To calculate these particles are assigned to bundles by their
3900 genParticleIndex then ranked and tagged as normal.
3901 @param responseCut: float min classifier response that considers two tracks to come from the same particle.
3902 If set to ``-1`` a cut value optimised to maximise the accuracy on a BBbar sample is used.
3903 Note 'cut' selector is binary 0/1.
3904 @param selectorType: string name of selector to use. The available options are 'cut' and 'mva'.
3905 It is strongly recommended to used the 'mva' selection. The 'cut' selection
3906 is based on BN1079 and is only calibrated for Belle data.
3907
3908 @param ptCut: Pre-selection cut on transverse momentum. Only tracks below that are considered as curler candidates.
3909
3910 @param expert_train: flag to set training mode if selector has a training mode (mva).
3911 @param expert_filename: set file name of produced training ntuple (mva).
3912 @param path: module is added to this path.
3913 """
3914
3915 import b2bii
3916 belle = b2bii.isB2BII()
3917
3918 if (not isinstance(particleLists, list)):
3919 particleLists = [particleLists] # in case user inputs a particle list as string
3920
3921 curlTagger = register_module('CurlTagger')
3922 curlTagger.set_name('CurlTagger_')
3923 curlTagger.param('particleLists', particleLists)
3924 curlTagger.param('belle', belle)
3925 curlTagger.param('mcTruth', mcTruth)
3926 curlTagger.param('responseCut', responseCut)
3927 if abs(responseCut + 1) < 1e-9:
3928 curlTagger.param('usePayloadCut', True)
3929 else:
3930 curlTagger.param('usePayloadCut', False)
3931
3932 curlTagger.param('selectorType', selectorType)
3933 curlTagger.param('ptCut', ptCut)
3934 curlTagger.param('train', expert_train)
3935 curlTagger.param('trainFilename', expert_filename)
3936
3937 path.add_module(curlTagger)
3938
3939
3940def applyChargedPidMVA(particleLists, path, trainingMode, chargeIndependent=False, binaryHypoPDGCodes=(0, 0)):
3941 """
3942 Use an MVA to perform particle identification for charged stable particles, using the `ChargedPidMVA` module.
3943
3944 The module decorates Particle objects in the input ParticleList(s) with variables
3945 containing the appropriate MVA score, which can be used to select candidates by placing a cut on it.
3946
3947 Note:
3948 The MVA algorithm used is a gradient boosted decision tree (**TMVA 4.3.0**, **ROOT 6.20/04**).
3949
3950 The module can perform either 'binary' PID between input S, B particle mass hypotheses according to the following scheme:
3951
3952 * e (11) vs. pi (211)
3953 * mu (13) vs. pi (211)
3954 * pi (211) vs. K (321)
3955 * K (321) vs. pi (211)
3956
3957 , or 'global' PID, namely "one-vs-others" separation. The latter exploits an MVA algorithm trained in multi-class mode,
3958 and it's the default behaviour. Currently, the multi-class training separates the following standard charged hypotheses:
3959
3960 - e (11), mu (13), pi (211), K (321)
3961
3962 Warning:
3963 In order to run the `ChargedPidMVA` and ensure the most up-to-date MVA training weights are applied,
3964 it is necessary to append the latest analysis global tag (GT) to the steering script.
3965
3966 Parameters:
3967 particleLists (list(str)): the input list of DecayStrings, where each selected (^) daughter should correspond to a
3968 standard charged ParticleList, e.g. ``['Lambda0:sig -> ^p+ ^pi-', 'J/psi:sig -> ^mu+ ^mu-']``.
3969 One can also directly pass a list of standard charged ParticleLists,
3970 e.g. ``['e+:my_electrons', 'pi+:my_pions']``.
3971 Note that charge-conjugated ParticleLists will automatically be included.
3972 path (basf2.Path): the module is added to this path.
3973 trainingMode (``Belle2.ChargedPidMVAWeights.ChargedPidMVATrainingMode``): enum identifier of the training mode.
3974 Needed to pick up the correct payload from the DB. Available choices:
3975
3976 * c_Classification=0
3977 * c_Multiclass=1
3978 * c_ECL_Classification=2
3979 * c_ECL_Multiclass=3
3980 * c_PSD_Classification=4
3981 * c_PSD_Multiclass=5
3982 * c_ECL_PSD_Classification=6
3983 * c_ECL_PSD_Multiclass=7
3984
3985 chargeIndependent (bool, ``optional``): use a BDT trained on a sample of inclusively charged particles.
3986 binaryHypoPDGCodes (tuple(int, int), ``optional``): the pdgIds of the signal, background mass hypothesis.
3987 Required only for binary PID mode.
3988 """
3989
3990 import b2bii
3991 if b2bii.isB2BII():
3992 B2ERROR("Charged PID via MVA is not available for Belle data.")
3993
3994 from ROOT import Belle2
3995
3996 TrainingMode = Belle2.ChargedPidMVAWeights.ChargedPidMVATrainingMode
3997 Const = Belle2.Const
3998
3999 plSet = set(particleLists)
4000
4001 # Map the training mode enum value to the actual name of the payload in the GT.
4002 payloadNames = {
4003 TrainingMode.c_Classification:
4004 {"mode": "Classification", "detector": "ALL"},
4005 TrainingMode.c_Multiclass:
4006 {"mode": "Multiclass", "detector": "ALL"},
4007 TrainingMode.c_ECL_Classification:
4008 {"mode": "ECL_Classification", "detector": "ECL"},
4009 TrainingMode.c_ECL_Multiclass:
4010 {"mode": "ECL_Multiclass", "detector": "ECL"},
4011 TrainingMode.c_PSD_Classification:
4012 {"mode": "PSD_Classification", "detector": "ALL"},
4013 TrainingMode.c_PSD_Multiclass:
4014 {"mode": "PSD_Multiclass", "detector": "ALL"},
4015 TrainingMode.c_ECL_PSD_Classification:
4016 {"mode": "ECL_PSD_Classification", "detector": "ECL"},
4017 TrainingMode.c_ECL_PSD_Multiclass:
4018 {"mode": "ECL_PSD_Multiclass", "detector": "ECL"},
4019 }
4020
4021 if payloadNames.get(trainingMode) is None:
4022 B2FATAL("The chosen training mode integer identifier:\n", trainingMode,
4023 "\nis not supported. Please choose among the following:\n",
4024 "\n".join(f"{key}:{val.get('mode')}" for key, val in sorted(payloadNames.items())))
4025
4026 mode = payloadNames.get(trainingMode).get("mode")
4027 detector = payloadNames.get(trainingMode).get("detector")
4028
4029 payloadName = f"ChargedPidMVAWeights_{mode}"
4030
4031 # Map pdgIds of std charged particles to name identifiers,
4032 # and binary bkg identifiers.
4033 stdChargedMap = {
4034 Const.electron.getPDGCode():
4035 {"pName": "e", "pFullName": "electron", "pNameBkg": "pi", "pdgIdBkg": Const.pion.getPDGCode()},
4036 Const.muon.getPDGCode():
4037 {"pName": "mu", "pFullName": "muon", "pNameBkg": "pi", "pdgIdBkg": Const.pion.getPDGCode()},
4038 Const.pion.getPDGCode():
4039 {"pName": "pi", "pFullName": "pion", "pNameBkg": "K", "pdgIdBkg": Const.kaon.getPDGCode()},
4040 Const.kaon.getPDGCode():
4041 {"pName": "K", "pFullName": "kaon", "pNameBkg": "pi", "pdgIdBkg": Const.pion.getPDGCode()},
4042 Const.proton.getPDGCode():
4043 {"pName": "p", "pFullName": "proton", "pNameBkg": "pi", "pdgIdBkg": Const.pion.getPDGCode()},
4044 Const.deuteron.getPDGCode():
4045 {"pName": "d", "pFullName": "deuteron", "pNameBkg": "pi", "pdgIdBkg": Const.pion.getPDGCode()},
4046 }
4047
4048 if binaryHypoPDGCodes == (0, 0):
4049
4050 # MULTI-CLASS training mode.
4051 chargedpid = register_module("ChargedPidMVAMulticlass")
4052 chargedpid.set_name(f"ChargedPidMVAMulticlass_{mode}")
4053
4054 else:
4055
4056 # BINARY training mode.
4057 # In binary mode, enforce check on input S, B hypotheses compatibility.
4058
4059 binaryOpts = [(pdgIdSig, info["pdgIdBkg"]) for pdgIdSig, info in stdChargedMap.items()]
4060
4061 if binaryHypoPDGCodes not in binaryOpts:
4062 B2FATAL("No charged pid MVA was trained to separate ", binaryHypoPDGCodes[0], " vs. ", binaryHypoPDGCodes[1],
4063 ". Please choose among the following pairs:\n",
4064 "\n".join(f"{opt[0]} vs. {opt[1]}" for opt in binaryOpts))
4065
4066 decayDescriptor = Belle2.DecayDescriptor()
4067 for name in plSet:
4068 if not decayDescriptor.init(name):
4069 raise ValueError(f"Invalid particle list {name} in applyChargedPidMVA!")
4070 msg = f"Input ParticleList: {name}"
4071 pdgs = [abs(decayDescriptor.getMother().getPDGCode())]
4072 daughter_pdgs = decayDescriptor.getSelectionPDGCodes()
4073 if len(daughter_pdgs) > 0:
4074 pdgs = daughter_pdgs
4075 for idaughter, pdg in enumerate(pdgs):
4076 if abs(pdg) not in binaryHypoPDGCodes:
4077 if daughter_pdgs:
4078 msg = f"Selected daughter {idaughter} in ParticleList: {name}"
4079 B2WARNING(
4080 f"{msg} (PDG={pdg}) is neither signal ({binaryHypoPDGCodes[0]}) nor background ({binaryHypoPDGCodes[1]}).")
4081
4082 chargedpid = register_module("ChargedPidMVA")
4083 chargedpid.set_name(f"ChargedPidMVA_{binaryHypoPDGCodes[0]}_vs_{binaryHypoPDGCodes[1]}_{mode}")
4084 chargedpid.param("sigHypoPDGCode", binaryHypoPDGCodes[0])
4085 chargedpid.param("bkgHypoPDGCode", binaryHypoPDGCodes[1])
4086
4087 chargedpid.param("particleLists", list(plSet))
4088 chargedpid.param("payloadName", payloadName)
4089 chargedpid.param("chargeIndependent", chargeIndependent)
4090
4091 # Ensure the module knows whether we are using ECL-only training mode.
4092 if detector == "ECL":
4093 chargedpid.param("useECLOnlyTraining", True)
4094
4095 path.add_module(chargedpid)
4096
4097
4098def calculateTrackIsolation(
4099 decay_string,
4100 path,
4101 *detectors,
4102 reference_list_name=None,
4103 vars_for_nearest_part=[],
4104 highest_prob_mass_for_ext=True,
4105 exclude_pid_det_weights=False):
4106 """
4107 Given an input decay string, compute variables that quantify track helix-based isolation of the charged
4108 stable particles in the input decay chain.
4109
4110 Note:
4111 An "isolation score" can be defined using the distance
4112 of each particle to its closest neighbour, defined as the segment connecting the two
4113 extrapolated track helices intersection points on a given cylindrical surface.
4114 The distance variables defined in the `VariableManager` is named `minET2ETDist`,
4115 the isolation scores are named `minET2ETIsoScore`, `minET2ETIsoScoreAsWeightedAvg`.
4116
4117 The definition of distance and the number of distances that are calculated per sub-detector is based on
4118 the following recipe:
4119
4120 * **CDC**: as the segmentation is very coarse along :math:`z`,
4121 the distance is defined as the cord length on the :math:`(\\rho=R, \\phi)` plane.
4122 A total of 9 distances are calculated: the cylindrical surfaces are defined at radiuses
4123 that correspond to the positions of the 9 CDC wire superlayers: :math:`R_{i}^{\\mathrm{CDC}}~(i \\in \\{0,...,8\\})`.
4124
4125 * **TOP**: as there is no segmentation along :math:`z`,
4126 the distance is defined as the cord length on the :math:`(\\rho=R, \\phi)` plane.
4127 Only one distance at the TOP entry radius :math:`R_{0}^{\\mathrm{TOP}}` is calculated.
4128
4129 * **ARICH**: as there is no segmentation along :math:`z`,
4130 the distance is defined as the distance on the :math:`(\\rho=R, \\phi)` plane at fixed :math:`z=Z`.
4131 Only one distance at the ARICH photon detector entry coordinate :math:`Z_{0}^{\\mathrm{ARICH}}` is calculated.
4132
4133 * **ECL**: the distance is defined on the :math:`(\\rho=R, \\phi, z)` surface in the barrel,
4134 on the :math:`(\\rho, \\phi, z=Z)` surface in the endcaps.
4135 Two distances are calculated: one at the ECL entry surface :math:`R_{0}^{\\mathrm{ECL}}` (barrel),
4136 :math:`Z_{0}^{\\mathrm{ECL}}` (endcaps), and one at :math:`R_{1}^{\\mathrm{ECL}}` (barrel),
4137 :math:`Z_{1}^{\\mathrm{ECL}}` (endcaps), corresponding roughly to the mid-point
4138 of the longitudinal size of the crystals.
4139
4140 * **KLM**: the distance is defined on the :math:`(\\rho=R, \\phi, z)` surface in the barrel,
4141 on the :math:`(\\rho, \\phi, z=Z)` surface in the endcaps.
4142 Only one distance at the KLM first strip entry surface :math:`R_{0}^{\\mathrm{KLM}}` (barrel),
4143 :math:`Z_{0}^{\\mathrm{KLM}}` (endcaps) is calculated.
4144
4145 Parameters:
4146 decay_string (str): name of the input decay string with selected charged stable daughters,
4147 for example: ``Lambda0:merged -> ^p+ ^pi-``.
4148 Alternatively, it can be a particle list for charged stable particles
4149 as defined in ``Const::chargedStableSet``, for example: ``mu+:all``.
4150 The charge-conjugate particle list will be also processed automatically.
4151 path (basf2.Path): path to which module(s) will be added.
4152 *detectors: detectors for which track isolation variables will be calculated.
4153 Choose among: ``{'CDC', 'TOP', 'ARICH', 'ECL', 'KLM'}``.
4154 reference_list_name (Optional[str]): name of the input charged stable particle list for the reference tracks.
4155 By default, the ``:all`` ParticleList of the same type
4156 of the selected particle in ``decay_string`` is used.
4157 The charge-conjugate particle list will be also processed automatically.
4158 vars_for_nearest_part (Optional[list(str)]): a list of variables to calculate for the nearest particle in the reference
4159 list at each detector surface. It uses the metavariable `minET2ETDistVar`.
4160 If unset, only the distances to the nearest neighbour
4161 per detector are calculated.
4162 highest_prob_mass_for_hex (Optional[bool]): if this option is set to True (default), the helix extrapolation
4163 for the particles will use the track fit result for the most
4164 probable mass hypothesis, namely, the one that gives the highest
4165 chi2Prob of the fit. Otherwise, it uses the mass hypothesis that
4166 corresponds to the particle lists PDG.
4167 exclude_pid_det_weights (Optional[bool]): if this option is set to False (default), the isolation score
4168 calculation will take into account the weight that each detector has on the PID
4169 for the particle species of interest.
4170
4171 Returns:
4172 dict(int, list(str)): a dictionary mapping the PDG of each reference particle list to its isolation variables.
4173
4174 """
4175
4176 import pdg
4177 from ROOT import Belle2, TDatabasePDG
4178
4179 decayDescriptor = Belle2.DecayDescriptor()
4180 if not decayDescriptor.init(decay_string):
4181 B2FATAL(f"Invalid particle list {decay_string} in calculateTrackIsolation!")
4182 no_reference_list_name = not reference_list_name
4183
4184 det_and_layers = {
4185 "CDC": list(range(9)),
4186 "TOP": [0],
4187 "ARICH": [0],
4188 "ECL": [0, 1],
4189 "KLM": [0],
4190 }
4191 if any(d not in det_and_layers for d in detectors):
4192 B2FATAL(
4193 "Your input detector list: ",
4194 detectors,
4195 " contains an invalid choice. Please select among: ",
4196 list(
4197 det_and_layers.keys()))
4198
4199 # The module allows only one daughter to be selected at a time,
4200 # that's why here we preprocess the input decay string.
4201 select_symbol = '^'
4202 processed_decay_strings = []
4203 if select_symbol in decay_string:
4204 splitted_ds = decay_string.split(select_symbol)
4205 for i in range(decay_string.count(select_symbol)):
4206 tmp = list(splitted_ds)
4207 tmp.insert(i+1, select_symbol)
4208 processed_decay_strings += [''.join(tmp)]
4209 else:
4210 processed_decay_strings += [decay_string]
4211
4212 reference_lists_to_vars = {}
4213
4214 for processed_dec in processed_decay_strings:
4215 if no_reference_list_name:
4216 decayDescriptor.init(processed_dec)
4217 selected_daughter_pdgs = decayDescriptor.getSelectionPDGCodes()
4218 if len(selected_daughter_pdgs) > 0:
4219 reference_list_name = f'{TDatabasePDG.Instance().GetParticle(abs(selected_daughter_pdgs[-1])).GetName()}:all'
4220 else:
4221 reference_list_name = f'{processed_dec.split(":")[0]}:all'
4222
4223 ref_pdg = pdg.from_name(reference_list_name.split(":")[0])
4224
4225 trackiso = path.add_module("TrackIsoCalculator",
4226 decayString=processed_dec,
4227 detectorNames=list(detectors),
4228 particleListReference=reference_list_name,
4229 useHighestProbMassForExt=highest_prob_mass_for_ext,
4230 excludePIDDetWeights=exclude_pid_det_weights)
4231 trackiso.set_name(f"TrackIsoCalculator_{'_'.join(detectors)}_{processed_dec}_VS_{reference_list_name}")
4232
4233 # Metavariables for the distances to the closest reference tracks at each detector surface.
4234 # Always calculate them.
4235 # Ensure the flag for the mass hypothesis of the fit is set.
4236 trackiso_vars = [
4237 f"minET2ETDist({d}, {d_layer}, {reference_list_name}, {int(highest_prob_mass_for_ext)})"
4238 for d in detectors for d_layer in det_and_layers[d]]
4239 # Track isolation score.
4240 trackiso_vars += [
4241 f"minET2ETIsoScore({reference_list_name}, {int(highest_prob_mass_for_ext)}, {', '.join(detectors)})",
4242 f"minET2ETIsoScoreAsWeightedAvg({reference_list_name}, {int(highest_prob_mass_for_ext)}, {', '.join(detectors)})",
4243 ]
4244 # Optionally, calculate the input variables for the nearest neighbour in the reference list.
4245 if vars_for_nearest_part:
4246 trackiso_vars.extend(
4247 [
4248 f"minET2ETDistVar({d}, {d_layer}, {reference_list_name}, {v})"
4249 for d in detectors for d_layer in det_and_layers[d] for v in vars_for_nearest_part
4250 ])
4251 trackiso_vars.sort()
4252
4253 reference_lists_to_vars[ref_pdg] = trackiso_vars
4254
4255 return reference_lists_to_vars
4256
4257
4258def calculateDistance(list_name, decay_string, mode='vertextrack', path=None):
4259 """
4260 Calculates distance between two vertices, distance of closest approach between a vertex and a track,\
4261 distance of closest approach between a vertex and btube. For track, this calculation ignores track curvature,\
4262 it's negligible for small distances.The user should use extraInfo(CalculatedDistance)\
4263 to get it. A full example steering file is at analysis/tests/test_DistanceCalculator.py
4264
4265 Example:
4266 .. code-block:: python
4267
4268 from modularAnalysis import calculateDistance
4269 calculateDistance('list_name', 'decay_string', "mode", path=user_path)
4270
4271 @param list_name name of the input ParticleList
4272 @param decay_string select particles between the distance of closest approach will be calculated
4273 @param mode Specifies how the distance is calculated
4274 vertextrack: calculate the distance of closest approach between a track and a\
4275 vertex, taking the first candidate as vertex, default
4276 trackvertex: calculate the distance of closest approach between a track and a\
4277 vertex, taking the first candidate as track
4278 2tracks: calculates the distance of closest approach between two tracks
4279 2vertices: calculates the distance between two vertices
4280 vertexbtube: calculates the distance of closest approach between a vertex and btube
4281 trackbtube: calculates the distance of closest approach between a track and btube
4282 @param path modules are added to this path
4283
4284 """
4285
4286 dist_mod = register_module('DistanceCalculator')
4287
4288 dist_mod.set_name('DistanceCalculator_' + list_name)
4289 dist_mod.param('listName', list_name)
4290 dist_mod.param('decayString', decay_string)
4291 dist_mod.param('mode', mode)
4292 path.add_module(dist_mod)
4293
4294
4295def addInclusiveDstarReconstruction(decayString, slowPionCut, DstarCut, path):
4296 """
4297 Adds the InclusiveDstarReconstruction module to the given path.
4298 This module creates a D* particle list by estimating the D* four momenta
4299 from slow pions, specified by a given cut. The D* energy is approximated
4300 as E(D*) = m(D*)/(m(D*) - m(D)) * E(pi). The absolute value of the D*
4301 momentum is calculated using the D* PDG mass and the direction is collinear
4302 to the slow pion direction. The charge of the given pion list has to be consistent
4303 with the D* charge
4304
4305 @param decayString Decay string, must be of form ``D* -> pi``
4306 @param slowPionCut Cut applied to the input pion list to identify slow pions
4307 @param DstarCut Cut applied to the output D* list
4308 @param path the module is added to this path
4309 """
4310
4311 incl_dstar = register_module("InclusiveDstarReconstruction")
4312 incl_dstar.param("decayString", decayString)
4313 incl_dstar.param("slowPionCut", slowPionCut)
4314 incl_dstar.param("DstarCut", DstarCut)
4315 path.add_module(incl_dstar)
4316
4317
4318def scaleError(outputListName, inputListName,
4319 scaleFactors=[1.149631, 1.085547, 1.151704, 1.096434, 1.086659],
4320 scaleFactorsNoPXD=[1.149631, 1.085547, 1.151704, 1.096434, 1.086659],
4321 d0Resolution=[0.00115328, 0.00134704],
4322 z0Resolution=[0.00124327, 0.0013272],
4323 d0MomThr=0.500000,
4324 z0MomThr=0.500000,
4325 path=None):
4326 """
4327 This module creates a new charged particle list.
4328 The helix errors of the new particles are scaled by constant factors.
4329 Two sets of five scale factors are defined for tracks with and without a PXD hit.
4330 The scale factors are in order of (d0, phi0, omega, z0, tanlambda).
4331 For tracks with a PXD hit, in order to avoid severe underestimation of d0 and z0 errors,
4332 lower limits (best resolution) can be set in a momentum-dependent form.
4333 This module is supposed to be used only for TDCPV analysis and for low-momentum (0-3 GeV/c) tracks in BBbar events.
4334 Details will be documented in a Belle II note, BELLE2-NOTE-PH-2021-038.
4335
4336 @param inputListName Name of input charged particle list to be scaled
4337 @param outputListName Name of output charged particle list with scaled error
4338 @param scaleFactors List of five constants to be multiplied to each of helix errors (for tracks with a PXD hit)
4339 @param scaleFactorsNoPXD List of five constants to be multiplied to each of helix errors (for tracks without a PXD hit)
4340 @param d0Resolution List of two parameters, (a [cm], b [cm/(GeV/c)]),
4341 defining d0 best resolution as sqrt{ a**2 + (b / (p*beta*sinTheta**1.5))**2 }
4342 @param z0Resolution List of two parameters, (a [cm], b [cm/(GeV/c)]),
4343 defining z0 best resolution as sqrt{ a**2 + (b / (p*beta*sinTheta**2.5))**2 }
4344 @param d0MomThr d0 best resolution is kept constant below this momentum
4345 @param z0MomThr z0 best resolution is kept constant below this momentum
4346
4347 """
4348
4349 scale_error = register_module("HelixErrorScaler")
4350 scale_error.set_name('ScaleError_' + inputListName)
4351 scale_error.param('inputListName', inputListName)
4352 scale_error.param('outputListName', outputListName)
4353 scale_error.param('scaleFactors_PXD', scaleFactors)
4354 scale_error.param('scaleFactors_noPXD', scaleFactorsNoPXD)
4355 scale_error.param('d0ResolutionParameters', d0Resolution)
4356 scale_error.param('z0ResolutionParameters', z0Resolution)
4357 scale_error.param('d0MomentumThreshold', d0MomThr)
4358 scale_error.param('z0MomentumThreshold', z0MomThr)
4359 path.add_module(scale_error)
4360
4361
4362def estimateAndAttachTrackFitResult(inputListName, path=None):
4363 """
4364 Create a TrackFitResult from the momentum of the Particle assuming it originates from the IP and make a relation between them.
4365 The covariance, detector hit information, and fit-related information (pValue, NDF) are assigned meaningless values. The input
4366 Particles must not have already Track or TrackFitResult and thus are supposed to be composite particles, recoil, dummy
4367 particles, and so on.
4368
4369
4370 .. warning:: Since the source type is not overwritten as Track, not all track-related variables are guaranteed to be available.
4371
4372
4373 @param inputListName Name of input ParticleList
4374 """
4375
4376 estimator = register_module("TrackFitResultEstimator")
4377 estimator.set_name("trackFitResultEstimator_" + inputListName)
4378 estimator.param("inputListName", inputListName)
4379 path.add_module(estimator)
4380
4381
4382def correctEnergyBias(inputListNames, tableName, path=None):
4383 """
4384 Scale energy of the particles according to the scaling factor.
4385 If the particle list contains composite particles, the energy of the daughters are scaled.
4386 Subsequently, the energy of the mother particle is updated as well.
4387
4388 Parameters:
4389 inputListNames (list(str)): input particle list names
4390 tableName : stored in localdb and created using ParticleWeightingLookUpCreator
4391 path (basf2.Path): module is added to this path
4392 """
4393
4394 import b2bii
4395 if b2bii.isB2BII():
4396 B2ERROR("The energy bias cannot be corrected with this tool for Belle data.")
4397
4398 correctenergybias = register_module('EnergyBiasCorrection')
4399 correctenergybias.param('particleLists', inputListNames)
4400 correctenergybias.param('tableName', tableName)
4401 path.add_module(correctenergybias)
4402
4403
4404def twoBodyISRPhotonCorrector(outputListName, inputListName, massiveParticle, path=None):
4405 """
4406 Sets photon kinematics to corrected values in two body decays with an ISR photon
4407 and a massive particle. The original photon kinematics are kept in the input
4408 particleList and can be accessed using the originalParticle() metavariable on the
4409 new list.
4410
4411 @param ouputListName new ParticleList filled with copied Particles
4412 @param inputListName input ParticleList with original Particles
4413 @param massiveParticle name or PDG code of massive particle participating in the two
4414 body decay with the ISR photon
4415 @param path modules are added to this path
4416 """
4417
4418 # set the corrected energy of the photon in a new list
4419 photon_energy_correction = register_module('TwoBodyISRPhotonCorrector')
4420 photon_energy_correction.set_name('TwoBodyISRPhotonCorrector_' + outputListName)
4421 photon_energy_correction.param('outputGammaList', outputListName)
4422 photon_energy_correction.param('inputGammaList', inputListName)
4423
4424 # prepare PDG code of massive particle
4425 if isinstance(massiveParticle, int):
4426 photon_energy_correction.param('massiveParticlePDGCode', massiveParticle)
4427 else:
4428 from ROOT import Belle2
4429 decayDescriptor = Belle2.DecayDescriptor()
4430 if not decayDescriptor.init(massiveParticle):
4431 raise ValueError("TwoBodyISRPhotonCorrector: value of massiveParticle must be" +
4432 " an int or valid decay string.")
4433 pdgCode = decayDescriptor.getMother().getPDGCode()
4434 photon_energy_correction.param('massiveParticlePDGCode', pdgCode)
4435
4436 path.add_module(photon_energy_correction)
4437
4438
4439def addPhotonEfficiencyRatioVariables(inputListNames, tableName, path=None):
4440 """
4441 Add photon Data/MC detection efficiency ratio weights to the specified particle list
4442
4443 Parameters:
4444 inputListNames (list(str)): input particle list names
4445 tableName : taken from database with appropriate name
4446 path (basf2.Path): module is added to this path
4447 """
4448
4449 import b2bii
4450 if b2bii.isB2BII():
4451 B2ERROR("For Belle data the photon data/MC detection efficiency ratio is not available with this tool.")
4452
4453 photon_efficiency_correction = register_module('PhotonEfficiencySystematics')
4454 photon_efficiency_correction.param('particleLists', inputListNames)
4455 photon_efficiency_correction.param('tableName', tableName)
4456 path.add_module(photon_efficiency_correction)
4457
4458
4459def addPi0VetoEfficiencySystematics(particleList, decayString, tableName, threshold, mode='standard', suffix='', path=None):
4460 """
4461 Add pi0 veto Data/MC efficiency ratio weights to the specified particle list
4462
4463 @param particleList the input ParticleList
4464 @param decayString specify hard photon to be performed pi0 veto (e.g. 'B+:sig -> rho+:sig ^gamma:hard')
4465 @param tableName table name corresponding to payload version (e.g. 'Pi0VetoEfficiencySystematics_Mar2022')
4466 @param threshold pi0 veto threshold (0.10, 0.11, ..., 0.99)
4467 @param mode choose one mode (same as writePi0EtaVeto) out of 'standard', 'tight', 'cluster' and 'both'
4468 @param suffix optional suffix to be appended to the usual extraInfo name
4469 @param path the module is added to this path
4470
4471 The following extraInfo are available related with the given particleList:
4472
4473 * Pi0VetoEfficiencySystematics_{mode}{suffix}_data_MC_ratio : weight of Data/MC for the veto efficiency
4474 * Pi0VetoEfficiencySystematics_{mode}{suffix}_data_MC_uncertainty_stat : the statistical uncertainty of the weight
4475 * Pi0VetoEfficiencySystematics_{mode}{suffix}_data_MC_uncertainty_sys : the systematic uncertainty of the weight
4476 * Pi0VetoEfficiencySystematics_{mode}{suffix}_data_MC_uncertainty_total : the total uncertainty of the weight
4477 * Pi0VetoEfficiencySystematics_{mode}{suffix}_threshold : threshold of the pi0 veto
4478 """
4479
4480 import b2bii
4481 if b2bii.isB2BII():
4482 B2ERROR("For Belle data the pi0 veto data/MC efficiency ratio weights are not available via this tool.")
4483
4484 pi0veto_efficiency_correction = register_module('Pi0VetoEfficiencySystematics')
4485 pi0veto_efficiency_correction.param('particleLists', particleList)
4486 pi0veto_efficiency_correction.param('decayString', decayString)
4487 pi0veto_efficiency_correction.param('tableName', tableName)
4488 pi0veto_efficiency_correction.param('threshold', threshold)
4489 pi0veto_efficiency_correction.param('mode', mode)
4490 pi0veto_efficiency_correction.param('suffix', suffix)
4491 path.add_module(pi0veto_efficiency_correction)
4492
4493
4494def getAnalysisGlobaltag(timeout=180) -> str:
4495 """
4496 Returns a string containing the name of the latest and recommended analysis globaltag.
4497
4498 Parameters:
4499 timeout: Seconds to wait for b2conditionsdb-recommend
4500 """
4501
4502 import b2bii
4503 if b2bii.isB2BII():
4504 B2ERROR("The getAnalysisGlobaltag() function cannot be used for Belle data.")
4505
4506 # b2conditionsdb-recommend relies on a different repository, so it's better to protect
4507 # this function against potential failures of check_output.
4508 try:
4509 tags = subprocess.check_output(
4510 ['b2conditionsdb-recommend', '--oneline'],
4511 timeout=timeout
4512 ).decode('UTF-8').rstrip().split(' ')
4513 analysis_tag = ''
4514 for tag in tags:
4515 if tag.startswith('analysis_tools'):
4516 analysis_tag = tag
4517 return analysis_tag
4518 # In case of issues with git, b2conditionsdb-recommend may take too much time.
4519 except subprocess.TimeoutExpired as te:
4520 B2FATAL(f'A {te} exception was raised during the call of getAnalysisGlobaltag(). '
4521 'The function took too much time to retrieve the requested information '
4522 'from the versioning repository.\n'
4523 'Please try to re-run your job. In case of persistent failures, there may '
4524 'be issues with the DESY collaborative services, so please contact the experts.')
4525 except subprocess.CalledProcessError as ce:
4526 B2FATAL(f'A {ce} exception was raised during the call of getAnalysisGlobaltag(). '
4527 'Please try to re-run your job. In case of persistent failures, please contact '
4528 'the experts.')
4529
4530
4531def getAnalysisGlobaltagB2BII() -> str:
4532 """
4533 Get recommended global tag for B2BII analysis.
4534 """
4535
4536 import b2bii
4537 if not b2bii.isB2BII():
4538 B2ERROR('The getAnalysisGlobaltagB2BII() function cannot be used for Belle II data.')
4539 from versioning import recommended_b2bii_analysis_global_tag
4540 return recommended_b2bii_analysis_global_tag()
4541
4542
4543def getECLKLID(particleList: str, variable='ECLKLID', path=None):
4544 """
4545 The function calculates the PID value for Klongs that are constructed from ECL cluster.
4546
4547 @param particleList the input ParticleList
4548 @param variable the variable name for Klong ID
4549 @param path modules are added to this path
4550 """
4551
4552 import b2bii
4553
4554 if b2bii.isB2BII():
4555 B2ERROR("The ECL variables based Klong Identification is only available for Belle II data.")
4556
4557 from variables import variables
4558 path.add_module('MVAExpert', listNames=particleList, extraInfoName='ECLKLID', identifier='ECLKLID')
4559
4560 variables.addAlias(variable, 'conditionalVariableSelector(isFromECL and PDG==130, extraInfo(ECLKLID), constant(NaN))')
4561
4562
4563def getNbarIDMVA(particleList: str, path=None):
4564 """
4565 This function can give a score to predict if it is a anti-n0.
4566 It is not used to predict n0.
4567 Currently, this can be used only for ECL cluster.
4568 output will be stored in extraInfo(nbarID); -1 means MVA invalid
4569
4570 @param particleList The input ParticleList name or a decay string which contains a full mother particle list name.
4571 Only one selected daughter is supported.
4572 @param path modules are added to this path
4573 """
4574 import b2bii
4575 from ROOT import Belle2
4576
4577 if b2bii.isB2BII():
4578 B2ERROR("The MVA-based anti-neutron PID is only available for Belle II data.")
4579
4580 from variables import variables
4581
4582 variables.addAlias('V1', 'clusterHasPulseShapeDiscrimination')
4583 variables.addAlias('V2', 'clusterE')
4584 variables.addAlias('V3', 'clusterLAT')
4585 variables.addAlias('V4', 'clusterE1E9')
4586 variables.addAlias('V5', 'clusterE9E21')
4587 variables.addAlias('V6', 'clusterZernikeMVA')
4588 variables.addAlias('V7', 'clusterAbsZernikeMoment40')
4589 variables.addAlias('V8', 'clusterAbsZernikeMoment51')
4590
4591 variables.addAlias(
4592 'nbarIDValid',
4593 'passesCut(V1 == 1 and V2 >= 0 and V3 >= 0 and V4 >= 0 and V5 >= 0 and V6 >= 0 and V7 >= 0 and V8 >= 0)')
4594 variables.addAlias('nbarIDmod', 'conditionalVariableSelector(nbarIDValid == 1, extraInfo(nbarIDFromMVA), constant(-1.0))')
4595
4596 path.add_module('MVAExpert', listNames=particleList, extraInfoName='nbarIDFromMVA', identifier='db_nbarIDECL')
4597 decayDescriptor = Belle2.DecayDescriptor()
4598 if not decayDescriptor.init(particleList):
4599 raise ValueError(f"Provided decay string is invalid: {particleList}")
4600 if decayDescriptor.getNDaughters() == 0:
4601 variablesToExtraInfo(particleList, {'nbarIDmod': 'nbarID'}, option=2, path=path)
4602 else:
4603 listname = decayDescriptor.getMother().getFullName()
4604 variablesToDaughterExtraInfo(listname, particleList, {'nbarIDmod': 'nbarID'}, option=2, path=path)
4605
4606
4607def reconstructDecayWithNeutralHadron(decayString, cut, allowGamma=False, allowAnyParticleSource=False, path=None, **kwargs):
4608 r"""
4609 Reconstructs decay with a long-lived neutral hadron e.g.
4610 :math:`B^0 \to J/\psi K_L^0`,
4611 :math:`B^0 \to p \bar{n} D^*(2010)^-`.
4612
4613 The calculation is done with IP constraint and mother mass constraint.
4614
4615 The decay string passed in must satisfy the following rules:
4616
4617 - The neutral hadron must be **selected** in the decay string with the
4618 caret (``^``) e.g. ``B0:sig -> J/psi:sig ^K_L0:sig``. (Note the caret
4619 next to the neutral hadron.)
4620 - There can only be **one neutral hadron in a decay**.
4621 - The neutral hadron has to be a direct daughter of its mother.
4622
4623 .. note:: This function forwards its arguments to `reconstructDecay`,
4624 so please check the documentation of `reconstructDecay` for all
4625 possible arguments.
4626
4627 @param decayString A decay string following the mentioned rules
4628 @param cut Cut to apply to the particle list
4629 @param allowGamma Whether allow the selected particle to be ``gamma``
4630 @param allowAnyParticleSource Whether allow the selected particle to be from any source.
4631 Should only be used when studying control sample.
4632 @param path The path to put in the module
4633 """
4634
4635 reconstructDecay(decayString, cut, path=path, **kwargs)
4636 module = register_module('NeutralHadron4MomentumCalculator')
4637 module.set_name('NeutralHadron4MomentumCalculator_' + decayString)
4638 module.param('decayString', decayString)
4639 module.param('allowGamma', allowGamma)
4640 module.param('allowAnyParticleSource', allowAnyParticleSource)
4641 path.add_module(module)
4642
4643
4644def updateMassHypothesis(particleList, pdg, writeOut=False, path=None):
4645 """
4646 Module to update the mass hypothesis of a given input particle list with the chosen PDG.
4647 A new particle list is created with updated mass hypothesis.
4648 The allowed mass hypotheses for both input and output are electrons, muons, pions, kaons and protons.
4649
4650 .. note:
4651 The new particle list is named after the input one, with the additional suffix ``_converted_from_OLDHYPOTHESIS``,
4652 e.g. ``e+:all`` converted to muons becomes ``mu+:all_converted_from_e``.
4653
4654 @param particleList The input particle list name
4655 @param pdg The PDG code for the new mass hypothesis, in [11, 13, 211, 321, 2212]
4656 @param writeOut Whether `RootOutput` module should save the new particle list
4657 @param path Modules are added to this path
4658 """
4659 mass_updater = register_module("ParticleMassHypothesesUpdater")
4660 mass_updater.set_name("ParticleMassHypothesesUpdater_" + particleList + "_to_" + str(pdg))
4661 mass_updater.param("particleList", particleList)
4662 mass_updater.param("writeOut", writeOut)
4663 mass_updater.param("pdgCode", pdg)
4664 path.add_module(mass_updater)
4665
4666
4667func_requiring_analysisGT = [
4668 correctTrackEnergy, scaleTrackMomenta, smearTrackMomenta, oldwritePi0EtaVeto, writePi0EtaVeto, lowEnergyPi0Identification,
4669 getBeamBackgroundProbability, getFakePhotonProbability, tagCurlTracks, applyChargedPidMVA, correctEnergyBias,
4670 addPhotonEfficiencyRatioVariables, addPi0VetoEfficiencySystematics, getNbarIDMVA, getECLKLID]
4671for _ in func_requiring_analysisGT:
4672 _.__doc__ += "\n .. note:: This function (optionally) requires a payload stored in the analysis GlobalTag. "\
4673 "Please append or prepend the latest one from `getAnalysisGlobaltag` or `getAnalysisGlobaltagB2BII`.\n"
4674
4675
4676if __name__ == '__main__':
4677 from basf2.utils import pretty_print_module
4678 pretty_print_module(__name__, "modularAnalysis")
isB2BII()
Definition b2bii.py:14
setB2BII()
Definition b2bii.py:21
tuple parse(str cut, verbose=False)
Definition b2parser.py:975
This class provides a set of constants for the framework.
Definition Const.h:34
The DecayDescriptor stores information about a decay tree or parts of a decay tree.
Describe one component of the Geometry.
Magnetic field map.
static DBStore & Instance()
Instance of a singleton DBStore.
Definition DBStore.cc:26
add_mdst_output(path, mc=True, filename='mdst.root', additionalBranches=[], dataDescription=None)
Definition mdst.py:38
from_name(name)
Definition pdg.py:63
add_udst_output(path, filename, particleLists=None, additionalBranches=None, dataDescription=None, mc=True)
Definition udst.py:27