correct Filter Bit4 counting

[u/mrichter/AliRoot.git] / CORRFW / AliCFTrackCutPid.cxx

diff --git a/CORRFW/AliCFTrackCutPid.cxx b/CORRFW/AliCFTrackCutPid.cxx
index cee3edd083091b649e2671f5726ea46d49c4d6ab..7d94b0efc4c43a599b96544322edb6f052abffd9 100644 (file)
--- a/CORRFW/AliCFTrackCutPid.cxx
+++ b/CORRFW/AliCFTrackCutPid.cxx
@@ -84,16 +84,18 @@ AliCFTrackCutPid::AliCFTrackCutPid() :
   fNbins(101),
   fDetRestr(-1),
   fiPartRestr(-1),
-  fDetProbRestr(1)
+  fDetProbRestr(1),
+  fProbThreshold(0.)
 { 
   //
   //Default constructor 
   //
   for(Int_t j=0; j< AliPID::kSPECIES; j++) {
     fPriors[j]=0.2;
+  }
+  for(Int_t j=0; j< AliPID::kSPECIES; j++) {
     fPriorsFunc[j]=0x0;
   }
-  
   for(Int_t jDet=0; jDet< kNdets; jDet++)  {
     fDets[jDet]=kFALSE;
     fDetsInAnd[jDet]=kFALSE;
@@ -119,16 +121,18 @@ AliCFTrackCutPid::AliCFTrackCutPid(const Char_t* name, const Char_t* title) :
   fNbins(101),
   fDetRestr(-1),
   fiPartRestr(-1),
-  fDetProbRestr(1)
+  fDetProbRestr(1),
+  fProbThreshold(0.)
 {
   //
   //Constructor
   // 
   for(Int_t j=0; j< AliPID::kSPECIES; j++) {
     fPriors[j]=0.2;
+  }
+  for(Int_t j=0; j< AliPID::kSPECIES; j++) {
     fPriorsFunc[j]=0x0;
   }
-  
   for(Int_t jDet=0; jDet< kNdets; jDet++)  {
     fDets[jDet]=kFALSE;
     fDetsInAnd[jDet]=kFALSE;
@@ -154,7 +158,8 @@ AliCFTrackCutPid::AliCFTrackCutPid(const AliCFTrackCutPid& c) :
   fNbins(c.fNbins),
   fDetRestr(c.fDetRestr),
   fiPartRestr(c.fiPartRestr),
-  fDetProbRestr(c.fDetProbRestr)
+  fDetProbRestr(c.fDetProbRestr),
+  fProbThreshold(c.fProbThreshold)
 {
   //
   //Copy constructor
@@ -169,6 +174,8 @@ AliCFTrackCutPid::AliCFTrackCutPid(const AliCFTrackCutPid& c) :
   }
   for(Int_t j=0; j< AliPID::kSPECIES; j++){
     fPriors[j]=c.fPriors[j];
+  }
+  for(Int_t j=0; j< AliPID::kSPECIES; j++){
     fPriorsFunc[j]=c.fPriorsFunc[j];
     fhCombResp[j]=c.fhCombResp[j];
     fhCombProb[j]=c.fhCombProb[j];
@@ -198,6 +205,7 @@ AliCFTrackCutPid& AliCFTrackCutPid::operator=(const AliCFTrackCutPid& c)
     this->fDetRestr=c.fDetRestr;
     this->fiPartRestr=c.fiPartRestr;
     this->fDetProbRestr=c.fDetProbRestr;
+    this->fProbThreshold=c.fProbThreshold;
   
     for(Int_t i=0; i< kNdets ; i++ )   {
       this->fDets[i]=c.fDets[i];
@@ -209,6 +217,8 @@ AliCFTrackCutPid& AliCFTrackCutPid::operator=(const AliCFTrackCutPid& c)
 
     for(Int_t j=0; j< AliPID::kSPECIES; j++){
       this->fPriors[j]=c.fPriors[j];
+    }
+    for(Int_t j=0; j< AliPID::kSPECIES; j++){
       this->fhCombResp[j]=c.fhCombResp[j];
       this->fhCombProb[j]=c.fhCombProb[j];
       this-> fPriorsFunc[j]=c.fPriorsFunc[j];
@@ -346,7 +356,7 @@ void AliCFTrackCutPid::SetPPriors(AliESDtrack *pTrk)
   
   for(Int_t i=0; i< AliPID::kSPECIES; i++) {
     if(pTrk->P()>fPriorsFunc[i]->GetXmin() && pTrk->P() < fPriorsFunc[i]->GetXmax()) fPriors[i]=fPriorsFunc[i]->Eval(pTrk->P());
-    else {AliInfo("the track momentum is not in the function range. Priors are equal") fPriors[i] = 0.2;}   
+    else {AliInfo("the track momentum is not in the function range. Priors are equal"); fPriors[i] = 0.2;}   
   }
 }   
 //______________________________________
@@ -480,9 +490,15 @@ Int_t AliCFTrackCutPid::Identify(Double_t pid[AliPID::kSPECIES]) const
     if(fIsQAOn) fhCombProb[iP]->Fill(probability[iP]);
   }
   
-  AliPID::EParticleType sel = getpid.GetMostProbable();
-  if(getpid.GetProbability(sel,fPriors)>fCut) iPart= (Int_t)sel;
-  AliDebug(2,Form("probabilities   : %f  %f  %f  %f  %f",probability[0],probability[1],probability[2],probability[3],probability[4]));
+
+  if (fProbThreshold > 0.) {
+    if (probability[fgParticleType] >= fProbThreshold) iPart=fgParticleType;
+  }
+  else {
+    AliPID::EParticleType sel = getpid.GetMostProbable();
+    if(getpid.GetProbability(sel,fPriors)>fCut) iPart= (Int_t)sel;
+    AliDebug(2,Form("probabilities   : %f  %f  %f  %f  %f",probability[0],probability[1],probability[2],probability[3],probability[4]));
+  }
 
   if(fCheckResponse && !Check(pid,iPart, fMinDiffResponse)) iPart=kCheckResp;
   if(fCheckSelection && !Check(probability,iPart,fMinDiffProbability)) iPart=kCheckProb;
@@ -510,7 +526,8 @@ Int_t AliCFTrackCutPid::IdentifyQA(const Double_t pid[AliPID::kSPECIES], Int_t i
     fhProb[idets][iP]->Fill(probability[iP]);
   }
   
-  AliPID toresp(pid,kTRUE); Double_t qapriors[5]={0.2,0.2,0.2,0.2,0.2};
+  AliPID toresp(pid,kTRUE); 
+  Double_t qapriors[10]={0.2,0.2,0.2,0.2,0.2,0,0,0,0,0};
   toresp.SetPriors(qapriors);
   for(Int_t iPr=0; iPr<AliPID::kSPECIES; iPr++) fhResp[idets][iPr]->Fill(toresp.GetProbability((AliPID::EParticleType)iPr));
   
@@ -532,6 +549,7 @@ Bool_t AliCFTrackCutPid::IsSelected(TObject *track){
   TString className(track->ClassName());
   if (className.CompareTo("AliESDtrack") == 0) {
   AliESDtrack *esdTrack = dynamic_cast<AliESDtrack*>(track); 
+  if (!esdTrack) return kFALSE;
   ULong_t status[kNdets+1]={0,0,0,0,0,0};
   Double_t pid[kNdets+1][AliPID::kSPECIES];
   TrackInfo(esdTrack,status,pid);
@@ -541,6 +559,7 @@ Bool_t AliCFTrackCutPid::IsSelected(TObject *track){
 
   if (className.CompareTo("AliAODTrack") == 0) {
   AliAODTrack *aodtrack = dynamic_cast<AliAODTrack *>(track);
+  if (!aodtrack) return kFALSE ;
   if(GetAODID(aodtrack) == fgParticleType) sel = kTRUE;
   }
 
@@ -664,35 +683,35 @@ void AliCFTrackCutPid::DefineHistograms()
   //
 
   if(fgIsAOD){
-     char *partic[AliPID::kSPECIES]={"electron","muon","pion","kaon","proton"}; 
+     const char *partic[AliPID::kSPECIES]={"electron","muon","pion","kaon","proton"}; 
 
      for(Int_t iPart =0; iPart < AliPID::kSPECIES; iPart++)
       {
-       fhCombResp[iPart] = new TH1F(Form("rCombPart%i",iPart),Form(" %s combined response  (AODTrack)  ",partic[iPart]),fNbins,fXmin,fXmax);
+       fhCombResp[iPart] = new TH1F(Form("%s_rCombPart%i",GetName(),iPart),Form(" %s combined response  (AODTrack)  ",partic[iPart]),fNbins,fXmin,fXmax);
        fhCombResp[iPart]->SetXTitle(Form(" %s combined response ",partic[iPart]));
        fhCombResp[iPart]->SetYTitle("entries");
-       AliDebug(1,Form(  "rCombPart%i is booked!!",iPart));
-       fhCombProb[iPart] = new TH1F(Form("pCombPart%i",iPart),Form("%s combined probability (AODTrack) ",partic[iPart]),fNbins,fXmin,fXmax);
-       fhCombProb[iPart]->SetXTitle(Form(" %s combined response ",partic[iPart]));
+       AliDebug(1,Form(  "%s is booked!!",fhCombResp[iPart]->GetName()));
+       fhCombProb[iPart] = new TH1F(Form("%s_pCombPart%i",GetName(),iPart),Form("%s combined probability (AODTrack) ",partic[iPart]),fNbins,fXmin,fXmax);
+       fhCombProb[iPart]->SetXTitle(Form(" %s combined probability ",partic[iPart]));
        fhCombProb[iPart]->SetYTitle("entries");
-       AliDebug(1,Form(  "rCombProb%i is booked!!",iPart));
+       AliDebug(1,Form(  "%s is booked!!",fhCombProb[iPart]->GetName()));
      }
    }
 
 
   else {
-   char *detect[kNdets]={"ITS","TPC","TRD","TOF","HMPID"};
-   char *partic[AliPID::kSPECIES]={"electron","muon","pion","kaon","proton"};
+   const char *detect[kNdets]={"ITS","TPC","TRD","TOF","HMPID"};
+   const char *partic[AliPID::kSPECIES]={"electron","muon","pion","kaon","proton"};
   
     for(Int_t iDet =0; iDet< kNdets; iDet++)
      {
        if(!fDets[iDet]) continue;
        for(Int_t iP =0; iP < AliPID::kSPECIES; iP++){
-       fhResp[iDet][iP] = new TH1F(Form("rDet%iPart%i",iDet,iP),Form("%s response for %s  ",detect[iDet],partic[iP]),fNbins,fXmin,fXmax);
+       fhResp[iDet][iP] = new TH1F(Form("%s_rDet%iPart%i",GetName(),iDet,iP),Form("%s response for %s  ",detect[iDet],partic[iP]),fNbins,fXmin,fXmax);
         fhResp[iDet][iP]->SetXTitle(Form(" %s response ",partic[iP]));
         fhResp[iDet][iP]->SetYTitle("entries");
-       fhProb[iDet][iP] = new TH1F(Form("pDet%iPart%i",iDet,iP),Form("%s calculated probability for %s",detect[iDet],partic[iP]),fNbins,fXmin,fXmax);
-        fhProb[iDet][iP]->SetXTitle(Form(" %s response ",partic[iP]));
+       fhProb[iDet][iP] = new TH1F(Form("%s_pDet%iPart%i",GetName(),iDet,iP),Form("%s calculated probability for %s",detect[iDet],partic[iP]),fNbins,fXmin,fXmax);
+        fhProb[iDet][iP]->SetXTitle(Form(" %s probability ",partic[iP]));
         fhProb[iDet][iP]->SetYTitle("entries");
        }
      }
@@ -702,14 +721,14 @@ void AliCFTrackCutPid::DefineHistograms()
      { 
       for(Int_t iPart =0; iPart < AliPID::kSPECIES; iPart++)
         {
-         fhCombResp[iPart] = new TH1F(Form("rCombPart%i",iPart),Form(" %s combined response    ",partic[iPart]),fNbins,fXmin,fXmax);
+         fhCombResp[iPart] = new TH1F(Form("%s_rCombPart%i",GetName(),iPart),Form(" %s combined response    ",partic[iPart]),fNbins,fXmin,fXmax);
           fhCombResp[iPart]->SetXTitle(Form(" %s response ",partic[iPart]));
           fhCombResp[iPart]->SetYTitle("entries");
-          AliDebug(1,Form(  "rCombPart%i is booked!!",iPart));
-         fhCombProb[iPart] = new TH1F(Form("pCombPart%i",iPart),Form("%s combined probability ",partic[iPart]),fNbins,fXmin,fXmax);
+          AliDebug(1,Form(  "%s is booked!!",fhCombResp[iPart]->GetName()));
+         fhCombProb[iPart] = new TH1F(Form("%s_pCombPart%i",GetName(),iPart),Form("%s combined probability ",partic[iPart]),fNbins,fXmin,fXmax);
           fhCombProb[iPart]->SetXTitle(Form(" %s response ",partic[iPart]));
           fhCombProb[iPart]->SetYTitle("entries");
-          AliDebug(1,Form(  "rCombProb%i is booked!!",iPart));
+          AliDebug(1,Form(  "%s is booked!!",fhCombProb[iPart]->GetName()));
         }
      }
    }