7 from AthenaCommon
import Logging
8 from ...decorators
import timed
9 from ...utility
import LHE, ProcessManager, SingleProcessThread
12 logger = Logging.logging.getLogger(
"PowhegControl")
15 @
timed(
"MadSpin post-processing")
17 """! Move output to correctly named file.
19 @param process External MadSpin process.
20 @param powheg_LHE_output Name of LHE file produced by PowhegBox.
22 @author James Robinson <james.robinson@cern.ch>
31 """! Construct MadSpin runcard.
33 @param input_LHE_events Input LHE file name.
34 @param process MadSpin process.
36 @author James Robinson <james.robinson@cern.ch>
39 logger.info(
"Constructing MadSpin runcard header")
40 pre_header =
"\n".
join([elem
for elem
in [LHE.opening_tag(input_LHE_events), LHE.comment_block(input_LHE_events)]
if elem])
41 header_contents = LHE.header_block(input_LHE_events).
replace(
"<header>",
"").
replace(
"</header>",
"").strip(
"\n")
42 post_header = LHE.init_block(input_LHE_events)
43 postamble = LHE.postamble(input_LHE_events)
46 with open(
"madspin_LHE_input.lhe",
"w")
as f_madspin_LHE:
47 f_madspin_LHE.write(
"{}\n".
format(pre_header))
48 f_madspin_LHE.write(
"<header>\n")
49 f_madspin_LHE.write(
"<mgversion>\n")
51 f_madspin_LHE.write(
"</mgversion>\n")
52 f_madspin_LHE.write(
"<mg5proccard>\n")
53 f_madspin_LHE.write(
"set group_subprocesses Auto\n")
54 f_madspin_LHE.write(
"set ignore_six_quark_processes False\n")
55 f_madspin_LHE.write(
"set loop_optimized_output True\n")
56 f_madspin_LHE.write(
"set gauge unitary\n")
57 f_madspin_LHE.write(
"set complex_mass_scheme False\n")
58 f_madspin_LHE.write(
"import model {}\n".
format(process.MadSpin_model))
59 if int(process.MadSpin_nFlavours) == 4:
60 f_madspin_LHE.write(
"define p = g u c d s u~ c~ d~ s~\n")
61 f_madspin_LHE.write(
"define j = g u c d s u~ c~ d~ s~\n")
62 elif int(process.MadSpin_nFlavours) == 5:
63 f_madspin_LHE.write(
"define p = g u c d s b u~ c~ d~ s~ b~\n")
64 f_madspin_LHE.write(
"define j = g u c d s b u~ c~ d~ s~ b~\n")
66 raise ValueError(
"'MadSpin_nFlavours' must be set to 4 or 5")
67 if process.MadSpin_taus_are_leptons:
68 f_madspin_LHE.write(
"define l+ = e+ mu+ ta+\n")
69 f_madspin_LHE.write(
"define l- = e- mu- ta-\n")
71 f_madspin_LHE.write(
"define l+ = e+ mu+\n")
72 f_madspin_LHE.write(
"define l- = e- mu-\n")
73 f_madspin_LHE.write(
"{}\n".
format(process.MadSpin_process))
74 f_madspin_LHE.write(
"output tchan\n")
75 f_madspin_LHE.write(
"</mg5proccard>\n")
76 f_madspin_LHE.write(
"<mgruncard>\n")
77 f_madspin_LHE.write(
"#0.01 = req_acc_FO ! needed for determining LO/NLO - see AGENE-1459\n")
78 f_madspin_LHE.write(
"{} = nevents\n".
format(LHE.event_counter(input_LHE_events)))
79 f_madspin_LHE.write(
"1 = lpp1 ! beam 1 type (0 = no PDF)\n")
80 f_madspin_LHE.write(
"1 = lpp2 ! beam 2 type (0 = no PDF)\n")
81 f_madspin_LHE.write(
"{} = ebeam1 ! beam 1 energy in GeV\n".
format(process.beam_energy))
82 f_madspin_LHE.write(
"{} = ebeam2 ! beam 2 energy in GeV\n".
format(process.beam_energy))
83 f_madspin_LHE.write(
"{} = bwcutoff\n".
format(process.bwcutoff))
84 f_madspin_LHE.write(
"</mgruncard>\n")
85 f_madspin_LHE.write(
"<slha>\n")
86 if process.MadSpin_paramcard !=
"":
87 with open(process.MadSpin_paramcard)
as file:
89 f_madspin_LHE.write(line)
91 f_madspin_LHE.write(
"######################################################################\n")
92 f_madspin_LHE.write(
"## PARAM_CARD AUTOMATICALY GENERATED BY MG5 ####\n")
93 f_madspin_LHE.write(
"######################################################################\n")
94 f_madspin_LHE.write(
"###################################\n")
95 f_madspin_LHE.write(
"## INFORMATION FOR LOOP\n")
96 f_madspin_LHE.write(
"###################################\n")
97 f_madspin_LHE.write(
"BLOCK LOOP #\n")
98 f_madspin_LHE.write(
" 1 {:e} # mu_r\n".
format(process.mass_Z))
99 f_madspin_LHE.write(
"###################################\n")
100 f_madspin_LHE.write(
"## INFORMATION FOR MASS\n")
101 f_madspin_LHE.write(
"###################################\n")
102 f_madspin_LHE.write(
"BLOCK MASS #\n")
103 f_madspin_LHE.write(
" 1 0.000000e+00 # d : 0.0\n")
104 f_madspin_LHE.write(
" 2 0.000000e+00 # u : 0.0\n")
105 f_madspin_LHE.write(
" 3 0.000000e+00 # s : 0.0\n")
106 f_madspin_LHE.write(
" 4 0.000000e+00 # c : 0.0\n")
107 f_madspin_LHE.write(
" 5 {:e} # mb\n".
format(process.mass_b))
108 f_madspin_LHE.write(
" 6 {:e} # mt\n".
format(process.mass_t))
109 f_madspin_LHE.write(
" 11 0.000000e+00 # e- : 0.0\n")
110 f_madspin_LHE.write(
" 12 0.000000e+00 # ve : 0.0\n")
111 f_madspin_LHE.write(
" 13 0.000000e+00 # mu- : 0.0\n")
112 f_madspin_LHE.write(
" 14 0.000000e+00 # vm : 0.0\n")
113 f_madspin_LHE.write(
" 15 {:e} # mta\n".
format(process.mass_tau))
114 f_madspin_LHE.write(
" 16 0.000000e+00 # vt : 0.0\n")
115 f_madspin_LHE.write(
" 21 0.000000e+00 # g : 0.0\n")
116 f_madspin_LHE.write(
" 22 0.000000e+00 # a : 0.0\n")
117 f_madspin_LHE.write(
" 23 {:e} # mz\n".
format(process.mass_Z))
118 f_madspin_LHE.write(
" 24 {:e} # w+\n".
format(process.mass_W))
119 f_madspin_LHE.write(
" 25 {:e} # mh\n".
format(process.mass_H))
120 f_madspin_LHE.write(
" 82 0.000000e+00 # gh : 0.0\n")
121 f_madspin_LHE.write(
"###################################\n")
122 f_madspin_LHE.write(
"## INFORMATION FOR SMINPUTS\n")
123 f_madspin_LHE.write(
"###################################\n")
124 f_madspin_LHE.write(
"BLOCK SMINPUTS #\n")
125 f_madspin_LHE.write(
" 1 {:e} # aewm1\n".
format(process.alphaem_inv))
126 f_madspin_LHE.write(
" 2 {:e} # gf\n".
format(process.G_F))
127 f_madspin_LHE.write(
" 3 {:e} # as\n".
format(process.alphaqcd))
128 f_madspin_LHE.write(
"###################################\n")
129 f_madspin_LHE.write(
"## INFORMATION FOR YUKAWA\n")
130 f_madspin_LHE.write(
"###################################\n")
131 f_madspin_LHE.write(
"BLOCK YUKAWA #\n")
132 f_madspin_LHE.write(
" 5 {:e} # ymb\n".
format(process.mass_b))
133 f_madspin_LHE.write(
" 6 {:e} # ymt\n".
format(process.mass_t))
134 f_madspin_LHE.write(
" 15 {:e} # ymtau\n".
format(process.mass_tau))
135 f_madspin_LHE.write(
"###################################\n")
136 f_madspin_LHE.write(
"## INFORMATION FOR QNUMBERS 82\n")
137 f_madspin_LHE.write(
"###################################\n")
138 f_madspin_LHE.write(
"BLOCK QNUMBERS 82 # gh\n")
139 f_madspin_LHE.write(
" 1 0 # 3 times electric charge\n")
140 f_madspin_LHE.write(
" 2 1 # number of spin states (2s+1)\n")
141 f_madspin_LHE.write(
" 3 8 # colour rep (1: singlet, 3: triplet, 8: octet)\n")
142 f_madspin_LHE.write(
" 4 1 # particle/antiparticle distinction (0=own anti)\n")
143 if "ckm" in process.MadSpin_model
or "CKM" in process.MadSpin_model:
144 logger.info(
"Adding Wolfenstein parameters in param_card for model {}".
format(process.MadSpin_model))
145 f_madspin_LHE.write(
"###################################\n")
146 f_madspin_LHE.write(
"INFORMATION FOR WOLFENSTEIN\n")
147 f_madspin_LHE.write(
"###################################\n")
148 f_madspin_LHE.write(
"Block wolfenstein\n")
149 f_madspin_LHE.write(
" 1 {:e} # lamWS\n".
format(process.Wolfenstein_lambda))
150 f_madspin_LHE.write(
" 2 {:e} # AWS\n".
format(process.Wolfenstein_A))
151 f_madspin_LHE.write(
" 3 {:e} # rhoWS\n".
format(process.Wolfenstein_rho))
152 f_madspin_LHE.write(
" 4 {:e} # etaWS\n".
format(process.Wolfenstein_eta))
153 f_madspin_LHE.write(
"#\n")
154 f_madspin_LHE.write(
"#*************************\n")
155 f_madspin_LHE.write(
"# Decay widths *\n")
156 f_madspin_LHE.write(
"#*************************\n")
157 f_madspin_LHE.write(
"#\n")
158 f_madspin_LHE.write(
"# PDG Width\n")
159 f_madspin_LHE.write(
"DECAY 1 0.000000e+00\n")
160 f_madspin_LHE.write(
"#\n")
161 f_madspin_LHE.write(
"# PDG Width\n")
162 f_madspin_LHE.write(
"DECAY 2 0.000000e+00\n")
163 f_madspin_LHE.write(
"#\n")
164 f_madspin_LHE.write(
"# PDG Width\n")
165 f_madspin_LHE.write(
"DECAY 3 0.000000e+00\n")
166 f_madspin_LHE.write(
"#\n")
167 f_madspin_LHE.write(
"# PDG Width\n")
168 f_madspin_LHE.write(
"DECAY 4 0.000000e+00\n")
169 f_madspin_LHE.write(
"#\n")
170 f_madspin_LHE.write(
"# PDG Width\n")
171 f_madspin_LHE.write(
"DECAY 5 0.000000e+00\n")
172 f_madspin_LHE.write(
"#\n")
173 if "ckm" in process.MadSpin_model
or "CKM" in process.MadSpin_model:
174 logger.info(
"Top decay modes not set in param_card for model {} - BRs will be calculated at LO by MadSpin".
format(process.MadSpin_model))
175 f_madspin_LHE.write(
"# PDG Width\n")
176 f_madspin_LHE.write(
"DECAY 6 {:e}\n".
format(process.width_t))
177 f_madspin_LHE.write(
"#\n")
180 BR_t_to_any = process.BR_t_to_Wb + process.BR_t_to_Ws + process.BR_t_to_Wd
181 f_madspin_LHE.write(
"# PDG Width\n")
182 f_madspin_LHE.write(
"DECAY 6 {:e}\n".
format(process.width_t))
183 f_madspin_LHE.write(
"# BR NDA ID1 ID2 ...\n")
184 f_madspin_LHE.write(
" {:e} 2 5 24 # 1.32\n".
format(process.BR_t_to_Wb / BR_t_to_any))
185 f_madspin_LHE.write(
" {:e} 2 3 24 # 1.32\n".
format(process.BR_t_to_Ws / BR_t_to_any))
186 f_madspin_LHE.write(
" {:e} 2 1 24 # 1.32\n".
format(process.BR_t_to_Wd / BR_t_to_any))
187 f_madspin_LHE.write(
"#\n")
188 f_madspin_LHE.write(
"# PDG Width\n")
189 f_madspin_LHE.write(
"DECAY 11 0.000000e+00\n")
190 f_madspin_LHE.write(
"#\n")
191 f_madspin_LHE.write(
"# PDG Width\n")
192 f_madspin_LHE.write(
"DECAY 12 0.000000e+00\n")
193 f_madspin_LHE.write(
"#\n")
194 f_madspin_LHE.write(
"# PDG Width\n")
195 f_madspin_LHE.write(
"DECAY 13 0.000000e+00\n")
196 f_madspin_LHE.write(
"#\n")
197 f_madspin_LHE.write(
"# PDG Width\n")
198 f_madspin_LHE.write(
"DECAY 14 0.000000e+00\n")
199 f_madspin_LHE.write(
"#\n")
200 f_madspin_LHE.write(
"# PDG Width\n")
201 f_madspin_LHE.write(
"DECAY 15 0.000000e+00\n")
202 f_madspin_LHE.write(
"#\n")
203 f_madspin_LHE.write(
"# PDG Width\n")
204 f_madspin_LHE.write(
"DECAY 16 0.000000e+00\n")
205 f_madspin_LHE.write(
"#\n")
206 f_madspin_LHE.write(
"# PDG Width\n")
207 f_madspin_LHE.write(
"DECAY 21 0.000000e+00\n")
208 f_madspin_LHE.write(
"#\n")
209 f_madspin_LHE.write(
"# PDG Width\n")
210 f_madspin_LHE.write(
"DECAY 22 0.000000e+00\n")
211 f_madspin_LHE.write(
"#\n")
212 f_madspin_LHE.write(
"# PDG Width\n")
213 f_madspin_LHE.write(
"DECAY 23 {:e}\n".
format(process.width_Z))
214 f_madspin_LHE.write(
"#\n")
215 if "ckm" in process.MadSpin_model
or "CKM" in process.MadSpin_model:
216 logger.info(
"W decay modes not set in param_card for model {} - BRs will be calculated at LO by MadSpin".
format(process.MadSpin_model))
217 f_madspin_LHE.write(
"# PDG Width\n")
218 f_madspin_LHE.write(
"DECAY 24 {}\n".
format(process.width_W))
219 f_madspin_LHE.write(
"#\n")
222 BR_W_to_any = process.BR_W_to_hadrons + process.BR_W_to_leptons
223 f_madspin_LHE.write(
"# PDG Width\n")
224 f_madspin_LHE.write(
"DECAY 24 {}\n".
format(process.width_W))
225 f_madspin_LHE.write(
"# BR NDA ID1 ID2 ...\n")
226 f_madspin_LHE.write(
" {:e} 2 -1 2\n".
format((process.BR_W_to_hadrons / 2.) / BR_W_to_any))
227 f_madspin_LHE.write(
" {:e} 2 -3 4\n".
format((process.BR_W_to_hadrons / 2.) / BR_W_to_any))
228 f_madspin_LHE.write(
" {:e} 2 -11 12\n".
format((process.BR_W_to_leptons / 3.) / BR_W_to_any))
229 f_madspin_LHE.write(
" {:e} 2 -13 14\n".
format((process.BR_W_to_leptons / 3.) / BR_W_to_any))
230 f_madspin_LHE.write(
" {:e} 2 -15 16\n".
format((process.BR_W_to_leptons / 3.) / BR_W_to_any))
231 f_madspin_LHE.write(
"#\n")
232 f_madspin_LHE.write(
"# PDG Width\n")
233 f_madspin_LHE.write(
"DECAY 25 {:e}\n".
format(process.width_H))
234 f_madspin_LHE.write(
"#\n")
235 f_madspin_LHE.write(
"# PDG Width\n")
236 f_madspin_LHE.write(
"DECAY 82 0.000000e+00\n")
237 f_madspin_LHE.write(
"</slha>\n")
238 f_madspin_LHE.write(
"{}\n".
format(header_contents))
239 f_madspin_LHE.write(
"</header>\n")
240 f_madspin_LHE.write(
"{}\n".
format(post_header))
241 for event
in LHE.event_iterator(input_LHE_events):
242 f_madspin_LHE.write(event)
243 f_madspin_LHE.write(postamble)
246 shutil.move(input_LHE_events,
"{}.undecayed".
format(input_LHE_events))
247 shutil.move(
"madspin_LHE_input.lhe", input_LHE_events)
250 with open(
"madspin_runcard.txt",
"w")
as f_madspin_runcard:
251 f_madspin_runcard.write(
"import {}\n".
format(input_LHE_events))
252 f_madspin_runcard.write(
"set spinmode {}\n".
format(process.MadSpin_mode))
253 if process.MadSpin_max_weight_ps_point>0:
254 f_madspin_runcard.write(
"set max_weight_ps_point {}\n".
format(
int(process.MadSpin_max_weight_ps_point)))
255 if process.MadSpin_Nevents_for_max_weight>0:
256 f_madspin_runcard.write(
"set Nevents_for_max_weight {}\n".
format(
int(process.MadSpin_Nevents_for_max_weight)))
257 for decay
in process.MadSpin_decays:
258 f_madspin_runcard.write(
"{0}\n".
format(decay))
259 f_madspin_runcard.write(
"launch\n")
260 f_madspin_runcard.write(
"quit\n")
263 @
timed(
"MadSpin executable")
265 """! Run MadSpin executable.
267 @author James Robinson <james.robinson@cern.ch>
269 if not os.path.isfile(executable):
270 raise OSError(
"MadSpin executable {} not found!".
format(executable))
271 logger.info(
"MadSpin executable: {}".
format(executable))
272 with open(
"madspin_runcard.txt",
"r")
as runcard_input:
273 processes = [SingleProcessThread([
'python', executable], stdin=runcard_input, ignore_output=[
"INFO:",
"MadSpin>"],
274 error_output=[
"Command \"launch\" interrupted with error:",
"MadSpinError"],
275 info_output=[
"generating the production square matrix element"],
276 warning_output=[
"stty: standard input: Invalid argument",
277 "stty: standard input: Inappropriate ioctl for device",
278 "no version information available",
279 "CRITICAL: Branching ratio larger than one for"])]
280 manager = ProcessManager(processes)
281 while manager.monitor():
286 """! Prepare MadSpin output.
288 @author James Robinson <james.robinson@cern.ch>
290 logger.info(
"Preparing MadSpin output")
293 subprocess.call(
"gunzip pwgevents_decayed.lhe.gz 2> /dev/null", shell=
True)
294 shutil.move(
"pwgevents_decayed.lhe",
"{}.decayed".
format(input_LHE_events))
297 powheg_LHE, madspin_LHE =
"{}.undecayed".
format(input_LHE_events),
"{}.decayed".
format(input_LHE_events)
298 pre_header =
"\n".
join([elem
for elem
in [LHE.opening_tag(powheg_LHE), LHE.comment_block(powheg_LHE)]
if elem])
299 header = LHE.header_block(madspin_LHE)
300 post_header = LHE.init_block(powheg_LHE)
301 postamble = LHE.postamble(powheg_LHE)
304 with open(
"pwgevents_with_MadSpin.lhe",
"w")
as f_combined_LHE:
305 f_combined_LHE.write(
"{}\n".
format(pre_header))
306 f_combined_LHE.write(
"{}\n".
format(header))
307 f_combined_LHE.write(
"{}\n".
format(post_header))
308 for event
in LHE.event_iterator(
"{}.decayed".
format(input_LHE_events)):
309 f_combined_LHE.write(event)
310 f_combined_LHE.write(postamble)
313 if os.path.isfile(input_LHE_events):
314 shutil.move(input_LHE_events,
"{}.undecayed_ready_for_madspin".
format(input_LHE_events))
315 shutil.move(
"pwgevents_with_MadSpin.lhe", input_LHE_events)
317 for LHE_tarball
in list(
set(glob.glob(
"*lhe*.gz") + glob.glob(
"*events*.gz"))):
318 logger.info(
"Cleaning up MadSpin tarball: {}".
format(LHE_tarball))
319 os.remove(LHE_tarball)