Implementing two different test pid procedures based on the raw dE/dx info and a...
authorpchrist <pchrist@f7af4fe6-9843-0410-8265-dc069ae4e863>
Wed, 8 Apr 2009 11:37:34 +0000 (11:37 +0000)
committerpchrist <pchrist@f7af4fe6-9843-0410-8265-dc069ae4e863>
Wed, 8 Apr 2009 11:37:34 +0000 (11:37 +0000)
PWG2/SPECTRA/AliProtonAnalysisBase.cxx
PWG2/SPECTRA/AliProtonAnalysisBase.h
PWG2/SPECTRA/AliProtonQAAnalysis.cxx
PWG2/SPECTRA/macros/configProtonAnalysis.C
PWG2/SPECTRA/macros/runProtonAnalysis.C
PWG2/SPECTRA/macros/runProtonAnalysisQA.C

index 23f6ae0..7ff7d3f 100644 (file)
@@ -55,7 +55,7 @@ AliProtonAnalysisBase::AliProtonAnalysisBase() :
   fMaxDCAXY(0), fMaxDCAXYTPC(0),
   fMaxDCAZ(0), fMaxDCAZTPC(0),
   fMaxDCA3D(0), fMaxDCA3DTPC(0),
-  fMaxConstrainChi2(0),
+  fMaxConstrainChi2(0), fMinTPCdEdxPoints(0),
   fMinTPCClustersFlag(kFALSE), fMinITSClustersFlag(kFALSE),
   fMaxChi2PerTPCClusterFlag(kFALSE), fMaxChi2PerITSClusterFlag(kFALSE),
   fMaxCov11Flag(kFALSE), fMaxCov22Flag(kFALSE), 
@@ -70,12 +70,18 @@ AliProtonAnalysisBase::AliProtonAnalysisBase() :
   fPointOnITSLayer1Flag(0), fPointOnITSLayer2Flag(0),
   fPointOnITSLayer3Flag(0), fPointOnITSLayer4Flag(0),
   fPointOnITSLayer5Flag(0), fPointOnITSLayer6Flag(0),
+  fMinTPCdEdxPointsFlag(kFALSE),
   fFunctionProbabilityFlag(kFALSE), 
+  fNSigma(0),
   fElectronFunction(0), fMuonFunction(0),
   fPionFunction(0), fKaonFunction(0), fProtonFunction(0),
   fDebugMode(kFALSE) {
   //Default constructor
   for(Int_t i = 0; i < 5; i++) fPartFrac[i] = 0.0;
+  for(Int_t i = 0; i < 5; i++) {
+    fdEdxMean[i] = 0.0;
+    fdEdxSigma[i] = 0.0;
+  }
 }
 
 //____________________________________________________________________//
@@ -392,6 +398,13 @@ Bool_t AliProtonAnalysisBase::IsAccepted(AliESDEvent *esd,
       return kFALSE;
       }
   }
+  if(fMinTPCdEdxPointsFlag) {
+    if(track->GetTPCsignalN() < fMinTPCdEdxPoints) {
+      if(fDebugMode)
+       Printf("IsAccepted: Track rejected because it has %d TPC points for the calculation of the energy loss (min. requested: %d)",track->GetTPCsignalN(),fMinTPCdEdxPoints);
+      return kFALSE;
+    }
+  }
   if(fITSRefitFlag) {
     if ((track->GetStatus() & AliESDtrack::kITSrefit) == 0) {
       if(fDebugMode)
@@ -606,7 +619,14 @@ TCanvas *AliProtonAnalysisBase::GetListOfCuts() {
   listOfCuts = "PID mode: "; 
   if(fProtonPIDMode == kBayesian) listOfCuts += "Bayesian PID";
   if(fProtonPIDMode == kRatio) listOfCuts += "Z = ln((dE/dx)_{exp.}/(dE/dx)_{theor.})"; 
-  if(fProtonPIDMode == kSigma) listOfCuts += "N__{#sigma} area"; 
+  if(fProtonPIDMode == kSigma1) {
+    listOfCuts += "N_{#sigma}(1) area: "; listOfCuts += fNSigma;
+    listOfCuts += " #sigma";
+  }
+  if(fProtonPIDMode == kSigma2) {
+    listOfCuts += "N_{#sigma}(2) area: "; listOfCuts += fNSigma;
+    listOfCuts += " #sigma";
+  }
   l.DrawLatex(0.1,0.58,listOfCuts.Data());
   listOfCuts = "Accepted vertex diamond: "; 
   l.DrawLatex(0.1,0.5,listOfCuts.Data());
@@ -754,7 +774,7 @@ TCanvas *AliProtonAnalysisBase::GetListOfCuts() {
 Bool_t AliProtonAnalysisBase::IsProton(AliESDtrack *track) {
   //Function that checks if a track is a proton
   Double_t probability[5];
-  Double_t gPt = 0.0, gP = 0.0;
+  Double_t gPt = 0.0, gP = 0.0, gEta = 0.0;
   Long64_t fParticleType = 0;
  
   //Bayesian approach for the PID
@@ -791,13 +811,80 @@ Bool_t AliProtonAnalysisBase::IsProton(AliESDtrack *track) {
     return kFALSE;
   }
   //Definition of an N-sigma area around the dE/dx vs P band
-  else if(fProtonPIDMode == kSigma) {
-    Printf("The kSigma mode is not implemented yet!!!");
-    return kFALSE;
+  else if(fProtonPIDMode == kSigma1) {
+    AliExternalTrackParam *tpcTrack = (AliExternalTrackParam *)track->GetTPCInnerParam();
+    if(tpcTrack) {
+      gPt = tpcTrack->Pt();
+      gP = tpcTrack->P();
+      gEta = tpcTrack->Eta();
+    }
+    //We start the P slices at >0.3GeV/c with a bining of 50MeV/c ==> Int_t(0.3001/0.05) = 6
+    Int_t nbinP = Int_t(gP/0.05) - 6;
+    Double_t tpcSignal = track->GetTPCsignal();
+    Double_t dEdxMean = fdEdxMean[nbinP];
+    Double_t dEdxSigma = fdEdxSigma[nbinP];
+    if((tpcSignal <= dEdxMean + fNSigma*dEdxSigma)&&
+       (tpcSignal <= dEdxMean + fNSigma*dEdxSigma))
+      return kTRUE;
+  }//kSigma1 PID method
+  //Another definition of an N-sigma area around the dE/dx vs P band
+  else if(fProtonPIDMode == kSigma2) {
+    AliExternalTrackParam *tpcTrack = (AliExternalTrackParam *)track->GetTPCInnerParam();
+    if(tpcTrack) {
+      gPt = tpcTrack->Pt();
+      gP = tpcTrack->P();
+      gEta = tpcTrack->Eta();
+    }
+    //We start the P slices at >0.3GeV/c with a bining of 50MeV/c ==> Int_t(0.3001/0.05) = 6
+    Int_t nbinP = Int_t(gP/0.05) - 6;
+    Double_t tpcSignal = track->GetTPCsignal();
+    Double_t dEdxTheory = Bethe(gP/9.38270000000000048e-01);
+    Double_t dEdxSigma = fdEdxSigma[nbinP];
+    Double_t nsigma = TMath::Abs(tpcSignal - dEdxTheory)/(tpcSignal*(dEdxSigma/TMath::Sqrt(track->GetTPCsignalN())));
+    if(nsigma <= fNSigma) 
+      return kTRUE;
   }
 
   return kFALSE;
 }
 
+//________________________________________________________________________
+void AliProtonAnalysisBase::SetdEdxBandInfo(const char *filename) {
+  // This function is used in case the kSigma1 or kSigma2 PID mode is selected
+  // It takes as an argument the name of the ascii file (for the time being) 
+  // that is generated as a prior process.
+  // This ascii file has three columns: The min. P value (bins of 50MeV/c) 
+  // the mean and the sigma of the dE/dx distributions for protons coming 
+  // from a gaussian fit.
+  ifstream in;
+  in.open(filename);
+
+  Double_t gPtMin = 0.0;
+  Int_t iCounter = 0;
+  while(in.good()) {
+    in >> gPtMin >> fdEdxMean[iCounter] >> fdEdxSigma[iCounter];
+    if(fDebugMode)
+      Printf("Momentum bin: %d - Min momentum: %lf - mean(dE/dx): %lf - sigma(dE/dx): %lf",iCounter+1,gPtMin,fdEdxMean[iCounter],fdEdxSigma[iCounter]);
+    iCounter += 1;
+  }
+}
+
+//________________________________________________________________________
+Double_t AliProtonAnalysisBase::Bethe(Double_t bg) {
+  // This is the Bethe-Bloch function normalised to 1 at the minimum
+  // We renormalize it based on the MC information
+  // WARNING: To be revised soon!!!
+  // This is just a temporary fix
+  Double_t normalization = 49.2;
+  Double_t bg2=bg*bg;
+  Double_t beta2 = bg2/(1.+ bg2);
+  Double_t bb = 8.62702e-2*(9.14550 - beta2 - TMath::Log(3.51000e-5 + 1./bg2))/beta2;
+  //
+  const Float_t kmeanCorrection =0.1;
+  Double_t meanCorrection =(1+(bb-1)*kmeanCorrection);
+  bb *= meanCorrection;
+
+  return normalization*bb; 
+}
 
 
index 8b1ee34..50ff059 100644 (file)
@@ -27,7 +27,7 @@ class AliProtonAnalysisBase : public TObject {
  public:
   enum TriggerMode { kMB1 = 0, kMB2, kSPDFASTOR };
   enum AnalysisMode { kInvalid = -1, kTPC = 0, kHybrid, kGlobal };
-  enum PIDMode { kBayesian = 0, kRatio, kSigma };
+  enum PIDMode { kBayesian = 0, kRatio, kSigma1, kSigma2 };
 
   AliProtonAnalysisBase();
   virtual ~AliProtonAnalysisBase();
@@ -202,6 +202,13 @@ class AliProtonAnalysisBase : public TObject {
   Bool_t  IsUsedMaxConstrainChi2() {return fMaxConstrainChi2Flag;}
   Double_t   GetMaxConstrainChi2() {return fMaxConstrainChi2;}
 
+  void    SetMinTPCdEdxPoints(Int_t mindEdxpoints) {
+    fMinTPCdEdxPoints = mindEdxpoints;
+    fMinTPCdEdxPointsFlag = kTRUE;
+  }
+  Bool_t  IsUsedMinTPCdEdxPoints() {return fMinTPCdEdxPointsFlag;}
+  Int_t   GetMinTPCdEdxPoints() {return fMinTPCdEdxPoints;}
+  
   void    SetITSRefit() {fITSRefitFlag = kTRUE;}
   Bool_t  IsUsedITSRefit() {return fITSRefitFlag;}
   void    SetTPCRefit() {fTPCRefitFlag = kTRUE;}
@@ -215,6 +222,9 @@ class AliProtonAnalysisBase : public TObject {
 
   //PID related functions
   Bool_t IsProton(AliESDtrack *track);
+  void SetNSigma(Int_t nsigma) {fNSigma = nsigma;}  
+  Int_t GetNSigma() {return fNSigma;}
+  void SetdEdxBandInfo(const char* filename);
   void SetPriorProbabilities(Double_t * const partFrac) {
     for(Int_t i = 0; i < AliPID::kSPECIESN; i++) fPartFrac[i] = partFrac[i];} 
   void SetPriorProbabilityFunctions(TF1 *const felectron, 
@@ -231,6 +241,7 @@ class AliProtonAnalysisBase : public TObject {
   }
   Bool_t IsPriorProbabilityFunctionUsed() {return fFunctionProbabilityFlag;}
   Double_t GetParticleFraction(Int_t i, Double_t p);
+  Double_t Bethe(Double_t bg);
 
   void SetDebugMode() {fDebugMode = kTRUE;}
   Bool_t GetDebugMode() {return fDebugMode;}
@@ -262,6 +273,7 @@ class AliProtonAnalysisBase : public TObject {
   Double_t fMaxDCAZ, fMaxDCAZTPC; //max DCA z
   Double_t fMaxDCA3D, fMaxDCA3DTPC; //max DCA 3D
   Double_t fMaxConstrainChi2; //max constrain chi2 - vertex
+  Int_t  fMinTPCdEdxPoints;//min number of TPC points used for the dE/dx
   Bool_t fMinTPCClustersFlag, fMinITSClustersFlag; //shows if this cut is used or not
   Bool_t fMaxChi2PerTPCClusterFlag, fMaxChi2PerITSClusterFlag; //shows if this cut is used or not
   Bool_t fMaxCov11Flag, fMaxCov22Flag, fMaxCov33Flag, fMaxCov44Flag, fMaxCov55Flag; //shows if this cut is used or not
@@ -276,9 +288,13 @@ class AliProtonAnalysisBase : public TObject {
   Bool_t fPointOnITSLayer1Flag, fPointOnITSLayer2Flag; //shows if this cut is used or not
   Bool_t fPointOnITSLayer3Flag, fPointOnITSLayer4Flag; //shows if this cut is used or not
   Bool_t fPointOnITSLayer5Flag, fPointOnITSLayer6Flag; //shows if this cut is used or not
-  
+  Bool_t fMinTPCdEdxPointsFlag; //shows if this cut is used or not
+
   //pid
   Bool_t fFunctionProbabilityFlag; //flag: kTRUE if functions used
+  Int_t fNSigma; //N-sigma cut in the dE/dx band
+  Double_t fdEdxMean[24]; //mean values of the dE/dx distributions for the proton band - P slices
+  Double_t fdEdxSigma[24]; //sigma values of the dE/dx distributions for the proton band - P slices
   Double_t fPartFrac[10]; //prior probabilities
   TF1  *fElectronFunction; //momentum dependence of the prior probs
   TF1  *fMuonFunction; //momentum dependence of the prior probs
index 8eea2f8..c2ce2fd 100644 (file)
@@ -342,6 +342,13 @@ void AliProtonQAAnalysis::FillQA(AliStack *stack,
        else if(track->HasPointOnITSLayer(5))
          ((TH1F *)(fQAPrimaryProtonsAcceptedList->At(25)))->Fill(0);
       }//point on SSD2
+      if(fProtonAnalysisBase->IsUsedMinTPCdEdxPoints()) {
+       if(track->GetTPCsignalN() < fProtonAnalysisBase->GetMinTPCdEdxPoints()) {
+         ((TH1F *)(fQAPrimaryProtonsRejectedList->At(26)))->Fill(track->GetTPCsignalN());
+       }
+       if(track->GetTPCsignalN() >= fProtonAnalysisBase->GetMinTPCdEdxPoints())
+         ((TH1F *)(fQAPrimaryProtonsAcceptedList->At(26)))->Fill(track->GetTPCsignalN());
+      }//number of TPC points for the dE/dx
     }//primary particle cut
 
     //Secondaries
@@ -533,6 +540,13 @@ void AliProtonQAAnalysis::FillQA(AliStack *stack,
        else if(track->HasPointOnITSLayer(5))
          ((TH1F *)(fQASecondaryProtonsAcceptedList->At(25)))->Fill(0);
       }//point on SSD2
+      if(fProtonAnalysisBase->IsUsedMinTPCdEdxPoints()) {
+       if(track->GetTPCsignalN() < fProtonAnalysisBase->GetMinTPCdEdxPoints()) {
+         ((TH1F *)(fQASecondaryProtonsRejectedList->At(26)))->Fill(track->GetTPCsignalN());
+       }
+       if(track->GetTPCsignalN() >= fProtonAnalysisBase->GetMinTPCdEdxPoints())
+         ((TH1F *)(fQASecondaryProtonsAcceptedList->At(26)))->Fill(track->GetTPCsignalN());
+      }//number of TPC points for the dE/dx
     }//secondary particle cut
   }//protons
 
@@ -726,6 +740,13 @@ void AliProtonQAAnalysis::FillQA(AliStack *stack,
        else if(track->HasPointOnITSLayer(5))
          ((TH1F *)(fQAPrimaryAntiProtonsAcceptedList->At(25)))->Fill(0);
       }//point on SSD2
+      if(fProtonAnalysisBase->IsUsedMinTPCdEdxPoints()) {
+       if(track->GetTPCsignalN() < fProtonAnalysisBase->GetMinTPCdEdxPoints()) {
+         ((TH1F *)(fQAPrimaryAntiProtonsRejectedList->At(26)))->Fill(track->GetTPCsignalN());
+       }
+       if(track->GetTPCsignalN() >= fProtonAnalysisBase->GetMinTPCdEdxPoints())
+         ((TH1F *)(fQAPrimaryAntiProtonsAcceptedList->At(26)))->Fill(track->GetTPCsignalN());
+      }//number of TPC points for the dE/dx
     }//primary particle cut
 
     //Secondaries
@@ -915,6 +936,13 @@ void AliProtonQAAnalysis::FillQA(AliStack *stack,
        else if(track->HasPointOnITSLayer(5))
          ((TH1F *)(fQASecondaryAntiProtonsAcceptedList->At(25)))->Fill(0);
       }//point on SSD2
+      if(fProtonAnalysisBase->IsUsedMinTPCdEdxPoints()) {
+       if(track->GetTPCsignalN() < fProtonAnalysisBase->GetMinTPCdEdxPoints()) {
+         ((TH1F *)(fQASecondaryAntiProtonsRejectedList->At(26)))->Fill(track->GetTPCsignalN());
+       }
+       if(track->GetTPCsignalN() >= fProtonAnalysisBase->GetMinTPCdEdxPoints())
+         ((TH1F *)(fQASecondaryAntiProtonsAcceptedList->At(26)))->Fill(track->GetTPCsignalN());
+      }//number of TPC points for the dE/dx
     }//secondary particle cut
   }//antiprotons
 }
@@ -1899,6 +1927,8 @@ void AliProtonQAAnalysis::InitQA() {
   TH1F *gPrimaryProtonsPointOnITSLayer6Pass = new TH1F("gPrimaryProtonsPointOnITSLayer6Pass",
                                             "",10,-1,1);
   fQAPrimaryProtonsAcceptedList->Add(gPrimaryProtonsPointOnITSLayer6Pass);
+  TH1F *gPrimaryProtonsNumberOfTPCdEdxPointsPass = new TH1F("gPrimaryProtonsNumberOfTPCdEdxPointsPass","",100,0,200);
+  fQAPrimaryProtonsAcceptedList->Add(gPrimaryProtonsNumberOfTPCdEdxPointsPass);
 
   //Rejected primary protons
   /*gDirectory->cd("../");
@@ -1999,6 +2029,8 @@ void AliProtonQAAnalysis::InitQA() {
   TH1F *gPrimaryProtonsPointOnITSLayer6Reject = new TH1F("gPrimaryProtonsPointOnITSLayer6Reject",
                                             "",10,-1,1);
   fQAPrimaryProtonsRejectedList->Add(gPrimaryProtonsPointOnITSLayer6Reject);
+  TH1F *gPrimaryProtonsNumberOfTPCdEdxPointsReject = new TH1F("gPrimaryProtonsNumberOfTPCdEdxPointsReject","",100,0,200);
+  fQAPrimaryProtonsRejectedList->Add(gPrimaryProtonsNumberOfTPCdEdxPointsReject);
 
   //________________________________________________________________//
   /*gDirectory->cd("../../");
@@ -2103,6 +2135,8 @@ void AliProtonQAAnalysis::InitQA() {
   TH1F *gSecondaryProtonsPointOnITSLayer6Pass = new TH1F("gSecondaryProtonsPointOnITSLayer6Pass",
                                                         "",10,-1,1);
   fQASecondaryProtonsAcceptedList->Add(gSecondaryProtonsPointOnITSLayer6Pass);
+  TH1F *gSecondaryProtonsNumberOfTPCdEdxPointsPass = new TH1F("gSecondaryProtonsNumberOfTPCdEdxPointsPass","",100,0,200);
+  fQASecondaryProtonsAcceptedList->Add(gSecondaryProtonsNumberOfTPCdEdxPointsPass);
 
   //Rejected secondary protons
   /*gDirectory->cd("../");
@@ -2203,7 +2237,8 @@ void AliProtonQAAnalysis::InitQA() {
   TH1F *gSecondaryProtonsPointOnITSLayer6Reject = new TH1F("gSecondaryProtonsPointOnITSLayer6Reject",
                                                           "",10,-1,1);
   fQASecondaryProtonsRejectedList->Add(gSecondaryProtonsPointOnITSLayer6Reject);
-  
+  TH1F *gSecondaryProtonsNumberOfTPCdEdxPointsReject = new TH1F("gSecondaryProtonsNumberOfTPCdEdxPointsReject","",100,0,200);
+  fQASecondaryProtonsRejectedList->Add(gSecondaryProtonsNumberOfTPCdEdxPointsReject);  
 
   /*gDirectory->cd("../../../");
 
@@ -2312,6 +2347,8 @@ void AliProtonQAAnalysis::InitQA() {
   TH1F *gPrimaryAntiProtonsPointOnITSLayer6Pass = new TH1F("gPrimaryAntiProtonsPointOnITSLayer6Pass",
                                                           "",10,-1,1);
   fQAPrimaryAntiProtonsAcceptedList->Add(gPrimaryAntiProtonsPointOnITSLayer6Pass);
+  TH1F *gPrimaryAntiProtonsNumberOfTPCdEdxPointsPass = new TH1F("gPrimaryAntiProtonsNumberOfTPCdEdxPointsPass","",100,0,200);
+  fQAPrimaryAntiProtonsAcceptedList->Add(gPrimaryAntiProtonsNumberOfTPCdEdxPointsPass);
   
   //Rejected primary antiprotons
   /*gDirectory->cd("../");
@@ -2412,6 +2449,8 @@ void AliProtonQAAnalysis::InitQA() {
   TH1F *gPrimaryAntiProtonsPointOnITSLayer6Reject = new TH1F("gPrimaryAntiProtonsPointOnITSLayer6Reject",
                                                             "",10,-1,1);
   fQAPrimaryAntiProtonsRejectedList->Add(gPrimaryAntiProtonsPointOnITSLayer6Reject);
+  TH1F *gPrimaryAntiProtonsNumberOfTPCdEdxPointsReject = new TH1F("gPrimaryAntiProtonsNumberOfTPCdEdxPointsReject","",100,0,200);
+  fQAPrimaryAntiProtonsRejectedList->Add(gPrimaryAntiProtonsNumberOfTPCdEdxPointsReject);
   
   //________________________________________________________________//
   /*gDirectory->cd("../../");
@@ -2516,7 +2555,9 @@ void AliProtonQAAnalysis::InitQA() {
   TH1F *gSecondaryAntiProtonsPointOnITSLayer6Pass = new TH1F("gSecondaryAntiProtonsPointOnITSLayer6Pass",
                                                             "",10,-1,1);
   fQASecondaryAntiProtonsAcceptedList->Add(gSecondaryAntiProtonsPointOnITSLayer6Pass);
-  
+  TH1F *gSecondaryAntiProtonsNumberOfTPCdEdxPointsPass = new TH1F("gSecondaryAntiProtonsNumberOfTPCdEdxPointsPass","",100,0,200);
+  fQASecondaryAntiProtonsAcceptedList->Add(gSecondaryAntiProtonsNumberOfTPCdEdxPointsPass);
+
   //Rejected secondary antiprotons
   /*gDirectory->cd("../");
   TDirectory *dirAntiProtonsSecondaryRejected = gDirectory->mkdir("Rejected");
@@ -2616,6 +2657,8 @@ void AliProtonQAAnalysis::InitQA() {
   TH1F *gSecondaryAntiProtonsPointOnITSLayer6Reject = new TH1F("gSecondaryAntiProtonsPointOnITSLayer6Reject",
                                                             "",10,-1,1);
   fQASecondaryAntiProtonsRejectedList->Add(gSecondaryAntiProtonsPointOnITSLayer6Reject);
+  TH1F *gSecondaryAntiProtonsNumberOfTPCdEdxPointsReject = new TH1F("gSecondaryAntiProtonsNumberOfTPCdEdxPointsReject","",100,0,200);
+  fQASecondaryAntiProtonsRejectedList->Add(gSecondaryAntiProtonsNumberOfTPCdEdxPointsReject);
 }
 
 //____________________________________________________________________//
index 6142da2..1e2193a 100644 (file)
@@ -121,7 +121,7 @@ AliProtonAnalysisBase *GetProtonAnalysisBaseObject(const char* analysisLevel = "
     case "Bayesian":
       baseAnalysis->SetPIDMode(AliProtonAnalysisBase::kBayesian);
       //Momentum dependent priors
-      /*TFile *f = TFile::Open("PriorProb/PriorProbabilities.root ");
+      /*TFile *f = TFile::Open("$ALICE_ROOT/PWG2/data/PriorProbabilities.root ");
        TF1 *fitElectrons = (TF1 *)f->Get("fitElectrons");
        TF1 *fitMuons = (TF1 *)f->Get("fitMuons");
        TF1 *fitPions = (TF1 *)f->Get("fitPions");
@@ -140,8 +140,15 @@ AliProtonAnalysisBase *GetProtonAnalysisBaseObject(const char* analysisLevel = "
     case "Ratio":
       baseAnalysis->SetPIDMode(AliProtonAnalysisBase::kRatio);
       break;
-    case "Sigma":
-      baseAnalysis->SetPIDMode(AliProtonAnalysisBase::kSigma);
+    case "Sigma1":
+      baseAnalysis->SetPIDMode(AliProtonAnalysisBase::kSigma1);
+      baseAnalysis->SetNSigma(3);
+      baseAnalysis->SetdEdxBandInfo("$ALICE_ROOT/PWG2/data/protonsdEdxInfo.dat");
+      break;
+    case "Sigma2":
+      baseAnalysis->SetPIDMode(AliProtonAnalysisBase::kSigma2);
+      baseAnalysis->SetNSigma(3);
+      baseAnalysis->SetdEdxBandInfo("$ALICE_ROOT/PWG2/data/protonsdEdxInfo.dat");
       break;
     default:
       break;
index dfa0598..ae0363f 100644 (file)
@@ -15,9 +15,10 @@ void runProtonAnalysis(const char* esdAnalysisType = "Hybrid",
   //       name and .  
   //Analysis mode can be: "MC", "ESD", "AOD"
   //ESD analysis type can be one of the three: "TPC", "Hybrid", "Global"
-  //PID mode can be one of the three: "Bayesian" (standard Bayesian approach) 
+  //PID mode can be one of the four: "Bayesian" (standard Bayesian approach) 
   //   "Ratio" (ratio of measured over expected/theoretical dE/dx a la STAR) 
-  //   "Sigma" (N-sigma area around the fitted dE/dx vs P band)
+  //   "Sigma1" (N-sigma area around the fitted dE/dx vs P band)
+  //   "Sigma2" (same as previous but taking into account the No of TPC points)
   TStopwatch timer;
   timer.Start();
   
index 0f8e902..16aac54 100644 (file)
@@ -15,9 +15,10 @@ void runProtonAnalysisQA(const char* analysisType = "Hybrid",
   //       name and .  
   //Analysis mode can be: "MC", "ESD", "AOD"
   //ESD analysis type can be one of the three: "TPC", "Hybrid", "Global"
-  //PID mode can be one of the three: "Bayesian" (standard Bayesian approach) 
+  //PID mode can be one of the four: "Bayesian" (standard Bayesian approach) 
   //   "Ratio" (ratio of measured over expected/theoretical dE/dx a la STAR) 
-  //   "Sigma" (N-sigma area around the fitted dE/dx vs P band)
+  //   "Sigma1" (N-sigma area around the fitted dE/dx vs P band)
+  //   "Sigma2" (same as previous but taking into account the No of TPC points)
   TStopwatch timer;
   timer.Start();