ViewVC Help
View File | Revision Log | Show Annotations | Root Listing
root/cvsroot/UserCode/VHbb/python/HistoMaker.py
Revision: 1.8
Committed: Fri Oct 12 09:09:01 2012 UTC (12 years, 7 months ago) by peller
Content type: text/x-python
Branch: MAIN
CVS Tags: hcpPreApp
Changes since 1.7: +10 -10 lines
Log Message:
BDT fix in HistoMaker

File Contents

# Content
1 from samplesclass import sample
2 from printcolor import printc
3 import pickle
4 import ROOT
5 from ROOT import TFile, TTree
6 import ROOT
7 from array import array
8 from BetterConfigParser import BetterConfigParser
9 import sys,os
10
11 class HistoMaker:
12 def __init__(self, path, config, region, optionsList,rescale=1,which_weightF='weightF'):
13 self.path = path
14 self.config = config
15 self.optionsList = optionsList
16 self.rescale = rescale
17 self.which_weightF=which_weightF
18 self.region = region
19 self.lumi=0.
20
21 def getScale(self,job,subsample=-1,MC_rescale_factor=1):
22 anaTag=self.config.get('Analysis','tag')
23 input = TFile.Open(self.path+'/'+job.getpath())
24 CountWithPU = input.Get("CountWithPU")
25 CountWithPU2011B = input.Get("CountWithPU2011B")
26 #print lumi*xsecs[i]/hist.GetBinContent(1)
27 if subsample>-1:
28 xsec=float(job.xsec[subsample])
29 sf=float(job.sf[subsample])
30 else:
31 xsec=float(job.xsec)
32 sf=float(job.sf)
33 theScale = 1.
34 if anaTag == '7TeV':
35 theScale = float(self.lumi)*xsec*sf/(0.46502*CountWithPU.GetBinContent(1)+0.53498*CountWithPU2011B.GetBinContent(1))*MC_rescale_factor/float(job.split)
36 elif anaTag == '8TeV':
37 theScale = float(self.lumi)*xsec*sf/(CountWithPU.GetBinContent(1))*MC_rescale_factor/float(job.split)
38 return theScale
39
40
41 def getHistoFromTree(self,job,subsample=-1):
42 if self.lumi == 0: raise Exception("You're trying to plot with no lumi")
43
44 hTreeList=[]
45 groupList=[]
46
47 #get the conversion rate in case of BDT plots
48 TrainFlag = eval(self.config.get('Analysis','TrainFlag'))
49 BDT_add_cut='EventForTraining == 0'
50
51
52 plot_path = self.config.get('Directories','plotpath')
53 addOverFlow=eval(self.config.get('Plot_general','addOverFlow'))
54
55 scratchDir = os.environ["TMPDIR"]
56 #scratchDir = '/shome/peller/'
57 # define treeCut
58 if job.type != 'DATA':
59 if type(self.region)==str:
60 cutcut=self.config.get('Cuts',self.region)
61 elif type(self.region)==list:
62 #replace vars with other vars in the cutstring (used in DC writer)
63 cutcut=self.config.get('Cuts',self.region[0])
64 cutcut=cutcut.replace(self.region[1],self.region[2])
65 #print cutcut
66 if subsample>-1:
67 treeCut='%s & %s'%(cutcut,job.subcuts[subsample])
68 else:
69 treeCut='%s'%(cutcut)
70 elif job.type == 'DATA':
71 cutcut=self.config.get('Cuts',self.region)
72 treeCut='%s'%(cutcut)
73
74 # get and skim the Trees
75 output=TFile.Open(scratchDir+'/tmp_plotCache_%s_%s.root'%(self.region,job.identifier),'recreate')
76 input = TFile.Open(self.path+'/'+job.getpath(),'read')
77 Tree = input.Get(job.tree)
78 output.cd()
79 CuttedTree=Tree.CopyTree(treeCut)
80
81 # get all Histos at once
82 weightF=self.config.get('Weights',self.which_weightF)
83 for options in self.optionsList:
84 if subsample>-1:
85 name=job.subnames[subsample]
86 group=job.group[subsample]
87 else:
88 name=job.name
89 group=job.group
90 treeVar=options[0]
91 name=options[1]
92 nBins=int(options[3])
93 xMin=float(options[4])
94 xMax=float(options[5])
95
96 #options
97
98 if job.type != 'DATA':
99 if CuttedTree.GetEntries():
100
101 if 'RTight' in treeVar or 'RMed' in treeVar: drawoption = '(%s)*(%s)'%(weightF,BDT_add_cut)
102 else: drawoption = '%s'%(weightF)
103 output.cd()
104 CuttedTree.Draw('%s>>%s(%s,%s,%s)' %(treeVar,name,nBins,xMin,xMax), drawoption, "goff,e")
105 full=True
106 else:
107 full=False
108 elif job.type == 'DATA':
109 if options[11] == 'blind':
110 output.cd()
111 if treeVar == 'H.mass':
112 CuttedTree.Draw('%s>>%s(%s,%s,%s)' %(treeVar,name,nBins,xMin,xMax),treeVar+'<80. || '+treeVar + '>150.' , "goff,e")
113 else:
114 CuttedTree.Draw('%s>>%s(%s,%s,%s)' %(treeVar,name,nBins,xMin,xMax),treeVar+'<0', "goff,e")
115
116 else:
117 output.cd()
118 CuttedTree.Draw('%s>>%s(%s,%s,%s)' %(treeVar,name,nBins,xMin,xMax),'1', "goff,e")
119 full = True
120 if full:
121 hTree = ROOT.gDirectory.Get(name)
122 else:
123 output.cd()
124 hTree = ROOT.TH1F('%s'%name,'%s'%name,nBins,xMin,xMax)
125 hTree.Sumw2()
126 if job.type != 'DATA':
127 if 'RTight' in treeVar or 'RMed' in treeVar:
128 if TrainFlag:
129 MC_rescale_factor=2.
130 print 'I RESCALE BY 2.0'
131 else: MC_rescale_factor = 1.
132 ScaleFactor = self.getScale(job,subsample,MC_rescale_factor)
133 else: ScaleFactor = self.getScale(job,subsample)
134 if ScaleFactor != 0:
135 hTree.Scale(ScaleFactor)
136 #print '\t-->import %s\t Integral: %s'%(job.name,hTree.Integral())
137 if addOverFlow:
138 uFlow = hTree.GetBinContent(0)+hTree.GetBinContent(1)
139 oFlow = hTree.GetBinContent(hTree.GetNbinsX()+1)+hTree.GetBinContent(hTree.GetNbinsX())
140 uFlowErr = ROOT.TMath.Sqrt(ROOT.TMath.Power(hTree.GetBinError(0),2)+ROOT.TMath.Power(hTree.GetBinError(1),2))
141 oFlowErr = ROOT.TMath.Sqrt(ROOT.TMath.Power(hTree.GetBinError(hTree.GetNbinsX()),2)+ROOT.TMath.Power(hTree.GetBinError(hTree.GetNbinsX()+1),2))
142 hTree.SetBinContent(1,uFlow)
143 hTree.SetBinContent(hTree.GetNbinsX(),oFlow)
144 hTree.SetBinError(1,uFlowErr)
145 hTree.SetBinError(hTree.GetNbinsX(),oFlowErr)
146 hTree.SetDirectory(0)
147 input.Close()
148 hTreeList.append(hTree)
149 groupList.append(group)
150
151 return hTreeList, groupList
152
153
154 ######################
155 def orderandadd(histos,typs,setup):
156 #ORDER AND ADD TOGETHER
157 ordnung=[]
158 ordnungtyp=[]
159 num=[0]*len(setup)
160 for i in range(0,len(setup)):
161 for j in range(0,len(histos)):
162 if typs[j] in setup[i]:
163 num[i]+=1
164 ordnung.append(histos[j])
165 ordnungtyp.append(typs[j])
166 del histos
167 del typs
168 histos=ordnung
169 typs=ordnungtyp
170 print typs
171 for k in range(0,len(num)):
172 for m in range(0,num[k]):
173 if m > 0:
174 #add
175 histos[k].Add(histos[k+1],1)
176 printc('magenta','','\t--> added %s to %s'%(typs[k],typs[k+1]))
177 del histos[k+1]
178 del typs[k+1]
179 del histos[len(setup):]
180 del typs[len(setup):]
181 return histos, typs