ATLAS Offline Software
Loading...
Searching...
No Matches
dqt_zlumi_alleff_HIST.py
Go to the documentation of this file.
1#!/usr/bin/env python
2# Copyright (C) 2002-2017 CERN for the benefit of the ATLAS collaboration
3import sys, os, glob
4import ROOT
5
6ROOT.gStyle.SetOptStat(0)
7
8#ACCEPTANCE = 3.173927e-01
9ACCEPTANCE = 3.323224e-01
10
11import argparse
12parser = argparse.ArgumentParser()
13parser.add_argument('infile', type=str, help='input HIST file')
14parser.add_argument('--out', type=str, help='output ROOT file')
15parser.add_argument('--plotdir', type=str, help='Directory to dump plots',
16 default='plots')
17parser.add_argument('--debug', action='store_true', help='Be verbose in output')
18parser.add_argument('--mode', type=str, help='Zee or Zmumu')
19
20
21args = parser.parse_args()
22
23infilename = args.infile
24infile = ROOT.TFile.Open(infilename, 'READ')
25
26runmode = args.mode
27print('Running in', runmode, 'mode')
28
29
30runname = None
31for k in infile.GetListOfKeys():
32 if k.GetName().startswith('run_'):
33 runname = k.GetName()
34 break
35if not runname:
36 print('Cannot find run directory in input file')
37 sys.exit(1)
38else:
39 print('Found runname', runname)
40
41lbdirs = []
42for k in infile.Get(runname).GetListOfKeys():
43 if k.GetName().startswith('lb_'):
44 lbdirs.append(k.GetName())
45
46print('Now to dump')
47lbnums = sorted([int(_[3:]) for _ in lbdirs])
48
49effcyt = ROOT.TH1F('effcyt', 'Trigger efficiency', lbnums[-1]-lbnums[0]+1, lbnums[0]-0.5,
50 lbnums[-1]+0.5)
51effcyr = ROOT.TH1F('effcyr', 'Loose muon reco efficiency', lbnums[-1]-lbnums[0]+1, lbnums[0]-0.5,
52 lbnums[-1]+0.5)
53effcya = ROOT.TH1F('effcya', 'Combined acc x efficiency', lbnums[-1]-lbnums[0]+1, lbnums[0]-0.5,
54 lbnums[-1]+0.5)
55effcydir = ROOT.TH1F('effcydir', 'Direct acc x efficiency', lbnums[-1]-lbnums[0]+1, lbnums[0]-0.5,
56 lbnums[-1]+0.5)
57
58from array import array
59fout = ROOT.TFile(args.out if args.out else '%s_all.root' % runname[4:], 'RECREATE')
60o_run = array('I', [0])
61o_lb = array('I', [0])
62o_lbwhen = array('d', [0., 0.])
63o_z_one = array('f', [0.])
64o_z_two = array('f', [0.])
65o_trigeff = array('f', [0.])
66o_trigeffstat = array('f', [0.])
67o_recoeff = array('f', [0.])
68o_recoeffstat = array('f', [0.])
69o_alleff = array('f', [0.])
70o_alleffstat = array('f', [0.])
71o_ae = array('f', [0.])
72o_aestat = array('f', [0.])
73tl = ROOT.TTree( 'lumitree', 'Luminosity tree' )
74tl.Branch('run', o_run, 'run/i')
75tl.Branch('lb', o_lb, 'lb/i')
76tl.Branch('lbwhen', o_lbwhen, 'lbwhen[2]/D')
77tl.Branch('z_one', o_z_one, 'z_one/F')
78tl.Branch('z_two', o_z_two, 'z_two/F')
79tl.Branch('trigeff', o_trigeff, 'trigeff/F')
80tl.Branch('trigeffstat', o_trigeffstat, 'trigeffstat/F')
81tl.Branch('recoeff', o_recoeff, 'recoeff/F')
82tl.Branch('recoeffstat', o_recoeffstat, 'recoeffstat/F')
83tl.Branch('alleff', o_alleff, 'alleff/F')
84tl.Branch('alleffstat', o_alleffstat, 'alleffstat/F')
85tl.Branch('ae', o_ae, 'ae/F')
86tl.Branch('aestat', o_aestat, 'aestat/F')
87
88
89from DQUtils import fetch_iovs
90#rset=set(_[0] for _ in rlb)
91#print(rset)
92lblb = fetch_iovs("LBLB", runs=int(runname[4:])).by_run
93for lb in sorted(lbdirs):
94 if runmode == "Zee":
95 h = infile.Get('%s/%s/GLOBAL/DQTGlobalWZFinder/m_eltrigtp_matches' % (runname, lb))
96 hmo = infile.Get('%s/%s/GLOBAL/DQTGlobalWZFinder/m_ele_tight_good_os' % (runname, lb))
97 hms = infile.Get('%s/%s/GLOBAL/DQTGlobalWZFinder/m_ele_tight_good_ss' % (runname, lb))
98 hno = infile.Get('%s/%s/GLOBAL/DQTGlobalWZFinder/m_ele_tight_bad_os' % (runname, lb))
99 hns = infile.Get('%s/%s/GLOBAL/DQTGlobalWZFinder/m_ele_tight_bad_ss' % (runname, lb))
100 if runmode == "Zmumu":
101 h = infile.Get('%s/%s/GLOBAL/DQTGlobalWZFinder/m_mutrigtp_matches' % (runname, lb))
102 hmo = infile.Get('%s/%s/GLOBAL/DQTGlobalWZFinder/m_muloosetp_match_os' % (runname, lb))
103 hms = infile.Get('%s/%s/GLOBAL/DQTGlobalWZFinder/m_muloosetp_match_ss' % (runname, lb))
104 hno = infile.Get('%s/%s/GLOBAL/DQTGlobalWZFinder/m_muloosetp_nomatch_os' % (runname, lb))
105 hns = infile.Get('%s/%s/GLOBAL/DQTGlobalWZFinder/m_muloosetp_nomatch_ss' % (runname, lb))
106 lbnum = int(lb[3:])
107 yld = (h[2], h[3])
108 ylderr = (h.GetBinError(2), h.GetBinError(3))
109 #print(yld, ylderr)
110 A, B = yld
111 o_z_one[0], o_z_two[0] = yld
112 if B == 0: continue
113 eff = 1./(float(A)/B/2.+1.)
114 inverrsq = ((1/2./B)*ylderr[0])**2+((A/2./B**2)*ylderr[1])**2
115 o_trigeff[0] = eff
116 o_trigeffstat[0] = (inverrsq**.5)*(eff**2)
117 o_run[0], o_lb[0] = int(runname[4:]), lbnum
118 try:
119 iov = lblb[int(runname[4:])][lbnum-1]
120 o_lbwhen[0], o_lbwhen[1] = iov.StartTime/1e9, iov.EndTime/1e9
121 except Exception as e:
122 o_lbwhen[0], o_lbwhen[1] = 0, 0
123 effcyt.SetBinContent(lbnum-lbnums[0]+1, eff)
124 effcyt.SetBinError(lbnum-lbnums[0]+1, o_trigeffstat[0])
125
126 def extract(histogram):
127 dbl = ROOT.Double()
128 rv1 = histogram.IntegralAndError(21, 30, dbl)
129 return (rv1, float(dbl))
130 matchos, matchoserr = extract(hmo)
131 matchss, matchsserr = extract(hms)
132 nomatchos, nomatchoserr = extract(hno)
133 nomatchss, nomatchsserr = extract(hns)
134 if args.debug:
135 print(lb)
136 print(' ->', matchos, matchoserr)
137 print(' ->', matchss, matchsserr)
138 print(' ->', nomatchos, nomatchoserr)
139 print(' ->', nomatchss, nomatchsserr)
140 A = float(matchos-matchss)
141 Aerr = (matchoserr**2+matchsserr**2)**.5
142 B = float(nomatchos-nomatchss)
143 Berr = (nomatchoserr**2+nomatchsserr**2)**.5
144 if Berr == 0: Berr = 1.
145 if A == 0 or B/A == -1:
146 eff = 1.
147 inverrsq = 1.
148 else:
149 eff = 1./(1.+B/A)
150 inverrsq = ((-B/A**2)*Aerr)**2+((1./A)*Berr)**2
151 o_recoeff[0] = eff
152 o_recoeffstat[0] = (inverrsq**.5)*(eff**2)
153 effcyr.SetBinContent(lbnum-lbnums[0]+1, eff)
154 effcyr.SetBinError(lbnum-lbnums[0]+1, o_recoeffstat[0])
155
156 o_ae[0] = ACCEPTANCE*(1-(1-o_trigeff[0])**2)*(o_recoeff[0])**2
157 o_aestat[0] = ACCEPTANCE*((o_recoeff[0]**2*2*(1-o_trigeff[0])*o_trigeffstat[0])**2+(2*o_recoeff[0]*(1-(1-o_trigeff[0])**2)*o_recoeffstat[0])**2)**.5
158 o_alleff[0] = (1-(1-o_trigeff[0])**2)*(o_recoeff[0])**2
159 o_alleffstat[0] = ((o_recoeff[0]**2*2*(1-o_trigeff[0])*o_trigeffstat[0])**2+(2*o_recoeff[0]*(1-(1-o_trigeff[0])**2)*o_recoeffstat[0])**2)**.5
160 effcya.SetBinContent(lbnum-lbnums[0]+1, o_ae[0])
161 effcya.SetBinError(lbnum-lbnums[0]+1, o_aestat[0])
162
163
164 tl.Fill()
165tl.Write()
166print('Done')
167
168c1 = ROOT.TCanvas()
169effcya.SetMarkerStyle(21)
170effcya.SetMarkerColor(ROOT.kBlue)
171effcya.GetYaxis().SetRangeUser(0.25,0.31)
172effcya.Draw('PE')
173c1.Print(os.path.join(args.plotdir, '%s_combined_efficiency.eps' % runname[4:]))
174fout.WriteTObject(effcya)
175c1.Clear()
176effcyt.SetMarkerStyle(21)
177effcyt.SetMarkerColor(ROOT.kBlue)
178effcyt.GetYaxis().SetRangeUser(0.66,0.86)
179effcyt.Draw('PE')
180c1.Print(os.path.join(args.plotdir, '%s_trigger_efficiency.eps' % runname[4:]))
181fout.WriteTObject(effcyt)
182c1.Clear()
183effcyr.SetMarkerStyle(21)
184effcyr.SetMarkerColor(ROOT.kBlue)
185effcyr.GetYaxis().SetRangeUser(0.9,1.0)
186effcyr.Draw('PE')
187c1.Print(os.path.join(args.plotdir, '%s_reco_efficiency.eps' % runname[4:]))
188fout.WriteTObject(effcyr)
189fout.Close()
190
191if runmode == "Zee":
192 sumweights = infile.Get('%s/GLOBAL/DQTDataFlow/m_sumweights' % runname)
193 ctr = infile.Get('%s/GLOBAL/DQTGlobalWZFinder/m_Z_Counter_el' % runname)
194if runmode == "Zmumu":
195 sumweights = infile.Get('%s/GLOBAL/DQTDataFlow/m_sumweights' % runname)
196 ctr = infile.Get('%s/GLOBAL/DQTGlobalWZFinder/m_Z_Counter_mu' % runname)
197if sumweights:
198 for ibin in xrange(1,sumweights.GetNbinsX()+1):
199 o_lb[0] = int(sumweights.GetBinCenter(ibin))
200 ctrbin = ctr.FindBin(o_lb[0])
201 print(ibin, o_lb[0], sumweights[ibin], ctr[ctrbin])
202 if sumweights[ibin] == 0: continue
203 p = ctr[ctrbin]/sumweights[ibin]
204 o_alleff[0]=p
205 try:
206 o_alleffstat[0]=(p*(1-p))**.5*(sumweights.GetBinError(ibin)/sumweights[ibin])
207 except ValueError:
208 o_alleffstat[0]=(sumweights.GetBinError(ibin)/sumweights[ibin])
209 effcydir.SetBinContent(effcydir.FindBin(o_lb[0]), p)
210 effcydir.SetBinError(effcydir.FindBin(o_lb[0]), o_alleffstat[0])
211
212 effcya.GetYaxis().SetRangeUser(0.27,0.31)
213 effcya.Draw('PE')
214 effcydir.SetMarkerStyle(20)
215 effcydir.SetMarkerColor(ROOT.kRed)
216 effcydir.Draw('SAME,PE')
217 leg=ROOT.TLegend(0.65, 0.7, 0.89, 0.89)
218 leg.AddEntry(effcya, 'Predicted A#epsilon', 'PE')
219 leg.AddEntry(effcydir, 'Actual A#epsilon', 'PE')
220 leg.Draw()
221 c1.Print(os.path.join(args.plotdir, '%s_tp_comparison.eps' % runname[4:]))
222
223 effcyrat=effcydir.Clone()
224 effcyrat.Divide(effcya)
225 effcyrat.SetTitle('MC Correction Factor')
226 effcyrat.SetXTitle('<#mu>')
227 effcyrat.Draw('PE')
228 effcyrat.Fit('pol1')
229 c1.Print(os.path.join(args.plotdir, '%s_tp_correction.eps' % runname[4:]))
print
void print(char *figname, TCanvas *c1)
Definition TRTCalib_StrawStatusPlots.cxx:26
array
STL class.
xrange
void xrange(TH1 *h, bool symmetric)
Definition computils.cxx:517
dqt_zlumi_alleff_HIST.extract
extract(histogram)
Definition dqt_zlumi_alleff_HIST.py:126
Generated on for ATLAS Offline Software by doxygen 1.14.0