ViewVC Help
View File | Revision Log | Show Annotations | Root Listing
root/cvsroot/UserCode/OSUT3Analysis/Configuration/scripts/makePlots.py
Revision: 1.44
Committed: Fri Jun 21 13:06:43 2013 UTC (11 years, 10 months ago) by wulsin
Content type: text/x-python
Branch: MAIN
Changes since 1.43: +19 -6 lines
Log Message:
Add options --ylog, --ymin, --ymax, to specify y-axis range

File Contents

# User Rev Content
1 lantonel 1.1 #!/usr/bin/env python
2     import sys
3     import os
4 ahart 1.6 import re
5 lantonel 1.31 from math import *
6 lantonel 1.1 from array import *
7 lantonel 1.3 from decimal import *
8 lantonel 1.17 from optparse import OptionParser
9 lantonel 1.2 from OSUT3Analysis.Configuration.configurationOptions import *
10 lantonel 1.7 from OSUT3Analysis.Configuration.processingUtilities import *
11 lantonel 1.1
12 biliu 1.35
13 lantonel 1.36
14     ### parse the command-line options
15    
16 lantonel 1.17 parser = OptionParser()
17 lantonel 1.7 parser = set_commandline_arguments(parser)
18 wulsin 1.44 parser.add_option("--ylog", action="store_true", dest="setLogY", default=False,
19     help="Set logarithmic scale on vertical axis on all plots")
20     parser.add_option("--ymin", dest="setYMin",
21     help="Minimum of y axis")
22     parser.add_option("--ymax", dest="setYMax",
23     help="Maximum of y axis")
24    
25 lantonel 1.17 (arguments, args) = parser.parse_args()
26 lantonel 1.1
27 wulsin 1.44
28 lantonel 1.16 if arguments.localConfig:
29 lantonel 1.1 sys.path.append(os.getcwd())
30 lantonel 1.16 exec("from " + arguments.localConfig.rstrip('.py') + " import *")
31 lantonel 1.1
32 lantonel 1.25 #### deal with conflicting arguments
33 lantonel 1.16 if arguments.normalizeToData and arguments.normalizeToUnitArea:
34 lantonel 1.13 print "Conflicting normalizations requsted, will normalize to unit area"
35 lantonel 1.16 arguments.normalizeToData = False
36     if arguments.normalizeToData and arguments.noStack:
37 lantonel 1.14 print "You have asked to scale non-stacked backgrounds to data. This is a very strange request. Will normalize to unit area instead"
38 lantonel 1.16 arguments.normalizeToData = False
39     arguments.normalizeToUnitArea = True
40 lantonel 1.25 if arguments.makeRatioPlots and arguments.makeDiffPlots:
41     print "You have requested both ratio and difference plots. Will make just ratio plots instead"
42     arguments.makeRatioPlots = False
43 lantonel 1.13
44 lantonel 1.36
45     from ROOT import TFile, gROOT, gStyle, gDirectory, TStyle, THStack, TH1F, TCanvas, TString, TLegend, TLegendEntry, THStack, TIter, TKey, TPaveLabel, gPad
46    
47    
48     ### setting ROOT options so our plots will look awesome and everyone will love us
49 lantonel 1.1
50     gROOT.SetBatch()
51     gStyle.SetOptStat(0)
52     gStyle.SetCanvasBorderMode(0)
53     gStyle.SetPadBorderMode(0)
54     gStyle.SetPadColor(0)
55     gStyle.SetCanvasColor(0)
56     gStyle.SetTextFont(42)
57 lantonel 1.36 gStyle.SetCanvasDefH(600)
58     gStyle.SetCanvasDefW(600)
59     gStyle.SetCanvasDefX(0)
60     gStyle.SetCanvasDefY(0)
61     gStyle.SetPadTopMargin(0.05)
62     gStyle.SetPadBottomMargin(0.13)
63     gStyle.SetPadLeftMargin(0.15)
64     gStyle.SetPadRightMargin(0.05)
65     gStyle.SetTitleColor(1, "XYZ")
66     gStyle.SetTitleFont(42, "XYZ")
67     gStyle.SetTitleSize(0.05, "XYZ")
68     gStyle.SetTitleXOffset(0.95)
69     gStyle.SetTitleYOffset(1.25)
70 lantonel 1.40 gStyle.SetTextAlign(12)
71 lantonel 1.36 gStyle.SetLabelColor(1, "XYZ")
72     gStyle.SetLabelFont(42, "XYZ")
73     gStyle.SetLabelOffset(0.007, "XYZ")
74     gStyle.SetLabelSize(0.05, "XYZ")
75     gStyle.SetAxisColor(1, "XYZ")
76     gStyle.SetStripDecimals(True)
77     gStyle.SetTickLength(0.03, "XYZ")
78     gStyle.SetNdivisions(510, "XYZ")
79     gStyle.SetPadTickX(1)
80     gStyle.SetPadTickY(1)
81 lantonel 1.3 gROOT.ForceStyle()
82 lantonel 1.1
83 biliu 1.35
84 lantonel 1.36 #set the text for the luminosity label
85     if(intLumi < 1000.):
86     LumiText = "L_{int} = " + str(intLumi) + " pb^{-1}"
87     LumiText = "L_{int} = " + str.format('{0:.1f}', LumiInPb) + " pb^{-1}"
88     else:
89     LumiInFb = intLumi/1000.
90     LumiText = "L_{int} = " + str.format('{0:.1f}', LumiInFb) + " fb^{-1}"
91 lantonel 1.1
92 lantonel 1.36
93     #bestest place for lumi. label, in top left corner
94     topLeft_x_left = 0.1375839
95     topLeft_x_right = 0.4580537
96     topLeft_y_bottom = 0.8479021
97     topLeft_y_top = 0.9475524
98     topLeft_y_offset = 0.035
99    
100    
101     ##########################################################################################################################################
102     ##########################################################################################################################################
103     ##########################################################################################################################################
104 lantonel 1.1
105 lantonel 1.37 # some fancy-ass code from Andrzej Zuranski to merge bins in the ratio plot until the error goes below some threshold
106     def ratioHistogram( dataHist, mcHist, relErrMax=0.10):
107    
108     def groupR(group):
109     Data,MC = [float(sum(hist.GetBinContent(i) for i in group)) for hist in [dataHist,mcHist]]
110     return (Data-MC)/MC if MC else 0
111    
112     def groupErr(group):
113     Data,MC = [float(sum(hist.GetBinContent(i) for i in group)) for hist in [dataHist,mcHist]]
114     dataErr2,mcErr2 = [sum(hist.GetBinError(i)**2 for i in group) for hist in [dataHist,mcHist]]
115     return abs(math.sqrt( (dataErr2+mcErr2)/(Data-MC)**2 + mcErr2/MC**2 ) * (Data-MC)/MC) if Data and MC else 0
116    
117     def regroup(groups):
118     err,iG = max( (groupErr(g),groups.index(g)) for g in groups )
119     if err < relErrMax or len(groups)<3 : return groups
120     iH = max( [iG-1,iG+1], key = lambda i: groupErr(groups[i]) if 0<=i<len(groups) else -1 )
121     iLo,iHi = sorted([iG,iH])
122     return regroup(groups[:iLo] + [groups[iLo]+groups[iHi]] + groups[iHi+1:])
123    
124 lantonel 1.41 #don't rebin the histograms of the number of a given object (except for the pileup ones)
125 wulsin 1.43 if ((dataHist.GetName().find("num") is not -1 and dataHist.GetName().find("Primaryvertexs") is -1) or
126     dataHist.GetName().find("CutFlow") is not -1 or
127     dataHist.GetName().find("GenMatch") is not -1):
128 lantonel 1.41 ratio = dataHist.Clone()
129     ratio.Add(mcHist,-1)
130     ratio.Divide(mcHist)
131     ratio.SetTitle("")
132     else:
133     groups = regroup( [(i,) for i in range(1,1+dataHist.GetNbinsX())] )
134     ratio = TH1F("ratio","",len(groups), array('d', [dataHist.GetBinLowEdge(min(g)) for g in groups ] + [dataHist.GetXaxis().GetBinUpEdge(dataHist.GetNbinsX())]) )
135     for i,g in enumerate(groups) :
136     ratio.SetBinContent(i+1,groupR(g))
137     ratio.SetBinError(i+1,groupErr(g))
138    
139     ratio.GetYaxis().SetTitle("#frac{Data-MC}{MC}")
140 lantonel 1.37 ratio.SetLineColor(1)
141     ratio.SetLineWidth(2)
142     return ratio
143    
144     ##########################################################################################################################################
145     ##########################################################################################################################################
146     ##########################################################################################################################################
147    
148    
149 lantonel 1.36 def MakeOneDHist(pathToDir,histogramName):
150 biliu 1.35
151 lantonel 1.36
152     numBgMCSamples = 0
153     numDataSamples = 0
154     numSignalSamples = 0
155    
156     Stack = THStack("stack",histogramName)
157 lantonel 1.1
158 lantonel 1.36 LumiLabel = TPaveLabel(topLeft_x_left,topLeft_y_bottom,topLeft_x_right,topLeft_y_top,LumiText,"NDC")
159     LumiLabel.SetBorderSize(0)
160     LumiLabel.SetFillColor(0)
161     LumiLabel.SetFillStyle(0)
162    
163     NormLabel = TPaveLabel()
164     NormLabel.SetDrawOption("NDC")
165     NormLabel.SetX1NDC(topLeft_x_left)
166     NormLabel.SetX2NDC(topLeft_x_right)
167    
168     NormLabel.SetBorderSize(0)
169     NormLabel.SetFillColor(0)
170     NormLabel.SetFillStyle(0)
171    
172     NormText = ""
173     if arguments.normalizeToUnitArea:
174     NormText = "Scaled to unit area"
175     elif arguments.normalizeToData:
176     NormText = "MC scaled to data"
177     NormLabel.SetLabel(NormText)
178 lantonel 1.1
179 lantonel 1.36
180     BgMCLegend = TLegend()
181     BgTitle = BgMCLegend.AddEntry(0, "Data & Bkgd. MC", "H")
182     BgTitle.SetTextAlign(22)
183     BgTitle.SetTextFont(62)
184     BgMCLegend.SetBorderSize(0)
185     BgMCLegend.SetFillColor(0)
186     BgMCLegend.SetFillStyle(0)
187     SignalMCLegend = TLegend()
188     SignalTitle = SignalMCLegend.AddEntry(0, "Signal MC", "H")
189     SignalTitle.SetTextAlign(22)
190     SignalTitle.SetTextFont(62)
191     SignalMCLegend.SetBorderSize(0)
192     SignalMCLegend.SetFillColor(0)
193     SignalMCLegend.SetFillStyle(0)
194    
195     outputFile.cd(pathToDir)
196     Canvas = TCanvas(histogramName)
197     BgMCHistograms = []
198     BgMCLegendEntries = []
199     SignalMCHistograms = []
200     SignalMCLegendEntries = []
201     DataHistograms = []
202     DataLegendEntries = []
203    
204    
205     backgroundIntegral = 0
206     dataIntegral = 0
207     scaleFactor = 1
208    
209     for sample in processed_datasets: # loop over different samples as listed in configurationOptions.py
210     dataset_file = "%s/%s.root" % (condor_dir,sample)
211     inputFile = TFile(dataset_file)
212     Histogram = inputFile.Get(pathToDir+"/"+histogramName).Clone()
213     Histogram.SetDirectory(0)
214     inputFile.Close()
215     if arguments.rebinFactor:
216     RebinFactor = int(arguments.rebinFactor)
217     #don't rebin histograms which will have less than 5 bins or any gen-matching histograms
218 wulsin 1.38 if Histogram.GetNbinsX() >= RebinFactor*5 and Histogram.GetTitle().find("GenMatch") is -1:
219 lantonel 1.36 Histogram.Rebin(RebinFactor)
220    
221     xAxisLabel = Histogram.GetXaxis().GetTitle()
222     histoTitle = Histogram.GetTitle()
223    
224     legLabel = labels[sample]
225     if (arguments.printYields):
226     yieldHist = Histogram.Integral()
227     legLabel = legLabel + " (%.1f)" % yieldHist
228 lantonel 1.10
229 lantonel 1.36 if( types[sample] == "bgMC"):
230 lantonel 1.10
231 lantonel 1.36 numBgMCSamples += 1
232     backgroundIntegral += Histogram.Integral()
233 lantonel 1.32
234 lantonel 1.36 Histogram.SetLineStyle(1)
235     if(arguments.noStack):
236     Histogram.SetFillStyle(0)
237     Histogram.SetLineColor(colors[sample])
238     Histogram.SetLineWidth(2)
239     else:
240     Histogram.SetFillStyle(1001)
241     Histogram.SetFillColor(colors[sample])
242     Histogram.SetLineColor(1)
243     Histogram.SetLineWidth(1)
244 wulsin 1.29
245 lantonel 1.36 BgMCLegendEntries.append(legLabel)
246     BgMCHistograms.append(Histogram)
247 wulsin 1.29
248 lantonel 1.10
249 lantonel 1.36 elif( types[sample] == "signalMC"):
250    
251     numSignalSamples += 1
252    
253     Histogram.SetFillStyle(0)
254     Histogram.SetLineColor(colors[sample])
255     Histogram.SetLineStyle(1)
256     Histogram.SetLineWidth(2)
257     if(arguments.normalizeToUnitArea and Histogram.Integral() > 0):
258     Histogram.Scale(1./Histogram.Integral())
259    
260     SignalMCLegendEntries.append(legLabel)
261     SignalMCHistograms.append(Histogram)
262 lantonel 1.10
263 lantonel 1.36 elif( types[sample] == "data"):
264 lantonel 1.11
265 lantonel 1.36 numDataSamples += 1
266     dataIntegral += Histogram.Integral()
267    
268     Histogram.SetFillStyle(0)
269     Histogram.SetLineColor(colors[sample])
270     Histogram.SetLineStyle(1)
271     Histogram.SetLineWidth(2)
272     if(arguments.normalizeToUnitArea and Histogram.Integral() > 0):
273     Histogram.Scale(1./Histogram.Integral())
274 lantonel 1.32
275 lantonel 1.36 DataLegendEntries.append(legLabel)
276     DataHistograms.append(Histogram)
277 lantonel 1.10
278 lantonel 1.36 #scaling histograms as per user's specifications
279     if dataIntegral > 0 and backgroundIntegral > 0:
280     scaleFactor = dataIntegral/backgroundIntegral
281     for bgMCHist in BgMCHistograms:
282     if arguments.normalizeToData:
283     bgMCHist.Scale(scaleFactor)
284    
285     if arguments.normalizeToUnitArea and not arguments.noStack and backgroundIntegral > 0:
286     bgMCHist.Scale(1./backgroundIntegral)
287     elif arguments.normalizeToUnitArea and arguments.noStack and bgMCHist.Integral() > 0:
288     bgMCHist.Scale(1./bgMCHist.Integral())
289 lantonel 1.10
290 lantonel 1.36 if not arguments.noStack:
291     Stack.Add(bgMCHist)
292 lantonel 1.3
293 lantonel 1.14
294    
295 lantonel 1.36 ### formatting data histograms and adding to legend
296     legendIndex = 0
297     for Histogram in DataHistograms:
298 lantonel 1.41 BgMCLegend.AddEntry(Histogram,DataLegendEntries[legendIndex],"LEP")
299 lantonel 1.36 legendIndex = legendIndex+1
300    
301    
302     ### creating the histogram to represent the statistical errors on the stack
303     if numBgMCSamples is not 0 and not arguments.noStack:
304     ErrorHisto = BgMCHistograms[0].Clone("errors")
305     ErrorHisto.SetFillStyle(3001)
306     ErrorHisto.SetFillColor(13)
307     ErrorHisto.SetLineWidth(0)
308 lantonel 1.41 BgMCLegend.AddEntry(ErrorHisto,"Stat. Errors","F")
309 lantonel 1.36 for Histogram in BgMCHistograms:
310     if Histogram is not BgMCHistograms[0]:
311     ErrorHisto.Add(Histogram)
312    
313    
314     ### formatting bgMC histograms and adding to legend
315     legendIndex = numBgMCSamples-1
316     for Histogram in reversed(BgMCHistograms):
317     if(arguments.noStack):
318 lantonel 1.41 BgMCLegend.AddEntry(Histogram,BgMCLegendEntries[legendIndex],"L")
319 lantonel 1.36 else:
320 lantonel 1.41 BgMCLegend.AddEntry(Histogram,BgMCLegendEntries[legendIndex],"F")
321 lantonel 1.36 legendIndex = legendIndex-1
322    
323    
324     ### formatting signalMC histograms and adding to legend
325     legendIndex = 0
326     for Histogram in SignalMCHistograms:
327 lantonel 1.41 SignalMCLegend.AddEntry(Histogram,SignalMCLegendEntries[legendIndex],"L")
328 lantonel 1.36 legendIndex = legendIndex+1
329    
330    
331     ### finding the maximum value of anything going on the canvas, so we know how to set the y-axis
332     finalMax = 0
333     if numBgMCSamples is not 0 and not arguments.noStack:
334     finalMax = ErrorHisto.GetMaximum() + ErrorHisto.GetBinError(ErrorHisto.GetMaximumBin())
335     else:
336     for bgMCHist in BgMCHistograms:
337     if(bgMCHist.GetMaximum() > finalMax):
338     finalMax = bgMCHist.GetMaximum()
339     for signalMCHist in SignalMCHistograms:
340     if(signalMCHist.GetMaximum() > finalMax):
341     finalMax = signalMCHist.GetMaximum()
342     for dataHist in DataHistograms:
343     if(dataHist.GetMaximum() > finalMax):
344     finalMax = dataHist.GetMaximum() + dataHist.GetBinError(dataHist.GetMaximumBin())
345     finalMax = 1.15*finalMax
346 wulsin 1.44 if arguments.setYMax:
347     finalMax = float(arguments.setYMax)
348 lantonel 1.15
349    
350 lantonel 1.36 ### Drawing histograms to canvas
351 lantonel 1.32
352 lantonel 1.36 outputFile.cd(pathToDir)
353    
354     makeRatioPlots = arguments.makeRatioPlots
355     makeDiffPlots = arguments.makeDiffPlots
356 lantonel 1.32
357 lantonel 1.36 if numBgMCSamples is 0 or numDataSamples is not 1:
358     makeRatioPlots = False
359     makeDiffPlots = False
360     if makeRatioPlots or makeDiffPlots:
361     Canvas.SetFillStyle(0)
362     Canvas.Divide(1,2)
363     Canvas.cd(1)
364     gPad.SetPad(0.01,0.25,0.99,0.99)
365     gPad.SetMargin(0.1,0.05,0.02,0.07)
366     gPad.SetFillStyle(0)
367     gPad.Update()
368     gPad.Draw()
369 wulsin 1.44 if arguments.setLogY:
370     gPad.SetLogy()
371 lantonel 1.36 Canvas.cd(2)
372     gPad.SetPad(0.01,0.01,0.99,0.25)
373     #format: gPad.SetMargin(l,r,b,t)
374     gPad.SetMargin(0.1,0.05,0.4,0.02)
375     gPad.SetFillStyle(0)
376     gPad.SetGridy(1)
377     gPad.Update()
378     gPad.Draw()
379    
380     Canvas.cd(1)
381 lantonel 1.32
382 lantonel 1.36 if numBgMCSamples is not 0: # the first thing to draw to the canvas is a bgMC sample
383    
384 wulsin 1.44 yAxisMin = 0.0001
385     if arguments.setYMin:
386     yAxisMin = float(arguments.setYMin)
387 lantonel 1.36 if not arguments.noStack: # draw unstacked background samples
388     Stack.SetTitle(histoTitle)
389     Stack.Draw("HIST")
390     Stack.GetXaxis().SetTitle(xAxisLabel)
391     Stack.SetMaximum(finalMax)
392 wulsin 1.44 Stack.SetMinimum(yAxisMin)
393 lantonel 1.36 if makeRatioPlots or makeDiffPlots:
394     Stack.GetHistogram().GetXaxis().SetLabelSize(0)
395     #draw shaded error bands
396     ErrorHisto.Draw("A E2 SAME")
397    
398     else: #draw the stacked backgrounds
399     BgMCHistograms[0].SetTitle(histoTitle)
400     BgMCHistograms[0].Draw("HIST")
401     BgMCHistograms[0].GetXaxis().SetTitle(xAxisLabel)
402     BgMCHistograms[0].SetMaximum(finalMax)
403 wulsin 1.44 BgMCHistograms[0].SetMinimum(yAxisMin)
404 lantonel 1.36 for bgMCHist in BgMCHistograms:
405     bgMCHist.Draw("A HIST SAME")
406 lantonel 1.33
407 lantonel 1.36 for signalMCHist in SignalMCHistograms:
408     signalMCHist.Draw("A HIST SAME")
409     for dataHist in DataHistograms:
410     dataHist.Draw("A E SAME")
411 lantonel 1.32
412 lantonel 1.36
413     elif numSignalSamples is not 0: # the first thing to draw to the canvas is a signalMC sample
414     SignalMCHistograms[0].SetTitle(histoTitle)
415     SignalMCHistograms[0].Draw("HIST")
416     SignalMCHistograms[0].GetXaxis().SetTitle(xAxisLabel)
417     SignalMCHistograms[0].SetMaximum(finalMax)
418 wulsin 1.44 SignalMCHistograms[0].SetMinimum(yAxisMin)
419 lantonel 1.36
420     for signalMCHist in SignalMCHistograms:
421     if(signalMCHist is not SignalMCHistograms[0]):
422     signalMCHist.Draw("A HIST SAME")
423     for dataHist in DataHistograms:
424     dataHist.Draw("A E SAME")
425    
426    
427     elif(numDataSamples is not 0): # the first thing to draw to the canvas is a data sample
428     DataHistograms[0].SetTitle(histoTitle)
429     DataHistograms[0].Draw("E")
430     DataHistograms[0].GetXaxis().SetTitle(xAxisLabel)
431     DataHistograms[0].SetMaximum(finalMax)
432 wulsin 1.44 DataHistograms[0].SetMinimum(yAxisMin)
433 lantonel 1.36 for dataHist in DataHistograms:
434     if(dataHist is not DataHistograms[0]):
435     dataHist.Draw("A E SAME")
436 lantonel 1.32
437 lantonel 1.1
438 lantonel 1.3
439 lantonel 1.36 #legend coordinates, empirically determined :-)
440     x_left = 0.6761745
441     x_right = 0.9328859
442     x_width = x_right - x_left
443     y_max = 0.9335664
444     entry_height = 0.05
445    
446     if(numBgMCSamples is not 0 or numDataSamples is not 0): #then draw the data & bgMC legend
447    
448     numExtraEntries = 1 # count the legend title
449     BgMCLegend.SetX1NDC(x_left)
450     if numBgMCSamples > 0:
451     numExtraEntries = numExtraEntries + 1 # count the stat. errors entry
452 lantonel 1.10
453 lantonel 1.36 BgMCLegend.SetY1NDC(y_max-entry_height*(numExtraEntries+numBgMCSamples+numDataSamples))
454     BgMCLegend.SetX2NDC(x_right)
455     BgMCLegend.SetY2NDC(y_max)
456     BgMCLegend.Draw()
457    
458     if(numSignalSamples is not 0): #then draw the signalMC legend to the left of the other one
459     SignalMCLegend.SetX1NDC(x_left-x_width)
460     SignalMCLegend.SetY1NDC(y_max-entry_height*(1+numSignalSamples)) # add one for the title
461     SignalMCLegend.SetX2NDC(x_left)
462     SignalMCLegend.SetY2NDC(y_max)
463     SignalMCLegend.Draw()
464    
465     elif numSignalSamples is not 0: #draw the signalMC legend in the upper right corner
466     SignalMCLegend.SetX1NDC(x_left)
467     SignalMCLegend.SetY1NDC(y_max-entry_height*(1+numSignalSamples)) # add one for the title
468     SignalMCLegend.SetX2NDC(x_right)
469     SignalMCLegend.SetY2NDC(y_max)
470     SignalMCLegend.Draw()
471    
472    
473     if not arguments.normalizeToUnitArea or numDataSamples > 0: #don't draw the lumi label if there's no data and it's scaled to unit area
474     LumiLabel.Draw()
475     if arguments.normalizeToUnitArea or arguments.normalizeToData:
476     #move the normalization label down before drawing if we drew the lumi. label
477     NormLabel.SetY1NDC(topLeft_y_bottom-topLeft_y_offset)
478     NormLabel.SetY2NDC(topLeft_y_top-topLeft_y_offset)
479     NormLabel.Draw()
480    
481     elif arguments.normalizeToUnitArea or arguments.normalizeToData:
482     NormLabel.SetY1NDC(topLeft_y_bottom)
483     NormLabel.SetY2NDC(topLeft_y_top)
484     NormLabel.Draw()
485    
486    
487     if makeRatioPlots or makeDiffPlots:
488     Canvas.cd(2)
489     BgSum = Stack.GetStack().Last()
490     if makeRatioPlots:
491 lantonel 1.41 Comparison = ratioHistogram(DataHistograms[0],BgSum)
492 lantonel 1.36 elif makeDiffPlots:
493 wulsin 1.42 Comparison = DataHistograms[0].Clone("diff")
494 lantonel 1.41 Comparison.Add(BgSum,-1)
495     Comparison.SetTitle("")
496 lantonel 1.36 Comparison.GetYaxis().SetTitle("Data-MC")
497 lantonel 1.41 Comparison.GetXaxis().SetTitle(xAxisLabel)
498 lantonel 1.36 Comparison.GetYaxis().CenterTitle()
499     Comparison.GetYaxis().SetTitleSize(0.1)
500 lantonel 1.41 Comparison.GetYaxis().SetTitleOffset(0.4)
501 lantonel 1.36 Comparison.GetXaxis().SetTitleSize(0.15)
502     Comparison.GetYaxis().SetLabelSize(0.1)
503     Comparison.GetXaxis().SetLabelSize(0.15)
504 lantonel 1.41 if makeRatioPlots:
505     RatioYRange = 1.15
506     if arguments.ratioYRange:
507     RatioYRange = float(arguments.ratioYRange)
508 wulsin 1.39 Comparison.GetYaxis().SetRangeUser(-1*RatioYRange, RatioYRange)
509 lantonel 1.36 elif makeDiffPlots:
510     YMax = Comparison.GetMaximum()
511     YMin = Comparison.GetMinimum()
512     if YMax <= 0 and YMin <= 0:
513     Comparison.GetYaxis().SetRangeUser(-1.2*YMin,0)
514     elif YMax >= 0 and YMin >= 0:
515     Comparison.GetYaxis().SetRangeUser(0,1.2*YMax)
516     else: #axis crosses y=0
517     if abs(YMax) > abs(YMin):
518     Comparison.GetYaxis().SetRangeUser(-1.2*YMax,1.2*YMax)
519 lantonel 1.32 else:
520 lantonel 1.36 Comparison.GetYaxis().SetRangeUser(-1.2*YMin,1.2*YMin)
521    
522     Comparison.GetYaxis().SetNdivisions(205)
523     Comparison.Draw()
524 lantonel 1.10
525 lantonel 1.36 Canvas.Write()
526     if arguments.savePDFs:
527     pathToDirString = pathToDir.replace(' ','_').replace('<','lt').replace('>','gt').replace('(','').replace(')','').replace('=','eq')
528     Canvas.SaveAs(condor_dir+"/stacked_histograms_pdfs/"+pathToDirString+"/"+histogramName+".pdf")
529 lantonel 1.15
530    
531 lantonel 1.36 ##########################################################################################################################################
532     ##########################################################################################################################################
533     ##########################################################################################################################################
534 lantonel 1.1
535 lantonel 1.36 def MakeTwoDHist(pathToDir,histogramName):
536 lantonel 1.10 numBgMCSamples = 0
537     numDataSamples = 0
538     numSignalSamples = 0
539    
540     LumiLabel = TPaveLabel(0.1,0.8,0.34,0.9,LumiText,"NDC")
541     LumiLabel.SetBorderSize(0)
542     LumiLabel.SetFillColor(0)
543     LumiLabel.SetFillStyle(0)
544    
545 ahart 1.30 BgMCLegend = TLegend(0.76,0.65,0.99,0.9)
546     BgMCLegend.AddEntry (0, "Data & Bkgd. MC", "H").SetTextFont (62)
547 lantonel 1.10 BgMCLegend.SetBorderSize(0)
548     BgMCLegend.SetFillColor(0)
549     BgMCLegend.SetFillStyle(0)
550 ahart 1.30 SignalMCLegend = TLegend(0.76,0.135,0.99,0.377)
551     SignalMCLegend.AddEntry (0, "Signal MC", "H").SetTextFont (62)
552 lantonel 1.10 SignalMCLegend.SetBorderSize(0)
553     SignalMCLegend.SetFillColor(0)
554     SignalMCLegend.SetFillStyle(0)
555    
556 lantonel 1.36 outputFile.cd(pathToDir)
557 lantonel 1.10 Canvas = TCanvas(histogramName)
558     Canvas.SetRightMargin(0.2413793);
559     BgMCHistograms = []
560     SignalMCHistograms = []
561     DataHistograms = []
562    
563     for sample in processed_datasets: # loop over different samples as listed in configurationOptions.py
564 ahart 1.28 dataset_file = "%s/%s.root" % (condor_dir,sample)
565 lantonel 1.10 inputFile = TFile(dataset_file)
566 lantonel 1.36 Histogram = inputFile.Get(pathToDir+"/"+histogramName).Clone()
567 lantonel 1.10 Histogram.SetDirectory(0)
568 lantonel 1.25 RebinFactor = int(arguments.rebinFactor)
569     if arguments.rebinFactor and Histogram.GetNbinsX() >= RebinFactor*10 and Histogram.GetNbinsY() >= RebinFactor*10:
570     Histogram.Rebin2D(RebinFactor)
571 lantonel 1.10 inputFile.Close()
572     xAxisLabel = Histogram.GetXaxis().GetTitle()
573     yAxisLabel = Histogram.GetYaxis().GetTitle()
574     histoTitle = Histogram.GetTitle()
575    
576     if( types[sample] == "bgMC"):
577 lantonel 1.4
578 lantonel 1.10 numBgMCSamples += 1
579     Histogram.SetMarkerColor(colors[sample])
580     Histogram.SetFillColor(colors[sample])
581 ahart 1.30 BgMCLegend.AddEntry(Histogram,labels[sample],"F").SetTextFont (42)
582 lantonel 1.10 BgMCHistograms.append(Histogram)
583    
584     elif( types[sample] == "signalMC"):
585 lantonel 1.3
586 lantonel 1.10 numSignalSamples += 1
587     Histogram.SetMarkerColor(colors[sample])
588     Histogram.SetFillColor(colors[sample])
589 ahart 1.30 SignalMCLegend.AddEntry(Histogram,labels[sample],"F").SetTextFont (42)
590 lantonel 1.10 SignalMCHistograms.append(Histogram)
591    
592     elif( types[sample] == "data"):
593    
594     numDataSamples += 1
595     Histogram.SetMarkerColor(colors[sample])
596 wulsin 1.29 Histogram.SetFillColor(colors[sample])
597 ahart 1.30 BgMCLegend.AddEntry(Histogram,labels[sample],"F").SetTextFont (42)
598 lantonel 1.10 DataHistograms.append(Histogram)
599    
600    
601 lantonel 1.36 outputFile.cd(pathToDir)
602 lantonel 1.10
603     if(numBgMCSamples is not 0):
604     BgMCHistograms[0].SetTitle(histoTitle)
605     BgMCHistograms[0].GetXaxis().SetTitle(xAxisLabel)
606     BgMCHistograms[0].GetYaxis().SetTitle(yAxisLabel)
607     BgMCHistograms[0].Draw()
608     for signalMCHist in SignalMCHistograms:
609     signalMCHist.Draw("SAME")
610     for dataHist in DataHistograms:
611     dataHist.Draw("SAME")
612 lantonel 1.3
613 lantonel 1.10 elif(numSignalSamples is not 0):
614     SignalMCHistograms[0].SetTitle(histoTitle)
615     SignalMCHistograms[0].Draw()
616     SignalMCHistograms[0].GetXaxis().SetTitle(xAxisLabel)
617     SignalMCHistograms[0].GetYaxis().SetTitle(yAxisLabel)
618     for signalMCHist in SignalMCHistograms:
619     if(signalMCHist is not SignalMCHistograms[0]):
620     signalMCHist.Draw("SAME")
621     for dataHist in DataHistograms:
622     dataHist.Draw("SAME")
623    
624     elif(numDataSamples is not 0):
625     DataHistograms[0].SetTitle(histoTitle)
626     DataHistograms[0].GetXaxis().SetTitle(xAxisLabel)
627     DataHistograms[0].GetYaxis().SetTitle(yAxisLabel)
628     DataHistograms[0].Draw()
629     for dataHist in DataHistograms:
630     if(dataHist is not DataHistograms[0]):
631     dataHist.Draw("SAME")
632 lantonel 1.3
633 lantonel 1.8
634 lantonel 1.10 if(numBgMCSamples is not 0 or numDataSamples is not 0):
635     BgMCLegend.Draw()
636     if(numSignalSamples is not 0):
637     SignalMCLegend.Draw()
638 lantonel 1.25 if not arguments.normalizeToUnitArea or numDataSamples > 0:
639     LumiLabel.Draw()
640 lantonel 1.4
641 lantonel 1.10 Canvas.Write()
642 biliu 1.35
643    
644    
645    
646 lantonel 1.36 ##########################################################################################################################################
647     ##########################################################################################################################################
648     ##########################################################################################################################################
649    
650     processed_datasets = []
651    
652     condor_dir = set_condor_output_dir(arguments)
653    
654     #### check which input datasets have valid output files
655     for sample in datasets:
656     fileName = condor_dir + "/" + sample + ".root"
657     if not os.path.exists(fileName):
658     continue
659     testFile = TFile(fileName)
660     if testFile.IsZombie() or not testFile.GetNkeys():
661     continue
662     processed_datasets.append(sample)
663    
664     if len(processed_datasets) is 0:
665     sys.exit("No datasets have been processed")
666    
667    
668     #### make output file
669     outputFileName = "stacked_histograms.root"
670     if arguments.outputFileName:
671     outputFileName = arguments.outputFileName
672    
673     outputFile = TFile(condor_dir + "/" + outputFileName, "RECREATE")
674    
675    
676    
677     #### use the first input file as a template and make stacked versions of all its histograms
678     inputFile = TFile(condor_dir + "/" + processed_datasets[0] + ".root")
679     inputFile.cd()
680     outputFile.cd()
681    
682     if arguments.savePDFs:
683     os.system("rm -r %s/stacked_histograms_pdfs" % (condor_dir))
684     os.system("mkdir %s/stacked_histograms_pdfs" % (condor_dir))
685    
686    
687     #get root directory in the first layer, generally "OSUAnalysis"
688     for key in inputFile.GetListOfKeys():
689     if (key.GetClassName() != "TDirectoryFile"):
690     continue
691     rootDirectory = key.GetName()
692     outputFile.mkdir(rootDirectory)
693     if arguments.savePDFs:
694     os.system("mkdir %s/stacked_histograms_pdfs/%s" % (condor_dir,rootDirectory.replace(' ','_').replace('<','lt').replace('>','gt').replace('(','').replace(')','').replace('=','eq')))
695    
696     #cd to root directory and look for histograms
697     inputFile.cd(rootDirectory)
698     for key2 in gDirectory.GetListOfKeys():
699    
700     if re.match ('TH1', key2.GetClassName()): # found a 1-D histogram
701     MakeOneDHist(rootDirectory,key2.GetName())
702     elif re.match ('TH2', key2.GetClassName()) and arguments.draw2DPlots: # found a 2-D histogram
703     MakeTwoDHist(rootDirectory,key2.GetName())
704    
705     elif (key2.GetClassName() == "TDirectoryFile"): # found a directory, cd there and look for histograms
706     level2Directory = rootDirectory+"/"+key2.GetName()
707    
708     #make a corresponding directory in the output file
709     outputFile.cd(rootDirectory)
710     gDirectory.mkdir(key2.GetName())
711     if arguments.savePDFs:
712     os.system("mkdir %s/stacked_histograms_pdfs/%s" % (condor_dir,level2Directory.replace(' ','_').replace('<','lt').replace('>','gt').replace('(','').replace(')','').replace('=','eq')))
713    
714     #####################################################
715     ### This layer is typically the "channels" layer ###
716     #####################################################
717    
718     inputFile.cd(level2Directory)
719     for key3 in gDirectory.GetListOfKeys():
720     if re.match ('TH1', key3.GetClassName()): # found a 1-D histogram
721     MakeOneDHist(level2Directory,key3.GetName())
722     elif re.match ('TH2', key3.GetClassName()) and arguments.draw2DPlots: # found a 2-D histogram
723     MakeTwoDHist(level2Directory,key3.GetName())
724    
725     elif (key3.GetClassName() == "TDirectoryFile"): # found a directory, cd there and look for histograms
726     level3Directory = level2Directory+"/"+key3.GetName()
727    
728     #make a corresponding directory in the output file
729     outputFile.cd(level2Directory)
730     gDirectory.mkdir(key3.GetName())
731     if arguments.savePDFs:
732     os.system("mkdir %s/stacked_histograms_pdfs/%s" % (condor_dir,level3Directory.replace(' ','_').replace('<','lt').replace('>','gt').replace('(','').replace(')','').replace('=','eq')))
733    
734     #################################################
735     ### This layer is typically the "cuts" layer ###
736     #################################################
737    
738     inputFile.cd(level3Directory)
739     for key3 in gDirectory.GetListOfKeys():
740     if re.match ('TH1', key3.GetClassName()): # found a 1-D histogram
741     MakeOneDHist(level3Directory,key3.GetName())
742     elif re.match ('TH2', key3.GetClassName()) and arguments.draw2DPlots: # found a 2-D histogram
743     MakeTwoDHist(level3Directory,key3.GetName())
744    
745 lantonel 1.1
746     outputFile.Close()