Analysis Software
Documentation for sPHENIX simulation software
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
INTTVtxZ.py
Go to the documentation of this file. Or view the newest version in sPHENIX GitHub for file INTTVtxZ.py
1 #! /usr/bin/env python
2 from optparse import OptionParser
3 import sys
4 import os
5 import datetime
6 from array import *
7 from ROOT import *
8 import numpy as np
9 import math
10 import glob
11 from plotUtil import Draw_2Dhist
12 
13 gROOT.LoadMacro('./sPHENIXStyle/sPhenixStyle.C')
14 gROOT.ProcessLine('SetsPhenixStyle()')
15 gROOT.SetBatch(True)
16 
17 TickSize = 0.03
18 AxisTitleSize = 0.05
19 AxisLabelSize = 0.04
20 LeftMargin = 0.15
21 RightMargin = 0.08
22 TopMargin = 0.08
23 BottomMargin = 0.13
24 
25 
26 def Draw_1Dhist_fitGaussian(hist, norm1, logy, ymaxscale, XaxisName, Ytitle_unit, outname):
27  hist.Sumw2()
28  binwidth = hist.GetXaxis().GetBinWidth(1)
29  c = TCanvas('c', 'c', 800, 700)
30  if norm1:
31  hist.Scale(1. / hist.Integral(-1, -1))
32  if logy:
33  c.SetLogy()
34  c.cd()
35  gPad.SetRightMargin(RightMargin)
36  gPad.SetTopMargin(TopMargin)
37  gPad.SetLeftMargin(LeftMargin)
38  gPad.SetBottomMargin(BottomMargin)
39  if norm1:
40  if Ytitle_unit == '':
41  hist.GetYaxis().SetTitle(
42  'Normalized entries / ({:g})'.format(binwidth))
43  else:
44  hist.GetYaxis().SetTitle(
45  'Normalized entries / ({:g} {unit})'.format(binwidth, unit=Ytitle_unit))
46  else:
47  if Ytitle_unit == '':
48  hist.GetYaxis().SetTitle('Entries / ({:g})'.format(binwidth))
49  else:
50  hist.GetYaxis().SetTitle(
51  'Entries / ({:g} {unit})'.format(binwidth, unit=Ytitle_unit))
52 
53  # hist.GetXaxis().SetRangeUser(hist.GetBinLowEdge(1)-binwidth, hist.GetBinLowEdge(hist.GetNbinsX())+2*binwidth)
54  if logy:
55  hist.GetYaxis().SetRangeUser(hist.GetMinimum(0)*0.5, (hist.GetMaximum()) * ymaxscale)
56  else:
57  hist.GetYaxis().SetRangeUser(0., (hist.GetMaximum()) * ymaxscale)
58  hist.GetXaxis().SetTitle(XaxisName)
59  hist.GetXaxis().SetTickSize(TickSize)
60  hist.GetXaxis().SetTitleSize(AxisTitleSize)
61  hist.GetXaxis().SetLabelSize(AxisLabelSize)
62  hist.GetYaxis().SetTickSize(TickSize)
63  hist.GetYaxis().SetTitleSize(AxisTitleSize)
64  hist.GetYaxis().SetLabelSize(AxisLabelSize)
65  hist.GetXaxis().SetTitleOffset(1.1)
66  hist.GetYaxis().SetTitleOffset(1.35)
67  hist.SetLineColor(1)
68  hist.SetLineWidth(2)
69  hist.Draw('hist')
70  f1 = TF1('f1', 'gaus', hist.GetXaxis().GetXmin(), hist.GetXaxis().GetXmax())
71  f1.SetParameter(1,hist.GetMean())
72  f1.SetParLimits(1,-1,1)
73  f1.SetParameter(2,hist.GetRMS())
74  f1.SetParLimits(2,0,100)
75  f1.SetLineColorAlpha(TColor.GetColor('#F54748'), 0.9)
76  hist.Fit(f1, 'B')
77  f1.Draw('same')
78  leg = TLegend((1-RightMargin)-0.45, (1-TopMargin)-0.13,
79  (1-RightMargin)-0.1, (1-TopMargin)-0.03)
80  leg.SetTextSize(0.040)
81  leg.SetFillStyle(0)
82  leg.AddEntry('', '#it{#bf{sPHENIX}} Simulation', '')
83  leg.AddEntry('','Au+Au #sqrt{s_{NN}}=200 GeV','')
84  leg.Draw()
85  leg2 = TLegend(0.54, 0.65, 0.88, 0.78)
86  leg2.SetTextSize(0.033)
87  leg2.SetFillStyle(0)
88  leg2.AddEntry(f1,'Gaussian fit','l')
89  leg2.AddEntry('', '#mu={0:.2e}#pm{1:.2e}'.format(f1.GetParameter(1), f1.GetParError(1)), '')
90  leg2.AddEntry('', '#sigma={0:.2e}#pm{1:.2e}'.format(f1.GetParameter(2), f1.GetParError(2)), '')
91  leg2.Draw()
92  c.RedrawAxis()
93  c.Draw()
94  c.SaveAs(outname+'.pdf')
95  c.SaveAs(outname+'.png')
96  if(c):
97  c.Close()
98  gSystem.ProcessEvents()
99  del c
100  c = 0
101  return f1.GetParameter(1), f1.GetParError(1), f1.GetParameter(2), f1.GetParError(2)
102 
103 
104 def Draw_HistGraph(h, xaxistitle, yaxistitle, yaxisrange, xaxisbinlabel, outname):
105  c = TCanvas('c1', 'c1', 800, 700)
106  gPad.SetRightMargin(0.08)
107  c.cd()
108  h.GetXaxis().SetTitle(xaxistitle)
109  h.GetYaxis().SetTitle(yaxistitle)
110  h.SetLineWidth(2)
111  h.SetMarkerSize(1.5)
112  h.GetXaxis().SetTitleOffset(1.2)
113  h.GetXaxis().SetNdivisions(10, 0, 0, kTRUE)
114  h.GetYaxis().SetTitleOffset(1.3)
115  h.GetYaxis().SetRangeUser(yaxisrange[0], yaxisrange[1])
116  # Set bin labels
117  for i in range(len(xaxisbinlabel)):
118  h.GetXaxis().SetBinLabel(i+1, xaxisbinlabel[i])
119 
120  h.Draw('PEX0')
121 
122  leg = TLegend(1 - RightMargin - 0.45, 1-TopMargin-0.15, 1 - RightMargin - 0.08, 1-TopMargin-0.035)
123  leg.SetTextSize(0.040)
124  leg.SetFillStyle(0)
125  leg.AddEntry('', '#it{#bf{sPHENIX}} Simulation', '')
126  leg.AddEntry('','Au+Au #sqrt{s_{NN}}=200 GeV','')
127  leg.Draw()
128  c.Draw()
129  c.SaveAs(outname+'.pdf')
130  c.SaveAs(outname+'.png')
131  if(c):
132  c.Close()
133  gSystem.ProcessEvents()
134  del c
135  c = 0
136 
137 def Draw_1DEff(EffErr, logx, XaxisName, axesrange, outname):
138  c = TCanvas('c', 'c', 800, 700)
139  if logx:
140  c.SetLogx()
141  c.cd()
142  gPad.SetRightMargin(RightMargin)
143  gPad.SetTopMargin(TopMargin)
144  gPad.SetLeftMargin(LeftMargin)
145  gPad.SetBottomMargin(BottomMargin)
146  EffErr.GetXaxis().SetTitle(XaxisName)
147  EffErr.GetYaxis().SetTitle('Efficiency')
148  EffErr.GetXaxis().SetRangeUser(axesrange[0], axesrange[1])
149  # EffErr.GetXaxis().SetMoreLogLabels()
150  EffErr.GetYaxis().SetRangeUser(axesrange[2], axesrange[3])
151  EffErr.SetMarkerColor(1)
152  EffErr.SetMarkerSize(1.5)
153  EffErr.SetMarkerStyle(20)
154  EffErr.GetXaxis().SetTitleOffset(1.2)
155  # EffErr.GetYaxis().SetTitleOffset(ytitleoffset)
156  # EffErr.GetXaxis().SetNdivisions(10, 0, 0, kTRUE)
157  EffErr.Draw('AP')
158  leg = TLegend(LeftMargin+0.03, 1-TopMargin-0.15, LeftMargin+0.1, 1-TopMargin-0.035)
159  leg.SetTextSize(0.040)
160  leg.SetFillStyle(0)
161  leg.AddEntry('', '#it{#bf{sPHENIX}} Simulation', '')
162  leg.AddEntry('','Au+Au #sqrt{s_{NN}}=200 GeV','')
163  leg.Draw()
164  c.Draw()
165  c.SaveAs(outname+'.pdf')
166  c.SaveAs(outname+'.png')
167  if(c):
168  c.Close()
169  gSystem.ProcessEvents()
170  del c
171  c = 0
172 
173 
174 def Draw_1DEffComp(leff, lcolor, logx, XaxisName, NdivisionArg, legtext, axesrange, outname):
175  c = TCanvas('c', 'c', 800, 700)
176  if logx:
177  c.SetLogx()
178  c.cd()
179  gPad.SetRightMargin(RightMargin)
180  gPad.SetTopMargin(TopMargin)
181  gPad.SetLeftMargin(LeftMargin)
182  gPad.SetBottomMargin(BottomMargin)
183  for i,eff in enumerate(leff):
184  if i == 0:
185  leff[i].GetXaxis().SetTitle(XaxisName)
186  leff[i].GetXaxis().SetTitleOffset(1.2)
187  leff[i].GetYaxis().SetTitle('Efficiency')
188  leff[i].GetXaxis().SetRangeUser(axesrange[0], axesrange[1])
189  leff[i].GetYaxis().SetRangeUser(axesrange[2], axesrange[3])
190  leff[i].SetLineColor(TColor.GetColor(lcolor[i]))
191  leff[i].SetMarkerColor(TColor.GetColor(lcolor[i]))
192  leff[i].SetMarkerSize(1.3)
193  leff[i].SetMarkerStyle(20)
194  if len(NdivisionArg) != 0:
195  leff[i].GetXaxis().SetNdivisions(NdivisionArg[0], NdivisionArg[1], NdivisionArg[2], NdivisionArg[3])
196  leff[i].Draw('AP')
197  else:
198  leff[i].SetLineColor(TColor.GetColor(lcolor[i]))
199  leff[i].SetMarkerColor(TColor.GetColor(lcolor[i]))
200  leff[i].SetMarkerSize(1.3)
201  leff[i].SetMarkerStyle(20)
202  leff[i].Draw('Psame')
203 
204  leg = TLegend(LeftMargin, (1-TopMargin)-0.15,
205  LeftMargin+0.25, (1-TopMargin)-0.03)
206  leg.SetTextSize(0.04)
207  leg.SetFillStyle(0)
208  leg.AddEntry('', '#it{#bf{sPHENIX}} Simulation', '')
209  leg.AddEntry('','Au+Au #sqrt{s_{NN}}=200 GeV','')
210  leg.Draw()
211  leg1 = TLegend((1-RightMargin)-0.43, BottomMargin+0.2,
212  (1-RightMargin)-0.03, BottomMargin+0.05)
213  leg1.SetTextSize(0.035)
214  leg1.SetFillStyle(0)
215  for i, eff in enumerate(leff):
216  leg1.AddEntry(eff, legtext[i], "pe")
217  leg1.Draw()
218  c.Draw()
219  c.SaveAs(outname+'.pdf')
220  c.SaveAs(outname+'.png')
221  if(c):
222  c.Close()
223  gSystem.ProcessEvents()
224  del c
225  c = 0
226 
227 def Draw_2Dhist_eff(hist, logx, logy, logz, norm1, rmargin, XaxisName, YaxisName, drawopt, outname):
228  gStyle.SetPalette(kThermometer)
229  gStyle.SetPaintTextFormat('1.3g')
230  c = TCanvas('c', 'c', 800, 700)
231  if logx:
232  c.SetLogx()
233  if logy:
234  c.SetLogy()
235  if logz:
236  c.SetLogz()
237  c.cd()
238  gPad.SetRightMargin(rmargin)
239  gPad.SetTopMargin(TopMargin)
240  gPad.SetLeftMargin(LeftMargin)
241  gPad.SetBottomMargin(BottomMargin)
242  if norm1:
243  hist.Scale(1. / hist.Integral(-1, -1, -1, -1))
244  hist.GetXaxis().SetTitle(XaxisName)
245  hist.GetYaxis().SetTitle(YaxisName)
246  hist.GetZaxis().SetTitle('Efficiency')
247  hist.GetXaxis().SetTickSize(TickSize)
248  hist.GetYaxis().SetTickSize(TickSize)
249  hist.GetXaxis().SetTitleSize(AxisTitleSize)
250  hist.GetYaxis().SetTitleSize(AxisTitleSize)
251  hist.GetXaxis().SetLabelSize(AxisLabelSize)
252  hist.GetYaxis().SetLabelSize(AxisLabelSize)
253  hist.GetXaxis().SetTitleOffset(1.1)
254  hist.GetYaxis().SetTitleOffset(1.3)
255  hist.GetZaxis().SetTitleOffset(1.2)
256  hist.GetZaxis().SetLabelSize(AxisLabelSize)
257  hist.SetContour(1000)
258  hist.SetMarkerSize(0.9)
259  hist.Draw(drawopt)
260  c.RedrawAxis()
261  c.Draw()
262  c.SaveAs(outname+'.pdf')
263  c.SaveAs(outname+'.png')
264  if(c):
265  c.Close()
266  gSystem.ProcessEvents()
267  del c
268  c = 0
269 
270 
271 
272 if __name__ == '__main__':
273  parser = OptionParser(usage="usage: %prog ver [options -h]")
274  parser.add_option("-i", "--infiledir", dest="infiledir", default='/sphenix/user/hjheng/TrackletAna/minitree/INTT/VtxEvtMap_ana382_zvtx-20cm_dummyAlignParams', help="Input file directory")
275  parser.add_option("-o", "--outdirprefix", dest="outdirprefix", default='ana382_zvtx-20cm_dummyAlignParams', help="Output file directory")
276  parser.add_option("-p", "--dphi", dest="dphi", default=0.3, help="dphi cut [degree]")
277  parser.add_option("-d", "--dca", dest="dca", default=0.15, help="dca cut [cm]")
278 
279  (opt, args) = parser.parse_args()
280  print('opt: {}'.format(opt))
281 
282  infiledir = opt.infiledir
283  outdirprefix = opt.outdirprefix
284  dphi = float(opt.dphi)
285  dca = float(opt.dca)
286  dphitext = '{:.3f}'.format(dphi).replace('.', 'p')
287  dcadtext = '{:.3f}'.format(dca).replace('.', 'p')
288 
289  plotpath = './RecoPV_ana/'
290  os.makedirs('{}/{}'.format(plotpath,outdirprefix), exist_ok=True)
291 
292  # Input file
293  fname = '{}/dphi{}deg_dca{}cm/INTTVtxZ.root'.format(infiledir, dphitext, dcadtext)
294 
295  # of clusters in inner layer, percentile
296  df = ROOT.RDataFrame('minitree', fname)
297  nparr_NClusInner = df.AsNumpy(columns=['NClusLayer1'])
298  NclusInner_percentile = []
299  Binedge_NclusInner_percentile = [0]
300  NclusInner_percentile_cut = [0]
301  NpercentileDiv = 20
302  for i in range(NpercentileDiv-1):
303  # print('percentile={}%, Nhits={}'.format((i+1)*5, np.percentile(np_NhitsL1['NClusLayer1'], (i+1)*5)))
304  NclusInner_percentile.append(np.percentile(nparr_NClusInner['NClusLayer1'], (i+1)*5))
305  Binedge_NclusInner_percentile.append(np.percentile(nparr_NClusInner['NClusLayer1'], (i+1)*5))
306  if i % 2 == 1:
307  NclusInner_percentile_cut.append(np.percentile(nparr_NClusInner['NClusLayer1'], (i+1)*5))
308  NclusInner_percentile_cut.append(10000)
309  Binedge_NclusInner_percentile.append(10000)
310 
311  # Prepare histograms
312  l_hM_DiffVtxZ = []
313  for i in range(10):
314  if i < 7:
315  l_hM_DiffVtxZ.append(TH1F('hM_DiffVtxZ_Cent{}'.format(5*(2*i+1)), 'hM_DiffVtxZ_Cent{}'.format(5*(2*i+1)), 100, -2, 2))
316  elif i >= 7:
317  l_hM_DiffVtxZ.append(TH1F('hM_DiffVtxZ_Cent{}'.format(5*(2*i+1)), 'hM_DiffVtxZ_Cent{}'.format(5*(2*i+1)), 100, -5, 5))
318  else:
319  print('ERROR: i={}'.format(i))
320  sys.exit()
321 
322  hM_NClusInner_NtruthVtx1 = TH1F('hM_NClusInner_NtruthVtx1', 'hM_NClusInner_NtruthVtx1', NpercentileDiv, np.asarray(Binedge_NclusInner_percentile))
323  hM_NClusInner_NtruthVtx1_absdiffle5cm = TH1F('hM_NClusInner_NtruthVtx1_absdiffle5cm', 'hM_NClusInner_NtruthVtx1_absdiffle5cm', NpercentileDiv, np.asarray(Binedge_NclusInner_percentile))
324  hM_NClusInner_NtruthVtx1_absdiffle2cm = TH1F('hM_NClusInner_NtruthVtx1_absdiffle2cm', 'hM_NClusInner_NtruthVtx1_absdiffle2cm', NpercentileDiv, np.asarray(Binedge_NclusInner_percentile))
325  hM_NClusInner_NtruthVtx1_absdiffle1cm = TH1F('hM_NClusInner_NtruthVtx1_absdiffle1cm', 'hM_NClusInner_NtruthVtx1_absdiffle1cm', NpercentileDiv, np.asarray(Binedge_NclusInner_percentile))
326  hM_Centrality_NtruthVtx1 = TH1F('hM_Centrality_NtruthVtx1', 'hM_Centrality_NtruthVtx1', 10, 0, 10)
327  hM_Centrality_NtruthVtx1_absdiffle5cm = TH1F('hM_Centrality_NtruthVtx1_absdiffle5cm', 'hM_Centrality_NtruthVtx1_absdiffle5cm', 10, 0, 10)
328  hM_Centrality_NtruthVtx1_absdiffle2cm = TH1F('hM_Centrality_NtruthVtx1_absdiffle2cm', 'hM_Centrality_NtruthVtx1_absdiffle2cm', 10, 0, 10)
329  hM_Centrality_NtruthVtx1_absdiffle1cm = TH1F('hM_Centrality_NtruthVtx1_absdiffle1cm', 'hM_Centrality_NtruthVtx1_absdiffle1cm', 10, 0, 10)
330  hM_NClusInner_TruthPVz_NtruthVtx1 = TH2F('hM_NClusInner_TruthPVz_NtruthVtx1', 'hM_NClusInner_TruthPVz_NtruthVtx1', NpercentileDiv, np.asarray(Binedge_NclusInner_percentile), 20, -45, 5)
331  hM_NClusInner_TruthPVz_NtruthVtx1_absdiffle5cm = TH2F('hM_NClusInner_TruthPVz_NtruthVtx1_absdiffle5cm', 'hM_NClusInner_TruthPVz_NtruthVtx1_absdiffle5cm', NpercentileDiv, np.asarray(Binedge_NclusInner_percentile), 20, -45, 5)
332  hM_NClusInner_TruthPVz_NtruthVtx1_absdiffle2cm = TH2F('hM_NClusInner_TruthPVz_NtruthVtx1_absdiffle2cm', 'hM_NClusInner_TruthPVz_NtruthVtx1_absdiffle2cm', NpercentileDiv, np.asarray(Binedge_NclusInner_percentile), 20, -45, 5)
333  hM_NClusInner_TruthPVz_NtruthVtx1_absdiffle1cm = TH2F('hM_NClusInner_TruthPVz_NtruthVtx1_absdiffle1cm', 'hM_NClusInner_TruthPVz_NtruthVtx1_absdiffle1cm', NpercentileDiv, np.asarray(Binedge_NclusInner_percentile), 20, -45, 5)
334  hM_Centrality_TruthPVz_NtruthVtx1 = TH2F('hM_Centrality_TruthPVz_NtruthVtx1', 'hM_Centrality_TruthPVz_NtruthVtx1', 10, 0, 10, 20, -45, 5)
335  hM_Centrality_TruthPVz_NtruthVtx1_absdiffle5cm = TH2F('hM_Centrality_TruthPVz_NtruthVtx1_absdiffle5cm', 'hM_Centrality_TruthPVz_NtruthVtx1_absdiffle5cm', 10, 0, 10, 20, -45, 5)
336  hM_Centrality_TruthPVz_NtruthVtx1_absdiffle2cm = TH2F('hM_Centrality_TruthPVz_NtruthVtx1_absdiffle2cm', 'hM_Centrality_TruthPVz_NtruthVtx1_absdiffle2cm', 10, 0, 10, 20, -45, 5)
337  hM_Centrality_TruthPVz_NtruthVtx1_absdiffle1cm = TH2F('hM_Centrality_TruthPVz_NtruthVtx1_absdiffle1cm', 'hM_Centrality_TruthPVz_NtruthVtx1_absdiffle1cm', 10, 0, 10, 20, -45, 5)
338 
339  # Loop events and fill histograms
340  f = TFile(fname, 'r')
341  tree = f.Get('minitree')
342  for idx in range(tree.GetEntries()):
343  tree.GetEntry(idx)
344  if tree.NTruthVtx == 1:
345  idx = int(((tree.Centrality_mbd / 5) - 1) / 2)
346  l_hM_DiffVtxZ[idx].Fill(tree.PV_z - tree.TruthPV_z)
347  hM_NClusInner_NtruthVtx1.Fill(tree.NClusLayer1)
348  if abs(tree.PV_z - tree.TruthPV_z) <= 5:
349  hM_NClusInner_NtruthVtx1_absdiffle5cm.Fill(tree.NClusLayer1)
350  if abs(tree.PV_z - tree.TruthPV_z) <= 2:
351  hM_NClusInner_NtruthVtx1_absdiffle2cm.Fill(tree.NClusLayer1)
352  if abs(tree.PV_z - tree.TruthPV_z) <= 1:
353  hM_NClusInner_NtruthVtx1_absdiffle1cm.Fill(tree.NClusLayer1)
354 
355  centbin = int(((tree.Centrality_mbd / 5) - 1) / 2)
356  hM_Centrality_NtruthVtx1.Fill(centbin)
357  if abs(tree.PV_z - tree.TruthPV_z) <= 5:
358  hM_Centrality_NtruthVtx1_absdiffle5cm.Fill(centbin)
359  if abs(tree.PV_z - tree.TruthPV_z) <= 2:
360  hM_Centrality_NtruthVtx1_absdiffle2cm.Fill(centbin)
361  if abs(tree.PV_z - tree.TruthPV_z) <= 1:
362  hM_Centrality_NtruthVtx1_absdiffle1cm.Fill(centbin)
363 
364  hM_NClusInner_TruthPVz_NtruthVtx1.Fill(tree.NClusLayer1, tree.TruthPV_z)
365  if abs(tree.PV_z - tree.TruthPV_z) <= 5:
366  hM_NClusInner_TruthPVz_NtruthVtx1_absdiffle5cm.Fill(tree.NClusLayer1, tree.TruthPV_z)
367  if abs(tree.PV_z - tree.TruthPV_z) <= 2:
368  hM_NClusInner_TruthPVz_NtruthVtx1_absdiffle2cm.Fill(tree.NClusLayer1, tree.TruthPV_z)
369  if abs(tree.PV_z - tree.TruthPV_z) <= 1:
370  hM_NClusInner_TruthPVz_NtruthVtx1_absdiffle1cm.Fill(tree.NClusLayer1, tree.TruthPV_z)
371 
372  hM_Centrality_TruthPVz_NtruthVtx1.Fill(centbin, tree.TruthPV_z)
373  if abs(tree.PV_z - tree.TruthPV_z) <= 5:
374  hM_Centrality_TruthPVz_NtruthVtx1_absdiffle5cm.Fill(centbin, tree.TruthPV_z)
375  if abs(tree.PV_z - tree.TruthPV_z) <= 2:
376  hM_Centrality_TruthPVz_NtruthVtx1_absdiffle2cm.Fill(centbin, tree.TruthPV_z)
377  if abs(tree.PV_z - tree.TruthPV_z) <= 1:
378  hM_Centrality_TruthPVz_NtruthVtx1_absdiffle1cm.Fill(centbin, tree.TruthPV_z)
379 
380 
381  l_res_vtxZ = []
382  l_reserr_vtxZ = []
383  h_resolution_cent = TH1F('h_resolution_cent', 'h_resolution_cent', 10, 0, 10)
384  h_bias_cent = TH1F('h_bias_cent', 'h_bias_cent', 10, 0, 10)
385  for i, hist in enumerate(l_hM_DiffVtxZ):
386  mean, meanerr, sigma, sigmaerr = Draw_1Dhist_fitGaussian(hist, False, True, 50, '#Deltaz(PV_{Truth}, PV_{Reco}) [cm]', 'cm', '{}/{}/DiffVtxZ_Cent{}'.format(plotpath,outdirprefix,5*(2*i+1)))
387  l_res_vtxZ.append(sigma)
388  l_reserr_vtxZ.append(sigmaerr)
389  h_bias_cent.SetBinContent(i+1, mean)
390  h_bias_cent.SetBinError(i+1, meanerr)
391  h_resolution_cent.SetBinContent(i+1, sigma)
392  h_resolution_cent.SetBinError(i+1, sigmaerr)
393 
394  Draw_HistGraph(h_bias_cent, 'Centrality', '#mu_{#Deltaz}^{Gaussian} [cm]', [-0.35,0.35], ['{}-{}%'.format(10*i, 10*(i+1)) for i in range(len(l_reserr_vtxZ))], '{}/{}/VtxZBias_Cent'.format(plotpath,outdirprefix))
395  Draw_HistGraph(h_resolution_cent, 'Centrality', '#sigma_{#Deltaz}^{Gaussian} [cm]', [0, h_resolution_cent.GetMaximum()*1.5], ['{}-{}%'.format(10*i, 10*(i+1)) for i in range(len(l_reserr_vtxZ))], '{}/{}/VtxZReolution_Cent'.format(plotpath,outdirprefix))
396 
397  # Efficiency
398  err_RecoPVEff_NClusInner_absdiffle5cm = TGraphAsymmErrors()
399  err_RecoPVEff_NClusInner_absdiffle5cm.Divide(hM_NClusInner_NtruthVtx1_absdiffle5cm, hM_NClusInner_NtruthVtx1)
400  err_RecoPVEff_NClusInner_absdiffle2cm = TGraphAsymmErrors()
401  err_RecoPVEff_NClusInner_absdiffle2cm.Divide(hM_NClusInner_NtruthVtx1_absdiffle2cm, hM_NClusInner_NtruthVtx1)
402  err_RecoPVEff_NClusInner_absdiffle1cm = TGraphAsymmErrors()
403  err_RecoPVEff_NClusInner_absdiffle1cm.Divide(hM_NClusInner_NtruthVtx1_absdiffle1cm, hM_NClusInner_NtruthVtx1)
404  list_eff_NClusInner = [err_RecoPVEff_NClusInner_absdiffle5cm, err_RecoPVEff_NClusInner_absdiffle2cm, err_RecoPVEff_NClusInner_absdiffle1cm]
405  list_color = ['#03001C', '#035397', '#9B0000']
406  list_leg = ['|\Deltaz(PV_{Reco},PV_{Truth})|\leq5cm', '|\Deltaz(PV_{Reco},PV_{Truth})|\leq2cm', '|\Deltaz(PV_{Reco},PV_{Truth})|\leq1cm']
407  Draw_1DEffComp(list_eff_NClusInner, list_color, True, 'Number of clusters (inner)', [], list_leg, [0,10100,0,1.4], '{}/{}/RecoPVEff_NClusInner_EffComp'.format(plotpath,outdirprefix))
408 
409  err_RecoPVEff_Centrality_absdiffle5cm = TGraphAsymmErrors()
410  err_RecoPVEff_Centrality_absdiffle5cm.Divide(hM_Centrality_NtruthVtx1_absdiffle5cm, hM_Centrality_NtruthVtx1)
411  err_RecoPVEff_Centrality_absdiffle5cm.GetXaxis().SetMaxDigits(2)
412  err_RecoPVEff_Centrality_absdiffle2cm = TGraphAsymmErrors()
413  err_RecoPVEff_Centrality_absdiffle2cm.Divide(hM_Centrality_NtruthVtx1_absdiffle2cm, hM_Centrality_NtruthVtx1)
414  err_RecoPVEff_Centrality_absdiffle1cm = TGraphAsymmErrors()
415  err_RecoPVEff_Centrality_absdiffle1cm.Divide(hM_Centrality_NtruthVtx1_absdiffle1cm, hM_Centrality_NtruthVtx1)
416  list_eff_Centrality = [err_RecoPVEff_Centrality_absdiffle5cm, err_RecoPVEff_Centrality_absdiffle2cm, err_RecoPVEff_Centrality_absdiffle1cm]
417  Draw_1DEffComp(list_eff_Centrality, list_color, False, 'Centrality bin', [-10,0,0,kFALSE], list_leg, [0,10,0,1.4], '{}/{}/RecoPVEff_CentralityBin_EffComp'.format(plotpath,outdirprefix))
418 
419  # 2D efficiency
420  eff_RecoPVEff_NClusInner_TruthPVz_absdiffle5cm = hM_NClusInner_TruthPVz_NtruthVtx1_absdiffle5cm.Clone()
421  eff_RecoPVEff_NClusInner_TruthPVz_absdiffle5cm.Divide(hM_NClusInner_TruthPVz_NtruthVtx1)
422  eff_RecoPVEff_NClusInner_TruthPVz_absdiffle2cm = hM_NClusInner_TruthPVz_NtruthVtx1_absdiffle2cm.Clone()
423  eff_RecoPVEff_NClusInner_TruthPVz_absdiffle2cm.Divide(hM_NClusInner_TruthPVz_NtruthVtx1)
424  eff_RecoPVEff_NClusInner_TruthPVz_absdiffle1cm = hM_NClusInner_TruthPVz_NtruthVtx1_absdiffle1cm.Clone()
425  eff_RecoPVEff_NClusInner_TruthPVz_absdiffle1cm.Divide(hM_NClusInner_TruthPVz_NtruthVtx1)
426  # Draw_2Dhist_eff(hist, logx, logy, logz, norm1, rmargin, XaxisName, YaxisName, drawopt, outname)
427  Draw_2Dhist_eff(eff_RecoPVEff_NClusInner_TruthPVz_absdiffle5cm, True, False, False, False, 0.18, 'Number of clusters (inner)', 'Truth PV_{z} [cm]', 'colz text45', '{}/{}/RecoPVEff_NClusInner_TruthPVz_absdiffle5cm'.format(plotpath,outdirprefix))
428  Draw_2Dhist_eff(eff_RecoPVEff_NClusInner_TruthPVz_absdiffle2cm, True, False, False, False, 0.18, 'Number of clusters (inner)', 'Truth PV_{z} [cm]', 'colz text45', '{}/{}/RecoPVEff_NClusInner_TruthPVz_absdiffle2cm'.format(plotpath,outdirprefix))
429  Draw_2Dhist_eff(eff_RecoPVEff_NClusInner_TruthPVz_absdiffle1cm, True, False, False, False, 0.18, 'Number of clusters (inner)', 'Truth PV_{z} [cm]', 'colz text45', '{}/{}/RecoPVEff_NClusInner_TruthPVz_absdiffle1cm'.format(plotpath,outdirprefix))
430 
431  eff_RecoPVEff_Centrality_TruthPVz_absdiffle5cm = hM_Centrality_TruthPVz_NtruthVtx1_absdiffle5cm.Clone()
432  eff_RecoPVEff_Centrality_TruthPVz_absdiffle5cm.Divide(hM_Centrality_TruthPVz_NtruthVtx1)
433  eff_RecoPVEff_Centrality_TruthPVz_absdiffle2cm = hM_Centrality_TruthPVz_NtruthVtx1_absdiffle2cm.Clone()
434  eff_RecoPVEff_Centrality_TruthPVz_absdiffle2cm.Divide(hM_Centrality_TruthPVz_NtruthVtx1)
435  eff_RecoPVEff_Centrality_TruthPVz_absdiffle1cm = hM_Centrality_TruthPVz_NtruthVtx1_absdiffle1cm.Clone()
436  eff_RecoPVEff_Centrality_TruthPVz_absdiffle1cm.Divide(hM_Centrality_TruthPVz_NtruthVtx1)
437  Draw_2Dhist_eff(eff_RecoPVEff_Centrality_TruthPVz_absdiffle5cm, False, False, False, False, 0.18, 'Centrality bin', 'Truth PV_{z} [cm]', 'colz text45', '{}/{}/RecoPVEff_Centrality_TruthPVz_absdiffle5cm'.format(plotpath,outdirprefix))
438  Draw_2Dhist_eff(eff_RecoPVEff_Centrality_TruthPVz_absdiffle2cm, False, False, False, False, 0.18, 'Centrality bin', 'Truth PV_{z} [cm]', 'colz text45', '{}/{}/RecoPVEff_Centrality_TruthPVz_absdiffle2cm'.format(plotpath,outdirprefix))
439  Draw_2Dhist_eff(eff_RecoPVEff_Centrality_TruthPVz_absdiffle1cm, False, False, False, False, 0.18, 'Centrality bin', 'Truth PV_{z} [cm]', 'colz text45', '{}/{}/RecoPVEff_Centrality_TruthPVz_absdiffle1cm'.format(plotpath,outdirprefix))
440