ViewVC Help
View File | Revision Log | Show Annotations | Root Listing
root/cvsroot/UserCode/OSUT3Analysis/Configuration/scripts/makePlots.py
Revision: 1.17
Committed: Wed Mar 13 09:27:11 2013 UTC (12 years, 1 month ago) by lantonel
Content type: text/x-python
Branch: MAIN
Changes since 1.16: +3 -4 lines
Log Message:
switched back to optparse, since CMSSW_5_X_X has an older version of python

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