Belle II Software  release-05-02-19
PXDPositionEstimation Class Reference
Inheritance diagram for PXDPositionEstimation:
Collaboration diagram for PXDPositionEstimation:

Public Member Functions

def initialize (self)
 
def event (self)
 
def terminate (self)
 

Public Attributes

 nclusters
 Counter for all clusters.
 
 nfound_shapes
 Counter for cluster where shape likelyhood was found in payload.
 
 nfound_offset
 Counter for clusters where position correction was found in payload.
 
 rfile
 Output file to store all plots.
 
 hist_map_momentum
 Histograms for true particle momenta.
 
 hist_map_theta_u
 Histograms for true particle angle thetaU.
 
 hist_map_theta_v
 Histograms for true particle angle thetaV.
 
 hist_map_clustercharge
 Histograms for cluster charge related to particle.
 
 hist_map_residual_u
 Histograms for u residuals.
 
 hist_map_residual_v
 Histograms for v residuals.
 
 hist_map_residual_v_special
 Histograms for v residuals for smaller thetaV range.
 
 hist_map_residual_pull_u
 Histograms for u residual pulls.
 
 hist_map_residual_pull_v
 Histograms for v residual pulls.
 
 binlimits
 ThetaV angle ranges for v residuals.
 

Detailed Description

Histogram the difference position residuals and pulls between clusters and truehits.

Definition at line 21 of file test_cluster_position_estimator.py.

Member Function Documentation

◆ event()

def event (   self)
Fill the residual and pull histograms

Definition at line 97 of file test_cluster_position_estimator.py.

97  def event(self):
98  """Fill the residual and pull histograms"""
99 
100  # Get truehits
101  truehits = Belle2.PyStoreArray("PXDTrueHits")
102 
103  for truehit in truehits:
104  if isinstance(truehit, Belle2.PXDTrueHit):
105  sensor_info = Belle2.VXD.GeoCache.get(truehit.getSensorID())
106  clusters = truehit.getRelationsFrom("PXDClusters")
107 
108  # now check if we find a cluster
109  for j, cls in enumerate(clusters):
110  # we ignore all clusters where less then 100 electrons come from
111  # our truehit
112  if clusters.weight(j) < 100:
113  continue
114 
115  mom = truehit.getMomentum()
116  tu = mom[0] / mom[2]
117  tv = mom[1] / mom[2]
118  thetaU = math.atan(tu) * 180 / math.pi
119  thetaV = math.atan(tv) * 180 / math.pi
120 
121  # Only look at primary particles
122  for mcp in truehit.getRelationsFrom("MCParticles"):
123  if not mcp.hasStatus(1):
124  reject = True
125 
126  # Get instance of position estimator
128  clusterkind = cls.getKind()
129 
130  # Clusterkinds 0,1,2,3 refer to all cases which can currently
131  # be corrected. Cases where a cluster pixel touches a sensor
132  # edge or contains pixel with varying vPitch are excluded here.
133  if clusterkind <= 3 and mom.Mag() > 0.02:
134 
135  self.nclusters += 1
136 
137  # Fill momentum and angles for clusterkind
138  self.hist_map_momentum[clusterkind].Fill(mom.Mag())
139  self.hist_map_theta_u[clusterkind].Fill(thetaU)
140  self.hist_map_theta_v[clusterkind].Fill(thetaV)
141  self.hist_map_clustercharge[clusterkind].Fill(cls.getCharge())
142 
143  # Fill clusterkind=4 for all PXD sensors
144  self.hist_map_momentum[4].Fill(mom.Mag())
145  self.hist_map_theta_u[4].Fill(thetaU)
146  self.hist_map_theta_v[4].Fill(thetaV)
147  self.hist_map_clustercharge[4].Fill(cls.getCharge())
148 
149  # Fill the histograms (mode=2)
150  mode = 2
151  pull_u = (truehit.getU() - cls.getU()) / cls.getUSigma()
152  pull_v = (truehit.getV() - cls.getV()) / cls.getVSigma()
153 
154  self.hist_map_residual_u[(clusterkind, mode)].Fill(truehit.getU() - cls.getU())
155  self.hist_map_residual_v[(clusterkind, mode)].Fill(truehit.getV() - cls.getV())
156  self.hist_map_residual_pull_u[(clusterkind, mode)].Fill(pull_u)
157  self.hist_map_residual_pull_v[(clusterkind, mode)].Fill(pull_v)
158 
159  if thetaV >= self.binlimits[0][0] and thetaV < self.binlimits[0][1]:
160  self.hist_map_residual_v_special[(clusterkind, mode, 0)].Fill(truehit.getV() - cls.getV())
161  elif thetaV >= self.binlimits[1][0] and thetaV < self.binlimits[1][1]:
162  self.hist_map_residual_v_special[(clusterkind, mode, 1)].Fill(truehit.getV() - cls.getV())
163  else:
164  self.hist_map_residual_v_special[(clusterkind, mode, 2)].Fill(truehit.getV() - cls.getV())
165 
166  shape_likelyhood = PositionEstimator.getShapeLikelyhood(cls, tu, tv)
167  if shape_likelyhood > 0:
168  self.nfound_shapes += 1
169 
170  offset = PositionEstimator.getClusterOffset(cls, tu, tv)
171  if offset:
172  # Now, we can safely querry the correction
173  self.nfound_offset += 1
174 
175  # We need to explicitely add a shift to the offsets
176  # This is not needed when working with PXDRecoHits
177  shiftU = sensor_info.getUCellPosition(cls.getUStart())
178  shiftV = sensor_info.getVCellPosition(cls.getVStart())
179 
180  # Fill the histograms (mode=0)
181  mode = 0
182  pull_u = (truehit.getU() - shiftU - offset.getU()) / (math.sqrt(offset.getUSigma2()))
183  pull_v = (truehit.getV() - shiftV - offset.getV()) / (math.sqrt(offset.getVSigma2()))
184 
185  self.hist_map_residual_u[(clusterkind, mode)].Fill(truehit.getU() - shiftU - offset.getU())
186  self.hist_map_residual_v[(clusterkind, mode)].Fill(truehit.getV() - shiftV - offset.getV())
187  self.hist_map_residual_pull_u[(clusterkind, mode)].Fill(pull_u)
188  self.hist_map_residual_pull_v[(clusterkind, mode)].Fill(pull_v)
189 
190  if thetaV >= self.binlimits[0][0] and thetaV < self.binlimits[0][1]:
191  self.hist_map_residual_v_special[(clusterkind, mode, 0)].Fill(
192  truehit.getV() - shiftV - offset.getV())
193  elif thetaV >= self.binlimits[1][0] and thetaV < self.binlimits[1][1]:
194  self.hist_map_residual_v_special[(clusterkind, mode, 1)].Fill(
195  truehit.getV() - shiftV - offset.getV())
196  else:
197  self.hist_map_residual_v_special[(clusterkind, mode, 2)].Fill(
198  truehit.getV() - shiftV - offset.getV())
199 
200  # Fill the histograms (mode=1)
201  mode = 1
202  self.hist_map_residual_u[(clusterkind, mode)].Fill(truehit.getU() - shiftU - offset.getU())
203  self.hist_map_residual_v[(clusterkind, mode)].Fill(truehit.getV() - shiftV - offset.getV())
204  self.hist_map_residual_pull_u[(clusterkind, mode)].Fill(pull_u)
205  self.hist_map_residual_pull_v[(clusterkind, mode)].Fill(pull_v)
206 
207  if thetaV >= self.binlimits[0][0] and thetaV < self.binlimits[0][1]:
208  self.hist_map_residual_v_special[(clusterkind, mode, 0)].Fill(
209  truehit.getV() - shiftV - offset.getV())
210  elif thetaV >= self.binlimits[1][0] and thetaV < self.binlimits[1][1]:
211  self.hist_map_residual_v_special[(clusterkind, mode, 1)].Fill(
212  truehit.getV() - shiftV - offset.getV())
213  else:
214  self.hist_map_residual_v_special[(clusterkind, mode, 2)].Fill(
215  truehit.getV() - shiftV - offset.getV())
216 
217  else:
218 
219  # Fill the histograms (mode=1)
220  mode = 1
221  pull_u = (truehit.getU() - cls.getU()) / cls.getUSigma()
222  pull_v = (truehit.getV() - cls.getV()) / cls.getVSigma()
223 
224  self.hist_map_residual_u[(clusterkind, mode)].Fill(truehit.getU() - cls.getU())
225  self.hist_map_residual_v[(clusterkind, mode)].Fill(truehit.getV() - cls.getV())
226  self.hist_map_residual_pull_u[(clusterkind, mode)].Fill(pull_u)
227  self.hist_map_residual_pull_v[(clusterkind, mode)].Fill(pull_v)
228 
229  if thetaV >= self.binlimits[0][0] and thetaV < self.binlimits[0][1]:
230  self.hist_map_residual_v_special[(clusterkind, mode, 0)].Fill(truehit.getV() - cls.getV())
231  elif thetaV >= self.binlimits[1][0] and thetaV < self.binlimits[1][1]:
232  self.hist_map_residual_v_special[(clusterkind, mode, 1)].Fill(truehit.getV() - cls.getV())
233  else:
234  self.hist_map_residual_v_special[(clusterkind, mode, 2)].Fill(truehit.getV() - cls.getV())
235 

◆ initialize()

def initialize (   self)
Create histograms for pulls and residuals

Definition at line 26 of file test_cluster_position_estimator.py.

◆ terminate()

def terminate (   self)
Format and write all histograms and plot them

Definition at line 236 of file test_cluster_position_estimator.py.


The documentation for this class was generated from the following file:
Belle2::PXDTrueHit
Class PXDTrueHit - Records of tracks that either enter or leave the sensitive volume.
Definition: PXDTrueHit.h:35
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
ClusterEfficiency.ClusterEfficiency.event
def event(self)
Definition: ClusterEfficiency.py:146
Belle2::PXD::PXDClusterPositionEstimator::getInstance
static PXDClusterPositionEstimator & getInstance()
Main (and only) way to access the PXDClusterPositionEstimator.
Definition: PXDClusterPositionEstimator.cc:51
Belle2::PyStoreArray
a (simplified) python wrapper for StoreArray.
Definition: PyStoreArray.h:58