1 //---------------------------------------------------------------------------------
2 // The AliKFParticle class
4 // @author S.Gorbunov, I.Kisel
8 // Class to reconstruct and store the decayed particle parameters.
9 // The method is described in CBM-SOFT note 2007-003,
10 // ``Reconstruction of decayed particles based on the Kalman filter'',
11 // http://www.gsi.de/documents/DOC-2007-May-14-1.pdf
13 // This class is ALICE interface to general mathematics in AliKFParticleBase
15 // -= Copyright © ALICE HLT Group =-
16 //_________________________________________________________________________________
18 #ifndef ALIKFPARTICLE_H
19 #define ALIKFPARTICLE_H
21 #include "AliKFParticleBase.h"
27 class AliKFParticle :public AliKFParticleBase
36 //* Set magnetic field for all particles
38 static void SetField( Double_t Bz );
40 //* Constructor (empty)
42 AliKFParticle():AliKFParticleBase(){ ; }
44 //* Destructor (empty)
48 //* Construction of mother particle by its 2-3-4 daughters
50 AliKFParticle( const AliKFParticle &d1, const AliKFParticle &d2 );
52 AliKFParticle( const AliKFParticle &d1, const AliKFParticle &d2,
53 const AliKFParticle &d3 );
55 AliKFParticle( const AliKFParticle &d1, const AliKFParticle &d2,
56 const AliKFParticle &d3, const AliKFParticle &d4 );
58 //* Initialisation from "cartesian" coordinates ( X Y Z Px Py Pz )
59 //* Parameters, covariance matrix, charge and PID hypothesis should be provided
61 void Create( const Double_t Param[], const Double_t Cov[], Int_t Charge, Int_t PID );
63 //* Initialisation from ALICE track, PID hypothesis shoould be provided
65 AliKFParticle( const AliVTrack &track, Int_t PID );
67 //* Initialisation from VVertex
69 AliKFParticle( const AliVVertex &vertex );
71 //* Initialise covariance matrix and set current parameters to 0.0
75 //* Set decay vertex parameters for linearisation
77 void SetVtxGuess( Double_t x, Double_t y, Double_t z );
85 Double_t GetX () const ; //* x of current position
86 Double_t GetY () const ; //* y of current position
87 Double_t GetZ () const ; //* z of current position
88 Double_t GetPx () const ; //* x-compoment of 3-momentum
89 Double_t GetPy () const ; //* y-compoment of 3-momentum
90 Double_t GetPz () const ; //* z-compoment of 3-momentum
91 Double_t GetE () const ; //* energy
92 Double_t GetS () const ; //* decay length / momentum
93 Int_t GetQ () const ; //* charge
94 Double_t GetChi2 () const ; //* chi^2
95 Int_t GetNDF () const ; //* Number of Degrees of Freedom
97 Double_t GetParameter ( int i ) const ;
98 Double_t GetCovariance( int i ) const ;
99 Double_t GetCovariance( int i, int j ) const ;
101 //* Accessors with calculations, value returned w/o error flag
103 Double_t GetP () const; //* momentum
104 Double_t GetPt () const; //* transverse momentum
105 Double_t GetEta () const; //* pseudorapidity
106 Double_t GetPhi () const; //* phi
107 Double_t GetMomentum () const; //* momentum (same as GetP() )
108 Double_t GetMass () const; //* mass
109 Double_t GetDecayLength () const; //* decay length
110 Double_t GetLifeTime () const; //* life time
111 Double_t GetR () const; //* distance to the origin
113 //* Accessors to estimated errors
115 Double_t GetErrX () const ; //* x of current position
116 Double_t GetErrY () const ; //* y of current position
117 Double_t GetErrZ () const ; //* z of current position
118 Double_t GetErrPx () const ; //* x-compoment of 3-momentum
119 Double_t GetErrPy () const ; //* y-compoment of 3-momentum
120 Double_t GetErrPz () const ; //* z-compoment of 3-momentum
121 Double_t GetErrE () const ; //* energy
122 Double_t GetErrS () const ; //* decay length / momentum
123 Double_t GetErrP () const ; //* momentum
124 Double_t GetErrPt () const ; //* transverse momentum
125 Double_t GetErrEta () const ; //* pseudorapidity
126 Double_t GetErrPhi () const ; //* phi
127 Double_t GetErrMomentum () const ; //* momentum
128 Double_t GetErrMass () const ; //* mass
129 Double_t GetErrDecayLength () const ; //* decay length
130 Double_t GetErrLifeTime () const ; //* life time
131 Double_t GetErrR () const ; //* distance to the origin
133 //* Accessors with calculations( &value, &estimated sigma )
134 //* error flag returned (0 means no error during calculations)
136 int GetP ( Double_t &P, Double_t &SigmaP ) const ; //* momentum
137 int GetPt ( Double_t &Pt, Double_t &SigmaPt ) const ; //* transverse momentum
138 int GetEta ( Double_t &Eta, Double_t &SigmaEta ) const ; //* pseudorapidity
139 int GetPhi ( Double_t &Phi, Double_t &SigmaPhi ) const ; //* phi
140 int GetMomentum ( Double_t &P, Double_t &SigmaP ) const ; //* momentum
141 int GetMass ( Double_t &M, Double_t &SigmaM ) const ; //* mass
142 int GetDecayLength ( Double_t &L, Double_t &SigmaL ) const ; //* decay length
143 int GetLifeTime ( Double_t &T, Double_t &SigmaT ) const ; //* life time
144 int GetR ( Double_t &R, Double_t &SigmaR ) const ; //* R
163 Double_t & Parameter ( int i ) ;
164 Double_t & Covariance( int i ) ;
165 Double_t & Covariance( int i, int j ) ;
166 Double_t * Parameters () ;
167 Double_t * CovarianceMatrix() ;
170 //* CONSTRUCTION OF THE PARTICLE BY ITS DAUGHTERS AND MOTHER
171 //* USING THE KALMAN FILTER METHOD
175 //* Add daughter to the particle
177 void AddDaughter( const AliKFParticle &Daughter );
179 //* Add daughter via += operator: ex.{ D0; D0+=Pion; D0+= Kaon; }
181 void operator +=( const AliKFParticle &Daughter );
183 //* Set production vertex
185 void SetProductionVertex( const AliKFParticle &Vtx );
187 //* Set mass constraint
189 void SetMassConstraint( Double_t Mass, Double_t SigmaMass = 0 );
191 //* Set no decay length for resonances
193 void SetNoDecayLength();
195 //* Everything in one go
197 void Construct( const AliKFParticle *vDaughters[], int NDaughters,
198 const AliKFParticle *ProdVtx=0, Double_t Mass=-1, Bool_t IsConstrained=0 );
203 //* ( main transportation parameter is S = SignedPath/Momentum )
204 //* ( parameters of decay & production vertices are stored locally )
207 //* Transport the particle to its decay vertex
209 void TransportToDecayVertex();
211 //* Transport the particle to its production vertex
213 void TransportToProductionVertex();
215 //* Transport the particle close to xyz[] point
217 void TransportToPoint( const Double_t xyz[] );
219 //* Transport the particle close to VVertex
221 void TransportToVertex( const AliVVertex &v );
223 //* Transport the particle close to another particle p
225 void TransportToParticle( const AliKFParticle &p );
227 //* Transport the particle on dS parameter (SignedPath/Momentum)
229 void TransportToDS( Double_t dS );
231 //* Get dS to a certain space point
233 Double_t GetDStoPoint( const Double_t xyz[] ) const ;
235 //* Get dS to other particle p (dSp for particle p also returned)
237 void GetDStoParticle( const AliKFParticle &p,
238 Double_t &DS, Double_t &DSp ) const ;
240 //* Get dS to other particle p in XY-plane
242 void GetDStoParticleXY( const AliKFParticleBase &p,
243 Double_t &DS, Double_t &DSp ) const ;
250 //* Calculate distance from another object [cm]
252 Double_t GetDistanceFromVertex( const Double_t vtx[] ) const ;
253 Double_t GetDistanceFromVertex( const AliKFParticle &Vtx ) const ;
254 Double_t GetDistanceFromVertex( const AliVVertex &Vtx ) const ;
255 Double_t GetDistanceFromParticle( const AliKFParticle &p ) const ;
257 //* Calculate sqrt(Chi2/ndf) deviation from another object
258 //* ( v = [xyz]-vertex, Cv=[Cxx,Cxy,Cyy,Cxz,Cyz,Czz]-covariance matrix )
260 Double_t GetDeviationFromVertex( const Double_t v[], const Double_t Cv[]=0 ) const ;
261 Double_t GetDeviationFromVertex( const AliKFParticle &Vtx ) const ;
262 Double_t GetDeviationFromVertex( const AliVVertex &Vtx ) const ;
263 Double_t GetDeviationFromParticle( const AliKFParticle &p ) const ;
265 //* Calculate distance from another object [cm] in XY-plane
267 Double_t GetDistanceFromVertexXY( const Double_t vtx[] ) const ;
268 Double_t GetDistanceFromVertexXY( const AliKFParticle &Vtx ) const ;
269 Double_t GetDistanceFromVertexXY( const AliVVertex &Vtx ) const ;
270 Double_t GetDistanceFromParticleXY( const AliKFParticle &p ) const ;
272 //* Calculate sqrt(Chi2/ndf) deviation from another object in XY plane
273 //* ( v = [xyz]-vertex, Cv=[Cxx,Cxy,Cyy,Cxz,Cyz,Czz]-covariance matrix )
275 Double_t GetDeviationFromVertexXY( const Double_t v[], const Double_t Cv[]=0 ) const ;
276 Double_t GetDeviationFromVertexXY( const AliKFParticle &Vtx ) const ;
277 Double_t GetDeviationFromVertexXY( const AliVVertex &Vtx ) const ;
278 Double_t GetDeviationFromParticleXY( const AliKFParticle &p ) const ;
280 //* Calculate opennig angle between two particles
282 Double_t GetAngle ( const AliKFParticle &p ) const ;
283 Double_t GetAngleXY( const AliKFParticle &p ) const ;
284 Double_t GetAngleRZ( const AliKFParticle &p ) const ;
286 //* Subtract the particle from the vertex
288 void SubtractFromVertex( AliKFParticle &v ) const ;
296 //* Method to access ALICE field
298 static Double_t GetFieldAlice();
300 //* Other methods required by the abstract AliKFParticleBase class
302 void GetFieldValue( const Double_t xyz[], Double_t B[] ) const ;
303 void GetDStoParticle( const AliKFParticleBase &p, Double_t &DS, Double_t &DSp )const ;
304 void Transport( Double_t dS, Double_t P[], Double_t C[] ) const ;
305 static void GetExternalTrackParam( const AliKFParticleBase &p, Double_t &X, Double_t &Alpha, Double_t P[5] ) ;
307 //void GetDStoParticleALICE( const AliKFParticleBase &p, Double_t &DS, Double_t &DS1 ) const;
312 static Double_t fgBz; //* Bz compoment of the magnetic field
314 ClassDef( AliKFParticle, 1 );
320 //---------------------------------------------------------------------
322 // Inline implementation of the AliKFParticle methods
324 //---------------------------------------------------------------------
327 inline void AliKFParticle::SetField( Double_t Bz )
333 inline AliKFParticle::AliKFParticle( const AliKFParticle &d1,
334 const AliKFParticle &d2 )
336 AliKFParticle mother;
342 inline AliKFParticle::AliKFParticle( const AliKFParticle &d1,
343 const AliKFParticle &d2,
344 const AliKFParticle &d3 )
346 AliKFParticle mother;
353 inline AliKFParticle::AliKFParticle( const AliKFParticle &d1,
354 const AliKFParticle &d2,
355 const AliKFParticle &d3,
356 const AliKFParticle &d4 )
358 AliKFParticle mother;
367 inline void AliKFParticle::Initialize()
369 AliKFParticleBase::Initialize();
372 inline void AliKFParticle::SetVtxGuess( Double_t x, Double_t y, Double_t z )
374 AliKFParticleBase::SetVtxGuess(x,y,z);
377 inline Double_t AliKFParticle::GetX () const
379 return AliKFParticleBase::GetX();
382 inline Double_t AliKFParticle::GetY () const
384 return AliKFParticleBase::GetY();
387 inline Double_t AliKFParticle::GetZ () const
389 return AliKFParticleBase::GetZ();
392 inline Double_t AliKFParticle::GetPx () const
394 return AliKFParticleBase::GetPx();
397 inline Double_t AliKFParticle::GetPy () const
399 return AliKFParticleBase::GetPy();
402 inline Double_t AliKFParticle::GetPz () const
404 return AliKFParticleBase::GetPz();
407 inline Double_t AliKFParticle::GetE () const
409 return AliKFParticleBase::GetE();
412 inline Double_t AliKFParticle::GetS () const
414 return AliKFParticleBase::GetS();
417 inline Int_t AliKFParticle::GetQ () const
419 return AliKFParticleBase::GetQ();
422 inline Double_t AliKFParticle::GetChi2 () const
424 return AliKFParticleBase::GetChi2();
427 inline Int_t AliKFParticle::GetNDF () const
429 return AliKFParticleBase::GetNDF();
432 inline Double_t AliKFParticle::GetParameter ( int i ) const
434 return AliKFParticleBase::GetParameter(i);
437 inline Double_t AliKFParticle::GetCovariance( int i ) const
439 return AliKFParticleBase::GetCovariance(i);
442 inline Double_t AliKFParticle::GetCovariance( int i, int j ) const
444 return AliKFParticleBase::GetCovariance(i,j);
448 inline Double_t AliKFParticle::GetP () const
451 if( AliKFParticleBase::GetMomentum( par, err ) ) return 0;
455 inline Double_t AliKFParticle::GetPt () const
458 if( AliKFParticleBase::GetPt( par, err ) ) return 0;
462 inline Double_t AliKFParticle::GetEta () const
465 if( AliKFParticleBase::GetEta( par, err ) ) return 0;
469 inline Double_t AliKFParticle::GetPhi () const
472 if( AliKFParticleBase::GetPhi( par, err ) ) return 0;
476 inline Double_t AliKFParticle::GetMomentum () const
479 if( AliKFParticleBase::GetMomentum( par, err ) ) return 0;
483 inline Double_t AliKFParticle::GetMass () const
486 if( AliKFParticleBase::GetMass( par, err ) ) return 0;
490 inline Double_t AliKFParticle::GetDecayLength () const
493 if( AliKFParticleBase::GetDecayLength( par, err ) ) return 0;
497 inline Double_t AliKFParticle::GetLifeTime () const
500 if( AliKFParticleBase::GetLifeTime( par, err ) ) return 0;
504 inline Double_t AliKFParticle::GetR () const
507 if( AliKFParticleBase::GetR( par, err ) ) return 0;
511 inline Double_t AliKFParticle::GetErrX () const
513 return TMath::Sqrt(TMath::Abs( GetCovariance(0,0) ));
516 inline Double_t AliKFParticle::GetErrY () const
518 return TMath::Sqrt(TMath::Abs( GetCovariance(1,1) ));
521 inline Double_t AliKFParticle::GetErrZ () const
523 return TMath::Sqrt(TMath::Abs( GetCovariance(2,2) ));
526 inline Double_t AliKFParticle::GetErrPx () const
528 return TMath::Sqrt(TMath::Abs( GetCovariance(3,3) ));
531 inline Double_t AliKFParticle::GetErrPy () const
533 return TMath::Sqrt(TMath::Abs( GetCovariance(4,4) ));
536 inline Double_t AliKFParticle::GetErrPz () const
538 return TMath::Sqrt(TMath::Abs( GetCovariance(5,5) ));
541 inline Double_t AliKFParticle::GetErrE () const
543 return TMath::Sqrt(TMath::Abs( GetCovariance(6,6) ));
546 inline Double_t AliKFParticle::GetErrS () const
548 return TMath::Sqrt(TMath::Abs( GetCovariance(7,7) ));
551 inline Double_t AliKFParticle::GetErrP () const
554 if( AliKFParticleBase::GetMomentum( par, err ) ) return 1.e10;
558 inline Double_t AliKFParticle::GetErrPt () const
561 if( AliKFParticleBase::GetPt( par, err ) ) return 1.e10;
565 inline Double_t AliKFParticle::GetErrEta () const
568 if( AliKFParticleBase::GetEta( par, err ) ) return 1.e10;
572 inline Double_t AliKFParticle::GetErrPhi () const
575 if( AliKFParticleBase::GetPhi( par, err ) ) return 1.e10;
579 inline Double_t AliKFParticle::GetErrMomentum () const
582 if( AliKFParticleBase::GetMomentum( par, err ) ) return 1.e10;
586 inline Double_t AliKFParticle::GetErrMass () const
589 if( AliKFParticleBase::GetMass( par, err ) ) return 1.e10;
593 inline Double_t AliKFParticle::GetErrDecayLength () const
596 if( AliKFParticleBase::GetDecayLength( par, err ) ) return 1.e10;
600 inline Double_t AliKFParticle::GetErrLifeTime () const
603 if( AliKFParticleBase::GetLifeTime( par, err ) ) return 1.e10;
607 inline Double_t AliKFParticle::GetErrR () const
610 if( AliKFParticleBase::GetR( par, err ) ) return 1.e10;
615 inline int AliKFParticle::GetP( Double_t &P, Double_t &SigmaP ) const
617 return AliKFParticleBase::GetMomentum( P, SigmaP );
620 inline int AliKFParticle::GetPt( Double_t &Pt, Double_t &SigmaPt ) const
622 return AliKFParticleBase::GetPt( Pt, SigmaPt );
625 inline int AliKFParticle::GetEta( Double_t &Eta, Double_t &SigmaEta ) const
627 return AliKFParticleBase::GetEta( Eta, SigmaEta );
630 inline int AliKFParticle::GetPhi( Double_t &Phi, Double_t &SigmaPhi ) const
632 return AliKFParticleBase::GetPhi( Phi, SigmaPhi );
635 inline int AliKFParticle::GetMomentum( Double_t &P, Double_t &SigmaP ) const
637 return AliKFParticleBase::GetMomentum( P, SigmaP );
640 inline int AliKFParticle::GetMass( Double_t &M, Double_t &SigmaM ) const
642 return AliKFParticleBase::GetMass( M, SigmaM );
645 inline int AliKFParticle::GetDecayLength( Double_t &L, Double_t &SigmaL ) const
647 return AliKFParticleBase::GetDecayLength( L, SigmaL );
650 inline int AliKFParticle::GetLifeTime( Double_t &T, Double_t &SigmaT ) const
652 return AliKFParticleBase::GetLifeTime( T, SigmaT );
655 inline int AliKFParticle::GetR( Double_t &R, Double_t &SigmaR ) const
657 return AliKFParticleBase::GetR( R, SigmaR );
660 inline Double_t & AliKFParticle::X()
662 return AliKFParticleBase::X();
665 inline Double_t & AliKFParticle::Y()
667 return AliKFParticleBase::Y();
670 inline Double_t & AliKFParticle::Z()
672 return AliKFParticleBase::Z();
675 inline Double_t & AliKFParticle::Px()
677 return AliKFParticleBase::Px();
680 inline Double_t & AliKFParticle::Py()
682 return AliKFParticleBase::Py();
685 inline Double_t & AliKFParticle::Pz()
687 return AliKFParticleBase::Pz();
690 inline Double_t & AliKFParticle::E()
692 return AliKFParticleBase::E();
695 inline Double_t & AliKFParticle::S()
697 return AliKFParticleBase::S();
700 inline Int_t & AliKFParticle::Q()
702 return AliKFParticleBase::Q();
705 inline Double_t & AliKFParticle::Chi2()
707 return AliKFParticleBase::Chi2();
710 inline Int_t & AliKFParticle::NDF()
712 return AliKFParticleBase::NDF();
715 inline Double_t & AliKFParticle::Parameter ( int i )
717 return AliKFParticleBase::Parameter(i);
720 inline Double_t & AliKFParticle::Covariance( int i )
722 return AliKFParticleBase::Covariance(i);
725 inline Double_t & AliKFParticle::Covariance( int i, int j )
727 return AliKFParticleBase::Covariance(i,j);
730 inline Double_t * AliKFParticle::Parameters ()
735 inline Double_t * AliKFParticle::CovarianceMatrix()
741 inline void AliKFParticle::operator +=( const AliKFParticle &Daughter )
743 AliKFParticleBase::operator +=( Daughter );
747 inline void AliKFParticle::AddDaughter( const AliKFParticle &Daughter )
749 AliKFParticleBase::AddDaughter( Daughter );
752 inline void AliKFParticle::SetProductionVertex( const AliKFParticle &Vtx )
754 AliKFParticleBase::SetProductionVertex( Vtx );
757 inline void AliKFParticle::SetMassConstraint( Double_t Mass, Double_t SigmaMass )
759 AliKFParticleBase::SetMassConstraint( Mass, SigmaMass );
762 inline void AliKFParticle::SetNoDecayLength()
764 AliKFParticleBase::SetNoDecayLength();
767 inline void AliKFParticle::Construct( const AliKFParticle *vDaughters[], int NDaughters,
768 const AliKFParticle *ProdVtx, Double_t Mass, Bool_t IsConstrained )
770 AliKFParticleBase::Construct( ( const AliKFParticleBase**)vDaughters, NDaughters,
771 ( const AliKFParticleBase*)ProdVtx, Mass, IsConstrained );
774 inline void AliKFParticle::TransportToDecayVertex()
776 AliKFParticleBase::TransportToDecayVertex();
779 inline void AliKFParticle::TransportToProductionVertex()
781 AliKFParticleBase::TransportToProductionVertex();
784 inline void AliKFParticle::TransportToPoint( const Double_t xyz[] )
786 TransportToDS( GetDStoPoint(xyz) );
789 inline void AliKFParticle::TransportToVertex( const AliVVertex &v )
791 TransportToPoint( AliKFParticle(v).fP );
794 inline void AliKFParticle::TransportToParticle( const AliKFParticle &p )
797 GetDStoParticle( p, dS, dSp );
801 inline void AliKFParticle::TransportToDS( Double_t dS )
803 AliKFParticleBase::TransportToDS( dS );
806 inline Double_t AliKFParticle::GetDStoPoint( const Double_t xyz[] ) const
808 return AliKFParticleBase::GetDStoPointBz( GetFieldAlice(), xyz );
812 inline void AliKFParticle::GetDStoParticle( const AliKFParticle &p,
813 Double_t &DS, Double_t &DSp ) const
815 GetDStoParticleXY( p, DS, DSp );
819 inline Double_t AliKFParticle::GetDistanceFromVertex( const Double_t vtx[] ) const
821 return AliKFParticleBase::GetDistanceFromVertex( vtx );
824 inline Double_t AliKFParticle::GetDeviationFromVertex( const Double_t v[],
825 const Double_t Cv[] ) const
827 return AliKFParticleBase::GetDeviationFromVertex( v, Cv);
830 inline Double_t AliKFParticle::GetDistanceFromVertex( const AliKFParticle &Vtx ) const
832 return AliKFParticleBase::GetDistanceFromVertex( Vtx );
835 inline Double_t AliKFParticle::GetDeviationFromVertex( const AliKFParticle &Vtx ) const
837 return AliKFParticleBase::GetDeviationFromVertex( Vtx );
840 inline Double_t AliKFParticle::GetDistanceFromVertex( const AliVVertex &Vtx ) const
842 return GetDistanceFromVertex( AliKFParticle(Vtx) );
845 inline Double_t AliKFParticle::GetDeviationFromVertex( const AliVVertex &Vtx ) const
847 return GetDeviationFromVertex( AliKFParticle(Vtx) );
850 inline Double_t AliKFParticle::GetDistanceFromParticle( const AliKFParticle &p ) const
852 return AliKFParticleBase::GetDistanceFromParticle( p );
855 inline Double_t AliKFParticle::GetDeviationFromParticle( const AliKFParticle &p ) const
857 return AliKFParticleBase::GetDeviationFromParticle( p );
860 inline void AliKFParticle::SubtractFromVertex( AliKFParticle &v ) const
862 AliKFParticleBase::SubtractFromVertex( v.fP, v.fC, v.fChi2, v.fNDF);
865 inline Double_t AliKFParticle::GetFieldAlice()
870 inline void AliKFParticle::GetFieldValue( const Double_t * /*xyz*/, Double_t B[] ) const
873 B[2] = GetFieldAlice();
876 inline void AliKFParticle::GetDStoParticle( const AliKFParticleBase &p,
877 Double_t &DS, Double_t &DSp )const
879 GetDStoParticleXY( p, DS, DSp );
882 inline void AliKFParticle::GetDStoParticleXY( const AliKFParticleBase &p,
883 Double_t &DS, Double_t &DSp ) const
885 AliKFParticleBase::GetDStoParticleBz( GetFieldAlice(), p, DS, DSp ) ;
886 //GetDStoParticleALICE( p, DS, DSp ) ;
889 inline void AliKFParticle::Transport( Double_t dS, Double_t P[], Double_t C[] ) const
891 AliKFParticleBase::TransportBz( GetFieldAlice(), dS, P, C );