X-Git-Url: http://git.uio.no/git/?a=blobdiff_plain;f=STEER%2FAliESDtrack.cxx;h=23e442bfe14ce4d1f01a8e8d4dcef7c8c7731ce7;hb=30bda8ef7270a200dfe2dcf55badbedb6d36a9ac;hp=75a9d6faf74b159f0a6e9f7385b913199f05a8cd;hpb=00dce61aa16f19ed5226054a59cf5fa8fc9c4532;p=u%2Fmrichter%2FAliRoot.git

diff --git a/STEER/AliESDtrack.cxx b/STEER/AliESDtrack.cxx
index 75a9d6faf74..23e442bfe14 100644
--- a/STEER/AliESDtrack.cxx
+++ b/STEER/AliESDtrack.cxx
@@ -21,6 +21,7 @@
 //-----------------------------------------------------------------
 
 #include <TMath.h>
+#include <TParticle.h>
 
 #include "AliESDVertex.h"
 #include "AliESDtrack.h"
@@ -73,9 +74,9 @@ AliESDtrack::AliESDtrack() :
   fOp(0),
   fITSchi2(0),
   fITSncls(0),
+  fITSClusterMap(0),
   fITSsignal(0),
   fITSLabel(0),
-  fITSFakeRatio(0),
   fTPCchi2(0),
   fTPCncls(0),
   fTPCnclsF(0),
@@ -96,17 +97,17 @@ AliESDtrack::AliESDtrack() :
   fTOFCalChannel(-1),
   fTOFsignal(-1),
   fTOFsignalToT(0),
-  fRICHchi2(1e10),
-  fRICHncls(0),
-  fRICHindex(0),
-  fRICHsignal(-1),
-  fRICHtheta(-1),
-  fRICHphi(-1),
-  fRICHdx(-1),
-  fRICHdy(-1),
-  fRICHmipX(-1),
-  fRICHmipY(-1),
-  fPoints(0),
+  fHMPIDchi2(1e10),
+  fHMPIDqn(-1),
+  fHMPIDcluIdx(-1),
+  fHMPIDsignal(-1),
+  fHMPIDtrkTheta(-1),
+  fHMPIDtrkPhi(-1),
+  fHMPIDtrkX(-1),
+  fHMPIDtrkY(-1),
+  fHMPIDmipX(-1),
+  fHMPIDmipY(-1),
+  fEMCALindex(kEMCALNoMatch),
   fFriendTrack(new AliESDfriendTrack())
 {
   //
@@ -120,7 +121,7 @@ AliESDtrack::AliESDtrack() :
     fTPCr[i]=1.;
     fTRDr[i]=1.;
     fTOFr[i]=1.;
-    fRICHr[i]=1.;
+    fHMPIDr[i]=1.;
   }
   
   for (i=0; i<3; i++)   { fKinkIndexes[i]=0;}
@@ -152,9 +153,9 @@ AliESDtrack::AliESDtrack(const AliESDtrack& track):
   fOp(0),
   fITSchi2(track.fITSchi2),
   fITSncls(track.fITSncls),
+  fITSClusterMap(track.fITSClusterMap),
   fITSsignal(track.fITSsignal),
   fITSLabel(track.fITSLabel),
-  fITSFakeRatio(track.fITSFakeRatio),
   fTPCchi2(track.fTPCchi2),
   fTPCncls(track.fTPCncls),
   fTPCnclsF(track.fTPCnclsF),
@@ -175,18 +176,18 @@ AliESDtrack::AliESDtrack(const AliESDtrack& track):
   fTOFCalChannel(track.fTOFCalChannel),
   fTOFsignal(track.fTOFsignal),
   fTOFsignalToT(track.fTOFsignalToT),
-  fRICHchi2(track.fRICHchi2),
-  fRICHncls(track.fRICHncls),
-  fRICHindex(track.fRICHindex),
-  fRICHsignal(track.fRICHsignal),
-  fRICHtheta(track.fRICHtheta),
-  fRICHphi(track.fRICHphi),
-  fRICHdx(track.fRICHdx),
-  fRICHdy(track.fRICHdy),
-  fRICHmipX(track.fRICHmipX),
-  fRICHmipY(track.fRICHmipY),
-  fPoints(0),
-  fFriendTrack(new AliESDfriendTrack(*(track.fFriendTrack)))
+  fHMPIDchi2(track.fHMPIDchi2),
+  fHMPIDqn(track.fHMPIDqn),
+  fHMPIDcluIdx(track.fHMPIDcluIdx),
+  fHMPIDsignal(track.fHMPIDsignal),
+  fHMPIDtrkTheta(track.fHMPIDtrkTheta),
+  fHMPIDtrkPhi(track.fHMPIDtrkPhi),
+  fHMPIDtrkX(track.fHMPIDtrkX),
+  fHMPIDtrkY(track.fHMPIDtrkY),
+  fHMPIDmipX(track.fHMPIDmipX),
+  fHMPIDmipY(track.fHMPIDmipY),
+  fEMCALindex(track.fEMCALindex),
+  fFriendTrack(0)
 {
   //
   //copy constructor
@@ -211,13 +212,182 @@ AliESDtrack::AliESDtrack(const AliESDtrack& track):
   for (Int_t i=0;i<AliPID::kSPECIES;i++) fTOFr[i]=track.fTOFr[i];
   for (Int_t i=0;i<3;i++) fTOFLabel[i]=track.fTOFLabel[i];
   for (Int_t i=0;i<10;i++) fTOFInfo[i]=track.fTOFInfo[i];
-  for (Int_t i=0;i<AliPID::kSPECIES;i++) fRICHr[i]=track.fRICHr[i];
+  for (Int_t i=0;i<AliPID::kSPECIES;i++) fHMPIDr[i]=track.fHMPIDr[i];
 
   if (track.fCp) fCp=new AliExternalTrackParam(*track.fCp);
   if (track.fIp) fIp=new AliExternalTrackParam(*track.fIp);
   if (track.fOp) fOp=new AliExternalTrackParam(*track.fOp);
 
-  if (track.fPoints) fPoints=new AliTrackPointArray(*(track.fPoints));
+  if (track.fFriendTrack) fFriendTrack=new AliESDfriendTrack(*(track.fFriendTrack));
+}
+
+//_______________________________________________________________________
+AliESDtrack::AliESDtrack(TParticle * part) : 
+  AliExternalTrackParam(),
+  fFlags(0),
+  fLabel(0),
+  fID(0),
+  fTrackLength(0),
+  fD(0),fZ(0),
+  fCdd(0),fCdz(0),fCzz(0),
+  fStopVertex(0),
+  fCp(0),
+  fCchi2(1e10),
+  fIp(0),
+  fOp(0),
+  fITSchi2(0),
+  fITSncls(0),
+  fITSClusterMap(0),
+  fITSsignal(0),
+  fITSLabel(0),
+  fTPCchi2(0),
+  fTPCncls(0),
+  fTPCnclsF(0),
+  fTPCClusterMap(159),//number of padrows
+  fTPCsignal(0),
+  fTPCsignalN(0),
+  fTPCsignalS(0),
+  fTPCLabel(0),
+  fTRDchi2(0),
+  fTRDncls(0),
+  fTRDncls0(0),
+  fTRDsignal(0),
+  fTRDLabel(0),
+  fTRDQuality(0),
+  fTRDBudget(0),
+  fTOFchi2(0),
+  fTOFindex(0),
+  fTOFCalChannel(-1),
+  fTOFsignal(-1),
+  fTOFsignalToT(0),
+  fHMPIDchi2(1e10),
+  fHMPIDqn(-1),
+  fHMPIDcluIdx(-1),
+  fHMPIDsignal(-1),
+  fHMPIDtrkTheta(-1),
+  fHMPIDtrkPhi(-1),
+  fHMPIDtrkX(-1),
+  fHMPIDtrkY(-1),
+  fHMPIDmipX(-1),
+  fHMPIDmipY(-1),
+  fEMCALindex(kEMCALNoMatch),
+  fFriendTrack(0)
+{
+  //
+  // ESD track from TParticle
+  //
+
+  // Reset all the arrays
+  Int_t i, j;
+  for (i=0; i<AliPID::kSPECIES; i++) {
+    fTrackTime[i]=0.;
+    fR[i]=0.;
+    fITSr[i]=0.;
+    fTPCr[i]=0.;
+    fTRDr[i]=0.;
+    fTOFr[i]=0.;
+    fHMPIDr[i]=0.;
+  }
+  
+  for (i=0; i<3; i++)   { fKinkIndexes[i]=0;}
+  for (i=0; i<3; i++)   { fV0Indexes[i]=-1;}
+  for (i=0;i<kNPlane;i++) {
+    for (j=0;j<kNSlice;j++) {
+      fTRDsignals[i][j]=0.; 
+    }
+    fTRDTimBin[i]=-1;
+  }
+  for (i=0;i<4;i++) {fTPCPoints[i]=-1;}
+  for (i=0;i<3;i++) {fTOFLabel[i]=-1;}
+  for (i=0;i<10;i++) {fTOFInfo[i]=-1;}
+
+  // Calculate the AliExternalTrackParam content
+
+  Double_t xref;
+  Double_t alpha;
+  Double_t param[5];
+  Double_t covar[15];
+
+  // Calculate alpha: the rotation angle of the corresponding local system (TPC sector)
+  alpha = part->Phi()*180./TMath::Pi();
+  if (alpha<0) alpha+= 360.;
+  if (alpha>360) alpha -= 360.;
+
+  Int_t sector = (Int_t)(alpha/20.);
+  alpha = 10. + 20.*sector;
+  alpha /= 180;
+  alpha *= TMath::Pi();
+
+  // Covariance matrix: no errors, the parameters are exact
+  for (Int_t i=0; i<15; i++) covar[i]=0.;
+
+  // Get the vertex of origin and the momentum
+  TVector3 ver(part->Vx(),part->Vy(),part->Vz());
+  TVector3 mom(part->Px(),part->Py(),part->Pz());
+
+  // Rotate to the local coordinate system (TPC sector)
+  ver.RotateZ(-alpha);
+  mom.RotateZ(-alpha);
+
+  // X of the referense plane
+  xref = ver.X();
+
+  Int_t pdgCode = part->GetPdgCode();
+
+  Double_t charge = 
+    TDatabasePDG::Instance()->GetParticle(pdgCode)->Charge();
+
+  param[0] = ver.Y();
+  param[1] = ver.Z();
+  param[2] = TMath::Sin(mom.Phi());
+  param[3] = mom.Pz()/mom.Pt();
+  param[4] = TMath::Sign(1/mom.Pt(),charge);
+
+  // Set AliExternalTrackParam
+  Set(xref, alpha, param, covar);
+
+  // Set the PID
+  Int_t indexPID = 99;
+
+  switch (TMath::Abs(pdgCode)) {
+
+  case  11: // electron
+    indexPID = 0;
+    break;
+
+  case 13: // muon
+    indexPID = 1;
+    break;
+
+  case 211: // pion
+    indexPID = 2;
+    break;
+
+  case 321: // kaon
+    indexPID = 3;
+    break;
+
+  case 2212: // proton
+    indexPID = 4;
+    break;
+
+  default:
+    break;
+  }
+
+  // If the particle is not e,mu,pi,K or p the PID probabilities are set to 0
+  if (indexPID < AliPID::kSPECIES) {
+    fR[indexPID]=1.;
+    fITSr[indexPID]=1.;
+    fTPCr[indexPID]=1.;
+    fTRDr[indexPID]=1.;
+    fTOFr[indexPID]=1.;
+    fHMPIDr[indexPID]=1.;
+
+  }
+  // AliESD track label
+  SetLabel(part->GetUniqueID());
+
 }
 
 //_______________________________________________________________________
@@ -230,7 +400,6 @@ AliESDtrack::~AliESDtrack(){
   delete fOp;
   delete fCp; 
   delete fFriendTrack;
-  delete fPoints;
 }
 
 void AliESDtrack::AddCalibObject(TObject * object){
@@ -259,11 +428,7 @@ void AliESDtrack::MakeMiniESDtrack(){
   // fD: Impact parameter in XY-plane
   // fZ: Impact parameter in Z 
   // fR[AliPID::kSPECIES]: combined "detector response probability"
-  // Running track parameters
-  // fRalpha: track rotation angle
-  // fRx: X-coordinate of the track reference plane 
-  // fRp[5]: external track parameters  
-  // fRc[15]: external cov. matrix of the track parameters
+  // Running track parameters in the base class (AliExternalTrackParam)
   
   fTrackLength = 0;
   for (Int_t i=0;i<AliPID::kSPECIES;i++) fTrackTime[i] = 0;
@@ -282,10 +447,10 @@ void AliESDtrack::MakeMiniESDtrack(){
   // Reset ITS track related information
   fITSchi2 = 0;
   fITSncls = 0;       
+  fITSClusterMap=0;
   fITSsignal = 0;     
   for (Int_t i=0;i<AliPID::kSPECIES;i++) fITSr[i]=0; 
   fITSLabel = 0;       
-  fITSFakeRatio = 0;   
 
   // Reset TPC related track information
   fTPCchi2 = 0;       
@@ -327,19 +492,21 @@ void AliESDtrack::MakeMiniESDtrack(){
   for (Int_t i=0;i<3;i++) fTOFLabel[i] = 0;
   for (Int_t i=0;i<10;i++) fTOFInfo[i] = 0;
 
-  // Reset RICH related track information
-  fRICHchi2 = 0;     
-  fRICHncls = 0;     
-  fRICHindex = 0;     
-  fRICHsignal = 0;     
-  for (Int_t i=0;i<AliPID::kSPECIES;i++) fRICHr[i] = 0;
-  fRICHtheta = 0;     
-  fRICHphi = 0;      
-  fRICHdx = 0;     
-  fRICHdy = 0;      
+  // Reset HMPID related track information
+  fHMPIDchi2 = 0;     
+  fHMPIDqn = -1;     
+  fHMPIDcluIdx = -1;     
+  fHMPIDsignal = 0;     
+  for (Int_t i=0;i<AliPID::kSPECIES;i++) fHMPIDr[i] = 0;
+  fHMPIDtrkTheta = -1;     
+  fHMPIDtrkPhi = -1;      
+  fHMPIDtrkX = -1;     
+  fHMPIDtrkY = -1;      
+  fHMPIDmipX = -1;
+  fHMPIDmipY = -1;
+  fEMCALindex = kEMCALNoMatch;
 
   delete fFriendTrack; fFriendTrack = 0;
-  delete fPoints; fPoints = 0;
 } 
 //_______________________________________________________________________
 Double_t AliESDtrack::GetMass() const {
@@ -382,30 +549,36 @@ Bool_t AliESDtrack::UpdateTrackParams(const AliKalmanTrack *t, ULong_t flags){
     SetIntegratedLength(t->GetIntegratedLength());
   }
 
-  Set(*t);
+  Set(t->GetX(),t->GetAlpha(),t->GetParameter(),t->GetCovariance());
   
   switch (flags) {
     
   case kITSin: case kITSout: case kITSrefit:
-    index=fFriendTrack->GetITSindices(); 
-    for (Int_t i=0;i<AliESDfriendTrack::kMaxITScluster;i++) 
-         index[i]=t->GetClusterIndex(i);
     fITSncls=t->GetNumberOfClusters();
+    index=fFriendTrack->GetITSindices(); 
+    for (Int_t i=0;i<AliESDfriendTrack::kMaxITScluster;i++) {
+        index[i]=t->GetClusterIndex(i);
+	if (i<fITSncls) {
+           Int_t l=(index[i] & 0xf0000000) >> 28;
+           SETBIT(fITSClusterMap,l);                 
+        }
+    }
     fITSchi2=t->GetChi2();
     fITSsignal=t->GetPIDsignal();
     fITSLabel = t->GetLabel();
-    fITSFakeRatio = t->GetFakeRatio();
     break;
     
   case kTPCin: case kTPCrefit:
     fTPCLabel = t->GetLabel();
     if (!fIp) fIp=new AliExternalTrackParam(*t);
-    else fIp->Set(*t);
+    else 
+      fIp->Set(t->GetX(),t->GetAlpha(),t->GetParameter(),t->GetCovariance());
   case kTPCout:
     index=fFriendTrack->GetTPCindices(); 
     if (flags & kTPCout){
       if (!fOp) fOp=new AliExternalTrackParam(*t);
-      else fOp->Set(*t);
+      else 
+        fOp->Set(t->GetX(),t->GetAlpha(),t->GetParameter(),t->GetCovariance());
     }
     fTPCncls=t->GetNumberOfClusters();    
     fTPCchi2=t->GetChi2();
@@ -474,7 +647,8 @@ Bool_t AliESDtrack::UpdateTrackParams(const AliKalmanTrack *t, ULong_t flags){
     break;
   case kTRDbackup:
     if (!fOp) fOp=new AliExternalTrackParam(*t);
-    else fOp->Set(*t);
+    else 
+      fOp->Set(t->GetX(),t->GetAlpha(),t->GetParameter(),t->GetCovariance());
     fTRDncls0 = t->GetNumberOfClusters(); 
     break;
   case kTOFin: 
@@ -793,16 +967,16 @@ void AliESDtrack::SetTOFInfo(Float_t*info) {
 
 
 //_______________________________________________________________________
-void AliESDtrack::SetRICHpid(const Double_t *p) {  
-  // Sets the probability of each particle type (in RICH)
-  SetPIDValues(fRICHr,p,AliPID::kSPECIES);
-  SetStatus(AliESDtrack::kRICHpid);
+void AliESDtrack::SetHMPIDpid(const Double_t *p) {  
+  // Sets the probability of each particle type (in HMPID)
+  SetPIDValues(fHMPIDr,p,AliPID::kSPECIES);
+  SetStatus(AliESDtrack::kHMPIDpid);
 }
 
 //_______________________________________________________________________
-void AliESDtrack::GetRICHpid(Double_t *p) const {
-  // Gets probabilities of each particle type (in RICH)
-  for (Int_t i=0; i<AliPID::kSPECIES; i++) p[i]=fRICHr[i];
+void AliESDtrack::GetHMPIDpid(Double_t *p) const {
+  // Gets probabilities of each particle type (in HMPID)
+  for (Int_t i=0; i<AliPID::kSPECIES; i++) p[i]=fHMPIDr[i];
 }
 
 
@@ -835,6 +1009,9 @@ Bool_t AliESDtrack::RelateToVertex
   //    In the case of success, the returned value is kTRUE
   //    (otherwise, it's kFALSE)
   //  
+
+  if (!vtx) return kFALSE;
+
   Double_t alpha=GetAlpha();
   Double_t sn=TMath::Sin(alpha), cs=TMath::Cos(alpha);
   Double_t x=GetX(), y=GetParameter()[0], snp=GetParameter()[2];
@@ -864,9 +1041,14 @@ Bool_t AliESDtrack::RelateToVertex
   fZ = GetParameter()[1] - zv;
   
   Double_t cov[6]; vtx->GetCovMatrix(cov);
-  fCdd = GetCovariance()[0] + cov[2];      // neglecting non-diagonals
-  fCdz = GetCovariance()[1];               //     in the vertex's    
-  fCzz = GetCovariance()[2] + cov[5];      //    covariance matrix
+
+  //***** Improvements by A.Dainese
+  alpha=GetAlpha(); sn=TMath::Sin(alpha); cs=TMath::Cos(alpha);
+  Double_t s2ylocvtx = cov[0]*sn*sn + cov[2]*cs*cs - 2.*cov[1]*cs*sn;
+  fCdd = GetCovariance()[0] + s2ylocvtx;   // neglecting correlations
+  fCdz = GetCovariance()[1];               // between (x,y) and z
+  fCzz = GetCovariance()[2] + cov[5];      // in vertex's covariance matrix
+  //*****
 
   {//Try to constrain 
     Double_t p[2]={yv,zv}, c[3]={cov[2],0.,cov[5]};
@@ -885,15 +1067,6 @@ Bool_t AliESDtrack::RelateToVertex
   return kTRUE;
 }
 
-void AliESDtrack::SetTrackPointArray(AliTrackPointArray *points) { 
-    fPoints=points;
-    //fFriendTrack->SetTrackPointArray(points); 
-}
-const AliTrackPointArray *AliESDtrack::GetTrackPointArray() const { 
-    return fPoints;
-  //return fFriendTrack->GetTrackPointArray(); 
-}
-
 //_______________________________________________________________________
 void AliESDtrack::Print(Option_t *) const {
   // Prints info on the track
@@ -929,11 +1102,60 @@ void AliESDtrack::Print(Option_t *) const {
       printf("%f, ", p[index]) ;
     printf("\n           signal = %f\n", GetTOFsignal()) ;
   }
-  if( IsOn(kRICHpid) ){
-    printf("From RICH: ") ; 
-    GetRICHpid(p) ; 
+  if( IsOn(kHMPIDpid) ){
+    printf("From HMPID: ") ; 
+    GetHMPIDpid(p) ; 
     for(index = 0 ; index < AliPID::kSPECIES; index++) 
       printf("%f, ", p[index]) ;
-    printf("\n           signal = %f\n", GetRICHsignal()) ;
+    printf("\n           signal = %f\n", GetHMPIDsignal()) ;
   }
 } 
+
+Bool_t AliESDtrack::PropagateTo(Double_t xToGo, Double_t b, Double_t mass, 
+Double_t maxStep, Bool_t rotateTo, Double_t maxSnp){
+  //----------------------------------------------------------------
+  //
+  // MI's function
+  //
+  // Propagates this track to the plane X=xk (cm) 
+  // in the magnetic field "b" (kG),
+  // the correction for the material is included
+  //
+  // mass     - mass used in propagation - used for energy loss correction
+  // maxStep  - maximal step for propagation
+  //----------------------------------------------------------------
+  const Double_t kEpsilon = 0.00001;
+  Double_t xpos     = GetX();
+  Double_t dir      = (xpos<xToGo) ? 1.:-1.;
+  //
+  while ( (xToGo-xpos)*dir > kEpsilon){
+    Double_t step = dir*TMath::Min(TMath::Abs(xToGo-xpos), maxStep);
+    Double_t x    = xpos+step;
+    Double_t xyz0[3],xyz1[3],param[7];
+    GetXYZ(xyz0);   //starting global position
+    if (!GetXYZAt(x,b,xyz1)) return kFALSE;   // no prolongation
+    xyz1[2]+=kEpsilon; // waiting for bug correction in geo
+    AliKalmanTrack::MeanMaterialBudget(xyz0,xyz1,param);	
+    if (TMath::Abs(GetSnpAt(x,b)) >= maxSnp) return kFALSE;
+    if (!AliExternalTrackParam::PropagateTo(x,b))  return kFALSE;
+
+    Double_t rho=param[0],x0=param[1],distance=param[4];
+    Double_t d=distance*rho/x0;
+
+    if (!CorrectForMaterial(d,x0,mass)) return kFALSE;
+    if (rotateTo){
+      if (TMath::Abs(GetSnp()) >= maxSnp) return kFALSE;
+      GetXYZ(xyz0);   // global position
+      Double_t alphan = TMath::ATan2(xyz0[1], xyz0[0]); 
+      //
+      Double_t ca=TMath::Cos(alphan-GetAlpha()), 
+               sa=TMath::Sin(alphan-GetAlpha());
+      Double_t sf=GetSnp(), cf=TMath::Sqrt(1.- sf*sf);
+      Double_t sinNew =  sf*ca - cf*sa;
+      if (TMath::Abs(sinNew) >= maxSnp) return kFALSE;
+      if (!Rotate(alphan)) return kFALSE;
+    }
+    xpos = GetX();
+  }
+  return kTRUE;
+}