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