Modifications needed to use PID framework based mass during tracking and
authorshahoian <ruben.shahoyan@cern.ch>
Fri, 17 Jan 2014 16:52:18 +0000 (17:52 +0100)
committershahoian <ruben.shahoyan@cern.ch>
Fri, 17 Jan 2014 16:52:18 +0000 (17:52 +0100)
to include heavy nuclei into tracking.
The ESDtrack PID probability arrays are made transient (kept for backward compatibility
only)

20 files changed:
HMPID/AliHMPIDtrack.cxx
ITS/AliITStrackV2.cxx
ITS/AliITStrackerMI.cxx
ITS/UPGRADE/AliITSUTrackHyp.cxx
ITS/UPGRADE/AliITSUTrackerGlo.cxx
STEER/ESD/AliESDpid.cxx
STEER/ESD/AliESDpid.h
STEER/ESD/AliESDtrack.cxx
STEER/ESD/AliESDtrack.h
STEER/ESD/AliKalmanTrack.cxx
STEER/ESD/AliTrackerBase.cxx
STEER/STEER/AliReconstruction.cxx
STEER/STEER/AliReconstruction.h
STEER/STEERBase/AliExternalTrackParam.cxx
STEER/STEERBase/AliPID.cxx
STEER/STEERBase/AliPID.h
TOF/AliTOFtrack.cxx
TPC/Rec/AliTPCtrack.cxx
TRD/AliTRDtrackV1.cxx
TRD/AliTRDtrackerV1.cxx

index 451a83c..1f74210 100644 (file)
@@ -41,7 +41,7 @@ AliHMPIDtrack::AliHMPIDtrack(const AliESDtrack& t):AliKalmanTrack()
   //
   SetLabel(t.GetLabel());
   SetChi2(0.);
-  SetMass(t.GetMass());
+  SetMass(t.GetMassForTracking());
 
   Set(t.GetX(),t.GetAlpha(),t.GetParameter(),t.GetCovariance());
 
index ae71822..86a2fe5 100644 (file)
@@ -66,7 +66,7 @@ AliITStrackV2::AliITStrackV2(AliESDtrack& t,Bool_t c):
   Set(par->GetX(),par->GetAlpha(),par->GetParameter(),par->GetCovariance());
 
   SetLabel(t.GetITSLabel());
-  SetMass(t.GetMass());
+  SetMass(t.GetMassForTracking());
   SetNumberOfClusters(t.GetITSclusters(fIndex));
 
   if (t.GetStatus()&AliESDtrack::kTIME) {
@@ -456,6 +456,7 @@ Bool_t AliITStrackV2::Improve(Double_t x0,Double_t xyz[3],Double_t ers[3]) {
   Double_t dx = x - xv, dy = par[0] - yv, dz = par[1] - zv;
   Double_t r2=dx*dx + dy*dy;
   Double_t p2=(1.+ GetTgl()*GetTgl())/(GetSigned1Pt()*GetSigned1Pt());
+  if (GetMass()<0) p2 *= 4; // q=2
   Double_t beta2=p2/(p2 + GetMass()*GetMass());
   x0*=TMath::Sqrt((1.+ GetTgl()*GetTgl())/(1.- GetSnp()*GetSnp()));
   Double_t theta2=14.1*14.1/(beta2*p2*1e6)*x0;
@@ -690,6 +691,7 @@ Bool_t AliITStrackV2::ImproveKalman(Double_t xyz[3],Double_t ers[3], const Doubl
   double ms44t = p34*p34;
   //
   double p2=(1.+ par[3]*par[3])/(par[4]*par[4]);
+  if (GetMass()<0) p2 *= 4; // q=2
   double beta2 = p2/(p2+GetMass()*GetMass());
   double theta2t = 14.1*14.1/(beta2*p2*1e6) * (1. + par[3]*par[3]);
   //
index e687742..15cc909 100644 (file)
@@ -563,7 +563,6 @@ Int_t AliITStrackerMI::Clusters2Tracks(AliESDEvent *event) {
   // temporary
   Int_t noesd = 0;
   {/* Read ESD tracks */
-    Double_t pimass = TDatabasePDG::Instance()->GetParticle(211)->Mass();
     Int_t nentr=event->GetNumberOfTracks();
     noesd=nentr;
     //    Info("Clusters2Tracks", "Number of ESD tracks: %d\n", nentr);
@@ -580,9 +579,6 @@ Int_t AliITStrackerMI::Clusters2Tracks(AliESDEvent *event) {
       t->GetDZ(GetX(),GetY(),GetZ(),t->GetDP());              //I.B.
       Double_t vdist = TMath::Sqrt(t->GetD(0)*t->GetD(0)+t->GetD(1)*t->GetD(1));
 
-
-      // look at the ESD mass hypothesys !
-      if (t->GetMass()<0.9*pimass) t->SetMass(pimass); 
       // write expected q
       t->SetExpQ(TMath::Max(0.8*t->GetESDtrack()->GetTPCsignal(),30.));
 
index fcbf4ce..28b4bd1 100644 (file)
@@ -81,7 +81,7 @@ AliITSUTrackHyp::AliITSUTrackHyp(const AliESDtrack &src)
 {
   // copy c-tor. from ESD track: take only kinematics, mass and time integral
   AliExternalTrackParam::operator=((const AliExternalTrackParam&)src);
-  SetMass(src.GetMass());
+  SetMass(src.GetMassForTracking());
   if (src.IsOn(AliESDtrack::kTIME)) {
     StartTimeIntegral();
     SetIntegratedLength(src.GetIntegratedLength());
@@ -130,7 +130,7 @@ AliITSUTrackHyp &AliITSUTrackHyp::operator=(const AliESDtrack &src)
 {
   // copy oparator from ESD track: take only kinematics, mass and time integral
   AliExternalTrackParam::operator=((const AliExternalTrackParam&)src);
-  SetMass(src.GetMass());
+  SetMass(src.GetMassForTracking());
   if (src.IsOn(AliESDtrack::kTIME)) {
     StartTimeIntegral();
     SetIntegratedLength(src.GetIntegratedLength());
index c42c13c..fc843fc 100644 (file)
@@ -263,7 +263,7 @@ Int_t AliITSUTrackerGlo::Clusters2Tracks(AliESDEvent *esdEv)
       AliLog::SetClassDebugLevel("AliITSUTrackerGlo",dbg ? 10:0);
       */
 #ifdef _ITSU_DEBUG_
-      AliDebug(1,Form("Processing track %d(esd%d) | M=%.3f Pt=%.3f | MCLabel: %d",itr,trID,fCurrESDtrack->GetMass(kTRUE),fCurrESDtrack->Pt(),fCurrESDtrMClb));//RS
+      AliDebug(1,Form("Processing track %d(esd%d) | M=%.3f Pt=%.3f | MCLabel: %d",itr,trID,fCurrESDtrack->GetMassForTracking(kTRUE),fCurrESDtrack->Pt(),fCurrESDtrMClb));//RS
 #endif
       FindTrack(fCurrESDtrack, trID);
     }   
@@ -680,8 +680,7 @@ AliITSUTrackHyp* AliITSUTrackerGlo::InitHypothesis(AliESDtrack *esdTr, Int_t esd
   //
   fCountProlongationTrials++;
   //
-  fCurrMass = esdTr->GetMass();
-  if (fCurrMass<kPionMass*0.9) fCurrMass = kPionMass; // don't trust to mu, e identification from TPCin
+  fCurrMass = esdTr->GetMassForTracking();
   //
   hyp = new AliITSUTrackHyp(fNLrActive);
   hyp->SetESDTrack(esdTr);
index 89f6d92..66020a6 100644 (file)
@@ -467,3 +467,29 @@ Float_t AliESDpid::GetNumberOfSigmasTOFold(const AliVParticle *track, AliPID::EP
   Double_t expTime = fTOFResponse.GetExpectedSignal(vtrack,type);
   return (tofTime - fTOFResponse.GetStartTime(vtrack->P()) - expTime)/fTOFResponse.GetExpectedSigma(vtrack->P(),expTime,AliPID::ParticleMassZ(type));
 }
+
+//_________________________________________________________________________
+void AliESDpid::SetMassForTracking(AliESDtrack *esdtr) const
+{
+  // assign mass for tracking
+  //
+  int pid = AliPID::kPion; // this should be substituted by real most probable TPC pid (e,mu -> pion) or poin if no PID possible
+  //
+  if (pid<AliPID::kPion || pid>AliPID::kSPECIESC-1) pid = AliPID::kPion;
+  //
+  esdtr->SetMassForTracking( AliPID::ParticleCharge(pid)==1 ? AliPID::ParticleMass(pid) : -AliPID::ParticleMass(pid));
+  //
+}
+
+
+//_________________________________________________________________________
+void AliESDpid::MakePIDForTracking(AliESDEvent *event) const
+{
+  // assign masses using for tracking
+  Int_t nTrk=event->GetNumberOfTracks();
+  for (Int_t iTrk=nTrk; iTrk--;) {  
+    AliESDtrack *track = event->GetTrack(iTrk);
+    SetMassForTracking(track);
+  }
+
+}
index 2c16383..6e1d527 100644 (file)
@@ -32,6 +32,8 @@ AliESDpid& operator=(const AliESDpid& a){if (this==&a) return *this; AliPIDRespo
   virtual ~AliESDpid() {}
   
   Int_t MakePID(AliESDEvent *event, Bool_t TPCOnly = kFALSE, Float_t timeZeroTOF=9999) const;
+  void  MakePIDForTracking(AliESDEvent *event) const;
+
   void MakeTPCPID(AliESDtrack *track) const;
   void MakeITSPID(AliESDtrack *track) const;
   void MakeTOFPID(AliESDtrack *track, Float_t /*timeZeroTOF*/) const;
@@ -40,6 +42,8 @@ AliESDpid& operator=(const AliESDpid& a){if (this==&a) return *this; AliPIDRespo
   void MakeTRDPID(AliESDtrack *track) const;
   void CombinePID(AliESDtrack *track) const;
 
+  void SetMassForTracking(AliESDtrack *track) const;
+
 //   Float_t NumberOfSigmasTOF(const AliVParticle *track, AliPID::EParticleType type) const {return AliPIDResponse::NumberOfSigmasTOF(track,type);}
 //   Float_t GetNumberOfSigmasTOF(const AliVParticle *track, AliPID::EParticleType type, const Float_t timeZeroTOF) const;
   
index 0e228f7..da01af8 100644 (file)
@@ -185,6 +185,7 @@ AliESDtrack::AliESDtrack() :
   fHMPIDqn(0),
   fHMPIDcluIdx(-1),
   fCaloIndex(kEMCALNoMatch),
+  fMassForTracking(0.13957),
   fHMPIDtrkTheta(0),
   fHMPIDtrkPhi(0),
   fHMPIDsignal(0),
@@ -259,14 +260,7 @@ AliESDtrack::AliESDtrack() :
 
   Int_t i;
   for (i=kNITSchi2Std;i--;) fITSchi2Std[i] = 0;
-  for (i=0; i<AliPID::kSPECIES; i++) {
-    fR[i]=0.;
-    fITSr[i]=0.;
-    fTPCr[i]=0.;
-    fTRDr[i]=0.;
-    fTOFr[i]=0.;
-    fHMPIDr[i]=0.;
-  }
+  for (i=0; i<AliPID::kSPECIES; i++) fR[i]=fITSr[i]=fTPCr[i]=fTRDr[i]=fTOFr[i]=fHMPIDr[i]=0.;
   
   for (i=0; i<3; i++)   { fKinkIndexes[i]=0;}
   for (i=0; i<3; i++)   { fV0Indexes[i]=0;}
@@ -305,6 +299,7 @@ AliESDtrack::AliESDtrack(const AliESDtrack& track):
   fHMPIDqn(track.fHMPIDqn),
   fHMPIDcluIdx(track.fHMPIDcluIdx),
   fCaloIndex(track.fCaloIndex),
+  fMassForTracking(track.fMassForTracking),
   fHMPIDtrkTheta(track.fHMPIDtrkTheta),
   fHMPIDtrkPhi(track.fHMPIDtrkPhi),
   fHMPIDsignal(track.fHMPIDsignal),
@@ -454,6 +449,7 @@ AliESDtrack::AliESDtrack(const AliVTrack *track) :
   fHMPIDqn(0),
   fHMPIDcluIdx(-1),
   fCaloIndex(kEMCALNoMatch),
+  fMassForTracking(0.13957),
   fHMPIDtrkTheta(0),
   fHMPIDtrkPhi(0),
   fHMPIDsignal(0),
@@ -536,14 +532,7 @@ AliESDtrack::AliESDtrack(const AliVTrack *track) :
   // Reset all the arrays
   Int_t i;
   for (i=kNITSchi2Std;i--;) fITSchi2Std[i] = 0;
-  for (i=0; i<AliPID::kSPECIES; i++) {
-    fR[i]=0.;
-    fITSr[i]=0.;
-    fTPCr[i]=0.;
-    fTRDr[i]=0.;
-    fTOFr[i]=0.;
-    fHMPIDr[i]=0.;
-  }
+  for (i=0; i<AliPID::kSPECIES; i++) fR[i]=fITSr[i]=fTPCr[i]=fTRDr[i]=fTOFr[i]=fHMPIDr[i]=0.;
   
   for (i=0; i<3; i++)   { fKinkIndexes[i]=0;}
   for (i=0; i<3; i++)   { fV0Indexes[i]=-1;}
@@ -647,6 +636,7 @@ AliESDtrack::AliESDtrack(TParticle * part) :
   fHMPIDqn(0),
   fHMPIDcluIdx(-1),
   fCaloIndex(kEMCALNoMatch),
+  fMassForTracking(0.13957),
   fHMPIDtrkTheta(0),
   fHMPIDtrkPhi(0),
   fHMPIDsignal(0),
@@ -721,14 +711,7 @@ AliESDtrack::AliESDtrack(TParticle * part) :
   // Reset all the arrays
   Int_t i;
   for (i=kNITSchi2Std;i--;) fITSchi2Std[i] = 0;
-  for (i=0; i<AliPID::kSPECIES; i++) {
-    fR[i]=0.;
-    fITSr[i]=0.;
-    fTPCr[i]=0.;
-    fTRDr[i]=0.;
-    fTOFr[i]=0.;
-    fHMPIDr[i]=0.;
-  }
+  for (i=0; i<AliPID::kSPECIES; i++) fR[i]=fITSr[i]=fTPCr[i]=fTRDr[i]=fTOFr[i]=fHMPIDr[i]=0.;
   
   for (i=0; i<3; i++)   { fKinkIndexes[i]=0;}
   for (i=0; i<3; i++)   { fV0Indexes[i]=-1;}
@@ -787,34 +770,10 @@ AliESDtrack::AliESDtrack(TParticle * part) :
 
   // Set the PID
   Int_t indexPID = 99;
+  if (pdgCode<0) pdgCode = -pdgCode;
+  for (i=0;i<AliPID::kSPECIES;i++) if (pdgCode==AliPID::ParticleCode(i)) {indexPID = i; break;}
 
-  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 the particle is not valid charged the PID probabilities are set to 0
   if (indexPID < AliPID::kSPECIES) {
     fR[indexPID]=1.;
     fITSr[indexPID]=1.;
@@ -979,6 +938,8 @@ AliESDtrack &AliESDtrack::operator=(const AliESDtrack &source){
     fTOFr[i]  = source.fTOFr[i];
     fHMPIDr[i] = source.fHMPIDr[i];
   }
+  
+  fMassForTracking = source.fMassForTracking;
 
   fHMPIDtrkTheta = source.fHMPIDtrkTheta;
   fHMPIDtrkPhi   = source.fHMPIDtrkPhi;
@@ -1317,6 +1278,7 @@ void AliESDtrack::MakeMiniESDtrack(){
   fHMPIDcluIdx = -1;     
   fHMPIDsignal = 0;     
   for (Int_t i=0;i<AliPID::kSPECIES;i++) fHMPIDr[i] = 0;
+  fMassForTracking = 0.13957;
   fHMPIDtrkTheta = 0;     
   fHMPIDtrkPhi = 0;      
   fHMPIDtrkX = 0;     
index 97b86e4..3446802 100644 (file)
@@ -91,6 +91,8 @@ public:
   Int_t    GetPID(Bool_t tpcOnly=kFALSE)  const;
   Int_t    GetTOFBunchCrossing(Double_t b=0, Bool_t pidTPConly=kTRUE) const;
   Double_t GetMass(Bool_t tpcOnly=kFALSE) const {return AliPID::ParticleMass(GetPID(tpcOnly));}
+  Double_t GetMassForTracking() const {return fMassForTracking;}
+  void     SetMassForTracking(Double_t m) {fMassForTracking = m;}
   Double_t M() const;
   Double_t E() const;
   Double_t Y() const;
@@ -470,12 +472,13 @@ protected:
   Int_t     fKinkIndexes[3]; // array of indexes of posible kink candidates 
   Int_t     fV0Indexes[3];   // array of indexes of posible kink candidates 
 
-  Double32_t   fR[AliPID::kSPECIES]; //[0.,0.,8] combined "detector response probability"
-  Double32_t   fITSr[AliPID::kSPECIES]; //[0.,0.,8] "detector response probabilities" (for the PID)
-  Double32_t   fTPCr[AliPID::kSPECIES]; //[0.,0.,8] "detector response probabilities" (for the PID)
-  Double32_t   fTRDr[AliPID::kSPECIES]; //[0.,0.,8] "detector response probabilities" (for the PID)  
-  Double32_t   fTOFr[AliPID::kSPECIES]; //[0.,0.,8] "detector response probabilities" (for the PID)
-  Double32_t   fHMPIDr[AliPID::kSPECIES];//[0.,0.,8] "detector response probabilities" (for the PID)
+  Double32_t   fR[AliPID::kSPECIES]; //! [0.,0.,8] combined "detector response probability"
+  Double32_t   fITSr[AliPID::kSPECIES]; //! [0.,0.,8] "detector response probabilities" (for the PID)
+  Double32_t   fTPCr[AliPID::kSPECIES]; //! [0.,0.,8] "detector response probabilities" (for the PID)
+  Double32_t   fTRDr[AliPID::kSPECIES]; //! [0.,0.,8] "detector response probabilities" (for the PID)  
+  Double32_t   fTOFr[AliPID::kSPECIES]; //! [0.,0.,8] "detector response probabilities" (for the PID)
+  Double32_t   fHMPIDr[AliPID::kSPECIES];//! [0.,0.,8] "detector response probabilities" (for the PID)
+  Double32_t   fMassForTracking;         // mass used for tracking
 
   Double32_t fHMPIDtrkTheta;//[-2*pi,2*pi,16] theta of the track extrapolated to the HMPID, LORS
   // how much of this is needed?
@@ -579,7 +582,7 @@ protected:
   static bool fgkOnlineMode; //! indicate the online mode to skip some of the functionality
 
   AliESDtrack & operator=(const AliESDtrack & );
-  ClassDef(AliESDtrack,69)  //ESDtrack 
+  ClassDef(AliESDtrack,70)  //ESDtrack 
 };
 
 
index 70aef79..4b47a35 100644 (file)
@@ -124,21 +124,18 @@ void AliKalmanTrack:: AddTimeStep(Double_t length)
   fIntegratedLength += length;
 
   Double_t xr, param[5];
-  Double_t pt, tgl;
   
   GetExternalParameters(xr, param);
-  pt =  1/param[4] ;
-  tgl = param[3];
+  double tgl = param[3];
 
-  Double_t p = TMath::Abs(pt * TMath::Sqrt(1+tgl*tgl));
+  Double_t p2inv = param[4]*param[4]/(1+tgl*tgl);
 
   //  if (length > 100) return;
 
   for (Int_t i=0; i<AliPID::kSPECIESC; i++) {
     
-    Double_t mass = AliPID::ParticleMass(i);
-    Double_t pCharge = p * AliPID::ParticleCharge(i);
-    Double_t correction = TMath::Sqrt( pCharge * pCharge + mass * mass ) / pCharge;
+    Double_t massz = AliPID::ParticleMassZ(i);
+    Double_t correction = TMath::Sqrt( 1. + massz*massz*p2inv ); // 1/beta
     Double_t time = length * correction / kcc;
 
     fIntegratedTime[i] += time;
index 3073494..c7811aa 100644 (file)
@@ -274,7 +274,7 @@ AliTrackerBase::PropagateTrackTo(AliExternalTrackParam *track, Double_t xToGo,
   // Propagates the track to the plane X=xk (cm) using the magnetic field map 
   // and correcting for the crossed material.
   //
-  // mass     - mass used in propagation - used for energy loss correction
+  // mass     - mass used in propagation - used for energy loss correction (if <0 then q = 2)
   // maxStep  - maximal step for propagation
   //
   //  Origin: Marian Ivanov,  Marian.Ivanov@cern.ch
@@ -343,7 +343,7 @@ AliTrackerBase::PropagateTrackToBxByBz(AliExternalTrackParam *track,
   // taking into account all the three components of the magnetic field 
   // and correcting for the crossed material.
   //
-  // mass     - mass used in propagation - used for energy loss correction
+  // mass     - mass used in propagation - used for energy loss correction (if <0 then q=2)
   // maxStep  - maximal step for propagation
   //
   //  Origin: Marian Ivanov,  Marian.Ivanov@cern.ch
index d1913e1..d55d7c2 100644 (file)
@@ -300,6 +300,8 @@ AliReconstruction::AliReconstruction(const char* gAliceFilename) :
   fWriteQAExpertData(kTRUE), 
   fRunPlaneEff(kFALSE),
 
+  fESDpid(NULL),
+
   fesd(NULL),
   fhltesd(NULL),
   fesdf(NULL),
@@ -431,6 +433,8 @@ AliReconstruction::AliReconstruction(const AliReconstruction& rec) :
   fWriteQAExpertData(rec.fWriteQAExpertData), 
   fRunPlaneEff(rec.fRunPlaneEff),
 
+  fESDpid(NULL),
+
   fesd(NULL),
   fhltesd(NULL),
   fesdf(NULL),
@@ -610,6 +614,7 @@ AliReconstruction& AliReconstruction::operator = (const AliReconstruction& rec)
   fWriteQAExpertData           = rec.fWriteQAExpertData;
   fRunPlaneEff                 = rec.fRunPlaneEff;
   for (int i=2;i--;) for (int j=2;j--;) fBeamInt[i][j] = rec.fBeamInt[i][j];
+  fESDpid  = NULL;
   fesd     = NULL;
   fhltesd  = NULL;
   fesdf    = NULL;
@@ -1826,6 +1831,9 @@ void AliReconstruction::SlaveBegin(TTree*)
   gSystem->GetProcInfo(&procInfo);
   AliInfo(Form("Current memory usage %ld %ld", procInfo.fMemResident, procInfo.fMemVirtual));
   
+  // PID
+  fESDpid = new AliESDpid();
+
   //QA
   //Initialize the QA and start of cycle 
   if (fRunQA || fRunGlobalQA) 
@@ -1915,8 +1923,6 @@ Bool_t AliReconstruction::ProcessEvent(Int_t iEvent)
 
   AliCodeTimerAuto("",0);
 
-  AliESDpid pid;
-
   AliSysInfo::AddStamp(Form("StartEv_%d",iEvent), 0,0,iEvent);
 
   if (iEvent >= fRunLoader->GetNumberOfEvents()) {
@@ -1961,7 +1967,7 @@ Bool_t AliReconstruction::ProcessEvent(Int_t iEvent)
       if (reconstructor && fRecoParam.GetDetRecoParamArray(iDet)) {
         const AliDetectorRecoParam *par = fRecoParam.GetDetRecoParam(iDet);
         reconstructor->SetRecoParam(par);
-       reconstructor->GetPidSettings(&pid);
+       reconstructor->GetPidSettings(fESDpid);
        reconstructor->SetEventInfo(&fEventInfo);
         if (fRunQA) {
           AliQAManager::QAManager()->SetEventInfo(&fEventInfo) ;
@@ -2122,7 +2128,7 @@ Bool_t AliReconstruction::ProcessEvent(Int_t iEvent)
 
     // barrel tracking
     if (!fRunTracking.IsNull()) {
-      if (!RunTracking(fesd,pid)) {
+      if (!RunTracking(fesd,*fESDpid)) {
        if (fStopOnError) {CleanUp(); return kFALSE;}
       }
     }
@@ -2166,7 +2172,7 @@ Bool_t AliReconstruction::ProcessEvent(Int_t iEvent)
       ok = kFALSE;
       if (tpcTrack)
        ok = AliTracker::
-         PropagateTrackToBxByBz(tpcTrack,kRadius,track->GetMass(),kMaxStep,kFALSE);
+         PropagateTrackToBxByBz(tpcTrack,kRadius,track->GetMassForTracking(),kMaxStep,kFALSE);
 
       if (ok) {
        Int_t n=trkArray.GetEntriesFast();
@@ -2178,7 +2184,7 @@ Bool_t AliReconstruction::ProcessEvent(Int_t iEvent)
       if (track->IsOn(AliESDtrack::kITSrefit)) continue;
 
       AliTracker::
-         PropagateTrackToBxByBz(track,kRadius,track->GetMass(),kMaxStep,kFALSE);
+         PropagateTrackToBxByBz(track,kRadius,track->GetMassForTracking(),kMaxStep,kFALSE);
       Double_t x[3]; track->GetXYZ(x);
       Double_t b[3]; AliTracker::GetBxByBz(x,b);
       track->RelateToVertexBxByBz(fesd->GetPrimaryVertexSPD(), b, kVeryBig);
@@ -2303,10 +2309,10 @@ Bool_t AliReconstruction::ProcessEvent(Int_t iEvent)
 
     // AdC+FN
     if (fReconstructor[3])
-      GetReconstructor(3)->FillEventTimeWithTOF(fesd,&pid);
+      GetReconstructor(3)->FillEventTimeWithTOF(fesd,fESDpid);
 
     // combined PID
-    pid.MakePID(fesd);
+    //    fESDpid->MakePID(fesd);
 
     if (fFillTriggerESD) {
       if (!FillTriggerESD(fesd)) {
@@ -2972,7 +2978,7 @@ Bool_t AliReconstruction::RunTracking(AliESDEvent*& esd,AliESDpid &PID)
     // preliminary PID in TPC needed by the ITS tracker
     if (iDet == 1) {
       GetReconstructor(1)->FillESD((TTree*)NULL, (TTree*)NULL, esd);
-      PID.MakePID(esd,kTRUE);
+      PID.MakePIDForTracking(esd);
       AliSysInfo::AddStamp(Form("MakePID0%s_%d",fgkDetectorName[iDet],eventNr), iDet,4,eventNr);
     } 
   }
@@ -3027,7 +3033,7 @@ Bool_t AliReconstruction::RunTracking(AliESDEvent*& esd,AliESDpid &PID)
     if (iDet == 1) {
       //GetReconstructor(1)->FillESD((TTree*)NULL, (TTree*)NULL, esd);
       //AliESDpid::MakePID(esd);
-      PID.MakePID(esd,kTRUE);
+      PID.MakePIDForTracking(esd);
       AliSysInfo::AddStamp(Form("MakePID1%s_%d",fgkDetectorName[iDet],eventNr), iDet,4,eventNr);
     }
 
@@ -3559,6 +3565,9 @@ void AliReconstruction::CleanUp()
   delete fParentRawReader;
   fParentRawReader=NULL;
 
+  delete fESDpid;
+  fESDpid = NULL;
+
   if (ffile) {
     ffile->Close();
     delete ffile;
index 6946a67..f10ef53 100644 (file)
@@ -363,6 +363,9 @@ private:
   // Plane Efficiency Evaluation
   Bool_t         fRunPlaneEff ;      // Evaluate Plane Efficiency
 
+  // PID
+  AliESDpid*     fESDpid;                      // PID object
+
   // New members needed in order to split Run method
   // into InitRun,RunEvent,FinishRun methods
   AliESDEvent*         fesd;        //! Pointer to the ESD event object
@@ -407,7 +410,7 @@ private:
   Int_t                fMaxVMEM;        //  max VMEM memory, MB
   static const char*   fgkStopEvFName;  //  filename for stop.event stamp
   //
-  ClassDef(AliReconstruction, 46)      // class for running the reconstruction
+  ClassDef(AliReconstruction, 47)      // class for running the reconstruction
 };
 
 #endif
index 0f5bca4..950f1d5 100644 (file)
@@ -418,7 +418,7 @@ Bool_t AliExternalTrackParam::CorrectForMeanMaterialdEdx
   // "xTimesRho" - is the product length*density (g/cm^2).
   //     It should be passed as negative when propagating tracks 
   //     from the intreaction point to the outside of the central barrel. 
-  // "mass" - the mass of this particle (GeV/c^2).
+  // "mass" - the mass of this particle (GeV/c^2). Negative mass means charge=2 particle
   // "dEdx" - mean enery loss (GeV/(g/cm^2)
   // "anglecorr" - switch for the angular correction
   //------------------------------------------------------------------
@@ -439,6 +439,7 @@ Bool_t AliExternalTrackParam::CorrectForMeanMaterialdEdx
   } 
 
   Double_t p=GetP();
+  if (mass<0) p += p; // q=2 particle 
   Double_t p2=p*p;
   Double_t beta2=p2/(p2 + mass*mass);
 
@@ -454,6 +455,7 @@ Bool_t AliExternalTrackParam::CorrectForMeanMaterialdEdx
       double lt = 1+0.038*TMath::Log(TMath::Abs(xOverX0));
       if (lt>0) theta2 *= lt*lt;
     }
+    if (mass<0) theta2 *= 4; // q=2 particle
     if(theta2>TMath::Pi()*TMath::Pi()) return kFALSE;
     cC22 = theta2*((1.-fP2)*(1.+fP2))*(1. + fP3*fP3);
     cC33 = theta2*(1. + fP3*fP3)*(1. + fP3*fP3);
@@ -467,7 +469,6 @@ Bool_t AliExternalTrackParam::CorrectForMeanMaterialdEdx
      Double_t dE=dEdx*xTimesRho;
      Double_t e=TMath::Sqrt(p2 + mass*mass);
      if ( TMath::Abs(dE) > 0.3*e ) return kFALSE; //30% energy loss is too much!
-     //cP4 = (1.- e/p2*dE);
      if ( (1.+ dE/p2*(dE + 2*e)) < 0. ) return kFALSE;
      cP4 = 1./TMath::Sqrt(1.+ dE/p2*(dE + 2*e));  //A precise formula by Ruben !
      if (TMath::Abs(fP4*cP4)>100.) return kFALSE; //Do not track below 10 MeV/c
@@ -502,16 +503,20 @@ Bool_t AliExternalTrackParam::CorrectForMeanMaterial
   // "xTimesRho" - is the product length*density (g/cm^2). 
   //     It should be passed as negative when propagating tracks 
   //     from the intreaction point to the outside of the central barrel. 
-  // "mass" - the mass of this particle (GeV/c^2).
+  // "mass" - the mass of this particle (GeV/c^2). mass<0 means charge=2
   // "anglecorr" - switch for the angular correction
   // "Bethe" - function calculating the energy loss (GeV/(g/cm^2)) 
   //------------------------------------------------------------------
   Double_t bg=GetP()/mass;
-  if (bg<kAlmost0) {
-    AliDebug(2,Form("Non-positive beta*gamma = %e, mass = %e",bg,mass));
-    return kFALSE;
+  if (mass<0) {
+    if (mass<-990) {
+      AliDebug(2,Form("Mass %f corresponds to unknown PID particle = %e",mass));
+      return kFALSE;
+    }
+    bg = -2*bg;
   }
   Double_t dEdx=Bethe(bg);
+  if (mass<0) dEdx *= 4;
 
   return CorrectForMeanMaterialdEdx(xOverX0,xTimesRho,mass,dEdx,anglecorr);
 }
@@ -531,7 +536,7 @@ Bool_t AliExternalTrackParam::CorrectForMeanMaterialZA
   // "xTimesRho" - is the product length*density (g/cm^2). 
   //     It should be passed as negative when propagating tracks 
   //     from the intreaction point to the outside of the central barrel. 
-  // "mass" - the mass of this particle (GeV/c^2).
+  // "mass" - the mass of this particle (GeV/c^2). mass<0 means charge=2 particle
   // "density"  - mean density (g/cm^3)
   // "zOverA"   - mean Z/A
   // "exEnergy" - mean exitation energy (GeV)
@@ -544,8 +549,15 @@ Bool_t AliExternalTrackParam::CorrectForMeanMaterialZA
   //------------------------------------------------------------------
 
   Double_t bg=GetP()/mass;
+  if (mass<0) {
+    if (mass<-990) {
+      AliDebug(2,Form("Mass %f corresponds to unknown PID particle = %e",mass));
+      return kFALSE;
+    }
+    bg = -2*bg;
+  }
   Double_t dEdx=BetheBlochGeant(bg,density,jp1,jp2,exEnergy,zOverA);
-
+  if (mass<0) dEdx *= 4;
   return CorrectForMeanMaterialdEdx(xOverX0,xTimesRho,mass,dEdx,anglecorr);
 }
 
index 96d2030..ecc4180 100644 (file)
@@ -174,6 +174,8 @@ const Int_t AliPID::fgkParticleCode[AliPID::kSPECIESCN+1] = {
   */
 };
 
+Char_t AliPID::fgkParticleCharge[AliPID::kSPECIESCN+1] = { 0,0,0,0,0,0,0,0,0,0,0,0,0,0,0 };
+
 Double_t AliPID::fgPrior[kSPECIESCN] = {0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0};
 
 
@@ -269,15 +271,21 @@ AliPID& AliPID::operator = (const AliPID& pid)
 void AliPID::Init() 
 {
   //
-  // Initialise the masses
+  // Initialise the masses, charges
   //
   // Initialise only once... 
   if(!fgkParticleMass[0]) {
     AliPDG::AddParticlesToPdgDataBase();
     for (Int_t i = 0; i < kSPECIESC; i++) {
       fgkParticleMass[i] = M(i);
-      if (i == kHe3 || i == kAlpha) fgkParticleMassZ[i] = M(i)/2.;
-      else fgkParticleMassZ[i]=M(i);
+      if (i == kHe3 || i == kAlpha) {
+       fgkParticleMassZ[i] = M(i)/2.;
+       fgkParticleCharge[i] = 2;
+      }
+      else {
+       fgkParticleMassZ[i]=M(i);
+       fgkParticleCharge[i]=1;
+      }
     }
   }
 }
index c4b6875..760a736 100644 (file)
@@ -44,6 +44,10 @@ class AliPID : public TObject {
     kUnknown = 14
   };
   
+  static Int_t         ParticleCharge(Int_t iType) {
+     if(!fgkParticleMass[0]) Init(); 
+     return fgkParticleCharge[iType];
+  }
   static Float_t       ParticleMass(Int_t iType) {
      if(!fgkParticleMass[0]) Init(); 
      return fgkParticleMass[iType];
@@ -52,11 +56,6 @@ class AliPID : public TObject {
      if(!fgkParticleMass[0]) Init(); 
      return fgkParticleMassZ[iType];
   }
-  static Int_t         ParticleCharge(Int_t iType){
-    if(!fgkParticleMass[0]) Init();
-    if (iType<0||iType>=kSPECIESC) return 0;
-    return TMath::Nint(fgkParticleMass[iType]/fgkParticleMassZ[iType]);
-  }
   static const char*   ParticleName(Int_t iType) 
     {return fgkParticleName[iType];};
   static const char*   ParticleShortName(Int_t iType) 
@@ -100,12 +99,13 @@ class AliPID : public TObject {
 
   static /*const*/ Float_t fgkParticleMass[kSPECIESCN+1];  // particle masses
   static /*const*/ Float_t fgkParticleMassZ[kSPECIESCN+1]; // particle masses/charge
+  static /*const*/ Char_t  fgkParticleCharge[kSPECIESCN+1]; // particle charge (in e units!)
   static const char*   fgkParticleName[kSPECIESCN+1];      // particle names
   static const char*   fgkParticleShortName[kSPECIESCN+1]; // particle names
   static const char*   fgkParticleLatexName[kSPECIESCN+1]; // particle names
   static const Int_t   fgkParticleCode[kSPECIESCN+1];      // particle codes
 
-  ClassDef(AliPID, 3)                                      // particle id probability densities
+  ClassDef(AliPID, 4)                                      // particle id probability densities
 };
 
 
index 64bd058..eee16aa 100644 (file)
@@ -66,7 +66,7 @@ AliTOFtrack::AliTOFtrack(const AliESDtrack& t) :
   //
   SetLabel(t.GetLabel());
   SetChi2(0.);
-  SetMass(t.GetMass());
+  SetMass(t.GetMassForTracking());
 
   Set(t.GetX(),t.GetAlpha(),t.GetParameter(),t.GetCovariance());
 
index 5cd0ec6..5398c2b 100644 (file)
@@ -142,7 +142,7 @@ AliTPCtrack::AliTPCtrack(const AliESDtrack& t, TTreeSRedirector *pcstream) :
   const Double_t kmaxC[4]={10,10,0.1,0.1};  // cuts on the rms /fP0,fP1,fP2,fP3
   SetNumberOfClusters(t.GetTPCclusters(fIndex));
   SetLabel(t.GetLabel());
-  SetMass(t.GetMass());
+  SetMass(t.GetMassForTracking());
   for (Int_t i=0; i<4;i++) fPoints[i]=0.;
   for (Int_t i=0; i<12;i++) fKinkPoint[i]=0.;
   for (Int_t i=0; i<3;i++) fKinkIndexes[i]=0;
index a13ddaa..aceb181 100644 (file)
@@ -132,7 +132,7 @@ AliTRDtrackV1::AliTRDtrackV1(const AliESDtrack &t) : AliKalmanTrack()
   SetLabel(t.GetLabel());
   SetChi2(0.0);
 
-  SetMass(t.GetMass(kTRUE));
+  SetMass(t.GetMassForTracking());
   AliKalmanTrack::SetNumberOfClusters(t.GetTRDncls()); 
   Int_t ti[]={-1, -1, -1, -1, -1, -1}; t.GetTRDtracklets(&ti[0]);
   for(int ip=0; ip<kNplane; ip++){ 
@@ -669,6 +669,7 @@ Bool_t AliTRDtrackV1::PropagateTo(Double_t xk, Double_t xx0, Double_t xrho)
 
     // Energy losses
     Double_t p2    = (1.0 + GetTgl()*GetTgl()) / (GetSigned1Pt()*GetSigned1Pt());
+    if (fMass<0) p2 *= 4; // q=2
     Double_t beta2 = p2 / (p2 + fMass*fMass);
     if ((beta2 < 1.0e-10) || 
         ((5940.0 * beta2/(1.0 - beta2 + 1.0e-10) - beta2) < 0.0)) {
@@ -678,6 +679,7 @@ Bool_t AliTRDtrackV1::PropagateTo(Double_t xk, Double_t xx0, Double_t xrho)
     Double_t dE    = 0.153e-3 / beta2 
                    * (TMath::Log(5940.0 * beta2/(1.0 - beta2 + 1.0e-10)) - beta2)
                    * xrho;
+    if (fMass<0) dE *= 4; // q=2
     fBudget[0] += xrho;
 
     /*
index 0cbc79c..4ab9e5d 100644 (file)
@@ -1867,7 +1867,7 @@ Int_t AliTRDtrackerV1::PropagateToX(AliTRDtrackV1 &t, Double_t xToGo, Double_t m
 
 /*    // Correct for mean material budget
     Double_t dEdx(0.),
-             bg(t.GetP()/mass);
+             bg(TMath::Abs(t.GetP()/mass));
     if(AliLog::GetDebugLevel("TRD", "AliTRDtrackerV1")>=3){
       const char *pn[] = {"rho", "x/X0", "<A>", "<Z>", "L", "<Z/A>", "Nb"};
       printf("D-AliTRDtrackerV1::PropagateTo(): x[%6.2f] bg[%6.2f]\n", xpos, bg);