Analysis Software
Documentation for sPHENIX simulation software
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
genfit::MplTrackRep Class Reference

#include <GenFit/blob/master/trackReps/include/MplTrackRep.h>

+ Inheritance diagram for genfit::MplTrackRep:
+ Collaboration diagram for genfit::MplTrackRep:

Public Member Functions

 MplTrackRep ()
 
 MplTrackRep (int pdgCode, float magCharge, char propDir=0)
 
 ~MplTrackRep ()
 
double RKPropagate (M1x7 &state7, M7x7 *jacobian, M1x3 &SA, double S, bool varField=true, bool calcOnlyLastRowOfJ=false) const override
 The actual Runge Kutta propagation.
 
double getCharge (const StateOnPlane &state) const override
 Get the (fitted) charge of a state. This is not always equal the pdg charge (e.g. if the charge sign was flipped during the fit).
 
- Public Member Functions inherited from genfit::RKTrackRep
 RKTrackRep ()
 
 RKTrackRep (int pdgCode, char propDir=0)
 
virtual ~RKTrackRep ()
 
virtual AbsTrackRepclone () const override
 Clone the trackRep.
 
virtual double extrapolateToPlane (StateOnPlane &state, const SharedPlanePtr &plane, bool stopAtBoundary=false, bool calcJacobianNoise=false) const override
 Extrapolates the state to plane, and returns the extrapolation length and, via reference, the extrapolated state.
 
virtual double extrapolateToLine (StateOnPlane &state, const TVector3 &linePoint, const TVector3 &lineDirection, bool stopAtBoundary=false, bool calcJacobianNoise=false) const override
 Extrapolates the state to the POCA to a line, and returns the extrapolation length and, via reference, the extrapolated state.
 
virtual double extrapolateToPoint (StateOnPlane &state, const TVector3 &point, bool stopAtBoundary=false, bool calcJacobianNoise=false) const override
 Extrapolates the state to the POCA to a point, and returns the extrapolation length and, via reference, the extrapolated state.
 
virtual double extrapolateToPoint (StateOnPlane &state, const TVector3 &point, const TMatrixDSym &G, bool stopAtBoundary=false, bool calcJacobianNoise=false) const override
 Extrapolates the state to the POCA to a point in the metric of G, and returns the extrapolation length and, via reference, the extrapolated state.
 
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 override
 Extrapolates the state to the cylinder surface, and returns the extrapolation length and, via reference, the extrapolated state.
 
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 override
 Extrapolates the state to the cone surface, and returns the extrapolation length and, via reference, the extrapolated state.
 
virtual double extrapolateToSphere (StateOnPlane &state, double radius, const TVector3 &point=TVector3(0., 0., 0.), bool stopAtBoundary=false, bool calcJacobianNoise=false) const override
 Extrapolates the state to the sphere surface, and returns the extrapolation length and, via reference, the extrapolated state.
 
virtual double extrapolateBy (StateOnPlane &state, double step, bool stopAtBoundary=false, bool calcJacobianNoise=false) const override
 Extrapolates the state by step (cm) and returns the extrapolation length and, via reference, the extrapolated state.
 
unsigned int getDim () const override
 Get the dimension of the state vector used by the track representation.
 
virtual TVector3 getPos (const StateOnPlane &state) const override
 Get the cartesian position of a state.
 
virtual TVector3 getMom (const StateOnPlane &state) const override
 Get the cartesian momentum vector of a state.
 
virtual void getPosMom (const StateOnPlane &state, TVector3 &pos, TVector3 &mom) const override
 Get cartesian position and momentum vector of a state.
 
virtual double getMomMag (const StateOnPlane &state) const override
 get the magnitude of the momentum in GeV.
 
virtual double getMomVar (const MeasuredStateOnPlane &state) const override
 get the variance of the absolute value of the momentum .
 
virtual TMatrixDSym get6DCov (const MeasuredStateOnPlane &state) const override
 Get the 6D covariance.
 
virtual void getPosMomCov (const MeasuredStateOnPlane &state, TVector3 &pos, TVector3 &mom, TMatrixDSym &cov) const override
 Translates MeasuredStateOnPlane into 3D position, momentum and 6x6 covariance.
 
virtual double getQop (const StateOnPlane &state) const override
 Get charge over momentum.
 
double getSpu (const StateOnPlane &state) const
 
double getTime (const StateOnPlane &state) const override
 Get the time corresponding to the StateOnPlane. Extrapolation.
 
virtual void getForwardJacobianAndNoise (TMatrixD &jacobian, TMatrixDSym &noise, TVectorD &deltaState) const override
 Get the jacobian and noise matrix of the last extrapolation.
 
virtual void getBackwardJacobianAndNoise (TMatrixD &jacobian, TMatrixDSym &noise, TVectorD &deltaState) const override
 Get the jacobian and noise matrix of the last extrapolation if it would have been done in opposite direction.
 
std::vector< genfit::MatStepgetSteps () const override
 Get stepsizes and material properties of crossed materials of the last extrapolation.
 
virtual double getRadiationLenght () const override
 Get the accumulated X/X0 (path / radiation length) of the material crossed in the last extrapolation.
 
virtual void setPosMom (StateOnPlane &state, const TVector3 &pos, const TVector3 &mom) const override
 Set position and momentum of state.
 
virtual void setPosMom (StateOnPlane &state, const TVectorD &state6) const override
 Set position and momentum of state.
 
virtual void setPosMomErr (MeasuredStateOnPlane &state, const TVector3 &pos, const TVector3 &mom, const TVector3 &posErr, const TVector3 &momErr) const override
 Set position and momentum and error of state.
 
virtual void setPosMomCov (MeasuredStateOnPlane &state, const TVector3 &pos, const TVector3 &mom, const TMatrixDSym &cov6x6) const override
 Set position, momentum and covariance of state.
 
virtual void setPosMomCov (MeasuredStateOnPlane &state, const TVectorD &state6, const TMatrixDSym &cov6x6) const override
 Set position, momentum and covariance of state.
 
virtual void setChargeSign (StateOnPlane &state, double charge) const override
 Set the sign of the charge according to charge.
 
virtual void setQop (StateOnPlane &state, double qop) const override
 Set charge/momentum.
 
void setSpu (StateOnPlane &state, double spu) const
 
void setTime (StateOnPlane &state, double time) const override
 Set time at which the state was defined.
 
virtual bool isSameType (const AbsTrackRep *other) override
 check if other is of same type (e.g. RKTrackRep).
 
virtual bool isSame (const AbsTrackRep *other) override
 check if other is of same type (e.g. RKTrackRep) and has same pdg code.
 
- Public Member Functions inherited from genfit::AbsTrackRep
 AbsTrackRep ()
 
 AbsTrackRep (int pdgCode, char propDir=0)
 
virtual ~AbsTrackRep ()
 
virtual double extrapolateToLine (StateOnPlane &state, const TVector3 &point1, const TVector3 &point2, TVector3 &poca, TVector3 &dirInPoca, TVector3 &poca_onwire, bool stopAtBoundary=false, bool calcJacobianNoise=false) const
 Resembles the interface of GFAbsTrackRep in old versions of genfit.
 
double extrapolateToMeasurement (StateOnPlane &state, const AbsMeasurement *measurement, bool stopAtBoundary=false, bool calcJacobianNoise=false) const
 extrapolate to an AbsMeasurement
 
TVector3 getDir (const StateOnPlane &state) const
 Get the direction vector of a state.
 
void getPosDir (const StateOnPlane &state, TVector3 &pos, TVector3 &dir) const
 Get cartesian position and direction vector of a state.
 
virtual TVectorD get6DState (const StateOnPlane &state) const
 Get the 6D state vector (x, y, z, p_x, p_y, p_z).
 
virtual void get6DStateCov (const MeasuredStateOnPlane &state, TVectorD &stateVec, TMatrixDSym &cov) const
 Translates MeasuredStateOnPlane into 6D state vector (x, y, z, p_x, p_y, p_z) and 6x6 covariance.
 
int getPDG () const
 Get the pdg code.
 
double getPDGCharge () const
 Get the charge of the particle of the pdg code.
 
double getMass (const StateOnPlane &state) const
 Get tha particle mass in GeV/c^2.
 
char getPropDir () const
 Get propagation direction. (-1, 0, 1) -> (backward, auto, forward).
 
void calcJacobianNumerically (const genfit::StateOnPlane &origState, const genfit::SharedPlanePtr destPlane, TMatrixD &jacobian) const
 Calculate Jacobian of transportation numerically. Slow but accurate. Can be used to validate (semi)analytic calculations.
 
bool switchPDGSign ()
 try to multiply pdg code with -1. (Switch from particle to anti-particle and vice versa).
 
void setPropDir (int dir)
 Set propagation direction. (-1, 0, 1) -> (backward, auto, forward).
 
void switchPropDir ()
 Switch propagation direction. Has no effect if propDir_ is set to 0.
 
virtual void setDebugLvl (unsigned int lvl=1)
 
virtual void Print (const Option_t *="") const
 

Private Attributes

const double m_magCharge
 
const double m_mass
 

Additional Inherited Members

- Protected Member Functions inherited from genfit::AbsTrackRep
 AbsTrackRep (const AbsTrackRep &)
 protect from calling copy c'tor from outside the class. Use clone() if you want a copy!
 
AbsTrackRepoperator= (const AbsTrackRep &)
 protect from calling assignment operator from outside the class. Use clone() instead!
 
- Protected Attributes inherited from genfit::RKTrackRep
StateOnPlane lastStartState_
 
StateOnPlane lastEndState_
 state where the last extrapolation has started
 

Detailed Description

Monopole track representation. It is a minimal modification of RKTrackRep with a different equations of motion for magnetic charges.

In the current implementation the states on plane are 5-d: u, v, u', v', q/p except that q in this case is magnetic, and the monopole has no electic charge.

Definition at line 33 of file MplTrackRep.h.

View newest version in sPHENIX GitHub at line 33 of file MplTrackRep.h

Constructor & Destructor Documentation

genfit::MplTrackRep::MplTrackRep ( )
inline

Definition at line 37 of file MplTrackRep.h.

View newest version in sPHENIX GitHub at line 37 of file MplTrackRep.h

MplTrackRep::MplTrackRep ( int  pdgCode,
float  magCharge,
char  propDir = 0 
)

Definition at line 33 of file MplTrackRep.cc.

View newest version in sPHENIX GitHub at line 33 of file MplTrackRep.cc

MplTrackRep::~MplTrackRep ( )

Definition at line 40 of file MplTrackRep.cc.

View newest version in sPHENIX GitHub at line 40 of file MplTrackRep.cc

Member Function Documentation

double MplTrackRep::getCharge ( const StateOnPlane state) const
overridevirtual

Get the (fitted) charge of a state. This is not always equal the pdg charge (e.g. if the charge sign was flipped during the fit).

Reimplemented from genfit::RKTrackRep.

Definition at line 44 of file MplTrackRep.cc.

View newest version in sPHENIX GitHub at line 44 of file MplTrackRep.cc

References genfit::AbsTrackRep::getPDGCharge(), genfit::StateOnPlane::getState(), m_magCharge, and genfit::Exception::setFatal().

+ Here is the call graph for this function:

double MplTrackRep::RKPropagate ( M1x7 state7,
M7x7 jacobian,
M1x3 SA,
double  S,
bool  varField = true,
bool  calcOnlyLastRowOfJ = false 
) const
overridevirtual

The actual Runge Kutta propagation.

propagate state7 with step S. Fills SA (Start directions derivatives dA/S). This is a single Runge-Kutta step. If jacobian is nullptr, only the state is propagated, otherwise also the 7x7 jacobian is calculated. If varField is false, the magnetic field will only be evaluated at the starting position. The return value is an estimation on how good the extrapolation is, and it is usually fine if it is > 1. It gives a suggestion how you must scale S so that the quality will be sufficient.

Reimplemented from genfit::RKTrackRep.

Definition at line 61 of file MplTrackRep.cc.

View newest version in sPHENIX GitHub at line 61 of file MplTrackRep.cc

References A, D1, D2, D4, genfit::AbsTrackRep::debugLvl_, genfit::debugOut, Acts::UnitConstants::e, F1, F2, F3, genfit::FieldManager::getFieldVal(), genfit::FieldManager::getInstance(), i, Acts::UnitConstants::J, m_magCharge, m_mass, Acts::IntegrationTest::R, physmon_track_finding_ttbar::r, S(), sign(), and start.

+ Here is the call graph for this function:

Member Data Documentation

const double genfit::MplTrackRep::m_magCharge
private

Definition at line 54 of file MplTrackRep.h.

View newest version in sPHENIX GitHub at line 54 of file MplTrackRep.h

Referenced by getCharge(), and RKPropagate().

const double genfit::MplTrackRep::m_mass
private

Definition at line 55 of file MplTrackRep.h.

View newest version in sPHENIX GitHub at line 55 of file MplTrackRep.h

Referenced by RKPropagate().


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