//--------------------------------------------------------------------------------- // The AliKFParticle class // . // @author S.Gorbunov, I.Kisel // @version 1.0 // @since 13.05.07 // // Class to reconstruct and store the decayed particle parameters. // The method is described in CBM-SOFT note 2007-003, // ``Reconstruction of decayed particles based on the Kalman filter'', // http://www.gsi.de/documents/DOC-2007-May-14-1.pdf // // This class is ALICE interface to general mathematics in AliKFParticleBase // // -= Copyright © ALICE HLT Group =- //_________________________________________________________________________________ #ifndef ALIKFPARTICLE_H #define ALIKFPARTICLE_H #include "AliKFParticleBase.h" #include "TMath.h" class KFPTrack; class KFPVertex; class AliKFParticle :public AliKFParticleBase { public: //* //* INITIALIZATION //* //* Set magnetic field for all particles static void SetField( double Bz ); //* Constructor (empty) AliKFParticle():AliKFParticleBase(){}; //* Destructor (empty) ~AliKFParticle(){}; //* Construction of mother particle by its 2-3-4 daughters AliKFParticle( const AliKFParticle &d1, const AliKFParticle &d2 ); AliKFParticle( const AliKFParticle &d1, const AliKFParticle &d2, const AliKFParticle &d3 ); AliKFParticle( const AliKFParticle &d1, const AliKFParticle &d2, const AliKFParticle &d3, const AliKFParticle &d4 ); //* Initialisation from "cartesian" coordinates ( X Y Z Px Py Pz ) //* Parameters, covariance matrix, charge and PID hypothesis should be provided void Create( const double Param[], const double Cov[], int Charge, int PID ); //* Initialisation from ALICE track, PID hypothesis shoould be provided AliKFParticle( const KFPTrack &track, int PID ); //* Initialisation from VVertex AliKFParticle( const KFPVertex &vertex ); //* Initialise covariance matrix and set current parameters to 0.0 void Initialize(); //* Set decay vertex parameters for linearisation void SetVtxGuess( double x, double y, double z ); //* //* ACCESSORS //* //* Simple accessors double GetX () const ; //* x of current position double GetY () const ; //* y of current position double GetZ () const ; //* z of current position double GetPx () const ; //* x-compoment of 3-momentum double GetPy () const ; //* y-compoment of 3-momentum double GetPz () const ; //* z-compoment of 3-momentum double GetE () const ; //* energy double GetS () const ; //* decay length / momentum int GetQ () const ; //* charge double GetChi2 () const ; //* chi^2 int GetNDF () const ; //* Number of Degrees of Freedom double GetParameter ( int i ) const ; double GetCovariance( int i ) const ; double GetCovariance( int i, int j ) const ; //* Accessors with calculations, value returned w/o error flag double GetP () const; //* momentum double GetPt () const; //* transverse momentum double GetEta () const; //* pseudorapidity double GetPhi () const; //* phi double GetMomentum () const; //* momentum (same as GetP() ) double GetMass () const; //* mass double GetDecayLength () const; //* decay length double GetLifeTime () const; //* life time double GetR () const; //* distance to the origin //* Accessors to estimated errors double GetErrX () const ; //* x of current position double GetErrY () const ; //* y of current position double GetErrZ () const ; //* z of current position double GetErrPx () const ; //* x-compoment of 3-momentum double GetErrPy () const ; //* y-compoment of 3-momentum double GetErrPz () const ; //* z-compoment of 3-momentum double GetErrE () const ; //* energy double GetErrS () const ; //* decay length / momentum double GetErrP () const ; //* momentum double GetErrPt () const ; //* transverse momentum double GetErrEta () const ; //* pseudorapidity double GetErrPhi () const ; //* phi double GetErrMomentum () const ; //* momentum double GetErrMass () const ; //* mass double GetErrDecayLength () const ; //* decay length double GetErrLifeTime () const ; //* life time double GetErrR () const ; //* distance to the origin //* Accessors with calculations( &value, &estimated sigma ) //* error flag returned (0 means no error during calculations) int GetP ( double &P, double &SigmaP ) const ; //* momentum int GetPt ( double &Pt, double &SigmaPt ) const ; //* transverse momentum int GetEta ( double &Eta, double &SigmaEta ) const ; //* pseudorapidity int GetPhi ( double &Phi, double &SigmaPhi ) const ; //* phi int GetMomentum ( double &P, double &SigmaP ) const ; //* momentum int GetMass ( double &M, double &SigmaM ) const ; //* mass int GetDecayLength ( double &L, double &SigmaL ) const ; //* decay length int GetLifeTime ( double &T, double &SigmaT ) const ; //* life time int GetR ( double &R, double &SigmaR ) const ; //* R //* //* MODIFIERS //* double & X () ; double & Y () ; double & Z () ; double & Px () ; double & Py () ; double & Pz () ; double & E () ; double & S () ; int & Q () ; double & Chi2 () ; int & NDF () ; double & Parameter ( int i ) ; double & Covariance( int i ) ; double & Covariance( int i, int j ) ; double * Parameters () ; double * CovarianceMatrix() ; //* //* CONSTRUCTION OF THE PARTICLE BY ITS DAUGHTERS AND MOTHER //* USING THE KALMAN FILTER METHOD //* //* Add daughter to the particle void AddDaughter( const AliKFParticle &Daughter ); //* Add daughter via += operator: ex.{ D0; D0+=Pion; D0+= Kaon; } void operator +=( const AliKFParticle &Daughter ); //* Set production vertex void SetProductionVertex( const AliKFParticle &Vtx ); //* Set mass constraint void SetMassConstraint( double Mass, double SigmaMass = 0 ); //* Set no decay length for resonances void SetNoDecayLength(); //* Everything in one go void Construct( const AliKFParticle *vDaughters[], int NDaughters, const AliKFParticle *ProdVtx=0, double Mass=-1, bool IsConstrained=0 ); //* //* TRANSPORT //* //* ( main transportation parameter is S = SignedPath/Momentum ) //* ( parameters of decay & production vertices are stored locally ) //* //* Transport the particle to its decay vertex void TransportToDecayVertex(); //* Transport the particle to its production vertex void TransportToProductionVertex(); //* Transport the particle close to xyz[] point void TransportToPoint( const double xyz[] ); //* Transport the particle close to VVertex void TransportToVertex( const KFPVertex &v ); //* Transport the particle close to another particle p void TransportToParticle( const AliKFParticle &p ); //* Transport the particle on dS parameter (SignedPath/Momentum) void TransportToDS( double dS ); //* Get dS to a certain space point double GetDStoPoint( const double xyz[] ) const ; //* Get dS to other particle p (dSp for particle p also returned) void GetDStoParticle( const AliKFParticle &p, double &DS, double &DSp ) const ; //* Get dS to other particle p in XY-plane void GetDStoParticleXY( const AliKFParticleBase &p, double &DS, double &DSp ) const ; //* //* OTHER UTILITIES //* //* Calculate distance from another object [cm] double GetDistanceFromVertex( const double vtx[] ) const ; double GetDistanceFromVertex( const AliKFParticle &Vtx ) const ; double GetDistanceFromVertex( const KFPVertex &Vtx ) const ; double GetDistanceFromParticle( const AliKFParticle &p ) const ; //* Calculate sqrt(Chi2/ndf) deviation from another object //* ( v = [xyz]-vertex, Cv=[Cxx,Cxy,Cyy,Cxz,Cyz,Czz]-covariance matrix ) double GetDeviationFromVertex( const double v[], const double Cv[]=0 ) const ; double GetDeviationFromVertex( const AliKFParticle &Vtx ) const ; double GetDeviationFromVertex( const KFPVertex &Vtx ) const ; double GetDeviationFromParticle( const AliKFParticle &p ) const ; //* Calculate distance from another object [cm] in XY-plane double GetDistanceFromVertexXY( const double vtx[] ) const ; double GetDistanceFromVertexXY( const AliKFParticle &Vtx ) const ; double GetDistanceFromVertexXY( const KFPVertex &Vtx ) const ; double GetDistanceFromParticleXY( const AliKFParticle &p ) const ; //* Calculate sqrt(Chi2/ndf) deviation from another object in XY plane //* ( v = [xyz]-vertex, Cv=[Cxx,Cxy,Cyy,Cxz,Cyz,Czz]-covariance matrix ) double GetDeviationFromVertexXY( const double v[], const double Cv[]=0 ) const ; double GetDeviationFromVertexXY( const AliKFParticle &Vtx ) const ; double GetDeviationFromVertexXY( const KFPVertex &Vtx ) const ; double GetDeviationFromParticleXY( const AliKFParticle &p ) const ; //* Calculate opennig angle between two particles double GetAngle ( const AliKFParticle &p ) const ; double GetAngleXY( const AliKFParticle &p ) const ; double GetAngleRZ( const AliKFParticle &p ) const ; //* Subtract the particle from the vertex void SubtractFromVertex( AliKFParticle &v ) const ; protected: //* //* INTERNAL STUFF //* //* Method to access ALICE field static double GetFieldAlice(); //* Other methods required by the abstract AliKFParticleBase class void GetFieldValue( const double xyz[], double B[] ) const ; void GetDStoParticle( const AliKFParticleBase &p, double &DS, double &DSp )const ; void Transport( double dS, double P[], double C[] ) const ; static void GetExternalTrackParam( const AliKFParticleBase &p, double &X, double &Alpha, double P[5] ) ; //void GetDStoParticleALICE( const AliKFParticleBase &p, double &DS, double &DS1 ) const; private: static double fgBz; //* Bz compoment of the magnetic field }; //--------------------------------------------------------------------- // // Inline implementation of the AliKFParticle methods // //--------------------------------------------------------------------- inline void AliKFParticle::SetField( double Bz ) { fgBz = -Bz;//!!! } inline AliKFParticle::AliKFParticle( const AliKFParticle &d1, const AliKFParticle &d2 ) { AliKFParticle mother; mother+= d1; mother+= d2; *this = mother; } inline AliKFParticle::AliKFParticle( const AliKFParticle &d1, const AliKFParticle &d2, const AliKFParticle &d3 ) { AliKFParticle mother; mother+= d1; mother+= d2; mother+= d3; *this = mother; } inline AliKFParticle::AliKFParticle( const AliKFParticle &d1, const AliKFParticle &d2, const AliKFParticle &d3, const AliKFParticle &d4 ) { AliKFParticle mother; mother+= d1; mother+= d2; mother+= d3; mother+= d4; *this = mother; } inline void AliKFParticle::Initialize() { AliKFParticleBase::Initialize(); } inline void AliKFParticle::SetVtxGuess( double x, double y, double z ) { AliKFParticleBase::SetVtxGuess(x,y,z); } inline double AliKFParticle::GetX () const { return AliKFParticleBase::GetX(); } inline double AliKFParticle::GetY () const { return AliKFParticleBase::GetY(); } inline double AliKFParticle::GetZ () const { return AliKFParticleBase::GetZ(); } inline double AliKFParticle::GetPx () const { return AliKFParticleBase::GetPx(); } inline double AliKFParticle::GetPy () const { return AliKFParticleBase::GetPy(); } inline double AliKFParticle::GetPz () const { return AliKFParticleBase::GetPz(); } inline double AliKFParticle::GetE () const { return AliKFParticleBase::GetE(); } inline double AliKFParticle::GetS () const { return AliKFParticleBase::GetS(); } inline int AliKFParticle::GetQ () const { return AliKFParticleBase::GetQ(); } inline double AliKFParticle::GetChi2 () const { return AliKFParticleBase::GetChi2(); } inline int AliKFParticle::GetNDF () const { return AliKFParticleBase::GetNDF(); } inline double AliKFParticle::GetParameter ( int i ) const { return AliKFParticleBase::GetParameter(i); } inline double AliKFParticle::GetCovariance( int i ) const { return AliKFParticleBase::GetCovariance(i); } inline double AliKFParticle::GetCovariance( int i, int j ) const { return AliKFParticleBase::GetCovariance(i,j); } inline double AliKFParticle::GetP () const { double par, err; if( AliKFParticleBase::GetMomentum( par, err ) ) return 0; else return par; } inline double AliKFParticle::GetPt () const { double par, err; if( AliKFParticleBase::GetPt( par, err ) ) return 0; else return par; } inline double AliKFParticle::GetEta () const { double par, err; if( AliKFParticleBase::GetEta( par, err ) ) return 0; else return par; } inline double AliKFParticle::GetPhi () const { double par, err; if( AliKFParticleBase::GetPhi( par, err ) ) return 0; else return par; } inline double AliKFParticle::GetMomentum () const { double par, err; if( AliKFParticleBase::GetMomentum( par, err ) ) return 0; else return par; } inline double AliKFParticle::GetMass () const { double par, err; if( AliKFParticleBase::GetMass( par, err ) ) return 0; else return par; } inline double AliKFParticle::GetDecayLength () const { double par, err; if( AliKFParticleBase::GetDecayLength( par, err ) ) return 0; else return par; } inline double AliKFParticle::GetLifeTime () const { double par, err; if( AliKFParticleBase::GetLifeTime( par, err ) ) return 0; else return par; } inline double AliKFParticle::GetR () const { double par, err; if( AliKFParticleBase::GetR( par, err ) ) return 0; else return par; } inline double AliKFParticle::GetErrX () const { return TMath::Sqrt(TMath::Abs( GetCovariance(0,0) )); } inline double AliKFParticle::GetErrY () const { return TMath::Sqrt(TMath::Abs( GetCovariance(1,1) )); } inline double AliKFParticle::GetErrZ () const { return TMath::Sqrt(TMath::Abs( GetCovariance(2,2) )); } inline double AliKFParticle::GetErrPx () const { return TMath::Sqrt(TMath::Abs( GetCovariance(3,3) )); } inline double AliKFParticle::GetErrPy () const { return TMath::Sqrt(TMath::Abs( GetCovariance(4,4) )); } inline double AliKFParticle::GetErrPz () const { return TMath::Sqrt(TMath::Abs( GetCovariance(5,5) )); } inline double AliKFParticle::GetErrE () const { return TMath::Sqrt(TMath::Abs( GetCovariance(6,6) )); } inline double AliKFParticle::GetErrS () const { return TMath::Sqrt(TMath::Abs( GetCovariance(7,7) )); } inline double AliKFParticle::GetErrP () const { double par, err; if( AliKFParticleBase::GetMomentum( par, err ) ) return 1.e10; else return err; } inline double AliKFParticle::GetErrPt () const { double par, err; if( AliKFParticleBase::GetPt( par, err ) ) return 1.e10; else return err; } inline double AliKFParticle::GetErrEta () const { double par, err; if( AliKFParticleBase::GetEta( par, err ) ) return 1.e10; else return err; } inline double AliKFParticle::GetErrPhi () const { double par, err; if( AliKFParticleBase::GetPhi( par, err ) ) return 1.e10; else return err; } inline double AliKFParticle::GetErrMomentum () const { double par, err; if( AliKFParticleBase::GetMomentum( par, err ) ) return 1.e10; else return err; } inline double AliKFParticle::GetErrMass () const { double par, err; if( AliKFParticleBase::GetMass( par, err ) ) return 1.e10; else return err; } inline double AliKFParticle::GetErrDecayLength () const { double par, err; if( AliKFParticleBase::GetDecayLength( par, err ) ) return 1.e10; else return err; } inline double AliKFParticle::GetErrLifeTime () const { double par, err; if( AliKFParticleBase::GetLifeTime( par, err ) ) return 1.e10; else return err; } inline double AliKFParticle::GetErrR () const { double par, err; if( AliKFParticleBase::GetR( par, err ) ) return 1.e10; else return err; } inline int AliKFParticle::GetP( double &P, double &SigmaP ) const { return AliKFParticleBase::GetMomentum( P, SigmaP ); } inline int AliKFParticle::GetPt( double &Pt, double &SigmaPt ) const { return AliKFParticleBase::GetPt( Pt, SigmaPt ); } inline int AliKFParticle::GetEta( double &Eta, double &SigmaEta ) const { return AliKFParticleBase::GetEta( Eta, SigmaEta ); } inline int AliKFParticle::GetPhi( double &Phi, double &SigmaPhi ) const { return AliKFParticleBase::GetPhi( Phi, SigmaPhi ); } inline int AliKFParticle::GetMomentum( double &P, double &SigmaP ) const { return AliKFParticleBase::GetMomentum( P, SigmaP ); } inline int AliKFParticle::GetMass( double &M, double &SigmaM ) const { return AliKFParticleBase::GetMass( M, SigmaM ); } inline int AliKFParticle::GetDecayLength( double &L, double &SigmaL ) const { return AliKFParticleBase::GetDecayLength( L, SigmaL ); } inline int AliKFParticle::GetLifeTime( double &T, double &SigmaT ) const { return AliKFParticleBase::GetLifeTime( T, SigmaT ); } inline int AliKFParticle::GetR( double &R, double &SigmaR ) const { return AliKFParticleBase::GetR( R, SigmaR ); } inline double & AliKFParticle::X() { return AliKFParticleBase::X(); } inline double & AliKFParticle::Y() { return AliKFParticleBase::Y(); } inline double & AliKFParticle::Z() { return AliKFParticleBase::Z(); } inline double & AliKFParticle::Px() { return AliKFParticleBase::Px(); } inline double & AliKFParticle::Py() { return AliKFParticleBase::Py(); } inline double & AliKFParticle::Pz() { return AliKFParticleBase::Pz(); } inline double & AliKFParticle::E() { return AliKFParticleBase::E(); } inline double & AliKFParticle::S() { return AliKFParticleBase::S(); } inline int & AliKFParticle::Q() { return AliKFParticleBase::Q(); } inline double & AliKFParticle::Chi2() { return AliKFParticleBase::Chi2(); } inline int & AliKFParticle::NDF() { return AliKFParticleBase::NDF(); } inline double & AliKFParticle::Parameter ( int i ) { return AliKFParticleBase::Parameter(i); } inline double & AliKFParticle::Covariance( int i ) { return AliKFParticleBase::Covariance(i); } inline double & AliKFParticle::Covariance( int i, int j ) { return AliKFParticleBase::Covariance(i,j); } inline double * AliKFParticle::Parameters () { return fP; } inline double * AliKFParticle::CovarianceMatrix() { return fC; } inline void AliKFParticle::operator +=( const AliKFParticle &Daughter ) { AliKFParticleBase::operator +=( Daughter ); } inline void AliKFParticle::AddDaughter( const AliKFParticle &Daughter ) { AliKFParticleBase::AddDaughter( Daughter ); } inline void AliKFParticle::SetProductionVertex( const AliKFParticle &Vtx ) { AliKFParticleBase::SetProductionVertex( Vtx ); } inline void AliKFParticle::SetMassConstraint( double Mass, double SigmaMass ) { AliKFParticleBase::SetMassConstraint( Mass, SigmaMass ); } inline void AliKFParticle::SetNoDecayLength() { AliKFParticleBase::SetNoDecayLength(); } inline void AliKFParticle::Construct( const AliKFParticle *vDaughters[], int NDaughters, const AliKFParticle *ProdVtx, double Mass, bool IsConstrained ) { AliKFParticleBase::Construct( ( const AliKFParticleBase**)vDaughters, NDaughters, ( const AliKFParticleBase*)ProdVtx, Mass, IsConstrained ); } inline void AliKFParticle::TransportToDecayVertex() { AliKFParticleBase::TransportToDecayVertex(); } inline void AliKFParticle::TransportToProductionVertex() { AliKFParticleBase::TransportToProductionVertex(); } inline void AliKFParticle::TransportToPoint( const double xyz[] ) { TransportToDS( GetDStoPoint(xyz) ); } inline void AliKFParticle::TransportToVertex( const KFPVertex &v ) { TransportToPoint( AliKFParticle(v).fP ); } inline void AliKFParticle::TransportToParticle( const AliKFParticle &p ) { double dS, dSp; GetDStoParticle( p, dS, dSp ); TransportToDS( dS ); } inline void AliKFParticle::TransportToDS( double dS ) { AliKFParticleBase::TransportToDS( dS ); } inline double AliKFParticle::GetDStoPoint( const double xyz[] ) const { return AliKFParticleBase::GetDStoPointBz( GetFieldAlice(), xyz ); } inline void AliKFParticle::GetDStoParticle( const AliKFParticle &p, double &DS, double &DSp ) const { GetDStoParticleXY( p, DS, DSp ); } inline double AliKFParticle::GetDistanceFromVertex( const double vtx[] ) const { return AliKFParticleBase::GetDistanceFromVertex( vtx ); } inline double AliKFParticle::GetDeviationFromVertex( const double v[], const double Cv[] ) const { return AliKFParticleBase::GetDeviationFromVertex( v, Cv); } inline double AliKFParticle::GetDistanceFromVertex( const AliKFParticle &Vtx ) const { return AliKFParticleBase::GetDistanceFromVertex( Vtx ); } inline double AliKFParticle::GetDeviationFromVertex( const AliKFParticle &Vtx ) const { return AliKFParticleBase::GetDeviationFromVertex( Vtx ); } inline double AliKFParticle::GetDistanceFromVertex( const KFPVertex &Vtx ) const { return GetDistanceFromVertex( AliKFParticle(Vtx) ); } inline double AliKFParticle::GetDeviationFromVertex( const KFPVertex &Vtx ) const { return GetDeviationFromVertex( AliKFParticle(Vtx) ); } inline double AliKFParticle::GetDistanceFromParticle( const AliKFParticle &p ) const { return AliKFParticleBase::GetDistanceFromParticle( p ); } inline double AliKFParticle::GetDeviationFromParticle( const AliKFParticle &p ) const { return AliKFParticleBase::GetDeviationFromParticle( p ); } inline void AliKFParticle::SubtractFromVertex( AliKFParticle &v ) const { AliKFParticleBase::SubtractFromVertex( v.fP, v.fC, v.fChi2, v.fNDF); } inline double AliKFParticle::GetFieldAlice() { return fgBz; } inline void AliKFParticle::GetFieldValue( const double * /*xyz*/, double B[] ) const { B[0] = B[1] = 0; B[2] = GetFieldAlice(); } inline void AliKFParticle::GetDStoParticle( const AliKFParticleBase &p, double &DS, double &DSp )const { GetDStoParticleXY( p, DS, DSp ); } inline void AliKFParticle::GetDStoParticleXY( const AliKFParticleBase &p, double &DS, double &DSp ) const { AliKFParticleBase::GetDStoParticleBz( GetFieldAlice(), p, DS, DSp ) ; //GetDStoParticleALICE( p, DS, DSp ) ; } inline void AliKFParticle::Transport( double dS, double P[], double C[] ) const { AliKFParticleBase::TransportBz( GetFieldAlice(), dS, P, C ); } #endif