![]() |
GeneralBrokenLines
Rev: 2.2.0
|
Point on trajectory. More...
#include <GblPoint.h>
Public Member Functions | |
| GblPoint (const Matrix5d &aJacobian) | |
| 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... | |
| template<typename Precision , typename std::enable_if<(Precision::ColsAtCompileTime==2)>::type * = nullptr> | |
| void | addScatterer (const Eigen::Vector2d &aResiduals, const Eigen::MatrixBase< Precision > &aPrecision) |
| Add a (thin) scatterer to a point. More... | |
| template<typename Precision , typename std::enable_if<(Precision::ColsAtCompileTime==1)>::type * = nullptr> | |
| void | addScatterer (const Eigen::Vector2d &aResiduals, const Eigen::MatrixBase< Precision > &aPrecision) |
| Add a (thin) 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 | hasMeasurement () const |
| Check for measurement at a point. More... | |
| double | getMeasPrecMin () const |
| get precision cutoff. More... | |
| void | getMeasurement (Matrix5d &aProjection, Vector5d &aResiduals, Vector5d &aPrecision) const |
| Retrieve measurement of a point. More... | |
| void | getMeasTransformation (Eigen::MatrixXd &aTransformation) const |
| Get measurement transformation (from diagonalization). More... | |
| bool | hasScatterer () const |
| Check for scatterer at a point. More... | |
| void | getScatterer (Eigen::Matrix2d &aTransformation, Eigen::Vector2d &aResiduals, Eigen::Vector2d &aPrecision) const |
| Retrieve scatterer of a point. More... | |
| void | getScatTransformation (Eigen::Matrix2d &aTransformation) const |
| Get scatterer transformation (from diagonalization). More... | |
| unsigned int | getNumLocals () const |
| Retrieve number of local derivatives from a point. More... | |
| const Eigen::MatrixXd & | getLocalDerivatives () const |
| Retrieve local derivatives from a point. More... | |
| unsigned int | getNumGlobals () const |
| Retrieve number of global derivatives from a point. More... | |
| void | getGlobalLabels (std::vector< int > &aLabels) const |
| Retrieve global derivatives labels from a point. More... | |
| void | getGlobalDerivatives (Eigen::MatrixXd &aDerivatives) const |
| Retrieve global derivatives from a point. More... | |
| void | getGlobalLabelsAndDerivatives (unsigned int aRow, std::vector< int > &aLabels, std::vector< double > &aDerivatives) const |
| Retrieve global derivatives from a point for a single row. 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... | |
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 | 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... | |
| Matrix5d | p2pJacobian |
| Point-to-point jacobian from previous point. More... | |
| Matrix5d | prevJacobian |
| Jacobian to previous scatterer (or first measurement) More... | |
| Matrix5d | nextJacobian |
| Jacobian to next scatterer (or last measurement) More... | |
| unsigned int | measDim |
| Dimension of measurement (1-5), 0 indicates absence of measurement. More... | |
| double | measPrecMin |
| Minimal measurement precision (for usage) More... | |
| Matrix5d | measProjection |
| Projection from measurement to local system. More... | |
| Vector5d | measResiduals |
| Measurement residuals. More... | |
| Vector5d | measPrecision |
| Measurement precision (diagonal of inverse covariance matrix) More... | |
| bool | transFlag |
| Transformation exists? More... | |
| Eigen::Matrix< double, Eigen::Dynamic, Eigen::Dynamic, Eigen::ColMajor, 5, 5 > | measTransformation |
| Transformation of diagonalization (of meas. precision matrix) More... | |
| bool | scatFlag |
| Scatterer present? More... | |
| Eigen::Matrix2d | scatTransformation |
| Transformation of diagonalization (of scat. precision matrix) More... | |
| Eigen::Vector2d | scatResiduals |
| Scattering residuals (initial kinks if iterating) More... | |
| Eigen::Vector2d | scatPrecision |
| Scattering precision (diagonal of inverse covariance matrix) More... | |
| Eigen::MatrixXd | localDerivatives |
| Derivatives of measurement vs additional local (fit) parameters. More... | |
| std::vector< int > | globalLabels |
| Labels of global (MP-II) derivatives. More... | |
| Eigen::MatrixXd | globalDerivatives |
| Derivatives of measurement vs additional global (MP-II) parameters. 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 69 of file GblPoint.h.
| gbl::GblPoint::GblPoint | ( | const Matrix5d & | aJacobian | ) |
Create a point.
Create point on (initial) trajectory. Needs transformation jacobian from previous point.
| [in] | aJacobian | Transformation jacobian from previous point |
Definition at line 41 of file GblPoint.cpp.
|
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 403 of file GblPoint.h.
References globalDerivatives, globalLabels, measDim, measTransformation, and transFlag.
| 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 391 of file GblPoint.h.
References localDerivatives, measDim, measTransformation, and transFlag.
| 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 290 of file GblPoint.h.
| 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 332 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 460 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 445 of file GblPoint.cpp.
References prevJacobian.
| void gbl::GblPoint::addScatterer | ( | const Eigen::Vector2d & | aResiduals, |
| const Eigen::MatrixBase< Precision > & | 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 |
| Precision | Precision matrix or vector (with diagonal) |
| [in] | aResiduals | Scatterer residuals |
| [in] | aPrecision | Scatterer precision (full matrix) |
Definition at line 366 of file GblPoint.h.
| void gbl::GblPoint::addScatterer | ( | const Eigen::Vector2d & | aResiduals, |
| const Eigen::MatrixBase< Precision > & | 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 |
| Precision | Precision matrix or vector (with diagonal) |
| [in] | aResiduals | Scatterer residuals |
| [in] | aPrecision | Scatterer precision (vector with diagonal) |
| 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 474 of file GblPoint.cpp.
References nextJacobian, and prevJacobian.
Referenced by gbl::GblTrajectory::getFitToKinkJacobian(), and gbl::GblTrajectory::getFitToLocalJacobian().
| void gbl::GblPoint::getGlobalDerivatives | ( | Eigen::MatrixXd & | aDerivatives | ) | const |
Retrieve global derivatives from a point.
| [out] | aDerivatives | Global derivatives |
Definition at line 390 of file GblPoint.cpp.
References globalDerivatives.
| void gbl::GblPoint::getGlobalLabels | ( | std::vector< int > & | aLabels | ) | const |
Retrieve global derivatives labels from a point.
| [out] | aLabels | Global labels |
Definition at line 382 of file GblPoint.cpp.
References globalLabels.
| void gbl::GblPoint::getGlobalLabelsAndDerivatives | ( | unsigned int | aRow, |
| std::vector< int > & | aLabels, | ||
| std::vector< double > & | aDerivatives | ||
| ) | const |
Retrieve global derivatives from a point for a single row.
| [in] | aRow | Row number |
| [out] | aLabels | Global labels |
| [out] | aDerivatives | Global derivatives |
Definition at line 400 of file GblPoint.cpp.
References globalDerivatives, and globalLabels.
| unsigned int gbl::GblPoint::getLabel | ( | ) | const |
| const MatrixXd & gbl::GblPoint::getLocalDerivatives | ( | ) | const |
Retrieve local derivatives from a point.
Definition at line 340 of file GblPoint.cpp.
References localDerivatives.
| double gbl::GblPoint::getMeasPrecMin | ( | ) | const |
get precision cutoff.
Definition at line 194 of file GblPoint.cpp.
References measPrecMin.
| void gbl::GblPoint::getMeasTransformation | ( | Eigen::MatrixXd & | aTransformation | ) | const |
Get measurement transformation (from diagonalization).
| [out] | aTransformation | Transformation matrix |
Definition at line 216 of file GblPoint.cpp.
References measDim, measTransformation, and transFlag.
| void gbl::GblPoint::getMeasurement | ( | Matrix5d & | aProjection, |
| Vector5d & | aResiduals, | ||
| Vector5d & | aPrecision | ||
| ) | const |
Retrieve measurement of a point.
| [out] | aProjection | Projection from (diagonalized) measurement to local system |
| [out] | aResiduals | Measurement residuals |
| [out] | aPrecision | Measurement precision (diagonal) |
Definition at line 204 of file GblPoint.cpp.
References measDim, measPrecision, measProjection, and measResiduals.
| unsigned int gbl::GblPoint::getNumGlobals | ( | ) | const |
Retrieve number of global derivatives from a point.
Definition at line 374 of file GblPoint.cpp.
References globalDerivatives.
| unsigned int gbl::GblPoint::getNumLocals | ( | ) | const |
Retrieve number of local derivatives from a point.
Definition at line 335 of file GblPoint.cpp.
References localDerivatives.
| int gbl::GblPoint::getOffset | ( | ) | const |
Retrieve offset for point.
Definition at line 432 of file GblPoint.cpp.
References theOffset.
Referenced by gbl::GblTrajectory::getFitToKinkJacobian(), and gbl::GblTrajectory::getFitToLocalJacobian().
| const Matrix5d & gbl::GblPoint::getP2pJacobian | ( | ) | const |
Retrieve point-to-(previous)point jacobian.
Definition at line 437 of file GblPoint.cpp.
References p2pJacobian.
Referenced by gbl::GblTrajectory::calcJacobians().
| void gbl::GblPoint::getScatterer | ( | Eigen::Matrix2d & | aTransformation, |
| Eigen::Vector2d & | aResiduals, | ||
| Eigen::Vector2d & | 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 290 of file GblPoint.cpp.
References scatPrecision, scatResiduals, and scatTransformation.
| void gbl::GblPoint::getScatTransformation | ( | Eigen::Matrix2d & | aTransformation | ) | const |
Get scatterer transformation (from diagonalization).
| [out] | aTransformation | Transformation matrix |
Definition at line 301 of file GblPoint.cpp.
References scatFlag, and scatTransformation.
| unsigned int gbl::GblPoint::hasMeasurement | ( | ) | const |
Check for measurement at a point.
Get dimension of measurement (0 = none).
Definition at line 186 of file GblPoint.cpp.
References measDim.
| bool gbl::GblPoint::hasScatterer | ( | ) | const |
| void gbl::GblPoint::printPoint | ( | unsigned int | level = 0 | ) | const |
Print GblPoint.
| [in] | level | print level (0: minimum, >0: more) |
Definition at line 506 of file GblPoint.cpp.
References globalDerivatives, globalLabels, localDerivatives, measDim, measPrecision, measPrecMin, measProjection, measResiduals, nextJacobian, p2pJacobian, prevJacobian, scatFlag, scatPrecision, scatResiduals, theLabel, theOffset, and transFlag.
|
private |
Define label of point (by GBLTrajectory constructor)
| [in] | aLabel | Label identifying point |
Definition at line 414 of file GblPoint.cpp.
References theLabel.
|
private |
Define offset for point (by GBLTrajectory constructor)
| [in] | anOffset | Offset number |
Definition at line 427 of file GblPoint.cpp.
References theOffset.
|
friend |
Definition at line 259 of file GblPoint.h.
|
private |
Derivatives of measurement vs additional global (MP-II) parameters.
Definition at line 285 of file GblPoint.h.
Referenced by addGlobals(), getGlobalDerivatives(), getGlobalLabelsAndDerivatives(), getNumGlobals(), and printPoint().
|
private |
Labels of global (MP-II) derivatives.
Definition at line 284 of file GblPoint.h.
Referenced by addGlobals(), getGlobalLabels(), getGlobalLabelsAndDerivatives(), and printPoint().
|
private |
Derivatives of measurement vs additional local (fit) parameters.
Definition at line 283 of file GblPoint.h.
Referenced by addLocals(), getLocalDerivatives(), getNumLocals(), and printPoint().
|
private |
Dimension of measurement (1-5), 0 indicates absence of measurement.
Definition at line 270 of file GblPoint.h.
Referenced by addGlobals(), addLocals(), getMeasTransformation(), getMeasurement(), hasMeasurement(), and printPoint().
|
private |
Measurement precision (diagonal of inverse covariance matrix)
Definition at line 275 of file GblPoint.h.
Referenced by getMeasurement(), and printPoint().
|
private |
Minimal measurement precision (for usage)
Definition at line 271 of file GblPoint.h.
Referenced by getMeasPrecMin(), and printPoint().
|
private |
Projection from measurement to local system.
Definition at line 272 of file GblPoint.h.
Referenced by getMeasurement(), and printPoint().
|
private |
Measurement residuals.
Definition at line 274 of file GblPoint.h.
Referenced by getMeasurement(), and printPoint().
|
private |
Transformation of diagonalization (of meas. precision matrix)
Definition at line 278 of file GblPoint.h.
Referenced by addGlobals(), addLocals(), and getMeasTransformation().
|
private |
Jacobian to next scatterer (or last measurement)
Definition at line 269 of file GblPoint.h.
Referenced by addNextJacobian(), getDerivatives(), and printPoint().
|
private |
Point-to-point jacobian from previous point.
Definition at line 267 of file GblPoint.h.
Referenced by getP2pJacobian(), and printPoint().
|
private |
Jacobian to previous scatterer (or first measurement)
Definition at line 268 of file GblPoint.h.
Referenced by addPrevJacobian(), getDerivatives(), and printPoint().
|
private |
Scatterer present?
Definition at line 279 of file GblPoint.h.
Referenced by getScatTransformation(), hasScatterer(), and printPoint().
|
private |
Scattering precision (diagonal of inverse covariance matrix)
Definition at line 282 of file GblPoint.h.
Referenced by getScatterer(), and printPoint().
|
private |
Scattering residuals (initial kinks if iterating)
Definition at line 281 of file GblPoint.h.
Referenced by getScatterer(), and printPoint().
|
private |
Transformation of diagonalization (of scat. precision matrix)
Definition at line 280 of file GblPoint.h.
Referenced by getScatterer(), and getScatTransformation().
|
private |
Label identifying point.
Definition at line 265 of file GblPoint.h.
Referenced by getLabel(), printPoint(), and setLabel().
|
private |
Offset number at point if not negative (else interpolation needed)
Definition at line 266 of file GblPoint.h.
Referenced by getOffset(), printPoint(), and setOffset().
|
private |
Transformation exists?
Definition at line 276 of file GblPoint.h.
Referenced by addGlobals(), addLocals(), getMeasTransformation(), and printPoint().