ViewVC Help
View File | Revision Log | Show Annotations | Root Listing
root/cvsroot/UserCode/OSUT3Analysis/Configuration/scripts/makePlots.py
Revision: 1.26
Committed: Wed Mar 27 14:01:36 2013 UTC (12 years, 1 month ago) by lantonel
Content type: text/x-python
Branch: MAIN
Changes since 1.25: +1 -1 lines
Log Message:
turned off 2D histgram plotting by default, can be enabled with --2D option

File Contents

# Content
1 #!/usr/bin/env python
2 import sys
3 import os
4 import re
5 from array import *
6 from decimal import *
7 from optparse import OptionParser
8 from OSUT3Analysis.Configuration.configurationOptions import *
9 from OSUT3Analysis.Configuration.processingUtilities import *
10
11 parser = OptionParser()
12 parser = set_commandline_arguments(parser)
13 (arguments, args) = parser.parse_args()
14
15 if arguments.localConfig:
16 sys.path.append(os.getcwd())
17 exec("from " + arguments.localConfig.rstrip('.py') + " import *")
18
19
20 outputFileName = "stacked_histograms.root"
21 if arguments.outputFileName:
22 outputFileName = arguments.outputFileName
23
24 condor_dir = set_condor_output_dir(arguments)
25
26
27
28 #### deal with conflicting arguments
29 if arguments.normalizeToData and arguments.normalizeToUnitArea:
30 print "Conflicting normalizations requsted, will normalize to unit area"
31 arguments.normalizeToData = False
32 if arguments.normalizeToData and arguments.noStack:
33 print "You have asked to scale non-stacked backgrounds to data. This is a very strange request. Will normalize to unit area instead"
34 arguments.normalizeToData = False
35 arguments.normalizeToUnitArea = True
36 if arguments.makeRatioPlots and arguments.makeDiffPlots:
37 print "You have requested both ratio and difference plots. Will make just ratio plots instead"
38 arguments.makeRatioPlots = False
39
40 from ROOT import TFile, gROOT, gStyle, gDirectory, TStyle, THStack, TH1F, TCanvas, TString, TLegend, TArrow, THStack, TIter, TKey, TPaveLabel, gPad
41
42 gROOT.SetBatch()
43 gStyle.SetOptStat(0)
44 gStyle.SetCanvasBorderMode(0)
45 gStyle.SetPadBorderMode(0)
46 gStyle.SetPadColor(0)
47 gStyle.SetCanvasColor(0)
48 gStyle.SetTextFont(42)
49 gROOT.ForceStyle()
50 outputFile = TFile(condor_dir + "/" + outputFileName, "RECREATE")
51
52 channels = []
53 processed_datasets = []
54
55 #### check which input datasets have valid output files
56 for sample in datasets:
57 fileName = condor_dir + "/" + sample + ".root"
58 if not os.path.exists(fileName):
59 continue
60 testFile = TFile(fileName)
61 if testFile.IsZombie() or not testFile.GetNkeys():
62 continue
63 processed_datasets.append(sample)
64
65 if len(processed_datasets) is 0:
66 sys.exit("No datasets have been processed")
67
68 #### open first input file and re-make its directory structure in the output file
69 testFile = TFile(condor_dir + "/" + processed_datasets[0] + ".root")
70 testFile.cd()
71 for key in testFile.GetListOfKeys():
72 if (key.GetClassName() != "TDirectoryFile"):
73 continue
74 outputFile.cd()
75 outputFile.mkdir(key.GetName())
76 rootDirectory = key.GetName()
77
78 testFile.cd(key.GetName())
79 for key2 in gDirectory.GetListOfKeys():
80 if (key2.GetClassName() != "TDirectoryFile"):
81 continue
82 outputFile.cd(key.GetName())
83 gDirectory.mkdir(key2.GetName())
84 channels.append(key2.GetName())
85
86
87 weight = intLumi / 10000.0
88 for dataset in processed_datasets:
89 dataset_file = "%s/%s.root" % (condor_dir,dataset)
90 fin = TFile (dataset_file)
91 flags = fin.Get ("flags")
92 noWeights = flags and flags.GetBinContent (1)
93 fin.Close ()
94 if types[dataset] != "data" and not noWeights:
95 os.system("mergeTFileServiceHistograms -i %s -o %s -w %g" % (dataset_file, dataset_file + "_tmp", weight))
96 else:
97 os.system("mergeTFileServiceHistograms -i %s -o %s -w %g" % (dataset_file, dataset_file + "_tmp", 1.0))
98
99 for channel in channels: # loop over final states, which each have their own directory
100
101 testFile.cd(rootDirectory+"/"+channel)
102
103 for key in gDirectory.GetListOfKeys(): # loop over histograms in the current directory
104 histogramName = key.GetName()
105
106 if re.match ('TH1', key.GetClassName()): # plot a 1-D histogram
107
108 numBgMCSamples = 0
109 numDataSamples = 0
110 numSignalSamples = 0
111
112 Stack = THStack("stack",histogramName)
113
114 if(intLumi < 1000.):
115 LumiText = "L_{int} = " + str(intLumi) + " pb^{-1}"
116 else:
117 getcontext().prec = 2
118 LumiInFb = intLumi/1000.
119 LumiText = "L_{int} = " + str(LumiInFb) + " fb^{-1}"
120
121 LumiLabel = TPaveLabel(0.1,0.8,0.34,0.9,LumiText,"NDC")
122 LumiLabel.SetBorderSize(0)
123 LumiLabel.SetFillColor(0)
124 LumiLabel.SetFillStyle(0)
125
126 BgMCLegend = TLegend(0.70,0.65,0.94,0.89, "Data & Bkgd. MC")
127 BgMCLegend.SetBorderSize(0)
128 BgMCLegend.SetFillColor(0)
129 BgMCLegend.SetFillStyle(0)
130 SignalMCLegend = TLegend(0.45,0.65,0.70,0.89,"Signal MC")
131 SignalMCLegend.SetBorderSize(0)
132 SignalMCLegend.SetFillColor(0)
133 SignalMCLegend.SetFillStyle(0)
134
135 outputFile.cd(rootDirectory+"/"+channel)
136 Canvas = TCanvas(histogramName)
137 BgMCHistograms = []
138 SignalMCHistograms = []
139 DataHistograms = []
140
141 backgroundIntegral = 0
142 dataIntegral = 0
143 scaleFactor = 1
144
145 for sample in processed_datasets: # loop over different samples as listed in configurationOptions.py
146 dataset_file = "%s/%s.root_tmp" % (condor_dir,sample)
147 inputFile = TFile(dataset_file)
148 Histogram = inputFile.Get(rootDirectory+"/"+channel+"/"+histogramName).Clone()
149 Histogram.SetDirectory(0)
150 RebinFactor = int(arguments.rebinFactor)
151 if arguments.rebinFactor and Histogram.GetNbinsX() >= RebinFactor*10:
152 Histogram.Rebin(RebinFactor)
153 inputFile.Close()
154 xAxisLabel = Histogram.GetXaxis().GetTitle()
155 histoTitle = Histogram.GetTitle()
156
157 if( types[sample] == "bgMC"):
158
159 numBgMCSamples += 1
160
161 if(arguments.noStack):
162 Histogram.SetFillStyle(0)
163 Histogram.SetLineColor(colors[sample])
164 Histogram.SetLineWidth(2)
165 BgMCLegend.AddEntry(Histogram,labels[sample],"L")
166 else:
167 Histogram.SetFillStyle(1001)
168 Histogram.SetFillColor(colors[sample])
169 Histogram.SetLineColor(1)
170 Histogram.SetLineWidth(1)
171 BgMCLegend.AddEntry(Histogram,labels[sample],"F")
172 Histogram.SetLineStyle(1)
173
174 backgroundIntegral += Histogram.Integral()
175
176 BgMCHistograms.append(Histogram)
177
178 elif( types[sample] == "signalMC"):
179
180 numSignalSamples += 1
181
182 Histogram.SetFillStyle(0)
183 Histogram.SetLineColor(colors[sample])
184 Histogram.SetLineStyle(1)
185 Histogram.SetLineWidth(2)
186 if(arguments.normalizeToUnitArea and Histogram.Integral() > 0):
187 Histogram.Scale(1./Histogram.Integral())
188 SignalMCLegend.AddEntry(Histogram,labels[sample],"L")
189 SignalMCHistograms.append(Histogram)
190
191 elif( types[sample] == "data"):
192
193 numDataSamples += 1
194
195 Histogram.SetFillStyle(0)
196 Histogram.SetLineColor(colors[sample])
197 Histogram.SetLineStyle(1)
198 Histogram.SetLineWidth(2)
199 if(arguments.normalizeToUnitArea and Histogram.Integral() > 0):
200 Histogram.Scale(1./Histogram.Integral())
201
202 dataIntegral += Histogram.Integral()
203
204 BgMCLegend.AddEntry(Histogram,labels[sample],"LEP")
205 DataHistograms.append(Histogram)
206
207 if dataIntegral > 0 and backgroundIntegral > 0:
208 scaleFactor = dataIntegral/backgroundIntegral
209 for bgMCHist in BgMCHistograms:
210 if arguments.normalizeToData:
211 bgMCHist.Scale(scaleFactor)
212 if arguments.normalizeToUnitArea and not arguments.noStack and backgroundIntegral > 0:
213 bgMCHist.Scale(1./backgroundIntegral)
214 elif arguments.normalizeToUnitArea and arguments.noStack and bgMCHist.Integral() > 0:
215 bgMCHist.Scale(1./bgMCHist.Integral())
216 if not arguments.noStack:
217 Stack.Add(bgMCHist)
218
219
220
221 finalMax = 0
222 if not arguments.noStack:
223 finalMax = Stack.GetMaximum()
224 else:
225 for bgMCHist in BgMCHistograms:
226 if(bgMCHist.GetMaximum() > finalMax):
227 finalMax = bgMCHist.GetMaximum()
228 for signalMCHist in SignalMCHistograms:
229 if(signalMCHist.GetMaximum() > finalMax):
230 finalMax = signalMCHist.GetMaximum()
231 for dataHist in DataHistograms:
232 if(dataHist.GetMaximum() > finalMax):
233 finalMax = dataHist.GetMaximum()
234
235
236
237 if len(DataHistograms) is 1:
238 dataIntegral += DataHistograms[0].Integral()
239
240
241 ### Drawing histograms to canvas
242
243
244 outputFile.cd(rootDirectory+"/"+channel)
245
246 makeRatioPlots = arguments.makeRatioPlots
247 makeDiffPlots = arguments.makeDiffPlots
248
249 if numBgMCSamples is 0 or numDataSamples is not 1:
250 makeRatioPlots = False
251 makeDiffPlots = False
252 if makeRatioPlots or makeDiffPlots:
253 Canvas.SetFillStyle(0)
254 Canvas.Divide(1,2)
255 Canvas.cd(1)
256 gPad.SetPad(0.01,0.25,0.99,0.99)
257 gPad.SetMargin(0.1,0.05,0.02,0.07)
258 gPad.SetFillStyle(0)
259 gPad.Update()
260 gPad.Draw()
261 Canvas.cd(2)
262 gPad.SetPad(0.01,0.01,0.99,0.25)
263 #format: gPad.SetMargin(l,r,b,t)
264 gPad.SetMargin(0.1,0.05,0.4,0.02)
265 gPad.SetFillStyle(0)
266 gPad.SetGridy(1)
267 gPad.Update()
268 gPad.Draw()
269
270 Canvas.cd(1)
271
272 if(numBgMCSamples is not 0):
273 if not arguments.noStack:
274 Stack.SetTitle(histoTitle)
275 Stack.Draw("HIST")
276 Stack.GetXaxis().SetTitle(xAxisLabel)
277 Stack.SetMaximum(1.1*finalMax)
278 Stack.SetMinimum(0.0001)
279 if makeRatioPlots or makeDiffPlots:
280 Stack.GetHistogram().GetXaxis().SetLabelSize(0)
281 else:
282 BgMCHistograms[0].SetTitle(histoTitle)
283 BgMCHistograms[0].Draw("HIST")
284 BgMCHistograms[0].GetXaxis().SetTitle(xAxisLabel)
285 BgMCHistograms[0].SetMaximum(1.1*finalMax)
286 BgMCHistograms[0].SetMinimum(0.0001)
287 for bgMCHist in BgMCHistograms:
288 bgMCHist.Draw("HIST SAME")
289 for signalMCHist in SignalMCHistograms:
290 signalMCHist.Draw("HIST SAME")
291 for dataHist in DataHistograms:
292 dataHist.Draw("E SAME")
293
294 elif(numSignalSamples is not 0):
295 SignalMCHistograms[0].SetTitle(histoTitle)
296 SignalMCHistograms[0].Draw("HIST")
297 SignalMCHistograms[0].GetXaxis().SetTitle(xAxisLabel)
298 SignalMCHistograms[0].SetMaximum(1.1*finalMax)
299 SignalMCHistograms[0].SetMinimum(0.0001)
300 for signalMCHist in SignalMCHistograms:
301 if(signalMCHist is not SignalMCHistograms[0]):
302 signalMCHist.Draw("HIST SAME")
303 for dataHist in DataHistograms:
304 dataHist.Draw("E SAME")
305
306 elif(numDataSamples is not 0):
307 DataHistograms[0].SetTitle(histoTitle)
308 DataHistograms[0].Draw("E")
309 DataHistograms[0].GetXaxis().SetTitle(xAxisLabel)
310 DataHistograms[0].SetMaximum(1.1*finalMax)
311 DataHistograms[0].SetMinimum(0.0001)
312 for dataHist in DataHistograms:
313 if(dataHist is not DataHistograms[0]):
314 dataHist.Draw("E SAME")
315
316
317 if(numBgMCSamples is not 0 or numDataSamples is not 0):
318 BgMCLegend.Draw()
319 if(numSignalSamples is not 0):
320 SignalMCLegend.Draw()
321
322 if not arguments.normalizeToUnitArea or numDataSamples > 0:
323 LumiLabel.Draw()
324 if arguments.normalizeToData and numBgMCSamples > 0 and numDataSamples > 0:
325 NormLabel = TPaveLabel(0.1,0.75,0.35,0.85,"MC scaled to data","NDC")
326 NormLabel.SetBorderSize(0)
327 NormLabel.SetFillColor(0)
328 NormLabel.SetFillStyle(0)
329 NormLabel.Draw()
330 elif arguments.normalizeToUnitArea:
331 NormLabel = TPaveLabel(0.1,0.75,0.35,0.85,"Scaled to unit area","NDC")
332 NormLabel.SetBorderSize(0)
333 NormLabel.SetFillColor(0)
334 NormLabel.SetFillStyle(0)
335 NormLabel.Draw()
336
337
338 if makeRatioPlots or makeDiffPlots:
339 Canvas.cd(2)
340 BgSum = Stack.GetStack().Last()
341 Comparison = DataHistograms[0].Clone()
342 Comparison.Add(BgSum,-1)
343 if not makeDiffPlots:
344 Comparison.Divide(BgSum)
345 Comparison.SetTitle("")
346 Comparison.GetXaxis().SetTitle(xAxisLabel)
347 if makeRatioPlots:
348 Comparison.GetYaxis().SetTitle("#frac{Data-MC}{MC}")
349 elif makeDiffPlots:
350 Comparison.GetYaxis().SetTitle("Data-MC")
351 Comparison.GetYaxis().CenterTitle()
352 Comparison.GetYaxis().SetTitleSize(0.1)
353 Comparison.GetYaxis().SetTitleOffset(0.35)
354 Comparison.GetXaxis().SetTitleSize(0.15)
355 Comparison.GetYaxis().SetLabelSize(0.1)
356 Comparison.GetXaxis().SetLabelSize(0.15)
357 if makeRatioPlots:
358 Comparison.GetYaxis().SetRangeUser(-1,1)
359 elif makeDiffPlots:
360 YMax = Comparison.GetMaximum()
361 YMin = Comparison.GetMinimum()
362 if YMax <= 0 and YMin <= 0:
363 Comparison.GetYaxis().SetRangeUser(-1.2*YMin,0)
364 elif YMax >= 0 and YMin >= 0:
365 Comparison.GetYaxis().SetRangeUser(0,1.2*YMax)
366 else: #axis crosses y=0
367 if abs(YMax) > abs(YMin):
368 Comparison.GetYaxis().SetRangeUser(-1.2*YMax,1.2*YMax)
369 else:
370 Comparison.GetYaxis().SetRangeUser(-1.2*YMin,1.2*YMin)
371
372 Comparison.GetYaxis().SetNdivisions(205)
373 Comparison.Draw()
374
375 Canvas.Write()
376
377
378 if re.match ('TH2', key.GetClassName()) and arguments.draw2DPlots: # plot a 2-D histogram
379
380 numBgMCSamples = 0
381 numDataSamples = 0
382 numSignalSamples = 0
383
384 if(intLumi < 1000.):
385 LumiText = "L_{int} = " + str(intLumi) + " pb^{-1}"
386 else:
387 getcontext().prec = 2
388 LumiInFb = intLumi/1000.
389 LumiText = "L_{int} = " + str(LumiInFb) + " fb^{-1}"
390
391 LumiLabel = TPaveLabel(0.1,0.8,0.34,0.9,LumiText,"NDC")
392 LumiLabel.SetBorderSize(0)
393 LumiLabel.SetFillColor(0)
394 LumiLabel.SetFillStyle(0)
395
396 BgMCLegend = TLegend(0.76,0.65,0.99,0.9, "Data & Bkgd. MC")
397 BgMCLegend.SetBorderSize(0)
398 BgMCLegend.SetFillColor(0)
399 BgMCLegend.SetFillStyle(0)
400 SignalMCLegend = TLegend(0.76,0.135,0.99,0.377,"Signal MC")
401 SignalMCLegend.SetBorderSize(0)
402 SignalMCLegend.SetFillColor(0)
403 SignalMCLegend.SetFillStyle(0)
404
405 outputFile.cd(rootDirectory+"/"+channel)
406 Canvas = TCanvas(histogramName)
407 Canvas.SetRightMargin(0.2413793);
408 BgMCHistograms = []
409 SignalMCHistograms = []
410 DataHistograms = []
411
412 for sample in processed_datasets: # loop over different samples as listed in configurationOptions.py
413 dataset_file = "%s/%s.root_tmp" % (condor_dir,sample)
414 inputFile = TFile(dataset_file)
415 Histogram = inputFile.Get(rootDirectory+"/"+channel+"/"+histogramName).Clone()
416 Histogram.SetDirectory(0)
417 RebinFactor = int(arguments.rebinFactor)
418 if arguments.rebinFactor and Histogram.GetNbinsX() >= RebinFactor*10 and Histogram.GetNbinsY() >= RebinFactor*10:
419 Histogram.Rebin2D(RebinFactor)
420 inputFile.Close()
421 xAxisLabel = Histogram.GetXaxis().GetTitle()
422 yAxisLabel = Histogram.GetYaxis().GetTitle()
423 histoTitle = Histogram.GetTitle()
424
425 if( types[sample] == "bgMC"):
426
427 numBgMCSamples += 1
428 Histogram.SetMarkerColor(colors[sample])
429 Histogram.SetFillColor(colors[sample])
430 BgMCLegend.AddEntry(Histogram,labels[sample],"F")
431 BgMCHistograms.append(Histogram)
432
433 elif( types[sample] == "signalMC"):
434
435 numSignalSamples += 1
436 Histogram.SetMarkerColor(colors[sample])
437 Histogram.SetFillColor(colors[sample])
438 SignalMCLegend.AddEntry(Histogram,labels[sample],"F")
439 SignalMCHistograms.append(Histogram)
440
441 elif( types[sample] == "data"):
442
443 numDataSamples += 1
444 Histogram.SetMarkerColor(colors[sample])
445 Histogram.SetFillColor(colors[sample])
446 BgMCLegend.AddEntry(Histogram,labels[sample],"F")
447 DataHistograms.append(Histogram)
448
449
450 outputFile.cd(rootDirectory+"/"+channel)
451
452 if(numBgMCSamples is not 0):
453 BgMCHistograms[0].SetTitle(histoTitle)
454 BgMCHistograms[0].GetXaxis().SetTitle(xAxisLabel)
455 BgMCHistograms[0].GetYaxis().SetTitle(yAxisLabel)
456 BgMCHistograms[0].Draw()
457 for signalMCHist in SignalMCHistograms:
458 signalMCHist.Draw("SAME")
459 for dataHist in DataHistograms:
460 dataHist.Draw("SAME")
461
462 elif(numSignalSamples is not 0):
463 SignalMCHistograms[0].SetTitle(histoTitle)
464 SignalMCHistograms[0].Draw()
465 SignalMCHistograms[0].GetXaxis().SetTitle(xAxisLabel)
466 SignalMCHistograms[0].GetYaxis().SetTitle(yAxisLabel)
467 for signalMCHist in SignalMCHistograms:
468 if(signalMCHist is not SignalMCHistograms[0]):
469 signalMCHist.Draw("SAME")
470 for dataHist in DataHistograms:
471 dataHist.Draw("SAME")
472
473 elif(numDataSamples is not 0):
474 DataHistograms[0].SetTitle(histoTitle)
475 DataHistograms[0].GetXaxis().SetTitle(xAxisLabel)
476 DataHistograms[0].GetYaxis().SetTitle(yAxisLabel)
477 DataHistograms[0].Draw()
478 for dataHist in DataHistograms:
479 if(dataHist is not DataHistograms[0]):
480 dataHist.Draw("SAME")
481
482
483 if(numBgMCSamples is not 0 or numDataSamples is not 0):
484 BgMCLegend.Draw()
485 if(numSignalSamples is not 0):
486 SignalMCLegend.Draw()
487 if not arguments.normalizeToUnitArea or numDataSamples > 0:
488 LumiLabel.Draw()
489
490 Canvas.Write()
491
492
493
494 for dataset in processed_datasets:
495 dataset_file = "%s/%s.root_tmp" % (condor_dir,dataset)
496 os.remove(dataset_file)
497
498 outputFile.Close()