ViewVC Help
View File | Revision Log | Show Annotations | Root Listing
root/cvsroot/UserCode/VHbb/python/evaluateMVA.py
(Generate patch)

Comparing UserCode/VHbb/python/evaluateMVA.py (file contents):
Revision 1.17 by peller, Tue Nov 27 09:00:27 2012 UTC vs.
Revision 1.18 by nmohr, Wed Nov 28 13:04:59 2012 UTC

# Line 2 | Line 2
2   import sys
3   import os
4   import ROOT
5 from ROOT import TFile
5   from array import array
6   from math import sqrt
7   from copy import copy
# Line 11 | Line 10 | import warnings
10   warnings.filterwarnings( action='ignore', category=RuntimeWarning, message='creating converter.*' )
11   from optparse import OptionParser
12   from BetterConfigParser import BetterConfigParser
14 from samplesclass import sample
15 from mvainfos import mvainfo
13   import pickle
17 from progbar import progbar
18 from printcolor import printc
14  
15   #CONFIGURE
16 + ROOT.gROOT.SetBatch(True)
17   print 'hello'
18   #load config
19   #os.mkdir(path+'/sys')
# Line 36 | Line 32 | parser.add_option("-S", "--samples", des
32   parser.add_option("-C", "--config", dest="config", default=[], action="append",
33                        help="configuration file")
34   (opts, args) = parser.parse_args(argv)
35 +
36 + from samplesclass import sample
37 + from mvainfos import mvainfo
38 + from progbar import progbar
39 + from printcolor import printc
40 +
41   if opts.config =="":
42          opts.config = "config"
43   config = BetterConfigParser()
# Line 49 | Line 51 | MVASubdir=config.get('Directories','MVAd
51   samplesinfo=config.get('Directories','samplesinfo')
52  
53   #systematics
52 systematics=config.get('systematics','systematics')
53 systematics=systematics.split(' ')
54  
55 #TreeVar Array
56 #MVA_Vars={}
57 #for systematic in systematics:
58 #    MVA_Vars[systematic]=config.get('treeVars',systematic)
59 #    MVA_Vars[systematic]=MVA_Vars[systematic].split(' ')
55  
61 ######################
62 #Evaluate multi: Must Have same treeVars!!!
63
64 #OUTpath=opts.outpath
65 #INpath=opts.inpath
56   INpath = config.get('Directories','MVAin')
57   OUTpath = config.get('Directories','MVAout')
58  
# Line 77 | Line 67 | namelist=namelistIN.split(',')
67   doinfo=bool(int(opts.update))
68  
69   MVAlist=arglist.split(',')
70 < MVAdir=config.get('Directories','vhbbpath')
70 >
71  
72   #CONFIG
73   #factory
74   factoryname=config.get('factory','factoryname')
75   #MVA
86 #MVAnames=[]
87 #for MVA in MVAlist:
88 #    print MVA
89 #    MVAnames.append(config.get(MVA,'MVAname'))
90 #print Wdir+'/weights/'+factoryname+'_'+MVAname+'.info'
91 #MVAinfofiles=[]
76   MVAinfos=[]
77 + MVAdir=config.get('Directories','vhbbpath')
78   for MVAname in MVAlist:
79      MVAinfofile = open(MVAdir+'/data/'+factoryname+'_'+MVAname+'.info','r')
80      MVAinfos.append(pickle.load(MVAinfofile))
81      MVAinfofile.close()
82      
98 treeVarSet=MVAinfos[0].varset
99 #variables
100 #TreeVar Array
101 MVA_Vars={}
102 for systematic in systematics:
103    MVA_Vars[systematic]=config.get(treeVarSet,systematic)
104    MVA_Vars[systematic]=MVA_Vars[systematic].split(' ')
105 #Spectators:
106 #spectators=config.get(treeVarSet,'spectators')
107 #spectators=spectators.split(' ')
108 #progbar quatsch
83   longe=40
84   #Workdir
85   workdir=ROOT.gDirectory.GetPath()
112 #os.mkdir(Apath+'/MVAout')
113
114 #Book TMVA readers: MVAlist=["MMCC_bla","CC5050_bla"]
115 readers=[]
116 for MVA in MVAlist:
117    readers.append(ROOT.TMVA.Reader("!Color:!Silent"))
118
119 #define variables and specatators
120 MVA_var_buffer = []
121 MVA_var_buffer4 = []
122 for i in range(len( MVA_Vars['Nominal'])):
123    MVA_var_buffer.append(array( 'f', [ 0 ] ))
124    for reader in readers:
125        reader.AddVariable( MVA_Vars['Nominal'][i],MVA_var_buffer[i])
126 #MVA_spectator_buffer = []
127 #for i in range(len(spectators)):
128 #    MVA_spectator_buffer.append(array( 'f', [ 0 ] ))
129 #    for reader in readers:
130 #        reader.AddSpectator(spectators[i],MVA_spectator_buffer[i])
131 #Load raeder
132 for i in range(0,len(readers)):
133    readers[i].BookMVA(MVAinfos[i].MVAname,MVAdir+'/data/'+MVAinfos[i].getweightfile())
134 #--> Now the MVA is booked
86  
87   #Apply samples
88   infofile = open(samplesinfo,'r')
89   Ainfo = pickle.load(infofile)
90   infofile.close()
91  
92 +
93 + class MvaEvaluater:
94 +    def __init__(self, config, MVAinfo):
95 +        self.varset = MVAinfo.varset
96 +        #Define reader
97 +        self.reader = ROOT.TMVA.Reader("!Color:!Silent")
98 +        MVAdir=config.get('Directories','vhbbpath')
99 +        self.systematics=config.get('systematics','systematics').split(' ')
100 +        self.MVA_Vars={}
101 +        self.MVAname = MVAinfo.MVAname
102 +        for systematic in self.systematics:
103 +            self.MVA_Vars[systematic]=config.get(self.varset,systematic)
104 +            self.MVA_Vars[systematic]=self.MVA_Vars[systematic].split(' ')
105 +        #define variables and specatators
106 +        self.MVA_var_buffer = []
107 +        for i in range(len( self.MVA_Vars['Nominal'])):
108 +            self.MVA_var_buffer.append(array( 'f', [ 0 ] ))
109 +            self.reader.AddVariable( self.MVA_Vars['Nominal'][i],self.MVA_var_buffer[i])
110 +        self.reader.BookMVA(MVAinfo.MVAname,MVAdir+'/data/'+MVAinfo.getweightfile())
111 +        #--> Now the MVA is booked
112 +
113 +    def setBranches(self,tree,job):
114 +        #Set formulas for all vars
115 +        self.MVA_formulas={}
116 +        for systematic in self.systematics:
117 +            if job.type == 'DATA' and not systematic == 'Nominal': continue
118 +            self.MVA_formulas[systematic]=[]
119 +            for j in range(len( self.MVA_Vars['Nominal'])):
120 +                self.MVA_formulas[systematic].append(ROOT.TTreeFormula("MVA_formula%s_%s"%(j,systematic),self.MVA_Vars[systematic][j],tree))
121 +
122 +    def evaluate(self,MVAbranches,job):
123 +        #Evaluate all vars and fill the branches
124 +        for systematic in self.systematics:
125 +            for j in range(len( self.MVA_Vars['Nominal'])):
126 +                if job.type == 'DATA' and not systematic == 'Nominal': continue
127 +                self.MVA_var_buffer[j][0] = self.MVA_formulas[systematic][j].EvalInstance()                
128 +            MVAbranches[self.systematics.index(systematic)] = self.reader.EvaluateMVA(self.MVAname)
129 +
130 +
131 + theMVAs = []
132 + for mva in MVAinfos:
133 +    theMVAs.append(MvaEvaluater(config,mva))
134 +
135 +
136   #eval
137   for job in Ainfo:
138      if eval(job.active):
139          if job.name in namelist:
140              #get trees:
141              print INpath+'/'+job.prefix+job.identifier+'.root'
142 <            input = TFile.Open(INpath+'/'+job.prefix+job.identifier+'.root','read')
142 >            input = ROOT.TFile.Open(INpath+'/'+job.prefix+job.identifier+'.root','read')
143              print OUTpath+'/'+job.prefix+job.identifier+'.root'
144 <            outfile = TFile.Open(OUTpath+'/'+job.prefix+job.identifier+'.root','recreate')
144 >            outfile = ROOT.TFile.Open(OUTpath+'/'+job.prefix+job.identifier+'.root','recreate')
145              input.cd()
146              obj = ROOT.TObject
147              for key in ROOT.gDirectory.GetListOfKeys():
# Line 162 | Line 157 | for job in Ainfo:
157              nEntries = tree.GetEntries()
158              outfile.cd()
159              newtree = tree.CloneTree(0)
160 <            #input.Close()
160 >            
161  
162 <            #MCs:
163 <            if job.type != 'DATA':
164 <                MVA_formulas={}
165 <                MVA_formulas4={}
166 <                for systematic in systematics:
167 <                    #print '\t\t - ' + systematic
168 <                    MVA_formulas[systematic]=[]
169 <                    MVA_formulas4[systematic]=[]
175 <                    #create TTreeFormulas
176 <                    for j in range(len( MVA_Vars['Nominal'])):
177 <                        MVA_formulas[systematic].append(ROOT.TTreeFormula("MVA_formula%s_%s"%(j,systematic),MVA_Vars[systematic][j],tree))
178 <                        MVA_formulas4[systematic].append(ROOT.TTreeFormula("MVA_formula4%s_%s"%(j,systematic),MVA_Vars['Nominal'][j]+'+('+MVA_Vars[systematic][j]+'-'+MVA_Vars['Nominal'][j]+')*4',tree))#HERE change
179 <                outfile.cd()
180 <                #Setup Branches
181 <                MVAbranches=[]
182 <                MVAbranches4=[]
183 <                for i in range(0,len(readers)):
184 <                    MVAbranches.append(array('f',[0]*11))
185 <                    MVAbranches4.append(array('f',[0]*11))
186 <                    newtree.Branch(MVAinfos[i].MVAname,MVAbranches[i],'nominal:JER_up:JER_down:JES_up:JES_down:beff_up:beff_down:bmis_up:bmis_down:beff1_up:beff1_down/F')
187 <                    newtree.Branch(MVAinfos[i].MVAname+'_4',MVAbranches4[i],'nominal:JER_up:JER_down:JES_up:JES_down:beff_up:beff_down:bmis_up:bmis_down:beff1_up:beff1_down/F')
188 <                print '\n--> ' + job.name +':'
189 <                #progbar setup
190 <                if nEntries >= longe:
191 <                    step=int(nEntries/longe)
192 <                    long=longe
193 <                else:
194 <                    long=nEntries
195 <                    step = 1
196 <                bar=progbar(long)
197 <                #Fill event by event:
198 <                for entry in range(0,nEntries):
199 <                    if entry % step == 0:
200 <                        bar.move()
201 <                    #load entry
202 <                    tree.GetEntry(entry)
203 <                    for systematic in systematics:
204 <                        for j in range(len( MVA_Vars['Nominal'])):
205 <                            MVA_var_buffer[j][0] = MVA_formulas[systematic][j].EvalInstance()
206 <                            
207 <                        for j in range(0,len(readers)):
208 <                            MVAbranches[j][systematics.index(systematic)] = readers[j].EvaluateMVA(MVAinfos[j].MVAname)
209 <                            
210 <                        for j in range(len( MVA_Vars['Nominal'])):
211 <                            MVA_var_buffer[j][0] = MVA_formulas4[systematic][j].EvalInstance()
212 <                            
213 <                        for j in range(0,len(readers)):
214 <                            MVAbranches4[j][systematics.index(systematic)] = readers[j].EvaluateMVA(MVAinfos[j].MVAname)
215 <                    #Fill:
216 <                    newtree.Fill()
217 <                newtree.AutoSave()
218 <                outfile.Close()
219 <                
220 <            #DATA:
221 <            if job.type == 'DATA':
222 <                #MVA Formulas
223 <                MVA_formulas_Nominal = []
224 <                #create TTreeFormulas
225 <                for j in range(len( MVA_Vars['Nominal'])):
226 <                    MVA_formulas_Nominal.append(ROOT.TTreeFormula("MVA_formula%s_Nominal"%j, MVA_Vars['Nominal'][j],tree))
227 <                outfile.cd()
228 <                MVAbranches=[]
229 <                for i in range(0,len(readers)):
162 >            #Set branch adress for all vars
163 >            for i in range(0,len(theMVAs)):
164 >                theMVAs[i].setBranches(tree,job)
165 >            outfile.cd()
166 >            #Setup Branches
167 >            MVAbranches=[]
168 >            for i in range(0,len(theMVAs)):
169 >                if job.type == 'Data':
170                      MVAbranches.append(array('f',[0]))
171                      newtree.Branch(MVAinfos[i].MVAname,MVAbranches[i],'nominal/F')
232                    newtree.Branch(MVAinfos[i].MVAname+'_4',MVAbranches[i],'nominal/F')
233                #progbar          
234                print '\n--> ' + job.name +':'
235                if nEntries >= longe:
236                    step=int(nEntries/longe)
237                    long=longe
172                  else:
173 <                    long=nEntries
174 <                    step = 1
175 <                bar=progbar(long)
176 <                #Fill event by event:
177 <                for entry in range(0,nEntries):
178 <                    if entry % step == 0:
179 <                        bar.move()
180 <                    #load entry
181 <                    tree.GetEntry(entry)
182 <                    #nominal:
183 <                    for j in range(len( MVA_Vars['Nominal'])):
184 <                            MVA_var_buffer[j][0] = MVA_formulas_Nominal[j].EvalInstance()
173 >                    MVAbranches.append(array('f',[0]*11))
174 >                    newtree.Branch(theMVAs[i].MVAname,MVAbranches[i],'nominal:JER_up:JER_down:JES_up:JES_down:beff_up:beff_down:bmis_up:bmis_down:beff1_up:beff1_down/F')
175 >                MVA_formulas_Nominal = []
176 >            print '\n--> ' + job.name +':'
177 >            #progbar setup
178 >            if nEntries >= longe:
179 >                step=int(nEntries/longe)
180 >                long=longe
181 >            else:
182 >                long=nEntries
183 >                step = 1
184 >            bar=progbar(long)
185 >            #Fill event by event:
186 >            for entry in range(0,nEntries):
187 >                if entry % step == 0:
188 >                    bar.move()
189 >                #load entry
190 >                tree.GetEntry(entry)
191                              
192 <                    for j in range(0,len(readers)):
193 <                        MVAbranches[j][0]= readers[j].EvaluateMVA(MVAinfos[j].MVAname)
194 <                    newtree.Fill()
195 <                newtree.AutoSave()
196 <                outfile.Close()
197 <
192 >                for i in range(0,len(theMVAs)):
193 >                    theMVAs[i].evaluate(MVAbranches[i],job)
194 >                #Fill:
195 >                newtree.Fill()
196 >            newtree.AutoSave()
197 >            outfile.Close()
198 >                
199   print '\n'
200  
201   #Update Info:

Diff Legend

Removed lines
+ Added lines
< Changed lines
> Changed lines