ATLAS Offline Software
Hto4lPowhegMerge.py
Go to the documentation of this file.
1 # Copyright (C) 2002-2022 CERN for the benefit of the ATLAS collaboration
2 
3 import os, subprocess, tarfile
4 from AthenaCommon import Logging
5 from PowhegControl.utility import LHE
6 # from ...utility import LHE
7 # from xml.etree import ElementTree
8 import shutil
9 
10 
14 
15  __run_directory = os.environ['PATH']
16 
17 
18  __logger = Logging.logging.getLogger('Hto4lPowhegMerger')
19 
20 
21  _merger_executable = 'mergeHto4l4f.exe'
22 
23  def __init__( self, runArgs=None, opts=None ) :
24 
25 
26  self.__output_events_file_name = 'Hto4lPowhegMergedOTF._1.events'
27  # self.__output_events_file_name = 'events.lhe'
28 
29 
30 
31  self.__output_hto4l4e_file_name = 'Hto4lOTF4e._1.events'
32  self.__output_hto4l4mu_file_name = 'Hto4lOTF4mu._1.events'
33  self.__output_hto4l2e2mu_file_name = 'Hto4lOTF2e2mu._1.events'
34  self.__random_seed = 0
35 
36 
37  if runArgs is None :
38  self.logger.warning( 'No run arguments found! Using defaults.' )
39  else :
40  # Read values from runArgs
41  # if hasattr(runArgs,'outputTXTFile') :
42  # self.logger.info( 'outputTXTFile {0} {1}'.format (runArgs.outputTXTFile, type(runArgs.outputTXTFile) ) )
43  # self.__output_events_file_name = runArgs.outputTXTFile
44  # self.output_events_file_name = runArgs.outputTXTFile.split('.tar.gz')[0]+'.events'
45  # Set inputGeneratorFile to match output events file. Otherwise Generate_trf check will fail.
46 
47  runArgs.inputGeneratorFile = self.output_events_file_name
48  # runArgs.inputGeneratorFile = "Hto4lOTF4mu._1.events"
49 
50 
51  def setUpInput(self) :
52 
53  self.logger.info( 'entering setUpInput')
54 
55  # we require four input files, 1 for 4e, 1 for 4mu, and 2 for 2e2mu. The two for 2e2mu will be merged
56  myinputfiles = self.input_powheg_file_name
57  genInputFiles = myinputfiles.split(',')
58  inputsPath = os.path.dirname(genInputFiles[0])
59  numberOfFiles = len(genInputFiles)
60  # if there is a single file, make a symlink. If multiple files, merge them into one output eventsFile
61  if numberOfFiles > 0:
62  allFiles = []
63  for file in genInputFiles:
64  file_1 = file
65  # untar as needed
66  if tarfile.is_tarfile(file):
67  tar = tarfile.open(file)
68  tar.extractall(inputsPath)
69  tar.close()
70  file_1 = file.replace("tar.gz.1","events")
71  self.logger.info( 'Extracted tar file, and renaming {0} to {1}'.format ( file, file_1 ) )
72  pass
73 
74  # The only input format where merging is permitted is LHE
75  allFiles.append(file_1)
76  with open(file_1, 'r') as f:
77  first_line = f.readline()
78  self.logger.info( 'first_line {0}'.format ( first_line ) )
79  if(not ("LesHouche" in first_line)):
80  raise RuntimeError("%s is NOT a LesHouche file" % file)
81  pass
82  pass
83 
84  self.logger.info( 'Found files: nf {0}'.format( len(allFiles) ) )
85  for file in allFiles:
86  self.logger.info( ' {0}'.format( file ) )
87  pass
88 
89  # Merge files for 4e, 4mu and 2e2mu
90  self.__input_powheg_to_hto4l_file_name = "merged_hto4l_events.lhe"
92  self.logger.info( 'Merged: {0} to {1}'.format ( allFiles, self.input_powheg_to_hto4l_file_name ) )
93  else:
94  self.__input_powheg_to_hto4l_file_name = myinputfiles
95  self.logger.info( 'Using single input file: {0}'.format ( self.input_powheg_to_hto4l_file_name ) )
96 
97 
98  def merge(self) :
99 
100  self.logger.info( 'Starting merge' )
101 
102  hto4lLHE4e = self.output_hto4l4e_file_name
103  hto4lLHE4mu = self.output_hto4l4mu_file_name
104  hto4lLHE2e2mu = self.output_hto4l2e2mu_file_name
105 # powhegLHE_output = self.output_events_file_name
106 
107  self.logger.info( 'Input files: {0} {1} {2}'.format( hto4lLHE4e, hto4lLHE4mu, hto4lLHE2e2mu ) )
108 
109  allFiles = []
110  allFiles += [ hto4lLHE4e ]
111  allFiles += [ hto4lLHE4mu ]
112  allFiles += [ hto4lLHE2e2mu ]
113  self.merge_lhe_files(allFiles, self.output_events_file_name)
114 
115  self.logger.info( 'Merged to output file: {0}'.format( self.output_events_file_name ) )
116 
117  # Now fix the output file for negative weights
119 
120  self.logger.info( 'Fixed weights in file: {}'.format( self.output_events_file_name ) )
121 
122 
123  return
124 
125  def reweight_for_negative_weights(self, powheg_LHE_output):
126  """Post-process the LHE file to update the weights for a negative Hto4l weight.
127  We do two things:
128  1) check if the event weight (XWGTUP) is -1 (can be +/-1). If so, remember this so that the
129  weights in the <rwgt> block can be corrected by multiplying by -1
130  2) After correcting the weights in the <rwgt> block, overwrite XWGTUP by the weight at id=0.
131  Not completely sure if this step is necessary, but is done for Prophecy4f as well.
132  Note: we do NOT save the original XWGTUP as an extra weight (i.e. +/-1). Normally, the sign of
133  any weight should indicate this. The one exception is when BOTH PowHeg and Hto4l has a negative
134  weight and so all weights are positive. This is pretty rare.
135 
136  #@param powheg_LHE_output Name of LHE file produced by merge of Hto4l files.
137 
138  @author RD Schaffer <r.d.schaffer@cern.ch>
139  #"""
140 
141  self.logger.info( 'Starting reweight_for_negative_weights' )
142 
143  # Get opening and closing strings
144  preamble = LHE.preamble(powheg_LHE_output)
145  postamble = LHE.postamble(powheg_LHE_output)
146 
147  # open new file for updates
148  powheg_LHE_updated = "{}.updated".format(powheg_LHE_output)
149  with open(powheg_LHE_updated, "w") as f_output:
150  f_output.write("{}".format(preamble))
151  # for each event, check for an -1 Hto4l weight in XWGTUP. If it is -1, correct the <rwgt>
152  # weights (-1 * wgt), and save the id=0 weight in XWGTUP.
153  for input_event in LHE.event_iterator(powheg_LHE_output):
154 
155  output_event = ""
156  has_neg_wgt = False
157 
158  # First check for negative weights
159  input_lines = input_event.splitlines(True)
160 
161  try: # interpret 2nd line as a general event info line
162  NUP, IDPRUP, XWGTUP, SCALUP, AQEDUP, AQCDUP = input_lines[1].split()
163  has_neg_wgt = float(XWGTUP) < 0.
164 
165  self.logger.info( 'has_neg_wgt {}, XWGTUP {}'.format(has_neg_wgt, XWGTUP))
166 
167  except ValueError: # this is not a general event info line
168  self.logger.warning( 'could not get first line of event')
169  pass
170 
171  # Now correct weights
172  output_event = ""
173  if has_neg_wgt:
174  for input_line in input_event.splitlines(True):
175  self.logger.info( 'input line {}'.format(input_line))
176  output_line = None
177  if input_line.find("<wgt") >= 0:
178  weights = input_line.split()
179  self.logger.info( 'weightsinput {}'.format(weights))
180  try:
181  wgt_tag, wid, wgt, wgt_tag_end = input_line.split()
182  # We 'multiply' the two weights - if wgt is negative, it now becomes positive.
183  # Otherwise, set it to negative
184  if wgt[0] == '-':
185  output_line = "{} {} {} {} \n".format(wgt_tag, wid, wgt[1:], wgt_tag_end)
186  else:
187  output_line = "{} {} -{} {} \n".format(wgt_tag, wid, wgt, wgt_tag_end)
188  except Exception:
189  pass
190  self.logger.info( 'output line {}'.format(output_line))
191 
192  # save event lines
193  output_event += output_line if output_line is not None else input_line
194  else:
195  output_event = input_event
196 
197  # Overwrite (XWGTUP)
198  output_event_updated = LHE.update_XWGTUP_with_reweighted_nominal(output_event)
199  f_output.write(output_event_updated)
200 
201  f_output.write(postamble)
202  # f_output.close()
203 
204  # Make a backup of the original events
205  shutil.move(powheg_LHE_output, "{}.lhe_nominal_weight_updater_backup".format(powheg_LHE_output))
206  shutil.move(powheg_LHE_updated, powheg_LHE_output)
207 
208  self.logger.info( 'wrote file {}'.format(powheg_LHE_output))
209 
210 
211  @property
213  return self.__output_events_file_name
214 
215 
216  @property
218  return self.__input_powheg_file_name
219 
220 
221  @input_powheg_file_name.setter
222  def input_powheg_file_name(self, value) :
224 
225 
226  @property
229 
230 
231  @property
233  return self.__output_hto4l4e_file_name
234 
235 
236  @property
238  return self.__output_hto4l4mu_file_name
239 
240 
241  @property
244 
245 
246  @output_hto4l4e_file_name.setter
247  def output_hto4l4e_file_name(self, value) :
248  self.__output_hto4l4e_file_name = value
249 
250 
251  @output_hto4l4mu_file_name.setter
252  def output_hto4l4mu_file_name(self, value) :
253  self.__output_hto4l4mu_file_name = value
254 
255 
256  @output_hto4l2e2mu_file_name.setter
257  def output_hto4l2e2mu_file_name(self, value) :
258  self.__output_hto4l2e2mu_file_name = value
259 
260  @property
261  def random_seed(self) :
262  return self.__random_seed
263 
264  @random_seed.setter
265  def random_seed(self, value) :
266  self.__random_seed = value
267 
268 
269  @property
270  def logger(self) :
271  return self.__logger
272 
273  # This function merges a list of input LHE file to make one outputFile. The header is taken from the first
274  # file, but the number of events is updated to equal the total number of events in all the input files
275  def merge_lhe_files(self, listOfFiles, outputFile):
276  if(os.path.exists(outputFile)):
277  self.logger.info( 'outputFile {0} already exists. Will rename to {1}.OLD'.format ( outputFile, outputFile ) )
278  os.rename(outputFile,outputFile+".OLD")
279  output = open(outputFile,'w')
280  holdHeader = ""
281  nevents=0
282  for file in listOfFiles:
283  cmd = "grep /event "+file+" | wc -l"
284  nevents+=int(subprocess.check_output(cmd,stderr=subprocess.STDOUT,shell=True))
285 
286  for file in listOfFiles:
287  inHeader = True
288  header = ""
289  self.logger.info( '*** Starting to merge file {0}'.format ( file ) )
290  for line in open(file,"r"):
291 
296  if("<event" in line and inHeader):
297  inHeader = False
298  if(len(holdHeader)<1):
299  holdHeader = header
300  output.write(header)
301  output.write(line)
302 
304  elif(not inHeader and not ("</LesHouchesEvents>" in line)):
305  output.write(line)
306  if(inHeader):
307 
308  if("nevents" in line):
309 
310  tmp = line.split("=")
311  line = line.replace(tmp[0],str(nevents))
312  elif("numevts" in line):
313 
314  tmp = line.split(" ")
315  nnn = str(nevents)
316  line = line.replace(tmp[1],nnn)
317  header+=line
318  output.write("</LesHouchesEvents>\n")
319  output.close()
320 
python.Hto4lPowhegMerge.Hto4lPowhegMerge.reweight_for_negative_weights
def reweight_for_negative_weights(self, powheg_LHE_output)
Definition: Hto4lPowhegMerge.py:125
python.Hto4lPowhegMerge.Hto4lPowhegMerge.__output_hto4l2e2mu_file_name
__output_hto4l2e2mu_file_name
Definition: Hto4lPowhegMerge.py:33
vtune_athena.format
format
Definition: vtune_athena.py:14
python.Hto4lPowhegMerge.Hto4lPowhegMerge.merge
def merge(self)
Initialise runcard with generic options.
Definition: Hto4lPowhegMerge.py:98
python.Hto4lPowhegMerge.Hto4lPowhegMerge.output_hto4l4e_file_name
def output_hto4l4e_file_name(self)
Get output Hto4l file name for 4e.
Definition: Hto4lPowhegMerge.py:232
python.Hto4lPowhegMerge.Hto4lPowhegMerge.__input_powheg_file_name
__input_powheg_file_name
Definition: Hto4lPowhegMerge.py:223
python.Hto4lPowhegMerge.Hto4lPowhegMerge.input_powheg_file_name
def input_powheg_file_name(self)
Get input Powheg file name.
Definition: Hto4lPowhegMerge.py:217
python.Hto4lPowhegMerge.Hto4lPowhegMerge.output_hto4l2e2mu_file_name
def output_hto4l2e2mu_file_name(self)
Get output Hto4l file name for 2e2mu.
Definition: Hto4lPowhegMerge.py:242
python.Hto4lPowhegMerge.Hto4lPowhegMerge.__output_hto4l4e_file_name
__output_hto4l4e_file_name
Using default output names from PowhegConfig_base and Hto4lConfig.
Definition: Hto4lPowhegMerge.py:31
python.Hto4lPowhegMerge.Hto4lPowhegMerge.merge_lhe_files
def merge_lhe_files(self, listOfFiles, outputFile)
Definition: Hto4lPowhegMerge.py:275
python.Hto4lPowhegMerge.Hto4lPowhegMerge.setUpInput
def setUpInput(self)
Initialise runcard with generic options.
Definition: Hto4lPowhegMerge.py:51
python.Hto4lPowhegMerge.Hto4lPowhegMerge.input_powheg_to_hto4l_file_name
def input_powheg_to_hto4l_file_name(self)
Get input Hto4l file name.
Definition: Hto4lPowhegMerge.py:227
python.Hto4lPowhegMerge.Hto4lPowhegMerge.__random_seed
__random_seed
Definition: Hto4lPowhegMerge.py:34
python.Hto4lPowhegMerge.Hto4lPowhegMerge.__init__
def __init__(self, runArgs=None, opts=None)
Definition: Hto4lPowhegMerge.py:23
python.Hto4lPowhegMerge.Hto4lPowhegMerge
Base class for configurable objects in the jobOptions.
Definition: Hto4lPowhegMerge.py:13
python.Hto4lPowhegMerge.Hto4lPowhegMerge.output_events_file_name
def output_events_file_name(self)
Get output file name.
Definition: Hto4lPowhegMerge.py:212
python.Hto4lPowhegMerge.Hto4lPowhegMerge.logger
def logger(self)
Get handle to logger.
Definition: Hto4lPowhegMerge.py:270
python.Hto4lPowhegMerge.Hto4lPowhegMerge.__logger
__logger
Setup athena-compatible logger.
Definition: Hto4lPowhegMerge.py:18
Trk::open
@ open
Definition: BinningType.h:40
python.Hto4lPowhegMerge.Hto4lPowhegMerge.output_hto4l4mu_file_name
def output_hto4l4mu_file_name(self)
Get output Hto4l file name for 4mu.
Definition: Hto4lPowhegMerge.py:237
python.CaloAddPedShiftConfig.int
int
Definition: CaloAddPedShiftConfig.py:45
if
if(febId1==febId2)
Definition: LArRodBlockPhysicsV0.cxx:567
pickleTool.object
object
Definition: pickleTool.py:29
str
Definition: BTagTrackIpAccessor.cxx:11
python.Hto4lPowhegMerge.Hto4lPowhegMerge.__input_powheg_to_hto4l_file_name
__input_powheg_to_hto4l_file_name
Definition: Hto4lPowhegMerge.py:90
python.Hto4lPowhegMerge.Hto4lPowhegMerge.random_seed
def random_seed(self)
Definition: Hto4lPowhegMerge.py:261
python.ParticleTypeUtil.info
def info
Definition: ParticleTypeUtil.py:87
python.Hto4lPowhegMerge.Hto4lPowhegMerge.__output_events_file_name
__output_events_file_name
This needs to be set so that Generate_trf finds an appropriate file format for showering.
Definition: Hto4lPowhegMerge.py:26
Trk::split
@ split
Definition: LayerMaterialProperties.h:38
python.Hto4lPowhegMerge.Hto4lPowhegMerge.__output_hto4l4mu_file_name
__output_hto4l4mu_file_name
Definition: Hto4lPowhegMerge.py:32
python.LArMinBiasAlgConfig.float
float
Definition: LArMinBiasAlgConfig.py:65