Belle II Software  release-05-01-25
cdc_tracking_validation.py
1 # -*- coding: utf-8 -*-
2 # """
3 # CDC tracking validation
4 # """
5 
6 import basf2
7 from prompt import ValidationSettings
8 import sys
9 import os
10 
11 
12 settings = ValidationSettings(name='CDC Tracking',
13  description=__doc__,
14  download_files=[],
15  expert_config=None)
16 
17 
18 def findLastIteration(job_path, algorithm):
19  root = f"{job_path}/{algorithm}/"
20  list_dirs = [item for item in os.listdir(root) if item != 'outputdb']
21  list_dirs.sort(reverse=True)
22  print(f"The following iterations are available for {algorithm} algorithm: {list_dirs}")
23  return list_dirs[0]
24 
25 
26 def run_validation(job_path, input_data_path, requested_iov, expert_config, **kwargs):
27  # job_path will be replaced with path/to/calibration_results
28  # input_data_path will be replaced with path/to/data_path used for calibration, e.g. /group/belle2/dataprod/Data/PromptSkim/
29 
30  # will be included in the expert_config in the future
31  file_extension = 'png'
32  algo_tz = 'tz2'
33  algo_tw = 'tw0'
34  algo_sr = 'sr0'
35  algo_xt = 'xt0'
36 
37  import os
38  import ROOT
39 
40  ROOT.gROOT.SetBatch(True)
41  ROOT.gStyle.SetOptStat(0)
42 
43  plot_directory = "plots"
44  if not os.path.exists(plot_directory):
45  os.makedirs(plot_directory)
46 
47 
48  print("**** T0 validation plots ****")
49 
50  lastIt = findLastIteration(job_path, algo_tz)
51  histT0_tz = f'{job_path}/{algo_tz}/{lastIt}/algorithm_output/histT0_{algo_tz}.root'
52  f_histT0_tz = ROOT.TFile(histT0_tz)
53  print(f"Plots from {algo_tz}/{lastIt}/algorithm_output/histT0_{algo_tz}.root")
54 
55  # plot total histograms of the fitted means and sigmas for each wire
56  hm_All = f_histT0_tz.Get("hm_All")
57  can1 = ROOT.TCanvas()
58  hm_All.GetXaxis().SetRangeUser(-1, 0)
59  hm_All.GetXaxis().SetTitle(r"<\Delta t> [ns]")
60  hm_All.Draw("pe")
61  can1.Draw()
62  can1.SaveAs(f"{plot_directory}/{algo_tz}_hm_All.{file_extension}")
63 
64  hs_All = f_histT0_tz.Get("hs_All")
65  can2 = ROOT.TCanvas()
66  hs_All.GetXaxis().SetRangeUser(2, 8)
67  hs_All.GetXaxis().SetTitle(r"\sigma(\Delta t) [ns]")
68  hs_All.Draw("pe")
69  can2.Draw()
70  can2.SaveAs(f"{plot_directory}/{algo_tz}_hs_All.{file_extension}")
71 
72  # plot total histogram of the fitted means for each wire + fit
73  hTotal = f_histT0_tz.Get("hTotal")
74  can3 = ROOT.TCanvas()
75  hTotal.GetXaxis().SetTitle(r"<\Delta t> [ns]")
76  hTotal.Draw("pe")
77  can3.Draw()
78  can3.SaveAs(f"{plot_directory}/{algo_tz}_hTotal.{file_extension}")
79 
80  # plot histograms of the fitted means per layer
81  gr1 = [f_histT0_tz.Get(f'DeltaT0/lay{i}') for i in range(56)]
82  cs = [ROOT.TCanvas(f'cs{i}', '', 1200, 750) for i in range(7)]
83  for i in range(7):
84  cs[i].Divide(4, 2)
85  for j in range(7):
86  pad = [cs[j].GetPrimitive(f'cs{j}_{i + 1}') for i in range(8)]
87  cs[0].cd(j + 1)
88  for i in range(8):
89  pad[i].cd()
90  pad[i].SetGrid()
91  gid = 8 * j + i
92  if gr1[gid]:
93  gr1[gid].Draw("AP")
94  for i in range(7):
95  cs[i].Draw()
96  cs[i].SaveAs(f"{plot_directory}/{algo_tz}_layers-{i}.{file_extension}")
97 
98 
99  print("**** TW validation plots ****")
100 
101  lastIt = findLastIteration(job_path, algo_tw)
102  histTW_tw = f'{job_path}/{algo_tw}/{lastIt}/algorithm_output/histTW_{algo_tw}.root'
103  f_histTW_tw = ROOT.TFile(histTW_tw)
104  print(f"Plots from {algo_tw}/{lastIt}/algorithm_output/histTW_{algo_tw}.root")
105 
106  # plot histograms for each board
107  rangeBorad = range(1, 301)
108  board_1D = [f_histTW_tw.Get(f'h1D/board_{boardID}_1') for boardID in rangeBorad]
109  can = [ROOT.TCanvas(f'c{c}', f'c{c}', 2000, 1500) for c in range(12)]
110  for c in range(12):
111  ni, nj = 5, 5
112  can[c].Divide(ni, nj)
113  for j in range(ni * nj):
114  can[c].cd(j + 1)
115  boardNumber = j + c * (ni * nj)
116  if board_1D[boardNumber]:
117  board_1D[boardNumber].SetMarkerStyle(2)
118  board_1D[boardNumber].SetTitle(f'Board number {boardNumber}')
119  board_1D[boardNumber].GetXaxis().SetTitle("ADC count")
120  board_1D[boardNumber].GetYaxis().SetTitle(r"$\Delta t$ [ns]")
121  board_1D[boardNumber].Draw()
122  can[c].Draw()
123  can[c].SaveAs(f"{plot_directory}/{algo_tw}_boards-{c}.{file_extension}")
124 
125 
126  print("**** sigma res validation plots ****")
127 
128  lastIt = findLastIteration(job_path, algo_sr)
129  histsr_sr = f'{job_path}/{algo_sr}/{lastIt}/algorithm_output/histSigma_{algo_sr}.root'
130  f_histsr_sr = ROOT.TFile(histsr_sr)
131  print(f"Plots from {algo_sr}/{lastIt}/algorithm_output/histSigma_{algo_sr}.root")
132 
133  alpha = 4
134  theta = 3
135  for LR in [0, 1]:
136  histograms = [f_histsr_sr.Get(f'lay_{ilay}/sigma2_lay{ilay}_lr{LR}_al{alpha}_th{theta}') for ilay in range(56)]
137  ncols = 5
138  count_h = 0
139  for h in histograms:
140  if h:
141  count_h += 1
142  div, mod = count_h // ncols, count_h % ncols
143  nrows = div
144  if mod != 0:
145  nrows += 1
146 
147  print(f"number of valid histograms = {count_h} => canvas layout = ({ncols}, {nrows})")
148 
149  c2 = ROOT.TCanvas('c2', '', ncols * 700, nrows * 400) # ncols * 700, nrows * 400)
150  c2.Divide(ncols, nrows)
151  j = 0
152  for h in histograms:
153  if h:
154  j = j + 1
155  c2.cd(j)
156  h.GetXaxis().SetTitle("drift lenght [cm]")
157  h.GetYaxis().SetTitle("#sigma_{r}^{2} = #sigma_{u}.#sigma_{d}")
158  h.Draw("AP")
159  h.Draw()
160  c2.Draw()
161  c2.SaveAs(f"{plot_directory}/{algo_sr}_lr{LR}_al{alpha}_th{theta}.{file_extension}")
162 
163 
164  print("**** XT validation plots ****")
165 
166  lastIt = findLastIteration(job_path, algo_xt)
167  histXT_xt = f'{job_path}/{algo_xt}/{lastIt}/algorithm_output/histXT_{algo_xt}.root'
168  f_histXT_xt = ROOT.TFile(histXT_xt)
169  print(f"Plots from {algo_xt}/{lastIt}/algorithm_output/histXT_{algo_xt}.root")
170 
171  # get histograms for all layers for a specific incident angles
172  thetaID = 1
173  alphaID = 3
174  for LR in [0, 1]:
175  histograms = [f_histXT_xt.Get(f'lay_{layer}/m_hProf{layer}_{LR}_{alphaID}_{thetaID}') for layer in range(56)]
176  # get the number of plots in the canvas
177  ncols = 5
178  count_h = 0
179  for h in histograms:
180  if h:
181  count_h += 1
182  div, mod = count_h // ncols, count_h % ncols
183  nrows = div
184  if mod != 0:
185  nrows += 1
186 
187  c2 = ROOT.TCanvas('c2', '', ncols * 900, nrows * 600)
188  c2.Divide(ncols, nrows)
189  j = 0
190  for histo in histograms:
191  if histo:
192  j = j + 1
193  c2.cd(j)
194  histo.GetYaxis().SetRangeUser(-1.2, 1.2)
195  histo.SetMarkerStyle(2)
196  histo.Draw()
197  c2.Draw()
198  c2.SaveAs(f"{plot_directory}/{algo_xt}_lr{LR}_al{alphaID}_th{thetaID}.{file_extension}")
199 
200 
201 if __name__ == "__main__":
202  run_validation(*sys.argv[1:])