ATLAS Offline Software
Loading...
Searching...
No Matches
SFGenUtils.py
Go to the documentation of this file.
1# Copyright (C) 2002-2022 CERN for the benefit of the ATLAS collaboration
2
3
4import subprocess, os, shlex, re
5
6from AthenaCommon import Logging
7
8
9logger = Logging.logging.getLogger("SFGen_i")
10
11# This module defines the SFGen Utilities in the Athena framework
12
14
15# Init definition
16
17 def __init__(self, runArgs):
18 self.sfgenpath = os.environ['SFGENPATH']
19
20 #SFGen specific variables for input.DAT, see writeInputDAT function for more elaboration
21
22 self.rts = 13000 #collision energy (GeV)
23 if hasattr(runArgs,"ecmEnergy"):
24 self.rts = runArgs.ecmEnergy
25
26 self.coll = False
27 self.PDFname = 'SF_MSHT20qed_nnlo'
28 self.PDFmember = 0
29 self.proc = 1 #Please consult SFGen Manual https://sfgen.hepforge.org/
30 self.outtag = 'out'
31 self.diff = 'tot' #Interaction: elastic ('el'), single ('sd','sda','sdb') and double ('dd') dissociation.
32 self.SFerror = False
33 self.mixed = False
34 self.Zinit = False
35 self.subt = False
36 self.lep1 = 'mu+'
37 self.lep2 = 'mu+'
38 self.kmu = 1.
39 self.ncall = 50000
40 self.itmx = 10
41 self.prec = 1
42 self.ncall1 = 10000
43 self.inccall = 10000
44 self.itend = 1000
45
46 self.iseed = 13
47 if hasattr(runArgs,"randomSeed"):
48 self.iseed = runArgs.randomSeed
49
50 self.genunw = False
51
52 self.nev = "10"
53 if hasattr(runArgs,"maxEvents"):
54 self.nev = runArgs.maxEvents
55
56 self.erec = 'lhe' #output file type
57 self.ymin = -2.4 # Minimum dilepton rapidity
58 self.ymax = 2.4 # Maximum dilepton rapidity
59 self.mmin = 10. # Minimum dilepton mass
60 self.mmax = 1000. # Maximum dilepton mass
61
62 self.gencuts = True # Generate cuts below
63
64 self.ptamin = 20. # Minimum pT of outgoing object a
65 self.ptbmin = 20. # Minimum pT of outgoing object b
66 self.etaamin = -2.4 # Minimum eta of outgoing object a
67 self.etaamax = 2.4 # Maximum eta of outgoing object a
68 self.etabmin = -2.4 # Minimum eta of outgoing object b
69 self.etabmax = 2.4 # Maximum eta of outgoing object b
70 self.ptllmin = 0. # Minimun pT of dilepton pair
71
72# Writing of the input.DAT file
73
74 def toFortran(self):
75
76 def fortDouble(x):
77 return str(x)+"d0"
78 def fortInt(x):
79 return str(x)
80 def fortBool(x):
81 return '.true.' if x else '.false.'
82 def fortStr(x):
83 return "'{}'".format(x)
84
85 conf = ""
86 conf+="***********************************************************************************\n"
87 conf+="***********************************************************************************\n"
88 conf+="***********************************************************************************\n"
89 conf+=fortDouble(self.rts) + " ! [rts] : CMS collision energy (GeV) \n"
90 conf+="***********************************************************************************\n"
91 conf+="***********************************************************************************\n"
92 conf+=fortInt(self.proc) + " ! [proc] : Process number \n"
93 conf+=fortStr(self.outtag) + " ! [outtag] : for output file \n"
94 conf+=fortStr(self.diff) + " ! [diff] : elastic ('el'), single/double dissociation ('sd'/'dd') \n"
95 conf+=fortBool(self.SFerror) + " ! [SFerror] : Include error from SF input - increases run time \n"
96 conf+=fortBool(self.mixed) + " ! [mixed] : include mixed gam/Z + q diagrams \n"
97 conf+=fortBool(self.Zinit) + " ! [Zinit] : include Z bosons in initial state \n"
98 conf+=fortBool(self.subt) + " ! [subt] : calculate *only* (positive) subtraction term \n"
99 conf+=fortStr(self.lep1) + " ! [lep1] : for lepton-lepton scattering (proc=4) \n"
100 conf+=fortStr(self.lep2) + " ! [lep2] : for lepton-lepton scattering (proc=4) \n"
101 conf+=fortDouble(self.kmu) + " ! [kmu] : = mu(f,r)/mu0 \n"
102 conf+="***********************************************************************************\n"
103 conf+="************************** To run in collinear mode *****************************\n"
104 conf+="***********************************************************************************\n"
105 conf+=fortBool(self.coll) + " ! [coll] : use collinear approach + lhapdf \n"
106 conf+=fortStr(self.PDFname) + " ! [PDFname] : PDF set \n"
107 conf+=fortInt(self.PDFmember) + " ! [PDFmember] : PDF member \n"
108 conf+="***********************************************************************************\n"
109 conf+="*************Integration parameters************************************************\n"
110 conf+="***********************************************************************************\n"
111 conf+=fortInt(self.ncall) + " ! [ncall] : Number of calls for preconditioning \n"
112 conf+=fortInt(self.itmx) + " ! [itmx] : Number of iterations for preconditioning \n"
113 conf+=fortDouble(self.prec) + " ! [prec] : Relative accuracy (in %) in main run \n"
114 conf+=fortInt(self.ncall1) + " ! [ncall1] : Number of calls in first iteration \n"
115 conf+=fortInt(self.inccall) + " ! [inccall] : Number of increase calls per iteration \n"
116 conf+=fortInt(self.itend) + " ! [itend] : Maximum number of iterations \n"
117 conf+=fortInt(self.iseed) + " ! [iseed] : Random number seed (integer > 0) \n"
118 conf+="***********************************************************************************\n"
119 conf+="********************Unweighted events**********************************************\n"
120 conf+="***********************************************************************************\n"
121 conf+=fortBool(self.genunw) + " ! [genunw] : Generate unweighted events \n"
122 conf+=fortInt(int(self.nev)) + " ! [nev] : Number of events (preferably controlled by maxEvents option in Gen_tf command) \n"
123 conf+=fortStr(self.erec) + " ! [erec] : Event record format ('hepmc','lhe','hepevt') \n"
124 conf+="***********************************************************************************\n"
125 conf+="******************* general cuts ************************************************\n"
126 conf+="***********************************************************************************\n"
127 conf+=fortDouble(self.ymin) + " ! [ymin] : Minimum dilepton rapidity \n"
128 conf+=fortDouble(self.ymax) + " ! [ymax] : Maximum dilepton rapidity \n"
129 conf+=fortDouble(self.mmin) + " ! [mmin] : Minimum dilepton mass \n"
130 conf+=fortDouble(self.mmax) + " ! [mmax] : Maximum dilepton mass \n"
131 conf+=fortBool(self.gencuts) + " ! [gencuts] : Generate cuts below \n"
132 conf+="***********************************************************************************\n"
133 conf+="********** 2 body final states : p(a) + p(b) **************************************\n"
134 conf+="***********************************************************************************\n"
135 conf+=fortDouble(self.ptamin) + " ! [ptamin] \n"
136 conf+=fortDouble(self.ptbmin) + " ! [ptbmin] \n"
137 conf+=fortDouble(self.etaamin) + " ! [etaamin] \n"
138 conf+=fortDouble(self.etaamax) + " ! [etaamax] \n"
139 conf+=fortDouble(self.etabmin) + " ! [etabmin] \n"
140 conf+=fortDouble(self.etabmax) + " ! [etabmax] \n"
141 conf+=fortDouble(self.ptllmin) + " ! [ptllmin] \n"
142 conf+="***********************************************************************************\n"
143 conf+="***********************************************************************************\n"
144
145 return conf
146
147 def outputLHEFile(self):
148 return "evrecs/evrec"+self.outtag+".dat"
149
150
152
153 with open("input.DAT", "w") as outF:
154 outF.write(Init.toFortran())
155
156 return
157
158
159def run_command(command, stdin = None ):
160 """
161 Run a command and print output continuously
162 """
163 process = subprocess.Popen(shlex.split(command), stdout=subprocess.PIPE, stdin=stdin)
164 while True:
165 output = process.stdout.readline().decode("utf-8")
166 if output == '' and process.poll() is not None:
167 break
168 if output:
169 # remove ANSI escape formatting characters
170 reaesc = re.compile(r'(\x9B|\x1B\[)[0-?]*[ -\/]*[@-~]')
171 text = reaesc.sub('', output.strip())
172 logger.info(text)
173
174 rc = process.poll()
175 return rc
176
177# SFGen initialization
178
179def SFGenInitialize(Init, stdin=None):
180
181 logger.info("Starting SFGen Initialization")
182
183 if not os.path.exists('evrecs'):
184 os.makedirs('evrecs')
185 if not os.path.exists('outputs'):
186 os.makedirs('outputs')
187
188
189 try:
190 inputDAT = open('input.DAT')
191
192 except IOError:
193 raise Exception("Problem with file IO; potentially input.DAT not created correctly")
194 else:
195
196 try:
197 rc = run_command(Init.sfgenpath+"/bin/SFGen", inputDAT)
198
199 except OSError:
200 raise Exception("File not found")
201
202 except Exception:
203 raise Exception("Non-OSError or IOError in execution block")
204
205 if rc:
206 raise Exception('Unexpected error in sfgen execution in SFGenInitialize')
207
208 return
209
210# SFGen execution
211
212def SFGenExecute(Init):
213
214 logger.info("Starting SFGen Itself")
215
216 if not os.path.exists('evrecs'):
217 os.makedirs('evrecs')
218 if not os.path.exists('outputs'):
219 os.makedirs('outputs')
220
221
222 try:
223 inputDAT = open('input.DAT')
224 except IOError:
225 raise Exception ("Problem with IO; potentially input.DAT not created correctly")
226 else:
227
228 try:
229
230 rc = run_command(Init.sfgenpath+'/bin/SFGen', stdin=inputDAT)
231
232 except OSError:
233 raise Exception("SFGen executable or file not found")
234
235 except Exception:
236 raise Exception("Non-OSError or IOError in SFGen execution block")
237
238 if rc:
239 raise Exception('Unexpected error in sfgen execution in SFGenExecute')
240
241 return
242
243# SFGen running
244
245def SFGenRun(Init, genSeq):
246
247 # dump the job configuration for fortran code
248 print(Init.toFortran())
249
250 # attach config to genSequence for later using in showering
251 genSeq.SFGenConfig = Init
252
253 writeInputDAT(Init)
254 SFGenInitialize(Init)
255 SFGenExecute(Init)
256
257 return
void print(char *figname, TCanvas *c1)
__init__(self, runArgs)
Definition SFGenUtils.py:17
SFGenInitialize(Init, stdin=None)
SFGenRun(Init, genSeq)
run_command(command, stdin=None)
SFGenExecute(Init)
writeInputDAT(Init)