Belle II Software  release-06-02-00
ReadTSIMInfo.py
1 #!/usr/bin/env python3
2 # -*- coding: utf-8 -*-
3 
4 
11 
12 import basf2 as b2
13 import ROOT
14 import math
15 from ROOT import Belle2
16 import modularAnalysis as ma
17 from array import array
18 import sys # get argv
19 
20 from heapq import nlargest
21 from effCalculation import EffCalculation
22 
23 PI = 3.1415926
24 Fac = 180.0 / PI
25 
26 
27 argvs = sys.argv
28 argn = len(argvs)
29 if argn != 3:
30  sys.exit("ztsim02.py> # of arg is strange. Exit.")
31 
32 
33 def Count_mcpart(parts):
34  npart = [0, 0]
35  dec_cdc = [17, 150]
36  dec_ecl = [14.1, 155]
37  for part in parts:
38  pdglist = [11, 13, 2212, 321, 211]
39  pdg = math.fabs(part.getPDG())
40  charge = part.getCharge()
41  theta = part.getMomentum().Theta() * Fac
42  if theta > dec_cdc[0] and theta < dec_cdc[1]:
43  if pdg in pdglist:
44  npart[0] += 1
45 
46  if math.fabs(pdg) == 22:
47  if theta > dec_ecl[0] and theta < dec_ecl[1]:
48  npart[1] += 1
49  return npart
50 
51 
52 def Count_part_inDetector(parts):
53  npart = [0, 0]
54  for part in parts:
55  pdg = part.getPDG()
56  charge = part.getCharge()
57  indec = part.hasSeenInDetector(Belle2.Const.DetectorSet(Belle2.Const.ECL))
58  if indec:
59  if math.fabs(charge) > 0:
60  npart[0] += 1
61  elif math.fabs(pdg) == 22:
62  npart[1] += 1
63  return npart
64 
65 
66 def Theta_MCTruth(part):
67  theta = part.getMomentum().Theta() * Fac
68  return theta
69 
70 
71 def NeutralCluster(clusters):
72  neutral_cluster = []
73  for cluster in clusters:
74  if cluster.getRelatedFrom('TRGNNTracks'):
75  continue
76  else:
77  neutral_cluster.append(cluster)
78  return neutral_cluster
79 
80 
81 def Vec_Cluster(cluster, CMS=False):
82  x = cluster.getPositionX()
83  y = cluster.getPositionY()
84  z = cluster.getPositionZ()
85  e = cluster.getEnergyDep()
86  vec = ROOT.TVector3(x, y, z)
87  vec_unit = vec.Unit()
88  vec_unit *= e
89  v_mom = ROOT.TLorentzVector(vec_unit, e)
90  if CMS:
91  new_vec = trans.rotateLabToCms() * v_mom
92  else:
93  new_vec = v_mom
94  new_theta = new_vec.Theta() * Fac
95  new_phi = new_vec.Phi() * Fac
96  if(new_phi < 0):
97  new_phi += 2 * PI
98  new_e = new_vec.E()
99  NVC = [new_e, new_theta, new_phi]
100  return NVC
101 
102 
103 def MatchTrk_Info(matchtrk, clusters):
104  RV = [0., 0., 0., 0., 0., 0.]
105  RV[0] = matchtrk.getPt()
106  RV[1] = matchtrk.getPz()
107  RV[2] = (PI - matchtrk.getTheta()) * Fac
108  RV[3] = matchtrk.getPhi() * Fac
109  clusterId = matchtrk.getECLClusterId()
110  for cluster in clusters:
111  clusterId2 = cluster.getClusterId()
112  if clusterId == clusterId2:
113  vec_ = Vec_Cluster(cluster, False)
114  vec_cms = Vec_Cluster(cluster, False)
115  RV[4] = vec_[0]
116  RV[5] = vec_cms[0]
117  return RV
118 
119 
120 def Max_Cluster(neuclus, CMS, n):
121  E = []
122  for neuclu in neuclus:
123  NV = Vec_Cluster(neuclu, CMS)
124  E.append(NV)
125  if E and len(E) >= (n + 1):
126  seqen = nlargest(n + 1, E, key=lambda item: item[0])
127  return seqen[n]
128  else:
129  return [-1, -999, -999]
130 
131 
132 def Etot_Cluster(eclcluster):
133  Elist = []
134  for clu in eclcluster:
135  e = clu.getEnergyDep()
136  Elist.append(e)
137  return sum(Elist)
138 
139 
140 def Max_DeltPhi_trk(trk_2d_list):
141  Delt = []
142  for i, trk in enumerate(trk_2d_list):
143  phi1 = (trk.getPhi0()) * Fac
144  for j, trk2 in enumerate(trk_2d_list):
145  if j <= i:
146  continue
147  else:
148  phi2 = (trk2.getPhi0()) * Fac
149  Delt_phi = math.fabs(phi1 - phi2)
150  Delt.append(Delt_phi)
151  return max(Delt or [-100.0])
152 
153 
154 def BhabhaVeto1(matchtrks):
155  Delt = []
156  if len(matchtrks) <= 1:
157  return Delt
158  for i, matchtrk1 in enumerate(matchtrks):
159  if i == len(matchtrks) - 1:
160  continue
161  cdctrk1 = matchtrk1.getRelatedTo('TRGNNTracks')
162  cluster1 = matchtrk1.getRelatedTo('TRGECLClusters')
163  if cdctrk1 and cluster1:
164  tanLam1 = cdctrk1.getTanLambda()
165  theta1 = math.acos(tanLam1 / math.sqrt(1. + tanLam1 * tanLam1)) * Fac
166  phi1 = (cdctrk1.getPhi0()) * Fac
167  if phi1 < 0.:
168  phi1 += 2 * PI
169  e1 = cluster1.getEnergyDep()
170  for j, matchtrk2 in enumerate(matchtrks):
171  if j <= i:
172  continue
173  else:
174  cdctrk2 = matchtrk2.getRelatedTo('TRGNNTracks')
175  cluster2 = matchtrk2.getRelatedTo('TRGECLClusters')
176  if cdctrk2 and cluster2:
177  tanLam2 = cdctrk2.getTanLambda()
178  theta2 = math.acos(tanLam2 / math.sqrt(1. + tanLam2 * tanLam2)) * Fac
179  phi2 = (cdctrk2.getPhi0()) * Fac
180  if phi2 < 0.:
181  phi2 += 2 * PI
182  e2 = cluster2.getEnergyDep()
183  Delt_theta = theta1 + theta2 - 180
184  Delt_phi = math.fabs(phi1 - phi2) - 180
185  DeltArray = [Delt_theta, Delt_phi, e1, e2, e1 + e2]
186  Delt.append(DeltArray)
187  if len(Delt) >= 1:
188  return Delt
189  else:
190  return [-999., -999., -999., -999., -999]
191 
192 
193 def SBhabhaVeto(matchtrks):
194  Elist = []
195  if len(matchtrks) <= 1:
196  return Elist
197  for i, matchtrk1 in enumerate(matchtrks):
198  cdctrk1 = matchtrk1.getRelatedTo('TRGNNTracks')
199  cluster1 = matchtrk1.getRelatedTo('TRGECLClusters')
200  if cdctrk1 and cluster1:
201  e1 = cluster1.getEnergyDep()
202  Elist.append([e1])
203  if len(Elist) >= 1:
204  return Elist
205  else:
206  return [-999.]
207 
208 
209 def eclBhabhaVeto(eclclusters):
210  neclBhabha = []
211  for i, cluster1 in enumerate(eclclusters):
212  if i == len(eclclusters) - 1:
213  continue
214  e1 = cluster1.getEnergyDep()
215  x1 = cluster1.getPositionX()
216  y1 = cluster1.getPositionY()
217  z1 = cluster1.getPositionZ()
218  vec1 = ROOT.TVector3(x1, y1, z1)
219  theta1 = vec1.Theta()
220  phi1 = vec1.Phi()
221  if(phi1 < 0):
222  phi1 += 2 * PI
223  for j, cluster2 in enumerate(eclclusters):
224  e2 = cluster2.getEnergyDep()
225  if j <= i:
226  continue
227  x2 = cluster2.getPositionX()
228  y2 = cluster2.getPositionY()
229  z2 = cluster2.getPositionZ()
230  vec2 = ROOT.TVector3(x2, y2, z2)
231  theta2 = vec2.Theta()
232  phi2 = vec2.Phi()
233  if phi2 < 0.0:
234  phi2 += 2 * PI
235  delt_theta = math.fabs((theta1 + theta2) * Fac - 180.0)
236  delt_phi = math.fabs(phi1 - phi2) * Fac - 180.0
237  if theta1 < theta2:
238  efr = e1
239  ebr = e2
240  else:
241  efr = e2
242  ebr = e1
243  bhacl = [delt_theta, delt_phi, efr, ebr, efr + ebr]
244  neclBhabha.append(bhacl)
245  return neclBhabha
246 
247 
248 def Max_DeltPhi_cluster(eclclusters):
249  eclphi_col = []
250  for i, cluster1 in enumerate(eclclusters):
251  e1 = cluster1.getEnergyDep()
252  if e1 < 0.1:
253  continue
254  x1 = cluster1.getPositionX()
255  y1 = cluster1.getPositionY()
256  z1 = cluster1.getPositionZ()
257  vec1 = ROOT.TVector3(x1, y1, z1)
258  phi1 = vec1.Phi()
259  if(phi1 < 0):
260  phi1 += 2 * PI
261  for j, cluster2 in enumerate(eclclusters):
262  e2 = cluster2.getEnergyDep()
263  if j <= i:
264  continue
265  if e2 < 0.1:
266  continue
267  x2 = cluster2.getPositionX()
268  y2 = cluster2.getPositionY()
269  z2 = cluster2.getPositionZ()
270  vec2 = ROOT.TVector3(x2, y2, z2)
271  phi2 = vec2.Phi()
272  if phi2 < 0.0:
273  phi2 += 2 * PI
274  delt_phi = math.fabs(phi1 - phi2) * 180.0 / 3.1415926
275  eclphi_col.append(delt_phi)
276  return max(eclphi_col or [-1.0])
277 
278 
279 def Time_Window(clusters, eventtime):
280  new_clusters = []
281  energy = []
282  for ev in eventtime:
283  event_time = ev.m_eventtiming
284  event_tot = ev.m_etot
285  energy.append([event_time, event_tot])
286  if energy:
287  tmp = max(energy, key=lambda item: item[1])
288  event_time = tmp[0]
289  else:
290  event_time = -999999.
291  for cluster in clusters:
292  ctime_ori = cluster.getTimeAve()
293  ctime = ctime_ori - event_time
294  if math.fabs(ctime) < 100:
295  new_clusters.append(cluster)
296  return new_clusters
297 
298 
299 def Cluster_Threshold(clusters, threshold, CMS):
300  new_clusters = []
301  for cluster in clusters:
302  eng = 0
303  if CMS:
304  newv = Vec_Cluster(cluster, CMS)
305  eng = newv[0]
306  else:
307  eng = cluster.getEnergyDep()
308  if eng > threshold:
309  new_clusters.append(cluster)
310  return new_clusters
311 
312 
313 def Back_to_Back(clusters1, clusters2):
314  npai = 0
315  for cluster1 in clusters1:
316  cid1 = cluster1.getClusterId()
317  vec1 = Vec_Cluster(cluster1, False)
318  theta1 = vec1[1]
319  phi1 = vec1[2]
320  for cluster2 in clusters2:
321  cid2 = cluster2.getClusterId()
322  if cid1 == cid2:
323  continue
324  vec2 = Vec_Cluster(cluster2, False)
325  theta2 = vec2[1]
326  phi2 = vec2[2]
327  delttheta = math.fabs(theta1 + theta2 - 180)
328  deltphi = math.fabs(math.fabs(phi1 - phi2) - 180)
329  if delttheta < 100 and deltphi < 100:
330  npai += 1
331  return npai
332 
333 
334 def Time_Cluster(clusters, eventtime):
335  time_list = []
336  energy = []
337  event_time = 0.
338  for ev in eventtime:
339  event_time = ev.m_eventtiming
340  event_tot = ev.m_etot
341  energy.append([event_time, event_tot])
342  if energy:
343  tmp = max(energy, key=lambda item: item[1])
344  event_time = tmp[0]
345  for cluster in clusters:
346  ctime_ori = cluster.getTimeAve()
347  ctime = ctime_ori - event_time
348  time_list.append(ctime)
349  return time_list
350 
351 
352 def PrintBranchDef():
353  print('\n')
354  print('ntrk_2dfinder: # 2d finder track')
355  print('ntrk_2dfitter: # 2d fitter track')
356  print('ntrk_3dfitter: # 3d fitter track')
357  print('ntrk_NN : # Neuro network track')
358  print('ntrk_2Dmatch : # 2d finder track w/ associated ecl cluster')
359  print('ntrk_3Dmatch : # 3d finder track w/ associated ecl cluster')
360  print('ntrk_klm : # KLM track')
361  print('nhit_klm : # KLM hit', '\n')
362 
363  print('ncluster: # ecl cluster')
364  print('ncluster_1000b: # ecl cluster with threshold >1.0GeV, exclude TC ID 1,2, 17')
365  print('ncluster_2000e: # ecl cluster with threshold >2.0GeV in TC ID 1, 17')
366  print('max_cluster[3]: [energy, theta, phi] of the largest energetic ecl cluster')
367  print('smax_cluster[3]: [energy, theta, phi] of the secondary energetic ecl cluster')
368  print('ncluster_neutral: # ecl cluster w/o associated cdc track')
369  print('max_cluster_neutral[3]: [energy, theta, phi] of the largest energetic ecl neutral cluster')
370  print('smax_cluster_neutral[3]: [energy, theta, phi] of the secondary energetic ecl neutral cluster')
371  print('time_cluster: the cluster timing obtained with cluster.timing-event.timing', '\n')
372 
373  print('nbbc: # back to back cluster pairs')
374  print('nbbtc: # back to back track and cluster pairs')
375  print('bhabha: bhabha veto logic, 1: bhabha, 0: non bhabha')
376  print('sbhabha: bhabha with single track veto logic, 1: bhabha, 0: non bhabha')
377  print('eclbhabha: eclbhabha veto logic, 1: bhabha, 0: non bhabha')
378  print('bhabha_var[5]: variables used in bhabha logic, [Delt_theta, Delt_phi, E1, E2, E1+E2]')
379  print(' For two tracks: Delt_theta: theta1+theta2-180, Delt_phi: |phi1-phi2|-180')
380  print(' E1, E2 are the ecl clusters energy accociated with the two tracks')
381  print('eclbhabha_var[5]: variables used in eclbhabha logic, [Delt_theta, Delt_phi, E1, E2, E1+E2]')
382  print(' For two eclclusers: Delt_theta: theta1+theta2-180, Delt_phi: |phi1-phi2|-180')
383  print(' E1, E2 are the ecl clusters energy')
384  print('\n')
385 
386 
387 class CreateLogics(b2.Module):
388  """This module is to calculate some variables which are useful for the trigger development"""
389 
390  file = ROOT.TFile(argvs[2], 'recreate')
391 
392  tgrl = ROOT.TTree('tgrl', 'tree with GRL_Logic')
393 
394  ntrk_2dfinder_t = array('i', [-1])
395 
396  ntrk_2dfitter_t = array('i', [-1])
397 
398  ntrk_3dfitter_t = array('i', [-1])
399 
400  ntrk_NN_t = array('i', [-1])
401 
402  ntrk_2Dmatch_t = array('i', [-1])
403 
404  ntrk_3Dmatch_t = array('i', [-1])
405 
406  max_deltphi_2dfinder_t = array('f', [0.0])
407 
408  cpair_t = array('i', 8 * [0])
409 
410  ncomp_clu = 3
411 
412  etot_t = array('f', [0.0])
413 
414  ncluster_1000b_t = array('i', [-1])
415 
416  ncluster_2000e_t = array('i', [-1])
417 
418  ncluster_t = array('i', [-1])
419 
420  ncluster_neutral_t = array('i', [-1])
421 
422  max_cluster_neutral_t = array('f', ncomp_clu * [0.0])
423 
424  max_cms_cluster_neutral_t = array('f', ncomp_clu * [0.0])
425 
426  max_cluster_t = array('f', ncomp_clu * [0.0])
427 
428  max_cms_cluster_t = array('f', ncomp_clu * [0.0])
429 
430  smax_cluster_neutral_t = array('f', ncomp_clu * [0.0])
431 
432  smax_cms_cluster_neutral_t = array('f', ncomp_clu * [0.0])
433 
434  smax_cluster_t = array('f', ncomp_clu * [0.0])
435 
436  smax_cms_cluster_t = array('f', ncomp_clu * [0.0])
437 
438  max_deltphi_cluster_t = array('f', [0.0])
439 
440  time_cluster_t = array('f', 100 * [-99999.])
441 
442  ntrk_klm_t = array('i', [-1])
443 
444  nhit_klm_t = array('i', [-1])
445 
446  npair_tc_t = array('i', [-1])
447 
448  npair_cc_t = array('i', [-1])
449 
450  nbha_var = 5
451 
452  bhabha_t = array('i', [0])
453 
454  sbhabha_t = array('i', [0])
455 
456  eclbhabha_t = array('i', [0])
457 
458  bhabha_var_t = array('f', nbha_var * [0.0])
459 
460  eclbhabha_var_t = array('f', nbha_var * [0.0])
461 
462  # track
463  tgrl.Branch('ntrk_2dfinder', ntrk_2dfinder_t, 'ntrk_2dfinder/i')
464  tgrl.Branch('ntrk_2dfitter', ntrk_2dfitter_t, 'ntrk_2dfitter/i')
465  tgrl.Branch('ntrk_3dfitter', ntrk_3dfitter_t, 'ntrk_3dfitter/i')
466  tgrl.Branch('ntrk_NN', ntrk_NN_t, 'ntrk_NN/i')
467  tgrl.Branch('ntrk_2Dmatch', ntrk_2Dmatch_t, 'ntrk_2Dmatch/i')
468  tgrl.Branch('ntrk_3Dmatch', ntrk_3Dmatch_t, 'ntrk_3Dmatch/i')
469  tgrl.Branch('npair_tc', npair_tc_t, 'npair_tc/i')
470  tgrl.Branch('npair_cc', npair_cc_t, 'npair_cc/i')
471  tgrl.Branch('max_deltphi_2dfinder', max_deltphi_2dfinder_t, 'max_deltphi_2dfinder/f')
472 
473  # cluster
474  tgrl.Branch('etot', etot_t, 'etot/f')
475  tgrl.Branch('ncluster', ncluster_t, 'ncluster/i')
476  tgrl.Branch('ncluster_1000b', ncluster_1000b_t, 'ncluster_1000b/i')
477  tgrl.Branch('ncluster_2000e', ncluster_2000e_t, 'ncluster_2000e/i')
478  tgrl.Branch('ncluster_neutral', ncluster_neutral_t, 'ncluster_neutral/i')
479  tgrl.Branch('max_cluster_neutral', max_cluster_neutral_t, 'max_cluster_neutral[3]/f')
480  # tgrl.Branch('max_cms_cluster_neutral', max_cms_cluster_neutral_t, 'max_cms_cluster_neutral[3]/f')
481  tgrl.Branch('max_cluster', max_cluster_t, 'max_cluster[3]/f')
482  # tgrl.Branch('max_cms_cluster', max_cms_cluster_t, 'max_cms_cluster[3]/f')
483  tgrl.Branch('smax_cluster_neutral', smax_cluster_neutral_t, 'smax_cluster_neutral[3]/f')
484  # tgrl.Branch('smax_cms_cluster_neutral', smax_cms_cluster_neutral_t, 'smax_cms_cluster_neutral[3]/f')
485  tgrl.Branch('smax_cluster', smax_cluster_t, 'smax_cluster[3]/f')
486  # tgrl.Branch('smax_cms_cluster', smax_cms_cluster_t, 'smax_cms_cluster[3]/f')
487  tgrl.Branch('max_deltphi_cluster', max_deltphi_cluster_t, 'max_deltphi_cluster/f')
488  tgrl.Branch('time_cluster', time_cluster_t, 'time_cluster[100]/f')
489 
490  # klm
491  tgrl.Branch('ntrk_klm', ntrk_klm_t, 'ntrk_klm/i')
492  tgrl.Branch('nhit_klm', nhit_klm_t, 'nhit_klm/i')
493 
494  # bhabha
495  tgrl.Branch('bhabha', bhabha_t, 'bhabha/i')
496  tgrl.Branch('sbhabha', sbhabha_t, 'sbhabha/i')
497  tgrl.Branch('eclbhabha', eclbhabha_t, 'eclbhabha/i')
498  tgrl.Branch('bhabha_var', bhabha_var_t, 'bhabha_var[5]/f')
499  tgrl.Branch('eclbhabha_var', eclbhabha_var_t, 'eclbhabha_var[5]/f')
500 
501  tgrl.Branch('cpair', cpair_t, 'cpair[8]/i')
502 
503  def event(self):
504  """
505  MCPart = Belle2.PyStoreArray('MCParticles')
506  if len(MCPart)>=4:
507  self.eplus_t[0]=Theta_MCTruth(MCPart[2])
508  self.eminus_t[0]=Theta_MCTruth(MCPart[3])
509  if self.eplus_t[0] > self.eminus_t[0]:
510  self.efrd_t[0]=self.eminus_t[0]
511  self.ebkd_t[0]=self.eplus_t[0]
512  else:
513  self.efrd_t[0]=self.eplus_t[0]
514  self.ebkd_t[0]=self.eminus_t[0]
515  count_part = Count_part_inDetector(MCPart)
516  self.n_par_t[0]=count_part[0]
517  self.n_par_t[1]=count_part[1]
518  count_mcpart = Count_mcpart(MCPart)
519  self.n_mcpar_t[0]=count_mcpart[0]
520  self.n_mcpar_t[1]=count_mcpart[1]
521  """
522  trk_2d_finder = Belle2.PyStoreArray('TRG2DFinderTracks')
523  self.ntrk_2dfinder_tntrk_2dfinder_t[0] = len(trk_2d_finder)
524  trk_2d_fitter = Belle2.PyStoreArray('TRG2DFitterTracks')
525  self.ntrk_2dfitter_tntrk_2dfitter_t[0] = len(trk_2d_fitter)
526  trk_3d_fitter = Belle2.PyStoreArray('TRG3DFitterTracks')
527  self.ntrk_3dfitter_tntrk_3dfitter_t[0] = len(trk_3d_fitter)
528  tracks = Belle2.PyStoreArray('TRGNNTracks')
529  self.ntrk_NN_tntrk_NN_t[0] = len(tracks)
530 
531  # clusters
532  clusters_original = Belle2.PyStoreArray('TRGECLClusters')
533  eventtime = Belle2.PyStoreArray('TRGECLTrgs')
534  # get clusters list in time window [-100,100]
535  clusters = Time_Window(clusters_original, eventtime)
536  self.ncluster_tncluster_t[0] = len(clusters)
537  # get the cluster timing before time window requirement
538  time_clu_list = Time_Cluster(clusters_original, eventtime)
539  for i in range(len(time_clu_list)):
540  if i < 100:
541  self.time_cluster_ttime_cluster_t[i] = time_clu_list[i]
542 
543  clusters_100 = Cluster_Threshold(clusters, 0.1, False)
544  clusters_100_cms = Cluster_Threshold(clusters, 0.1, False)
545  clusters_300_cms = Cluster_Threshold(clusters, 0.3, False)
546  clusters_400_cms = Cluster_Threshold(clusters, 0.4, False)
547  clusters_500_cms = Cluster_Threshold(clusters, 0.5, False)
548  clusters_700_cms = Cluster_Threshold(clusters, 0.7, False)
549  clusters_1000_cms = Cluster_Threshold(clusters, 1.0, False)
550  clusters_2000_cms = Cluster_Threshold(clusters, 2.0, False)
551  clusters_2500_cms = Cluster_Threshold(clusters, 2.5, False)
552  self.cpair_tcpair_t[0] = Back_to_Back(clusters, clusters)
553  self.cpair_tcpair_t[1] = Back_to_Back(clusters_300_cms, clusters)
554  self.cpair_tcpair_t[2] = Back_to_Back(clusters_400_cms, clusters)
555  self.cpair_tcpair_t[3] = Back_to_Back(clusters_500_cms, clusters)
556  self.cpair_tcpair_t[4] = Back_to_Back(clusters_700_cms, clusters)
557  self.cpair_tcpair_t[5] = Back_to_Back(clusters_1000_cms, clusters)
558  self.cpair_tcpair_t[6] = Back_to_Back(clusters_2000_cms, clusters)
559  self.cpair_tcpair_t[7] = Back_to_Back(clusters_2500_cms, clusters)
560 
561  klmtrkcol = Belle2.PyStoreArray('TRGKLMTracks')
562  self.ntrk_klm_tntrk_klm_t[0] = len(klmtrkcol)
563  klmhitcol = Belle2.PyStoreArray('TRGKLMHits')
564  self.nhit_klm_tnhit_klm_t[0] = len(klmhitcol)
565 
566  matchlist = Belle2.PyStoreArray('TRG3DMatchTracks')
567  self.ntrk_3Dmatch_tntrk_3Dmatch_t[0] = len(matchlist)
568 
569  trginfo = Belle2.PyStoreObj('TRGGRLObjects')
570  self.npair_tc_tnpair_tc_t[0] = trginfo.getNbbTrkCluster()
571  self.npair_cc_tnpair_cc_t[0] = trginfo.getNbbCluster()
572  self.ncluster_1000b_tncluster_1000b_t[0] = trginfo.getNhighcluster2()
573  self.ncluster_2000e_tncluster_2000e_t[0] = trginfo.getNhighcluster4()
574  self.max_deltphi_2dfinder_tmax_deltphi_2dfinder_t[0] = Max_DeltPhi_trk(trk_2d_finder)
575 
576  neutral_clusters = NeutralCluster(clusters)
577  self.ncluster_neutral_tncluster_neutral_t[0] = len(neutral_clusters)
578  max_cluster_neu = Max_Cluster(neutral_clusters, False, 0)
579  max_cluster = Max_Cluster(clusters, False, 0)
580  # max_cms_cluster_neu = Max_Cluster(neutral_clusters, False, 0)
581  # max_cms_cluster = Max_Cluster(clusters, False, 0)
582  smax_cluster_neu = Max_Cluster(neutral_clusters, False, 1)
583  smax_cluster = Max_Cluster(clusters, False, 1)
584  # smax_cms_cluster_neu = Max_Cluster(neutral_clusters, False, 1)
585  # smax_cms_cluster = Max_Cluster(clusters, False, 1)
586  for i in range(self.ncomp_cluncomp_clu):
587  self.max_cluster_neutral_tmax_cluster_neutral_t[i] = max_cluster_neu[i]
588  # self.max_cms_cluster_neutral_t[i] = max_cms_cluster_neu[i]
589  self.max_cluster_tmax_cluster_t[i] = max_cluster[i]
590  # self.max_cms_cluster_t[i] = max_cms_cluster[i]
591  self.smax_cluster_neutral_tsmax_cluster_neutral_t[i] = smax_cluster_neu[i]
592  # self.smax_cms_cluster_neutral_t[i] = smax_cms_cluster_neu[i]
593  self.smax_cluster_tsmax_cluster_t[i] = smax_cluster[i]
594  # self.smax_cms_cluster_t[i] = smax_cms_cluster[i]
595 
596  self.etot_tetot_t[0] = Etot_Cluster(clusters)
597 
598  # bhabha
599  bhabhaveto_1 = BhabhaVeto1(matchlist)
600  bha_logic = 0
601  for bha in bhabhaveto_1:
602  if math.fabs(bha[0]) < 50 and math.fabs(bha[0]) > 10 and math.fabs(bha[1]) < 20:
603  if bha[2] > 2.0 and bha[3] > 2.0 and bha[4] > 6.0:
604  if bha[2] > 3.0 or bha[3] > 3.0:
605  if len(trk_2d_finder) == 2:
606  bha_logic = 1
607  self.bhabha_tbhabha_t[0] = bha_logic
608 
609  if len(bhabhaveto_1) >= 1:
610  for i in range(self.nbha_varnbha_var):
611  self.bhabha_var_tbhabha_var_t[i] = bhabhaveto_1[0][i]
612 
613  # eclbhabha
614  eclbhabhaveto = eclBhabhaVeto(clusters)
615  eclbha_logic = 0
616  for eclbha in eclbhabhaveto:
617  if math.fabs(eclbha[0]) < 50 and math.fabs(eclbha[1]) < 50:
618  if eclbha[2] > 2.0 and eclbha[3] > 2.0 and eclbha[4] > 6.0:
619  if eclbha[2] > 3.0 or eclbha[3] > 3.0:
620  eclbha_logic = 1
621  self.eclbhabha_teclbhabha_t[0] = eclbha_logic
622 
623  # sbhahba
624  sbha_logic = 0
625  if len(trk_2d_finder) == 1:
626  if eclbha_logic == 1:
627  ecol = SBhabhaVeto(matchlist)
628  for i, etr in ecol:
629  if etr[i] > 1.0:
630  sbha_logic = 1
631  self.sbhabha_tsbhabha_t[0] = sbha_logic
632 
633  self.max_deltphi_cluster_tmax_deltphi_cluster_t[0] = Max_DeltPhi_cluster(clusters)
634  if len(eclbhabhaveto) >= 1:
635  for i in range(self.nbha_varnbha_var):
636  self.eclbhabha_var_teclbhabha_var_t[i] = eclbhabhaveto[0][i]
637  self.tgrltgrl.Fill()
638 
639  def terminate(self):
640  """Write and close the file"""
641  self.filefile.cd()
642  self.filefile.Write()
643  self.filefile.Close()
644  PrintBranchDef()
645 
646 
647 if __name__ == "__main__":
648  main = b2.create_path()
649  ma.inputMdst('default', argvs[1], main)
650  main.add_module(CreateLogics())
651  EffCalculation(main)
652  b2.process(main)
The DetectorSet class for sets of detector IDs in the form of EDetector values.
Definition: Const.h:71
a (simplified) python wrapper for StoreArray.
Definition: PyStoreArray.h:56
a (simplified) python wrapper for StoreObjPtr.
Definition: PyStoreObj.h:67
cpair_t
#cluster pairs with different energy threshold
eclbhabha_var_t
variables used in eclbhabha logic, [Delt_theta, Delt_phi, E1, E2, E1+E2]
ncluster_1000b_t
#ecl cluster with threshold >1.0GeV, exclude TC ID 1,2, 17
max_cluster_neutral_t
[energy, theta, phi] of the largest energetic ecl neutral cluster
ntrk_3dfitter_t
#3d fitter tracks
eclbhabha_t
eclbhabha veto logic, 1: bhabha, 0: non bhabha
etot_t
the total deposited cluster energy in ecl
smax_cluster_neutral_t
[energy, theta, phi] of the secondary energetic ecl neutral cluster
int ncomp_clu
#array components
smax_cluster_t
[energy, theta, phi] of the secondary energetic ecl cluster
tgrl
the tree in the output file
ntrk_2dfinder_t
#2d finder tracks
ntrk_3Dmatch_t
#3d matched tracks
time_cluster_t
the timing of clusters
bhabha_t
bhabha veto logic, 1: bhabha, 0: non bhabha
ntrk_2dfitter_t
#2d fitter tracks
max_cluster_t
[energy, theta, phi] of the largest energetic ecl cluster
int nbha_var
#array components
sbhabha_t
bhabha with single track veto logic, 1: bhabha, 0: non bhabha
max_deltphi_cluster_t
max delt phi angle between two clusters
bhabha_var_t
variables used in bhabha logic, [Delt_theta, Delt_phi, E1, E2, E1+E2]
max_deltphi_2dfinder_t
#max phi angle between two 2d finder tracks