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