ViewVC Help
View File | Revision Log | Show Annotations | Root Listing
root/cvsroot/UserCode/bm409/autoPlotter.py
(Generate patch)

Comparing UserCode/bm409/autoPlotter.py (file contents):
Revision 1.4 by bm409, Mon Aug 9 12:01:32 2010 UTC vs.
Revision 1.7 by bm409, Fri Aug 20 12:05:42 2010 UTC

# Line 1 | Line 1
1 + #!/usr/bin/env python
2 +
3 + '''Created by Bryn Mathias
4 + bryn.mathias@cern.ch
5 + '''
6 +
7   import ROOT as Root
8 + import math
9   from time import strftime
10   import os, commands
11 < Root.gROOT.SetBatch(True)
11 > Root.gROOT.SetBatch(True) # suppress the creation of canvases on the screen.. much much faster if over a remote connection
12   Root.gROOT.SetStyle("Plain") #To set plain bkgds for slides
13   Root.gStyle.SetTitleBorderSize(0)
14   Root.gStyle.SetCanvasBorderMode(0)
# Line 25 | Line 32 | Root.gStyle.SetLabelColor(1, "XYZ");
32   Root.gStyle.SetLabelFont(42, "XYZ");
33   Root.gStyle.SetLabelOffset(0.007, "XYZ");
34   Root.gStyle.SetLabelSize(0.05, "XYZ");
35 < intlumi = 0.2683 #inv pico barns
35 > Root.gStyle.SetHatchesLineWidth(3)
36 >
37 >
38 > #Please enter the integrated lumi for the plots here:
39 > intlumi = 0.838 #inv pico barns
40 > print "The Integrated Luminosity your plots are being scaled to is: ", intlumi , "pb^{-1}"
41 > #Set your output file here
42 > outputfile = "../pdfs/Calo/"
43 > print "you are outputting to: " ,outputfile
44 > #SetWhere your results from your analysis are
45 > resultsDir = " ../results/"
46  
47 < #close list is a global variable that has the file name appended to it so at the end of each hist loop the open files are closed.
47 > #close list is a global variable that has the file name appended to it so at the end of each hist loop the open files are closed. It is cleared each time a new hist is called.
48  
49   def GetHist(DataSetName,col,norm,Legend):
50 <    a = Root.TFile.Open(DataSetName)
51 <    closeList.append(a)
52 <    b = a.Get(DirKeys[dir].GetTitle())
53 <    Hist = b.Get(hist)
54 <    leg.AddEntry(Hist,Legend,"LP")
50 >    a = Root.TFile.Open(DataSetName) #open the file
51 >    closeList.append(a) # append the file to the close list
52 >    b = a.Get(DirKeys[dir].GetTitle()) #open the directory in the root file
53 >    Hist = b.Get(hist) # get your histogram by name
54 >    if Legend != 0:
55 >      leg.AddEntry(Hist,Legend,"LP") # add a legend entry
56      Hist.SetLineWidth(3)
57 <    Hist.SetLineColor(col)
57 >    Hist.SetLineColor(col) #set colour
58      if norm != 0:
59 <       Hist.Scale(intlumi/100.)
59 >       Hist.Scale(intlumi/100.) #if not data normilse to the data by lumi, MC is by default weighted to 100pb-1, if you have changed this change here!
60      return Hist
61  
62  
# Line 64 | Line 82 | def HistogramMaxX(H):
82  
83  
84  
85 <
85 > def Systematics(H1,H2,H3):
86 >  Standard = H1.Clone()
87 >  UpperError = H2.Clone()
88 >  LowerError = H3.Clone()
89 >  for bin in range(1,Standard.GetNbinsX()):
90 >    StandardUpper = math.sqrt((Standard.GetBinError(bin))**2 + ((UpperError.GetBinContent(bin) - LowerError.GetBinContent(bin))/2)**2)
91 >    Standard.SetBinError(bin, StandardUpper)
92 >  return Standard
93 >
94 >
95 > def PassingCutNumbers(Hist, lowerBound):
96 >  lowbin = Hist.FindBin(lowerBound)
97 >  errorVal = Root.Double(0)
98 >  passingCut = Hist.IntegralAndError(lowbin, Hist.GetNbinsX(), errorVal)
99 >  print passingCut, errorVal
100 >  textLine = "Sample =" + DirKeys[dir].GetTitle() + hist +  " , Number passing cut of" + "" + " is "
101 >  CutNumbers.write(textLine)
102 > #+ str(passingCut) + " pm "
103  
104  
105   time = strftime("%Y_%m_%d")
71 resultsDir = " ../results/"
106   print time
107   mergedPlots = Root.TFile("StandardPlots.root", "RECREATE")
108  
109  
110 < #List of files from which you wish to plot
111 < RootFileList = ["../results/AK5Calo_JetMETTau_Run2010A_Compleate.root","../results/AK5Calo_ZJets-madgraph.root","../results/AK5Calo_Zinvisible_jets.root","../results/AK5Calo_WJets-madgraph.root","../results/AK5Calo_ttbarTauola.root","../results/AK5Calo_LM0.root","../results/AK5Calo_TotalBackgrounds.root","../results/AK5Calo_QCD_AllPtBins_7TeV_Pythia6.root"]
110 > #A file to open and scan for histograms inside directories
111 > RootFileList = ["../results/AK5Calo_JetMET_Compleate.root"]
112  
113 < Webpage = open("../pdfs/HadronicPlots"+time+".html",'w')
113 > CutNumbers = open(outputfile+"CutTable.txt",'w')
114 >
115 > Webpage = open(outputfile+"StandardPlots"+time+".html",'w')
116 > #put your name in the contact email field
117   header = '''
118   <body>
119   <CENTER>
# Line 85 | Line 122 | header = '''
122   <table border="1" bordercolor="#FFCC00" style="background-color:#FFFFCC" width="1200" cellpadding="3" cellspacing="3">
123    <tr>
124      <td>DiJet</td>
125 <    <td>>=3 Jet</td>
125 >    <td>=3 Jet</td>
126 >    <td>>=4 Jet</td>
127 >    <td>>=2 Jet</td>
128 >
129    </tr>'''
130   footer = '''  </table>
131    </CENTER>
# Line 93 | Line 133 | footer = '''  </table>
133   Webpage.write(header)
134  
135  
96 print RootFileList
136  
137 < temp = Root.TFile.Open(RootFileList[1])
137 > temp = Root.TFile.Open(RootFileList[0])
138   DirKeys = temp.GetListOfKeys()
139  
140   HistKeys = [ (dir.ReadObj()).GetListOfKeys() for dir in DirKeys]
141   HistNames = [ [k.GetName() for k in D] for D in HistKeys]
142 <
142 > #Loops though all the histograms that have been read from the first input file, this is done by histogram name
143  
144   closeList = []
145   for dir in range(0,len(DirKeys)):
146    for hist in HistNames[dir]:
147      if "/" in hist : continue
148 <    leg = Root.TLegend(0.6, 0.6, 0.8, 0.8)
148 >    leg = Root.TLegend(0.7, 0.6, 0.89, 0.85)
149      leg.SetShadowColor(0)
150 +    leg.SetBorderSize(0)
151 +    leg.SetFillStyle(4100)
152      leg.SetFillColor(0)
153      leg.SetLineColor(0)
154      if "PFMetVsMHT" in hist: continue #This is the 2d Histogram in the file, dont want to plot this!
155      col=2
156      fileno = 1
157 <    MaxY = 1.0
158 <    MaxX = 1.0
159 <    MinX = 1000.
157 >    MaxY = 1.0#These three are the low terms for the axis and are reset for each histogram see below
158 >    MaxX = 1.0# to see the numbers set by the functions that are defined in the beginning of the program
159 >    MinX = 1000.#
160      c1 = Root.TCanvas("canvas"+hist,"canname"+hist,1200,1200)
161 <    # c1.Divide(1,2)#,0.4,0.4)
162 <    mainPad = Root.TPad("","",0.01,0.21,0.99,0.99);
161 >    #Make two pads, one small for the ratio plot and one large for the plot its self
162 >    mainPad = Root.TPad("","",0.01,0.25,0.99,0.99);
163      mainPad.SetNumber(1);
164      mainPad.Draw();
165 <    ratioPad = Root.TPad("","",0.01,0.01,0.99,0.2);
165 >    ratioPad = Root.TPad("","",0.01,0.01,0.99,0.25);
166      ratioPad.SetNumber(2);
167      ratioPad.Draw();
168  
169      c1.cd(1)
170 <    Data = GetHist("../results/AK5Calo_JetMETTau_Run2010A_Compleate.root",1,0,"Data")
170 >    #make your histograms form each file, read in the files you want below
171 >    # GetHist("RootFile",Colour, scale to lumi 0=false anything else = true, Legend entry)
172 >    #NB the order in which you book the histos is the order in which they appear in the legend
173 >    Data = GetHist("../results/AK5Calo_JetMET_Compleate.root",1,0,"Data")
174 >    QCD = GetHist("../results/AK5Calo_QCD_AllPtBins_7TeV_Pythia6.root",Root.kPink+4,1,"QCD")
175      Data.SetMarkerStyle(20)
176      Data.SetLineWidth(3)
177      Data.SetLineColor(1)
178      Data.SetFillColor(1)
179 <    ZJets  = GetHist("../results/AK5Calo_ZJets-madgraph.root",Root.kBlue,1,"Z+Jets (Madgraph)")
180 <    Zinv = GetHist("../results/AK5Calo_Zinvisible_jets.root",Root.kGray,1,"Z->Invisiable")
181 <    WJets = GetHist("../results/AK5Calo_WJets-madgraph.root",Root.kGreen,1,"W+Jets")
182 <    ttbar = GetHist("../results/AK5Calo_ttbarTauola.root",Root.kAzure,1,"TTBar")
183 <    LM0 = GetHist("../results/AK5Calo_LM0.root",Root.kRed,1,"LM0")
184 <    QCD = GetHist("../results/AK5Calo_QCD_AllPtBins_7TeV_Pythia6.root",Root.kGray,1,"QCD")
140 <    # Total = GetHist("../results/AK5Calo_TotalBackgrounds.root",Root.kGreen-3,1,"Standard Model Backgrounds")
179 >    ZJets  = GetHist("../results/AK5Calo_ZJets-madgraph.root",Root.kTeal-7,1,"Z+Jets (Madgraph)")
180 >    Zinv = GetHist("../results/AK5Calo_Zinvisible_jets.root",91,1,"Z->#nu#nu + Jets")
181 >    WJets = GetHist("../results/AK5Calo_WJets-madgraph.root",Root.kPink+7,1,"W+Jets")
182 >    ttbar = GetHist("../results/AK5Calo_ttbarTauola.root",Root.kBlue+1,1,"TTBar")
183 >    LM0 = GetHist("../results/AK5Calo_LM0.root",2,1,"LM0")
184 >
185  
186 +    # ScaledUp = GetHist("../results/AK5Calo_Scale1.1_QCD_Pythia6.root",Root.kTeal-7,1,0)
187 +    # ScaledDown = GetHist("../results/AK5Calo_Scale0.9_QCD_Pythia6.root",Root.kTeal-7,1,0)
188 +    #Make a histogram of the sum of all the SM backgrounds
189      Total = QCD.Clone()
190      Total.Add(Zinv)
191      Total.Add(ZJets)
192      Total.Add(WJets)
193      Total.Add(ttbar)
194 <    Total.SetLineWidth(2)
195 <    Total.SetLineColor(Root.kAzure-4)
196 <    leg.AddEntry(Total,"Standard Model Backgrounds","LPf")
194 >    if "AlphaT_all" in hist : PassingCutNumbers(LM0, 0.55)
195 >    #Some funkey stuff from Henning for drawing the total background with some interesting lines around it
196 >    hcen=Total.Clone()
197 >    herr=Total.Clone()
198 >    # Total = Systematics(Total,ScaledUp,ScaledDown)
199 >    herr.SetLineColor(Root.kTeal+3)
200 >    herr.SetMarkerColor(Root.kAzure+6)
201 >    herr.SetFillColor(Root.kAzure+6)
202 >
203 >    herr.SetLineWidth(3)
204 >    Total.SetLineWidth(3)
205 >
206 >    Total.SetFillColor(Root.kAzure+2)
207 >    Total.SetLineColor(Root.kAzure+2)
208 >    Total.SetFillStyle(3245)
209 >    hcen.SetFillColor(0)
210 >    hcen.SetMarkerColor(Root.kTeal+3)
211 >    hcen.SetLineColor(Root.kTeal+3)
212 >    hcen.SetLineWidth(3)
213  
214 +    leg.AddEntry(Total,"Standard Model","LPf")
215  
216 +    #Defind the ranges of the histogram for the two highest histograms ie the data and the total
217      if HistogramMinX(Total) < MinX :
218         MinX = HistogramMinX(Total)
219      if HistogramMaxX(Total) > MaxX :
# Line 164 | Line 229 | for dir in range(0,len(DirKeys)):
229         MaxY = HistogramMaxY(Data)
230  
231  
232 <    Total.SetFillColor(Root.kAzure+10)
233 <    # Total.SetFillStyle(3001)
234 <    Total.GetXaxis().SetRangeUser(MinX,MaxX)
235 <    Total.GetYaxis().SetRangeUser(0.00001,MaxY*1.2)
236 <    Total.Draw("9HIST")
232 >    herr.GetXaxis().SetRangeUser(MinX,MaxX)
233 >    herr.GetYaxis().SetRangeUser(0.00001,MaxY*5.)
234 >    #Draw order is important!
235 >    herr.Draw("hist")
236 >    hcen.Draw("histsame")
237 >    Total.Draw("9E2same")
238      ZJets.Draw("9Sameh")
239      Zinv.Draw("9SAMEh")
240      WJets.Draw("9SAMEh")
241      ttbar.Draw("9SAMEh")
242      QCD.Draw("9SAMEh")
243      LM0.Draw("9SAMEh")
244 +    #Draw Data last to get the points above everything else
245      Data.Draw("9SAMEP")
246  
247 <    title = Root.TLatex(0.2,0.95,(DirKeys[dir].GetTitle())+"_"+hist[0:-4])
247 >    #Draw Text for histogram titile and for lumi
248 >    prelim = Root.TLatex(0.5,0.75,"CMS preliminary 2010");
249 >    lumi = Root.TLatex(0.4,.85,"#scale[0.8]{#int L dt = " + str(intlumi) + "pb^{-1}, #sqrt{s} = 7 TeV}");
250 >    prelim.SetNDC()
251 >    lumi.SetNDC()
252 >    # prelim.Draw()
253 >    lumi.Draw()
254 >    title = Root.TLatex(0.05,0.95,(DirKeys[dir].GetTitle())+"_"+hist[0:-4])
255      title.SetNDC()
256      title.Draw()
257 +    lumi.Draw("Same")
258 +
259  
260 +    #Ratio plot is between the sum of the background and the data - clone the two relevant histograms
261      RatioBottom = Total.Clone()
262      RatioTop = Data.Clone()
263 +    RatioTop.GetYaxis().SetTitle("data / sim")
264      RatioTop.Divide(RatioBottom)
265      c1.cd(1).Update()
266      # c1.cd(1).SetLogy()
267      leg.Draw()
190
268      # c1.cd(1).Update()
269      c1.cd(2)
193    RatioTop.GetXaxis().SetRangeUser(MinX,MaxX)
270  
271 +    RatioTop.SetTitleSize(0.1, "XYZ");
272 +    RatioTop.SetTitleOffset(0.55, "X");
273 +    RatioTop.SetTitleOffset(0.3, "Y");
274 +    RatioTop.SetLabelSize(0.06,"XY")
275 +
276 +    # RatioTop.SetTitleXOffset(0.9);
277 +    # RatioTop.SetTitleYOffset(0.9);
278 +    RatioTop.GetXaxis().SetRangeUser(MinX,MaxX)
279      RatioTop.GetYaxis().SetRangeUser(0.,2.0)
280      RatioTop.Draw()
281 +    #Draw a line though the perfectly matching point
282      unity = Root.TLine();
283      unity.SetLineWidth(2);
284      # unity.SetLineStyle(Root.kDashed);
# Line 202 | Line 287 | for dir in range(0,len(DirKeys)):
287  
288  
289  
290 <
290 >    #Do not save all of the histograms produced by the analysis framework, only the ones containing all objects in your selection, if this is not what you want to do then remove the if "all in hist" line
291      if "all" in hist :
292  
293 <      c1.SaveAs("../pdfs/"+(DirKeys[dir].GetTitle())+"_"+hist+".png")
293 >      # c1.SaveAs(outputfile+(DirKeys[dir].GetTitle())+"_"+hist+".png")
294        if (DirKeys[dir].GetTitle())[0] != "n":
295  
296          Images ='''<tr>
297 <          <td><a href="'''+(DirKeys[dir].GetTitle())+"_"+hist+'.png\"><img src=\"'+(DirKeys[dir].GetTitle())+"_"+hist+'.png\" width =\"598\" height=\"286\" /></a></td> \n' + "<td><a href=\""+"n"+(DirKeys[dir].GetTitle())+"_"+hist+'.png\"><img src=\"'+"n"+(DirKeys[dir].GetTitle())+"_"+hist+'.png\" width =\"598\" height=\"286\" /></a></td> <tr> \n' + "<td><a href=\""+(DirKeys[dir].GetTitle())+"_Log_"+hist+'.png\"><img src=\"'+(DirKeys[dir].GetTitle())+"_Log_"+hist+'.png\" width =\"598\" height=\"286\" /></a></td>\n' "<td><a href=\""+"n"+(DirKeys[dir].GetTitle())+"_Log_"+hist+'.png\"><img src=\"'+"n"+(DirKeys[dir].GetTitle())+"_Log_"+hist+'.png\" width =\"598\" height=\"286\" /></a></td>\n'
297 >          <td><a href="'''+(DirKeys[dir].GetTitle())+"_Log_"+hist+'.png\"><img src=\"'+(DirKeys[dir].GetTitle())+"_Log_"+hist+'.png\" width =\"286\" height=\"286\" /></a></td>\n'+  "<td><a href=\""+"n"+(DirKeys[dir].GetTitle())+"_Log_"+hist+'.png\"><img src=\"'+"n"+(DirKeys[dir].GetTitle())+"_Log_"+hist+'.png\" width =\"286\" height=\"286\" /></a></td>\n'
298          Webpage.write(Images)
299  
300        # mergedPlots.writeTObject(c1)
# Line 219 | Line 304 | for dir in range(0,len(DirKeys)):
304        # c1.cd(2)
305        # RatioTop.Draw()
306  
307 <      c1.SaveAs("../pdfs/"+(DirKeys[dir].GetTitle())+"_Log_"+hist+".png")
307 >      c1.SaveAs(outputfile+(DirKeys[dir].GetTitle())+"_Log_"+hist+".png")
308 >      c1.SaveAs(outputfile+(DirKeys[dir].GetTitle())+"_Log_"+hist+".pdf")
309  
310      if "AlphaT_0" in hist :
311        # c1.cd(2)
312        # RatioTop.Draw()
313 <      c1.SaveAs("../pdfs/"+(DirKeys[dir].GetTitle())+"_"+hist+".png")
313 >      # c1.SaveAs(outputfile+(DirKeys[dir].GetTitle())+"_"+hist+".png")
314        if (DirKeys[dir].GetTitle())[0] != "n":
315  
316          Images ='''<tr>
317 <          <td><a href="'''+(DirKeys[dir].GetTitle())+"_"+hist+'.png\"><img src=\"'+(DirKeys[dir].GetTitle())+"_"+hist+'.png\" width =\"598\" height=\"286\" /></a></td>' + "<td><a href=\""+"n"+(DirKeys[dir].GetTitle())+"_"+hist+'.png\"><img src=\"'+"n"+(DirKeys[dir].GetTitle())+"_"+hist+'.png\" width =\"598\" height=\"286\" /></a></td> <tr>' + "<td><a href=\""+(DirKeys[dir].GetTitle())+"_Log_"+hist+'.png\"><img src=\"'+(DirKeys[dir].GetTitle())+"_Log_"+hist+'.png\" width =\"598\" height=\"286\" /></a></td>'"<td><a href=\""+"n"+(DirKeys[dir].GetTitle())+"_Log_"+hist+'.png\"><img src=\"'+"n"+(DirKeys[dir].GetTitle())+"_Log_"+hist+'.png\" width =\"598\" height=\"286\" /></a></td>'
317 >          <td><a href="'''+(DirKeys[dir].GetTitle())+"_Log_"+hist+'.png\"><img src=\"'+(DirKeys[dir].GetTitle())+"_Log_"+hist+'.png\" width =\"286\" height=\"286\" /></a></td>'"<td><a href=\""+"n"+(DirKeys[dir].GetTitle())+"_Log_"+hist+'.png\"><img src=\"'+"n"+(DirKeys[dir].GetTitle())+"_Log_"+hist+'.png\" width =\"286\" height=\"286\" /></a></td>'
318          Webpage.write(Images)
319        # mergedPlots.writeTObject(c1)
320        c1.cd(1).SetLogy()
# Line 236 | Line 322 | for dir in range(0,len(DirKeys)):
322        # RatioTop.Draw()
323  
324        c1.Update()
325 <      c1.SaveAs("../pdfs/"+(DirKeys[dir].GetTitle())+"_Log_"+hist+".png")
325 >      c1.SaveAs(outputfile+(DirKeys[dir].GetTitle())+"_Log_"+hist+".png")
326 >      c1.SaveAs(outputfile+(DirKeys[dir].GetTitle())+"_Log_"+hist+".pdf")
327 >
328        c1.cd(1).SetLogy(0)
329  
330      for a in closeList :

Diff Legend

Removed lines
+ Added lines
< Changed lines
> Changed lines