ATLAS Offline Software
yearwise_efficiency.py
Go to the documentation of this file.
1 #!/usr/bin/env python
2 # Copyright (C) 2002-2024 CERN for the benefit of the ATLAS collaboration
3 
4 """
5 Plot trigger and reconstruction efficiencies over entire data-periods.
6 """
7 
8 import numpy as np
9 import pandas as pd
10 import ROOT as R
11 import python_tools as pt
12 import ZLumiScripts.tools.zlumi_mc_cf as dq_cf
13 import math
14 from array import array
15 import time
16 import argparse
17 
18 parser = argparse.ArgumentParser()
19 parser.add_argument('--year', type=str, help='15-18, all for full Run-2')
20 parser.add_argument('--channel', type=str, help='Zee or Zmumu')
21 parser.add_argument('--indir', type=str, help='Input directory for CSV files')
22 parser.add_argument('--outdir', type=str, help='Output directory for plots')
23 parser.add_argument('--dir_2022', type=str, help='Input directory for 2022 data')
24 parser.add_argument('--dir_2023', type=str, help='Input directory for 2023 data')
25 
26 args = parser.parse_args()
27 year = args.year
28 channel = args.channel
29 indir = args.indir
30 outdir = args.outdir
31 dir_2022 = args.dir_2022
32 dir_2023 = args.dir_2023
33 
34 if year == "run3":
35  years = ["22", "23"]
36  out_tag = "_run3"
37  time_format = "%m/%y"
38  xtitle = 'Month / Year'
39  date_tag = "Run 3, #sqrt{s} = 13.6 TeV"
40  norm_type = "Run3"
41  if channel is not None:
42  xval = 0.30
43  yval = 0.33
44  else:
45  xval = 0.43
46  yval = 0.33
47  set_size = 1
48 else:
49  years = [year]
50  out_tag = year
51  time_format = "%d/%m"
52  ymin, ymax = 0.5, 1.1
53  xtitle = 'Date in 20' + year
54  date_tag = "Data 20" + year + ", #sqrt{s} = 13.6 TeV"
55  norm_type = "year"
56  xval = 0.235
57  yval = 0.86
58  set_size = 0
59 
60 def main():
61  plot_efficiency(channel, years)
62  plot_efficiency_comb(channel, years)
63 
64 def plot_efficiency_comb(channel, years):
65 
66  arr_date = []
67  arr_combeff = []
68  arr_comberr = []
69  run_num = []
70 
71  for year in years:
72  print("year = ", year)
73  grl = pt.get_grl(year)
74 
75  if year == "23":
76  maindir = args.indir + dir_2023
77  print("2023 grl = ", grl)
78 
79  elif year == "22":
80  maindir = args.indir + dir_2022
81  print("2022 grl = ", grl)
82 
83  for run in grl:
84 
85  print('Begin Run ', run, 'processing')
86 
87  dfz = pd.read_csv(maindir + "run_" + run + ".csv")
88  dfz_small = dfz
89  dfz_small['ZLumi'] = dfz_small[channel + 'Lumi']
90  dfz_small['ZLumiErr'] = dfz_small[channel + 'LumiErr']
91  dfz_small['CombEff'] = dfz_small[channel + 'EffComb']
92  dfz_small['CombErr'] = dfz_small[channel + 'ErrComb']
93  dfz_small['LBLive'] = dfz_small['LBLive']
94  dfz_small['OffMu'] = dfz_small['OffMu']
95  dfz_small = dfz_small.drop(dfz_small[dfz_small.ZLumi == 0].index)
96  dfz_small = dfz_small.drop(dfz_small[(dfz_small['LBLive']<10) | (dfz_small['PassGRL']==0)].index)
97 
98  # Cut out all runs shorter than 40 minutes
99  if dfz_small['LBLive'].sum()/60 < 40:
100  print("Skip Run", run, "because of live time", dfz_small['LBLive'].sum()/60, "min")
101  continue
102 
103  # Scale event-level efficiency with FMC
104  campaign = "mc23a"
105  dfz_small['CombEff'] *= dq_cf.correction(dfz_small['OffMu'], channel, campaign, int(run))
106  dfz_small['CombErr'] *= dq_cf.correction(dfz_small['OffMu'], channel, campaign, int(run))
107 
108  # Grab start of the run for plotting later on
109  run_start = dfz_small['LBStart'].iloc[0]
110  timestamp = time.gmtime(run_start)
111  timestamp = R.TDatime(timestamp[0], timestamp[1], timestamp[2], timestamp[3], timestamp[4], timestamp[5])
112  timestamp = timestamp.Convert()
113 
114  # Calculate average event-level efficiency
115  dfz_small['CombEff'] *= dfz_small['LBLive']
116  total_time = dfz_small['LBLive'].sum()
117  comb_eff_avg = dfz_small['CombEff'].sum()/total_time
118 
119  # Calculate average trigger efficiency error
120  dfz_small['CombErr'] *= dfz_small['LBLive']
121  dfz_small['CombErr'] *= dfz_small['CombErr']
122  total_time = dfz_small['LBLive'].sum()
123  comb_err_avg = math.sqrt(dfz_small['CombErr'].sum())/total_time
124 
125  arr_date.append(timestamp)
126  arr_combeff.append(comb_eff_avg)
127  arr_comberr.append(comb_err_avg)
128  run_num.append(run)
129 
130  arr_date = array('d', arr_date)
131 
132  arr_combeff = np.array(arr_combeff)
133  arr_comberr = np.array(arr_comberr)
134 
135  if channel == "Zee":
136  lep = "e"
137  channel_string = "Z #rightarrow ee"
138  ymin, ymax = 0.56, 0.74
139  elif channel == "Zmumu":
140  lep = "#mu"
141  channel_string = "Z #rightarrow #mu#mu"
142  ymin, ymax = 0.74, 0.80
143 
144  comb_graph = R.TGraphErrors(len(arr_date), arr_date, arr_combeff, R.nullptr,arr_comberr)
145  comb_graph.GetHistogram().SetYTitle("Efficiency")
146  comb_graph.GetHistogram().GetYaxis().SetRangeUser(ymin, ymax)
147  comb_graph.GetXaxis().SetTimeDisplay(2)
148  comb_graph.GetXaxis().SetNdivisions(9,R.kFALSE)
149  comb_graph.GetXaxis().SetTimeFormat(time_format)
150  comb_graph.GetXaxis().SetTimeOffset(0,"gmt")
151  comb_graph.SetMarkerSize(1)
152 
153  c1 = R.TCanvas()
154 
155  comb_graph.Draw("ap")
156 
157  if channel == "Zee":
158 
159  leg = R.TLegend(0.645, 0.4, 0.805, 0.6)
160  pt.drawAtlasLabel(0.2, ymax-0.06, "Internal")
161  if year in ['15', '16', '17', '18']:
162  pt.drawText(0.2, ymax-0.46, date_tag)
163  else:
164  pt.drawText(0.2, ymax-0.46, date_tag)
165  pt.drawText(0.2, ymax-0.52, channel_string + " counting")
166 
167  elif channel == "Zmumu":
168 
169  leg = R.TLegend(0.645, 0.4, 0.805, 0.6)
170  pt.drawAtlasLabel(0.2, ymax-0.4, "Internal")
171  if year in ['15', '16', '17', '18']:
172  pt.drawText(0.2, ymax-0.46, date_tag)
173  else:
174  pt.drawText(0.2, ymax-0.46, date_tag)
175  pt.drawText(0.2, ymax-0.52, channel_string + " counting")
176 
177  leg.SetBorderSize(0)
178  leg.SetTextSize(0.07)
179  leg.AddEntry(comb_graph, "#varepsilon_{event}^{single-"+lep+"}", "ep")
180 
181  leg.Draw()
182 
183  if channel == "Zee":
184  new_trig_line = R.TLine(1683743066.0, ymin, 1683743066.0, ymax)
185  new_trig_line.SetLineColor(R.kBlue)
186  new_trig_line.SetLineWidth(3)
187  new_trig_line.SetLineStyle(2)
188  new_trig_line.Draw("same")
189  R.gPad.Update()
190 
191  comb_graph.GetHistogram().SetXTitle("Date")
192  c1.SaveAs(outdir + "event_eff_v_time_"+channel+"_data"+out_tag+"_"+".eps")
193  c1.SaveAs(outdir + "event_eff_v_time_"+channel+"_data"+out_tag+"_"+".pdf")
194 
195 def plot_efficiency(channel, years):
196 
197  print("------------------------------------------")
198  print("Begin Yearwise Efficiency Plots vs Time")
199  print("------------------------------------------")
200 
201  arr_date = []
202  arr_trigeff = []
203  arr_trigerr = []
204  arr_recoeff = []
205  arr_recoerr = []
206  run_num = []
207 
208  for year in years:
209  print("year = ", year)
210  grl = pt.get_grl(year)
211 
212  if year == "23":
213 
214  maindir = args.indir + dir_2023
215  print("2023 grl = ", grl)
216 
217  elif year == "22":
218 
219  maindir = args.indir + dir_2022
220  print("2022 grl = ", grl)
221 
222  for run in grl:
223 
224  print('Begin Run ', run, 'processing')
225 
226  dfz = pd.read_csv(maindir + "run_" + run + ".csv")
227  dfz_small = dfz
228  dfz_small['ZLumi'] = dfz_small[channel + 'Lumi']
229  dfz_small['ZLumiErr'] = dfz_small[channel + 'LumiErr']
230  dfz_small['TrigEff'] = dfz_small[channel + 'EffTrig']
231  dfz_small['TrigErr'] = dfz_small[channel + 'ErrTrig']
232  dfz_small['RecoEff'] = dfz_small[channel + 'EffReco']
233  dfz_small['RecoErr'] = dfz_small[channel + 'ErrReco']
234  dfz_small['LBLive'] = dfz_small['LBLive']
235  dfz_small = dfz_small.drop(dfz_small[dfz_small.ZLumi == 0].index)
236  dfz_small = dfz_small.drop(dfz_small[(dfz_small['LBLive']<10) | (dfz_small['PassGRL']==0)].index)
237 
238  # Cut out all runs shorter than 40 minutes
239  if dfz_small['LBLive'].sum()/60 < 40:
240  print("Skip Run", run, "because of live time", dfz_small['LBLive'].sum()/60, "min")
241  continue
242 
243  # Grab start of the run for plotting later on
244  run_start = dfz_small['LBStart'].iloc[0]
245  timestamp = time.gmtime(run_start)
246  timestamp = R.TDatime(timestamp[0], timestamp[1], timestamp[2], timestamp[3], timestamp[4], timestamp[5])
247  timestamp = timestamp.Convert()
248 
249  # Calculate average trigger efficiency
250  dfz_small['TrigEff'] *= dfz_small['LBLive']
251  total_time = dfz_small['LBLive'].sum()
252  trig_eff_avg = dfz_small['TrigEff'].sum()/total_time
253 
254  # Calculate average reconstruction efficiency
255  dfz_small['RecoEff'] *= dfz_small['LBLive']
256  total_time = dfz_small['LBLive'].sum()
257  reco_eff_avg = dfz_small['RecoEff'].sum()/total_time
258 
259  # Calculate average trigger efficiency error
260  dfz_small['TrigErr'] *= dfz_small['LBLive']
261  dfz_small['TrigErr'] *= dfz_small['TrigErr']
262  total_time = dfz_small['LBLive'].sum()
263  trig_err_avg = math.sqrt(dfz_small['TrigErr'].sum())/total_time
264 
265  # Calculate average reconstruction efficiency error
266  dfz_small['RecoErr'] *= dfz_small['LBLive']
267  dfz_small['RecoErr'] *= dfz_small['RecoErr']
268  total_time = dfz_small['LBLive'].sum()
269  reco_err_avg = math.sqrt(dfz_small['RecoErr'].sum())/total_time
270 
271  arr_date.append(timestamp)
272  arr_trigeff.append(trig_eff_avg)
273  arr_trigerr.append(trig_err_avg)
274  arr_recoeff.append(reco_eff_avg)
275  arr_recoerr.append(reco_err_avg)
276  run_num.append(run)
277 
278  arr_date = array('d', arr_date)
279 
280  arr_trigeff = np.array(arr_trigeff)
281  arr_trigerr = np.array(arr_trigerr)
282  arr_recoeff = np.array(arr_recoeff)
283  arr_recoerr = np.array(arr_recoerr)
284 
285  if channel == "Zee":
286  lep = "e"
287  channel_string = "Z #rightarrow ee"
288  ymin, ymax = 0.64, 0.96
289  elif channel == "Zmumu":
290  lep = "#mu"
291  channel_string = "Z #rightarrow #mu#mu"
292  ymin, ymax = 0.64, 0.96
293 
294  trig_graph = R.TGraphErrors(len(arr_date), arr_date, arr_trigeff, R.nullptr,arr_trigerr)
295  trig_graph.GetHistogram().SetYTitle("Efficiency")
296  trig_graph.GetHistogram().GetYaxis().SetRangeUser(ymin, ymax)
297  trig_graph.GetXaxis().SetTimeDisplay(2)
298  trig_graph.GetXaxis().SetNdivisions(9,R.kFALSE)
299  trig_graph.GetXaxis().SetTimeFormat(time_format)
300  trig_graph.GetXaxis().SetTimeOffset(0,"gmt")
301  trig_graph.SetMarkerSize(1)
302 
303  reco_graph = R.TGraphErrors(len(arr_date), arr_date, arr_recoeff, R.nullptr,arr_recoerr)
304  reco_graph.GetHistogram().GetYaxis().SetRangeUser(ymin, ymax)
305  reco_graph.GetXaxis().SetTimeDisplay(2)
306  reco_graph.GetXaxis().SetNdivisions(9,R.kFALSE)
307  reco_graph.GetXaxis().SetTimeFormat(time_format)
308  reco_graph.GetXaxis().SetTimeOffset(0,"gmt")
309  reco_graph.SetMarkerSize(1)
310  reco_graph.SetMarkerStyle(21)
311  reco_graph.SetMarkerColor(R.kRed)
312  reco_graph.SetLineColor(R.kRed)
313 
314  c1 = R.TCanvas()
315 
316  trig_graph.Draw("ap")
317  reco_graph.Draw("p")
318 
319  if channel == "Zee":
320 
321  leg = R.TLegend(0.645, 0.2, 0.805, 0.4)
322  pt.drawAtlasLabel(0.2, ymax-0.64, "Internal")
323  if year in ['15', '16', '17', '18']:
324  pt.drawText(0.2, ymax-0.70, date_tag)
325  else:
326  pt.drawText(0.2, ymax-0.70, date_tag)
327  pt.drawText(0.2, ymax-0.76, channel_string + " counting")
328 
329  elif channel == "Zmumu":
330 
331  leg = R.TLegend(0.645, 0.45, 0.805, 0.65)
332  pt.drawAtlasLabel(0.2, ymax-0.36, "Internal")
333  if year in ['15', '16', '17', '18']:
334  pt.drawText(0.2, ymax-0.42, date_tag)
335  else:
336  pt.drawText(0.2, ymax-0.42, date_tag)
337  pt.drawText(0.2, ymax-0.48, channel_string + " counting")
338 
339  leg.SetBorderSize(0)
340  leg.SetTextSize(0.07)
341  leg.AddEntry(reco_graph, "#varepsilon_{reco}^{single-"+lep+"}", "ep")
342  leg.AddEntry(trig_graph, "#varepsilon_{trig}^{single-"+lep+"}", "ep")
343 
344  leg.Draw()
345 
346  if channel == "Zee":
347  new_trig_line = R.TLine(1683743066.0, ymin, 1683743066.0, ymax)
348  new_trig_line.SetLineColor(R.kBlue)
349  new_trig_line.SetLineWidth(3)
350  new_trig_line.SetLineStyle(2)
351  new_trig_line.Draw("same")
352  R.gPad.Update()
353 
354  trig_graph.GetHistogram().SetXTitle("Date")
355  c1.SaveAs(outdir + "eff_v_time_"+channel+"_data"+out_tag+"_"+".eps")
356  c1.SaveAs(outdir + "eff_v_time_"+channel+"_data"+out_tag+"_"+".pdf")
357 
358 if __name__ == "__main__":
359  pt.setAtlasStyle()
360  R.gROOT.SetBatch(R.kTRUE)
361  main()
CaloCellPos2Ntuple.int
int
Definition: CaloCellPos2Ntuple.py:24
plotting.yearwise_efficiency.plot_efficiency_comb
def plot_efficiency_comb(channel, years)
Definition: yearwise_efficiency.py:64
convertTimingResiduals.sum
sum
Definition: convertTimingResiduals.py:55
plotting.yearwise_efficiency.plot_efficiency
def plot_efficiency(channel, years)
Definition: yearwise_efficiency.py:195
plotting.yearwise_efficiency.main
def main()
Definition: yearwise_efficiency.py:60
array
Muon::print
std::string print(const MuPatSegment &)
Definition: MuonTrackSteering.cxx:28