ROOT logo
HYDRA - THE HADES ANALYSIS PACKAGE » (UNKNOWN) » HKalTrackState

class HKalTrackState: public TObject

_HADES_CLASS_DESCRIPTION


 The state of the particle track at a given site k can be described by a vector
 a_k containing a set of track parameters:
 1. x0: x position
 2. y0: y position
 3. tx: p_x / p_z
 4. ty: p_y / p_z
 5. qp: particle charage in elementary charges divided by
        magnitude of momentum in MeV/c

 The class also stores the matrices needed for the Kalman filter:
 1. propagator: F_k = @f_k(a_k) / @a_k
 2. projector : H_k = @h_k(a_k) / @a^(k-1)_k
 3. covariance
 4. process noise
 where f_k(a_k) is the propagator function that propagates the track to
 the next site and h_k(a_k) projects the track state onto a measurement vector.



Function Members (Methods)

public:
HKalTrackState()
HKalTrackState(const HKalTrackState&)
HKalTrackState(Kalman::kalFilterTypes stateType, Int_t measDim, Int_t stateDim)
HKalTrackState(Kalman::kalFilterTypes stateType, const TVectorD& sv, Int_t measDim)
virtual~HKalTrackState()
voidTObject::AbstractMethod(const char* method) const
virtual voidTObject::AppendPad(Option_t* option = "")
virtual voidTObject::Browse(TBrowser* b)
virtual voidcalcDir(TVector3& dir) const
static voidcalcDir(TVector3& dir, const TVectorD& sv)
virtual voidcalcMomVec(TVector3& dir) const
static voidcalcMomVec(TVector3& dir, const TVectorD& sv)
virtual Bool_tcalcPosAndDirAtPlane(TVector3& pos, TVector3& dir, const HKalPlane& plane) const
virtual Bool_tcalcPosAtPlane(TVector3& pos, const HKalPlane& plane) const
static Bool_tcalcPosAtPlane(TVector3& pos, const HKalPlane& plane, const TVectorD& sv)
static voidcalcStateVec(TVectorD& sv, Double_t qp, const TVector3& pos, const TVector3& dir)
static TClass*Class()
virtual const char*TObject::ClassName() const
virtual voidclear()
virtual voidTObject::Clear(Option_t* = "")
virtual TObject*TObject::Clone(const char* newname = "") const
virtual Int_tTObject::Compare(const TObject* obj) const
virtual voidTObject::Copy(TObject& object) const
virtual voidTObject::Delete(Option_t* option = "")MENU
virtual Int_tTObject::DistancetoPrimitive(Int_t px, Int_t py)
virtual voidTObject::Draw(Option_t* option = "")
virtual voidTObject::DrawClass() constMENU
virtual TObject*TObject::DrawClone(Option_t* option = "") constMENU
virtual voidTObject::Dump() constMENU
virtual voidTObject::Error(const char* method, const char* msgfmt) const
virtual voidTObject::Execute(const char* method, const char* params, Int_t* error = 0)
virtual voidTObject::Execute(TMethod* method, TObjArray* params, Int_t* error = 0)
virtual voidTObject::ExecuteEvent(Int_t event, Int_t px, Int_t py)
virtual voidTObject::Fatal(const char* method, const char* msgfmt) const
virtual TObject*TObject::FindObject(const char* name) const
virtual TObject*TObject::FindObject(const TObject* obj) const
virtual const TMatrixD&getCovMat() const
virtual Option_t*TObject::GetDrawOption() const
static Long_tTObject::GetDtorOnly()
virtual const char*TObject::GetIconName() const
virtual Int_tgetMeasDim() const
virtual const char*TObject::GetName() const
virtual char*TObject::GetObjectInfo(Int_t px, Int_t py) const
static Bool_tTObject::GetObjectStat()
virtual Option_t*TObject::GetOption() const
virtual const TMatrixD&getProcNoiseMat() const
virtual const TMatrixD&getProjMat() const
virtual const TMatrixD&getPropMat() const
virtual Int_tgetStateDim() const
virtual Double_tgetStateParam(Kalman::kalStateIdx par) const
virtual Int_tgetStateType() const
virtual const TVectorD&getStateVec() const
virtual const char*TObject::GetTitle() const
virtual UInt_tTObject::GetUniqueID() const
virtual Bool_tTObject::HandleTimer(TTimer* timer)
virtual ULong_tTObject::Hash() const
virtual voidTObject::Info(const char* method, const char* msgfmt) const
virtual Bool_tTObject::InheritsFrom(const char* classname) const
virtual Bool_tTObject::InheritsFrom(const TClass* cl) const
virtual voidTObject::Inspect() constMENU
voidTObject::InvertBit(UInt_t f)
virtual TClass*IsA() const
virtual Bool_tTObject::IsEqual(const TObject* obj) const
virtual Bool_tTObject::IsFolder() const
Bool_tTObject::IsOnHeap() const
virtual Bool_tTObject::IsSortable() const
Bool_tTObject::IsZombie() const
virtual voidTObject::ls(Option_t* option = "") const
voidTObject::MayNotUse(const char* method) const
virtual Bool_tTObject::Notify()
voidTObject::Obsolete(const char* method, const char* asOfVers, const char* removedFromVers) const
static voidTObject::operator delete(void* ptr)
static voidTObject::operator delete(void* ptr, void* vp)
static voidTObject::operator delete[](void* ptr)
static voidTObject::operator delete[](void* ptr, void* vp)
void*TObject::operator new(size_t sz)
void*TObject::operator new(size_t sz, void* vp)
void*TObject::operator new[](size_t sz)
void*TObject::operator new[](size_t sz, void* vp)
HKalTrackState&operator=(const HKalTrackState&)
virtual voidTObject::Paint(Option_t* option = "")
virtual voidTObject::Pop()
virtual voidprint(Option_t* opt = "") const
virtual voidTObject::Print(Option_t* option = "") const
virtual Int_tTObject::Read(const char* name)
virtual voidTObject::RecursiveRemove(TObject* obj)
voidTObject::ResetBit(UInt_t f)
virtual voidTObject::SaveAs(const char* filename = "", Option_t* option = "") constMENU
virtual voidTObject::SavePrimitive(ostream& out, Option_t* option = "")
voidTObject::SetBit(UInt_t f)
voidTObject::SetBit(UInt_t f, Bool_t set)
virtual voidsetCovMat(const TMatrixD& fCov)
virtual voidTObject::SetDrawOption(Option_t* option = "")MENU
static voidTObject::SetDtorOnly(void* obj)
static voidTObject::SetObjectStat(Bool_t stat)
virtual voidsetProcNoiseMat(const TMatrixD& fProc)
virtual voidsetProjMat(const TMatrixD& fProj)
virtual voidsetPropMat(const TMatrixD& fProp)
virtual voidsetStateVec(const TVectorD& sv)
virtual voidsetStateVec(Double_t qp, const TVector3& pos, const TVector3& dir)
virtual voidTObject::SetUniqueID(UInt_t uid)
virtual voidShowMembers(TMemberInspector&)
virtual voidStreamer(TBuffer&)
voidStreamerNVirtual(TBuffer& ClassDef_StreamerNVirtual_b)
virtual voidTObject::SysError(const char* method, const char* msgfmt) const
Bool_tTObject::TestBit(UInt_t f) const
Int_tTObject::TestBits(UInt_t f) const
virtual voidtransform(const TRotation& transMat, const HKalPlane& plane)
static voidtransformFromLayerToSector(TVectorD& svSec, const TVectorD& svLay, const HKalPlane& plane)
static voidtransformFromSectorToLayer(TVectorD& svLay, const TVectorD& svSec, const HKalPlane& plane)
virtual voidTObject::UseCurrentStyle()
virtual voidTObject::Warning(const char* method, const char* msgfmt) const
virtual Int_tTObject::Write(const char* name = 0, Int_t option = 0, Int_t bufsize = 0)
virtual Int_tTObject::Write(const char* name = 0, Int_t option = 0, Int_t bufsize = 0) const
protected:
virtual voidTObject::DoError(int level, const char* location, const char* fmt, va_list va) const
voidTObject::MakeZombie()

Data Members

private:
TMatrixDfCovariance! Covariance matrix.
TMatrixDfProcessNoise! Process noise matrix.
TMatrixDfProjector! Projector matrix.
TMatrixDfPropagator! Propagator matrix.
TVectorDstateVec! Vector with track state parameters.
Kalman::kalFilterTypestype! State is either for the prediction, filtering, smoothing or inverse filtering step in the Kalman filter.

Class Charts

Inheritance Inherited Members Includes Libraries
Class Charts

Function documentation

HKalTrackState(Kalman::kalFilterTypes stateType, Int_t measDim, Int_t stateDim)
 Create state with dummy state vector and matrices.
 stateType: state is either for the prediction, filtering, smoothing or inverse filtering step in the Kalman filter
 measDim:  dimension of measurement vector
 stateDim: dimension of state vector.
HKalTrackState(Kalman::kalFilterTypes stateType, const TVectorD& sv, Int_t measDim)
 Create a state with dummy matrices and sets the state vector equal to function
 parameter sv.
 stateType: state is either for the prediction, filtering, smoothing or inverse filtering step in the Kalman filter
 sv:        state vector.
 measDim:   dimension of measurement vector.
~HKalTrackState()
void calcDir(TVector3& dir) const
 Calculates a direction unit vector from this state.
 dir: track direction (return value).
void calcDir(TVector3& dir, const TVectorD& sv)
 Calculates a direction unit vector from a state vector given by function
 parameter sv.
 dir: track direction (return value)
 sv:  state vector to calculate direction from.
void calcMomVec(TVector3& dir) const
 Calculates the momentum vector from this state.
 dir: track momentum (return value).
void calcMomVec(TVector3& dir, const TVectorD& sv)
 Calculates momentum vector from a state vector given by function
 parameter sv.
 dir: track momentum (return value)
 sv:  state vector to calculate direction from.
Bool_t calcPosAtPlane(TVector3& pos, const HKalPlane& plane) const
 Calculate the track position as a three dimensional vector from this state.
 Only two coordinates are stored in the state vector. The z-coordinate is
 calculated by going in z-direction and find the intersection with a plane.
 pos:   track position (return value).
 plane: the plane the track is on.
Bool_t calcPosAtPlane(TVector3& pos, const HKalPlane& plane, const TVectorD& sv)
 Calculate the track position as a three dimensional vector from the state
 given by function parameter sv.
 Only two coordinates are stored in the state vector. The z-coordinate is
 calculated by going in z-direction and find the intersection with a plane.
 pos:   track position (return value).
 plane: the plane the track is on.
Bool_t calcPosAndDirAtPlane(TVector3& pos, TVector3& dir, const HKalPlane& plane) const
 Calculate the track position and direction as three dimensional vectors
 from this state.
 pos:   track position (return value)
 dir:   track direction (return value)
 plane: the plane the track is on
void calcStateVec(TVectorD& sv, Double_t qp, const TVector3& pos, const TVector3& dir)
 Return the state vector from a position, direction and momentum.
 sv:  the state vector (return value)
 qp:  particle charage in elementary charges divided by magnitude of momentum in MeV/c
 pos: track position
 dir: track direction.
void clear()
 Resets state vector to zeroes, matrices to unit matrix.
void print(Option_t* opt = "") const
 Prints information about the state.
 opt: Determines what information about the object will be printed.
 Capitalization does not matter.
 If opt contains:
 "S": print state vector
 "F": print propagator matrix
 "C": print covariance matrix
 "H": print projector matrix
 "Q": print process noise matrix
 If opt is empty all of the above will be printed.
void transform(const TRotation& transMat, const HKalPlane& plane)
 Transform the state vector of this object using a rotation matrix.
 transMat: the rotation matrix
 plane:    the plane the track is on. Needed to calculate the z-coordinate of the
           track position.
void transformFromLayerToSector(TVectorD& svSec, const TVectorD& svLay, const HKalPlane& plane)
 Transform a state vector given in virtual layer coordinates to sector coordinates.

 Output:
 svSec: Track state vector in sector coordinates.

 Input:
 svLay: Track state vector in virtual layer coordinates.
 plane: The plane defining the virtual layer system.

 State vector in sector coordinate system: (x,y,tx,ty,qp).
 State vector in layer  coordinate system: (u,v,tu,tv,qp).
 Layer coordinate system is defined by an origin vector O,
 two orthogonal Axis U and V on the layer and a normal vector W.
 All axis are given in sector coordinate system.

 Transformation from layer to sector:
 position in sector coordinates:
 P = O + u*U + v*V
 direction in sector coordinates:
 A = 1/sqrt(1 + tu^2 + tv^2) * (w + tu*U + tv*V)
 tx = A.x() / A.z()
 ty = A.y() / A.z()
void transformFromSectorToLayer(TVectorD& svLay, const TVectorD& svSec, const HKalPlane& plane)
 Transform a state vector given in virtual layer coordinates to sector coordinates.

 Output:
 svLay: Track state vector in virtual layer coordinates.

 Input:
 svSec: Track state vector in sector coordinates.
 plane: The plane defining the virtual layer system.

 State vector in sector coordinate system: (x,y,tx,ty,qp).
 State vector in layer  coordinate system: (u,v,tu,tv,qp).
 Layer coordinate system is defined by an origin vector O,
 two orthogonal Axis U and V on the layer and a normal vector W.
 All axis are given in sector coordinate system.

 Transformation from sector to layer:
 position on layer:
 u = (pos_sector - O) * U
 v = (pos_sector - O) * V
 direction on layer:
 tu = (A * U) / (A * W)
 tv = (A * V) / (A * W)
 with A = 1/sqrt(1 + tx^2 + ty^2) * (tx, ty, 1)^T the normalized direction vector in the sector coordinate system
void setCovMat(const TMatrixD& fCov)
 Replace the state's covariance matrix.
void setPropMat(const TMatrixD& fProp)
 Replace the state's propagator matrix.
void setProjMat(const TMatrixD& fProj)
 Replace the state's projector matrix.
void setProcNoiseMat(const TMatrixD& fProc)
 Replace the state's process noise matrix.
void setStateVec(const TVectorD& sv)
 Replace the state's track parameter vector.
HKalTrackState()
{}
HKalTrackState(Kalman::kalFilterTypes stateType, Int_t measDim, Int_t stateDim)
TMatrixD const& getCovMat() const
{ return fCovariance; }
TMatrixD const& getPropMat() const
{ return fPropagator; }
TMatrixD const& getProjMat() const
{ return fProjector; }
TMatrixD const& getProcNoiseMat() const
{ return fProcessNoise; }
TVectorD const& getStateVec() const
{ return stateVec; }
Int_t getStateType() const
{ return type; }
Double_t getStateParam(Kalman::kalStateIdx par) const
{ return stateVec(par); }
Int_t getMeasDim() const
{ return fProjector.GetNrows(); }
Int_t getStateDim() const
{ return stateVec.GetNrows(); }
void setStateVec(const TVectorD& sv)