CMS 3D CMS Logo

compare.py
Go to the documentation of this file.
1 from officialStyle import officialStyle
2 from array import array
3 from ROOT import gROOT, gStyle, TH1F, TH1D, TF1, TFile, TCanvas, TH2F, TLegend, TGraphAsymmErrors, Double, TLatex
4 import os, copy, sys
5 
6 gROOT.SetBatch(True)
7 officialStyle(gStyle)
8 gStyle.SetOptTitle(0)
9 #set_palette("color")
10 #gStyle.SetPaintTextFormat("2.0f")
11 
12 
13 argvs = sys.argv
14 argc = len(argvs)
15 
16 if argc != 2:
17  print 'Please specify the runtype : python tauPOGplot.py <ZTT, ZEE, ZMM, QCD>'
18  sys.exit(0)
19 
20 runtype = argvs[1]
21 print 'You selected', runtype
22 
23 
24 tlabel = 'Z #rightarrow #tau#tau'
25 xlabel = 'gen. tau p_{T}^{vis} (GeV)'
26 
27 if runtype == 'QCD':
28  tlabel = 'QCD'
29  xlabel = 'jet p_{T} (GeV)'
30 elif runtype == 'ZEE':
31  tlabel = 'Z #rightarrow ee'
32  xlabel = 'electron p_{T} (GeV)'
33 elif runtype == 'ZMM':
34  tlabel = 'Z #rightarrow #mu#mu'
35  xlabel = 'muon p_{T} (GeV)'
36 
37 
38 
39 
40 def ensureDir(directory):
41  if not os.path.exists(directory):
42  os.makedirs(directory)
43 
44 def save(canvas, name):
45  ensureDir('compare_' + runtype)
46  canvas.SaveAs(name.replace(' ','').replace('&&','')+'.pdf')
47  canvas.SaveAs(name.replace(' ','').replace('&&','')+'.gif')
48 
49 
50 def LegendSettings(leg, ncolumn):
51  leg.SetNColumns(ncolumn)
52  leg.SetBorderSize(0)
53  leg.SetFillColor(10)
54  leg.SetLineColor(0)
55  leg.SetFillStyle(0)
56  leg.SetTextSize(0.03)
57  leg.SetTextFont(42)
58 
59 
60 
61 def makeCompareVars(tree, var, sel, leglist, nbin, xmin, xmax, xtitle, ytitle, scale, name):
62 
63 # print leglist
64 
65  c = TCanvas()
66 
67  hists = []
68  col = [1,2,4,8,6]
69  ymax = 0
70 
71  for ii, isel in enumerate(sel):
72  hist = TH1F('h_' + str(ii), 'h_' + str(ii), nbin, xmin, xmax)
73  hist.GetXaxis().SetTitle(xtitle)
74  hist.GetYaxis().SetTitle(ytitle)
75  hist.GetYaxis().SetNdivisions(507)
76  hist.SetLineColor(col[ii])
77  hist.SetLineWidth(len(sel)+1-ii)
78  hist.SetLineStyle(1)
79  hist.SetMarkerSize(0)
80  hist.SetMinimum(0)
81  hist.Sumw2()
82 
83 # print hist.GetName(), var, isel
84  tree.Project(hist.GetName(), var, isel)
85  hist.Scale(1./hist.GetEntries())
86 
87  if ymax < hist.GetMaximum():
88  ymax = hist.GetMaximum()
89 
90  hists.append(hist)
91 
92 
93  leg = TLegend(0.6,0.65,0.91,0.9)
94  LegendSettings(leg,1)
95 
96  for ii, ihist in enumerate(hists):
97  ihist.SetMaximum(ymax*1.2)
98  ihist.SetMinimum(0.)
99 
100  if ii ==0:
101  ihist.Draw('h')
102  else:
103  ihist.Draw('hsame')
104 
105  if leglist[ii] != 'None':
106  leg.AddEntry(ihist, leglist[ii], "l")
107 
108 
109  if leglist[0]!='None':
110  leg.Draw()
111 
112 
113 # save(c, 'compare_' + runtype + '/compare_' + name)
114 
115 
116 
117 
118 
119 def overlay(hists, ytitle, header, addon):
120 
121  print 'number of histograms = ', len(hists)
122 
123  canvas = TCanvas()
124  leg = TLegend(0.2,0.7,0.5,0.9)
125  LegendSettings(leg, 1)
126 
127  col = [1,2,4,6,8,9,12]
128 
129  ymax = -1
130  ymin = 100
131 
132  for ii, hist in enumerate(hists):
133  hist.GetYaxis().SetTitle('efficiency')
134  hist.SetLineColor(col[ii])
135  hist.SetMarkerColor(col[ii])
136  hist.SetLineWidth(2)
137  hist.SetMarkerSize(1)
138 
139 
140  for ip in range(hist.GetN()):
141  x = Double(-1)
142  y = Double(-1)
143  hist.GetPoint(ip, x, y)
144 
145  if ymin > y:
146  ymin = y
147  if ymax < y:
148  ymax = y
149 
150 #
151 # if ymax < hist.GetMaximum():
152 # ymax = hist.GetMaximum()
153 # if ymin > hist.GetMinimum():
154 # ymin = hist.GetMinimum()
155 
156  if ii==0:
157  hist.Draw("Azp")
158  else:
159 
160  hist.Draw("pzsame")
161 
162 # print hist.GetName(), hist.GetTitle()
163  legname = hist.GetName()
164 
165  leg.AddEntry(hist, legname, 'lep')
166 
167 
168  for hist in hists:
169  hist.SetMaximum(ymax*2)
170 # hist.SetMinimum(ymin*0.5)
171 
172  leg.Draw()
173 
174  tex = TLatex( hists[-1].GetXaxis().GetXmin() + 0.01*(hists[-1].GetXaxis().GetXmax() - hists[-1].GetXaxis().GetXmin()), ymax*2.1, addon.replace('tau_',''))
175 
176  tex.SetTextFont(42)
177  tex.SetTextSize(0.03)
178  tex.Draw()
179 
180  tex2 = TLatex( hists[-1].GetXaxis().GetXmin() + 0.87*(hists[-1].GetXaxis().GetXmax() - hists[-1].GetXaxis().GetXmin()), ymax*2.1, tlabel)
181 
182  tex2.SetTextFont(42)
183  tex2.SetTextSize(0.03)
184  tex2.Draw()
185 
186 
187 
188  save(canvas, 'compare_' + runtype + '/' + header)
189 
190 
191 
192 def hoverlay(hists, xtitle, ytitle, name):
193 
194  c = TCanvas()
195 
196  ymax = 0
197  for hist in hists:
198  if ymax < hist.GetMaximum():
199  ymax = hist.GetMaximum()
200 
201 
202  leg = TLegend(0.6,0.65,0.91,0.9)
203  LegendSettings(leg,1)
204 
205  for ii, ihist in enumerate(hists):
206  ihist.SetMaximum(ymax*1.2)
207  ihist.SetMinimum(0.)
208  ihist.SetMarkerSize(0.)
209  ihist.GetXaxis().SetTitle(xtitle)
210  ihist.GetYaxis().SetTitle(ytitle)
211 
212  if ii ==0:
213  ihist.Draw('h')
214  else:
215  ihist.Draw('hsame')
216 
217  leg.AddEntry(ihist, ihist.GetName(), "l")
218 
219 
220  leg.Draw()
221 
222 
223  save(c, 'compare_' + runtype + '/hist_' + name)
224 
225 
226 
227 
228 def makeEffPlotsVars(tree, varx, vary, sel, nbinx, xmin, xmax, nbiny, ymin, ymax, xtitle, ytitle, leglabel = None, header='', addon='', option='pt', marker=20):
229 
230  binning = [20,30,40,50,60,70,80,100,150,200]
231 
232  c = TCanvas()
233 
234  if option=='pt':
235  _hist_ = TH1F('h_effp_' + addon, 'h_effp' + addon, len(binning)-1, array('d',binning))
236  _ahist_ = TH1F('ah_effp_' + addon, 'ah_effp' + addon, len(binning)-1, array('d',binning))
237  elif option=='eta':
238  _hist_ = TH1F('h_effp_' + addon, 'h_effp' + addon, nbinx, xmin, xmax)
239  _ahist_ = TH1F('ah_effp_' + addon, 'ah_effp' + addon, nbinx, xmin, xmax)
240  elif option=='nvtx':
241  _hist_ = TH1F('h_effp_' + addon, 'h_effp' + addon, len(vbinning)-1, array('d',vbinning))
242  _ahist_ = TH1F('ah_effp_' + addon, 'ah_effp' + addon, len(vbinning)-1, array('d',vbinning))
243 
244 
245  tree.Draw(varx + ' >> ' + _hist_.GetName(), sel)
246  tree.Draw(varx + ' >> ' + _ahist_.GetName(), sel + ' && ' + vary)
247 
248  g_efficiency = TGraphAsymmErrors()
249  g_efficiency.BayesDivide(_ahist_, _hist_)
250  g_efficiency.GetXaxis().SetTitle(xtitle)
251  g_efficiency.GetYaxis().SetTitle('efficiency')
252  g_efficiency.GetYaxis().SetNdivisions(507)
253  g_efficiency.SetLineWidth(3)
254  g_efficiency.SetName(header)
255  g_efficiency.SetMinimum(0.)
256  g_efficiency.GetYaxis().SetTitleOffset(1.3)
257  g_efficiency.SetMarkerStyle(marker)
258  g_efficiency.SetMarkerSize(1)
259  g_efficiency.Draw('ap')
260 
261 # save(c, 'plots/' + addon)
262 
263  return copy.deepcopy(g_efficiency)
264 
265 
266 
267 
268 
269 if __name__ == '__main__':
270 
271 
272  vardict = {
273  'againstMuonLoose3':{'var':'tau_againstMuonLoose3 > 0.5', 'nbin':2, 'min':-0.5, 'max':1.5, 'title':'againstMuonLoose3'},
274  'againstMuonTight3':{'var':'tau_againstMuonTight3 > 0.5', 'nbin':2, 'min':-0.5, 'max':1.5, 'title':'againstMuonTight3'},
275 
276  'againstElectronVLooseMVA5':{'var':'tau_againstElectronVLooseMVA5 > 0.5', 'nbin':2, 'min':-0.5, 'max':1.5, 'title':'againstElectronVLooseMVA5'},
277  'againstElectronLooseMVA5':{'var':'tau_againstElectronLooseMVA5 > 0.5', 'nbin':2, 'min':-0.5, 'max':1.5, 'title':'againstElectronLooseMVA5'},
278  'againstElectronMediumMVA5':{'var':'tau_againstElectronMediumMVA5 > 0.5', 'nbin':2, 'min':-0.5, 'max':1.5, 'title':'againstElectronMediumMVA5'},
279 
280  'againstElectronVLooseMVA6':{'var':'tau_againstElectronVLooseMVA6 > 0.5', 'nbin':2, 'min':-0.5, 'max':1.5, 'title':'againstElectronVLooseMVA6'},
281  'againstElectronLooseMVA6':{'var':'tau_againstElectronLooseMVA6 > 0.5', 'nbin':2, 'min':-0.5, 'max':1.5, 'title':'againstElectronLooseMVA6'},
282  'againstElectronMediumMVA6':{'var':'tau_againstElectronMediumMVA6 > 0.5', 'nbin':2, 'min':-0.5, 'max':1.5, 'title':'againstElectronMediumMVA6'},
283 
284  'byLoosePileupWeightedIsolation3Hits':{'var':'tau_byLoosePileupWeightedIsolation3Hits > 0.5 && tau_decayModeFindingOldDMs > 0.5', 'nbin':2, 'min':-0.5, 'max':1.5, 'title':'byLoosePileupWeightedIsolation3Hits'},
285  'byMediumPileupWeightedIsolation3Hits':{'var':'tau_byMediumPileupWeightedIsolation3Hits > 0.5 && tau_decayModeFindingOldDMs > 0.5', 'nbin':2, 'min':-0.5, 'max':1.5, 'title':'byMediumPileupWeightedIsolation3Hits'},
286  'byTightPileupWeightedIsolation3Hits':{'var':'tau_byTightPileupWeightedIsolation3Hits > 0.5 && tau_decayModeFindingOldDMs > 0.5', 'nbin':2, 'min':-0.5, 'max':1.5, 'title':'byTightPileupWeightedIsolation3Hits'},
287 
288  'byLooseCombinedIsolationDeltaBetaCorr3Hits':{'var':'tau_byLooseCombinedIsolationDeltaBetaCorr3Hits > 0.5 && tau_decayModeFindingOldDMs > 0.5', 'nbin':2, 'min':-0.5, 'max':1.5, 'title':'byLooseCombinedIsolationDeltaBetaCorr3Hits'},
289  'byMediumCombinedIsolationDeltaBetaCorr3Hits':{'var':'tau_byMediumCombinedIsolationDeltaBetaCorr3Hits > 0.5 && tau_decayModeFindingOldDMs > 0.5', 'nbin':2, 'min':-0.5, 'max':1.5, 'title':'byMediumCombinedIsolationDeltaBetaCorr3Hits'},
290  'byTightCombinedIsolationDeltaBetaCorr3Hits':{'var':'tau_byTightCombinedIsolationDeltaBetaCorr3Hits > 0.5 && tau_decayModeFindingOldDMs > 0.5', 'nbin':2, 'min':-0.5, 'max':1.5, 'title':'byTightCombinedIsolationDeltaBetaCorr3Hits'},
291 
292  'byLooseIsolationMVA3oldDMwLT':{'var':'tau_byLooseIsolationMVA3oldDMwLT > 0.5 && tau_decayModeFindingOldDMs > 0.5', 'nbin':2, 'min':-0.5, 'max':1.5, 'title':'byLooseIsolationMVA3oldDMwLT'},
293  'byMediumIsolationMVA3oldDMwLT':{'var':'tau_byMediumIsolationMVA3oldDMwLT > 0.5 && tau_decayModeFindingOldDMs > 0.5', 'nbin':2, 'min':-0.5, 'max':1.5, 'title':'byMediumIsolationMVA3oldDMwLT'},
294  'byTightIsolationMVA3oldDMwLT':{'var':'tau_byTightIsolationMVA3oldDMwLT > 0.5 && tau_decayModeFindingOldDMs > 0.5', 'nbin':2, 'min':-0.5, 'max':1.5, 'title':'byTightIsolationMVA3oldDMwLT'},
295 
296  'byLooseIsolationMVArun2v1DBoldDMwLT':{'var':'tau_byLooseIsolationMVArun2v1DBoldDMwLT > 0.5 && tau_decayModeFindingOldDMs > 0.5', 'nbin':2, 'min':-0.5, 'max':1.5, 'title':'byLooseIsolationMVArun2v1DBoldDMwLT'},
297  'byMediumIsolationMVArun2v1DBoldDMwLT':{'var':'tau_byMediumIsolationMVArun2v1DBoldDMwLT > 0.5 && tau_decayModeFindingOldDMs > 0.5', 'nbin':2, 'min':-0.5, 'max':1.5, 'title':'byMediumIsolationMVArun2v1DBoldDMwLT'},
298  'byTightIsolationMVArun2v1DBoldDMwLT':{'var':'tau_byTightIsolationMVArun2v1DBoldDMwLT > 0.5 && tau_decayModeFindingOldDMs > 0.5', 'nbin':2, 'min':-0.5, 'max':1.5, 'title':'byTightIsolationMVArun2v1DBoldDMwLT'},
299 
300  }
301 
302 
303  reco_cut = 'tau_pt > 20 && abs(tau_eta) < 2.3'
304  loose_id = 'tau_decayModeFindingOldDMs > 0.5 && tau_byLooseCombinedIsolationDeltaBetaCorr3Hits > 0.5'
305 
306  sampledict = {
307 # '7_6_0_pre7':{'file':'Myroot_7_6_0_pre7_' + runtype + '.root', 'col':2, 'marker':20, 'width':3},
308 # '7_6_0':{'file':'Myroot_7_6_0_' + runtype + '.root', 'col':1, 'marker':21, 'width':4},
309 # '7_6_1':{'file':'Myroot_7_6_1_' + runtype + '.root', 'col':4, 'marker':22, 'width':2},
310  '7_6_1_v3':{'file':'Myroot_7_6_1_v3_' + runtype + '.root', 'col':3, 'marker':23, 'width':1},
311  }
312 
313  for hname, hdict in sorted(vardict.iteritems()):
314 
315  hists = []
316 
317  for rel, rdict in sorted(sampledict.iteritems()):
318 
319  if rel.find('7_6_1')==-1 and (hname.find('MVA6')!=-1 or hname.find('MVArun2')!=-1): continue
320 
321  tfile = TFile(rdict['file'])
322  tree = tfile.Get('per_tau')
323 
324  num_sel = reco_cut
325  den_sel = '1'
326 
327  if hname.find('against')!=-1:
328  num_sel = '1'
329  den_sel = reco_cut + ' && ' + loose_id
330 
331  hists.append(makeEffPlotsVars(tree, 'tau_genpt', num_sel + '&&' + hdict['var'], den_sel, 30, 0, 300, hdict['nbin'], hdict['min'], hdict['max'], xlabel, hdict['title'], '', rel, rel, 'pt', rdict['marker']))
332 
333 
334 # if rel=='7_6_1' and (hname.find('MVA5')!=-1 or hname.find('IsolationMVA3')!=-1):
335 # xvar = hdict['var'].replace('IsolationMVA3', 'IsolationMVArun2v1DB').replace('MVA5','MVA6')
336 # print 'adding', xvar
337 #
338 # hists.append(makeEffPlotsVars(tree, 'tau_genpt', num_sel + '&&' + xvar, den_sel, 30, 0, 300, hdict['nbin'], hdict['min'], hdict['max'], xlabel, hdict['title'], '', rel + '(' + xvar.replace('tau_','').replace('> 0.5','').replace(' && decayModeFindingOldDMs ','') + ')', rel, 'pt', rdict['marker']))
339 
340 
341 
342  overlay(hists, hname, hname, hdict['title'])
343 
344 
345 
346 
347  hvardict = {
348  'tau_dm':{'var':'tau_dm', 'nbin':12, 'min':0., 'max':12, 'title':'decay Mode', 'sel':'1'},
349  'tau_mass_1prong':{'var':'tau_mass', 'nbin':30, 'min':0., 'max':2.5, 'title':'Tau mass, 1prong', 'sel':'tau_dm==0'},
350  'tau_mass_1prongp0':{'var':'tau_mass', 'nbin':30, 'min':0., 'max':2.5, 'title':'Tau mass, 1prong+#pi^{0}', 'sel':'tau_dm==1'},
351  'tau_mass_2prong':{'var':'tau_mass', 'nbin':30, 'min':0., 'max':2.5, 'title':'Tau mass, 2prong', 'sel':'(tau_dm==5 || tau_dm==6)'},
352  'tau_mass_3prong':{'var':'tau_mass', 'nbin':30, 'min':0., 'max':2.5, 'title':'Tau mass, 3prong (+#pi^{0})', 'sel':'(tau_dm==10 || tau_dm==11)'},
353  'pt_resolution_1prong':{'var':'(tau_genpt-tau_pt)/(tau_genpt)', 'nbin':30, 'min':-1., 'max':1., 'title':'pT resolution, 1prong', 'sel':'tau_dm==0'},
354  'pt_resolution_1prongp0':{'var':'(tau_genpt-tau_pt)/(tau_genpt)', 'nbin':30, 'min':-1., 'max':1., 'title':'pT resolution, 1prong+#pi^{0}', 'sel':'tau_dm==1'},
355  'pt_resolution_2prong':{'var':'(tau_genpt-tau_pt)/(tau_genpt)', 'nbin':30, 'min':-1., 'max':1., 'title':'pT resolution, 2prong', 'sel':'(tau_dm==5 || tau_dm==6)'},
356  'pt_resolution_3prong':{'var':'(tau_genpt-tau_pt)/(tau_genpt)', 'nbin':30, 'min':-1., 'max':1., 'title':'pT resolution, 3prong (+#pi^{0})', 'sel':'(tau_dm==10 || tau_dm==11)'},
357  }
358 
359 
360  for hname, hdict in sorted(hvardict.iteritems()):
361 
362  hists = []
363 
364  if runtype != 'ZTT' and hname.find('pt_resolution')!=-1: continue
365 
366 
367  for rel, rdict in sorted(sampledict.iteritems()):
368 
369  tfile = TFile(rdict['file'])
370  tree = tfile.Get('per_tau')
371 
372  hist = TH1F('h_' + hname + '_' + rel, 'h_' + hname + '_' + rel, hdict['nbin'], hdict['min'], hdict['max'])
373  hist.GetYaxis().SetNdivisions(507)
374  hist.SetLineColor(rdict['col'])
375  hist.SetLineWidth(rdict['width'])
376  hist.SetMinimum(0)
377  hist.SetName(rel)
378  hist.Sumw2()
379 
380  tree.Project(hist.GetName(), hdict['var'], hdict['sel'])
381  hist.Scale(1./hist.GetEntries())
382 
383  hists.append(copy.deepcopy(hist))
384 
385  hoverlay(hists, hdict['title'], 'a.u.', hname)
def replace(string, replacements)
def save(canvas, name)
Definition: compare.py:44
def LegendSettings(leg, ncolumn)
Definition: compare.py:50
def overlay(hists, ytitle, header, addon)
Definition: compare.py:119
def hoverlay(hists, xtitle, ytitle, name)
Definition: compare.py:192
def makeEffPlotsVars(tree, varx, vary, sel, nbinx, xmin, xmax, nbiny, ymin, ymax, xtitle, ytitle, leglabel=None, header='', addon='', option='pt', marker=20)
Definition: compare.py:228
def ensureDir(directory)
Definition: compare.py:40
def makeCompareVars(tree, var, sel, leglist, nbin, xmin, xmax, xtitle, ytitle, scale, name)
Definition: compare.py:61