Removing some temporary solutions for the track parameters at the primary vertex...
authorhristov <hristov@f7af4fe6-9843-0410-8265-dc069ae4e863>
Wed, 19 Nov 2003 14:04:43 +0000 (14:04 +0000)
committerhristov <hristov@f7af4fe6-9843-0410-8265-dc069ae4e863>
Wed, 19 Nov 2003 14:04:43 +0000 (14:04 +0000)
HBTAN/AliHBTReaderESD.cxx
STEER/AliESDtest.C
STEER/AliESDtrack.cxx
STEER/AliESDtrack.h
STEER/AliKalmanTrack.cxx
STEER/AliKalmanTrack.h

index 03f18847f9c161374c8af15d113fbc416ee755e8..42dc41f4593a2a135cd06efb848b5a03fee16594 100644 (file)
@@ -141,7 +141,8 @@ Int_t AliHBTReaderESD::ReadNext()
            continue;
          }
         
-        if (esdtrack->HasVertexParameters() == kFALSE)
+        //if (esdtrack->HasVertexParameters() == kFALSE) 
+        if ((esdtrack->GetStatus() & AliESDtrack::kITSrefit) == kFALSE)
          {
            if (AliHBTParticle::fgDebug > 2) 
              Info("ReadNext","Particle skipped: Data at vertex not available.");
@@ -156,8 +157,10 @@ Int_t AliHBTReaderESD::ReadNext()
          }
        
         esdtrack->GetESDpid(pidtable);
-        esdtrack->GetVertexPxPyPz(mom[0],mom[1],mom[2]);
-        esdtrack->GetVertexXYZ(pos[0],pos[1],pos[2]);
+        //esdtrack->GetVertexPxPyPz(mom[0],mom[1],mom[2]); 
+        esdtrack->GetPxPyPz(mom);
+        //esdtrack->GetVertexXYZ(pos[0],pos[1],pos[2]);
+        esdtrack->GetXYZ(pos);
         
         Double_t extx;
         Double_t extp[5];
index 8683b0355fd4dc731a7fa26937ff344c680834d1..da6ddf1dc6969590631d43f7074ea87dce105b31 100644 (file)
@@ -109,6 +109,8 @@ Int_t AliESDtest(Int_t nev=1) {
 
    //An instance of the ITS tracker
    AliITStrackerV2 itsTracker(geom);
+   {Double_t xyz[]={0.,0.,0.}, ers[]={0.005, 0.005, 0.010};
+   itsTracker.SetVertex(xyz,ers);}
    
    //An instance of the ITS PID maker
    Double_t parITS[]={34.,0.15,10.};
index 1de3543f9e372f1e3b103d95f077df37180c212c..60076dd448ea02e1bb01b7a83fd3cc76b29b0b4e 100644 (file)
@@ -39,13 +39,6 @@ fRx(0),
 fITSchi2(0),
 fITSncls(0),
 fITSsignal(0),
-fVertexX(0),
-fVertexY(0),
-fVertexZ(0),
-fVertexPx(0),
-fVertexPy(0),
-fVertexPz(0),
-fVertex(kFALSE),
 fTPCchi2(0),
 fTPCncls(0),
 fTPCsignal(0),
@@ -96,18 +89,37 @@ Bool_t AliESDtrack::UpdateTrackParams(AliKalmanTrack *t, ULong_t flags) {
   //
   // This function updates track's running parameters 
   //
+  SetStatus(flags);
+  fLabel=t->GetLabel();
+
+  if (t->IsStartedTimeIntegral()) {
+    SetStatus(kTIME);
+    Double_t times[10];t->GetIntegratedTimes(times); SetIntegratedTimes(times);
+    SetIntegratedLength(t->GetIntegratedLength());
+  }
+
+  fRalpha=t->GetAlpha();
+  t->GetExternalParameters(fRx,fRp);
+  t->GetExternalCovariance(fRc);
+
   switch (flags) {
     
-  case kITSin:
-  case kITSout: 
-  case kITSrefit:
+  case kITSin: case kITSout: case kITSrefit:
     fITSncls=t->GetNumberOfClusters();
     fITSchi2=t->GetChi2();
     for (Int_t i=0;i<fITSncls;i++) fITSindex[i]=t->GetClusterIndex(i);
     fITSsignal=t->GetPIDsignal();
     break;
     
-  case kTPCin: case kTPCout: case kTPCrefit:
+  case kTPCin: case kTPCrefit:
+    fIalpha=fRalpha;
+    fIx=fRx;
+    {
+      Int_t i;
+      for (i=0; i<5; i++) fIp[i]=fRp[i];
+      for (i=0; i<15;i++) fIc[i]=fRc[i];
+    }
+  case kTPCout:
     fTPCncls=t->GetNumberOfClusters();
     fTPCchi2=t->GetChi2();
     for (Int_t i=0;i<fTPCncls;i++) fTPCindex[i]=t->GetClusterIndex(i);
@@ -117,58 +129,19 @@ Bool_t AliESDtrack::UpdateTrackParams(AliKalmanTrack *t, ULong_t flags) {
     else if (mass<0.4) fR[2]=1.;    // the ITS reconstruction
     else fR[3]=1.;}                 //
     break;
+
   case kTRDin: case kTRDout: case kTRDrefit:
     fTRDncls=t->GetNumberOfClusters();
     fTRDchi2=t->GetChi2();
     for (Int_t i=0;i<fTRDncls;i++) fTRDindex[i]=t->GetClusterIndex(i);
     fTRDsignal=t->GetPIDsignal();
     break;
+
   default: 
     Error("UpdateTrackParams()","Wrong flag !\n");
     return kFALSE;
   }
 
-  SetStatus(flags);
-  fLabel=t->GetLabel();
-
-  if (t->IsStartedTimeIntegral()) {
-    SetStatus(kTIME);
-    Double_t times[10];t->GetIntegratedTimes(times); SetIntegratedTimes(times);
-    SetIntegratedLength(t->GetIntegratedLength());
-  }
-
-  fRalpha=t->GetAlpha();
-  t->GetExternalParameters(fRx,fRp);
-  t->GetExternalCovariance(fRc);
-  
-  if (flags == kITSin)
-   {
-     AliKalmanTrack *itstrack = t;
-     if (itstrack)
-      {
-        itstrack->PropagateTo(3.,0.0028,65.19);
-        itstrack->PropagateToVertex();
-        
-        Double_t ralpha=t->GetAlpha();
-        Double_t rx;      // X-coordinate of the track reference plane 
-        Double_t rp[5];   // external track parameters  
-        t->GetExternalParameters(rx,rp);
-   
-        Double_t phi=TMath::ASin(rp[2]) + ralpha;
-        Double_t pt=1./TMath::Abs(rp[4]);
-        Double_t r=TMath::Sqrt(rx*rx + rp[0]*rp[0]);
-        
-        fVertexX=r*TMath::Cos(phi); 
-        fVertexY=r*TMath::Sin(phi); 
-        fVertexZ=rp[1]; 
-        
-        fVertexPx = pt*TMath::Cos(phi); 
-        fVertexPy = pt*TMath::Sin(phi); 
-        fVertexPz = pt*rp[3]; 
-        fVertex = kTRUE;
-      }
-   }
-  
   return kTRUE;
 }
 
@@ -208,6 +181,26 @@ void AliESDtrack::GetXYZ(Double_t *xyz) const {
   xyz[0]=r*TMath::Cos(phi); xyz[1]=r*TMath::Sin(phi); xyz[2]=fRp[1]; 
 }
 
+void AliESDtrack::GetInnerPxPyPz(Double_t *p) const {
+  //---------------------------------------------------------------------
+  // This function returns the global track momentum components
+  // af the entrance of the TPC
+  //---------------------------------------------------------------------
+  Double_t phi=TMath::ASin(fIp[2]) + fIalpha;
+  Double_t pt=1./TMath::Abs(fIp[4]);
+  p[0]=pt*TMath::Cos(phi); p[1]=pt*TMath::Sin(phi); p[2]=pt*fIp[3]; 
+}
+
+void AliESDtrack::GetInnerXYZ(Double_t *xyz) const {
+  //---------------------------------------------------------------------
+  // This function returns the global track position
+  // af the entrance of the TPC
+  //---------------------------------------------------------------------
+  Double_t phi=TMath::ATan2(fIp[0],fIx) + fIalpha;
+  Double_t r=TMath::Sqrt(fIx*fIx + fIp[0]*fIp[0]);
+  xyz[0]=r*TMath::Cos(phi); xyz[1]=r*TMath::Sin(phi); xyz[2]=fIp[1]; 
+}
+
 //_______________________________________________________________________
 void AliESDtrack::GetExternalCovariance(Double_t c[15]) const {
   //---------------------------------------------------------------------
@@ -333,17 +326,3 @@ void AliESDtrack::GetESDpid(Double_t *p) const {
   for (Int_t i=0; i<kSPECIES; i++) p[i]=fR[i];
 }
 
-void AliESDtrack::GetVertexXYZ(Double_t& x,Double_t& y, Double_t&z) const
-{
-//returns track position in DCA to vertex  
-  x = fVertexX;
-  y = fVertexY;
-  z = fVertexZ;
-}
-void AliESDtrack::GetVertexPxPyPz(Double_t& px,Double_t& py, Double_t& pz) const
-{
-//returns track momentum in DCA to vertex  
-  px = fVertexPx;
-  py = fVertexPy;
-  pz = fVertexPz;
-}
index 798e5796736d4237c7740485081ac1abf4f30b28..79e1de9131a65dbbea871ebb27139a1b6815aa1f 100644 (file)
@@ -37,6 +37,9 @@ public:
   void GetPxPyPz(Double_t *p) const;
   void GetXYZ(Double_t *r) const;
 
+  void GetInnerPxPyPz(Double_t *p) const;
+  void GetInnerXYZ(Double_t *r) const;
+
   void SetTPCpid(const Double_t *p);
   void GetTPCpid(Double_t *p) const;
   Float_t GetTPCsignal() const {return fTPCsignal;}
@@ -61,10 +64,6 @@ public:
   UInt_t  GetTOFcluster() const {return fTOFindex;}
   void  SetTOFcluster(UInt_t index) {fTOFindex=index;}
   
-  void  GetVertexXYZ(Double_t& x,Double_t& y, Double_t&z) const;
-  void  GetVertexPxPyPz(Double_t& px,Double_t& py, Double_t& pz) const;
-  Bool_t  HasVertexParameters() const {return fVertex;}
-
   enum {
     kITSin=0x0001,kITSout=0x0002,kITSrefit=0x0004,kITSpid=0x0008,
     kTPCin=0x0010,kTPCout=0x0020,kTPCrefit=0x0040,kTPCpid=0x0080,
@@ -91,8 +90,8 @@ protected:
   Double_t fRp[5];   // external track parameters  
   Double_t fRc[15];  // external cov. matrix of the track parameters
 
-//Track parameters at the innermost measured point
-  //Double_t fIalpha,fIx,fIp[5],fIc[15];
+//Track parameters at the innermost measured point in the TPC
+  Double_t fIalpha,fIx,fIp[5],fIc[15];
 
 //Track parameters at the outermost measured point
   //Double_t fOalpha,fOx,fOp[5],fOc[15];
@@ -103,17 +102,7 @@ protected:
   UInt_t  fITSindex[6];    //! indices of the assigned ITS clusters
   Float_t fITSsignal;      // detector's PID signal
   Float_t fITSr[kSPECIES]; //! "detector response probabilities" (for the PID)
-  
-  Double_t fVertexX; // X coordinate of point of closest approach to the vertex
-  Double_t fVertexY; // Y coordinate of point of closest approach to the vertex
-  Double_t fVertexZ; // Z coordinate of point of closest approach to the vertex
-  
-  Double_t fVertexPx; // Px at point of closest approach to the vertex
-  Double_t fVertexPy; // Py at point of closest approach to the vertex
-  Double_t fVertexPz; // Pz at point of closest approach to the vertex
-  
-  Bool_t   fVertex; // TRUE if the track was prolongated to the vertex
-  
+
   // TPC related track information
   Float_t fTPCchi2;        // chi2 in the TPC
   Int_t   fTPCncls;        // number of clusters assigned in the TPC
index 7bbc416d423df42a5e8f0791a23bacd4a96fbf83..9f5092256d1288f65b3e77f283acd813a6813ace 100644 (file)
@@ -25,9 +25,6 @@
 #include "AliPDG.h"
 #include "TPDGCode.h"
 #include "TDatabasePDG.h"
-#include "AliRunLoader.h"
-#include "AliRun.h"
-#include "AliMagF.h"
 
 ClassImp(AliKalmanTrack)
 
@@ -44,8 +41,7 @@ AliKalmanTrack::AliKalmanTrack():
   // Default constructor
   //
     if (fgConvConst==0) {
-      Warning("AliKalmanTrack()", "The magnetic field has not been set!");
-      SetConvConst();
+      Fatal("AliKalmanTrack()", "The magnetic field has not been set!");
     }
     
     fStartTimeIntegral = kFALSE;
@@ -65,9 +61,8 @@ AliKalmanTrack::AliKalmanTrack(const AliKalmanTrack &t):
   // Copy constructor
   //
   if (fgConvConst==0) {
-    Warning("AliKalmanTrack(const AliKalmanTrack&)", 
+    Fatal("AliKalmanTrack(const AliKalmanTrack&)", 
            "The magnetic field has not been set!");
-    SetConvConst();
   }
 
   fStartTimeIntegral = t.fStartTimeIntegral;
@@ -77,24 +72,6 @@ AliKalmanTrack::AliKalmanTrack(const AliKalmanTrack &t):
     fIntegratedTime[i] = t.fIntegratedTime[i];
 }
 
-
-//_______________________________________________________________________
-void AliKalmanTrack::SetConvConst()
-{
-  // Sets the conversion constants for the magnetic field 
-  // (Momentum in GeV/c -> curvature in mm)
-  ::Info("SetConvConst()", "tryinig to get the magnetic field from the AliRun object..."); 
-  AliRunLoader* loader = AliRunLoader::GetRunLoader();
-  if (!loader) ::Fatal("SetConvConst()", "No run loader found"); 
-  if (!loader->GetAliRun()) loader->LoadgAlice();
-  AliRun* alirun = loader->GetAliRun();
-  if (!alirun) ::Fatal("SetConvConst()", "No AliRun object found");
-
-  Double_t field = alirun->Field()->SolenoidField();
-  SetConvConst(1000/0.299792458/field);
-  ::Info("SetConvConst()", "Magnetic field set to %f kGauss\n", field);
-}
-
 //_______________________________________________________________________
 Double_t AliKalmanTrack::GetX() const
 {
@@ -344,9 +321,6 @@ void AliKalmanTrack:: AddTimeStep(Double_t length)
     Double_t correction = TMath::Sqrt( pt*pt * (1 + tgl*tgl) + mass * mass ) / p;
     Double_t time = length * correction / kcc;
 
-    //cout << mass << "\t" << pt << "\t" << p << "\t" 
-    //     << correction << endl;
-
     fIntegratedTime[i] += time;
   }
 }
index 3339c3eb32264796dee1a1eeb739b5a1f2ee8e7f..d16aa2c632347e8dbb02da292e3cd8830aeff2f7 100644 (file)
@@ -8,8 +8,8 @@
 
 //-------------------------------------------------------------------------
 //                          Class AliKalmanTrack
-//
-//         Origin: Iouri Belikov, CERN, Jouri.Belikov@cern.ch 
+//      fixed the interface for the derived reconstructed track classes 
+//            Origin: Iouri Belikov, CERN, Jouri.Belikov@cern.ch 
 //-------------------------------------------------------------------------
 
 #include <TObject.h>
@@ -81,15 +81,12 @@ public:
   virtual Double_t P() const;
 
   virtual Double_t GetPredictedChi2(const AliCluster *) const {return 0.;}
-  virtual 
-    Int_t PropagateTo(Double_t /*xr*/, Double_t /*x0*/, Double_t /*rho*/) {return 0;}
-  virtual
-    Int_t PropagateToPrimVertex(Double_t /*x0*/,Double_t /*rho*/){return 0;}
-    
-  virtual Int_t Update(const AliCluster*, Double_t /*chi2*/, UInt_t) {return 0;}
+  virtual Int_t 
+  PropagateTo(Double_t /*xr*/, Double_t /*x0*/, Double_t /*rho*/) {return 0;}
+  virtual Int_t 
+  Update(const AliCluster*, Double_t /*chi2*/, UInt_t) {return 0;}
 
   static void SetConvConst(Double_t cc) {fgConvConst=cc;}
-  static void SetConvConst();
   static Double_t GetConvConst() {return fgConvConst;}
 
   static void SetMagneticField(Double_t f) {// f - Magnetic field in T
@@ -109,12 +106,6 @@ public:
   Double_t GetIntegratedLength() const {return fIntegratedLength;}
   void PrintTime() const;
 
-
-  //__________________
-  virtual Int_t PropagateToVertex(Double_t =0., Double_t =0.) {
-Fatal("PropagateToVertex","Not implemented!\n");return -99;}
-  
-
 protected:
   void SetChi2(Double_t chi2) {fChi2=chi2;} 
   void SetMass(Double_t mass) {fMass=mass;}