24 #ifndef genfit_RKTrackRep_h
25 #define genfit_RKTrackRep_h
45 std::fill(state7_.
begin(), state7_.
end(), 0);
58 std::fill(jac7_.
begin(), jac7_.
end(), 0);
59 std::fill(noise7_.
begin(), jac7_.
end(), 0);
85 bool stopAtBoundary =
false,
86 bool calcJacobianNoise =
false)
const;
91 const TVector3& linePoint,
92 const TVector3& lineDirection,
93 bool stopAtBoundary =
false,
94 bool calcJacobianNoise =
false)
const;
97 const TVector3& point,
98 bool stopAtBoundary =
false,
99 bool calcJacobianNoise =
false)
const {
100 return extrapToPoint(state, point, NULL, stopAtBoundary, calcJacobianNoise);
104 const TVector3& point,
105 const TMatrixDSym& G,
106 bool stopAtBoundary =
false,
107 bool calcJacobianNoise =
false)
const {
108 return extrapToPoint(state, point, &G, stopAtBoundary, calcJacobianNoise);
113 const TVector3& linePoint = TVector3(0.,0.,0.),
114 const TVector3& lineDirection = TVector3(0.,0.,1.),
115 bool stopAtBoundary =
false,
116 bool calcJacobianNoise =
false)
const;
121 const TVector3& linePoint = TVector3(0.,0.,0.),
122 const TVector3& lineDirection = TVector3(0.,0.,1.),
123 bool stopAtBoundary =
false,
124 bool calcJacobianNoise =
false)
const;
128 const TVector3& point = TVector3(0.,0.,0.),
129 bool stopAtBoundary =
false,
130 bool calcJacobianNoise =
false)
const;
134 bool stopAtBoundary =
false,
135 bool calcJacobianNoise =
false)
const;
159 std::vector<genfit::MatStep>
getSteps()
const;
188 bool varField =
true,
189 bool calcOnlyLastRowOfJ =
false)
const;
199 const TVector3& point,
200 const TMatrixDSym* G = NULL,
201 bool stopAtBoundary =
false,
202 bool calcJacobianNoise =
false)
const;
210 void calcJ_pM_5x7(
M5x7& J_pM,
const TVector3& U,
const TVector3& V,
const M1x3& pTilde,
double spu)
const;
219 void calcJ_Mp_7x5(
M7x5& J_Mp,
const TVector3& U,
const TVector3& V,
const TVector3& W,
const M1x3& A)
const;
222 const M1x7& destState7,
const DetPlane& destPlane)
const;
244 M1x7* J_MMT_unprojected_lastRow,
245 double& coveredDistance,
248 M7x7& noiseProjection,
250 bool onlyOneStep =
false,
251 bool calcOnlyLastRowOfJ =
false)
const;
256 const double& charge,
260 TVector3
pocaOnLine(
const TVector3& linePoint,
261 const TVector3& lineDirection,
262 const TVector3& point)
const;
280 bool fillExtrapSteps,
281 TMatrixDSym* cov =
nullptr,
282 bool onlyOneStep =
false,
283 bool stopAtBoundary =
false,
284 double maxStep = 1.E99)
const;
320 #endif // genfit_RKTrackRep_h
virtual double extrapolateToPoint(StateOnPlane &state, const TVector3 &point, const TMatrixDSym &G, bool stopAtBoundary=false, bool calcJacobianNoise=false) const
Extrapolates the state to the POCA to a point in the metric of G, and returns the extrapolation lengt...
void transformPM7(const MeasuredStateOnPlane &state, M7x7 &out7x7) const
void calcForwardJacobianAndNoise(const M1x7 &startState7, const DetPlane &startPlane, const M1x7 &destState7, const DetPlane &destPlane) const
Helper to store different limits on the stepsize for the RKTRackRep.
virtual TVector3 getPos(const StateOnPlane &state) const
Get the cartesian position of a state.
M7x7 noiseProjection_
noise matrix of the last extrapolation
TVector3 pocaOnLine(const TVector3 &linePoint, const TVector3 &lineDirection, const TVector3 &point) const
virtual double getQop(const StateOnPlane &state) const
Get charge over momentum.
boost::shared_ptr< genfit::DetPlane > SharedPlanePtr
Shared Pointer to a DetPlane.
virtual void getPosMomCov(const MeasuredStateOnPlane &state, TVector3 &pos, TVector3 &mom, TMatrixDSym &cov) const
Translates MeasuredStateOnPlane into 3D position, momentum and 6x6 covariance.
double momMag(const M1x7 &state7) const
virtual double extrapolateToPlane(StateOnPlane &state, const SharedPlanePtr &plane, bool stopAtBoundary=false, bool calcJacobianNoise=false) const
Extrapolates the state to plane, and returns the extrapolation length and, via reference, the extrapolated state.
double RKPropagate(M1x7 &state7, M7x7 *jacobian, M1x3 &SA, double S, bool varField=true, bool calcOnlyLastRowOfJ=false) const
The actual Runge Kutta propagation.
bool RKutta(const M1x4 &SU, const DetPlane &plane, double charge, double mass, M1x7 &state7, M7x7 *jacobianT, M1x7 *J_MMT_unprojected_lastRow, double &coveredDistance, double &flightTime, bool &checkJacProj, M7x7 &noiseProjection, StepLimits &limits, bool onlyOneStep=false, bool calcOnlyLastRowOfJ=false) const
Propagates the particle through the magnetic field.
virtual double extrapolateBy(StateOnPlane &state, double step, bool stopAtBoundary=false, bool calcJacobianNoise=false) const
Extrapolates the state by step (cm) and returns the extrapolation length and, via reference...
void checkCache(const StateOnPlane &state, const SharedPlanePtr *plane) const
const TVectorD & getState() const
virtual double extrapolateToSphere(StateOnPlane &state, double radius, const TVector3 &point=TVector3(0., 0., 0.), bool stopAtBoundary=false, bool calcJacobianNoise=false) const
Extrapolates the state to the sphere surface, and returns the extrapolation length and...
virtual void setPosMom(StateOnPlane &state, const TVector3 &pos, const TVector3 &mom) const
Set position and momentum of state.
void transformM6P(const M6x6 &in6x6, const M1x7 &state7, MeasuredStateOnPlane &state) const
void calcJ_pM_5x7(M5x7 &J_pM, const TVector3 &U, const TVector3 &V, const M1x3 &pTilde, double spu) const
int RKStepsFXStart_
RungeKutta steps made in the last extrapolation.
virtual void setQop(StateOnPlane &state, double qop) const
Set charge/momentum.
virtual double getRadiationLenght() const
Get the accumulated X/X0 (path / radiation length) of the material crossed in the last extrapolation...
TMatrixD fJacobian_
steps made in Extrap during last extrapolation
virtual double extrapolateToPoint(StateOnPlane &state, const TVector3 &point, bool stopAtBoundary=false, bool calcJacobianNoise=false) const
Extrapolates the state to the POCA to a point, and returns the extrapolation length and...
virtual void setPosMomErr(MeasuredStateOnPlane &state, const TVector3 &pos, const TVector3 &mom, const TVector3 &posErr, const TVector3 &momErr) const
Set position and momentum and error of state.
virtual AbsTrackRep * clone() const
Clone the trackRep.
void setSpu(StateOnPlane &state, double spu) const
void setTime(StateOnPlane &state, double time) const
Set time at which the state was defined.
StateOnPlane lastEndState_
state where the last extrapolation has started
virtual double extrapolateToLine(StateOnPlane &state, const TVector3 &linePoint, const TVector3 &lineDirection, bool stopAtBoundary=false, bool calcJacobianNoise=false) const
Extrapolates the state to the POCA to a line, and returns the extrapolation length and...
virtual bool isSameType(const AbsTrackRep *other)
check if other is of same type (e.g. RKTrackRep).
StateOnPlane with additional covariance matrix.
Simple struct containing MaterialProperties and stepsize in the material.
virtual double getCharge(const StateOnPlane &state) const
Get the (fitted) charge of a state. This is not always equal the pdg charge (e.g. if the charge sign ...
virtual TVector3 getMom(const StateOnPlane &state) const
Get the cartesian momentum vector of a state.
virtual TMatrixDSym get6DCov(const MeasuredStateOnPlane &state) const
Get the 6D covariance.
unsigned int cachePos_
use cached RKSteps_ for extrapolation
virtual void getPosMom(const StateOnPlane &state, TVector3 &pos, TVector3 &mom) const
Get cartesian position and momentum vector of a state.
virtual void getBackwardJacobianAndNoise(TMatrixD &jacobian, TMatrixDSym &noise, TVectorD &deltaState) const
Get the jacobian and noise matrix of the last extrapolation if it would have been done in opposite di...
Abstract base class for a track representation.
virtual void setPosMomCov(MeasuredStateOnPlane &state, const TVector3 &pos, const TVector3 &mom, const TMatrixDSym &cov6x6) const
Set position, momentum and covariance of state.
double estimateStep(const M1x7 &state7, const M1x4 &SU, const DetPlane &plane, const double &charge, double &relMomLoss, StepLimits &limits) const
double getSpu(const StateOnPlane &state) const
virtual double extrapolateToCone(StateOnPlane &state, double radius, const TVector3 &linePoint=TVector3(0., 0., 0.), const TVector3 &lineDirection=TVector3(0., 0., 1.), bool stopAtBoundary=false, bool calcJacobianNoise=false) const
Extrapolates the state to the cone surface, and returns the extrapolation length and, via reference, the extrapolated state.
void getState7(const StateOnPlane &state, M1x7 &state7) const
double getTime(const StateOnPlane &state) const
Get the time corresponding to the StateOnPlane. Extrapolation.
void calcJ_Mp_7x5(M7x5 &J_Mp, const TVector3 &U, const TVector3 &V, const TVector3 &W, const M1x3 &A) const
std::vector< genfit::MatStep > getSteps() const
Get stepsizes and material properties of crossed materials of the last extrapolation.
StateOnPlane lastStartState_
void getState5(StateOnPlane &state, const M1x7 &state7) const
unsigned int getDim() const
Get the dimension of the state vector used by the track representation.
AbsTrackRep with 5D track parameterization in plane coordinates: (q/p, u', v', u, v) ...
double Extrap(const DetPlane &startPlane, const DetPlane &destPlane, double charge, double mass, bool &isAtBoundary, M1x7 &state7, double &flightTime, bool fillExtrapSteps, TMatrixDSym *cov=nullptr, bool onlyOneStep=false, bool stopAtBoundary=false, double maxStep=1.E99) const
Handles propagation and material effects.
virtual double getMomVar(const MeasuredStateOnPlane &state) const
get the variance of the absolute value of the momentum .
virtual bool isSame(const AbsTrackRep *other)
check if other is of same type (e.g. RKTrackRep) and has same pdg code.
A state with arbitrary dimension defined in a DetPlane.
virtual double getMomMag(const StateOnPlane &state) const
get the magnitude of the momentum in GeV.
std::vector< RKStep > RKSteps_
state where the last extrapolation has ended
virtual double extrapolateToLine(StateOnPlane &state, const TVector3 &linePoint, const TVector3 &lineDirection, bool stopAtBoundary=false, bool calcJacobianNoise=false) const =0
Extrapolates the state to the POCA to a line, and returns the extrapolation length and...
virtual void getForwardJacobianAndNoise(TMatrixD &jacobian, TMatrixDSym &noise, TVectorD &deltaState) const
Get the jacobian and noise matrix of the last extrapolation.
virtual double extrapToPoint(StateOnPlane &state, const TVector3 &point, const TMatrixDSym *G=NULL, bool stopAtBoundary=false, bool calcJacobianNoise=false) const
virtual double extrapolateToCylinder(StateOnPlane &state, double radius, const TVector3 &linePoint=TVector3(0., 0., 0.), const TVector3 &lineDirection=TVector3(0., 0., 1.), bool stopAtBoundary=false, bool calcJacobianNoise=false) const
Extrapolates the state to the cylinder surface, and returns the extrapolation length and...
void transformPM6(const MeasuredStateOnPlane &state, M6x6 &out6x6) const
std::vector< ExtrapStep > ExtrapSteps_
void transformM7P(const M7x7 &in7x7, const M1x7 &state7, MeasuredStateOnPlane &state) const
virtual void setChargeSign(StateOnPlane &state, double charge) const
Set the sign of the charge according to charge.
Defines for I/O streams used for error and debug printing.