Belle II Software  release-08-00-04
TRGValidation.py
1 #!/usr/bin/env python3
2 
3 
10 
11 # The VariablesToHistogram module saves variables from the VariableManager
12 # to TH1F and TH2F here is an example of how to use it.
13 #
14 # For full documentation please refer to https://software.belle2.org
15 # Anything unclear? Ask questions at https://questions.belle2.org
16 
17 """
18 <header>
19  <input>TRGValidationGen.root</input>
20  <output>TRGValidation.root</output>
21  <contact>yinjh2012@korea.ac.kr</contact>
22  <description>makes validation plots for TRG</description>
23 </header>
24 
25 """
26 
27 import basf2
28 import modularAnalysis as ma # a shorthand for the analysis tools namespace
29 import stdCharged as stdc
30 import variables.utils as vu
31 
32 import ROOT
33 from ROOT import Belle2, TH1F, TFile, TNamed, TEfficiency
34 from math import pi as PI
35 
36 Fac = 180.0 / PI
37 inputBits = ["t3_0", "ty_0", "t2_0", "ts_0", "ta_0", "typ", "ehigh", "elow", "elum", "ecl_3dbha", "cdc_open90",
38  "clst_0", "clst_1", "clst_2", "clst_3", "klm_hit", "klm_0", "klm_1", "klm_2", "eklm_0", "eklm_1", "eklm_2"]
39 outputBits = ["fff", "ffy", "ffz", "fzo", "fso", "fyo", "ffb", "fsb", "ssb", "stt", "hie", "c4", "bha3d",
40  "lml1", "mu_b2b", "mu_eb2b", "eklm2", "beklm", "cdcklm1", "cdcklm2", "seklm1", "seklm2", "ecleklm1"]
41 moniInbits = []
42 moniOutbits = []
43 
44 inputBits_expert = [
45  't3_0',
46  't3_1',
47  't3_2',
48  't3_3',
49  'ty_0',
50  'ty_1',
51  'ty_2',
52  'ty_3',
53  't2_0',
54  't2_1',
55  't2_2',
56  't2_3',
57  'ts_0',
58  'ts_1',
59  'ts_2',
60  'ts_3',
61  'ta_0',
62  'ta_1',
63  'ta_2',
64  'ta_3',
65  'typ',
66  'typ4',
67  'typ5',
68  'cdc_open90',
69  'cdc_active',
70  'cdc_b2b3',
71  'cdc_b2b5',
72  'cdc_b2b7',
73  'cdc_b2b9',
74  'itsfb2b',
75  'ti\t 0',
76  'i2io',
77  'i2fo',
78  'f2f30',
79  's2f30',
80  's2s30',
81  's2s3',
82  's2s5',
83  's2so',
84  's2f3',
85  's2f5',
86  's2fo',
87  'fwd_s',
88  'bwd_s',
89  'track',
90  'trkflt',
91  'ehigh',
92  'elow',
93  'elum',
94  'ecl_bha',
95  'ecl_3dbha',
96  'bha_veto',
97  'bha_type_0',
98  'bha_type_1',
99  'bha_intrk',
100  'bha_theta_0',
101  'bha_theta_1',
102  'clst_0',
103  'clst_1',
104  'clst_2',
105  'clst_3',
106  'ecl_active',
107  'ecl_timing_fwd',
108  'ecl_timing_brl',
109  'ecl_timing_bwd',
110  'ecl_phys',
111  'ecl_oflo',
112  'ecl_lml_0',
113  'ecl_lml_1',
114  'ecl_lml_2',
115  'ecl_lml_3',
116  'ecl_lml_4',
117  'ecl_lml_5',
118  'ecl_lml_6',
119  'ecl_lml_7',
120  'ecl_lml_8',
121  'ecl_lml_9',
122  'ecl_lml_10',
123  'ecl_lml_12',
124  'ecl_lml_13',
125  'ecl_mumu',
126  'ecl_bhapur',
127  'ecl_bst',
128  'top_0',
129  'top_1',
130  'top_2',
131  'top_bb',
132  'top_active',
133  'klm_hit',
134  'klm_0',
135  'klm_1',
136  'klm_2',
137  'klmb2b',
138  'eklm_hit',
139  'eklm_0',
140  'eklm_1',
141  'eklm_2',
142  'eklmb2b',
143  'revo',
144  'her_kick',
145  'ler_kick',
146  'bha_delay',
147  'pseud_rand',
148  'plsin',
149  'poissonin',
150  'veto',
151  'injv',
152  'secl',
153  'iecl_0',
154  'iecl_1',
155  'samhem',
156  'opohem',
157  'd3',
158  'd5',
159  'd7',
160  'p3',
161  'p5',
162  'p7',
163  'typ6',
164  'cdcecl_0',
165  'cdcecl_1',
166  'cdcecl_2',
167  'cdcecl_3',
168  'c2gev_0',
169  'c2gev_1',
170  'c2gev_2',
171  'c2gev_3',
172  'cdctop_0',
173  'cdctop_1',
174  'cdctop_2',
175  'cdctop_3',
176  'cdcklm_0',
177  'cdcklm_1',
178  'seklm_0',
179  'seklm_1',
180  'ecleklm',
181  'ieklm',
182  'fwdsb',
183  'bwdsb',
184  'fwdnb',
185  'bwdnb',
186  'brlfb1',
187  'brlfb2',
188  'brlnb1',
189  'brlnb2',
190  'trkbha1',
191  'trkbha2',
192  'grlgg1',
193  'grlgg2',
194  'nimin0',
195  'nimin1',
196  'ecl_taub2b',
197  'ehigh1',
198  'ehigh2',
199  'ehigh3']
200 outputBits_expert = [
201  'fff',
202  'ffs',
203  'fss',
204  'sss',
205  'ffz',
206  'fzz',
207  'zzz',
208  'ffy',
209  'fyy',
210  'yyy',
211  'ff',
212  'fs',
213  'ss',
214  'fz',
215  'zz',
216  'fy',
217  'yy',
218  'ffo',
219  'fso',
220  'sso',
221  'fzo',
222  'fyo',
223  'ffb',
224  'fsb',
225  'ssb',
226  'fzb',
227  'fyb',
228  'aaa',
229  'aaao',
230  'aao',
231  'aab',
232  'aa',
233  'hie',
234  'lowe',
235  'lume',
236  'hade',
237  'c2',
238  'c3',
239  'c4',
240  'c5',
241  'bha3d',
242  'bhabha',
243  'bhabha_trk',
244  'bhabha_brl',
245  'bhabha_ecp',
246  'bhapur',
247  'eclmumu',
248  'bhauni',
249  'ecloflo',
250  'eclbst',
251  'g_high',
252  'g_c1',
253  'gg',
254  'eed',
255  'fed',
256  'fp',
257  'sp',
258  'zp',
259  'yp',
260  'd_5',
261  'shem',
262  'ohem',
263  'toptiming',
264  'ecltiming',
265  'cdctiming',
266  'cdcbb',
267  'mu_pair',
268  'mu_b2b',
269  'klmhit',
270  'mu_epair',
271  'mu_eb2b',
272  'eklmhit',
273  'revolution',
274  'random',
275  'bg',
276  'pls',
277  'poisson',
278  'vetout',
279  'f',
280  's',
281  'z',
282  'y',
283  'a',
284  'n1gev1',
285  'n1gev2',
286  'n1gev3',
287  'n1gev4',
288  'n2gev1',
289  'n2gev2',
290  'n2gev3',
291  'n2gev4',
292  'c2gev1',
293  'c2gev2',
294  'c2gev3',
295  'c2gev4',
296  'cdcecl1',
297  'cdcecl2',
298  'cdcecl3',
299  'cdcecl4',
300  'cdcklm1',
301  'cdcklm2',
302  'cdcklm3',
303  'cdcklm4',
304  'cdctop1',
305  'cdctop2',
306  'cdctop3',
307  'cdctop4',
308  'c1hie',
309  'c1lume',
310  'n1hie',
311  'n1lume',
312  'c3hie',
313  'c3lume',
314  'n3hie',
315  'n3lume',
316  'lml0',
317  'lml1',
318  'lml2',
319  'lml3',
320  'lml4',
321  'lml5',
322  'lml6',
323  'lml7',
324  'lml8',
325  'lml9',
326  'lml10',
327  'lml12',
328  'lml13',
329  'zzzv',
330  'yyyv',
331  'fffv',
332  'zzv',
333  'yyv',
334  'ffov',
335  'fffov',
336  'hiev',
337  'lumev',
338  'c4v',
339  'bhabhav',
340  'mu_pairv',
341  'bha3dv',
342  'bffo',
343  'bhie',
344  'ffoc',
345  'ffoc2',
346  'fffo',
347  'ff30',
348  'fs30',
349  'ss30',
350  'itsf_b2b',
351  'gggrl',
352  'ggtsf',
353  'ggbrl',
354  'bhabrl',
355  'bhamtc1',
356  'bhamtc2',
357  'bhaf',
358  'nim0',
359  'nima01',
360  'nimo01']
361 moniInbits_expert = []
362 moniOutbits_expert = []
363 
364 
365 class MakePlots(basf2.Module):
366  '''
367  Make validation histograms for trg ecl/cdc/klm
368  '''
369 
370  def set_descr(self, histogram, description, check):
371  '''
372  Sets description, check and contact to validation histogram.
373  :param h validation histogram
374  :param Descr description text
375  '''
376 
377  descr = TNamed("Description", description)
378  histogram.GetListOfFunctions().Add(descr)
379  Check = TNamed("Check", check)
380  histogram.GetListOfFunctions().Add(Check)
381  contact = TNamed("Contact", "yinjh2012@korea.ac.kr")
382  histogram.GetListOfFunctions().Add(contact)
383  Meta = TNamed("MetaOptions", "shifter")
384  histogram.GetListOfFunctions().Add(Meta)
385 
386  def set_style(self, histogram, xtitle, ytitle):
387  '''
388  Sets x-y titles, and sets histogram style.
389  :param xtitle X-axis title
390  :param xtitle Y-axis title
391  '''
392 
393  histogram.GetXaxis().SetTitle(xtitle)
394  histogram.GetYaxis().SetTitle(ytitle)
395  histogram.GetYaxis().SetTitleOffset(1.4)
396  histogram.GetXaxis().SetTitleSize(0.045)
397  histogram.GetYaxis().SetLabelSize(0.020)
398  histogram.SetLineColor(ROOT.kBlack)
399 
400  def get_relative(self, hist1, hist2, title, particle, trgbit):
401  '''
402  Get relative ratio between two hists.
403  :param hist1 numerator
404  :param hist2 denominator
405  :param title new histogram title
406  '''
407 
408  bin1 = hist1.GetXaxis().GetNbins()
409  bin2 = hist2.GetXaxis().GetNbins()
410  Xmin1 = hist1.GetXaxis().GetXmin()
411  Xmax1 = hist1.GetXaxis().GetXmax()
412  Xmin2 = hist2.GetXaxis().GetXmin()
413  Xmax2 = hist2.GetXaxis().GetXmax()
414  if bin1 != bin2 or Xmin1 != Xmin2 or Xmax1 != Xmax2:
415  print("warning: two histograms cannot be compared!!")
416  return 0
417 
418  teff = TEfficiency(hist1, hist2)
419  htitle = f'h_eff_{title}_{particle}_{trgbit}'
420  heff = TH1F(htitle, f"efficiency histogram of {particle} for {trgbit}", bin1, Xmin1, Xmax1)
421  for ibin in range(1, bin1+1):
422  binc0 = teff.GetEfficiency(ibin)
423  bine0 = (teff.GetEfficiencyErrorUp(ibin) + teff.GetEfficiencyErrorLow(ibin)) / 2.0
424  heff.SetBinContent(ibin, binc0)
425  heff.SetBinError(ibin, bine0)
426  heff.GetYaxis().SetRangeUser(0.0, 1.2)
427  return heff
428 
429  def initialize(self):
430 
431  print('start read')
432  self.tfiletfile = TFile("./TRGValidation.root", "update")
433  self.NeventNevent = 0
434 
435  m_dbinput = Belle2.PyDBObj("TRGGDLDBInputBits")
436  m_dbftdl = Belle2.PyDBObj("TRGGDLDBFTDLBits")
437  for i in range(320):
438  inbitname = m_dbinput.getinbitname(i)
439  outbitname = m_dbftdl.getoutbitname(i)
440  if inbitname in inputBits:
441  moniInbits.append({"name": inbitname, "index": i})
442  if outbitname in outputBits:
443  moniOutbits.append({"name": outbitname, "index": i})
444  if inbitname in inputBits_expert:
445  moniInbits_expert.append({"name": inbitname, "index": i})
446  if outbitname in outputBits_expert:
447  moniOutbits_expert.append({"name": outbitname, "index": i})
448 
449  mc_px = "MCParticles.m_momentum_x"
450  mc_py = "MCParticles.m_momentum_y"
451  trk2d_omega = "TRGCDC2DFinderTracks.m_omega"
452  trk2d_phi = "TRGCDC2DFinderTracks.m_phi0"
453 
454  ROOT.gStyle.SetOptStat(ROOT.kFALSE)
455 
456 
457  n_inbit_test = len(inputBits)
458  self.hist_inbithist_inbit = TH1F("hin", "trigger input bits", n_inbit_test, 0, n_inbit_test)
459  self.hist_inbithist_inbit.GetXaxis().SetRangeUser(0, n_inbit_test)
460  self.hist_inbithist_inbit.GetXaxis().SetLabelSize(0.04)
461  self.set_descrset_descr(self.hist_inbithist_inbit, "trigger input bits for shifters", "")
462  self.set_styleset_style(self.hist_inbithist_inbit, "", "Efficiency")
463 
464 
465  n_oubit_test = len(outputBits)
466  self.hist_outbithist_outbit = TH1F("hout", "trigger output bits", n_oubit_test, 0, n_oubit_test)
467  self.hist_outbithist_outbit.GetXaxis().SetRangeUser(0, n_oubit_test)
468  self.hist_outbithist_outbit.GetXaxis().SetLabelSize(0.04)
469  self.set_descrset_descr(self.hist_outbithist_outbit, "trigger ftdl bits for shifters", "")
470  self.set_styleset_style(self.hist_outbithist_outbit, "", "Efficiency")
471 
472  for i in range(n_inbit_test):
473  self.hist_inbithist_inbit.GetXaxis().SetBinLabel(
474  self.hist_inbithist_inbit.GetXaxis().FindBin(i + 0.5), inputBits[i])
475  for i in range(n_oubit_test):
476  self.hist_outbithist_outbit.GetXaxis().SetBinLabel(
477  self.hist_outbithist_outbit.GetXaxis().FindBin(i + 0.5), outputBits[i])
478 
479 
480  self.h_E_ECLh_E_ECL = TH1F("h_E_ECL", "ECL cluster energy [50 MeV]", 200, 0, 10)
481  self.set_descrset_descr(self.h_E_ECLh_E_ECL, "TRG ECL cluster energy", "")
482  self.set_styleset_style(self.h_E_ECLh_E_ECL, "ECL cluster energy (GeV)", "Events/(50 MeV)")
483 
484 
485  self.h_Esum_ECLh_Esum_ECL = TH1F("h_Esum_ECL", "sum of ECL cluster energy [50 MeV]", 200, 0, 10)
486  self.set_descrset_descr(self.h_Esum_ECLh_Esum_ECL, "sum of TRG ECL cluster energy", "")
487  self.set_styleset_style(self.h_Esum_ECLh_Esum_ECL, "sum of ECL cluster energy (GeV)", "Events/(50 MeV)")
488 
489 
490  self.h_theta_ECLh_theta_ECL = TH1F("h_theta_ECL", "TRG ECL cluster #theta [1.4 degrees]", 128, 0, 180)
491  self.set_descrset_descr(self.h_theta_ECLh_theta_ECL, "TRG ECL cluster polar angle", "unform in barrel")
492  self.set_styleset_style(self.h_theta_ECLh_theta_ECL, "TRG ECL cluster polar angle (degree)", "Events/(1.4 degree)")
493 
494 
495  self.h_thetaID_ECLh_thetaID_ECL = TH1F("h_thetaID_ECL", "ECL cluster TC ID", 610, 0, 610)
496  self.set_descrset_descr(self.h_thetaID_ECLh_thetaID_ECL, "TRG ECL cluster theta ID", "unform in barrel")
497  self.set_styleset_style(self.h_thetaID_ECLh_thetaID_ECL, "ECL cluster TC ID", "Events/(1)")
498 
499 
500  self.h_phi_ECLh_phi_ECL = TH1F("h_phi_ECL", "TRG ECL cluster phi [2.0 degrees]", 180, -180, 180)
501  self.set_descrset_descr(self.h_phi_ECLh_phi_ECL, "TRG ECL cluster phi distribution", "distribute unformly")
502  self.set_styleset_style(self.h_phi_ECLh_phi_ECL, "ECL cluster #phi (degree) ", "Events/(2.0 degrees)")
503 
504 
505  self.h_sector_BKLMh_sector_BKLM = TH1F("h_sector_BKLM", "BKLM TRG hit sector", 10, 0, 10)
506  self.set_descrset_descr(self.h_sector_BKLMh_sector_BKLM, "# of BKLM TRG sector", "peak at 0")
507  self.set_styleset_style(self.h_sector_BKLMh_sector_BKLM, "# of BKLM TRG sector", "Events/(1)")
508 
509 
510  self.h_sector_EKLMh_sector_EKLM = TH1F("h_sector_EKLM", "EKLM TRG hit sector", 10, 0, 10)
511  self.set_descrset_descr(self.h_sector_EKLMh_sector_EKLM, "# of EKLM TRG sector", "peak at 0")
512  self.set_styleset_style(self.h_sector_EKLMh_sector_EKLM, "# of EKLM TRG sector", "Events/(1)")
513 
514 
515  n_inbit_expert = len(inputBits_expert)
516  self.hist_inbit_experthist_inbit_expert = TH1F("hin_expert", "trigger input bits", n_inbit_expert, 0, n_inbit_expert)
517  self.hist_inbit_experthist_inbit_expert.GetXaxis().SetRangeUser(0, n_inbit_expert)
518  self.set_styleset_style(self.hist_inbit_experthist_inbit_expert, "", "Efficiency")
519 
520 
521  n_oubit_expert = len(outputBits_expert)
522  self.hist_outbit_experthist_outbit_expert = TH1F("hout_expert", "trigger output bits", n_oubit_expert, 0, n_oubit_expert)
523  self.hist_outbit_experthist_outbit_expert.GetXaxis().SetRangeUser(0, n_oubit_expert)
524  self.set_styleset_style(self.hist_outbit_experthist_outbit_expert, "", "Efficiency")
525 
526  for i in range(n_inbit_expert):
527  self.hist_inbit_experthist_inbit_expert.GetXaxis().SetBinLabel(
528  self.hist_inbit_experthist_inbit_expert.GetXaxis().FindBin(i + 0.5), inputBits_expert[i])
529  for i in range(n_oubit_expert):
530  self.hist_outbit_experthist_outbit_expert.GetXaxis().SetBinLabel(
531  self.hist_outbit_experthist_outbit_expert.GetXaxis().FindBin(i + 0.5), outputBits_expert[i])
532 
533  mc = "abs(MCParticles.m_pdg)==11&&MCParticles.m_status==11"
534  tree = ROOT.TChain("tree")
535  tree.Add("../TRGValidationGen.root")
536 
537 
538  self.d_wd_w = TH1F("d_w", "#Deltaw of CDC 2D finder, w = 0.00449/p_{t}", 50, -0.02, 0.02)
539  self.d_w_2d_w_2 = TH1F("d_w_2", "d_w_2", 50, -0.02, 0.02)
540  tree.Draw("({0} - 0.00449)/sqrt({1}*{1} + {2}*{2})>> d_w".format(trk2d_omega, mc_px, mc_py),
541  "MCParticles.m_pdg<0&&" + mc)
542  tree.Draw("({0} + 0.00449)/sqrt({1}*{1} + {2}*{2}) >> d_w_2".format(trk2d_omega, mc_px, mc_py),
543  "MCParticles.m_pdg>0 && " + mc)
544  self.d_wd_w.Add(self.d_w_2d_w_2)
545  self.set_descrset_descr(self.d_wd_w, "Comparison on w (=0.00449/pt) of a track between CDC 2D finder output and MC.",
546  "A clear peak at 0 with tail.")
547  self.set_styleset_style(self.d_wd_w, "#Deltaw", "Events/(0.08)")
548 
549 
550  self.d_phid_phi = TH1F("d_phi", "#Delta#phi of CDC 2D finder", 50, -0.5, 0.5)
551  tree.Draw("{0}-atan({1}/{2})>>d_phi".format(trk2d_phi, mc_py, mc_px),
552  "MCParticles.m_status==11&&abs(MCParticles.m_pdg)==11",
553  "fabs({0}-atan({1}/{2}))<{3}&&".format(trk2d_phi, mc_py, mc_px, PI) + mc)
554 
555  self.d_phi_2d_phi_2 = TH1F("d_phi_2", "d_phi_2", 50, -0.5, 0.5)
556  tree.Draw("{0}-atan({1}/{2})-{3}>>d_phi_2".format(trk2d_phi, mc_py, mc_px, PI),
557  "MCParticles.m_status==11&&abs(MCParticles.m_pdg)==11",
558  "{0}-atan({1}/{2})>={3} &&".format(trk2d_phi, mc_py, mc_px, PI) + mc)
559 
560  self.d_phi_3d_phi_3 = TH1F("d_phi_3", "d_phi_3", 50, -0.5, 0.5)
561  tree.Draw("{0}-atan({1}/{2})+{3}>>d_phi_2".format(trk2d_phi, mc_py, mc_px, PI),
562  "MCParticles.m_status==11&&abs(MCParticles.m_pdg)==11",
563  "{0}-atan({1}/{2})<=-{3}&&".format(trk2d_phi, mc_py, mc_px, PI) + mc)
564 
565  self.d_phid_phi.Add(self.d_phi_2d_phi_2)
566  self.d_phid_phi.Add(self.d_phi_3d_phi_3)
567  self.set_descrset_descr(self.d_phid_phi, "Comparison on phi_i of a track between CDC 2D finder output and MC.",
568  "A Gaussian peak at 0.")
569  self.set_styleset_style(self.d_phid_phi, "#Delta#phi [rad]", "Events/(0.02 rad)")
570 
571 
572  self.d_z0_3dd_z0_3d = TH1F("d_z0_3d", "#Deltaz0 of CDC 3D fitter", 60, -30, 30)
573  tree.Draw("TRGCDC3DFitterTracks.m_z0-MCParticles.m_productionVertex_z>>d_z0_3d", mc)
574  self.set_descrset_descr(self.d_z0_3dd_z0_3d, "Comparison on z0 of a track between CDC 2D fitter output and MC.",
575  "A Gaussian peak at 0 with small tail.")
576  self.set_styleset_style(self.d_z0_3dd_z0_3d, "#Deltaz0 [cm]", "Events/(1 cm)")
577 
578 
579  self.d_z0_nnd_z0_nn = TH1F("d_z0_nn", "#Deltaz0 of CDC Neuro", 60, -30, 30)
580  tree.Draw("TRGCDCNeuroTracks.m_z0-MCParticles.m_productionVertex_z>>d_z0_nn", mc)
581  self.set_descrset_descr(self.d_z0_nnd_z0_nn, "Comparison on z0 of a track between CDC Neuro output and MC.",
582  "A Gaussian peak at 0 with small tail.")
583  self.set_styleset_style(self.d_z0_nnd_z0_nn, "#Deltaz0 [cm]", "Events/(1 cm)")
584 
585 
586  self.d_E_ECLd_E_ECL = TH1F("d_E_ECL", "#DeltaE of ECL clustering", 100, -6, 6)
587  tree.Draw("TRGECLClusters.m_edep-MCParticles.m_energy>>d_E_ECL", mc)
588  self.set_descrset_descr(self.d_E_ECLd_E_ECL, "Comparison on deposit energy between ECL cluster output and MC.",
589  "A peak around -0.5 ~ 0 with a tail toward -6.")
590  self.set_styleset_style(self.d_E_ECLd_E_ECL, "#DeltaE [GeV]", "Events/(0.12 GeV)")
591 
592  def beginRun(self):
593  self.hist_inbithist_inbit.Reset()
594  self.hist_outbithist_outbit.Reset()
595  self.hist_inbit_experthist_inbit_expert.Reset()
596  self.hist_outbit_experthist_outbit_expert.Reset()
597  self.h_Esum_ECLh_Esum_ECL.Reset()
598  self.h_E_ECLh_E_ECL.Reset()
599  self.h_theta_ECLh_theta_ECL.Reset()
600  self.h_thetaID_ECLh_thetaID_ECL.Reset()
601  self.h_phi_ECLh_phi_ECL.Reset()
602 
603  def event(self):
604 
605  trg_summary = Belle2.PyStoreObj("TRGSummary")
606  clusters = Belle2.PyStoreArray("TRGECLClusters")
607  klmSummary = Belle2.PyStoreObj("KLMTrgSummary")
608 
609  for bit in moniInbits:
610  binIndex = inputBits.index(bit["name"])
611  bitIndex = bit["index"]
612  if trg_summary.testInput(bitIndex):
613  self.hist_inbithist_inbit.Fill(binIndex + 0.5)
614  for bit in moniOutbits:
615  binIndex = outputBits.index(bit["name"])
616  bitIndex = bit["index"]
617  if trg_summary.testFtdl(bitIndex):
618  self.hist_outbithist_outbit.Fill(binIndex + 0.5)
619 
620  for bit in moniInbits_expert:
621  binIndex = inputBits_expert.index(bit["name"])
622  bitIndex = bit["index"]
623  if trg_summary.testInput(bitIndex):
624  self.hist_inbit_experthist_inbit_expert.Fill(binIndex + 0.5)
625  for bit in moniOutbits_expert:
626  binIndex = outputBits_expert.index(bit["name"])
627  bitIndex = bit["index"]
628  if trg_summary.testFtdl(bitIndex):
629  self.hist_outbit_experthist_outbit_expert.Fill(binIndex + 0.5)
630 
631  etot = 0
632  for cluster in clusters:
633  x = cluster.getPositionX()
634  y = cluster.getPositionY()
635  z = cluster.getPositionZ()
636  e = cluster.getEnergyDep()
637  self.h_E_ECLh_E_ECL.Fill(e)
638  etot += e
639  vec = ROOT.TVector3(x, y, z)
640  self.h_theta_ECLh_theta_ECL.Fill(vec.Theta() * Fac)
641  self.h_thetaID_ECLh_thetaID_ECL.Fill(cluster.getMaxTCId())
642  self.h_phi_ECLh_phi_ECL.Fill(vec.Phi() * Fac)
643 
644  if etot > 0:
645  self.h_Esum_ECLh_Esum_ECL.Fill(etot)
646 
647  self.h_sector_BKLMh_sector_BKLM.Fill(klmSummary.getBKLM_n_trg_sectors())
648  self.h_sector_EKLMh_sector_EKLM.Fill(klmSummary.getEKLM_n_trg_sectors())
649 
650  self.NeventNevent = self.NeventNevent + 1
651 
652  def endRun(self):
653  print("end")
654 
655  def terminate(self):
656 
657  bits = ['ty_0', 'ehigh', 'clst_0', 'klm_0']
658  for part in ["e", "mu"]:
659  h_gen_p = self.tfiletfile.Get(f"{part}_all/p")
660  h_gen_E = self.tfiletfile.Get(f"{part}_all/clusterE")
661  h_gen_theta = self.tfiletfile.Get(f"{part}_all/theta")
662  h_gen_phi = self.tfiletfile.Get(f"{part}_all/phi")
663 
664  for bit in bits:
665  h_p = self.tfiletfile.Get(f"{part}_{bit}/p")
666  h_E = self.tfiletfile.Get(f"{part}_{bit}/clusterE")
667  h_theta = self.tfiletfile.Get(f"{part}_{bit}/theta")
668  h_phi = self.tfiletfile.Get(f"{part}_{bit}/phi")
669 
670  h_eff_p = self.get_relativeget_relative(h_p, h_gen_p, 'p', part, bit)
671  self.set_styleset_style(h_eff_p, "Momentum", "Efficiency")
672  if (bit == 'ty_0' or bit == 'klm_0') and part == 'mu':
673  self.set_descrset_descr(h_eff_p, "Momentum dependent efficiency for shifters", "")
674 
675  h_eff_E = self.get_relativeget_relative(h_E, h_gen_E, 'E', part, bit)
676  self.set_styleset_style(h_eff_E, "Energy", "Efficiency")
677  if (bit == 'ehigh' or bit == 'clst_0') and part == 'e':
678  self.set_descrset_descr(h_eff_E, "Energy dependent efficiency for shifters", "")
679 
680  h_eff_theta = self.get_relativeget_relative(h_theta, h_gen_theta, 'theta', part, bit)
681  self.set_styleset_style(h_eff_theta, "Polar angle", "Efficiency")
682 
683  h_eff_phi = self.get_relativeget_relative(h_phi, h_gen_phi, 'phi', part, bit)
684  self.set_styleset_style(h_eff_phi, "phi angle", "Efficiency")
685 
686  h_eff_p.Write()
687  h_eff_E.Write()
688  h_eff_theta.Write()
689  h_eff_phi.Write()
690 
691  self.hist_inbithist_inbit.Scale(1./self.NeventNevent)
692  self.hist_outbithist_outbit.Scale(1.0/self.NeventNevent)
693  self.hist_inbit_experthist_inbit_expert.Scale(1./self.NeventNevent)
694  self.hist_outbit_experthist_outbit_expert.Scale(1.0/self.NeventNevent)
695 
696  self.hist_inbithist_inbit.Write()
697  self.hist_outbithist_outbit.Write()
698  self.hist_inbit_experthist_inbit_expert.Write()
699  self.hist_outbit_experthist_outbit_expert.Write()
700  self.h_Esum_ECLh_Esum_ECL.Write()
701  self.h_E_ECLh_E_ECL.Write()
702  self.h_theta_ECLh_theta_ECL.Write()
703  self.h_thetaID_ECLh_thetaID_ECL.Write()
704  self.h_phi_ECLh_phi_ECL.Write()
705  self.h_sector_BKLMh_sector_BKLM.Write()
706  self.h_sector_EKLMh_sector_EKLM.Write()
707  self.d_wd_w.Write()
708  self.d_phid_phi.Write()
709  self.d_z0_3dd_z0_3d.Write()
710  self.d_z0_nnd_z0_nn.Write()
711  self.d_E_ECLd_E_ECL.Write()
712  self.tfiletfile.Close()
713 
714 
715 
716 
717 mypath = basf2.Path() # create a new path
718 
719 # add input data and ParticleLoader modules to the path
720 ma.inputMdst("../TRGValidationGen.root", path=mypath)
721 
722 
723 stdc.stdMu(listtype='all', path=mypath)
724 stdc.stdE(listtype='all', path=mypath)
725 
726 ma.matchMCTruth('mu+:all', path=mypath)
727 ma.matchMCTruth('e+:all', path=mypath)
728 
729 ma.cutAndCopyList('mu+:true', 'mu+:all', cut='isSignal>0 and mcPrimary>0', path=mypath)
730 ma.cutAndCopyList('e+:true', 'e+:all', cut='isSignal>0 and mcPrimary>0', path=mypath)
731 
732 
733 pvars = ["p", "theta", "phi"]
734 dvars = vu.create_aliases(pvars, 'daughter(0,{variable})', "emu")
735 
736 bits = ['ty_0', 'ehigh', 'clst_0', 'klm_0']
737 varlists = []
738 for bit in bits:
739  inb_cut = f"L1Input({bit})>0"
740  ma.cutAndCopyList(f'mu+:{bit}', 'mu+:true', cut=inb_cut, path=mypath)
741  ma.cutAndCopyList(f'e+:{bit}', 'e+:true', cut=inb_cut, path=mypath)
742 
743  ma.variablesToHistogram(decayString=f'mu+:{bit}',
744  variables=[('p', 40, 1.0, 3.0),
745  ('clusterE', 40, 1.0, 3.0),
746  ('theta', 40, 0.5, 2.5),
747  ('phi', 40, -3.141593, 3.141593)],
748  filename='TRGValidation.root',
749  directory=f'mu_{bit}',
750  path=mypath
751  )
752 
753  ma.variablesToHistogram(decayString=f'e+:{bit}',
754  variables=[('p', 40, 1.0, 3.0),
755  ('clusterE', 40, 1.0, 3.0),
756  ('theta', 40, 0.5, 2.5),
757  ('phi', 40, -3.141593, 3.141593)],
758  filename='TRGValidation.root',
759  directory=f'e_{bit}',
760  path=mypath
761  )
762 
763 
764 ma.variablesToHistogram(decayString='mu+:true',
765  variables=[('p', 40, 1.0, 3.0),
766  ('clusterE', 40, 1.0, 3.0),
767  ('theta', 40, 0.5, 2.5),
768  ('phi', 40, -3.141593, 3.141593)],
769  filename='TRGValidation.root',
770  directory='mu_all',
771  path=mypath
772  )
773 
774 ma.variablesToHistogram(decayString='e+:true',
775  variables=[('p', 40, 1.0, 3.0),
776  ('clusterE', 40, 1.0, 3.0),
777  ('theta', 40, 0.5, 2.5),
778  ('phi', 40, -3.141593, 3.141593)],
779  filename='TRGValidation.root',
780  directory='e_all',
781  path=mypath
782  )
783 
784 
785 # process the data
786 basf2.process(mypath)
787 
788 # INput
789 main = basf2.create_path()
790 root_input = basf2.register_module("RootInput")
791 root_input.param("inputFileName", "../TRGValidationGen.root")
792 main.add_module(root_input)
793 main.add_module(MakePlots())
794 main.add_module('Progress')
795 basf2.process(main)
Class to access a DBObjPtr from Python.
Definition: PyDBObj.h:50
A (simplified) python wrapper for StoreArray.
Definition: PyStoreArray.h:72
a (simplified) python wrapper for StoreObjPtr.
Definition: PyStoreObj.h:67
h_Esum_ECL
validation histogram
h_sector_BKLM
validation histogram
d_w
validation histogram
hist_outbit
validation histogram
h_E_ECL
validation histogram
hist_outbit_expert
validation histogram for experts
d_E_ECL
validation histogram
h_thetaID_ECL
validation histogram
def set_style(self, histogram, xtitle, ytitle)
hist_inbit_expert
validation histogram for experts
hist_inbit
validation histogram
d_z0_nn
validation histogram
h_theta_ECL
validation histogram
def set_descr(self, histogram, description, check)
def get_relative(self, hist1, hist2, title, particle, trgbit)
h_sector_EKLM
validation histogram
h_phi_ECL
validation histogram
d_z0_3d
validation histogram
d_phi
validation histogram