Belle II Software  release-05-02-19
SVDValidationTTreeCluster.py
1 #!/usr/bin/env python3
2 # -*- coding: utf-8 -*-
3 
4 """
5 <header>
6  <contact> SVD Software Group, svd-software@belle2.org </contact>
7  <description>
8  This module is used for the SVD validation.
9  It gets information about clusters and truehits, saving in a ttree in a ROOT file.
10  </description>
11 </header>
12 """
13 
14 import math
15 
16 from basf2 import *
17 
18 # Some ROOT tools
19 import ROOT
20 from ROOT import Belle2
21 from ROOT import gROOT, AddressOf
22 from ROOT import TVector3
23 
24 # Define a ROOT struct to hold output data in the TTree
25 gROOT.ProcessLine('struct EventDataCluster {\
26  int sensor_id;\
27  int layer;\
28  int ladder;\
29  int sensor;\
30  int sensor_type;\
31  int strip_dir;\
32  int cluster_truehits_number;\
33  float strip_pitch;\
34  float cluster_theta;\
35  float cluster_phi;\
36  float cluster_position;\
37  float cluster_positionSigma;\
38  float cluster_clsTime;\
39  float cluster_clsTimeSigma;\
40  float cluster_charge;\
41  float cluster_seedCharge;\
42  float cluster_size;\
43  float cluster_snr;\
44  float cluster_interstripPosition;\
45  float cluster_pull;\
46  float cluster_residual;\
47  float truehit_position;\
48  float truehit_interstripPosition;\
49  float truehit_deposEnergy;\
50  float truehit_lossmomentum;\
51  float truehit_time;\
52 };')
53 
54 
55 from ROOT import EventDataCluster
56 
57 
58 class SVDValidationTTreeCluster(Module):
59  '''class to produce the ttree for cluster validation'''
60 
61  def __init__(self):
62  """Initialize the module"""
63 
64  super(SVDValidationTTreeCluster, self).__init__()
65  self.file = ROOT.TFile('../SVDValidationTTreeCluster.root', 'recreate')
66  '''Output ROOT file'''
67  self.tree = ROOT.TTree('tree', 'Event data of SVD validation events')
68  '''TTrees for output data'''
69  self.data = EventDataCluster()
70  '''Instance of the EventData class'''
71 
72  # Declare tree branches
73  for key in EventDataCluster.__dict__:
74  if '__' not in key:
75  formstring = '/F'
76  if isinstance(self.data.__getattribute__(key), int):
77  formstring = '/I'
78  self.tree.Branch(key, AddressOf(self.data, key), key + formstring)
79 
80  def event(self):
81  """ Find clusters with a truehit and save needed information """
82  clusters = Belle2.PyStoreArray('SVDClusters')
83  for cluster in clusters:
84  cluster_truehits = cluster.getRelationsTo('SVDTrueHits') # SVDClustersToSVDTrueHits
85  cluster_TrueHit_Length = len(cluster_truehits)
86  if (cluster_TrueHit_Length == 0) or (cluster_TrueHit_Length != 1):
87  # Sensor identification
88  sensorID = cluster.getSensorID()
89  self.data.sensor_id = int(sensorID)
90  sensorNum = sensorID.getSensorNumber()
91  self.data.sensor = sensorNum
92  layerNum = sensorID.getLayerNumber()
93  self.data.layer = layerNum
94  if (layerNum == 3):
95  sensorType = 1
96  else:
97  if (sensorNum == 1):
98  sensorType = 0
99  else:
100  sensorType = 1
101  self.data.sensor_type = sensorType
102  ladderNum = sensorID.getLadderNumber()
103  self.data.ladder = ladderNum
104  #
105  self.data.strip_dir = -1
106  self.data.cluster_truehits_number = cluster_TrueHit_Length
107  # Fill tree
108  self.file.cd()
109  self.tree.Fill()
110  else:
111  # We want only clusters with exactly one associated TrueHit
112  for truehit in cluster_truehits:
113  self.data.cluster_truehits_number = cluster_TrueHit_Length
114  sensorInfo = Belle2.VXD.GeoCache.get(cluster.getSensorID())
115  # Sensor identification
116  sensorID = cluster.getSensorID()
117  self.data.sensor_id = int(sensorID)
118  sensorNum = sensorID.getSensorNumber()
119  self.data.sensor = sensorNum
120  layerNum = sensorID.getLayerNumber()
121  self.data.layer = layerNum
122  if (layerNum == 3):
123  sensorType = 1
124  else:
125  if (sensorNum == 1):
126  sensorType = 0
127  else:
128  sensorType = 1
129  self.data.sensor_type = sensorType
130  ladderNum = sensorID.getLadderNumber()
131  self.data.ladder = ladderNum
132  # Cluster information
133  self.data.cluster_clsTime = cluster.getClsTime()
134  self.data.cluster_clsTimeSigma = cluster.getClsTimeSigma()
135  self.data.cluster_charge = cluster.getCharge()
136  self.data.cluster_seedCharge = cluster.getSeedCharge()
137  self.data.cluster_size = cluster.getSize()
138  self.data.cluster_snr = cluster.getSNR()
139  cluster_position = cluster.getPosition()
140  if cluster.isUCluster():
141  cluster_position = cluster.getPosition(truehit.getV())
142  # Interstrip position calculations
143  if cluster.isUCluster():
144  strip_dir = 0
145  strip_pitch = sensorInfo.getUPitch(truehit.getV())
146  else:
147  strip_dir = 1
148  strip_pitch = sensorInfo.getVPitch(truehit.getU())
149  self.data.strip_dir = strip_dir
150  self.data.strip_pitch = strip_pitch
151  cluster_interstripPosition = cluster_position % strip_pitch / strip_pitch
152  self.data.cluster_interstripPosition = cluster_interstripPosition
153  # theta and phi definitions
154  if cluster.isUCluster():
155  uPos = cluster_position
156  vPos = 0
157  else:
158  uPos = 0
159  vPos = cluster_position
160  localPosition = TVector3(uPos, vPos, 0) # sensor center at (0, 0, 0)
161  globalPosition = sensorInfo.pointToGlobal(localPosition, True)
162  x = globalPosition[0]
163  y = globalPosition[1]
164  z = globalPosition[2]
165  # see https://d2comp.kek.jp/record/242?ln=en for the Belle II
166  # coordinate system and related variables
167  rho = math.sqrt(x * x + y * y)
168  r = math.sqrt(x * x + y * y + z * z)
169  # get theta as arccosine(z/r)
170  thetaRadians = math.acos(z / r)
171  theta = (thetaRadians * 180) / math.pi
172  # get phi as arccosine(x/rho)
173  phiRadians = math.acos(x / rho)
174  if y < 0:
175  phi = 360 - (phiRadians * 180) / math.pi
176  else:
177  phi = (phiRadians * 180) / math.pi
178  self.data.cluster_theta = theta
179  self.data.cluster_phi = phi
180  # Pull calculations
181  clusterPos = cluster_position
182  clusterPosSigma = cluster.getPositionSigma()
183  if cluster.isUCluster():
184  truehitPos = truehit.getU()
185  else:
186  truehitPos = truehit.getV()
187  cluster_residual = clusterPos - truehitPos
188  cluster_pull = cluster_residual / clusterPosSigma
189  self.data.cluster_position = clusterPos
190  self.data.cluster_positionSigma = clusterPosSigma
191  self.data.cluster_residual = cluster_residual
192  self.data.cluster_pull = cluster_pull
193  # Truehit information
194  self.data.truehit_position = truehitPos
195  truehit_interstripPosition = truehitPos % strip_pitch / strip_pitch
196  self.data.truehit_interstripPosition = truehit_interstripPosition
197  self.data.truehit_deposEnergy = truehit.getEnergyDep()
198  self.data.truehit_lossmomentum = truehit.getEntryMomentum().Mag() - truehit.getExitMomentum().Mag()
199  self.data.truehit_time = truehit.getGlobalTime()
200  # Fill tree
201  self.file.cd()
202  self.tree.Fill()
203 
204  def terminate(self):
205  """Close the output file. """
206  self.file.cd()
207  self.file.Write()
208  self.file.Close()
SVDValidationTTreeCluster.SVDValidationTTreeCluster.data
data
Definition: SVDValidationTTreeCluster.py:69
Belle2::VXD::GeoCache::get
static const SensorInfoBase & get(Belle2::VxdID id)
Return a reference to the SensorInfo of a given SensorID.
Definition: GeoCache.h:141
SVDValidationTTreeCluster.SVDValidationTTreeCluster
Definition: SVDValidationTTreeCluster.py:58
SVDValidationTTreeCluster.SVDValidationTTreeCluster.event
def event(self)
Definition: SVDValidationTTreeCluster.py:80
SVDValidationTTreeCluster.SVDValidationTTreeCluster.__init__
def __init__(self)
Definition: SVDValidationTTreeCluster.py:61
SVDValidationTTreeCluster.SVDValidationTTreeCluster.tree
tree
Definition: SVDValidationTTreeCluster.py:67
SVDValidationTTreeCluster.SVDValidationTTreeCluster.file
file
Definition: SVDValidationTTreeCluster.py:65
SVDValidationTTreeCluster.SVDValidationTTreeCluster.terminate
def terminate(self)
Definition: SVDValidationTTreeCluster.py:204
Belle2::PyStoreArray
a (simplified) python wrapper for StoreArray.
Definition: PyStoreArray.h:58