Belle II Software development
top.py
1#!/usr/bin/env python3
2
3
10
11'''
12Validation plots for TOP calibration.
13'''
14
15from prompt import ValidationSettings
16from ROOT import TFile, TH1F, TGraph, TCanvas, TLegend, gROOT, gStyle, PyConfig
17from basf2 import B2ERROR
18import sys
19
20
21settings = ValidationSettings(name='TOP post-tracking calibration',
22 description=__doc__,
23 download_files=[],
24 expert_config=None)
25
26
27def run_validation(job_path, input_data_path, requested_iov, expert_config):
28 '''
29 Makes validation plots
30 :job_path: path to top posttracking calibration output
31 :input_data_path: required argument but not used
32 :requested_iov: required argument but not used
33 :expert_config: required argument but not used
34 '''
35
36 # input/output file names
37
38 inputFileName = job_path + '/TOP_validation/0/algorithm_output/TOPCalValidation.root'
39 outputFileName = 'TOPValidation.pdf'
40
41 # configuration
42
43 PyConfig.IgnoreCommandLineOptions = True
44 gROOT.SetBatch(True)
45 gROOT.ForceStyle()
46 gStyle.SetOptStat(0)
47 gStyle.SetOptFit(1)
48 gStyle.SetHistMinimumZero(True)
49 gStyle.SetLabelSize(0.05, "X")
50 gStyle.SetLabelSize(0.05, "Y")
51 gStyle.SetTitleSize(0.06, "X")
52 gStyle.SetTitleSize(0.06, "Y")
53 gStyle.SetTitleOffset(0.7, "X")
54 gStyle.SetTitleOffset(0.8, "Y")
55 gStyle.SetGridColor(11)
56
57 # open input file and create the canvas
58
59 file_in = TFile.Open(inputFileName)
60 if (not file_in):
61 B2ERROR(inputFileName + ": file not found")
62 return
63 canvas = TCanvas("c1", "c1", 1000, 750)
64
65 # prepare and fit histograms for channel T0
66
67 channelT0s = []
68 residuals = []
69 residuals_all = TH1F("residuals_all", "ChannelT0 residuals; residuals [ns]", 200, -1.0, 1.0)
70 h_rms = TH1F("h_rms", "Width of residuals (rms); slot; r.m.s [ps]", 16, 0.5, 16.5)
71 h_core = TH1F("h_core", "Core width of residuals; slot; sigma [ps]", 16, 0.5, 16.5)
72 h_out = TH1F("h_out", "Outlayers (more than 3 sigma); slot; fraction [%]", 16, 0.5, 16.5)
73 for slot in range(1, 17):
74 name = "slot" + f'{slot:02d}'
75 h = file_in.Get("channelT0_" + name)
76 channelT0s.append(h)
77 r = TH1F("residuals_" + name, "ChannelT0, " + name + "; residuals [ns]", 200, -1.0, 1.0)
78 residuals.append(r)
79 if not h:
80 continue
81 for i in range(1, h.GetNbinsX() + 1):
82 if h.GetBinError(i) > 0:
83 r.Fill(h.GetBinContent(i))
84 residuals_all.Fill(h.GetBinContent(i))
85 r.Fit("gaus", "Q")
86 residuals_all.Fit("gaus", "Q")
87 h_rms.SetBinContent(slot, r.GetRMS() * 1000)
88 fun = r.GetFunction("gaus")
89 if not fun:
90 continue
91 fun.SetLineWidth(1)
92 h_core.SetBinContent(slot, fun.GetParameter(2) * 1000)
93 mean = fun.GetParameter(1)
94 wid = fun.GetParameter(2)
95 n_all = 0
96 n_in = 0
97 for i in range(1, h.GetNbinsX() + 1):
98 if h.GetBinError(i) == 0:
99 continue
100 n_all += 1
101 if abs(h.GetBinContent(i) - mean) < 3 * wid:
102 n_in += 1
103 if n_all > 0:
104 h_out.SetBinContent(slot, (1 - n_in / n_all) * 100)
105 residuals_all.Fit("gaus", "Q")
106
107 # open pdf file
108
109 canvas.Print(outputFileName + "[")
110
111 # make plots of channelT0 residuals vs. channel
112
113 canvas.Clear()
114 canvas.Divide(4, 4)
115 graphs = []
116 for i, h in enumerate(channelT0s):
117 canvas.cd(i + 1)
118 if h:
119 h.SetMinimum(-1.0)
120 h.SetMaximum(1.0)
121 h.Draw()
122 g_over = TGraph()
123 g_under = TGraph()
124 for ibin in range(1, h.GetNbinsX() + 1):
125 if h.GetBinContent(ibin) > 1.0:
126 g_over.SetPoint(g_over.GetN(), h.GetBinCenter(ibin), 0.95)
127 elif h.GetBinContent(ibin) < -1.0:
128 g_under.SetPoint(g_under.GetN(), h.GetBinCenter(ibin), -0.95)
129 if g_over.GetN() > 0:
130 g_over.SetMarkerStyle(22)
131 g_over.SetMarkerColor(2)
132 g_over.Draw("P")
133 graphs.append(g_over)
134 if g_under.GetN() > 0:
135 g_under.SetMarkerStyle(23)
136 g_under.SetMarkerColor(2)
137 g_under.Draw("P")
138 graphs.append(g_under)
139 canvas.Print(outputFileName)
140
141 # make plots of channelT0 residual distributions
142
143 canvas.Clear()
144 canvas.Divide(4, 4)
145 gStyle.SetOptFit(0)
146 gStyle.SetOptStat("oue")
147 W = gStyle.GetStatW()
148 H = gStyle.GetStatH()
149 gStyle.SetStatW(0.4)
150 gStyle.SetStatH(0.3)
151 for i, h in enumerate(residuals):
152 canvas.cd(i + 1)
153 h.Draw()
154 canvas.Print(outputFileName)
155 gStyle.SetStatW(W)
156 gStyle.SetStatH(H)
157
158 # make summary plots for channelT0
159
160 canvas.Clear()
161 canvas.Divide(2, 2)
162 gStyle.SetOptFit(111)
163 gStyle.SetOptStat(0)
164 histos = [residuals_all, h_rms, h_core, h_out]
165 for i, h in enumerate(histos):
166 canvas.cd(i + 1)
167 h.Draw()
168 canvas.Print(outputFileName)
169
170 # prepare commonT0 histogram
171
172 commonT0 = file_in.Get("commonT0")
173 if commonT0:
174 commonT0.Scale(1000)
175 commonT0.SetYTitle(commonT0.GetYaxis().GetTitle().replace('[ns]', '[ps]'))
176 commonT0.SetMinimum(-100)
177 commonT0.SetMaximum(100)
178 commonT0.Fit("pol0", "Q")
179
180 # make plots for commonT0
181
182 canvas.Clear()
183 canvas.Divide(1, 2)
184 canvas.GetPad(2).SetGrid()
185 gStyle.SetHistMinimumZero(False)
186 histos = [commonT0, file_in.Get("runIndex")]
187 for i, h in enumerate(histos):
188 canvas.cd(i + 1)
189 if h:
190 h.SetMarkerStyle(20)
191 h.SetMarkerSize(0.7)
192 h.Draw()
193 canvas.Print(outputFileName)
194 gStyle.SetHistMinimumZero(True)
195
196 # make plots for moduleT0
197
198 canvas.Clear()
199 canvas.Divide(4, 4)
200 histos = []
201 for slot in range(1, 17):
202 canvas.cd(slot)
203 h = file_in.Get("moduleT0_slot" + f'{slot:02d}')
204 if h:
205 h.Scale(1000)
206 h.SetYTitle(h.GetYaxis().GetTitle().replace('[ns]', '[ps]'))
207 h.SetMinimum(-100)
208 h.SetMaximum(100)
209 h.Fit("pol0", "Q")
210 histos.append(h)
211 canvas.Print(outputFileName)
212
213 # prepare moduleT0 summary histograms
214
215 h_offset = TH1F("h_offset", "ModuleT0: fitted offset; slot; p0 [ps]", 16, 0.5, 16.5)
216 h_chi = TH1F("h_chi", "ModuleT0: chi2/ndf; slot; #chi^{2}/ndf", 16, 0.5, 16.5)
217 h_res = TH1F("h_res", "ModuleT0: residuals; residual - p0 [ps]", 200, -100, 100)
218 h_pul = TH1F("h_pul", "ModuleT0: pulls; (residual - p0) / error", 200, -10, 10)
219 for i, h in enumerate(histos):
220 if not h:
221 continue
222 fun = h.GetFunction("pol0")
223 if not fun:
224 continue
225 p0 = fun.GetParameter(0)
226 err = fun.GetParError(0)
227 h_offset.SetBinContent(i + 1, p0)
228 h_offset.SetBinError(i + 1, err)
229 h_chi.SetBinContent(i + 1, fun.GetChisquare() / max(fun.GetNDF(), 1))
230 for ibin in range(1, h.GetNbinsX() + 1):
231 err = h.GetBinError(ibin)
232 if err > 0:
233 h_res.Fill(h.GetBinContent(ibin) - p0)
234 h_pul.Fill((h.GetBinContent(ibin) - p0) / err)
235
236 # plot moduleT0 summary histograms
237
238 canvas.Clear()
239 canvas.Divide(2, 2)
240 histos = [h_offset, h_chi, h_res, h_pul]
241 h_offset.SetStats(False)
242 h_chi.SetStats(False)
243 gStyle.SetOptStat("rme")
244 for i, h in enumerate(histos):
245 canvas.cd(i + 1)
246 h.Draw()
247 canvas.Print(outputFileName)
248 gStyle.SetOptStat(0)
249
250 # plot fractions of active and active-and-calibrated channels
251
252 canvas.Clear()
253 canvas.Divide(4, 4)
254 legends = []
255 for slot in range(1, 17):
256 canvas.cd(slot)
257 name = "slot" + f'{slot:02d}'
258 histos = [file_in.Get("numActive_" + name), file_in.Get("numActiveCalibrated_" + name)]
259 legend = TLegend(0.1, 0.1, 0.7, 0.3)
260 legends.append(legend)
261 opt = ''
262 for i, h in enumerate(histos):
263 if h:
264 h.SetMinimum(0.0)
265 h.SetMaximum(1.05)
266 text = h.GetTitle().split(',')[0]
267 h.SetTitle(name)
268 h.SetLineColor(i+1)
269 h.Draw(opt)
270 opt = 'same'
271 legend.AddEntry(h, text)
272 if legend.GetNRows() > 0:
273 legend.Draw("same")
274 canvas.Print(outputFileName)
275
276 # plot fractions of time base and channelT0 calibrated channels
277
278 canvas.Clear()
279 canvas.Divide(4, 4)
280 legends = []
281 for slot in range(1, 17):
282 canvas.cd(slot)
283 name = "slot" + f'{slot:02d}'
284 histos = [file_in.Get("numTBCalibrated_" + name), file_in.Get("numT0Calibrated_" + name)]
285 legend = TLegend(0.1, 0.1, 0.7, 0.3)
286 legends.append(legend)
287 opt = ''
288 for i, h in enumerate(histos):
289 if h:
290 h.SetMinimum(0.0)
291 h.SetMaximum(1.05)
292 text = h.GetTitle().split(',')[0]
293 h.SetTitle(name)
294 h.SetLineColor(i+1)
295 h.Draw(opt)
296 opt = 'same'
297 legend.AddEntry(h, text)
298 if legend.GetNRows() > 0:
299 legend.Draw("same")
300 canvas.Print(outputFileName)
301
302 # plot threshold efficiencies
303
304 canvas.Clear()
305 canvas.Divide(4, 4)
306 legends = []
307 for slot in range(1, 17):
308 canvas.cd(slot)
309 name = "slot" + f'{slot:02d}'
310 h = file_in.Get("thrEffi_slot" + f'{slot:02d}')
311 legend = TLegend(0.1, 0.1, 0.7, 0.3)
312 legends.append(legend)
313 if h:
314 h.SetMinimum(0.0)
315 h.SetMaximum(1.0)
316 text = h.GetTitle().split(',')[0]
317 h.SetTitle(name)
318 h.Draw()
319 legend.AddEntry(h, text)
320 legend.Draw("same")
321 canvas.Print(outputFileName)
322
323 # plot BS13d asic shifts
324
325 canvas.Clear()
326 canvas.Divide(2, 2)
327 for cb in range(4):
328 canvas.cd(cb + 1)
329 canvas.GetPad(cb + 1).SetGrid()
330 h = file_in.Get("asicShifts_" + str(cb))
331 if h:
332 h.SetMinimum(-50.0)
333 h.SetMaximum(50.0)
334 h.Draw()
335 canvas.Print(outputFileName)
336
337 # plot offsets
338
339 canvas.Clear()
340 canvas.Divide(2, 2)
341 offsets = [file_in.Get("svdOffset"), file_in.Get("cdcOffset")]
342 sigmas = [file_in.Get("svdSigma"), file_in.Get("cdcSigma")]
343 components = ['SVD', 'CDC']
344 colors = [4, 2]
345
346 canvas.cd(1)
347 hmin = 0
348 hmax = 0
349 for h in offsets:
350 if h:
351 hmin = min(hmin, h.GetMinimum())
352 hmax = max(hmax, h.GetMaximum())
353 hmin = round(hmin, 0) - 1
354 hmax = round(hmax, 0) + 1
355 for i, h in enumerate(offsets):
356 if h:
357 h.SetMinimum(hmin)
358 h.SetMaximum(hmax)
359 h.SetMarkerStyle(24)
360 h.SetMarkerColor(colors[i])
361 h.SetLineColor(colors[i])
362 first = True
363 legend1 = TLegend()
364 for i, h in enumerate(offsets):
365 if not h or not sigmas[i] or sigmas[i].Integral() == 0:
366 continue
367 if first:
368 h.Draw()
369 first = False
370 else:
371 h.Draw("same")
372 legend1.AddEntry(h, components[i])
373 legend1.Draw("same")
374
375 canvas.cd(2)
376 hmax = 0
377 for h in sigmas:
378 if h:
379 hmax = max(hmax, h.GetMaximum())
380 for i, h in enumerate(sigmas):
381 if h:
382 h.SetMinimum(0)
383 h.SetMaximum(hmax * 1.1)
384 h.SetMarkerStyle(24)
385 h.SetMarkerColor(colors[i])
386 h.SetLineColor(colors[i])
387 first = True
388 legend2 = TLegend()
389 for i, h in enumerate(sigmas):
390 if not h or h.Integral() == 0:
391 continue
392 if first:
393 h.Draw()
394 first = False
395 else:
396 h.Draw("same")
397 legend2.AddEntry(h, components[i])
398 legend2.Draw("same")
399
400 fillPatt = [file_in.Get("fillPatternOffset"), file_in.Get("fillPatternFract")]
401 for i, h in enumerate(fillPatt):
402 canvas.cd(3 + i)
403 if h:
404 if i == 1:
405 for k in range(h.GetNbinsX()):
406 h.SetBinError(k + 1, h.GetBinContent(k + 1) / 1000)
407 h.Draw()
408
409 canvas.Print(outputFileName)
410
411 # close pdf file
412
413 canvas.Print(outputFileName + "]")
414
415 # plot photon hits, use png format to shorten file size
416
417 canvas.Clear()
418 canvas.Divide(4, 4)
419 for slot in range(1, 17):
420 canvas.cd(slot)
421 h = file_in.Get("hits_slot" + f'{slot:02d}')
422 if h:
423 h.Draw("colz")
424 canvas.SaveAs(outputFileName.replace('.pdf', 'Hits.png'))
425
426
427if __name__ == "__main__":
428 run_validation(*sys.argv[1:])
run_validation
Definition: run_validation.py:1