ViewVC Help
View File | Revision Log | Show Annotations | Root Listing
root/cvsroot/UserCode/VHbb/python/tree_stack.py
Revision: 1.15
Committed: Wed Oct 3 12:03:49 2012 UTC (12 years, 7 months ago) by nmohr
Content type: text/x-python
Branch: MAIN
Changes since 1.14: +23 -12 lines
Log Message:
Plotting style

File Contents

# Content
1 #!/usr/bin/env python
2 from samplesclass import sample
3 from printcolor import printc
4 import pickle
5 import ROOT
6 from ROOT import TFile, TTree
7 import ROOT
8 from array import array
9 from BetterConfigParser import BetterConfigParser
10 import sys, os
11 from mvainfos import mvainfo
12 #from gethistofromtree import getHistoFromTree, orderandadd
13 from Ratio import getRatio
14 from optparse import OptionParser
15 from HistoMaker import HistoMaker, orderandadd
16 import TdrStyles
17
18 #CONFIGURE
19 argv = sys.argv
20 parser = OptionParser()
21 parser.add_option("-P", "--path", dest="path", default="",
22 help="path to samples")
23 parser.add_option("-R", "--reg", dest="region", default="",
24 help="region to plot")
25 parser.add_option("-C", "--config", dest="config", default=[], action="append",
26 help="configuration file")
27 (opts, args) = parser.parse_args(argv)
28 if opts.config =="":
29 opts.config = "config"
30 print opts.config
31 config = BetterConfigParser()
32 config.read(opts.config)
33 anaTag = config.get("Analysis","tag")
34
35 path = opts.path
36 region = opts.region
37
38 plotConfig = BetterConfigParser()
39 plotConfig.read('vhbbPlotDef.ini')
40
41 #get locations:
42 Wdir=config.get('Directories','Wdir')
43
44 section='Plot:%s'%region
45
46 Normalize = eval(config.get(section,'Normalize'))
47 log = eval(config.get(section,'log'))
48 blind = eval(config.get(section,'blind'))
49
50 infofile = open(path+'/samples.info','r')
51 info = pickle.load(infofile)
52 infofile.close()
53
54 #options = plot.split(',')
55
56 mass = config.get(section,'Signal')
57
58 vars = (config.get(section, 'vars')).split(',')
59
60 names = [plotConfig.get('plotDef:%s'%x,'relPath') for x in vars]
61 nBins = [eval(plotConfig.get('plotDef:%s'%x,'nBins')) for x in vars]
62 xMin = [eval(plotConfig.get('plotDef:%s'%x,'min')) for x in vars]
63 xMax = [eval(plotConfig.get('plotDef:%s'%x,'max')) for x in vars]
64 xAxis = [plotConfig.get('plotDef:%s'%x,'xAxis') for x in vars]
65
66 for p in range(0,len(names)):
67 if '<mass>' in names[p]:
68 newp= names[p].replace('<mass>',mass)
69 names[p]=newp
70 print names[p]
71
72 data = config.get(section,'Datas')
73 if config.has_option(section, 'Datacut'):
74 datacut=config.get(section, 'Datacut')
75 else:
76 datacut = region
77
78 options=[]
79
80 if blind: blindopt='blind'
81 else: blindopt = 'noblind'
82
83 for i in range(0,len(vars)):
84 options.append([names[i],'',xAxis[i],nBins[i],xMin[i],xMax[i],'%s_%s.pdf'%(region,vars[i]),region,datacut,mass,data,blindopt])
85
86
87 setup=config.get('Plot_general','setup')
88 setup=setup.split(',')
89
90 samples=config.get('Plot_general','samples')
91 samples=samples.split(',')
92
93 colorDict=eval(config.get('Plot_general','colorDict'))
94 #color=color.split(',')
95
96
97 weightF=config.get('Weights','weightF')
98 Group = eval(config.get('Plot_general','Group'))
99
100
101 #GETALL AT ONCE
102
103 Plotter=HistoMaker(path,config,region,options)
104
105 #print '\nProducing Plot of %s\n'%vars[v]
106 Lhistos = [[] for _ in range(0,len(vars))]
107 Ltyps = [[] for _ in range(0,len(vars))]
108 Ldatas = [[] for _ in range(0,len(vars))]
109 Ldatatyps = [[] for _ in range(0,len(vars))]
110 Ldatanames = [[] for _ in range(0,len(vars))]
111
112
113 #Find out Lumi:
114 for job in info:
115 if job.name in data: lumi_data=float(job.lumi)
116
117 Plotter.lumi=lumi_data
118
119 for job in info:
120 if eval(job.active):
121 if job.subsamples:
122 for subsample in range(0,len(job.subnames)):
123
124 if job.subnames[subsample] in samples:
125 hTempList, typList = Plotter.getHistoFromTree(job,subsample)
126 for v in range(0,len(vars)):
127 Lhistos[v].append(hTempList[v])
128 Ltyps[v].append(Group[job.subnames[subsample]])
129
130 else:
131 if job.name in samples:
132 #print job.getpath()
133 hTempList, typList = Plotter.getHistoFromTree(job)
134 for v in range(0,len(vars)):
135 Lhistos[v].append(hTempList[v])
136 Ltyps[v].append(Group[job.name])
137
138 elif job.name in data:
139 #print 'DATA'
140 hTemp, typ = Plotter.getHistoFromTree(job)
141 for v in range(0,len(vars)):
142 Ldatas[v].append(hTemp[v])
143 Ldatatyps[v].append(typ[v])
144 Ldatanames[v].append(job.name)
145
146
147 for v in range(0,len(vars)):
148
149 histos = Lhistos[v]
150 typs = Ltyps[v]
151 datas = Ldatas[v]
152 datatyps = Ldatatyps[v]
153 datanames= Ldatanames[v]
154
155 TdrStyles.tdrStyle()
156
157 c = ROOT.TCanvas(vars[v],'', 600, 600)
158 c.SetFillStyle(4000)
159 c.SetFrameFillStyle(1000)
160 c.SetFrameFillColor(0)
161
162 oben = ROOT.TPad('oben','oben',0,0.3 ,1.0,1.0)
163 oben.SetBottomMargin(0)
164 oben.SetFillStyle(4000)
165 oben.SetFrameFillStyle(1000)
166 oben.SetFrameFillColor(0)
167 unten = ROOT.TPad('unten','unten',0,0.0,1.0,0.3)
168 unten.SetTopMargin(0.)
169 unten.SetBottomMargin(0.35)
170 unten.SetFillStyle(4000)
171 unten.SetFrameFillStyle(1000)
172 unten.SetFrameFillColor(0)
173
174 oben.Draw()
175 unten.Draw()
176
177 oben.cd()
178
179 allStack = ROOT.THStack(vars[v],'')
180 l = ROOT.TLegend(0.75, 0.63, 0.88, 0.88)
181 MC_integral=0
182 MC_entries=0
183
184 for histo in histos:
185 MC_integral+=histo.Integral()
186 #MC_entries+=histo.GetEntries()
187 print "\033[1;32m\n\tMC integral = %s\033[1;m"%MC_integral
188
189 #ORDER AND ADD TOGETHER
190 #print typs
191 #print setup
192 histos, typs = orderandadd(histos,typs,setup)
193
194
195 k=len(histos)
196
197 for j in range(0,k):
198 #print histos[j].GetBinContent(1)
199 i=k-j-1
200 histos[i].SetFillColor(int(colorDict[setup[i]]))
201 histos[i].SetLineColor(1)
202 allStack.Add(histos[i])
203 l.AddEntry(histos[j],typs[j],'F')
204
205 d1 = ROOT.TH1F('noData','noData',nBins[v],xMin[v],xMax[v])
206 datatitle=''
207 for i in range(0,len(datas)):
208 d1.Add(datas[i],1)
209 if i ==0:
210 datatitle=datanames[i]
211 else:
212 datatitle=datatitle+ ' + '+datanames[i]
213 print "\033[1;32m\n\tDATA integral = %s\033[1;m"%d1.Integral()
214 flow = d1.GetEntries()-d1.Integral()
215 if flow > 0:
216 print "\033[1;31m\tU/O flow: %s\033[1;m"%flow
217 l.AddEntry(d1,datatitle,'PL')
218
219 if Normalize:
220 if MC_integral != 0: stackscale=d1.Integral()/MC_integral
221 stackhists=allStack.GetHists()
222 for blabla in stackhists:
223 if MC_integral != 0: blabla.Scale(stackscale)
224
225 allMC=ROOT.TH1F('allMC','allMC',nBins[v],xMin[v],xMax[v])
226 allMC.Sumw2()
227 for bin in range(0,nBins[v]):
228 allMC.SetBinContent(bin,allStack.GetStack().Last().GetBinContent(bin))
229 allMC.SetBinError(bin,allStack.GetStack().Last().GetBinError(bin))
230
231 allStack.SetTitle()
232 allStack.Draw("hist")
233 allStack.GetXaxis().SetTitle('')
234 yTitle = 'Entries'
235 if not '/' in yTitle:
236 yAppend = '%s' %(allStack.GetXaxis().GetBinWidth(1))
237 yTitle = '%s / %s' %(yTitle, yAppend)
238 allStack.GetYaxis().SetTitle(yTitle)
239 allStack.GetXaxis().SetRangeUser(xMin[v],xMax[v])
240 allStack.GetYaxis().SetRangeUser(0,20000)
241 theErrorGraph = ROOT.TGraphErrors(allMC)
242 theErrorGraph.SetFillColor(ROOT.kGray+3)
243 theErrorGraph.SetFillStyle(3013)
244 theErrorGraph.Draw('SAME2')
245 Ymax = max(allStack.GetMaximum(),d1.GetMaximum())*1.3
246 allStack.SetMaximum(Ymax)
247 allStack.SetMinimum(0.1)
248 c.Update()
249 if log:
250 ROOT.gPad.SetLogy()
251 ROOT.gPad.SetTicks(1,1)
252 #allStack.Draw("hist")
253 d1.Draw("E0same")
254 l.SetFillColor(0)
255 l.SetBorderSize(0)
256 l.Draw()
257
258
259
260 t = ROOT.TLatex()
261 t.SetNDC()
262 t.SetTextAlign(12)
263 t.SetTextSize(0.04)
264 t.DrawLatex(0.13,0.85,"CMS Preliminary")
265 t.SetTextSize(0.03)
266 t.DrawLatex(0.13,0.79,"#sqrt{s} = %s, L = %s fb^{-1}"%(anaTag,(float(lumi_data)/1000.)))
267
268 unten.cd()
269 ROOT.gPad.SetTicks(1,1)
270
271 ratio, error, ksScore, chiScore = getRatio(d1,allMC,xMin[v],xMax[v])
272 ksScore = allMC.KolmogorovTest( d1 )
273 chiScore = allMC.Chi2Test( d1 , "UWCHI2/NDF")
274 print ksScore
275 print chiScore
276 ratio.SetStats(0)
277 ratio.GetYaxis().SetRangeUser(0,2)
278 ratio.GetYaxis().SetNdivisions(502,0)
279 ratio.GetXaxis().SetTitle(xAxis[v])
280 ratio.Draw("E1")
281 ratio.SetTitle("")
282 m_one_line = ROOT.TLine(xMin[v],1,xMax[v],1)
283 m_one_line.SetLineStyle(7)
284 m_one_line.SetLineColor(4)
285 m_one_line.Draw("Same")
286
287 t = ROOT.TLatex()
288 t.SetNDC()
289 t.SetTextAlign(12)
290 t.SetTextSize(0.07)
291 t.DrawLatex(0.12,0.9,"K_{s}: %.2f"%(ksScore))
292 t.DrawLatex(0.12,0.4,"#chi_{#nu}^{2}: %.2f"%(chiScore))
293
294 name = '%s/%s' %(config.get('Directories','plotpath'),options[v][6])
295 c.Print(name)
296
297 os.system('rm %s/tmp_plotCache_%s*'%(config.get('Directories','plotpath'),region))
298 print 'i am done!\n'
299
300 sys.exit(0)