Belle II Software development
PXDRecoHit Class Reference

PXDRecoHit - an extended form of PXDCluster containing geometry information. More...

#include <PXDRecoHit.h>

Inheritance diagram for PXDRecoHit:
AlignablePXDRecoHit

Public Member Functions

 PXDRecoHit ()
 Default constructor for ROOT IO.
 
 PXDRecoHit (const PXDTrueHit *hit, const genfit::TrackCandHit *trackCandHit=NULL, float sigmaU=-1, float sigmaV=-1)
 Construct PXDRecoHit from a PXDTrueHit for Monte Carlo based tracking.
 
 PXDRecoHit (const PXDCluster *hit, float sigmaU, float sigmaV, float covUV)
 Construct PXDRecoHit from a PXD cluster.
 
 PXDRecoHit (const PXDCluster *hit, const genfit::TrackCandHit *trackCandHit=NULL)
 Construct PXDRecoHit from a PXD cluster This constructor is intended as a temporary solution for people who want to test the impact of realistic clusters right now using the current tracking stack and before the final error handling is implemented.
 
genfit::AbsMeasurement * clone () const override
 Creating a deep copy of this hit.
 
virtual std::vector< genfit::MeasurementOnPlane * > constructMeasurementsOnPlane (const genfit::StateOnPlane &state) const override
 Methods that actually interface to Genfit.
 
VxdID getSensorID () const
 Get the compact ID.
 
const PXDTrueHitgetTrueHit () const
 Get pointer to the TrueHit used when creating this RecoHit, can be NULL if created from something else.
 
const PXDClustergetCluster () const
 Get pointer to the Cluster used when creating this RecoHit, can be NULL if created from something else.
 
float getU () const
 Get u coordinate.
 
float getV () const
 Get v coordinate.
 
float getUVariance () const
 Get u coordinate variance.
 
float getVVariance () const
 Get v coordinate variance.
 
float getUVCov () const
 Get u-v error covariance.
 
float getEnergyDep () const
 Get deposited energy.
 
float getShapeLikelyhood (const genfit::StateOnPlane &state) const
 Get deposited energy error.
 
virtual const genfit::AbsHMatrix * constructHMatrix (const genfit::AbsTrackRep *) const override
 Construct the hessian matrix.
 

Private Types

enum  { HIT_DIMENSIONS = 2 }
 

Private Member Functions

void setDetectorPlane ()
 Set up Detector plane information.
 
TVectorD applyPlanarDeformation (TVectorD hitCoords, std::vector< double > planarParameters, const genfit::StateOnPlane &state) const
 Apply planar deformation of sensors.
 

Private Attributes

const PXDTrueHitm_trueHit
 Pointer to the TrueHit used when creating this object.
 
const PXDClusterm_cluster
 transient member (not written out during streaming)
 
float m_energyDep
 transient member (not written out during streaming)
 
unsigned short m_sensorID
 Unique sensor identifier.
 

Detailed Description

PXDRecoHit - an extended form of PXDCluster containing geometry information.

To create a list of PXDRecoHits for all PXDTrueHits belonging to one MCParticle do something like:

//Get the MCParticle in question
MCParticle* mcParticle = ...
//Assume some error on the position
float sigmaU = 10 * Unit::um;
float sigmaV = 15 * Unit::um;
//Iterate over the relation and create a list of hits
for(; it.first!=it.second; ++it.first){
hits.push_back(new PXDRecoHit(it.first->to, sigmaU, sigmaV));
}
A Class to store the Monte Carlo particle information.
Definition MCParticle.h:32
PXDRecoHit()
Default constructor for ROOT IO.
Definition PXDRecoHit.cc:27
Provides access to fast ( O(log n) ) bi-directional lookups on a specified relation.
range_from getElementsFrom(const FROM *from) const
Return a range of all elements pointing from the given object.
boost::iterator_range< iterator_from > range_from
Iterator range [first,second) of the from side.
static const double um
[micrometers]
Definition Unit.h:71
STL class.

Definition at line 53 of file PXDRecoHit.h.

Member Enumeration Documentation

◆ anonymous enum

anonymous enum
private
Enumerator
HIT_DIMENSIONS 

sensitive Dimensions of the Hit

Definition at line 138 of file PXDRecoHit.h.

138{ HIT_DIMENSIONS = 2 };

Constructor & Destructor Documentation

◆ PXDRecoHit() [1/4]

Default constructor for ROOT IO.

Definition at line 27 of file PXDRecoHit.cc.

27 :
28 genfit::PlanarMeasurement(HIT_DIMENSIONS), m_trueHit(0), m_cluster(0),
29 m_energyDep(0), m_sensorID(0)//, m_energyDepError(0)
30{}
float m_energyDep
transient member (not written out during streaming)
Definition PXDRecoHit.h:144
unsigned short m_sensorID
Unique sensor identifier.
Definition PXDRecoHit.h:147
const PXDTrueHit * m_trueHit
Pointer to the TrueHit used when creating this object.
Definition PXDRecoHit.h:141
const PXDCluster * m_cluster
transient member (not written out during streaming)
Definition PXDRecoHit.h:143
@ HIT_DIMENSIONS
sensitive Dimensions of the Hit
Definition PXDRecoHit.h:138

◆ PXDRecoHit() [2/4]

PXDRecoHit ( const PXDTrueHit * hit,
const genfit::TrackCandHit * trackCandHit = NULL,
float sigmaU = -1,
float sigmaV = -1 )
explicit

Construct PXDRecoHit from a PXDTrueHit for Monte Carlo based tracking.

This requires a valid random number generator to be initialized at gRandom. The Hit position will be smeared using a gaussian smearing with sigmaU and sigmaV along u and v respectively

If one of the errors is set <0, a default resolution will be assumed for both values by dividing the pixel size by sqrt(12).

Parameters
hitPXDTrueHit to use as base
trackCandHitactually not used, should be removed?
sigmaUError of the Hit along u
sigmaVError of the Hit along v

Definition at line 32 of file PXDRecoHit.cc.

32 :
33 genfit::PlanarMeasurement(HIT_DIMENSIONS), m_trueHit(hit), m_cluster(0),
34 m_energyDep(0), m_sensorID(0)//, m_energyDepError(0)
35{
36 if (!gRandom) B2FATAL("gRandom not initialized, please set up gRandom first");
37
38 // Set the sensor UID
39 m_sensorID = hit->getSensorID();
40
41 //If no error is given, estimate the error by dividing the pixel size by sqrt(12)
42 if (sigmaU < 0 || sigmaV < 0) {
43 const PXD::SensorInfo& geometry = dynamic_cast<const PXD::SensorInfo&>(VXD::GeoCache::getInstance().getSensorInfo(m_sensorID));
44 sigmaU = geometry.getUPitch(hit->getV()) / sqrt(12);
45 sigmaV = geometry.getVPitch(hit->getV()) / sqrt(12);
46 }
47
48 // Set positions
49 rawHitCoords_(0) = gRandom->Gaus(hit->getU(), sigmaU);
50 rawHitCoords_(1) = gRandom->Gaus(hit->getV(), sigmaV);
51 // Set the error covariance matrix
52 rawHitCov_(0, 0) = sigmaU * sigmaU;
53 rawHitCov_(0, 1) = 0;
54 rawHitCov_(1, 0) = 0;
55 rawHitCov_(1, 1) = sigmaV * sigmaV;
56 // Set physical parameters
57 m_energyDep = hit->getEnergyDep();
58 // Setup geometry information
60}
void setDetectorPlane()
Set up Detector plane information.
const SensorInfoBase & getSensorInfo(Belle2::VxdID id) const
Return a reference to the SensorInfo of a given SensorID.
Definition GeoCache.cc:67
static GeoCache & getInstance()
Return a reference to the singleton instance.
Definition GeoCache.cc:214
double sqrt(double a)
sqrt for double
Definition beamHelpers.h:28

◆ PXDRecoHit() [3/4]

PXDRecoHit ( const PXDCluster * hit,
float sigmaU,
float sigmaV,
float covUV )
explicit

Construct PXDRecoHit from a PXD cluster.

For users that want to supply their own errors on construction

Parameters
hitPXDCluster to use as base
sigmaUError of the Hit along u
sigmaVError of the Hit along v
covUVCovariance between u and v

Definition at line 62 of file PXDRecoHit.cc.

62 :
63 genfit::PlanarMeasurement(HIT_DIMENSIONS), m_trueHit(0), m_cluster(hit),
64 m_energyDep(0), m_sensorID(0)//, m_energyDepError(0)
65{
66 // Set the sensor UID
67 m_sensorID = hit->getSensorID();
68 // Set positions
69 rawHitCoords_(0) = hit->getU();
70 rawHitCoords_(1) = hit->getV();
71 // Set the error covariance matrix
72 rawHitCov_(0, 0) = sigmaU * sigmaU;
73 rawHitCov_(0, 1) = covUV;
74 rawHitCov_(1, 0) = covUV;
75 rawHitCov_(1, 1) = sigmaV * sigmaV;
76 // Set physical parameters
77 const PXD::SensorInfo& SensorInfo = dynamic_cast<const PXD::SensorInfo&>(VXD::GeoCache::getInstance().getSensorInfo(m_sensorID));
78 auto ADUToEnergy = PXD::PXDGainCalibrator::getInstance().getADUToEnergy(m_sensorID, SensorInfo.getUCellID(hit->getU()),
79 SensorInfo.getVCellID(hit->getV()));
80 m_energyDep = hit->getCharge() * ADUToEnergy;
81 //m_energyDepError = 0;
82 // Setup geometry information
84}
float getADUToEnergy(VxdID id, unsigned int uid, unsigned int vid) const
Get conversion factor from ADU to energy.
static PXDGainCalibrator & getInstance()
Main (and only) way to access the PXDGainCalibrator.
int getVCellID(double v, bool clamp=false) const
Return the corresponding pixel/strip ID of a given v coordinate.
int getUCellID(double u, double v=0, bool clamp=false) const
Return the corresponding pixel/strip ID of a given u coordinate.

◆ PXDRecoHit() [4/4]

PXDRecoHit ( const PXDCluster * hit,
const genfit::TrackCandHit * trackCandHit = NULL )
explicit

Construct PXDRecoHit from a PXD cluster This constructor is intended as a temporary solution for people who want to test the impact of realistic clusters right now using the current tracking stack and before the final error handling is implemented.

The pitch / sqrt(12) will be added as measurement error estimation. This is of course not the exact error of the cluster and one cannot expect perfect tracking results when using this constructor

Parameters
hitPXDCluster to use as base
trackCandHitactually not used, should be removed?

Definition at line 87 of file PXDRecoHit.cc.

87 :
88 genfit::PlanarMeasurement(HIT_DIMENSIONS), m_trueHit(0), m_cluster(hit),
89 m_energyDep(0), m_sensorID(0)//, m_energyDepError(0)
90{
91 // Set the sensor UID
92 m_sensorID = hit->getSensorID();
93 // Set positions
94 rawHitCoords_(0) = hit->getU();
95 rawHitCoords_(1) = hit->getV();
96 // Set the error covariance matrix
97 rawHitCov_(0, 0) = hit->getUSigma() * hit->getUSigma();
98 rawHitCov_(0, 1) = hit->getRho() * hit->getUSigma() * hit->getVSigma();
99 rawHitCov_(1, 0) = hit->getRho() * hit->getUSigma() * hit->getVSigma();
100 rawHitCov_(1, 1) = hit->getVSigma() * hit->getVSigma();
101 // Set physical parameters
102 const PXD::SensorInfo& SensorInfo = dynamic_cast<const PXD::SensorInfo&>(VXD::GeoCache::getInstance().getSensorInfo(m_sensorID));
103 auto ADUToEnergy = PXD::PXDGainCalibrator::getInstance().getADUToEnergy(m_sensorID, SensorInfo.getUCellID(hit->getU()),
104 SensorInfo.getVCellID(hit->getV()));
105 m_energyDep = hit->getCharge() * ADUToEnergy;
106 //m_energyDepError = 0;
107 // Setup geometry information
109}

Member Function Documentation

◆ applyPlanarDeformation()

TVectorD applyPlanarDeformation ( TVectorD hitCoords,
std::vector< double > planarParameters,
const genfit::StateOnPlane & state ) const
private

Apply planar deformation of sensors.

Definition at line 149 of file PXDRecoHit.cc.

151{
152 // Legendre parametrization of deformation
153 auto L1 = [](double x) {return x;};
154 auto L2 = [](double x) {return (3 * x * x - 1) / 2;};
155 auto L3 = [](double x) {return (5 * x * x * x - 3 * x) / 2;};
156 auto L4 = [](double x) {return (35 * x * x * x * x - 30 * x * x + 3) / 8;};
157
158 const PXD::SensorInfo& geometry = dynamic_cast<const PXD::SensorInfo&>(VXD::GeoCache::getInstance().getSensorInfo(m_sensorID));
159
160 double u = hitCoords[0];
161 double v = hitCoords[1];
162 double width = geometry.getWidth(v); // Width of sensor (U side)
163 double length = geometry.getLength(); // Length of sensor (V side)
164 u = u * 2 / width; // Legendre parametrization required U in (-1, 1)
165 v = v * 2 / length; // Legendre parametrization required V in (-1, 1)
166
167 /* Planar deformation using Legendre parametrization
168 w(u, v) = L_{31} * L2(u) + L_{32} * L1(u) * L1(v) + L_{33} * L2(v) +
169 L_{41} * L3(u) + L_{42} * L2(u) * L1(v) + L_{43} * L1(u) * L2(v) + L_{44} * L3(v) +
170 L_{51} * L4(u) + L_{52} * L3(u) * L1(v) + L_{53} * L2(u) * L2(v) + L_{54} * L1(u) * L3(v) + L_{55} * L4(v); */
171 double dw =
172 planarParameters[0] * L2(u) + planarParameters[1] * L1(u) * L1(v) + planarParameters[2] * L2(v) +
173 planarParameters[3] * L3(u) + planarParameters[4] * L2(u) * L1(v) + planarParameters[5] * L1(u) * L2(v) + planarParameters[6] * L3(
174 v) +
175 planarParameters[7] * L4(u) + planarParameters[8] * L3(u) * L1(v) + planarParameters[9] * L2(u) * L2(v) + planarParameters[10] * L1(
176 u) * L3(v) + planarParameters[11] * L4(v);
177
178 double du_dw = state.getState()[1]; // slope in U direction
179 double dv_dw = state.getState()[2]; // slope in V direction
180
181 u = u * width / 2; // from Legendre to Local parametrization
182 v = v * length / 2; // from Legendre to Local parametrization
183
184 TVectorD pos(2);
185
186 pos[0] = u + dw * du_dw;
187 pos[1] = v + dw * dv_dw;
188
189 return pos;
190}

◆ clone()

genfit::AbsMeasurement * clone ( ) const
override

Creating a deep copy of this hit.

Definition at line 111 of file PXDRecoHit.cc.

112{
113 return new PXDRecoHit(*this);
114}

◆ constructHMatrix()

virtual const genfit::AbsHMatrix * constructHMatrix ( const genfit::AbsTrackRep * ) const
inlineoverridevirtual

Construct the hessian matrix.

Definition at line 134 of file PXDRecoHit.h.

134{ return new genfit::HMatrixUV(); };

◆ constructMeasurementsOnPlane()

std::vector< genfit::MeasurementOnPlane * > constructMeasurementsOnPlane ( const genfit::StateOnPlane & state) const
overridevirtual

Methods that actually interface to Genfit.

Definition at line 192 of file PXDRecoHit.cc.

193{
194 // Track-based update only takes place when the RecoHit has an associated cluster
195 if (this->getCluster()) {
196 // Check if we can correct position coordinates based on track info
197 const TVectorD& state5 = state.getState();
198 auto offset = PXD::PXDClusterPositionEstimator::getInstance().getClusterOffset(*this->getCluster(), state5[1], state5[2]);
199
200 if (offset != nullptr) {
201 // Found a valid offset, lets apply it
202 const Belle2::VxdID& sensorID = (*this->getCluster()).getSensorID();
203 const Belle2::PXD::SensorInfo& Info = dynamic_cast<const Belle2::PXD::SensorInfo&>(VXD::GeoCache::getInstance().getSensorInfo(
204 sensorID));
205 double posU = Info.getUCellPosition((*this->getCluster()).getUStart());
206 double posV = Info.getVCellPosition((*this->getCluster()).getVStart());
207
208 TVectorD hitCoords(2);
209 hitCoords(0) = posU + offset->getU();
210 hitCoords(1) = posV + offset->getV();
211 TMatrixDSym hitCov(2);
212 hitCov(0, 0) = offset->getUSigma2();
213 hitCov(0, 1) = offset->getUVCovariance();
214 hitCov(1, 0) = offset->getUVCovariance();
215 hitCov(1, 1) = offset->getVSigma2();
216
217 // Apply planar deformation
218 TVectorD pos = applyPlanarDeformation(hitCoords, VXD::GeoCache::getInstance().getSensorInfo(m_sensorID).getSurfaceParameters(),
219 state);
220
221 return std::vector<genfit::MeasurementOnPlane*>(1, new genfit::MeasurementOnPlane(
222 pos, hitCov, state.getPlane(), state.getRep(), this->constructHMatrix(state.getRep())
223 ));
224 }
225 }
226
227 // Apply planar deformation
228 TVectorD pos = applyPlanarDeformation(rawHitCoords_, VXD::GeoCache::getInstance().getSensorInfo(m_sensorID).getSurfaceParameters(),
229 state);
230
231 // If we reach here, we can do no better than what we have
232 return std::vector<genfit::MeasurementOnPlane*>(1, new genfit::MeasurementOnPlane(
233 pos, rawHitCov_, state.getPlane(), state.getRep(), this->constructHMatrix(state.getRep())
234 ));
235}
TVectorD applyPlanarDeformation(TVectorD hitCoords, std::vector< double > planarParameters, const genfit::StateOnPlane &state) const
Apply planar deformation of sensors.
const PXDCluster * getCluster() const
Get pointer to the Cluster used when creating this RecoHit, can be NULL if created from something els...
Definition PXDRecoHit.h:110
VxdID getSensorID() const
Get the compact ID.
Definition PXDRecoHit.h:105
static PXDClusterPositionEstimator & getInstance()
Main (and only) way to access the PXDClusterPositionEstimator.
const PXDClusterOffsetPar * getClusterOffset(const PXDCluster &cluster, double tu, double tv) const
Return pointer to cluster offsets, can be nullptr.
double getVCellPosition(int vID) const
Return the position of a specific strip/pixel in v direction.
double getUCellPosition(int uID, int vID=-1) const
Return the position of a specific strip/pixel in u direction.

◆ getCluster()

const PXDCluster * getCluster ( ) const
inline

Get pointer to the Cluster used when creating this RecoHit, can be NULL if created from something else.

Definition at line 110 of file PXDRecoHit.h.

110{ return m_cluster; }

◆ getEnergyDep()

float getEnergyDep ( ) const
inline

Get deposited energy.

Definition at line 125 of file PXDRecoHit.h.

125{ return m_energyDep; }

◆ getSensorID()

VxdID getSensorID ( ) const
inline

Get the compact ID.

Definition at line 105 of file PXDRecoHit.h.

105{ return m_sensorID; }

◆ getShapeLikelyhood()

float getShapeLikelyhood ( const genfit::StateOnPlane & state) const

Get deposited energy error.

Get the likelihood that cluster shape is likely to be created from track state.

Definition at line 137 of file PXDRecoHit.cc.

138{
139 // We need an associated cluster
140 if (this->getCluster()) {
141 // Likelihood depends on the fitted incidence angles into the sensor
142 const TVectorD& state5 = state.getState();
143 return PXD::PXDClusterPositionEstimator::getInstance().getShapeLikelyhood(*this->getCluster(), state5[1], state5[2]);
144 }
145 // If we reach here, we can do no better than return zero
146 return 0;
147}
float getShapeLikelyhood(const PXDCluster &cluster, double tu, double tv) const
Return cluster shape likelihood.

◆ getTrueHit()

const PXDTrueHit * getTrueHit ( ) const
inline

Get pointer to the TrueHit used when creating this RecoHit, can be NULL if created from something else.

Definition at line 108 of file PXDRecoHit.h.

108{ return m_trueHit; }

◆ getU()

float getU ( ) const
inline

Get u coordinate.

Definition at line 113 of file PXDRecoHit.h.

113{ return rawHitCoords_(0); }

◆ getUVariance()

float getUVariance ( ) const
inline

Get u coordinate variance.

Definition at line 118 of file PXDRecoHit.h.

118{ return rawHitCov_(0, 0); }

◆ getUVCov()

float getUVCov ( ) const
inline

Get u-v error covariance.

Definition at line 122 of file PXDRecoHit.h.

122{ return rawHitCov_(0, 1); }

◆ getV()

float getV ( ) const
inline

Get v coordinate.

Definition at line 115 of file PXDRecoHit.h.

115{ return rawHitCoords_(1); }

◆ getVVariance()

float getVVariance ( ) const
inline

Get v coordinate variance.

Definition at line 120 of file PXDRecoHit.h.

120{ return rawHitCov_(1, 1); }

◆ setDetectorPlane()

void setDetectorPlane ( )
private

Set up Detector plane information.

Definition at line 117 of file PXDRecoHit.cc.

118{
119 // Construct a finite detector plane and set it.
120 const PXD::SensorInfo& geometry = dynamic_cast<const PXD::SensorInfo&>(VXD::GeoCache::getInstance().getSensorInfo(m_sensorID));
121
122 // Construct vectors o, u, v
123 ROOT::Math::XYZVector uLocal(1, 0, 0);
124 ROOT::Math::XYZVector vLocal(0, 1, 0);
125 ROOT::Math::XYZVector origin = geometry.pointToGlobal(ROOT::Math::XYZVector(0, 0, 0), true);
126 ROOT::Math::XYZVector uGlobal = geometry.vectorToGlobal(uLocal, true);
127 ROOT::Math::XYZVector vGlobal = geometry.vectorToGlobal(vLocal, true);
128
129 //Construct the detector plane
130 VXD::SensorPlane* finitePlane = new VXD::SensorPlane(m_sensorID, 20.0, 20.0);
131 genfit::SharedPlanePtr detPlane(new genfit::DetPlane(XYZToTVector(origin), XYZToTVector(uGlobal), XYZToTVector(vGlobal),
132 finitePlane));
133 setPlane(detPlane, m_sensorID);
134}
static constexpr auto XYZToTVector
Helper function to convert XYZVector to TVector3.
Definition VectorUtil.h:24

Member Data Documentation

◆ m_cluster

const PXDCluster* m_cluster
private

transient member (not written out during streaming)

Pointer to the Cluster used when creating this object

Definition at line 143 of file PXDRecoHit.h.

◆ m_energyDep

float m_energyDep
private

transient member (not written out during streaming)

deposited energy.

Definition at line 144 of file PXDRecoHit.h.

◆ m_sensorID

unsigned short m_sensorID
private

Unique sensor identifier.

Definition at line 147 of file PXDRecoHit.h.

◆ m_trueHit

const PXDTrueHit* m_trueHit
private

Pointer to the TrueHit used when creating this object.

Definition at line 141 of file PXDRecoHit.h.


The documentation for this class was generated from the following files: