CMS 3D CMS Logo

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