Belle II Software development
generators.py
1
8
9'''
10This module contains convenience functions to setup most commonly used physics
11generators correctly with their default settings. More information can be found
12in `BELLE2-NOTE-PH-2015-006`_
13
14.. _BELLE2-NOTE-PH-2015-006: https://docs.belle2.org/record/282
15'''
16
17import basf2 as b2
18import pdg
19
20
21def get_default_decayfile():
22 """Return the default DECAY.dec for Belle2"""
23 return b2.find_file("decfiles/dec/DECAY_BELLE2.DEC")
24
25
26def add_generator_preselection(
27 path,
28 emptypath,
29 nChargedMin=0,
30 nChargedMax=999,
31 MinChargedP=-1.0,
32 MinChargedPt=-1.0,
33 MinChargedTheta=0.0,
34 MaxChargedTheta=180.0,
35 nPhotonMin=0,
36 nPhotonMax=999,
37 MinPhotonEnergy=-1,
38 MinPhotonTheta=0.0,
39 MaxPhotonTheta=180.0,
40 applyInCMS=False,
41 stableParticles=False):
42 """
43 Adds generator preselection.
44 Should be added to the path after the generator.add_abc_generator but before simulation.add_simulation modules
45 It uses all particles from the event generator (i.e. primary, non-virtual, non-initial particles).
46 It checks if the required conditions are fulfilled.
47 If not, the events are given to the emptypath.
48 The main use case is a reduction of simulation time.
49 Note that you have to multiply the generated cross section by the retention fraction of the preselection.
50
51 Parameters:
52 path (basf2.Path): path where the generator should be added
53 emptypath (basf2.Path): path where the skipped events are given to
54 nChargedMin (int): minimum number of charged particles
55 nChargedMax (int): maximum number of charged particles
56 MinChargedP (float): minimum charged momentum [GeV]
57 MinChargedPt (float): minimum charged transverse momentum (pt) [GeV]
58 MinChargedTheta (float): minimum polar angle of charged particle [deg]
59 MaxChargedTheta (float): maximum polar angle of charged particle [deg]
60 nPhotonMin (int): minimum number of photons
61 nPhotonMax (int): maximum number of photons
62 MinPhotonEnergy (float): minimum photon energy [GeV]
63 MinPhotonTheta (float): minimum polar angle of photon [deg]
64 MaxPhotonTheta (float): maximum polar angle of photon [deg]
65 applyInCMS (bool): if true apply the P,Pt,theta, and energy cuts in the center of mass frame
66 stableParticles (bool): if true apply the selection criteria for stable particles in the generator
67 """
68
69 generatorpreselection = path.add_module('GeneratorPreselection',
70 nChargedMin=nChargedMin,
71 nChargedMax=nChargedMax,
72 MinChargedP=MinChargedP,
73 MinChargedPt=MinChargedPt,
74 MinChargedTheta=MinChargedTheta,
75 MaxChargedTheta=MaxChargedTheta,
76 nPhotonMin=nPhotonMin,
77 nPhotonMax=nPhotonMax,
78 MinPhotonEnergy=MinPhotonEnergy,
79 MinPhotonTheta=MinPhotonTheta,
80 MaxPhotonTheta=MaxPhotonTheta,
81 applyInCMS=applyInCMS,
82 stableParticles=stableParticles
83 )
84
85 # empty path for unwanted events
86 generatorpreselection.if_value('<11', emptypath)
87
88
89def add_aafh_generator(
90 path,
91 finalstate='',
92 preselection=False,
93 enableTauDecays=True,
94 minmass=0.5,
95 subweights=[],
96 maxsubweight=1,
97 maxfinalweight=3.0,
98 eventType=''):
99 """
100 Add the default two photon generator for four fermion final states
101
102 Parameters:
103 path (basf2.Path): path where the generator should be added
104 finalstate (str): either "e+e-e+e-", "e+e-mu+mu-", "e+e-tau+tau-", "mu+mu-mu+mu-" or "mu+mu-tau+tau-"
105 preselection (bool): if True, select events with at least one medium pt particle in the CDC acceptance
106 enableTauDecays (bool): if True, allow tau leptons to decay (using EvtGen)
107 minmass (float): minimum invariant mass
108 subweights (list(float)): list of four or eight values (first four are interpreted as WAP, rest as WBP)
109 which specify the relative weights for each of the four sub generators
110 maxsubweight (float): maximum expected subgenerator weight for rejection scheme
111 maxfinalweight (float): maximum expected final weight for rejection scheme
112 eventType (str) : event type information
113 """
114
115 if finalstate == 'e+e-e+e-':
116 aafh_mode = 5
117 if not subweights: # default subweights are for minmass=0.5
118 aafh_subgeneratorWeights = [1.0, 7.986e+01, 5.798e+04, 3.898e+05, 1.0, 1.664e+00, 2.812e+00, 7.321e-01]
119 else:
120 aafh_subgeneratorWeights = subweights
121 if abs(minmass - 0.5) > 0.01 and not subweights:
122 b2.B2WARNING("add_aafh_generator: non default invariant mass cut without updated subweights requested!")
123 elif finalstate == 'e+e-mu+mu-':
124 aafh_mode = 3
125 if not subweights: # default subweights are for minmass=0.5
126 aafh_subgeneratorWeights = [1.000e+00, 1.520e+01, 3.106e+03, 6.374e+03, 1.000e+00, 1.778e+00, 6.075e+00, 6.512e+00]
127 else:
128 aafh_subgeneratorWeights = subweights
129 if abs(minmass - 0.5) > 0.01 and not subweights:
130 b2.B2WARNING("add_aafh_generator: non default invariant mass cut without updated subweights requested!")
131 elif finalstate == 'e+e-tau+tau-':
132 aafh_mode = 4
133 minmass = 0
134 if not subweights:
135 aafh_subgeneratorWeights = [1.000e+00, 2.214e+00, 1.202e+01, 1.536e+01, 1.000e+00, 1.664e+00, 1.680e+01, 6.934e+00]
136 else:
137 aafh_subgeneratorWeights = subweights
138 if preselection:
139 b2.B2WARNING(
140 f"You requested a generator preselection for the final state {finalstate}: "
141 "please consider to remove it, since the cross section is small.")
142 elif finalstate == 'mu+mu-mu+mu-':
143 aafh_mode = 2
144 minmass = 0
145 maxsubweight = 3
146 if not subweights:
147 aafh_subgeneratorWeights = [0.000e+00, 0.000e+00, 1.000e+00, 3.726e+00, 1.000e+00, 1.778e+00, 1.000e+00, 1.094e+00]
148 else:
149 aafh_subgeneratorWeights = subweights
150 if preselection:
151 b2.B2WARNING(
152 f"You requested a generator preselection for the final state {finalstate}: "
153 "please consider to remove it, since the cross section is small.")
154 elif finalstate == 'mu+mu-tau+tau-':
155 aafh_mode = 1
156 minmass = 0
157 maxsubweight = 3
158 if not subweights:
159 aafh_subgeneratorWeights = [0.000e+00, 0.000e+00, 1.000e+00, 1.715e+00, 1.000e+00, 1.778e+00, 1.000e+00, 6.257e-01]
160 else:
161 aafh_subgeneratorWeights = subweights
162 if preselection:
163 b2.B2WARNING(
164 f"You requested a generator preselection for the final state {finalstate}: "
165 "please consider to remove it, since the cross section is small.")
166 elif finalstate == 'tau+tau-tau+tau-':
167 b2.B2FATAL(f"AAFH is not able to generate the {finalstate} final state. Please use KoralW instead.")
168 else:
169 b2.B2FATAL(f"add_aafh_generator final state not supported: {finalstate}")
170
171 aafh_maxSubgeneratorWeight = maxsubweight
172 aafh_maxFinalWeight = maxfinalweight
173
174 path.add_module(
175 'AafhInput',
176 mode=aafh_mode,
177 rejection=2,
178 maxSubgeneratorWeight=aafh_maxSubgeneratorWeight,
179 maxFinalWeight=aafh_maxFinalWeight,
180 subgeneratorWeights=aafh_subgeneratorWeights,
181 suppressionLimits=[1e100] * 4,
182 minMass=minmass,
183 eventType=eventType,
184 )
185
186 if preselection:
187 generator_emptypath = b2.create_path()
188 add_generator_preselection(
189 path=path,
190 emptypath=generator_emptypath,
191 nChargedMin=1,
192 MinChargedPt=0.1,
193 MinChargedTheta=17.0,
194 MaxChargedTheta=150.0)
195
196 if 'tau+tau-' in finalstate:
197 if enableTauDecays:
198 path.add_module('EvtGenDecay')
199 else:
200 b2.B2WARNING("The tau decays will not be generated.")
201
202
203def add_kkmc_generator(path, finalstate='', signalconfigfile='', useTauolaBelle=False, tauinputfile='', eventType=''):
204 """
205 Add the default muon pair and tau pair generator KKMC.
206 For tau decays, TauolaBelle and TauolaBelle2 are available.
207 Signal events can be produced setting a configuration file. Please notice that the configuration files for
208 TauolaBelle and TauolaBelle2 has a very different structure (see the examples below generators/examples).
209
210 Parameters:
211 path (basf2.Path): path where the generator should be added
212 finalstate(str): either "mu-mu+" or "tau-tau+"
213 signalconfigfile(str): File with configuration of the signal event to generate. It doesn't affect mu-mu+ decays.
214 useTauolaBelle(bool): If true, tau decay is driven by TauolaBelle. Otherwise TauolaBelle2 is used.
215 It doesn't affect mu-mu+ decays.
216 tauinputfile(str) : File to override KK2f_defaults. Only [sometimes] needed when tau decay is driven by TauolaBelle.
217 eventType (str) : event type information
218 """
219
220
221 kkmc_inputfile = b2.find_file('data/generators/kkmc/tauola_bbb.input.dat')
222
223
224 kkmc_logfile = 'kkmc_tautau.txt'
225
226
227 kkmc_config = b2.find_file('data/generators/kkmc/KK2f_defaults.dat')
228
229
230 kkmc_tauconfigfile = ''
231
232 if finalstate == 'tau+tau-':
233 b2.B2WARNING("add_kkmc_generator: please set finalstate as 'tau-tau+'. 'tau+tau-' will be deprecated in the future"
234 " for consistency in the configuration files.")
235 finalstate = 'tau-tau+'
236 if finalstate == 'mu+mu-':
237 b2.B2WARNING("add_kkmc_generator: please set finalstate as 'mu-mu+'. 'mu+mu-' will be deprecated in the future for"
238 " consistency in the configuration files.")
239 finalstate = 'mu-mu+'
240
241 if finalstate == 'tau-tau+':
242 if useTauolaBelle:
243 b2.B2INFO("Generating tau pair events with TauolaBelle")
244
245 kkmc_inputfile = b2.find_file('data/generators/kkmc/tau.input.dat')
246 kkmc_tauconfigfile = b2.find_file('data/generators/kkmc/tau_decaytable.dat')
247
248 if not signalconfigfile == '':
249 b2.B2INFO(f"Using config file defined by user: {signalconfigfile}")
250 if useTauolaBelle:
251 kkmc_tauconfigfile = b2.find_file(signalconfigfile)
252 else:
253 kkmc_inputfile = b2.find_file(signalconfigfile)
254
255 if not tauinputfile == '':
256 kkmc_inputfile = b2.find_file(tauinputfile)
257
258 elif finalstate == 'mu-mu+':
259 kkmc_inputfile = b2.find_file('data/generators/kkmc/mu.input.dat')
260 kkmc_logfile = 'kkmc_mumu.txt'
261
262 else:
263 b2.B2FATAL("add_kkmc_generator final state not supported: {}".format(finalstate))
264
265 # use KKMC to generate lepton pairs
266 path.add_module(
267 'KKGenInput',
268 tauinputFile=kkmc_inputfile,
269 KKdefaultFile=kkmc_config,
270 taudecaytableFile=kkmc_tauconfigfile,
271 kkmcoutputfilename=kkmc_logfile,
272 eventType=eventType,
273 )
274
275
276def add_evtgen_generator(path, finalstate='', signaldecfile=None, coherentMixing=True, parentParticle='Upsilon(4S)', eventType=''):
277 """
278 Add EvtGen for mixed and charged BB
279
280 Parameters:
281 path (basf2.Path): path where the generator should be added
282 finalstate (str): Either "charged" for generation of generic B+/B-, "mixed" for generic B0/anti-B0, or "signal" for
283 generation of user-defined signal mode
284 signaldecfile (str): path to decfile which defines the signal decay to be generated
285 (only needed if ``finalstate`` set to "signal")
286 coherentMixing: Either True or False. Switches on or off the coherent decay of the B0-B0bar pair.
287 It should always be True, unless you are generating Y(5,6S) -> BBar. In the latter case,
288 setting it False solves the internal limitation of Evtgen that allows to make a
289 coherent decay only starting from the Y(4S).
290 parentParticle (str): initial state (used only if it is not Upsilon(4S).
291 eventType (str) : event type information
292 """
293 evtgen_userdecfile = b2.find_file('data/generators/evtgen/charged.dec')
294
295 if parentParticle != 'Upsilon(3S)' and parentParticle != 'Upsilon(4S)'\
296 and parentParticle != 'Upsilon(5S)' and parentParticle != 'Upsilon(6S)':
297 b2.B2FATAL("add_evtgen_generator initial state not supported: {}".format(parentParticle))
298
299 if finalstate == 'charged':
300 pass
301 elif finalstate == 'mixed':
302 evtgen_userdecfile = b2.find_file('data/generators/evtgen/mixed.dec')
303 elif finalstate == 'signal':
304 evtgen_userdecfile = signaldecfile
305 else:
306 b2.B2FATAL("add_evtgen_generator final state not supported: {}".format(finalstate))
307
308 if signaldecfile and finalstate in ['charged', 'mixed']:
309 b2.B2WARNING("ignoring decfile: {}".format(signaldecfile))
310
311 # use EvtGen
312 if parentParticle == 'Upsilon(3S)':
313 if finalstate != 'signal':
314 b2.B2FATAL("add_evtgen_generator initial state {} is supported only with 'signal' final state".format(parentParticle))
315 if coherentMixing:
316 coherentMixing = False
317 b2.B2WARNING("add_evtgen_generator initial state {} has no BB mixing, now switching coherentMixing OFF"
318 .format(parentParticle))
319
320 if parentParticle == 'Upsilon(5S)':
321 if finalstate != 'signal':
322 b2.B2FATAL("add_evtgen_generator initial state {} is supported only with 'signal' final state".format(parentParticle))
323 if coherentMixing:
324 coherentMixing = False
325 b2.B2WARNING(
326 "add_evtgen_generator initial state {} is supported only with false coherentMixing, now switching it OFF"
327 .format(parentParticle))
328 pdg.load(b2.find_file('decfiles/dec/Y5S.pdl'))
329
330 if parentParticle == 'Upsilon(6S)':
331 if finalstate != 'signal':
332 b2.B2FATAL("add_evtgen_generator initial state {} is supported only with 'signal' final state".format(parentParticle))
333 if coherentMixing:
334 coherentMixing = False
335 b2.B2WARNING(
336 "add_evtgen_generator initial state {} is supported only with false coherentMixing, now switching it OFF"
337 .format(parentParticle))
338 pdg.load(b2.find_file('decfiles/dec/Y6S.pdl'))
339
340 path.add_module(
341 'EvtGenInput',
342 userDECFile=evtgen_userdecfile,
343 CoherentMixing=coherentMixing,
344 ParentParticle=parentParticle,
345 eventType=eventType,
346 )
347
348
349def add_continuum_generator(path, finalstate, userdecfile='', *, skip_on_failure=True, eventType=''):
350 """
351 Add the default continuum generators KKMC + PYTHIA including their default decfiles and PYTHIA settings
352
353 See Also:
354 `add_inclusive_continuum_generator()` to add continuum generation with preselected particles
355
356 Parameters:
357 path (basf2.Path): path where the generator should be added
358 finalstate (str): uubar, ddbar, ssbar, ccbar
359 userdecfile (str): EvtGen decfile used for particle decays
360 skip_on_failure (bool): If True stop event processing right after
361 fragmentation fails. Otherwise continue normally
362 eventType (str) : event type information
363 """
364
365
366 kkmc_inputfile = b2.find_file('data/generators/kkmc/uubar_nohadronization.input.dat')
367
368
369 kkmc_logfile = 'kkmc_uubar.txt'
370
371
372 pythia_config = b2.find_file('data/generators/modules/fragmentation/pythia_belle2.dat')
373
374
375 decay_user = b2.find_file('data/generators/modules/fragmentation/dec_belle2_qqbar.dec')
376 if userdecfile == '':
377 pass
378 else:
379 b2.B2INFO('Replacing default user decfile: {}'.format(userdecfile))
380 decay_user = userdecfile
381
382
383 kkmc_config = b2.find_file('data/generators/kkmc/KK2f_defaults.dat')
384
385
386 decay_file = get_default_decayfile()
387
388 if finalstate == 'uubar':
389 pass
390 elif finalstate == 'ddbar':
391 kkmc_inputfile = b2.find_file('data/generators/kkmc/ddbar_nohadronization.input.dat')
392 kkmc_logfile = 'kkmc_ddbar.txt'
393 elif finalstate == 'ssbar':
394 kkmc_inputfile = b2.find_file('data/generators/kkmc/ssbar_nohadronization.input.dat')
395 kkmc_logfile = 'kkmc_ssbar.txt'
396 elif finalstate == 'ccbar':
397 kkmc_inputfile = b2.find_file('data/generators/kkmc/ccbar_nohadronization.input.dat')
398 pythia_config = b2.find_file('data/generators/modules/fragmentation/pythia_belle2_charm.dat')
399 kkmc_logfile = 'kkmc_ccbar.txt'
400 else:
401 b2.B2FATAL("add_continuum_generator final state not supported: {}".format(finalstate))
402
403 # use KKMC to generate qqbar events (no fragmentation at this stage)
404 path.add_module(
405 'KKGenInput',
406 tauinputFile=kkmc_inputfile,
407 KKdefaultFile=kkmc_config,
408 taudecaytableFile='',
409 kkmcoutputfilename=kkmc_logfile,
410 eventType=eventType,
411 )
412
413 # add the fragmentation module to fragment the generated quarks into hadrons
414 # using PYTHIA8
415 fragmentation = path.add_module(
416 'Fragmentation',
417 ParameterFile=pythia_config,
418 ListPYTHIAEvent=0,
419 UseEvtGen=1,
420 DecFile=decay_file,
421 UserDecFile=decay_user,
422 )
423
424 if skip_on_failure:
425 # branch to an empty path if PYTHIA failed, this will change the number of events
426 # but the file meta data will contain the total number of generated events
427 generator_emptypath = b2.create_path()
428 fragmentation.if_value('<1', generator_emptypath)
429
430
431def add_inclusive_continuum_generator(path, finalstate, particles, userdecfile='',
432 *, include_conjugates=True, max_iterations=100000, eventType=''):
433 """
434 Add continuum generation but require at least one of the given particles be
435 present in the event.
436
437 For example to only generate ccbar events which contain a "D*+" or an
438 electron one could would use
439
440 >>> add_inclusive_continuum_generator(path, "ccbar", ["D*+", 11])
441
442 If you are unsure how the particles are named in Belle II please have a look
443 at the ``b2help-particles`` executable or the `pdg` python module.
444
445 See Also:
446 `add_continuum_generator()` to add continuum generation without preselection
447
448 Parameters:
449 path (basf2.Path): path to which the generator should be added
450 finalstate (str): uubar, ddbar, ssbar, ccbar
451 particles (list): A list of particle names or pdg codes. An event is
452 only accepted if at lease one of those particles appears in the event.
453 userdecfile (str): EvtGen decfile used for particle decays
454 include_conjugates (bool): If True (default) accept the event also if a
455 charge conjugate of the given particles is found
456 max_iterations (int): maximum tries per event to generate the requested
457 particle. If exceeded processing will be stopped with a
458 `FATAL <LogLevel.FATAL>` error so for rare particles one might need a
459 larger number.
460 eventType (str) : event type information
461 """
462 loop_path = b2.create_path()
463 # we might run this more than once so make sure we remove any particles
464 # before generating new ones
465 loop_path.add_module("PruneDataStore", keepMatchedEntries=False, matchEntries=["MCParticles", "EventExtraInfo"])
466 # add the generator but make sure it doesn't stop processing on
467 # fragmentation failure as is this currently not supported by do_while
468 add_continuum_generator(loop_path, finalstate, userdecfile, skip_on_failure=False, eventType=eventType)
469 # check for the particles we want
470 loop_path.add_module("InclusiveParticleChecker", particles=particles, includeConjugates=include_conjugates)
471 # Done, add this to the path and iterate it until we found our particle
472 path.do_while(loop_path, max_iterations=max_iterations)
473
474
475def add_bhwide_generator(path, minangle=0.5):
476 """
477 Add the high precision QED generator BHWIDE to the path. Settings are the default L1/HLT study settings
478 with a cross section of about 124000 nb (!)
479
480 Parameters:
481 path (basf2.Path): path where the generator should be added
482 minangle (float): minimum angle of the outgoing electron/positron in the CMS in degrees
483 """
484
485 if minangle < 0.0 or minangle > 180.0:
486 b2.B2FATAL("add_bhwide_generator minimum angle too small (<0.0) or too large (>180): {}".format(minangle))
487
488 bhwide = path.add_module("BHWideInput")
489 bhwide.param('ScatteringAngleRangePositron', [minangle, 180.0 - minangle])
490 bhwide.param('ScatteringAngleRangeElectron', [minangle, 180.0 - minangle])
491 bhwide.param('MaxAcollinearity', 180.0)
492 bhwide.param('MinEnergy', 0.10)
493 bhwide.param('VacuumPolarization', 'burkhardt')
494 bhwide.param('WeakCorrections', True)
495 bhwide.param('WtMax', 3.0)
496
497
498def add_babayaganlo_generator(
499 path,
500 finalstate='',
501 minenergy=0.01,
502 minangle=10.0,
503 fmax=-1.0,
504 generateInECLAcceptance=False,
505 eventType=''):
506 '''
507 Add the high precision QED generator BabaYaga@NLO to the path.
508
509 Parameters:
510 path (basf2.Path): path where the generator should be added.
511 finalstate (str): ee (e+e-), mm(mu+mu-) or gg (gammagamma).
512 minenergy (float): minimum particle (leptons for 'ee' or 'mm', photons for 'gg') energy in GeV.
513 minangle (float): angular range from minangle to 180-minangle for primary particles (in degrees).
514 fmax (float): maximum of differential cross section weight. This parameter should be set only by experts.
515 generateInECLAcceptance (bool): if True, the GeneratorPreselection module is used to select only events
516 with both the primary particles within the ECL acceptance.
517 eventType (str) : event type information
518 '''
519
520 if finalstate != 'ee' and finalstate != 'gg' and finalstate != 'mm':
521 b2.B2FATAL(f'add_babayaganlo_generator final state not supported: {finalstate}')
522
523 if not (fmax == -1.0):
524 b2.B2WARNING(f'The BabayagaNLOInput parameter "FMax" will be set to {fmax} instead to the default value (-1.0). '
525 'Please do not do this, unless you are extremely sure about this choice.')
526
527 path.add_module(
528 'BabayagaNLOInput',
529 eventType=eventType,
530 FinalState=finalstate,
531 ScatteringAngleRange=[minangle, 180.0 - minangle],
532 MinEnergy=minenergy,
533 FMax=fmax
534 )
535
536 if generateInECLAcceptance:
537 b2.B2INFO(f'The final state {finalstate} is preselected requiring both primary particles within the ECL acceptance.')
538 emptypath = b2.Path()
539 add_generator_preselection(path=path,
540 emptypath=emptypath,
541 applyInCMS=False)
542
543 if finalstate == 'ee' or finalstate == 'mm':
544 b2.set_module_parameters(path=path,
545 name='GeneratorPreselection',
546 nChargedMin=2,
547 MinChargedTheta=12.4,
548 MaxChargedTheta=155.1,)
549
550 elif finalstate == 'gg':
551 b2.set_module_parameters(path=path,
552 name='GeneratorPreselection',
553 nPhotonMin=2,
554 MinPhotonTheta=12.4,
555 MaxPhotonTheta=155.1)
556
557
558def add_phokhara_generator(path, finalstate='', eventType=''):
559 """
560 Add the high precision QED generator PHOKHARA to the path. Almost full
561 acceptance settings for photons and hadrons/muons.
562
563 Parameters:
564 path (basf2.Path): path where the generator should be added
565 finalstate (str): One of the following final states: "mu+mu-", "pi+pi-", "pi+pi-pi0", "pi+pi-pi+pi-" (or "2(pi+pi-)"),
566 "pi+pi-pi0pi0" or ("pi+pi-2pi0"), "pi+pi-eta", "K+K-", "K0K0bar", "ppbar", "n0n0bar" or "Lambda0Lambda0bar"
567 eventType (str) : event type information
568 """
569
570 if finalstate == 'mu+mu-':
571 path.add_module(
572 'PhokharaInput',
573 FinalState=0, # mu+mu-
574 LO=0, NLO=1, QED=0, # use full two loop corrections
575 MinInvMassHadrons=0.,
576 eventType=eventType,
577 ).set_name('PHOKHARA_mumuISR')
578
579 elif finalstate == 'pi+pi-':
580 path.add_module(
581 'PhokharaInput',
582 FinalState=1, # pi+pi-
583 LO=0, NLO=1, QED=0, # use full two loop corrections
584 MinInvMassHadrons=0.,
585 eventType=eventType,
586 ).set_name('PHOKHARA_pipiISR')
587
588 elif finalstate == 'pi+pi-pi0':
589 path.add_module(
590 'PhokharaInput',
591 FinalState=8, # pi+pi-pi0
592 LO=0, NLO=1, QED=0, # use full two loop corrections
593 MinInvMassHadrons=0.,
594 eventType=eventType,
595 ).set_name('PHOKHARA_pipipi0ISR')
596
597 elif finalstate == 'pi+pi-pi+pi-' or finalstate == '2(pi+pi-)':
598 path.add_module(
599 'PhokharaInput',
600 FinalState=3, # 2(pi+pi-)
601 LO=0, NLO=1, QED=0, # use full two loop corrections
602 MinInvMassHadrons=0.,
603 eventType=eventType,
604 ).set_name('PHOKHARA_pipipipiISR')
605
606 elif finalstate == 'pi+pi-pi0pi0' or finalstate == 'pi+pi-2pi0':
607 path.add_module(
608 'PhokharaInput',
609 FinalState=2, # pi+pi-2pi0
610 LO=0, NLO=1, QED=0, # use full two loop corrections
611 MinInvMassHadrons=0.,
612 eventType=eventType,
613 ).set_name('PHOKHARA_pipipi0pi0ISR')
614
615 elif finalstate == 'pi+pi-eta':
616 path.add_module(
617 'PhokharaInput',
618 FinalState=10, # pi+pi-eta
619 LO=0, NLO=1, QED=0, # use full two loop corrections
620 MinInvMassHadrons=0.,
621 eventType=eventType,
622 ).set_name('PHOKHARA_pipietaISR')
623
624 elif finalstate == 'K+K-':
625 path.add_module(
626 'PhokharaInput',
627 FinalState=6, # K+K-
628 LO=0, NLO=1, QED=0, # use full two loop corrections
629 MinInvMassHadrons=0.,
630 eventType=eventType,
631 ).set_name('PHOKHARA_KKISR')
632
633 elif finalstate == 'K0K0bar':
634 path.add_module(
635 'PhokharaInput',
636 FinalState=7, # K0K0bar
637 LO=0, NLO=1, QED=0, # use full two loop corrections
638 MinInvMassHadrons=0.,
639 eventType=eventType,
640 ).set_name('PHOKHARA_K0K0barISR')
641
642 elif finalstate == 'ppbar':
643 path.add_module(
644 'PhokharaInput',
645 FinalState=4, # ppbar
646 LO=0, NLO=1, QED=0, # use full two loop corrections
647 MinInvMassHadrons=0.,
648 eventType=eventType,
649 ).set_name('PHOKHARA_ppbarISR')
650
651 elif finalstate == 'n0n0bar':
652 path.add_module(
653 'PhokharaInput',
654 FinalState=5, # n0n0bar
655 LO=0, NLO=1, QED=0, # use full two loop corrections
656 MinInvMassHadrons=0.,
657 eventType=eventType,
658 ).set_name('PHOKHARA_n0n0barISR')
659
660 elif finalstate == 'Lambda0Lambda0bar':
661 b2.B2WARNING(
662 "The Lambda pair production and decays are implemented at the leading order only,"
663 " but as the expected number of events is modest, the accuracy of the code is"
664 " sufficient for the description of this process.")
665 path.add_module(
666 'PhokharaInput',
667 FinalState=9, # Lambda0Lambda0bar
668 LO=0, NLO=0, QED=0, # use only one loop corrections
669 MinInvMassHadrons=0.,
670 eventType=eventType,
671 ).set_name('PHOKHARA_Lambda0Lambda0barISR')
672
673 else:
674 b2.B2FATAL(f"add_phokhara_generator final state not supported: {finalstate}")
675
676
677def add_phokhara_evtgen_combination(
678 path, final_state_particles, user_decay_file,
679 beam_energy_spread=True, isr_events=False, min_inv_mass_vpho=0.0,
680 max_inv_mass_vpho=0.0, eventType=''):
681 """
682 Add combination of PHOKHARA and EvtGen to the path. Phokhara is
683 acting as ISR generator by generating e+ e- -> mu+ mu-, the muon pair is
684 then replaced by a virtual photon. Finally, the virtual photon is
685 decayed by EvtGen.
686
687 Parameters:
688 path (basf2.Path): Path where the generator should be added.
689 final_state_particles (list): List of final-state particles of
690 the virtual-photon decay. It is necessary to define the correct
691 mass threshold in PHOKHARA. For example, for the process
692 e+ e- -> J/psi eta_c, the list should be ['J/psi', 'eta_c'];
693 it does not depend on subsequent J/psi or eta_c decays.
694 user_decay_file (str): Name of EvtGen user decay file. The initial
695 particle must be the virtual photon (vpho).
696 beam_energy_spread (bool): Whether beam-energy spread should be
697 simulated
698 isr_events (bool): If true, then PHOKHARA is used with ph0 (LO in
699 basf2) = 0, i.e. generation of processes with one photon.
700 min_inv_mass_hadrons(float): Minimum invariant mass of the virtual
701 photon. This parameter is used only if isr_events is true,
702 otherwise the minimum mass is computed from the masses of
703 the final-state particles.
704 max_inv_mass_hadrons(float): Maximum invariant mass of the virtual
705 photon. This parameter is used only if isr_events is true,
706 otherwise the maximum mass is not restricted.
707 eventType (str) : event type information
708 """
709
710 import pdg
711
712 phokhara = path.add_module('PhokharaInput')
713 phokhara.param('eventType', eventType)
714
715 # Generate muons and replace them by a virtual photon.
716 phokhara.param('FinalState', 0)
717 phokhara.param('ReplaceMuonsByVirtualPhoton', True)
718
719 # Simulate beam-energy spread. This performs initialization for every
720 # event, and, thus, very slow, but usually necessary except for testing.
721 phokhara.param('BeamEnergySpread', beam_energy_spread)
722
723 # Soft photon cutoff.
724 phokhara.param('Epsilon', 0.0001)
725
726 # Maximum search.
727 phokhara.param('SearchMax', 5000)
728
729 # Number of generation attempts for event.
730 phokhara.param('nMaxTrials', 25000)
731
732 # Use NNLO.
733 if isr_events:
734 phokhara.param('LO', 0)
735 else:
736 phokhara.param('LO', 1)
737 phokhara.param('NLO', 1)
738
739 # Use ISR only.
740 phokhara.param('QED', 0)
741
742 # No interference.
743 phokhara.param('IFSNLO', 0)
744
745 # Vacuum polarization by Fred Jegerlehner.
746 phokhara.param('Alpha', 1)
747
748 # Do not include narrow resonances.
749 phokhara.param('NarrowRes', 0)
750
751 # Angular ragnes.
752 phokhara.param('ScatteringAngleRangePhoton', [0., 180.])
753 phokhara.param('ScatteringAngleRangeFinalStates', [0., 180.])
754
755 # Minimal invariant mass of the muons and tagged photon combination.
756 phokhara.param('MinInvMassHadronsGamma', 0.)
757
758 if isr_events:
759 # Minimum squared invariant mass of muons (final state).
760 phokhara.param('MinInvMassHadrons',
761 min_inv_mass_vpho * min_inv_mass_vpho)
762 # Maximum squared invariant mass of muons (final state).
763 phokhara.param('MaxInvMassHadrons',
764 max_inv_mass_vpho * max_inv_mass_vpho)
765 else:
766 # Minimum squared invariant mass of muons (final state).
767 # Force application of this cut.
768 mass = 0
769 for particle in final_state_particles:
770 p = pdg.get(particle)
771 mass = mass + p.Mass()
772 phokhara.param('MinInvMassHadrons', mass * mass)
773 phokhara.param('ForceMinInvMassHadronsCut', True)
774 # Maximum squared invariant mass of muons (final state).
775 # It is set to a large value, resulting in unrestricted mass.
776 phokhara.param('MaxInvMassHadrons', 200.0)
777
778 # Minimal photon energy.
779 phokhara.param('MinEnergyGamma', 0.01)
780
781 # EvtGen.
782 evtgen_decay = path.add_module('EvtGenDecay')
783 evtgen_decay.param('UserDecFile', user_decay_file)
784
785
786def add_koralw_generator(path, finalstate='', enableTauDecays=True, eventType=''):
787 """
788 Add KoralW generator for radiative four fermion final states (only four leptons final states are currently supported).
789
790 Parameters:
791 path (basf2.Path): path where the generator should be added
792 finalstate (str): either 'e+e-e+e-', 'e+e-mu+mu-', 'e+e-tau+tau-', 'mu+mu-mu+mu-', 'mu+mu-tau+tau-' or 'tau+tau-tau+tau-'
793 enableTauDecays (bool): if True, allow tau leptons to decay (using EvtGen)
794 eventType (str) : event type information
795 """
796
797 decayFile = ''
798 if finalstate == 'e+e-e+e-':
799 decayFile = b2.find_file('data/generators/koralw/KoralW_eeee.data')
800 elif finalstate == 'e+e-mu+mu-':
801 decayFile = b2.find_file('data/generators/koralw/KoralW_eeMuMu.data')
802 elif finalstate == 'e+e-tau+tau-':
803 decayFile = b2.find_file('data/generators/koralw/KoralW_eeTauTau.data')
804 elif finalstate == 'mu+mu-mu+mu-':
805 decayFile = b2.find_file('data/generators/koralw/KoralW_MuMuMuMu.data')
806 elif finalstate == 'mu+mu-tau+tau-':
807 decayFile = b2.find_file('data/generators/koralw/KoralW_MuMuTauTau.data')
808 elif finalstate == 'tau+tau-tau+tau-':
809 decayFile = b2.find_file('data/generators/koralw/KoralW_TauTauTauTau.data')
810 else:
811 b2.B2FATAL(f'add_koralw_generator final state not supported: {finalstate}')
812
813 path.add_module('KoralWInput',
814 UserDataFile=decayFile,
815 eventType=eventType)
816
817 if 'tau+tau-' in finalstate:
818 if enableTauDecays:
819 path.add_module('EvtGenDecay')
820 else:
821 b2.B2WARNING('The tau decays will not be generated.')
822
823
824def add_cosmics_generator(path, components=None,
825 global_box_size=None, accept_box=None, keep_box=None,
826 geometry_xml_file='geometry/Beast2_phase2.xml',
827 cosmics_data_dir='data/generators/modules/cryinput/',
828 setup_file='generators/scripts/cry.setup',
829 data_taking_period='gcr2017', top_in_counter=False):
830 """
831 Add the cosmics generator CRY with the default parameters to the path.
832
833 Warning:
834 Please remember to also change the reconstruction accordingly, if you
835 set "special" parameters here!
836
837 Parameters:
838 path (basf2.Path): path where the generator should be added
839 components (list(str)): list of geometry components to add in the
840 geometry module, or None for all components.
841 global_box_size (tuple(float, float, float)): sets global length, width
842 and height (in meters) in which to generate.
843 Default is ``[100, 100, 100]``
844 accept_box (tuple(float, float, float)): sets the size of the accept box in meter.
845 As a default it is set to ``[8.0, 8.0, 8.0]`` (the Belle II detector size).
846 keep_box (tuple(float, float, float)): sets the size of the keep box (keep box >= accept box).
847 geometry_xml_file (str): Name of the xml file to use for the geometry.
848 cosmics_data_dir (str): parameter CosmicDataDir for the cry module (absolute or relative to the basf2 repo).
849 setup_file (str): location of the cry.setup file (absolute or relative to the basf2 repo)
850 data_taking_period (str): The cosmics generation will be added using the
851 parameters, that where used in this period of data taking. The
852 periods can be found in ``cdc/cr/__init__.py``.
853 top_in_counter (bool): time of propagation from the hit point to the PMT in the trigger counter is subtracted
854 (assuming PMT is put at -z of the counter).
855 """
856
857 b2.B2FATAL('''The function "add_cosmics_generator()" is outdated and it is currently not working: please replace
858
859 add_cosmics_generator(path=path)
860
861with
862
863 path.add_module('CRYInput')
864
865in your steering file (the module parameter "acceptance" has to be set, see the module documentation).''')
866
867 import cdc.cr as cosmics_setup
868
869 if global_box_size is None:
870 global_box_size = [100, 100, 100]
871 if accept_box is None:
872 accept_box = [8, 8, 8]
873 if keep_box is None:
874 keep_box = [8, 8, 8]
875
876 cosmics_setup.set_cdc_cr_parameters(data_taking_period)
877
878 if cosmics_setup.cosmics_period == "201607":
879 b2.B2FATAL("The data taking period 201607 is very special (geometry setup, PMTs etc). This is not handled "
880 "by this script! Please ask the CDC group, if you want to simulate this.")
881
882 if 'Gearbox' not in path:
883 override = [("/Global/length", str(global_box_size[0]), "m"),
884 ("/Global/width", str(global_box_size[1]), "m"),
885 ("/Global/height", str(global_box_size[2]), "m")]
886
887 if cosmics_setup.globalPhi:
888 override += [("/DetectorComponent[@name='CDC']//GlobalPhiRotation", str(cosmics_setup.globalPhi), "deg")]
889
890 path.add_module('Gearbox', override=override, fileName=geometry_xml_file,)
891
892 # detector geometry
893 if 'Geometry' not in path:
894 geometry = path.add_module('Geometry')
895 if components:
896 geometry.param('components', components)
897
898 cry = path.add_module('CRYInput')
899
900 # cosmic data input
901 cry.param('CosmicDataDir', b2.find_file(cosmics_data_dir))
902
903 # user input file
904 cry.param('SetupFile', b2.find_file(setup_file))
905
906 # acceptance half-lengths - at least one particle has to enter that box to use that event
907 cry.param('acceptLength', accept_box[0])
908 cry.param('acceptWidth', accept_box[1])
909 cry.param('acceptHeight', accept_box[2])
910 cry.param('maxTrials', 100000)
911
912 # keep half-lengths - all particles that do not enter the box are removed (keep box >= accept box)
913 # default was 6.0
914 cry.param('keepLength', keep_box[0])
915 cry.param('keepWidth', keep_box[1])
916 cry.param('keepHeight', keep_box[2])
917
918 # minimal kinetic energy - all particles below that energy are ignored
919 cry.param('kineticEnergyThreshold', 0.01)
920
921 # TODO: I still do not fully understand, when the cosmics selector is needed and when not
922 if cosmics_setup.cosmics_period not in ["normal", "gcr2017"]:
923 # Selector module.
924 cosmics_selector = b2.register_module('CDCCosmicSelector',
925 lOfCounter=cosmics_setup.lengthOfCounter,
926 wOfCounter=cosmics_setup.widthOfCounter,
927 xOfCounter=cosmics_setup.triggerPos[0],
928 yOfCounter=cosmics_setup.triggerPos[1],
929 zOfCounter=cosmics_setup.triggerPos[2],
930 phiOfCounter=0.,
931 TOP=top_in_counter,
932 propSpeed=cosmics_setup.lightPropSpeed,
933 TOF=1,
934 cryGenerator=True
935 )
936
937 path.add_module(cosmics_selector)
938
939 empty_path = b2.create_path()
940 cosmics_selector.if_false(empty_path)
941
942
943def add_treps_generator(path, finalstate='', useDiscreteAndSortedW=False, eventType=''):
944 """
945 Add TREPS generator to produce hadronic two-photon processes.
946
947 Parameters:
948 path (basf2.Path): path where the generator should be added
949 finalstate(str): "e+e-pi+pi-", "e+e-K+K-" or "e+e-ppbar"
950 useDiscreteAndSortedW(bool): if True, wListTableFile is used for discrete and sorted W. evtNumList must be set proper value.
951 eventType (str) : event type information
952 """
953
954 if finalstate == 'e+e-pi+pi-':
955 parameterFile = b2.find_file('generators/treps/data/parameterFiles/treps_par_pipi.dat')
956 differentialCrossSectionFile = b2.find_file('generators/treps/data/differentialCrossSectionFiles/pipidcs.dat')
957 wListTableFile = b2.find_file('generators/treps/data/wListFiles/wlist_table_pipi.dat')
958 elif finalstate == 'e+e-K+K-':
959 parameterFile = b2.find_file('generators/treps/data/parameterFiles/treps_par_kk.dat')
960 differentialCrossSectionFile = b2.find_file('generators/treps/data/differentialCrossSectionFiles/kkdcs.dat')
961 wListTableFile = b2.find_file('generators/treps/data/wListFiles/wlist_table_kk.dat')
962 elif finalstate == 'e+e-ppbar':
963 parameterFile = b2.find_file('generators/treps/data/parameterFiles/treps_par_ppbar.dat')
964 differentialCrossSectionFile = b2.find_file('generators/treps/data/differentialCrossSectionFiles/ppbardcs.dat')
965 wListTableFile = b2.find_file('generators/treps/data/wListFiles/wlist_table_ppbar.dat')
966 else:
967 b2.B2FATAL("add_treps_generator final state not supported: {}".format(finalstate))
968
969 # use TREPS to generate two-photon events.
970 path.add_module(
971 'TrepsInput',
972 ParameterFile=parameterFile,
973 DifferentialCrossSectionFile=differentialCrossSectionFile,
974 WListTableFile=wListTableFile,
975 UseDiscreteAndSortedW=useDiscreteAndSortedW,
976 MaximalQ2=1.0,
977 MaximalAbsCosTheta=1.01,
978 ApplyCosThetaCutCharged=True,
979 MinimalTransverseMomentum=0,
980 ApplyTransverseMomentumCutCharged=True,
981 eventType=eventType,
982 )
get(name)
Definition pdg.py:48
load(filename)
Definition pdg.py:122