ATLAS Offline Software
Loading...
Searching...
No Matches
AFPSiDBasicKalmanToolTrack Class Reference

Class representing a reconstructed track with Kalman filter. More...

#include <AFPSiDBasicKalmanToolTrack.h>

Collaboration diagram for AFPSiDBasicKalmanToolTrack:

Public Member Functions

 AFPSiDBasicKalmanToolTrack (const xAOD::AFPSiHitsCluster *firstCluster, const xAOD::AFPSiHitsCluster *secondCluster, const CLHEP::HepMatrix &observationModel, const CLHEP::HepMatrix &observationNoise, const CLHEP::HepMatrix &covarianceMatrix)
 Main constructor that creates the track seed from two clusters.
CLHEP::HepMatrix transitionModel (const double firstPos, const double secondPos) const
 Create state transition matrix to calculate x and y position and slopes between to Z positions.
CLHEP::HepMatrix transitionModel (const double zPosition) const
 Create state transition matrix to calculate x and y position and slopes at given Z.
CLHEP::HepVector predictNextPoint (const double zPosition) const
 Predict x and y positions and slopes of the track at given Z.
CLHEP::HepMatrix predictNextPointCov (const double zPosition) const
 Predict covariance matrix of the track at given Z.
const xAOD::AFPSiHitsClusterfindNearestCluster (const std::vector< const xAOD::AFPSiHitsCluster * > &clusters, const float maxAllowedDistance) const
 Finds the cluster which is nearest to the track, returns nullptr if no is found.
void addCluster (const xAOD::AFPSiHitsCluster *cluster, const float htiMaxChi2)
 Adds a new cluster to the track and updates history and position.
int holes () const
void addHole ()
 Increase number of holes by 1. (see m_holes)
const std::list< const xAOD::AFPSiHitsCluster * > & clustersInTrack () const
 Clusters used to reconstruct the track.
void clearSmoothMatrices ()
 Clears matrices used for smoothing the track.
const std::list< CLHEP::HepVector > & positionAndSlopeSmooth () const
 Vectors of reconstructed true positions for each step after smoothing \(x_{k}\).
const std::list< CLHEP::HepVector > & positionAndSlopeHistory () const
 Vectors of reconstructed true positions for each step \(x_{k}\).
const std::list< double > & chi2Smooth () const
 History of chi2 for each step after smoothing.
double trkChi2NDFSmooth () const
 Chi2 per degrees of freedom (clusters) for the whole smoothed track.
void smooth ()
 Run smoothing algorithm.

Protected Attributes

std::list< float > m_zPositionHistory
 Vectors of reconstructed true positions for each step.
std::list< CLHEP::HepVector > m_positionAndSlopeHistory
 Vectors of reconstructed true positions for each step \(x_{k}\).
std::list< CLHEP::HepVector > m_positionAndSlopePredictedHistory
 Vectors of predicted true position of a track for each step.
std::list< CLHEP::HepMatrix > m_positionAndSlopeCovHistory
 Matrices of updated (a posteriori) estimate covariance of the position and slope for each step \(P_{k|k}\).
std::list< CLHEP::HepMatrix > m_positionAndSlopeCovPredictedHistory
 Matrices of predicted estimate covariance of the position and slope for each step.
std::list< CLHEP::HepMatrix > m_transitionHistory
 Matrices of transitionModel used at each step of track reconstruction.
std::list< double > m_chi2History
 History of chi2 for each step.
std::list< CLHEP::HepVector > m_positionAndSlopeSmooth
 Vectors of reconstructed true positions for each step after smoothing \(x_{k}\).
std::list< CLHEP::HepMatrix > m_positionAndSlopeSmoothCov
 Matrices of updated (a posteriori) estimate covariance of the position and slope for each step after smoothing \(P_{k|k}\).
std::list< double > m_chi2Smooth
 History of chi2 for each step after smoothing.
std::list< const xAOD::AFPSiHitsCluster * > m_clustersInTrack
 Clusters used to reconstruct the track.
const CLHEP::HepMatrix & m_observationModel
 Observation model, describing transformation from true to measured state.
const CLHEP::HepMatrix & m_observationNoise
 Matrix representing observation noise disturbing the measurement.
int m_holes
 Number of layers in which hit from a track is expected, but is not present.

Detailed Description

Class representing a reconstructed track with Kalman filter.

Definition at line 23 of file AFPSiDBasicKalmanToolTrack.h.

Constructor & Destructor Documentation

◆ AFPSiDBasicKalmanToolTrack()

AFPSiDBasicKalmanToolTrack::AFPSiDBasicKalmanToolTrack ( const xAOD::AFPSiHitsCluster * firstCluster,
const xAOD::AFPSiHitsCluster * secondCluster,
const CLHEP::HepMatrix & observationModel,
const CLHEP::HepMatrix & observationNoise,
const CLHEP::HepMatrix & covarianceMatrix )

Main constructor that creates the track seed from two clusters.

It adds the two clusters to the list of clusters used to reconstruct the track. It sets the position of the last point of the track to be equal to the position of the second cluster and calculates track slopes using first and second clusters. It adds first entry to m_chi2History and sets it to 0.

Parameters
firstClusterpointer to the first cluster of the track
secondClusterpointer to the second cluster of the track
observationModelmatrix representing observation model (see m_observationModel)
observationNoisesee m_observationNoise
covarianceMatrixinitial error covariance matrix (see m_positionAndSlopeCovHistory)

Definition at line 19 of file AFPSiDBasicKalmanToolTrack.cxx.

19 :
20 m_observationModel(observationModel),
21 m_observationNoise(observationNoise),
22 m_holes(0)
23{
24 m_positionAndSlopeHistory.emplace_back(4);
25
26 const double zDistance = secondCluster->zLocal() - firstCluster->zLocal(); // absolute value taken in case global z position is used, in that case the distance would be negative for stations on the A side
27
28 CLHEP::HepVector& positionAndSlope = m_positionAndSlopeHistory.back();
29 positionAndSlope[0] = secondCluster->xLocal(); // set track X position
30 positionAndSlope[1] = (secondCluster->xLocal() - firstCluster->xLocal())/zDistance; // set track X slope
31 positionAndSlope[2] = secondCluster->yLocal(); // set track Y position
32 positionAndSlope[3] = (secondCluster->yLocal() - firstCluster->yLocal())/zDistance; // set track Y slope
33
34 // set z position of the track
35 m_zPositionHistory.push_back(secondCluster->zLocal());
36
37 // add clusters to the list
38 m_clustersInTrack.push_back(firstCluster);
39 m_clustersInTrack.push_back(secondCluster);
40
41 // set the first covariance matrix
42 m_positionAndSlopeCovHistory.emplace_back(covarianceMatrix);
43
44 // set chi2 to 0 for two points
45 m_chi2History.emplace_back(0);
46}
std::list< const xAOD::AFPSiHitsCluster * > m_clustersInTrack
Clusters used to reconstruct the track.
std::list< CLHEP::HepMatrix > m_positionAndSlopeCovHistory
Matrices of updated (a posteriori) estimate covariance of the position and slope for each step .
int m_holes
Number of layers in which hit from a track is expected, but is not present.
const CLHEP::HepMatrix & m_observationNoise
Matrix representing observation noise disturbing the measurement.
std::list< CLHEP::HepVector > m_positionAndSlopeHistory
Vectors of reconstructed true positions for each step .
const CLHEP::HepMatrix & m_observationModel
Observation model, describing transformation from true to measured state.
std::list< float > m_zPositionHistory
Vectors of reconstructed true positions for each step.
std::list< double > m_chi2History
History of chi2 for each step.
float zLocal() const
Cluster position along Z axis in station local coordinate system.
float xLocal() const
Cluster position along X axis in station local coordinate system.
float yLocal() const
Cluster position along Y axis in station local coordinate system.

Member Function Documentation

◆ addCluster()

void AFPSiDBasicKalmanToolTrack::addCluster ( const xAOD::AFPSiHitsCluster * cluster,
const float htiMaxChi2 )

Adds a new cluster to the track and updates history and position.

Definition at line 48 of file AFPSiDBasicKalmanToolTrack.cxx.

49{
50 // prepare a vector with new observation point coordinates
51 CLHEP::HepVector newPoint (2);
52 newPoint[0] = cluster->xLocal();
53 newPoint[1] = cluster->yLocal();
54
55 // get the values for the predicted new point
56 const double newZ = cluster->zLocal();
57
58 const CLHEP::HepVector predictedTruePosition = predictNextPoint(newZ);
59 const CLHEP::HepMatrix predictedCovariance = predictNextPointCov(newZ);
60
61 // difference between measured and predicted (called innovation y_k and S_k)
62 const CLHEP::HepVector residualPredicted = newPoint - m_observationModel*predictedTruePosition;
63 const CLHEP::HepMatrix residualPredictedCov = m_observationNoise + m_observationModel*predictedCovariance*m_observationModel.T();
64
65 // optimal Kalman gain (Kk)
66 const CLHEP::HepMatrix kalmanGain = predictedCovariance*m_observationModel.T()*qr_inverse(residualPredictedCov);
67
68 // === calculated new values ===
69 // Updated (a posteriori) state estimate
70 const CLHEP::HepVector newPosition = predictedTruePosition + kalmanGain*residualPredicted;
71
72 // Updated (a posteriori) residuals
73 const CLHEP::HepVector newResidual = newPoint - m_observationModel*newPosition;
74 const CLHEP::HepMatrix newResidualCov = m_observationNoise - m_observationModel*kalmanGain*m_observationNoise;
75
76 // chi2 statistics
77 const CLHEP::HepVector chi2Vector = newResidual.T()*qr_inverse(newResidualCov)*newResidual;
78
79 // check if hit is good
80 if (chi2Vector[0] < clusterMaxChi2) {
81 // add cluster and all information
82 m_clustersInTrack.push_back(cluster);
83 m_zPositionHistory.push_back(newZ);
84 m_transitionHistory.emplace_back(transitionModel(newZ));
85 m_positionAndSlopePredictedHistory.push_back(predictedTruePosition);
86 m_positionAndSlopeCovPredictedHistory.push_back(predictedCovariance);
87 m_positionAndSlopeHistory.push_back(newPosition);
88 m_chi2History.emplace_back(chi2Vector[0]);
89 // Updated (a posteriori) estimate covariance
90 m_positionAndSlopeCovHistory.emplace_back(predictedCovariance - kalmanGain*m_observationModel*predictedCovariance);
91 }
92 else {
93 addHole();
94 }
95}
std::list< CLHEP::HepVector > m_positionAndSlopePredictedHistory
Vectors of predicted true position of a track for each step.
void addHole()
Increase number of holes by 1. (see m_holes)
std::list< CLHEP::HepMatrix > m_positionAndSlopeCovPredictedHistory
Matrices of predicted estimate covariance of the position and slope for each step.
CLHEP::HepVector predictNextPoint(const double zPosition) const
Predict x and y positions and slopes of the track at given Z.
CLHEP::HepMatrix predictNextPointCov(const double zPosition) const
Predict covariance matrix of the track at given Z.
std::list< CLHEP::HepMatrix > m_transitionHistory
Matrices of transitionModel used at each step of track reconstruction.
CLHEP::HepMatrix transitionModel(const double firstPos, const double secondPos) const
Create state transition matrix to calculate x and y position and slopes between to Z positions.

◆ addHole()

void AFPSiDBasicKalmanToolTrack::addHole ( )
inline

Increase number of holes by 1. (see m_holes)

Definition at line 97 of file AFPSiDBasicKalmanToolTrack.h.

97{++m_holes;}

◆ chi2Smooth()

const std::list< double > & AFPSiDBasicKalmanToolTrack::chi2Smooth ( ) const
inline

History of chi2 for each step after smoothing.

Definition at line 112 of file AFPSiDBasicKalmanToolTrack.h.

112{return m_chi2Smooth;}
std::list< double > m_chi2Smooth
History of chi2 for each step after smoothing.

◆ clearSmoothMatrices()

void AFPSiDBasicKalmanToolTrack::clearSmoothMatrices ( )

Clears matrices used for smoothing the track.

Definition at line 115 of file AFPSiDBasicKalmanToolTrack.cxx.

116{
119 m_chi2Smooth.clear();
120}
std::list< CLHEP::HepMatrix > m_positionAndSlopeSmoothCov
Matrices of updated (a posteriori) estimate covariance of the position and slope for each step after ...
std::list< CLHEP::HepVector > m_positionAndSlopeSmooth
Vectors of reconstructed true positions for each step after smoothing .

◆ clustersInTrack()

const std::list< const xAOD::AFPSiHitsCluster * > & AFPSiDBasicKalmanToolTrack::clustersInTrack ( ) const
inline

Clusters used to reconstruct the track.

Definition at line 100 of file AFPSiDBasicKalmanToolTrack.h.

100{return m_clustersInTrack;}

◆ findNearestCluster()

const xAOD::AFPSiHitsCluster * AFPSiDBasicKalmanToolTrack::findNearestCluster ( const std::vector< const xAOD::AFPSiHitsCluster * > & clusters,
const float maxAllowedDistance ) const

Finds the cluster which is nearest to the track, returns nullptr if no is found.

The method predicts what are track X and Y coordinates for the Z of a cluster. Next, loops over the clusters in the vector and calculates distance \(r = \sqrt{\Delta x^2 + \Delta y^2}\). Returns nullptr if there is no cluster which is closer than #maxAllowedDistance or if there are clusters with r < #maxAllowedDistance a pointer to the closest one.

Parameters
clustersvector of clusters from which the closest one is to be found
maxAllowedDistanceif no cluster is closer than this argument nullptr is returned

Definition at line 97 of file AFPSiDBasicKalmanToolTrack.cxx.

98{
99 double minDisntace = maxAllowedDistance;
100 const xAOD::AFPSiHitsCluster* nearestCluster = nullptr;
101 for (const xAOD::AFPSiHitsCluster* theCluster : clusters) {
102 const CLHEP::HepVector predictedPoint = predictNextPoint(theCluster->zLocal());
103 const double xDiff = predictedPoint[0] - theCluster->xLocal();
104 const double yDiff = predictedPoint[2] - theCluster->yLocal(); // index 2, because this is (x, dx, y, dy) vector
105 const double distance = std::sqrt(xDiff*xDiff + yDiff*yDiff);
106 if (distance < minDisntace) {
107 minDisntace = distance;
108 nearestCluster = theCluster;
109 }
110 }
111
112 return nearestCluster;
113}
float distance(const Amg::Vector3D &p1, const Amg::Vector3D &p2)
calculates the distance between two point in 3D space
AFPSiHitsCluster_v1 AFPSiHitsCluster

◆ holes()

int AFPSiDBasicKalmanToolTrack::holes ( ) const
inline

Definition at line 94 of file AFPSiDBasicKalmanToolTrack.h.

94{return m_holes;}

◆ positionAndSlopeHistory()

const std::list< CLHEP::HepVector > & AFPSiDBasicKalmanToolTrack::positionAndSlopeHistory ( ) const
inline

Vectors of reconstructed true positions for each step \(x_{k}\).

Each 4D vector containing x, y position of the track in a given step as well as slopes in these two directions \(\Delta x, \Delta y\). \(x_{k} = (x, \Delta x, y, \Delta y)\)

Definition at line 109 of file AFPSiDBasicKalmanToolTrack.h.

◆ positionAndSlopeSmooth()

const std::list< CLHEP::HepVector > & AFPSiDBasicKalmanToolTrack::positionAndSlopeSmooth ( ) const
inline

Vectors of reconstructed true positions for each step after smoothing \(x_{k}\).

Each 4D vector containing x, y position of the track in a given step as well as slopes in these two directions \(\Delta x, \Delta y\). \(x_{k} = (x, \Delta x, y, \Delta y)\)

Definition at line 106 of file AFPSiDBasicKalmanToolTrack.h.

◆ predictNextPoint()

CLHEP::HepVector AFPSiDBasicKalmanToolTrack::predictNextPoint ( const double zPosition) const
inline

Predict x and y positions and slopes of the track at given Z.

Definition at line 67 of file AFPSiDBasicKalmanToolTrack.h.

68 {return transitionModel(zPosition)*m_positionAndSlopeHistory.back();}

◆ predictNextPointCov()

CLHEP::HepMatrix AFPSiDBasicKalmanToolTrack::predictNextPointCov ( const double zPosition) const
inline

Predict covariance matrix of the track at given Z.

Definition at line 71 of file AFPSiDBasicKalmanToolTrack.h.

72 {CLHEP::HepMatrix transition = transitionModel(zPosition);
73 return transition*m_positionAndSlopeCovHistory.back()*transition.T();}

◆ smooth()

void AFPSiDBasicKalmanToolTrack::smooth ( )

Run smoothing algorithm.

The algorithm should probably just reconstruct track starting from the last hit, but the covaraince matrices are changed, so the result maybe slightly different. Algorithm is just copied from the previous implementation.

Definition at line 122 of file AFPSiDBasicKalmanToolTrack.cxx.

123{
125
126 // prepare iterators for going backwards through the lists
127 std::list< CLHEP::HepVector >::const_iterator positionIter = m_positionAndSlopeHistory.end();
128 std::list< CLHEP::HepVector >::const_iterator positionBegin = m_positionAndSlopeHistory.begin();
129 // if no position is reconstructed there is nothing to be done
130 if (positionIter != positionBegin) {
131 // prepare remaining iterators
132 std::list< CLHEP::HepMatrix >::const_iterator covIter = m_positionAndSlopeCovHistory.end();
133 std::list< const xAOD::AFPSiHitsCluster* >::const_iterator clustersIter = m_clustersInTrack.end();
134 std::list< CLHEP::HepVector >::const_iterator positionPredictIter = m_positionAndSlopePredictedHistory.end();
135 std::list< CLHEP::HepMatrix >::const_iterator covPredictIter = m_positionAndSlopeCovPredictedHistory.end();
136 std::list< CLHEP::HepMatrix >::const_iterator transitionIter = m_transitionHistory.end();
137
138 // === prepare first point ===
139 // iterator moved back by 1, because end() points after the last element
140 const CLHEP::HepVector& lastPosition = *(--positionIter);
141 const CLHEP::HepMatrix& lastCovariance = *(--covIter);
142 const xAOD::AFPSiHitsCluster* lastCluster = *(--clustersIter);
143
144 m_positionAndSlopeSmooth.push_back(lastPosition);
145 m_positionAndSlopeSmoothCov.push_back(lastCovariance);
146
147 // prepare a vector with new observation point coordinates
148 CLHEP::HepVector lastPoint (2);
149 lastPoint[0] = lastCluster->xLocal();
150 lastPoint[1] = lastCluster->yLocal();
151
152 const CLHEP::HepVector residualSmooth = lastPoint - m_observationModel*lastPosition;
153 const CLHEP::HepMatrix residualSmoothCov = m_observationNoise - m_observationModel*lastCovariance*m_observationModel.T();
154 const CLHEP::HepVector chi2Smooth = residualSmooth.T()*qr_inverse(residualSmoothCov)*residualSmooth;
155 m_chi2Smooth.push_back(chi2Smooth[0]);
156
157 // === loop over remaining clusters ===
158 while (positionIter != positionBegin) {
159 // --- read needed variables ---
160 const CLHEP::HepVector& position = *(--positionIter);
161 const CLHEP::HepMatrix& covariance = *(--covIter);
162 const xAOD::AFPSiHitsCluster* cluster = *(--clustersIter);
163
164 // shifted by 1 in history with respect to the position and covariance
165 const CLHEP::HepVector& prevPositionPred = *(--positionPredictIter);
166 const CLHEP::HepMatrix& prevCovPred = *(--covPredictIter);
167 const CLHEP::HepMatrix& prevTransition = *(--transitionIter);
168
169 // --- calculate the smoothed position and covariance ---
170 const CLHEP::HepMatrix Ai = qr_inverse(prevTransition); // when Qk == 0
171 // const CLHEP::HepMatrix Ai = covariance * prevTransition.T() * qr_inverse(prevCovPred);
172 const CLHEP::HepVector positionSmooth = position + Ai*( m_positionAndSlopeSmooth.back() - prevPositionPred );
173 const CLHEP::HepMatrix covarianceSmooth = covariance + Ai*( m_positionAndSlopeSmoothCov.back() - prevCovPred )*Ai.T();
174
175 // --- calculate smoothed chi2 ---
176 CLHEP::HepVector point (2);
177 point[0] = cluster->xLocal();
178 point[1] = cluster->yLocal();
179 const CLHEP::HepVector residualSmooth = point - m_observationModel*positionSmooth;
180 const CLHEP::HepMatrix residualCovSmooth = m_observationNoise - m_observationModel*covarianceSmooth*m_observationModel.T();
181 CLHEP::HepVector chi2Smooth = residualSmooth.T()*qr_inverse(residualCovSmooth)*residualSmooth;
182
183 // --- add results to the lists ---
184 m_positionAndSlopeSmooth.push_back(positionSmooth);
185 m_positionAndSlopeSmoothCov.push_back(covarianceSmooth);
186 m_chi2Smooth.push_back(chi2Smooth[0]);
187 } // end while over points/clusters
188
189
190 // === extrapolate to the first plane ===
191 clustersIter = m_clustersInTrack.begin();
192 const xAOD::AFPSiHitsCluster* firstCluster = *clustersIter;
193 const xAOD::AFPSiHitsCluster* secondCluster = *(++clustersIter);
194 const CLHEP::HepMatrix firstTransition = transitionModel(firstCluster->zLocal(), secondCluster->zLocal());
195
196 const CLHEP::HepMatrix Ai = qr_inverse(firstTransition); // when Qk == 0
197 const CLHEP::HepVector firstPosition = Ai*(m_positionAndSlopeSmooth.back());
198 const CLHEP::HepMatrix firstCovariance = m_positionAndSlopeCovHistory.front() + Ai*( m_positionAndSlopeSmoothCov.back() - m_positionAndSlopeCovHistory.front())*Ai.T();
199
200 // --- calculate smoothed chi2 ---
201 const CLHEP::HepVector firstResidual = m_observationModel*Ai*(m_positionAndSlopeHistory.front() - m_positionAndSlopeSmooth.back());
202 CLHEP::HepVector firstPoint (2);
203 firstPoint[0] = firstCluster->xLocal();
204 firstPoint[1] = firstCluster->yLocal();
205 // const CLHEP::HepVector firstResidual = firstPoint - m_observationModel*firstPosition;
206 const CLHEP::HepMatrix firstResidualCov = m_observationNoise - m_observationModel*firstCovariance*m_observationModel.T();
207 CLHEP::HepVector firstChi2 = firstResidual.T()*qr_inverse(firstResidualCov)*firstResidual;
208
209 // --- add results to the lists ---
210 m_positionAndSlopeSmooth.push_back(firstPosition);
211 m_positionAndSlopeSmoothCov.push_back(firstCovariance);
212 m_chi2Smooth.push_back(firstChi2[0]);
213 }
214}
void clearSmoothMatrices()
Clears matrices used for smoothing the track.
const std::list< double > & chi2Smooth() const
History of chi2 for each step after smoothing.

◆ transitionModel() [1/2]

CLHEP::HepMatrix AFPSiDBasicKalmanToolTrack::transitionModel ( const double firstPos,
const double secondPos ) const
inline

Create state transition matrix to calculate x and y position and slopes between to Z positions.

The matrix leaves slopes unchanged and adds to the previous positions slopes multiplied by the z distance i.e. secondPos - firstPos.

\[\left( \begin{array}{cccc} 1 & \Delta z & 0 & 0 0 & 1 & 0 & 0 0 & 0 & 1 & \Delta z 0 & 0 & 0 & 1 \end{array} \right) \]

Definition at line 187 of file AFPSiDBasicKalmanToolTrack.h.

188{
189 CLHEP::HepMatrix transitionModel (4, 4, 1); // create unit matrix (1 on diagonal, 0 offdiagonal)
190 const double zDistance = secondPos - firstPos;
191 transitionModel[0][1] = zDistance;// set linear transformation which sets x' = x + dx*dz
192 transitionModel[2][3] = zDistance; // set linear transformation which sets y' = y + dy*dz
193
194 return transitionModel;
195}

◆ transitionModel() [2/2]

CLHEP::HepMatrix AFPSiDBasicKalmanToolTrack::transitionModel ( const double zPosition) const
inline

Create state transition matrix to calculate x and y position and slopes at given Z.

The matrix for transition between the last point and given Z position is calculated.

Definition at line 63 of file AFPSiDBasicKalmanToolTrack.h.

64 {return transitionModel (m_zPositionHistory.back(), zPosition);}

◆ trkChi2NDFSmooth()

double AFPSiDBasicKalmanToolTrack::trkChi2NDFSmooth ( ) const

Chi2 per degrees of freedom (clusters) for the whole smoothed track.

Sum chi2 for each cluster and divide by the number of clusters.

Definition at line 216 of file AFPSiDBasicKalmanToolTrack.cxx.

217{
218 double chi2 = 0;
219 for (double val : m_chi2Smooth) {
220 chi2 += val;
221 }
222
223 return chi2/m_chi2Smooth.size();
224}
double chi2(TH1 *h0, TH1 *h1)

Member Data Documentation

◆ m_chi2History

std::list< double> AFPSiDBasicKalmanToolTrack::m_chi2History
protected

History of chi2 for each step.

Definition at line 156 of file AFPSiDBasicKalmanToolTrack.h.

◆ m_chi2Smooth

std::list< double > AFPSiDBasicKalmanToolTrack::m_chi2Smooth
protected

History of chi2 for each step after smoothing.

Definition at line 169 of file AFPSiDBasicKalmanToolTrack.h.

◆ m_clustersInTrack

std::list<const xAOD::AFPSiHitsCluster*> AFPSiDBasicKalmanToolTrack::m_clustersInTrack
protected

Clusters used to reconstruct the track.

Definition at line 172 of file AFPSiDBasicKalmanToolTrack.h.

◆ m_holes

int AFPSiDBasicKalmanToolTrack::m_holes
protected

Number of layers in which hit from a track is expected, but is not present.

Definition at line 184 of file AFPSiDBasicKalmanToolTrack.h.

◆ m_observationModel

const CLHEP::HepMatrix& AFPSiDBasicKalmanToolTrack::m_observationModel
protected

Observation model, describing transformation from true to measured state.

Definition at line 175 of file AFPSiDBasicKalmanToolTrack.h.

◆ m_observationNoise

const CLHEP::HepMatrix& AFPSiDBasicKalmanToolTrack::m_observationNoise
protected

Matrix representing observation noise disturbing the measurement.

Used to estimate covariance matrix of the track.

Definition at line 180 of file AFPSiDBasicKalmanToolTrack.h.

◆ m_positionAndSlopeCovHistory

std::list< CLHEP::HepMatrix > AFPSiDBasicKalmanToolTrack::m_positionAndSlopeCovHistory
protected

Matrices of updated (a posteriori) estimate covariance of the position and slope for each step \(P_{k|k}\).

Definition at line 147 of file AFPSiDBasicKalmanToolTrack.h.

◆ m_positionAndSlopeCovPredictedHistory

std::list< CLHEP::HepMatrix > AFPSiDBasicKalmanToolTrack::m_positionAndSlopeCovPredictedHistory
protected

Matrices of predicted estimate covariance of the position and slope for each step.

Definition at line 150 of file AFPSiDBasicKalmanToolTrack.h.

◆ m_positionAndSlopeHistory

std::list< CLHEP::HepVector > AFPSiDBasicKalmanToolTrack::m_positionAndSlopeHistory
protected

Vectors of reconstructed true positions for each step \(x_{k}\).

Each 4D vector containing x, y position of the track in a given step as well as slopes in these two directions \(\Delta x, \Delta y\). \(x_{k} = (x, \Delta x, y, \Delta y)\)

Definition at line 137 of file AFPSiDBasicKalmanToolTrack.h.

◆ m_positionAndSlopePredictedHistory

std::list< CLHEP::HepVector > AFPSiDBasicKalmanToolTrack::m_positionAndSlopePredictedHistory
protected

Vectors of predicted true position of a track for each step.

Each 4D vector containing x, y position of the track in a given step as well as slopes in these two directions \(\Delta x, \Delta y\). \(x_{k} = (x, \Delta x, y, \Delta y)\)

Definition at line 144 of file AFPSiDBasicKalmanToolTrack.h.

◆ m_positionAndSlopeSmooth

std::list< CLHEP::HepVector > AFPSiDBasicKalmanToolTrack::m_positionAndSlopeSmooth
protected

Vectors of reconstructed true positions for each step after smoothing \(x_{k}\).

Each 4D vector containing x, y position of the track in a given step as well as slopes in these two directions \(\Delta x, \Delta y\). \(x_{k} = (x, \Delta x, y, \Delta y)\)

Definition at line 163 of file AFPSiDBasicKalmanToolTrack.h.

◆ m_positionAndSlopeSmoothCov

std::list< CLHEP::HepMatrix > AFPSiDBasicKalmanToolTrack::m_positionAndSlopeSmoothCov
protected

Matrices of updated (a posteriori) estimate covariance of the position and slope for each step after smoothing \(P_{k|k}\).

Definition at line 166 of file AFPSiDBasicKalmanToolTrack.h.

◆ m_transitionHistory

std::list< CLHEP::HepMatrix > AFPSiDBasicKalmanToolTrack::m_transitionHistory
protected

Matrices of transitionModel used at each step of track reconstruction.

Definition at line 153 of file AFPSiDBasicKalmanToolTrack.h.

◆ m_zPositionHistory

std::list< float > AFPSiDBasicKalmanToolTrack::m_zPositionHistory
protected

Vectors of reconstructed true positions for each step.

Definition at line 130 of file AFPSiDBasicKalmanToolTrack.h.


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