~maddevelopers/mg5amcnlo/WWW5_caching

« back to all changes in this revision

Viewing changes to users/mardelcourt/PROC_407857/PROC_407857/bin/internal/combine_runs.py

  • Committer: John Doe
  • Date: 2013-03-25 20:27:02 UTC
  • Revision ID: john.doe@gmail.com-20130325202702-5sk3t1r8h33ca4p4
first clean version

Show diffs side-by-side

added added

removed removed

Lines of Context:
 
1
################################################################################
 
2
#
 
3
# Copyright (c) 2012 The MadGraph Development team and Contributors
 
4
#
 
5
# This file is a part of the MadGraph 5 project, an application which 
 
6
# automatically generates Feynman diagrams and matrix elements for arbitrary
 
7
# high-energy processes in the Standard Model and beyond.
 
8
#
 
9
# It is subject to the MadGraph license which should accompany this 
 
10
# distribution.
 
11
#
 
12
# For more information, please visit: http://madgraph.phys.ucl.ac.be
 
13
#
 
14
################################################################################
 
15
"""Program to combine results from channels that have been
 
16
     split into multiple jobs. Multi-job channels are identified
 
17
     by local file mjobs.dat in the channel directory.
 
18
"""
 
19
from __future__ import division
 
20
import math
 
21
import os
 
22
import re
 
23
import logging
 
24
 
 
25
try:
 
26
    import madgraph.various.sum_html as sum_html
 
27
    from madgraph import InvalidCmd, MadGraph5Error, MG5DIR
 
28
except ImportError:
 
29
    import internal.sum_html as sum_html
 
30
    from internal import InvalidCmd, MadGraph5Error
 
31
    
 
32
logger = logging.getLogger('madevent.combine_run') # -> stdout
 
33
 
 
34
#usefull shortcut
 
35
pjoin = os.path.join
 
36
 
 
37
   
 
38
def get_inc_file(path):
 
39
    """read the information of fortran inc files and returns
 
40
       the definition in a dictionary format.
 
41
       This catch PARAMETER (NAME = VALUE)"""
 
42
       
 
43
    pat = re.compile(r'''PARAMETER\s*\((?P<name>[_\w]*)\s*=\s*(?P<value>[\+\-\ded]*)\)''',
 
44
                     re.I)
 
45
        
 
46
    out = {}   
 
47
    for name, value in pat.findall(open(path).read()):
 
48
        orig_value = str(value)
 
49
        try:
 
50
            out[name.lower()] = float(value.replace('d','e'))
 
51
        except ValueError:
 
52
            out[name] = orig_value
 
53
    return out
 
54
 
 
55
class CombineRuns(object):
 
56
    
 
57
    def __init__(self, me_dir, subproc=None):
 
58
        
 
59
        self.me_dir = me_dir
 
60
        
 
61
        if not subproc:
 
62
            subproc = [l.strip() for l in open(pjoin(self.me_dir,'SubProcesses', 
 
63
                                                                 'subproc.mg'))]
 
64
        self.subproc = subproc
 
65
        maxpart = get_inc_file(pjoin(me_dir, 'Source', 'maxparticles.inc'))
 
66
        self.maxparticles = maxpart['max_particles']
 
67
    
 
68
    
 
69
        for procname in self.subproc:
 
70
            path = pjoin(self.me_dir,'SubProcesses', procname)
 
71
            channels = self.get_channels(path)
 
72
            for channel in channels:
 
73
                self.sum_multichannel(channel)
 
74
    
 
75
    def sum_multichannel(self, channel):
 
76
        """Looks in channel to see if there are multiple runs that
 
77
        need to be combined. If so combines them into single run"""
 
78
       
 
79
        alphabet = "abcdefghijklmnopqrstuvwxyz"
 
80
 
 
81
        if os.path.exists(pjoin(channel, 'multijob.dat')):
 
82
            njobs = int(open(pjoin(channel, 'multijob.dat')).read())
 
83
        else:
 
84
            return
 
85
        results = sum_html.Combine_results(channel)
 
86
        if njobs:
 
87
            logger.debug('find %s multijob in %s' % (njobs, channel))
 
88
        else:
 
89
            return
 
90
        for i in range(njobs):
 
91
            if channel.endswith(os.path.pathsep):
 
92
                path = channel[:-1] + alphabet[i] 
 
93
            else:
 
94
                path = channel + alphabet[i] 
 
95
            results.add_results(name=alphabet[i], 
 
96
                                filepath=pjoin(path, 'results.dat'))
 
97
        
 
98
        results.compute_average()
 
99
        if results.xsec:
 
100
            results.write_results_dat(pjoin(channel, 'results.dat'))
 
101
        else:
 
102
            return
 
103
        ### Adding information in the log file
 
104
        fsock = open(pjoin(channel, 'log.txt'), 'a')
 
105
        fsock.write('--------------------- Multi run with %s jobs. ---------------------\n'
 
106
                    % njobs)
 
107
        for r in results:
 
108
            fsock.write('job %s : %s %s %s\n' % (r.name, r.xsec, r.axsec, r.nunwgt))  
 
109
            
 
110
        #Now read in all of the events and write them
 
111
        #back out with the appropriate scaled weight
 
112
        fsock = open(pjoin(channel, 'events.lhe'), 'w')
 
113
        wgt = results.xsec / results.nunwgt
 
114
        tot_nevents=0
 
115
        for result in results:  
 
116
            i = result.name
 
117
            if channel.endswith(os.path.pathsep):
 
118
                path = channel[:-1] + i 
 
119
            else:
 
120
                path = channel + i 
 
121
            nw = self.copy_events(fsock, pjoin(path,'events.lhe'), wgt)
 
122
            #tot_events += nw
 
123
        #logger.debug("Combined %s events to %s " % (tot_events, channel))
 
124
 
 
125
 
 
126
    def copy_events(self, fsock, input, new_wgt):
 
127
        """ Copy events from separate runs into one file w/ appropriate wgts"""
 
128
        
 
129
        def get_fortran_str(nb):
 
130
            data = '%E' % nb
 
131
            nb, power = data.split('E')
 
132
            nb = float(nb) /10
 
133
            power = int(power) + 1
 
134
            return '%.7fE%+03i' %(nb,power)
 
135
        new_wgt = get_fortran_str(new_wgt)
 
136
        old_line = ""
 
137
        for line in open(input):
 
138
            if old_line.startswith("<event>"):
 
139
                data = line.split()
 
140
                if not len(data) == 6:
 
141
                    raise MadGraph5Error, "Line after <event> should have 6 entries"
 
142
                if float(data[2]) > 0:
 
143
                    sign = ''
 
144
                else:
 
145
                    sign = '-'
 
146
                line= ' %s  %s%s  %s\n' % ('   '.join(data[:2]), sign,
 
147
                                           new_wgt, '  '.join(data[3:]))
 
148
            fsock.write(line)
 
149
            old_line = line
 
150
           
 
151
    def get_channels(self, proc_path):
 
152
        """Opens file symfact.dat to determine all channels"""
 
153
        sympath = os.path.join(proc_path, 'symfact.dat')
 
154
        
 
155
        #ncode is number of digits needed for the bw coding
 
156
        
 
157
        ncode = int(math.log10(3)*(self.maxparticles-3))+1
 
158
        channels = []
 
159
        for line in open(sympath):
 
160
            try:
 
161
                xi, j = line.split()
 
162
            except Exception:
 
163
                break
 
164
            xi, j  = float(xi), int(j)
 
165
            
 
166
            if j > 0:
 
167
                k = int(xi) 
 
168
                npos = int(math.log10(k))+1
 
169
                #Write with correct number of digits
 
170
                if xi == k:
 
171
                    dirname = 'G%i' % k
 
172
                else:
 
173
                    dirname = 'G%.{0}f'.format(ncode) % xi
 
174
                channels.append(os.path.join(proc_path,dirname))
 
175
        return channels
 
176
    
 
177
        
 
178
        
 
179