Belle II Software  release-05-01-25
generateTwoPhoton.py
1 #!/usr/bin/env python
2 # -*- coding: utf-8 -*-
3 
4 # -------------------------------------------------------------------------------------
5 # BG simulation: two-photon (QED)
6 # usage:
7 # basf2 generateTwoPhoton.py equivTime_us num [sampleType phase outdir]
8 # arguments:
9 # equivTime_us equivalent SuperKEKB running time in micro-seconds
10 # num output file number
11 # sampleType one of: study, usual, PXD, ECL
12 # phase 2, 31 (= early phase 3) or 3
13 # outdir output directory path
14 # -------------------------------------------------------------------------------------
15 
16 from basf2 import *
17 import sys
18 import os
19 from background import add_output
20 
21 # read parameters
22 
23 argvs = sys.argv
24 argc = len(argvs)
25 
26 if argc == 3:
27  equivTime = argvs[1] # equivalent SuperKEKB running time in micro-seconds
28  num = argvs[2] # output file number
29  sampleType = 'usual' # study, usual, PXD, ECL
30  phase = 3 # phase number
31  outdir = 'output' # output directory path
32 elif argc == 4:
33  equivTime = argvs[1]
34  num = argvs[2]
35  sampleType = argvs[3]
36  phase = 3
37  outdir = 'output'
38 elif argc == 5:
39  equivTime = argvs[1]
40  num = argvs[2]
41  sampleType = argvs[3]
42  phase = int(argvs[4])
43  outdir = 'output'
44 elif argc == 6:
45  equivTime = argvs[1]
46  num = argvs[2]
47  sampleType = argvs[3]
48  phase = int(argvs[4])
49  outdir = argvs[5]
50 else:
51  print('usage:')
52  print(' basf2', argvs[0], 'equivTime_us num [sampleType phase outdir]')
53  print('arguments:')
54  print(' equivTime_us equivalent SuperKEKB running time in micro-seconds')
55  print(' num output file number')
56  print(' sampleType one of: study, usual, PXD, ECL')
57  print(' phase 2, 31 (=early phase 3) or 3')
58  print(' outdir output directory path')
59  sys.exit()
60 
61 # set parameters
62 
63 bgType = 'twoPhoton'
64 crossect = 7.28e6 # nb
65 
66 if phase == 3:
67  lumi = 800 # /nb/s
68 elif phase == 31:
69  lumi = 250 # /nb/s
70 elif phase == 2:
71  lumi = 20 # /nb/s
72 else:
73  print('phase ', phase, 'not supported')
74  sys.exit()
75 
76 realTime = float(equivTime) * 1000 # ns
77 numEvents = int(crossect * lumi * realTime * 1e-9)
78 fname = bgType + '_' + sampleType + '-phase' + str(phase) + '-' + num
79 outputFile = outdir + '/' + fname + '.root'
80 
81 if numEvents == 0:
82  B2ERROR('number of events is 0 -> increase equivTime_us')
83  sys.exit()
84 
85 # make output directory if doesn't exists
86 
87 if not os.path.exists(outdir):
88  os.makedirs(outdir)
89 
90 # log messages
91 
92 B2RESULT('Events to be generated: ' + str(numEvents) + ' - corresponds to ' + equivTime +
93  ' us of running at ' + str(lumi) + ' /nb/s (phase ' + str(phase) + ')')
94 B2RESULT('Output file: ' + outputFile)
95 
96 # Suppress messages and warnings during processing:
97 set_log_level(LogLevel.RESULT)
98 
99 # Create path
100 main = create_path()
101 
102 # Event info setter
103 eventinfosetter = register_module('EventInfoSetter')
104 eventinfosetter.param('evtNumList', [numEvents])
105 main.add_module(eventinfosetter)
106 
107 # Gearbox
108 gearbox = register_module('Gearbox')
109 if phase == 2:
110  gearbox.param('fileName', 'geometry/Beast2_phase2.xml')
111 elif phase == 31:
112  gearbox.param('fileName', 'geometry/Belle2_earlyPhase3.xml')
113 if sampleType == 'study':
114  gearbox.param('override', [
115  ("/DetectorComponent[@name='PXD']//ActiveChips", 'true', ''),
116  ("/DetectorComponent[@name='PXD']//SeeNeutrons", 'true', ''),
117  ("/DetectorComponent[@name='SVD']//ActiveChips", 'true', ''),
118  ("/DetectorComponent[@name='SVD']//SeeNeutrons", 'true', ''),
119  ("/DetectorComponent[@name='TOP']//BeamBackgroundStudy", '1', ''),
120  ("/DetectorComponent[@name='ARICH']//BeamBackgroundStudy", '1', ''),
121  ("/DetectorComponent[@name='ECL']//BeamBackgroundStudy", '1', ''),
122  ("/DetectorComponent[@name='KLM']//BeamBackgroundStudy", '1', ''),
123  ])
124 main.add_module(gearbox)
125 
126 # Generator (settings from Nakayama-san 15th campaign)
127 aafh = register_module('AafhInput')
128 aafh.param({
129  'mode': 5,
130  'rejection': 2,
131  'suppressionLimits': [1e100] * 4,
132  'minMass': 0.001,
133  'maxFinalWeight': 2.5,
134  'maxSubgeneratorWeight': 1.0,
135  'subgeneratorWeights': [1.000e+00, 2.216e+01, 3.301e+03, 6.606e+03, 1.000e+00,
136  1.675e+00, 5.948e+00, 6.513e+00],
137 })
138 main.add_module(aafh)
139 
140 # Geant geometry
141 geometry = register_module('Geometry')
142 geometry.param('useDB', False)
143 addComp = ["MagneticField3dQuadBeamline"]
144 # add beast detectors for early phase3
145 if phase == 31 and sampleType == 'study':
146  addComp.extend(["BEAMABORT", "MICROTPC", "CLAWS", "HE3TUBE"])
147 
148 geometry.param({"excludedComponents": ["MagneticField"],
149  "additionalComponents": addComp})
150 main.add_module(geometry)
151 
152 # Geant simulation
153 fullsim = register_module('FullSim')
154 if sampleType == 'study':
155  fullsim.param('PhysicsList', 'FTFP_BERT_HP')
156  fullsim.param('UICommandsAtIdle', ['/process/inactivate nKiller'])
157  fullsim.param('StoreAllSecondaries', True)
158  fullsim.param('SecondariesEnergyCut', 0.000001) # [MeV] need for CDC EB neutron flux
159 main.add_module(fullsim)
160 
161 # Show progress of processing
162 progress = register_module('Progress')
163 main.add_module(progress)
164 
165 # Output
166 if phase == 31:
167  phase = 3
168 add_output(main, bgType, realTime, sampleType, phase, fileName=outputFile)
169 
170 # Process events
171 process(main)
172 
173 # Print call statistics
174 print(statistics)