Excluding from the combined PID the cases of apparent mismatching and particle specie...
authorhristov <hristov@f7af4fe6-9843-0410-8265-dc069ae4e863>
Tue, 23 Oct 2007 09:35:29 +0000 (09:35 +0000)
committerhristov <hristov@f7af4fe6-9843-0410-8265-dc069ae4e863>
Tue, 23 Oct 2007 09:35:29 +0000 (09:35 +0000)
ITS/AliITSReconstructor.cxx
ITS/AliITSpidESD1.cxx
STEER/AliESDpid.cxx
TOF/AliTOFpidESD.cxx
TOF/AliTOFpidESD.h
TPC/AliTPCReconstructor.cxx
TPC/AliTPCpidESD.cxx

index 3d65b7c..a634b5d 100644 (file)
@@ -159,7 +159,7 @@ AliTracker* AliITSReconstructor::CreateTracker() const
   }
   else{
     Info("FillESD","ITS default PID\n");
-    Double_t parITS[] = {79.,0.13, 10.}; //IB: this is  "pp tuning"
+    Double_t parITS[] = {79.,0.13, 5.}; //IB: this is  "pp tuning"
     nc->fItsPID = new AliITSpidESD1(parITS);
   }
  
index b57d845..2e012f0 100755 (executable)
@@ -60,19 +60,31 @@ Int_t AliITSpidESD1::MakePID(AliESDEvent *event)
       if ((t->GetStatus()&AliESDtrack::kITSout)==0) continue;
     Double_t mom=t->GetP();
     Double_t dedx=t->GetITSsignal()/fMIP;
-    Int_t ns=AliPID::kSPECIES;
     Double_t p[10];
-    for (Int_t j=0; j<ns; j++) {
+    Bool_t mismatch=kTRUE, heavy=kTRUE;
+    for (Int_t j=0; j<AliPID::kSPECIES; j++) {
       Double_t mass=AliPID::ParticleMass(j);//GeV/c^2
       Double_t bethe=AliITSpidESD::Bethe(mom,mass);
       Double_t sigma=fRes*bethe;
       if (TMath::Abs(dedx-bethe) > fRange*sigma) {
        p[j]=TMath::Exp(-0.5*fRange*fRange)/sigma;
-        continue;
+      } else {
+        p[j]=TMath::Exp(-0.5*(dedx-bethe)*(dedx-bethe)/(sigma*sigma))/sigma;
+        mismatch=kFALSE;
       }
-      p[j]=TMath::Exp(-0.5*(dedx-bethe)*(dedx-bethe)/(sigma*sigma))/sigma;
+
+      // Check for particles heavier than (AliPID::kSPECIES - 1)
+      if (dedx < (bethe + fRange*sigma)) heavy=kFALSE;
+
     }
+
+    if (mismatch)
+       for (Int_t j=0; j<AliPID::kSPECIES; j++) p[j]=1/AliPID::kSPECIES;
+
     t->SetITSpid(p);
+
+    if (heavy) t->ResetStatus(AliESDtrack::kITSpid);
+
   }
   return 0;
 }
index 8a3daa9..616e775 100644 (file)
@@ -40,57 +40,41 @@ Int_t AliESDpid::MakePID(AliESDEvent *event)
   for (Int_t i=0; i<ntrk; i++) {
     Int_t ns=AliPID::kSPECIES;
     Double_t p[10]={1.,1.,1.,1.,1.,1.,1.,1.,1.,1.};
-    const Double_t keps=1e-13;
 
     AliESDtrack *t=event->GetTrack(i);
 
-    if ((t->GetStatus()&AliESDtrack::kITSpid )!=0) {
+    if (t->IsOn(AliESDtrack::kITSpid)) {
       Double_t d[10];
       t->GetITSpid(d);
-      Int_t j, ok=0;
-      for (j=0; j<ns; j++) if (d[j]>keps) ok=1;
-      if (ok) 
-      for (j=0; j<ns; j++) p[j]*=d[j];
+      for (Int_t j=0; j<ns; j++) p[j]*=d[j];
     }
 
-    if ((t->GetStatus()&AliESDtrack::kTPCpid )!=0) {
+    if (t->IsOn(AliESDtrack::kTPCpid)) {
       Double_t d[10];
       t->GetTPCpid(d);
-      Int_t j, ok=0;
-      for (j=0; j<ns; j++) if (d[j]>keps) ok=1;
-      if (ok) 
-      for (j=0; j<ns; j++) p[j]*=d[j];
+      for (Int_t j=0; j<ns; j++) p[j]*=d[j];
     }
 
-    if ((t->GetStatus()&AliESDtrack::kTRDpid )!=0) {
+    if (t->IsOn(AliESDtrack::kTRDpid)) {
       Double_t d[10];
       t->GetTRDpid(d);
-      Int_t j, ok=0;
-      for (j=0; j<ns; j++) if (d[j]>keps) ok=1;
-      if (ok) 
-      for (j=0; j<ns; j++) p[j]*=d[j];
+      for (Int_t j=0; j<ns; j++) p[j]*=d[j];
     }
 
-    if (t->GetP()>0.7) // accept the TOF only for the high momenta
-    if ((t->GetStatus()&AliESDtrack::kTOFpid )!=0) {
+    if (t->IsOn(AliESDtrack::kTOFpid)) {
       Double_t d[10];
       t->GetTOFpid(d);
-      Int_t j, ok=0;
-      for (j=0; j<ns; j++) if (d[j]>keps) ok=1;
-      if (ok) 
-      for (j=0; j<ns; j++) p[j]*=d[j];
+      for (Int_t j=0; j<ns; j++) p[j]*=d[j];
     }
 
-    if ((t->GetStatus()&AliESDtrack::kHMPIDpid )!=0) {
+    if (t->IsOn(AliESDtrack::kHMPIDpid)) {
       Double_t d[10];
       t->GetHMPIDpid(d);
-      Int_t j, ok=0;
-      for (j=0; j<ns; j++) if (d[j]>keps) ok=1;
-      if (ok) 
-      for (j=0; j<ns; j++) p[j]*=d[j];
+      for (Int_t j=0; j<ns; j++) p[j]*=d[j];
     }
 
     t->SetESDpid(p);
   }
+
   return 0;
 }
index 61acc8f..5b7acb6 100644 (file)
 ClassImp(AliTOFpidESD)
 
 //_________________________________________________________________________
-  AliTOFpidESD::AliTOFpidESD(): 
-  fN(-1),
-  fEventN(-1),
+AliTOFpidESD::AliTOFpidESD(): 
   fSigma(0),
-  fRange(0)
+  fRange(0),
+  fPmax(0)         // zero at 0.5 GeV/c for pp
 {
 }
 //_________________________________________________________________________
 AliTOFpidESD::AliTOFpidESD(Double_t *param):
-  fN(0),
-  fEventN(0),
-  fSigma(0),
-  fRange(0)
- {
+  fSigma(param[0]),
+  fRange(param[1]),
+  fPmax(0)          // zero at 0.5 GeV/c for pp
+{
   //
   //  The main constructor
   //
-  fSigma=param[0];
-  fRange=param[1];
+  //
+
+  //fPmax=TMath::Exp(-0.5*3*3)/fSigma; // ~3 sigma at 0.5 GeV/c for PbPb 
+}
+
+Double_t 
+AliTOFpidESD::GetMismatchProbability(Double_t p, Double_t mass) const {
+  //
+  // Returns the probability of mismatching 
+  // assuming 1/(p*beta)^2 scaling
+  //
+  const Double_t m=0.5;                   // "reference" momentum (GeV/c)
+
+  Double_t ref2=m*m*m*m/(m*m + mass*mass);// "reference" (p*beta)^2
+  Double_t p2beta2=p*p*p*p/(p*p + mass*mass);
 
+  return fPmax*ref2/p2beta2;
 }
 
 //_________________________________________________________________________
@@ -62,6 +74,7 @@ Int_t AliTOFpidESD::MakePID(AliESDEvent *event, Double_t timeZero)
 
   AliDebug(1,Form("TOF PID Parameters: Sigma (ps)= %f, Range= %f",fSigma,fRange));
   AliDebug(1,"++++++++++++++++++++++++++++++++++++++++++++++++++++++++++ \n");
+
   Int_t ntrk=event->GetNumberOfTracks();
   AliESDtrack **tracks=new AliESDtrack*[ntrk];
 
@@ -79,6 +92,7 @@ Int_t AliTOFpidESD::MakePID(AliESDEvent *event, Double_t timeZero)
     Double_t time[10]; t->GetIntegratedTimes(time);
     Double_t p[10];
     Double_t mom=t->GetP();
+    Bool_t mismatch=kTRUE, heavy=kTRUE;
     for (Int_t j=0; j<AliPID::kSPECIES; j++) {
       Double_t mass=AliPID::ParticleMass(j);
       Double_t dpp=0.01;      //mean relative pt resolution;
@@ -87,11 +101,25 @@ Int_t AliTOFpidESD::MakePID(AliESDEvent *event, Double_t timeZero)
       sigma=TMath::Sqrt(sigma*sigma + fSigma*fSigma);
       if (TMath::Abs(tof-time[j]) > fRange*sigma) {
        p[j]=TMath::Exp(-0.5*fRange*fRange)/sigma;
-        continue;
-      }
-      p[j]=TMath::Exp(-0.5*(tof-time[j])*(tof-time[j])/(sigma*sigma))/sigma;
+      } else 
+        p[j]=TMath::Exp(-0.5*(tof-time[j])*(tof-time[j])/(sigma*sigma))/sigma;
+
+      // Check the mismatching
+      Double_t pm=GetMismatchProbability(mom,mass);
+      if (p[j]>pm) mismatch=kFALSE;
+
+      // Check for particles heavier than (AliPID::kSPECIES - 1)
+      if (tof < (time[j] + fRange*sigma)) heavy=kFALSE;
+
     }
+
+    if (mismatch)
+       for (Int_t j=0; j<AliPID::kSPECIES; j++) p[j]=1/AliPID::kSPECIES;
+
     t->SetTOFpid(p);
+
+    if (heavy) t->ResetStatus(AliESDtrack::kTOFpid);    
+
   }
 
   delete[] tracks;
@@ -123,6 +151,7 @@ Int_t AliTOFpidESD::MakePID(AliESDEvent *event)
     Double_t time[10]; t->GetIntegratedTimes(time);
     Double_t p[10];
     Double_t mom=t->GetP();
+    Bool_t mismatch=kTRUE, heavy=kTRUE;
     for (Int_t j=0; j<AliPID::kSPECIES; j++) {
       Double_t mass=AliPID::ParticleMass(j);
       Double_t dpp=0.01;      //mean relative pt resolution;
@@ -131,11 +160,25 @@ Int_t AliTOFpidESD::MakePID(AliESDEvent *event)
       sigma=TMath::Sqrt(sigma*sigma + fSigma*fSigma);
       if (TMath::Abs(tof-time[j]) > fRange*sigma) {
        p[j]=TMath::Exp(-0.5*fRange*fRange)/sigma;
-        continue;
-      }
-      p[j]=TMath::Exp(-0.5*(tof-time[j])*(tof-time[j])/(sigma*sigma))/sigma;
+      } else
+        p[j]=TMath::Exp(-0.5*(tof-time[j])*(tof-time[j])/(sigma*sigma))/sigma;
+
+      // Check the mismatching
+      Double_t pm=GetMismatchProbability(mom,mass);
+      if (p[j]>pm) mismatch=kFALSE;
+
+      // Check for particles heavier than (AliPID::kSPECIES - 1)
+      if (tof < (time[j] + fRange*sigma)) heavy=kFALSE;
+
     }
+
+    if (mismatch)
+       for (Int_t j=0; j<AliPID::kSPECIES; j++) p[j]=1/AliPID::kSPECIES;
+
     t->SetTOFpid(p);
+
+    if (heavy) t->ResetStatus(AliESDtrack::kTOFpid);    
+
   }
 
   delete[] tracks;
index 281b365..58365a8 100644 (file)
@@ -16,25 +16,25 @@ class AliESDEvent;
 class AliTOFGeometry;
 
 class AliTOFpidESD : public TObject {
-enum {kMaxCluster=77777}; //maximal number of the TOF clusters
 public:
- AliTOFpidESD();
- AliTOFpidESD(Double_t *param);
+  AliTOFpidESD();
+  AliTOFpidESD(Double_t *param);
  ~AliTOFpidESD(){}
+  void     SetMaxMismatchProbability(Double_t p) {fPmax=p;}
+  Double_t GetMaxMismatchProbability() {return fPmax;}
 
   Int_t MakePID(AliESDEvent *event);
   Int_t MakePID(AliESDEvent *event, Double_t timeZero);
-  void  SetEventNumber(Int_t n) {fEventN=n;}
-  Int_t GetEventNumber() const {return fEventN;}
 
 private:
-  Int_t fN;               // number of the TOF clusters
-  Int_t fEventN;          // event number
+  Double_t GetMismatchProbability(Double_t p,Double_t mass) const;
+
   Double_t fSigma;        // intrinsic TOF resolution
   Double_t fRange;        // one particle type PID range (in sigmas)
+  Double_t fPmax;         // "maximal" probability of mismathing (at ~0.5 GeV/c)
 
-  ClassDef(AliTOFpidESD,1)   // TOF PID class
+  ClassDef(AliTOFpidESD,2)   // TOF PID class
 };
 
 #endif
index 780e09e..7add769 100644 (file)
@@ -108,7 +108,7 @@ void AliTPCReconstructor::FillESD(TTree */*digitsTree*/, TTree */*clustersTree*/
 {
 // make PID
 
-  Double_t parTPC[] = {47., 0.10, 10.};
+  Double_t parTPC[] = {47., 0.07, 5.};
   AliTPCpidESD tpcPID(parTPC);
   tpcPID.MakePID(esd);
 }
index 2909499..58f6638 100644 (file)
@@ -63,21 +63,35 @@ Int_t AliTPCpidESD::MakePID(AliESDEvent *event)
     AliESDtrack *t=event->GetTrack(i);
     if ((t->GetStatus()&AliESDtrack::kTPCin )==0)
       if ((t->GetStatus()&AliESDtrack::kTPCout)==0) continue;
-    Int_t ns=AliPID::kSPECIES;
     Double_t p[10];
-    for (Int_t j=0; j<ns; j++) {
+    Double_t mom=t->GetP();
+    const AliExternalTrackParam *in=t->GetInnerParam();
+    if (in) mom=in->GetP();
+    Double_t dedx=t->GetTPCsignal()/fMIP;
+    Bool_t mismatch=kTRUE, heavy=kTRUE;
+    for (Int_t j=0; j<AliPID::kSPECIES; j++) {
       Double_t mass=AliPID::ParticleMass(j);
-      Double_t mom=t->GetP();
-      Double_t dedx=t->GetTPCsignal()/fMIP;
       Double_t bethe=Bethe(mom/mass); 
       Double_t sigma=fRes*bethe;
       if (TMath::Abs(dedx-bethe) > fRange*sigma) {
        p[j]=TMath::Exp(-0.5*fRange*fRange)/sigma;
-        continue;
+      } else {
+        p[j]=TMath::Exp(-0.5*(dedx-bethe)*(dedx-bethe)/(sigma*sigma))/sigma;
+        mismatch=kFALSE;
       }
-      p[j]=TMath::Exp(-0.5*(dedx-bethe)*(dedx-bethe)/(sigma*sigma))/sigma;
+
+      // Check for particles heavier than (AliPID::kSPECIES - 1)
+      if (dedx < (bethe + fRange*sigma)) heavy=kFALSE;
+
     }
+
+    if (mismatch)
+       for (Int_t j=0; j<AliPID::kSPECIES; j++) p[j]=1/AliPID::kSPECIES;
+
     t->SetTPCpid(p);
+
+    if (heavy) t->ResetStatus(AliESDtrack::kTPCpid);
+
   }
   return 0;
 }