//--------------------------------------------------------------------------------- // The AliKFParticleBaseSIMD 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 describes general mathematics which is used by AliKFParticle class // // -= Copyright © ALICE HLT Group =- //_________________________________________________________________________________ #ifndef ALIKFPARTICLEBASESIMD_H #define ALIKFPARTICLEBASESIMD_H #include "TObject.h" #include "vectors/P4_F32vec4.h" class AliKFParticleBaseSIMD { public: //* //* ABSTRACT METHODS HAVE TO BE DEFINED IN USER CLASS //* //* Virtual method to access the magnetic field virtual void GetFieldValue(const fvec xyz[], fvec B[]) const = 0; //* Virtual methods needed for particle transportation //* One can use particular implementations for collider (only Bz component) //* geometry and for fixed-target (CBM-like) geometry which are provided below //* in TRANSPORT section //* Get dS to xyz[] space point virtual fvec GetDStoPoint( const fvec xyz[] ) const = 0; //* Get dS to other particle p (dSp for particle p also returned) virtual void GetDStoParticle( const AliKFParticleBaseSIMD &p, fvec &DS, fvec &DSp ) const = 0; //* Transport on dS value along trajectory, output to P,C virtual void Transport( fvec dS, fvec P[], fvec C[] ) const = 0; //* //* INITIALIZATION //* //* Constructor AliKFParticleBaseSIMD(); //* Destructor virtual ~AliKFParticleBaseSIMD() { ; } //* Initialisation from "cartesian" coordinates ( X Y Z Px Py Pz ) //* Parameters, covariance matrix, charge, and mass hypothesis should be provided void Initialize( const fvec Param[], const fvec Cov[], fvec Charge, fvec Mass ); //* Initialise covariance matrix and set current parameters to 0.0 void Initialize(); //* Set decay vertex parameters for linearisation void SetVtxGuess( fvec x, fvec y, fvec z ); //* //* ACCESSORS //* //* Simple accessors fvec GetX () const { return fP[0]; } fvec GetY () const { return fP[1]; } fvec GetZ () const { return fP[2]; } fvec GetPx () const { return fP[3]; } fvec GetPy () const { return fP[4]; } fvec GetPz () const { return fP[5]; } fvec GetE () const { return fP[6]; } fvec GetS () const { return fP[7]; } fvec GetQ () const { return fQ; } fvec GetChi2 () const { return fChi2; } Int_t GetNDF () const { return fNDF; } fvec GetParameter ( Int_t i ) const { return fP[i]; } fvec GetCovariance( Int_t i ) const { return fC[i]; } fvec GetCovariance( Int_t i, Int_t j ) const { return fC[IJ(i,j)]; } //* Accessors with calculations( &value, &estimated sigma ) //* error flag returned (0 means no error during calculations) fvec GetMomentum ( fvec &P, fvec &SigmaP ) const ; fvec GetPt ( fvec &Pt, fvec &SigmaPt ) const ; fvec GetEta ( fvec &Eta, fvec &SigmaEta ) const ; fvec GetPhi ( fvec &Phi, fvec &SigmaPhi ) const ; fvec GetMass ( fvec &M, fvec &SigmaM ) const ; fvec GetDecayLength ( fvec &L, fvec &SigmaL ) const ; fvec GetDecayLengthXY ( fvec &L, fvec &SigmaL ) const ; fvec GetLifeTime ( fvec &T, fvec &SigmaT ) const ; fvec GetR ( fvec &R, fvec &SigmaR ) const ; //* //* MODIFIERS //* fvec & X () { return fP[0]; } fvec & Y () { return fP[1]; } fvec & Z () { return fP[2]; } fvec & Px () { return fP[3]; } fvec & Py () { return fP[4]; } fvec & Pz () { return fP[5]; } fvec & E () { return fP[6]; } fvec & S () { return fP[7]; } fvec & Q () { return fQ; } fvec & Chi2 () { return fChi2; } Int_t & NDF () { return fNDF; } fvec & Parameter ( Int_t i ) { return fP[i]; } fvec & Covariance( Int_t i ) { return fC[i]; } fvec & Covariance( Int_t i, Int_t j ) { return fC[IJ(i,j)]; } //* //* CONSTRUCTION OF THE PARTICLE BY ITS DAUGHTERS AND MOTHER //* USING THE KALMAN FILTER METHOD //* //* Simple way to add daughter ex. D0+= Pion; void operator +=( const AliKFParticleBaseSIMD &Daughter ); //* Add daughter track to the particle void AddDaughter( const AliKFParticleBaseSIMD &Daughter ); void AddDaughterWithEnergyFit( const AliKFParticleBaseSIMD &Daughter ); void AddDaughterWithEnergyCalc( const AliKFParticleBaseSIMD &Daughter ); //* Set production vertex void SetProductionVertex( const AliKFParticleBaseSIMD &Vtx ); //* Set mass constraint void SetMassConstraint( fvec Mass, fvec SigmaMass = 0 ); //* Set no decay length for resonances void SetNoDecayLength(); //* Everything in one go void Construct( const AliKFParticleBaseSIMD *vDaughters[], Int_t NDaughters, const AliKFParticleBaseSIMD *ProdVtx=0, Float_t Mass=-1, Bool_t 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 on dS parameter (SignedPath/Momentum) void TransportToDS( fvec dS ); //* Particular extrapolators one can use fvec GetDStoPointBz( fvec Bz, const fvec xyz[] ) const; void GetDStoParticleBz( fvec Bz, const AliKFParticleBaseSIMD &p, fvec &dS, fvec &dS1 ) const ; // fvec GetDStoPointCBM( const fvec xyz[] ) const; void TransportBz( fvec Bz, fvec dS, fvec P[], fvec C[] ) const; void TransportCBM( fvec dS, fvec P[], fvec C[] ) const; //* //* OTHER UTILITIES //* //* Calculate distance from another object [cm] fvec GetDistanceFromVertex( const fvec vtx[] ) const; fvec GetDistanceFromVertex( const AliKFParticleBaseSIMD &Vtx ) const; fvec GetDistanceFromParticle( const AliKFParticleBaseSIMD &p ) const; //* Calculate sqrt(Chi2/ndf) deviation from vertex //* v = [xyz], Cv=[Cxx,Cxy,Cyy,Cxz,Cyz,Czz]-covariance matrix fvec GetDeviationFromVertex( const fvec v[], const fvec Cv[]=0 ) const; fvec GetDeviationFromVertex( const AliKFParticleBaseSIMD &Vtx ) const; fvec GetDeviationFromParticle( const AliKFParticleBaseSIMD &p ) const; //* Subtract the particle from the vertex void SubtractFromVertex( AliKFParticleBaseSIMD &Vtx ) const; //* Special method for creating gammas void ConstructGammaBz( const AliKFParticleBaseSIMD &daughter1, const AliKFParticleBaseSIMD &daughter2, fvec Bz ); protected: static Int_t IJ( Int_t i, Int_t j ){ return ( j<=i ) ? i*(i+1)/2+j :j*(j+1)/2+i; } fvec & Cij( Int_t i, Int_t j ){ return fC[IJ(i,j)]; } static fvec atan2SIMD(const fvec &y, const fvec &x); void Convert( bool ToProduction ); void TransportLine( fvec S, fvec P[], fvec C[] ) const ; fvec GetDStoPointLine( const fvec xyz[] ) const; static fvec InvertSym3( const fvec A[], fvec Ainv[] ); static void MultQSQt( const fvec Q[], const fvec S[], fvec SOut[] ); static fvec GetSCorrection( const fvec Part[], const fvec XYZ[] ); void GetMeasurement( const fvec XYZ[], fvec m[], fvec V[] ) const ; fvec fP[8]; //* Main particle parameters {X,Y,Z,Px,Py,Pz,E,S[=DecayLength/P]} fvec fC[36]; //* Low-triangle covariance matrix of fP fvec fQ; //* Particle charge Int_t fNDF; //* Number of degrees of freedom fvec fChi2; //* Chi^2 fvec fSFromDecay; //* Distance from decay vertex to current position Bool_t fAtProductionVertex; //* Flag shows that the particle error along //* its trajectory is taken from production vertex fvec fVtxGuess[3]; //* Guess for the position of the decay vertex //* ( used for linearisation of equations ) Bool_t fIsLinearized; //* Flag shows that the guess is present //* Determines the method for the particle construction. //* 0 - Energy considered as an independent veriable, fitted independently from momentum, without any constraints on mass //* 1 - Energy considered as a dependent variable, calculated from the momentum and mass hypothesis static const Int_t fConstructMethod = 0; }; inline fvec AliKFParticleBaseSIMD::atan2SIMD(const fvec &y, const fvec &x) { static const fvec ZERO = 0; static const fvec small = 1.e-10f; static const fvec minusOne = -1; static const fvec pi = fvec(3.1415926535); static const fvec pi_2 = pi/2.f; static const fvec pi_4 = pi/4.f; static const fvec tan_3pi_8 = fvec(2.414213562373095); static const fvec tan_pi_8 = fvec(0.414213562373095); const fvec xZero = fvec( fabs(x) < small ); const fvec yZero = fvec( fabs(y) < small ); const fvec xNeg = fvec( x < ZERO ); const fvec yNeg = fvec( y < ZERO ); const fvec absX = fabs(x); const fvec absY = fabs(y); fvec a = absY / absX; const fvec gt_tan_3pi_8 = fvec(tan_3pi_8 < a); const fvec gt_tan_pi_8 = fvec(tan_pi_8 < a) & (!gt_tan_3pi_8); fvec b = ZERO; b += fvec(gt_tan_3pi_8 & pi_2); b += fvec(gt_tan_pi_8 & pi_4); a = fvec(gt_tan_3pi_8 & (minusOne / a)) + fvec((!gt_tan_3pi_8) & a); a = fvec(gt_tan_pi_8 & ((absY - absX) / (absY + absX))) + fvec((!gt_tan_pi_8) & a); const fvec a2 = a * a; static const fvec C1 = 8.05374449538e-2, C2 = 1.38776856032E-1, C3 = 1.99777106478E-1, C4 = 3.33329491539E-1; b += (((C1 * a2 - C2) * a2 + C3) * a2 - C4) * a2 * a + a; const fvec neg = fvec(xNeg ^ yNeg); b = (neg & (-b)) + ((!neg) & b); const fvec xNegy = fvec(xNeg & (!yNeg)); b += (xNegy & pi); const fvec allNeg = fvec(xNeg & yNeg); b -= (allNeg & pi); const fvec allZero = fvec(xZero & yZero); b = ((!allZero) & b); const fvec xZeroyNeg = fvec(xZero & yNeg); b = -(xZeroyNeg & pi_2) + ((!xZeroyNeg) & b); return b; } #endif