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