~maddevelopers/mg5amcnlo/2.9.4

« back to all changes in this revision

Viewing changes to madgraph/various/combine_plots.py

pass to v2.0.0

Show diffs side-by-side

added added

removed removed

Lines of Context:
 
1
#!/usr/bin/env python
 
2
################################################################################
 
3
#
 
4
# Copyright (c) 2013 The MadGraph5_aMC@NLO Development team and Contributors
 
5
#
 
6
# This file is a part of the MadGraph5_aMC@NLO project, an application which 
 
7
# automatically generates Feynman diagrams and matrix elements for arbitrary
 
8
# high-energy processes in the Standard Model and beyond.
 
9
#
 
10
# It is subject to the MadGraph5_aMC@NLO license which should accompany this 
 
11
# distribution.
 
12
#
 
13
# For more information, visit madgraph.phys.ucl.ac.be and amcatnlo.web.cern.ch
 
14
#
 
15
################################################################################
 
16
import glob
 
17
import os
 
18
import re
 
19
import subprocess
 
20
 
 
21
pjoin = os.path.join
 
22
 
 
23
class histograms:
 
24
    """ x and y values """
 
25
    def __init__(self):
 
26
        pass
 
27
 
 
28
 
 
29
class one_plot:
 
30
    """ All information about one plot """
 
31
    def __init__(self):
 
32
        self.histo = {}
 
33
        self.top = ""
 
34
        self.end = ""
 
35
        self.max = 0
 
36
        self.max_file =''
 
37
        
 
38
    def add_comment(self, top_comment, end_comment=""):
 
39
        
 
40
        self.top = top_comment
 
41
        self.end = end_comment       
 
42
        
 
43
    
 
44
    def add_histo(self, string_values, tag):
 
45
        """Add the histogram in the data base (string format) """
 
46
        if tag in self.histo.values():
 
47
            print "Warning: skyping the histogram with tag " + tag
 
48
            print "         since it is already present in the data base"
 
49
            return
 
50
        self.histo[tag] = {}
 
51
        self.histo[tag]["values"] = string_values
 
52
        old_max = self.max
 
53
        for line in string_values.split('\n'):
 
54
            split = line.split()
 
55
            if len(split) ==3:
 
56
                self.max = max(self.max, float(line.split()[1]))
 
57
        if self.max != old_max:
 
58
            self.max_file = tag
 
59
            
 
60
    def get_histo(self, tag, norm):
 
61
        """return a string with the histogram values, and the normalization """
 
62
        if tag not in self.histo or self.histo[tag]["values"] == '':
 
63
            return ''
 
64
        histo = "SET ORDER X Y " + str(norm) + " \n"
 
65
        histo += self.histo[tag]["values"]
 
66
        return histo    
 
67
    
 
68
 
 
69
class load_data:
 
70
    """ Import all the plots from top drawer files """
 
71
    def __init__(self):
 
72
 
 
73
        self.plots = {}  # list of the tags associated with the plots
 
74
        self.normalization = 1.0
 
75
        self.order_plots = []
 
76
        self.cross = {}
 
77
 
 
78
    def import_data(self, file_name, file_nb, pretag=""):
 
79
        """ Read the file with the name file_name, and import the data for all the plots"""
 
80
 
 
81
        trappe = open(file_name, 'r')
 
82
        while 1:
 
83
            line = trappe.readline()
 
84
            if line == "": break
 
85
            if line.find("Cross-section") > -1:
 
86
                pos = line.find("is:")
 
87
                list = line[pos + 4:].split()
 
88
                self.cross1 = float(list[0])
 
89
            if line.find("NEW PLOT") > -1 or line.find("SET INTENSITY") > -1: 
 
90
                self.readplot(trappe, file_nb) 
 
91
 
 
92
    def print_plots(self, outputfile, file1, file2):
 
93
        """Write out histograms """
 
94
        trappe = open(outputfile, 'w')
 
95
        trappe.write("SET DEVICE POSTSCRIPT ORIENTATION=3 \n")
 
96
        trappe.write("set intensity 5 \n")
 
97
 
 
98
 
 
99
 
 
100
        for index, tag_plot in enumerate(self.order_plots):
 
101
            norm1 = 1 
 
102
            norm2 = 1
 
103
            if self.plots[tag_plot].max_file == file1:
 
104
                color1, color2 = 'WHITE','BLUE'
 
105
                histtype1, histtype2 = 'HIST SOLID\n', 'SET PATTERN .05 .07\n' 
 
106
            else:
 
107
                file1, file2 = file2, file1
 
108
                color1, color2 = 'BLUE', 'WHITE'
 
109
                histtype2, histtype1 = 'HIST SOLID\n', 'SET PATTERN .05 .07\n' 
 
110
            if self.plots[tag_plot].get_histo(file1, norm1) == '':
 
111
                continue
 
112
    
 
113
            trappe.write("NEW PLOT \n")
 
114
            top = self.plots[tag_plot].top
 
115
            if top.find('# particles') > -1:
 
116
                top = top.replace('LOG', 'LIN')
 
117
            trappe.write(top)
 
118
#            trappe.write("SET COLOR %s \n" % color2)
 
119
#            trappe.write("HIST PATTERNED %s \n" % color1)
 
120
 
 
121
                
 
122
            try:
 
123
                trappe.write(histtype1)                
 
124
                trappe.write(self.plots[tag_plot].get_histo(file1, norm1))
 
125
                trappe.write(histtype1)
 
126
                trappe.write("HIST PATTERNED %s \n" % color1)
 
127
#                trappe.write("SET COLOR %s \n" % color1)
 
128
            except:
 
129
                print "warning: cannot find histo in file 1 for tag " + tag_plot 
 
130
        
 
131
            try:
 
132
                trappe.write(self.plots[tag_plot].get_histo(file2, norm2))
 
133
                trappe.write(histtype2)
 
134
                trappe.write("HIST PATTERNED %s \n" % color2)
 
135
                trappe.write("SET COLOR WHITE \n")
 
136
            except Exception, error:
 
137
                print error
 
138
                print "warning: cannot find histo in file 2 for tag " + tag_plot 
 
139
                raise
 
140
            
 
141
            #trappe.write(self.plots[tag_plot].end)
 
142
    
 
143
            trappe.write("\n")
 
144
            trappe.write("\n")
 
145
        trappe.close()
 
146
 
 
147
    def readplot(self, trappe, file_nb):
 
148
        """ import the data of a single plot """
 
149
        top = ""
 
150
        end = ""
 
151
        histo = ""
 
152
        phase = 0
 
153
        plot_tag = ""
 
154
        newplot = 1
 
155
        while 1:
 
156
          line = trappe.readline()
 
157
          #print line
 
158
          if line == "": break
 
159
          if phase == 0:  # top of the histogram
 
160
             if line.find("TITLE TOP") > -1:
 
161
                index = line.find('''"''')
 
162
                if index < 0:
 
163
                  print "warning: unable to find the name of the plot in the title"
 
164
                  print "         skipping this plot (might not be a real plot)"
 
165
                  return
 
166
                else:
 
167
                  plot_tag = line[index + 1:]
 
168
                  plot_tag = plot_tag.replace('''"''', "")
 
169
                  plot_tag = plot_tag.replace("\n", "")
 
170
                  plot_tag = plot_tag.replace(" ", "")
 
171
             if line.find("SET LIMITS X") > -1:
 
172
                  tag = line.replace(".0 ", "")
 
173
                  tag = tag.replace(".00 ", "") 
 
174
                  tag = tag.replace(".000 ", "") 
 
175
                  tag = tag.replace(".0000 ", "") 
 
176
                  tag = tag.replace(".00000 ", "") 
 
177
                  tag = tag.replace(".0\n", "") 
 
178
                  tag = tag.replace(".00\n", "") 
 
179
                  tag = tag.replace(".000\n", "") 
 
180
                  tag = tag.replace(".0000\n", "") 
 
181
                  tag = tag.replace(".00000\n", "") 
 
182
                  tag = tag.replace(" ", "")
 
183
                  tag = tag.replace("\n", "")
 
184
                  tag = tag.replace("3.14160", "3.14159")
 
185
                  tag = tag.replace("9.42480", "9.42478")
 
186
                  tag = tag.replace("4.18880", "4.18879")
 
187
                  if plot_tag != "":
 
188
                    plot_tag += tag
 
189
                    if plot_tag not in self.plots:
 
190
                      self.plots[plot_tag] = one_plot()
 
191
                      self.order_plots.append(plot_tag)
 
192
                    else:
 
193
                      newplot = 0
 
194
                  else:
 
195
                    print "warning: unusual format, unable to extract the tag of the plot"
 
196
                    print "         skipping this plot (might not be a real plot)"
 
197
                    return
 
198
                  
 
199
             if line.find("SET ORDER") > -1 and plot_tag != "":
 
200
                line = trappe.readline()
 
201
                phase = 1
 
202
             else:
 
203
                top += line  # store the line 
 
204
          if phase == 1:  # histogram values
 
205
            if line.find("PLOT") < 0 and line.find("SET PATTERN") < 0 and line.find("HIST") < 0:
 
206
              histo += line  # store the line
 
207
            else:
 
208
              if line.find("HISTO") > -1: 
 
209
                histo_tag = file_nb
 
210
                self.plots[plot_tag].add_histo(histo, file_nb)
 
211
    
 
212
                if plot_tag.find("Weights") > -1:
 
213
                  pos = histo.find("1.0100")
 
214
                  if pos > -1:
 
215
                     list = histo[pos:].split()
 
216
                     self.cross[file_nb] = float(list[1])
 
217
                histo = ""
 
218
                phase = 2
 
219
    
 
220
    
 
221
          if phase == 2:
 
222
             if line.find("NEW PLOT") > -1:
 
223
              if (newplot):
 
224
                self.plots[plot_tag].add_comment(top, end_comment=end)
 
225
              self.readplot(trappe, file_nb)
 
226
             else:
 
227
                 if line != "":end += line
 
228
 
 
229
        if plot_tag and plot_tag in self.plots:
 
230
            self.plots[plot_tag].add_comment(top, end_comment=end)
 
231
 
 
232
def merge_all_plots(path1, path2, outputpath='/tmp', td='../../td/td', MA=None):
 
233
    """take a MA4 output and merge all the plots present in the HTML output"""
 
234
    
 
235
    #find which plots correspond to what
 
236
    pattern = re.compile(r'''TITLE TOP\s+\"\s*([^\"]*)\s*\"''')
 
237
    all_plot1 = {}
 
238
    for filepath in glob.glob(pjoin(path1, 'ma_*.top')):
 
239
        filename = os.path.basename(filepath)
 
240
        text = open(filepath).read()
 
241
        try:
 
242
            title = pattern.search(text).groups()[0].strip()
 
243
        except AttributeError:
 
244
            continue
 
245
        all_plot1[title] = filename
 
246
        
 
247
    # find the plot to add:
 
248
    for filepath in glob.glob(pjoin(path2, 'ma_*.top')):
 
249
        filename = os.path.basename(filepath)
 
250
        text = open(filepath).read()
 
251
        try:
 
252
            title = pattern.search(text).groups()[0].strip()
 
253
        except AttributeError:
 
254
            continue
 
255
        if title not in all_plot1:
 
256
            continue
 
257
        my_data = load_data()
 
258
        my_data.import_data(pjoin(path1, all_plot1[title]), 1)
 
259
        my_data.import_data(pjoin(path2, filename), 2)
 
260
        my_data.print_plots(pjoin(outputpath, filename), 1, 2)
 
261
        if 'DYLD_LIBRARY_PATH' not in os.environ:
 
262
            os.environ['DYLD_LIBRARY_PATH'] =  os.path.dirname(td)
 
263
        elif os.path.dirname(td) not in os.environ['DYLD_LIBRARY_PATH']:
 
264
            os.environ['DYLD_LIBRARY_PATH'] = '%s:%s' %( os.environ['DYLD_LIBRARY_PATH'], os.path.dirname(td))
 
265
        devnull = open(os.devnull,'w')
 
266
        subprocess.call([td, filename], cwd=outputpath, stdout=devnull)
 
267
        devnull.close()
 
268
        if MA:
 
269
            subprocess.call([pjoin(MA, 'epstosmth'),"--gsopt=\'-r60x60 -dGraphicsAlphaBits=4\'", 
 
270
                         "--gsdev=jpeg", filename.replace('.top', '.ps')], cwd=outputpath)
 
271
        #os.system('%s %s' % (td, path))      
 
272
        
 
273
        
 
274
 
 
275
if __name__ == "__main__":
 
276
 
 
277
    merge_all_plots('PROC_EWdim6_1/HTML/run_12/plots_parton_tag_1/',
 
278
                    'PROC_EWdim6_1/HTML/run_11_rw_15/plots_parton_reweight_tag_1/')
 
279