ViewVC Help
View File | Revision Log | Show Annotations | Root Listing
root/cvsroot/UserCode/VHbb/python/myutils/StackMaker.py
Revision: 1.12
Committed: Mon Mar 25 14:41:33 2013 UTC (12 years, 1 month ago) by nmohr
Content type: text/x-python
Branch: MAIN
CVS Tags: lhcp_UnblindFix, hcp_Unblind, lhcp_11April, LHCP_PreAppFixAfterFreeze, LHCP_PreAppFreeze
Changes since 1.11: +2 -2 lines
Log Message:
Fix

File Contents

# Content
1 import ROOT
2 ROOT.gROOT.SetBatch(True)
3 import sys,os
4 from BetterConfigParser import BetterConfigParser
5 import TdrStyles
6 from Ratio import getRatio
7 from HistoMaker import HistoMaker
8
9 class StackMaker:
10 def __init__(self, config, var,region,SignalRegion,setup=None):
11 section='Plot:%s'%region
12 self.var = var
13 self.SignalRegion=SignalRegion
14 self.normalize = eval(config.get(section,'Normalize'))
15 self.log = eval(config.get(section,'log'))
16 if config.has_option('plotDef:%s'%var,'log') and not self.log:
17 self.log = eval(config.get('plotDef:%s'%var,'log'))
18 self.blind = eval(config.get(section,'blind'))
19 #if self.blind: blindopt='True'
20 #else: blindopt = 'False'
21 if setup is None:
22 self.setup=config.get('Plot_general','setup')
23 if self.log:
24 self.setup=config.get('Plot_general','setupLog')
25 self.setup=self.setup.split(',')
26 else:
27 self.setup=setup
28 if not SignalRegion:
29 if 'ZH' in self.setup:
30 self.setup.remove('ZH')
31 if 'WH' in self.setup:
32 self.setup.remove('WH')
33 self.rebin = 1
34 if config.has_option(section,'rebin'):
35 self.rebin = eval(config.get(section,'rebin'))
36 if config.has_option(section,'nBins'):
37 self.nBins = int(eval(config.get(section,'nBins'))/self.rebin)
38 else:
39 self.nBins = int(eval(config.get('plotDef:%s'%var,'nBins'))/self.rebin)
40 print self.nBins
41 if config.has_option(section,'min'):
42 self.xMin = eval(config.get(section,'min'))
43 else:
44 self.xMin = eval(config.get('plotDef:%s'%var,'min'))
45 if config.has_option(section,'max'):
46 self.xMax = eval(config.get(section,'max'))
47 else:
48 self.xMax = eval(config.get('plotDef:%s'%var,'max'))
49 self.name = config.get('plotDef:%s'%var,'relPath')
50 self.mass = config.get(section,'Signal')
51 data = config.get(section,'Datas')
52 if '<mass>' in self.name:
53 self.name = self.name.replace('<mass>',self.mass)
54 print self.name
55 if config.has_option('Cuts',region):
56 cut = config.get('Cuts',region)
57 else:
58 cut = None
59 if config.has_option(section, 'Datacut'):
60 cut=config.get(section, 'Datacut')
61 if config.has_option(section, 'doFit'):
62 self.doFit=eval(config.get(section, 'doFit'))
63 else:
64 self.doFit = False
65
66 self.colorDict=eval(config.get('Plot_general','colorDict'))
67 self.typLegendDict=eval(config.get('Plot_general','typLegendDict'))
68 self.anaTag = config.get("Analysis","tag")
69 self.xAxis = config.get('plotDef:%s'%var,'xAxis')
70 self.options = {'var': self.name,'name':'','xAxis': self.xAxis, 'nBins': self.nBins, 'xMin': self.xMin, 'xMax': self.xMax,'pdfName': '%s_%s_%s.pdf'%(region,var,self.mass),'cut':cut,'mass': self.mass, 'data': data, 'blind': self.blind}
71 if config.has_option('Weights','weightF'):
72 self.options['weight'] = config.get('Weights','weightF')
73 else:
74 self.options['weight'] = None
75 self.plotDir = config.get('Directories','plotpath')
76 self.maxRatioUncert = 0.5
77 if self.SignalRegion:
78 self.maxRatioUncert = 1000.
79 self.config = config
80 self.datas = None
81 self.datatyps = None
82 self.overlay = None
83 self.lumi = None
84 self.histos = None
85 self.typs = None
86 self.AddErrors = None
87 print self.setup
88
89 @staticmethod
90 def myText(txt="CMS Preliminary",ndcX=0,ndcY=0,size=0.8):
91 ROOT.gPad.Update()
92 text = ROOT.TLatex()
93 text.SetNDC()
94 text.SetTextColor(ROOT.kBlack)
95 text.SetTextSize(text.GetTextSize()*size)
96 text.DrawLatex(ndcX,ndcY,txt)
97 return text
98
99 def doCompPlot(self,aStack,l):
100 c = ROOT.TCanvas(self.var+'Comp','',600,600)
101 c.SetFillStyle(4000)
102 c.SetFrameFillStyle(1000)
103 c.SetFrameFillColor(0)
104 k=len(self.histos)
105 l.Clear()
106 maximum = 0.
107 for j in range(0,k):
108 #print histos[j].GetBinContent(1)
109 i=k-j-1
110 self.histos[i].SetLineColor(int(self.colorDict[self.typs[i]]))
111 self.histos[i].SetFillColor(0)
112 self.histos[i].SetLineWidth(3)
113 if self.histos[i].Integral() > 0.:
114 self.histos[i].Scale(1./self.histos[i].Integral())
115 if self.histos[i].GetMaximum() > maximum:
116 maximum = self.histos[i].GetMaximum()
117 l.AddEntry(self.histos[j],self.typLegendDict[self.typs[j]],'l')
118 aStack.SetMinimum(0.)
119 aStack.SetMaximum(maximum*1.5)
120 aStack.GetXaxis().SetTitle(self.xAxis)
121 aStack.Draw('HISTNOSTACK')
122 if self.overlay:
123 if self.overlay.Integral() > 0.:
124 self.overlay.Scale(1./self.overlay.Integral())
125 self.overlay.Draw('hist,same')
126 l.AddEntry(self.overlay,self.typLegendDict['Overlay'],'L')
127 l.Draw()
128 name = '%s/comp_%s' %(self.plotDir,self.options['pdfName'])
129 c.Print(name)
130
131 def doPlot(self):
132 TdrStyles.tdrStyle()
133 histo_dict = HistoMaker.orderandadd([{self.typs[i]:self.histos[i]} for i in range(len(self.histos))],self.setup)
134 #sort
135 self.histos=[histo_dict[key] for key in self.setup]
136 self.typs=self.setup
137
138 c = ROOT.TCanvas(self.var,'', 600, 600)
139 c.SetFillStyle(4000)
140 c.SetFrameFillStyle(1000)
141 c.SetFrameFillColor(0)
142
143 oben = ROOT.TPad('oben','oben',0,0.3 ,1.0,1.0)
144 oben.SetBottomMargin(0)
145 oben.SetFillStyle(4000)
146 oben.SetFrameFillStyle(1000)
147 oben.SetFrameFillColor(0)
148 unten = ROOT.TPad('unten','unten',0,0.0,1.0,0.3)
149 unten.SetTopMargin(0.)
150 unten.SetBottomMargin(0.35)
151 unten.SetFillStyle(4000)
152 unten.SetFrameFillStyle(1000)
153 unten.SetFrameFillColor(0)
154
155 oben.Draw()
156 unten.Draw()
157
158 oben.cd()
159 allStack = ROOT.THStack(self.var,'')
160 l = ROOT.TLegend(0.63, 0.55,0.92,0.92)
161 l.SetLineWidth(2)
162 l.SetBorderSize(0)
163 l.SetFillColor(0)
164 l.SetFillStyle(4000)
165 l.SetTextFont(62)
166 l.SetTextSize(0.035)
167 MC_integral=0
168 MC_entries=0
169
170 for histo in self.histos:
171 MC_integral+=histo.Integral()
172 print "\033[1;32m\n\tMC integral = %s\033[1;m"%MC_integral
173
174 #ORDER AND ADD TOGETHER
175 #print typs
176 #print setup
177
178
179 if not 'DYc' in self.typs: self.typLegendDict.update({'DYlight':self.typLegendDict['DYlc']})
180 print self.typLegendDict
181
182 k=len(self.histos)
183
184 for j in range(0,k):
185 #print histos[j].GetBinContent(1)
186 i=k-j-1
187 self.histos[i].SetFillColor(int(self.colorDict[self.typs[i]]))
188 self.histos[i].SetLineColor(1)
189 allStack.Add(self.histos[i])
190
191 d1 = ROOT.TH1F('noData','noData',self.nBins,self.xMin,self.xMax)
192 datatitle='Data'
193 addFlag = ''
194 if 'Zee' in self.datanames and 'Zmm' in self.datanames:
195 addFlag = 'Z(l^{-}l^{+})H(b#bar{b})'
196 elif 'Zee' in self.datanames:
197 addFlag = 'Z(e^{-}e^{+})H(b#bar{b})'
198 elif 'Zmm' in self.datanames:
199 addFlag = 'Z(#mu^{-}#mu^{+})H(b#bar{b})'
200 elif 'Znn' in self.datanames:
201 addFlag = 'Z(#nu#nu)H(b#bar{b})'
202 elif 'Wmn' in self.datanames:
203 addFlag = 'W(#mu#nu)H(b#bar{b})'
204 elif 'Wen' in self.datanames:
205 addFlag = 'W(e#nu)H(b#bar{b})'
206 for i in range(0,len(self.datas)):
207 d1.Add(self.datas[i],1)
208 print "\033[1;32m\n\tDATA integral = %s\033[1;m"%d1.Integral()
209 flow = d1.GetEntries()-d1.Integral()
210 if flow > 0:
211 print "\033[1;31m\tU/O flow: %s\033[1;m"%flow
212
213 if self.overlay:
214 self.overlay.SetLineColor(2)
215 self.overlay.SetLineWidth(2)
216 self.overlay.SetFillColor(0)
217 self.overlay.SetFillStyle(4000)
218 self.overlay.SetNameTitle('Overlay','Overlay')
219
220 l.AddEntry(d1,datatitle,'P')
221 for j in range(0,k):
222 l.AddEntry(self.histos[j],self.typLegendDict[self.typs[j]],'F')
223 if self.overlay:
224 l.AddEntry(self.overlay,self.typLegendDict['Overlay'],'L')
225
226 if self.normalize:
227 if MC_integral != 0: stackscale=d1.Integral()/MC_integral
228 if self.overlay:
229 self.overlay.Scale(stackscale)
230 stackhists=allStack.GetHists()
231 for blabla in stackhists:
232 if MC_integral != 0: blabla.Scale(stackscale)
233
234 #if self.SignalRegion:
235 # allMC=allStack.GetStack().At(allStack.GetStack().GetLast()-1).Clone()
236 #else:
237 allMC=allStack.GetStack().Last().Clone()
238
239 allStack.SetTitle()
240 allStack.Draw("hist")
241 allStack.GetXaxis().SetTitle('')
242 yTitle = 'Entries'
243 if not '/' in yTitle:
244 yAppend = '%.2f' %(allStack.GetXaxis().GetBinWidth(1))
245 yTitle = '%s / %s' %(yTitle, yAppend)
246 allStack.GetYaxis().SetTitle(yTitle)
247 allStack.GetXaxis().SetRangeUser(self.xMin,self.xMax)
248 allStack.GetYaxis().SetRangeUser(0,20000)
249 theErrorGraph = ROOT.TGraphErrors(allMC)
250 theErrorGraph.SetFillColor(ROOT.kGray+3)
251 theErrorGraph.SetFillStyle(3013)
252 theErrorGraph.Draw('SAME2')
253 l.AddEntry(theErrorGraph,"MC uncert. (stat.)","fl")
254 Ymax = max(allStack.GetMaximum(),d1.GetMaximum())*1.7
255 if self.log:
256 allStack.SetMinimum(0.1)
257 Ymax = Ymax*ROOT.TMath.Power(10,1.2*(ROOT.TMath.Log(1.2*(Ymax/0.1))/ROOT.TMath.Log(10)))*(0.2*0.1)
258 ROOT.gPad.SetLogy()
259 allStack.SetMaximum(Ymax)
260 c.Update()
261 ROOT.gPad.SetTicks(1,1)
262 #allStack.Draw("hist")
263 l.SetFillColor(0)
264 l.SetBorderSize(0)
265
266 if self.overlay:
267 self.overlay.Draw('hist,same')
268 d1.Draw("E,same")
269 l.Draw()
270
271 tPrel = self.myText("CMS Preliminary",0.17,0.88,1.04)
272 tLumi = self.myText("#sqrt{s} = %s, L = %.1f fb^{-1}"%(self.anaTag,(float(self.lumi)/1000.)),0.17,0.83)
273 tAddFlag = self.myText(addFlag,0.17,0.78)
274
275 unten.cd()
276 ROOT.gPad.SetTicks(1,1)
277
278 l2 = ROOT.TLegend(0.5, 0.82,0.92,0.95)
279 l2.SetLineWidth(2)
280 l2.SetBorderSize(0)
281 l2.SetFillColor(0)
282 l2.SetFillStyle(4000)
283 l2.SetTextFont(62)
284 #l2.SetTextSize(0.035)
285 l2.SetNColumns(2)
286
287
288 ratio, error = getRatio(d1,allMC,self.xMin,self.xMax,"",self.maxRatioUncert)
289 ksScore = d1.KolmogorovTest( allMC )
290 chiScore = d1.Chi2Test( allMC , "UWCHI2/NDF")
291 print ksScore
292 print chiScore
293 ratio.SetStats(0)
294 ratio.GetXaxis().SetTitle(self.xAxis)
295 ratioError = ROOT.TGraphErrors(error)
296 ratioError.SetFillColor(ROOT.kGray+3)
297 ratioError.SetFillStyle(3013)
298 ratio.Draw("E1")
299 if self.doFit:
300 fitData = ROOT.TF1("fData", "gaus",0.7, 1.3)
301 fitMC = ROOT.TF1("fMC", "gaus",0.7, 1.3)
302 print 'Fit on data'
303 d1.Fit(fitData,"R")
304 print 'Fit on simulation'
305 allMC.Fit(fitMC,"R")
306
307
308 if not self.AddErrors == None:
309 self.AddErrors.SetFillColor(5)
310 self.AddErrors.SetFillStyle(1001)
311 self.AddErrors.Draw('SAME2')
312
313 l2.AddEntry(self.AddErrors,"MC uncert. (stat. + syst.)","f")
314
315 #ksScore = d1.KolmogorovTest( self.AddErrors )
316 #chiScore = d1.Chi2Test( self.AddErrors , "UWCHI2/NDF")
317
318
319 l2.AddEntry(ratioError,"MC uncert. (stat.)","f")
320
321 l2.Draw()
322
323 ratioError.Draw('SAME2')
324 ratio.Draw("E1SAME")
325 ratio.SetTitle("")
326 m_one_line = ROOT.TLine(self.xMin,1,self.xMax,1)
327 m_one_line.SetLineStyle(ROOT.kDashed)
328 m_one_line.Draw("Same")
329
330 if not self.blind:
331 tKsChi = self.myText("#chi_{#nu}^{2} = %.3f K_{s} = %.3f"%(chiScore,ksScore),0.17,0.9,1.5)
332 t0 = ROOT.TText()
333 t0.SetTextSize(ROOT.gStyle.GetLabelSize()*2.4)
334 t0.SetTextFont(ROOT.gStyle.GetLabelFont())
335 if not self.log:
336 t0.DrawTextNDC(0.1059,0.96, "0")
337 if not os.path.exists(self.plotDir):
338 os.makedirs(os.path.dirname(self.plotDir))
339 name = '%s/%s' %(self.plotDir,self.options['pdfName'])
340 c.Print(name)
341 self.doCompPlot(allStack,l)