GeneralBrokenLines V03-01-03
using EIGEN
|
Point on trajectory. More...
#include <GblPoint.h>
Public Member Functions | |
GblPoint (const Matrix5d &aJacobian, unsigned int numMeasReserve=0) | |
Create a point. More... | |
GblPoint (const GblPoint &)=default | |
GblPoint & | operator= (const GblPoint &)=default |
GblPoint (GblPoint &&)=default | |
GblPoint & | operator= (GblPoint &&)=default |
virtual | ~GblPoint () |
template<typename Projection , typename Residuals , typename Precision , typename std::enable_if<(Precision::ColsAtCompileTime !=1)>::type * = nullptr> | |
void | addMeasurement (const Eigen::MatrixBase< Projection > &aProjection, const Eigen::MatrixBase< Residuals > &aResiduals, const Eigen::MatrixBase< Precision > &aPrecision, double minPrecision=0.) |
Add a measurement to a point. More... | |
template<typename Projection , typename Residuals , typename Precision , typename std::enable_if<(Precision::ColsAtCompileTime==1)>::type * = nullptr> | |
void | addMeasurement (const Eigen::MatrixBase< Projection > &aProjection, const Eigen::MatrixBase< Residuals > &aResiduals, const Eigen::MatrixBase< Precision > &aPrecision, double minPrecision=0.) |
Add a measurement to a point. More... | |
template<typename Residuals , typename Precision , typename std::enable_if<(Precision::ColsAtCompileTime !=1)>::type * = nullptr> | |
void | addMeasurement (const Eigen::MatrixBase< Residuals > &aResiduals, const Eigen::MatrixBase< Precision > &aPrecision, double minPrecision=0.) |
Add a measurement to a point. More... | |
template<typename Residuals , typename Precision , typename std::enable_if<(Precision::ColsAtCompileTime==1)>::type * = nullptr> | |
void | addMeasurement (const Eigen::MatrixBase< Residuals > &aResiduals, const Eigen::MatrixBase< Precision > &aPrecision, double minPrecision=0.) |
Add a measurement to a point. More... | |
void | addScatterer (const Eigen::Vector2d &aResiduals, const Eigen::Matrix2d &aPrecision) |
Add a thin scatterer to a point. More... | |
void | addScatterer (const Eigen::Vector2d &aResiduals, const Eigen::Vector2d &aPrecision) |
Add a thin scatterer to a point. More... | |
void | addThickScatterer (const Eigen::Vector4d &aResiduals, const Eigen::Matrix4d &aPrecision) |
Add a thick scatterer to a point. More... | |
template<typename Derivative > | |
void | addLocals (const Eigen::MatrixBase< Derivative > &aDerivatives) |
Add local derivatives to a point. More... | |
template<typename Derivative > | |
void | addGlobals (const std::vector< int > &aLabels, const Eigen::MatrixBase< Derivative > &aDerivatives) |
Add global derivatives to a point. More... | |
unsigned int | numMeasurements () const |
Check for measurements at a point. More... | |
unsigned int | getScatDim () const |
Get scatterer dimension. More... | |
void | getScatterer (Eigen::Matrix4d &aTransformation, Eigen::Vector4d &aResiduals, Eigen::Vector4d &aPrecision) const |
Retrieve scatterer of a point. More... | |
void | getReducedScatterer (Eigen::Matrix4d &aTransformation, Eigen::Vector4d &aResiduals, Eigen::Vector4d &aPrecision) const |
Retrieve (reduced) scatterer of a point. More... | |
void | getScatTransformation (Eigen::MatrixXd &aTransformation) const |
Get scatterer transformation (from diagonalization). More... | |
unsigned int | getLabel () const |
Retrieve label of point. More... | |
int | getOffset () const |
Retrieve offset for point. More... | |
const Matrix5d & | getP2pJacobian () const |
Retrieve point-to-(previous)point jacobian. More... | |
void | getDerivatives (int aDirection, Eigen::Matrix2d &matW, Eigen::Matrix2d &matWJ, Eigen::Vector2d &vecWd) const |
Retrieve derivatives of local track model. More... | |
void | printPoint (unsigned int level=0) const |
Print GblPoint. More... | |
std::vector< GblMeasurement >::iterator | getMeasBegin () |
Get GblMeasurement iterator for begin. More... | |
std::vector< GblMeasurement >::iterator | getMeasEnd () |
Get GblMeasurement iterator for end. More... | |
void | getGlobalLabelsAndDerivatives (unsigned int aMeas, unsigned int aRow, std::vector< int > &aLabels, std::vector< double > &aDerivatives) const |
Retrieve global derivatives from a measurement at a point for a single row. More... | |
bool | isFirst () const |
Is first point on (sub)trajectory? More... | |
bool | isLast () const |
Is last point on (sub)trajectory? More... | |
Private Member Functions | |
void | setLabel (unsigned int aLabel) |
Define label of point (by GBLTrajectory constructor) More... | |
void | setOffset (int anOffset) |
Define offset for point (by GBLTrajectory constructor) More... | |
void | setType (int aType) |
Define type for point (by GBLTrajectory constructor) More... | |
void | addPrevJacobian (const Matrix5d &aJac) |
Define jacobian to previous scatterer (by GBLTrajectory constructor) More... | |
void | addNextJacobian (const Matrix5d &aJac) |
Define jacobian to next scatterer (by GBLTrajectory constructor) More... | |
Private Attributes | |
unsigned int | theLabel |
Label identifying point. More... | |
int | theOffset |
Offset number at point if not negative (else interpolation needed) More... | |
int | theType |
Type (-1: first, 0: inner, 1: last) More... | |
Matrix5d | p2pJacobian = Matrix5d::Zero(5,5) |
Point-to-point jacobian from previous point. More... | |
Matrix5d | prevJacobian = Matrix5d::Zero(5,5) |
Jacobian to previous scatterer (or first measurement) More... | |
Matrix5d | nextJacobian = Matrix5d::Zero(5,5) |
Jacobian to next scatterer (or last measurement) More... | |
unsigned int | scatDim |
Dimension of scatterer (0: none, 2: thin, 4:thick) More... | |
Eigen::Matrix4d | scatTransformation |
Transformation of diagonalization (of scat. precision matrix) More... | |
Eigen::Vector4d | scatResiduals |
Scattering residuals (initial kinks if iterating) More... | |
Eigen::Vector4d | scatPrecision |
Scattering precision (diagonal of inverse covariance matrix) More... | |
std::vector< GblMeasurement > | theMeasurements |
List of measurements at point. More... | |
Friends | |
class | GblTrajectory |
Point on trajectory.
User supplied point on (initial) trajectory.
Must have jacobian for propagation from previous point. May have:
Definition at line 62 of file GblPoint.h.
gbl::GblPoint::GblPoint | ( | const Matrix5d & | aJacobian, |
unsigned int | numMeasReserve = 0 |
||
) |
Create a point.
Create point on (initial) trajectory. Needs transformation jacobian from previous point.
[in] | aJacobian | Transformation jacobian from previous point |
[in] | numMeasReserve | number of measurements to reserve (space for) |
Definition at line 42 of file GblPoint.cpp.
References theMeasurements.
|
default |
|
default |
|
virtual |
Definition at line 65 of file GblPoint.cpp.
void gbl::GblPoint::addGlobals | ( | const std::vector< int > & | aLabels, |
const Eigen::MatrixBase< Derivative > & | aDerivatives | ||
) |
Add global derivatives to a point.
Point needs to have a measurement.
Derivative | Derivatives matrix |
[in] | aLabels | Global derivatives labels |
[in] | aDerivatives | Global derivatives (matrix) |
Definition at line 293 of file GblPoint.h.
References theMeasurements.
Referenced by exampleComposedGeo(), exampleComposedKin(), exampleDc(), exampleSit(), and GblPoint_addGlobals().
void gbl::GblPoint::addLocals | ( | const Eigen::MatrixBase< Derivative > & | aDerivatives | ) |
Add local derivatives to a point.
Point needs to have a measurement.
Derivative | Derivatives matrix |
[in] | aDerivatives | Local derivatives (matrix) |
Definition at line 287 of file GblPoint.h.
References theMeasurements.
Referenced by example3().
void gbl::GblPoint::addMeasurement | ( | const Eigen::MatrixBase< Projection > & | aProjection, |
const Eigen::MatrixBase< Residuals > & | aResiduals, | ||
const Eigen::MatrixBase< Precision > & | aPrecision, | ||
double | minPrecision = 0. |
||
) |
Add a measurement to a point.
Add measurement (in measurement system) with arbitrary precision (inverse covariance) matrix. Will be diagonalized. ((up to) 2D: position, 4D: slope+position, 5D: curvature+slope+position)
Projection | Projection matrix |
Residuals | Residuals vector |
Precision | Precision matrix or vector (with diagonal) |
[in] | aProjection | Projection from local to measurement system (derivative of measurement vs local parameters) |
[in] | aResiduals | Measurement residuals |
[in] | aPrecision | Measurement precision (matrix) |
[in] | minPrecision | Minimal precision to accept measurement |
Definition at line 238 of file GblPoint.h.
Referenced by example1(), example2(), example3(), example4(), exampleComposedGeo(), exampleComposedKin(), exampleDc(), exampleSit(), and GblPoint_addMeasurement2D().
void gbl::GblPoint::addMeasurement | ( | const Eigen::MatrixBase< Projection > & | aProjection, |
const Eigen::MatrixBase< Residuals > & | aResiduals, | ||
const Eigen::MatrixBase< Precision > & | aPrecision, | ||
double | minPrecision = 0. |
||
) |
Add a measurement to a point.
Add measurement (in measurement system) with diagonal precision (inverse covariance) matrix. ((up to) 2D: position, 4D: slope+position, 5D: curvature+slope+position)
Projection | Projection matrix |
Residuals | Residuals vector |
Precision | Precision matrix or vector (with diagonal) |
[in] | aProjection | Projection from local to measurement system (derivative of measurement vs local parameters) |
[in] | aResiduals | Measurement residuals |
[in] | aPrecision | Measurement precision (vector with diagonal) |
[in] | minPrecision | Minimal precision to accept measurement |
void gbl::GblPoint::addMeasurement | ( | const Eigen::MatrixBase< Residuals > & | aResiduals, |
const Eigen::MatrixBase< Precision > & | aPrecision, | ||
double | minPrecision = 0. |
||
) |
Add a measurement to a point.
Add measurement in local system with arbitrary precision (inverse covariance) matrix. Will be diagonalized. ((up to) 2D: position, 4D: slope+position, 5D: curvature+slope+position)
Residuals | Residuals vector |
Precision | Precision matrix or vector (with diagonal) |
[in] | aResiduals | Measurement residuals |
[in] | aPrecision | Measurement precision (matrix) |
[in] | minPrecision | Minimal precision to accept measurement |
Definition at line 268 of file GblPoint.h.
void gbl::GblPoint::addMeasurement | ( | const Eigen::MatrixBase< Residuals > & | aResiduals, |
const Eigen::MatrixBase< Precision > & | aPrecision, | ||
double | minPrecision = 0. |
||
) |
Add a measurement to a point.
Add measurement in local system with diagonal precision (inverse covariance) matrix. ((up to) 2D: position, 4D: slope+position, 5D: curvature+slope+position)
Residuals | Residuals vector |
Precision | Precision matrix or vector (with diagonal) |
[in] | aResiduals | Measurement residuals |
[in] | aPrecision | Measurement precision (vector with diagonal) |
[in] | minPrecision | Minimal precision to accept measurement |
|
private |
Define jacobian to next scatterer (by GBLTrajectory constructor)
[in] | aJac | Jacobian |
Definition at line 508 of file GblPoint.cpp.
References nextJacobian.
Referenced by gbl::GblTrajectory::calcJacobians().
|
private |
Define jacobian to previous scatterer (by GBLTrajectory constructor)
[in] | aJac | Jacobian |
Definition at line 493 of file GblPoint.cpp.
References prevJacobian.
void gbl::GblPoint::addScatterer | ( | const Eigen::Vector2d & | aResiduals, |
const Eigen::Matrix2d & | aPrecision | ||
) |
Add a thin scatterer to a point.
Add scatterer with arbitrary precision (inverse covariance) matrix. Will be diagonalized. Changes local track direction.
The precision matrix for the local slopes is defined by the angular scattering error theta_0 and the scalar products c_1, c_2 of the offset directions in the local frame with the track direction:
(1 - c_1*c_1 - c_2*c_2) | 1 - c_1*c_1 - c_1*c_2 | P = ----------------------- * | | theta_0*theta_0 | - c_1*c_2 1 - c_2*c_2 |
One offset is defined at the point. The comparison with the previous and next offsets allows to define a kink.
Precision | Precision matrix or vector (with diagonal) |
[in] | aResiduals | Scatterer residuals |
[in] | aPrecision | Scatterer precision (full matrix) |
Definition at line 159 of file GblPoint.cpp.
References scatDim, scatPrecision, scatResiduals, and scatTransformation.
Referenced by example1(), example2(), example3(), example4(), exampleComposedGeo(), exampleComposedKin(), exampleDc(), exampleSit(), and GblPoint_addScatterer().
void gbl::GblPoint::addScatterer | ( | const Eigen::Vector2d & | aResiduals, |
const Eigen::Vector2d & | aPrecision | ||
) |
Add a thin scatterer to a point.
Add scatterer with diagonal precision (inverse covariance) matrix. Changes local track direction.
The precision matrix for the local slopes is defined by the angular scattering error theta_0 and the scalar products c_1, c_2 of the offset directions in the local frame with the track direction:
(1 - c_1*c_1 - c_2*c_2) | 1 - c_1*c_1 - c_1*c_2 | P = ----------------------- * | | theta_0*theta_0 | - c_1*c_2 1 - c_2*c_2 |
For a diagonal precision matrix at least one of the scalar products c_1, c_2 must be zero.
One offset is defined at the point. The comparison with the previous and next offsets allows to define a kink.
Precision | Precision matrix or vector (with diagonal) |
[in] | aResiduals | Scatterer residuals |
[in] | aPrecision | Scatterer precision (vector with diagonal) |
Definition at line 194 of file GblPoint.cpp.
References scatDim, scatPrecision, scatResiduals, and scatTransformation.
void gbl::GblPoint::addThickScatterer | ( | const Eigen::Vector4d & | aResiduals, |
const Eigen::Matrix4d & | aPrecision | ||
) |
Add a thick scatterer to a point.
Add scatterer with arbitrary precision (inverse 4x4 covariance) matrix. Will be diagonalized. Local track direction and position change (correlated).
Combination of several scatterers extrapolated to point (yielding dense covariance matrix). E.g. sum of thin scatterers (covariance V_i, jacobian J_i for slopes and offsets):
P = [sum(J_i*V_i*J_i^t)]^-1
For each thin scatterer the covarinace matrix V_i is defined by the angular scattering error theta_0,i and the scalar products c_1,i, c_2,i of the offset directions in the local frame with the track direction:
| 1 - c_2,i^2 c_1,i*c_2,i 0 0 | | | theta_0,i^2 | c_1,i*c_2,i 1 - c_1,i^2 0 0 | V_i = ------------------------- * | | (1 - c_1,i^2 - c_2,i^2)^2 | 0 0 0 0 | | | | 0 0 0 0 |
Two offsets are defined at the point to describe a step (in the trajectory). The comparison with the previous and next offsets allows to define a kink.
Precision | Precision matrix |
[in] | aResiduals | Scatterer residuals |
[in] | aPrecision | Scatterer precision (full 4x4 matrix) |
Definition at line 232 of file GblPoint.cpp.
References scatDim, scatPrecision, scatResiduals, and scatTransformation.
Referenced by example4().
void gbl::GblPoint::getDerivatives | ( | int | aDirection, |
Eigen::Matrix2d & | matW, | ||
Eigen::Matrix2d & | matWJ, | ||
Eigen::Vector2d & | vecWd | ||
) | const |
Retrieve derivatives of local track model.
Linearized track model: F_u(q/p,u',u) = J*u + S*u' + d*q/p, W is inverse of S, negated for backward propagation.
[in] | aDirection | Propagation direction (>0 forward, else backward) |
[out] | matW | W |
[out] | matWJ | W*J |
[out] | vecWd | W*d |
std::overflow_error | : matrix S is singular. |
Definition at line 522 of file GblPoint.cpp.
References nextJacobian, and prevJacobian.
Referenced by gbl::GblTrajectory::getFitToKinkAndStepJacobian(), gbl::GblTrajectory::getFitToKinkJacobian(), and gbl::GblTrajectory::getFitToLocalJacobian().
void gbl::GblPoint::getGlobalLabelsAndDerivatives | ( | unsigned int | aMeas, |
unsigned int | aRow, | ||
std::vector< int > & | aLabels, | ||
std::vector< double > & | aDerivatives | ||
) | const |
Retrieve global derivatives from a measurement at a point for a single row.
[in] | aMeas | Measurement number |
[in] | aRow | Row number |
[out] | aLabels | Global labels |
[out] | aDerivatives | Global derivatives |
Definition at line 625 of file GblPoint.cpp.
References theMeasurements.
Referenced by GblPoint_getGlobalLabelsAndDerivatives().
unsigned int gbl::GblPoint::getLabel | ( | ) | const |
std::vector< GblMeasurement >::iterator gbl::GblPoint::getMeasBegin | ( | ) |
Get GblMeasurement iterator for begin.
Definition at line 609 of file GblPoint.cpp.
References theMeasurements.
Referenced by GblPoint_getGlobalLabelsAndDerivatives(), and GblPoint_getNumMeasurements().
std::vector< GblMeasurement >::iterator gbl::GblPoint::getMeasEnd | ( | ) |
Get GblMeasurement iterator for end.
Definition at line 614 of file GblPoint.cpp.
References theMeasurements.
Referenced by GblPoint_getGlobalLabelsAndDerivatives(), and GblPoint_getNumMeasurements().
int gbl::GblPoint::getOffset | ( | ) | const |
Retrieve offset for point.
Definition at line 472 of file GblPoint.cpp.
References theOffset.
Referenced by gbl::GblTrajectory::getFitToKinkAndStepJacobian(), gbl::GblTrajectory::getFitToKinkJacobian(), gbl::GblTrajectory::getFitToLocalJacobian(), and gbl::GblTrajectory::getFitToStepJacobian().
const Matrix5d & gbl::GblPoint::getP2pJacobian | ( | ) | const |
Retrieve point-to-(previous)point jacobian.
Definition at line 485 of file GblPoint.cpp.
References p2pJacobian.
Referenced by gbl::GblTrajectory::calcJacobians().
void gbl::GblPoint::getReducedScatterer | ( | Eigen::Matrix4d & | aTransformation, |
Eigen::Vector4d & | aResiduals, | ||
Eigen::Vector4d & | aPrecision | ||
) | const |
Retrieve (reduced) scatterer of a point.
Reduce to positional scatterer.
[out] | aTransformation | Scatterer transformation from diagonalization |
[out] | aResiduals | Scatterer residuals |
[out] | aPrecision | Scatterer precision (diagonal) |
Definition at line 391 of file GblPoint.cpp.
References scatPrecision, scatResiduals, and scatTransformation.
unsigned int gbl::GblPoint::getScatDim | ( | ) | const |
Get scatterer dimension.
Get dimension of scatterer (0 = none, 2 = thin, 4 = thick).
Definition at line 365 of file GblPoint.cpp.
References scatDim.
Referenced by gbl::GblTrajectory::getFitToLocalJacobian().
void gbl::GblPoint::getScatterer | ( | Eigen::Matrix4d & | aTransformation, |
Eigen::Vector4d & | aResiduals, | ||
Eigen::Vector4d & | aPrecision | ||
) | const |
Retrieve scatterer of a point.
[out] | aTransformation | Scatterer transformation from diagonalization |
[out] | aResiduals | Scatterer residuals |
[out] | aPrecision | Scatterer precision (diagonal) |
Definition at line 375 of file GblPoint.cpp.
References scatDim, scatPrecision, scatResiduals, and scatTransformation.
void gbl::GblPoint::getScatTransformation | ( | Eigen::MatrixXd & | aTransformation | ) | const |
Get scatterer transformation (from diagonalization).
[out] | aTransformation | Transformation matrix |
Definition at line 417 of file GblPoint.cpp.
References scatDim, and scatTransformation.
bool gbl::GblPoint::isFirst | ( | ) | const |
Is first point on (sub)trajectory?
Definition at line 633 of file GblPoint.cpp.
References theType.
Referenced by printPoint().
bool gbl::GblPoint::isLast | ( | ) | const |
Is last point on (sub)trajectory?
Definition at line 638 of file GblPoint.cpp.
References theType.
Referenced by printPoint().
unsigned int gbl::GblPoint::numMeasurements | ( | ) | const |
Check for measurements at a point.
Get number of measurement (0 = none).
Definition at line 135 of file GblPoint.cpp.
References theMeasurements.
void gbl::GblPoint::printPoint | ( | unsigned int | level = 0 | ) | const |
Print GblPoint.
[in] | level | print level (0: minimum, >0: more) |
Definition at line 554 of file GblPoint.cpp.
References isFirst(), isLast(), nextJacobian, p2pJacobian, prevJacobian, scatDim, scatPrecision, scatResiduals, theLabel, theMeasurements, and theOffset.
Referenced by GblPoint_printPoint().
|
private |
Define label of point (by GBLTrajectory constructor)
[in] | aLabel | Label identifying point |
Definition at line 454 of file GblPoint.cpp.
References theLabel.
|
private |
Define offset for point (by GBLTrajectory constructor)
[in] | anOffset | Offset number |
Definition at line 467 of file GblPoint.cpp.
References theOffset.
|
private |
Define type for point (by GBLTrajectory constructor)
[in] | aType | Type (-1: first, 0: inner, 1: last) |
Definition at line 480 of file GblPoint.cpp.
References theType.
|
friend |
Definition at line 214 of file GblPoint.h.
|
private |
Jacobian to next scatterer (or last measurement)
Definition at line 226 of file GblPoint.h.
Referenced by addNextJacobian(), getDerivatives(), and printPoint().
|
private |
Point-to-point jacobian from previous point.
Definition at line 224 of file GblPoint.h.
Referenced by getP2pJacobian(), and printPoint().
|
private |
Jacobian to previous scatterer (or first measurement)
Definition at line 225 of file GblPoint.h.
Referenced by addPrevJacobian(), getDerivatives(), and printPoint().
|
private |
Dimension of scatterer (0: none, 2: thin, 4:thick)
Definition at line 227 of file GblPoint.h.
Referenced by addScatterer(), addThickScatterer(), getScatDim(), getScatterer(), getScatTransformation(), and printPoint().
|
private |
Scattering precision (diagonal of inverse covariance matrix)
Definition at line 232 of file GblPoint.h.
Referenced by addScatterer(), addThickScatterer(), getReducedScatterer(), getScatterer(), and printPoint().
|
private |
Scattering residuals (initial kinks if iterating)
Definition at line 231 of file GblPoint.h.
Referenced by addScatterer(), addThickScatterer(), getReducedScatterer(), getScatterer(), and printPoint().
|
private |
Transformation of diagonalization (of scat. precision matrix)
Definition at line 230 of file GblPoint.h.
Referenced by addScatterer(), addThickScatterer(), getReducedScatterer(), getScatterer(), and getScatTransformation().
|
private |
Label identifying point.
Definition at line 221 of file GblPoint.h.
Referenced by getLabel(), printPoint(), and setLabel().
|
private |
List of measurements at point.
Definition at line 233 of file GblPoint.h.
Referenced by addGlobals(), addLocals(), GblPoint(), getGlobalLabelsAndDerivatives(), getMeasBegin(), getMeasEnd(), numMeasurements(), and printPoint().
|
private |
Offset number at point if not negative (else interpolation needed)
Definition at line 222 of file GblPoint.h.
Referenced by getOffset(), printPoint(), and setOffset().
|
private |
Type (-1: first, 0: inner, 1: last)
Definition at line 223 of file GblPoint.h.