Bug fix with matrix inversion in SetProductionVertex(); sign of GetDistanceToVertexXY...
authorbelikov <belikov@f7af4fe6-9843-0410-8265-dc069ae4e863>
Mon, 14 Jun 2010 20:05:49 +0000 (20:05 +0000)
committerbelikov <belikov@f7af4fe6-9843-0410-8265-dc069ae4e863>
Mon, 14 Jun 2010 20:05:49 +0000 (20:05 +0000)
STEER/AliKFParticle.cxx
STEER/AliKFParticle.h
STEER/AliKFParticleBase.cxx
STEER/AliKFParticleBase.h

index b65b2fc..11e1b26 100644 (file)
@@ -48,7 +48,6 @@ void AliKFParticle::Create( const Double_t Param[], const Double_t Cov[], Int_t
   AliKFParticleBase::Initialize( Param, C, Charge, mass );
 }
 
-
 AliKFParticle::AliKFParticle( const AliVTrack &track, Int_t PID )
 {
   // Constructor from ALICE track, PID hypothesis should be provided
@@ -74,7 +73,6 @@ AliKFParticle::AliKFParticle( const AliVVertex &vertex )
   fSFromDecay = 0;
 }
 
-
 void AliKFParticle::GetExternalTrackParam( const AliKFParticleBase &p, Double_t &X, Double_t &Alpha, Double_t P[5] ) 
 {
   // Conversion to AliExternalTrackParam parameterization
@@ -99,18 +97,80 @@ void AliKFParticle::GetExternalTrackParam( const AliKFParticleBase &p, Double_t
   P[4]= p.GetQ()*pti;
 }
 
+Bool_t AliKFParticle::GetDistanceFromVertexXY( const Double_t vtx[], const Double_t Cv[], Double_t &val, Double_t &err ) const
+{
+  //* Calculate DCA distance from vertex (transverse impact parameter) in XY
+  //* v = [xy], Cv=[Cxx,Cxy,Cyy ]-covariance matrix
+
+  Bool_t ret = 0;
+  
+  Double_t mP[8];
+  Double_t mC[36];
+  
+  Transport( GetDStoPoint(vtx), mP, mC );  
+
+  Double_t dx = mP[0] - vtx[0];
+  Double_t dy = mP[1] - vtx[1];
+  Double_t px = mP[3];
+  Double_t py = mP[4];
+  Double_t pt = TMath::Sqrt(px*px + py*py);
+  Double_t ex=0, ey=0;
+  if( pt<1.e-4 ){
+    ret = 1;
+    pt = 1.;
+    val = 1.e4;
+  } else{
+    ex = px/pt;
+    ey = ey/pt;
+    val = dy*ex - dx*ey;
+  }
 
+  Double_t h0 = -ey;
+  Double_t h1 = ex;
+  Double_t h3 = (dy*ey + dx*ex)*ey/pt;
+  Double_t h4 = -(dy*ey + dx*ex)*ex/pt;
+  
+  err = 
+    h0*(h0*GetCovariance(0,0) + h1*GetCovariance(0,1) + h3*GetCovariance(0,3) + h4*GetCovariance(0,4) ) +
+    h1*(h0*GetCovariance(1,0) + h1*GetCovariance(1,1) + h3*GetCovariance(1,3) + h4*GetCovariance(1,4) ) +
+    h3*(h0*GetCovariance(3,0) + h1*GetCovariance(3,1) + h3*GetCovariance(3,3) + h4*GetCovariance(3,4) ) +
+    h4*(h0*GetCovariance(4,0) + h1*GetCovariance(4,1) + h3*GetCovariance(4,3) + h4*GetCovariance(4,4) );
 
-Double_t AliKFParticle::GetDistanceFromVertexXY( const Double_t vtx[] ) const
+  if( Cv ){
+    err+= h0*(h0*Cv[0] + h1*Cv[1] ) + h1*(h0*Cv[1] + h1*Cv[2] ); 
+  }
+
+  err = TMath::Sqrt(TMath::Abs(err));
+
+  return ret;
+}
+
+Bool_t AliKFParticle::GetDistanceFromVertexXY( const Double_t vtx[], Double_t &val, Double_t &err ) const
+{
+  return GetDistanceFromVertexXY( vtx, 0, val, err );
+}
+
+
+Bool_t AliKFParticle::GetDistanceFromVertexXY( const AliKFParticle &Vtx, Double_t &val, Double_t &err ) const 
 {
   //* Calculate distance from vertex [cm] in XY-plane
 
-  Double_t mP[8], mC[36];
-  Transport( GetDStoPoint(vtx), mP, mC );
-  Double_t d[2]={ vtx[0]-mP[0], vtx[1]-mP[1] };
-  Double_t dist =  TMath::Sqrt( d[0]*d[0]+d[1]*d[1] );
-  Double_t sign = d[0]*mP[3] - d[1]*mP[4];  
-  return (sign>=0) ?dist :-dist;
+  return GetDistanceFromVertexXY( Vtx.fP, Vtx.fC, val, err );
+}
+
+Bool_t AliKFParticle::GetDistanceFromVertexXY( const AliVVertex &Vtx, Double_t &val, Double_t &err ) const 
+{
+  //* Calculate distance from vertex [cm] in XY-plane
+
+  return GetDistanceFromVertexXY( AliKFParticle(Vtx), val, err );
+}
+
+Double_t AliKFParticle::GetDistanceFromVertexXY( const Double_t vtx[] ) const
+{
+  //* Calculate distance from vertex [cm] in XY-plane
+  Double_t val, err;
+  GetDistanceFromVertexXY( vtx, 0, val, err );
+  return val;
 }
 
 Double_t AliKFParticle::GetDistanceFromVertexXY( const AliKFParticle &Vtx ) const 
@@ -165,39 +225,15 @@ Double_t AliKFParticle::GetDeviationFromParticleXY( const AliKFParticle &p ) con
 }
 
 
-Double_t AliKFParticle::GetDeviationFromVertexXY( const Double_t v[], const Double_t Cv[] ) const 
+Double_t AliKFParticle::GetDeviationFromVertexXY( const Double_t vtx[], const Double_t Cv[] ) const 
 {
   //* Calculate sqrt(Chi2/ndf) deviation from vertex
   //* v = [xyz], Cv=[Cxx,Cxy,Cyy,Cxz,Cyz,Czz]-covariance matrix
 
-  Double_t mP[8];
-  Double_t mC[36];
-  
-  Transport( GetDStoPoint(v), mP, mC );  
-
-  Double_t d[2]={ v[0]-mP[0], v[1]-mP[1] };
-
-  Double_t sigmaS = .1+10.*TMath::Sqrt( (d[0]*d[0]+d[1]*d[1] )/
-                                       (mP[3]*mP[3]+mP[4]*mP[4] )  );
-   
-  Double_t h[2] = { mP[3]*sigmaS, mP[4]*sigmaS };       
-  
-  Double_t mSi[3] = { mC[0] +h[0]*h[0], 
-                     mC[1] +h[1]*h[0], mC[2] +h[1]*h[1] };
-
-  if( Cv ){
-    mSi[0]+=Cv[0];
-    mSi[1]+=Cv[1];
-    mSi[2]+=Cv[2];
-  }
-  Double_t s = ( mSi[0]*mSi[2] - mSi[1]*mSi[1] );
-  s = ( s > 1.E-20 )  ?1./s :0;          
-  
-  Double_t mS[3] = { mSi[2], 
-                    -mSi[1], mSi[0] };      
-
-  return TMath::Sqrt( TMath::Abs(s*( ( mS[0]*d[0] + mS[1]*d[1] )*d[0]
-                                    +(mS[1]*d[0] + mS[2]*d[1] )*d[1] ))/1);
+  Double_t val, err;
+  Bool_t problem = GetDistanceFromVertexXY( vtx, Cv, val, err );
+  if( problem || err<1.e-20 ) return 1.e4;
+  else return val/err;
 }
 
 
@@ -218,7 +254,6 @@ Double_t AliKFParticle::GetDeviationFromVertexXY( const AliVVertex &Vtx ) const
   return GetDeviationFromVertexXY( v.fP, v.fC );
 }
 
-
 Double_t AliKFParticle::GetAngle  ( const AliKFParticle &p ) const 
 {
   //* Calculate the opening angle between two particles
index b2e3ec1..b3a39cc 100644 (file)
@@ -107,6 +107,7 @@ class AliKFParticle :public AliKFParticleBase
   Double_t GetMomentum    () const; //* momentum (same as GetP() )
   Double_t GetMass        () const; //* mass
   Double_t GetDecayLength () const; //* decay length
+  Double_t GetDecayLengthXY () const; //* decay length in XY
   Double_t GetLifeTime    () const; //* life time
   Double_t GetR           () const; //* distance to the origin
 
@@ -127,6 +128,7 @@ class AliKFParticle :public AliKFParticleBase
   Double_t GetErrMomentum    () const ; //* momentum
   Double_t GetErrMass        () const ; //* mass
   Double_t GetErrDecayLength () const ; //* decay length
+  Double_t GetErrDecayLengthXY () const ; //* decay length in XY
   Double_t GetErrLifeTime    () const ; //* life time
   Double_t GetErrR           () const ; //* distance to the origin
 
@@ -140,7 +142,8 @@ class AliKFParticle :public AliKFParticleBase
   int GetMomentum    ( Double_t &P, Double_t &SigmaP ) const ;   //* momentum
   int GetMass        ( Double_t &M, Double_t &SigmaM ) const ;   //* mass
   int GetDecayLength ( Double_t &L, Double_t &SigmaL ) const ;   //* decay length
-  int GetLifeTime    ( Double_t &T, Double_t &SigmaT ) const ;   //* life time
+  int GetDecayLengthXY ( Double_t &L, Double_t &SigmaL ) const ;   //* decay length in XY
+   int GetLifeTime    ( Double_t &T, Double_t &SigmaT ) const ;   //* life time
   int GetR           ( Double_t &R, Double_t &SigmaR ) const ; //* R
 
 
@@ -264,6 +267,11 @@ class AliKFParticle :public AliKFParticleBase
  
   //* Calculate distance from another object [cm] in XY-plane
 
+  Bool_t GetDistanceFromVertexXY( const Double_t vtx[], Double_t &val, Double_t &err ) const ;
+  Bool_t GetDistanceFromVertexXY( const Double_t vtx[], const Double_t Cv[], Double_t &val, Double_t &err ) const ;
+  Bool_t GetDistanceFromVertexXY( const AliKFParticle &Vtx, Double_t &val, Double_t &err ) const ;
+  Bool_t GetDistanceFromVertexXY( const AliVVertex &Vtx, Double_t &val, Double_t &err ) const ;
+
   Double_t GetDistanceFromVertexXY( const Double_t vtx[] ) const ;
   Double_t GetDistanceFromVertexXY( const AliKFParticle &Vtx ) const ;
   Double_t GetDistanceFromVertexXY( const AliVVertex &Vtx ) const ;
@@ -499,6 +507,13 @@ inline Double_t AliKFParticle::GetDecayLength () const
   else return par;
 }
 
+inline Double_t AliKFParticle::GetDecayLengthXY () const
+{
+  Double_t par, err;
+  if( AliKFParticleBase::GetDecayLengthXY( par, err ) ) return 0;
+  else return par;
+}
+
 inline Double_t AliKFParticle::GetLifeTime    () const
 {
   Double_t par, err;
@@ -602,6 +617,13 @@ inline Double_t AliKFParticle::GetErrDecayLength () const
   else return err;
 }
 
+inline Double_t AliKFParticle::GetErrDecayLengthXY () const
+{
+  Double_t par, err;
+  if( AliKFParticleBase::GetDecayLengthXY( par, err ) ) return 1.e10;
+  else return err;
+}
+
 inline Double_t AliKFParticle::GetErrLifeTime    () const
 {
   Double_t par, err;
@@ -652,6 +674,11 @@ inline int AliKFParticle::GetDecayLength( Double_t &L, Double_t &SigmaL ) const
   return AliKFParticleBase::GetDecayLength( L, SigmaL );
 }
 
+inline int AliKFParticle::GetDecayLengthXY( Double_t &L, Double_t &SigmaL ) const 
+{
+  return AliKFParticleBase::GetDecayLengthXY( L, SigmaL );
+}
+
 inline int AliKFParticle::GetLifeTime( Double_t &T, Double_t &SigmaT ) const 
 {
   return AliKFParticleBase::GetLifeTime( T, SigmaT );
@@ -851,7 +878,7 @@ inline Double_t AliKFParticle::GetDeviationFromVertex( const AliVVertex &Vtx ) c
 {
   return GetDeviationFromVertex( AliKFParticle(Vtx) );
 }
-  
 inline Double_t AliKFParticle::GetDistanceFromParticle( const AliKFParticle &p ) const 
 {
   return AliKFParticleBase::GetDistanceFromParticle( p );
@@ -903,4 +930,3 @@ inline void AliKFParticle::ConstructGamma( const AliKFParticle &daughter1,
 }
 
 #endif 
-
index 47e80bd..72cd889 100644 (file)
@@ -259,6 +259,28 @@ Int_t AliKFParticleBase::GetDecayLength( Double_t &l, Double_t &error ) const
   return 1;
 }
 
+Int_t AliKFParticleBase::GetDecayLengthXY( Double_t &l, Double_t &error ) const 
+{
+  //* Calculate particle decay length in XY projection [cm]
+
+  Double_t x = fP[3];
+  Double_t y = fP[4];
+  Double_t t = fP[7];
+  Double_t x2 = x*x;
+  Double_t y2 = y*y;
+  Double_t pt2 = x2+y2;
+  l = t*TMath::Sqrt(pt2);
+  if( pt2>1.e-4){
+    error = pt2*fC[35] + t*t/pt2*(x2*fC[9]+y2*fC[14] + 2*x*y*fC[13] )
+      + 2*t*(x*fC[31]+y*fC[32]);
+    error = TMath::Sqrt(TMath::Abs(error));
+    return 0;
+  }
+  error = 1.e20;
+  return 1;
+}
+
+
 Int_t AliKFParticleBase::GetLifeTime( Double_t &tauC, Double_t &error ) const 
 {
   //* Calculate particle decay time [s]
@@ -532,7 +554,7 @@ void AliKFParticleBase::SetProductionVertex( const AliKFParticleBase &Vtx )
 
   Double_t mAi[6];
 
-  InvertSym3( mAi, fC );
+  InvertSym3( fC, mAi );
 
   Double_t mB[5][3];
 
@@ -2087,3 +2109,4 @@ void AliKFParticleBase::MultQSQt( const Double_t Q[], const Double_t S[], Double
 
 // 72-charachters line to define the printer border
 //3456789012345678901234567890123456789012345678901234567890123456789012
+
index b6396b9..818087a 100644 (file)
@@ -110,6 +110,7 @@ class AliKFParticleBase :public TObject {
   Int_t GetPhi         ( Double_t &Phi, Double_t &SigmaPhi ) const ;
   Int_t GetMass        ( Double_t &M, Double_t &SigmaM ) const ;
   Int_t GetDecayLength ( Double_t &L, Double_t &SigmaL ) const ;
+  Int_t GetDecayLengthXY ( Double_t &L, Double_t &SigmaL ) const ;
   Int_t GetLifeTime    ( Double_t &T, Double_t &SigmaT ) const ;
   Int_t GetR           ( Double_t &R, Double_t &SigmaR ) const ;
 
@@ -268,4 +269,3 @@ class AliKFParticleBase :public TObject {
 };
 
 #endif 
-