ViewVC Help
View File | Revision Log | Show Annotations | Root Listing
root/cvsroot/UserCode/VHbb/python/myutils/StackMaker.py
Revision: 1.11
Committed: Mon Mar 25 13:53:02 2013 UTC (12 years, 1 month ago) by nmohr
Content type: text/x-python
Branch: MAIN
Changes since 1.10: +5 -1 lines
Log Message:
Fix key error

File Contents

# User Rev Content
1 peller 1.1 import ROOT
2 nmohr 1.3 ROOT.gROOT.SetBatch(True)
3 peller 1.1 import sys,os
4     from BetterConfigParser import BetterConfigParser
5     import TdrStyles
6     from Ratio import getRatio
7 peller 1.4 from HistoMaker import HistoMaker
8 peller 1.1
9     class StackMaker:
10 peller 1.7 def __init__(self, config, var,region,SignalRegion,setup=None):
11 peller 1.1 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 peller 1.8 #if self.blind: blindopt='True'
20     #else: blindopt = 'False'
21 peller 1.7 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 nmohr 1.11 if not SignalRegion:
29     if 'ZH' in setup:
30     self.setup.remove('ZH')
31     if 'WH' in setup:
32     self.setup.remove('WH')
33 peller 1.1 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 peller 1.7 if config.has_option('Cuts',region):
56     cut = config.get('Cuts',region)
57     else:
58     cut = None
59 peller 1.1 if config.has_option(section, 'Datacut'):
60 nmohr 1.2 cut=config.get(section, 'Datacut')
61 nmohr 1.10 if config.has_option(section, 'doFit'):
62     self.doFit=eval(config.get(section, 'doFit'))
63     else:
64     self.doFit = False
65 nmohr 1.2
66 peller 1.1 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 peller 1.9 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 peller 1.7 if config.has_option('Weights','weightF'):
72     self.options['weight'] = config.get('Weights','weightF')
73     else:
74     self.options['weight'] = None
75 peller 1.1 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 nmohr 1.2 @staticmethod
90     def myText(txt="CMS Preliminary",ndcX=0,ndcY=0,size=0.8):
91 peller 1.1 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 nmohr 1.6 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 peller 1.1
131     def doPlot(self):
132     TdrStyles.tdrStyle()
133 peller 1.4 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 peller 1.1
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 nmohr 1.10 tLumi = self.myText("#sqrt{s} = %s, L = %.1f fb^{-1}"%(self.anaTag,(float(self.lumi)/1000.)),0.17,0.83)
273 peller 1.1 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 nmohr 1.10 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 peller 1.1
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 nmohr 1.2 name = '%s/%s' %(self.plotDir,self.options['pdfName'])
340 peller 1.1 c.Print(name)
341 nmohr 1.6 self.doCompPlot(allStack,l)