32def plot_channel(channel):
33 dfz = pd.read_csv(infilename, delimiter=',')
34 if dfz.empty:
35 print(
'No data in file', infilename,
', exiting.')
36 return
37
38 run_number = dfz.RunNum[0]
39 lhc_fill = dfz.FillNum[0]
40
41 if args.t0:
42 rfo = R.TFile.Open(os.path.join(outdir, 'zlumi.root'), 'RECREATE')
43
44 Lumi = channel+'Lumi'
45 LumiErr = channel+'LumiErr'
46 LumiString = pt.plotlabel[channel]
47
48
49 dfz0 = dfz.copy()
50 if args.dropzero:
51 dfz = dfz.drop(dfz[(dfz[Lumi] == 0)].index)
52 dfz = dfz.drop(dfz[(dfz[
'LBLive']<pt.lblivetimecut) | (dfz[
'PassGRL']==0)].index)
53 if len(dfz) == 0:
54 print(
"No valid LBs found. Exiting")
55 return
56
57 dfz['OffDelLumi'] = dfz['OffLumi']*dfz['LBFull']
58
59
60
61 for entry in [Lumi, LumiErr, 'OffLumi']:
62 dfz[entry] *= dfz['LBLive']
63 normalisation = dfz[Lumi].sum() / dfz['OffLumi'].sum()
64 dfz[LumiErr] *= dfz[LumiErr]
65 if args.usemu:
66 dfz['OffMu'] = dfz['OffMu'].astype(int)
67 dfz = dfz.groupby(['OffMu']).sum()
68 else:
69 dfz['LBNum'] = (dfz['LBNum']//20)*20
70 dfz = dfz.groupby(['LBNum']).sum()
71 dfz[LumiErr] = np.sqrt(dfz[LumiErr])
72 for entry in [Lumi, LumiErr, 'OffLumi']:
73 dfz[entry] /= dfz['LBLive']
74
75 print(f
"{Lumi}/OffLumi mean ratio = {normalisation:.3f}")
76
77 if args.t0:
78 translist = []
79 timeformat = '%y/%m/%d %H:%M:%S'
80 import time
81 for k in range(len(dfz)):
82 thisrow = dfz.iloc[k]
83 lb = thisrow.name
84 lbselector = (lb <= dfz0['LBNum']) & (dfz0['LBNum'] < lb+20)
85 translist.append({'FillNum' : dfz0[lbselector]['FillNum'].iloc[0],
86 'beginTime' : time.strftime(timeformat,
87 time.gmtime(dfz0[lbselector]['LBStart'].iloc[0])),
88 'endTime': time.strftime(timeformat,
89 time.gmtime(dfz0[(dfz0['LBNum']//20)*20 == lb]['LBEnd'].iloc[-1])),
90 'ZmumuRate': thisrow['ZmumuRate']/thisrow['LBLive'],
91 'instDelLumi': thisrow['OffDelLumi']/thisrow['LBFull']/1e3,
92 'delLumi': thisrow['OffLumi']/1e3,
93 'ZDel': thisrow['ZmumuRate']
94 })
95 pdn = pd.DataFrame.from_records(translist)
96 pdn.to_csv(os.path.join(args.outdir, 'zrate.csv'), header=False, index=False)
97
98 if args.usemu:
99 leg = R.TLegend(0.6, 0.33, 0.75, 0.5)
100 xtitle = "<#mu>"
101 outfile = channel+"ATLAS_vs_mu"
102 else:
103 leg = R.TLegend(0.62, 0.75, 0.75, 0.9)
104 xtitle = "Luminosity Block Number"
105 outfile = channel+"ATLAS_vs_lb"
106
107 if args.absolute:
108 outfile += "_abs"
109 normalisation = 1.
110
111 arr_bins =
array(
'd', dfz.index)
112 arr_zlumi =
array(
'd', dfz[Lumi] / normalisation)
113 arr_zlumi_err =
array(
'd', dfz[LumiErr] / normalisation)
114 arr_rat =
array(
'd', dfz[Lumi] / dfz[
'OffLumi'] / normalisation)
115 arr_rat_err =
array(
'd', dfz[LumiErr] / dfz[
'OffLumi'] / normalisation)
116 arr_olumi =
array(
'd', dfz[
'OffLumi'])
117
118 hz = R.TGraphErrors(len(arr_bins), arr_bins, arr_zlumi, R.nullptr, arr_zlumi_err)
119 ho = R.TGraphErrors(len(arr_bins), arr_bins, arr_olumi, R.nullptr, R.nullptr)
120 hr = R.TGraphErrors(len(arr_bins), arr_bins, arr_rat, R.nullptr, arr_rat_err)
121
122 ho.SetLineColor(R.kAzure)
123 ho.SetLineWidth(3)
124
125 c1 = R.TCanvas()
126
127 pad1 = R.TPad("pad1", "pad1", 0, 0, 1, 1)
128 pad1.SetBottomMargin(0.3)
129 pad1.SetFillStyle(4000)
130 pad1.Draw()
131 pad1.cd()
132 pad1.RedrawAxis()
133 hz.GetXaxis().SetLabelSize(0)
134
135 xmin = hz.GetXaxis().GetXmin()
136 xmax = hz.GetXaxis().GetXmax()
137
138 hz.SetMarkerStyle(4)
139 hz.Draw('ap')
140 ho.Draw("same L")
141 hz.Draw('same p')
142 hz.GetYaxis().SetTitle("Luminosity [10^{33} cm^{-2}s^{-1}]")
143 hz.GetXaxis().SetTitle(xtitle)
144 hz.GetXaxis().SetTitleOffset(0.8)
145 ymax = hz.GetHistogram().GetMaximum()
146 ymin = hz.GetHistogram().GetMinimum()
147 if not args.usemu:
148 hz.GetYaxis().SetRangeUser(ymin*0.8, ymax*1.4)
149
150 leg.SetBorderSize(0)
151 leg.SetTextSize(0.055)
152 if args.absolute:
153 leg.AddEntry(hz, "L_{"+LumiString+"}", "ep")
154 else:
155 leg.AddEntry(hz, "L_{"+LumiString+"}^{normalised to L_{ATLAS}^{fill}}", "ep")
156 leg.AddEntry(ho, "L_{ATLAS}", "l")
157 leg.Draw()
158
159 pad2 = R.TPad("pad2", "pad2", 0, 0, 1, 1)
160 pad2.SetTopMargin(0.72)
161 pad2.SetBottomMargin(0.1)
162 pad2.SetFillStyle(4000)
163 pad2.Draw()
164 pad2.cd()
165
166 hr.Draw("ap0")
167 hr.GetYaxis().SetTitle("Ratio")
168 hr.GetXaxis().SetTitle(xtitle)
169 hr.GetXaxis().SetTitleOffset(0.865)
170
171 ymin, ymax = 0.95, 1.05
172 median = np.median(arr_rat)
173 spread = np.percentile(abs(arr_rat - median), 68)
174 if spread > 0.03: ymin, ymax = 0.9, 1.1
175 if args.absolute and args.t0: ymin, ymax = 0.9, 1.1
176 hr.GetYaxis().SetRangeUser(ymin, ymax)
177
178 line0 = R.TLine(xmin, 1.0, xmax, 1.0)
179 line0.SetLineStyle(2)
180 line0.Draw()
181
182 line1 = R.TLine(xmin, (ymin+1.)/2., xmax, (ymin+1.)/2.)
183 line1.SetLineStyle(2)
184 line1.Draw()
185
186 line2 = R.TLine(xmin, (ymax+1.)/2., xmax, (ymax+1.)/2.)
187 line2.SetLineStyle(2)
188 line2.Draw()
189
190 hr.GetYaxis().SetNdivisions(3)
191
192 pt.drawAtlasLabel(0.2, 0.88, "Internal")
193 if not args.t0:
194 pt.drawText(0.2, 0.82, pt.get_yearsqrtstxt(run_number), size=22)
195 pt.drawText(0.2, 0.77, "LHC Fill " + str(lhc_fill), size=22)
196 pt.drawText(0.2, 0.72, LumiString + " counting", size=22)
197
198 line4 = pt.make_bands(arr_bins, spread, median)
199 line4.Draw("same 3")
200 hr.Draw("same ep0")
201
202 if not args.t0:
203 c1.SaveAs(os.path.join(outdir, outfile + ".pdf"))
204 else:
205 rfo.WriteTObject(c1, outfile)
206 if channel == 'Zmumu':
207 rfo.WriteTObject(hz, 'z_lumi')
208 rfo.WriteTObject(hr, 'z_lumi_ratio')
209 c1.Clear()
210
211
212
213 c2 = R.TCanvas()
214 hr.GetXaxis().SetTitleOffset(1)
215 hr.Draw('ap0')
216
217 hr.Fit('pol0', 'q0')
218 hr.GetFunction('pol0').SetLineColor(R.kRed)
219
220 hr.GetYaxis().SetTitle("Ratio")
221 hr.GetXaxis().SetTitle(xtitle)
222 hr.GetYaxis().SetRangeUser(0.9, 1.1)
223 hr.GetYaxis().SetNdivisions()
224
225 mean = hr.GetFunction("pol0").GetParameter(0)
226 line0 = pt.make_bands(arr_bins, spread, median)
227 line0.Draw("same 3")
228 hr.Draw("same ep0")
229 hr.GetFunction('pol0').Draw('same l')
230
231 pt.drawAtlasLabel(0.2, 0.88, "Internal")
232 if not args.t0:
233 pt.drawText(0.2, 0.82, pt.get_yearsqrtstxt(run_number), size=22)
234 pt.drawText(0.2, 0.77, "LHC Fill " + str(lhc_fill), size=22)
235 pt.drawText(0.2, 0.72, LumiString + " counting", size=22)
236
237 leg = R.TLegend(0.17, 0.2, 0.90, 0.3)
238 leg.SetBorderSize(0)
239 leg.SetTextSize(0.05)
240 leg.SetNColumns(3)
241 leg.AddEntry(hr, "L_{"+LumiString+"}/L_{ATLAS}", "ep")
242 leg.AddEntry(hr.GetFunction("pol0"), f"Mean = {mean:.3f}", "l")
243 leg.AddEntry(line0, "68% band", "f")
244 leg.Draw()
245
246 if not args.t0:
247 c2.SaveAs(os.path.join(outdir, outfile+"_ratio.pdf"))
248 else:
249 rfo.WriteTObject(c2, f'{outfile}_ratio')
250
251 return
252
void print(char *figname, TCanvas *c1)
1.14.0