Check for PID added

[u/mrichter/AliRoot.git] / HBTAN / AliHBTReaderESD.cxx

diff --git a/HBTAN/AliHBTReaderESD.cxx b/HBTAN/AliHBTReaderESD.cxx
index d7a3915c1dc3db6ed4761b67db773a2b4067e755..717a5426093d257bf4c416d6e0cc1a10b6292231 100644 (file)
--- a/HBTAN/AliHBTReaderESD.cxx
+++ b/HBTAN/AliHBTReaderESD.cxx
@@ -41,13 +41,17 @@ AliHBTReaderESD::AliHBTReaderESD(const Char_t* esdfilename, const Char_t* galfil
  fRunLoader(0x0),
  fKeyIterator(0x0),
  fReadParticles(kFALSE),
+ fCheckParticlePID(kFALSE),
  fNTrackPoints(0),
  fdR(0.0),
  fClusterMap(kFALSE),
+ fMustTPC(kFALSE),
  fNTPCClustMin(0),
  fNTPCClustMax(150),
  fTPCChi2PerClustMin(0.0),
  fTPCChi2PerClustMax(10e5),
+ fChi2Min(0.0),
+ fChi2Max(10e5),
  fC00Min(0.0),
  fC00Max(10e5),
  fC11Min(0.0),
@@ -84,13 +88,17 @@ AliHBTReaderESD::AliHBTReaderESD(TObjArray* dirs,const Char_t* esdfilename, cons
  fRunLoader(0x0),
  fKeyIterator(0x0),
  fReadParticles(kFALSE),
+ fCheckParticlePID(kFALSE),
  fNTrackPoints(0),
  fdR(0.0),
  fClusterMap(kFALSE),
+ fMustTPC(kFALSE),
  fNTPCClustMin(0),
  fNTPCClustMax(150),
  fTPCChi2PerClustMin(0.0),
  fTPCChi2PerClustMax(10e5),
+ fChi2Min(0.0),
+ fChi2Max(10e5),
  fC00Min(0.0),
  fC00Max(10e5),
  fC11Min(0.0),
@@ -301,6 +309,15 @@ Int_t AliHBTReaderESD::ReadESD(AliESD* esd)
         continue;
       }
 
+     if (fMustTPC)
+      {
+       if ((esdtrack->GetStatus() & AliESDtrack::kTPCin) == kFALSE)
+        {
+          if (AliHBTParticle::GetDebug() > 2) 
+            Info("ReadNext","Particle skipped: Was not reconstructed in TPC.");
+          continue;
+        }
+      }     
      if ((esdtrack->GetStatus() & AliESDtrack::kESDpid) == kFALSE) 
       {
         if (AliHBTParticle::GetDebug() > 2) 
@@ -319,15 +336,13 @@ Int_t AliHBTReaderESD::ReadESD(AliESD* esd)
         continue;
       } 
      esdtrack->GetESDpid(pidtable);
-     //esdtrack->GetVertexPxPyPz(mom[0],mom[1],mom[2]); 
      esdtrack->GetConstrainedPxPyPz(mom);
-     //esdtrack->GetVertexXYZ(pos[0],pos[1],pos[2]);
      esdtrack->GetConstrainedXYZ(pos);
      pos[0] -= vertexpos[0];//we are interested only in relative position to Primary vertex at this point
      pos[1] -= vertexpos[1];
      pos[2] -= vertexpos[2];
 
-     Int_t charge = (extp[4] > 0)?-1:1;//if curvature=-charg/Pt is positive charge is negative
+     Int_t charge = (extp[4] > 0)?1:-1;//if curvature=charg/Pt is positive charge is positive
 
      //Particle from kinematics
      AliHBTParticle* particle = 0;
@@ -341,11 +356,22 @@ Int_t AliHBTReaderESD::ReadESD(AliESD* esd)
            Error("ReadNext","Can not find track with such label.");
            continue;
          }
+        if (fCheckParticlePID)
+         {
+           if(Pass(p->GetPdgCode())) 
+            {
+              if ( AliHBTParticle::GetDebug() > 5 )
+                Info("ReadNext","Simulated Particle PID (%d) did not pass the cut.",p->GetPdgCode());
+              continue; //check if we are intersted with particles of this type 
+            }
+         }
 //           if(p->GetPdgCode()<0) charge = -1;
         particle = new AliHBTParticle(*p,i);
 
       }
-
+      
+     if(CheckTrack(esdtrack)) continue;
+      
      //Here we apply Bayes' formula
      Double_t rc=0.;
      for (Int_t s=0; s<AliESDtrack::kSPECIES; s++) rc+=concentr[s]*pidtable[s];
@@ -451,6 +477,8 @@ Int_t AliHBTReaderESD::ReadESD(AliESD* esd)
          {
            Info("ReadNext","\n\nAdding Particle with incarnation %d",pdgcode);
            track->Print();
+           if (particle) particle->Print();
+           Info("ReadNext","\n----------------------------------------------\n");
          }
       }//for (Int_t s = 0; s<kNSpecies; s++)
 
@@ -563,6 +591,50 @@ Int_t AliHBTReaderESD::GetSpeciesPdgCode(ESpecies spec)//skowron
   return 0;
 }
 /********************************************************************/
+Bool_t AliHBTReaderESD::CheckTrack(AliESDtrack* t) const
+{
+  //Performs check of the track
+  
+  if ( (t->GetConstrainedChi2() < fChi2Min) || (t->GetConstrainedChi2() > fChi2Max) ) return kTRUE;
+  
+  if ( (t->GetTPCclusters(0x0) < fNTPCClustMin) || (t->GetTPCclusters(0x0) > fNTPCClustMax) ) return kTRUE;
+
+  if (t->GetTPCclusters(0x0) > 0)
+   {
+     Float_t chisqpercl = t->GetTPCchi2()/((Double_t)t->GetTPCclusters(0x0));
+     if ( (chisqpercl < fTPCChi2PerClustMin) || (chisqpercl > fTPCChi2PerClustMax) ) return kTRUE;
+   }
+
+  Double_t cc[15];
+  t->GetConstrainedExternalCovariance(cc);
+
+  if ( (cc[0]  < fC00Min) || (cc[0]  > fC00Max) ) return kTRUE;
+  if ( (cc[2]  < fC11Min) || (cc[2]  > fC11Max) ) return kTRUE;
+  if ( (cc[5]  < fC22Min) || (cc[5]  > fC22Max) ) return kTRUE;
+  if ( (cc[9]  < fC33Min) || (cc[9]  > fC33Max) ) return kTRUE;
+  if ( (cc[14] < fC44Min) || (cc[14] > fC44Max) ) return kTRUE;
+
+
+  t->GetInnerExternalCovariance(cc);
+
+  if ( (cc[0]  < fTPCC00Min) || (cc[0]  > fTPCC00Max) ) return kTRUE;
+  if ( (cc[2]  < fTPCC11Min) || (cc[2]  > fTPCC11Max) ) return kTRUE;
+  if ( (cc[5]  < fTPCC22Min) || (cc[5]  > fTPCC22Max) ) return kTRUE;
+  if ( (cc[9]  < fTPCC33Min) || (cc[9]  > fTPCC33Max) ) return kTRUE;
+  if ( (cc[14] < fTPCC44Min) || (cc[14] > fTPCC44Max) ) return kTRUE;
+
+  return kFALSE;
+
+}
+/********************************************************************/
+
+void AliHBTReaderESD::SetChi2Range(Float_t min, Float_t max)
+{
+  //sets range of Chi2 per Cluster
+  fChi2Min = min;
+  fChi2Max = max;
+}
+/********************************************************************/
 
 void AliHBTReaderESD::SetTPCNClustersRange(Int_t min,Int_t max)
 {