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