]> git.uio.no Git - u/mrichter/AliRoot.git/blobdiff - PWGHF/vertexingHF/AliRDHFCutsD0toKpi.cxx
Merge branch 'feature-movesplit'
[u/mrichter/AliRoot.git] / PWGHF / vertexingHF / AliRDHFCutsD0toKpi.cxx
index 340bdebe8cd5b6786dc47e1722e7a37fb9c9044c..7d98652a56655eff93278e66391da24f517477dc 100644 (file)
@@ -41,19 +41,22 @@ using std::endl;
 ClassImp(AliRDHFCutsD0toKpi)
 
 //--------------------------------------------------------------------------
-AliRDHFCutsD0toKpi::AliRDHFCutsD0toKpi(const char* name) : 
-AliRDHFCuts(name),
-fUseSpecialCuts(kFALSE),
-fLowPt(kTRUE),
-fDefaultPID(kFALSE),
-fUseKF(kFALSE),
-fPtLowPID(2.),
-fPtMaxSpecialCuts(9999.),
-fmaxPtrackForPID(9999),
-fCombPID(kFALSE),
-fWeightsPositive(new Double_t[AliPID::kSPECIES]),
-fWeightsNegative(new Double_t[AliPID::kSPECIES]),
-fBayesianStrategy(0)
+AliRDHFCutsD0toKpi::AliRDHFCutsD0toKpi(const char* name) :
+   AliRDHFCuts(name),
+   fUseSpecialCuts(kFALSE),
+   fLowPt(kTRUE),
+   fDefaultPID(kFALSE),
+   fUseKF(kFALSE),
+   fPtLowPID(2.),
+   fPtMaxSpecialCuts(9999.),
+   fmaxPtrackForPID(9999),
+   fCombPID(kFALSE),
+   fnSpecies(AliPID::kSPECIES),
+   fWeightsPositive(0),
+   fWeightsNegative(0),
+   fProbThreshold(0.5),
+   fBayesianStrategy(0),
+   fBayesianCondition(0)
 {
   //
   // Default Constructor
@@ -98,6 +101,14 @@ fBayesianStrategy(0)
   Float_t limits[2]={0,999999999.};
   SetPtBins(2,limits);
 
+  fWeightsNegative = new Double_t[AliPID::kSPECIES];
+  fWeightsPositive = new Double_t[AliPID::kSPECIES];
+  
+  for (Int_t i = 0; i<AliPID::kSPECIES; i++) {
+    fWeightsPositive[i] = 0.;
+    fWeightsNegative[i] = 0.;
+  }
 }
 //--------------------------------------------------------------------------
 AliRDHFCutsD0toKpi::AliRDHFCutsD0toKpi(const AliRDHFCutsD0toKpi &source) :
@@ -110,9 +121,12 @@ AliRDHFCutsD0toKpi::AliRDHFCutsD0toKpi(const AliRDHFCutsD0toKpi &source) :
   fPtMaxSpecialCuts(source.fPtMaxSpecialCuts),
   fmaxPtrackForPID(source.fmaxPtrackForPID),
   fCombPID(source.fCombPID),
+  fnSpecies(source.fnSpecies),
   fWeightsPositive(source.fWeightsPositive),
   fWeightsNegative(source.fWeightsNegative),
-  fBayesianStrategy(source.fBayesianStrategy)
+  fProbThreshold(source.fProbThreshold),
+  fBayesianStrategy(source.fBayesianStrategy),
+  fBayesianCondition(source.fBayesianCondition)
 {
   //
   // Copy constructor
@@ -136,9 +150,12 @@ AliRDHFCutsD0toKpi &AliRDHFCutsD0toKpi::operator=(const AliRDHFCutsD0toKpi &sour
   fPtMaxSpecialCuts=source.fPtMaxSpecialCuts;
   fmaxPtrackForPID=source.fmaxPtrackForPID;
   fCombPID = source.fCombPID;
+  fnSpecies = source.fnSpecies;
   fWeightsPositive = source.fWeightsPositive;
   fWeightsNegative = source.fWeightsNegative;
+  fProbThreshold = source.fProbThreshold;
   fBayesianStrategy = source.fBayesianStrategy;
+  fBayesianCondition = source.fBayesianCondition;
   return *this;
 }
 
@@ -405,31 +422,11 @@ Int_t AliRDHFCutsD0toKpi::IsSelected(TObject* obj,Int_t selectionLevel,AliAODEve
       fIsSelectedPID=returnvaluePID;
       if(!returnvaluePID) return 0;
    } else {
-      switch (fBayesianStrategy) {
-      case kBayesSimple: {
-       returnvaluePID = IsSelectedSimpleBayesianPID(d);
-       break;
-      }
-      case kBayesWeightNoFilter: {
-       CalculateBayesianWeights(d);
-       returnvaluePID = 3;
-       break;
-      }
-      case(kBayesWeight): {
-       returnvaluePID = IsSelectedCombPID(d);
-       break;
-      }
-      case(kBayesMomentum): {
-       returnvaluePID = IsSelectedCombPID(d);
-       break;
-      }
-      default: {
-       returnvaluePID = IsSelectedCombPID(d);
-       break;
-      }
+      returnvaluePID = IsSelectedCombPID(d);
+      if(!returnvaluePID) return 0;
       }
     }
-  }
+  
 
 
   Int_t returnvalueComb=CombineSelectionLevels(3,returnvalueCuts,returnvaluePID);
@@ -755,6 +752,7 @@ Int_t AliRDHFCutsD0toKpi::IsSelectedPID(AliAODRecoDecayHF* d)
   
   Bool_t checkPIDInfo[2]={kTRUE,kTRUE};
   Double_t sigma_tmp[3]={fPidHF->GetSigma(0),fPidHF->GetSigma(1),fPidHF->GetSigma(2)};
+  Bool_t isTOFused=fPidHF->GetTOF(),isCompat=fPidHF->GetCompat();
   for(Int_t daught=0;daught<2;daught++){
     //Loop con prongs
     AliAODTrack *aodtrack=(AliAODTrack*)d->GetDaughter(daught);
@@ -777,11 +775,11 @@ Int_t AliRDHFCutsD0toKpi::IsSelectedPID(AliAODRecoDecayHF* d)
       }else{
        fPidHF->SetTOF(kFALSE);
        combinedPID[daught][1]=fPidHF->MakeRawPid(aodtrack,2);
-       fPidHF->SetTOF(kTRUE);
-       fPidHF->SetCompat(kTRUE);
+       if(isTOFused)fPidHF->SetTOF(kTRUE);
+       if(isCompat)fPidHF->SetCompat(kTRUE);
       }
     }
-
+  
     if(combinedPID[daught][0]<=-1&&combinedPID[daught][1]<=-1){ // if not a K- and not a pi- both D0 and D0bar excluded
       isD0D0barPID[0]=0;
       isD0D0barPID[1]=0;
@@ -1234,7 +1232,8 @@ void AliRDHFCutsD0toKpi::SetStandardCutsPP2010() {
   esdTrackCuts->SetPtRange(0.3,1.e10);
   
   AddTrackCuts(esdTrackCuts);
-
+  delete esdTrackCuts;
+  esdTrackCuts=NULL;
   
   const Int_t nptbins =14;
   const Double_t ptmax = 9999.;
@@ -1329,6 +1328,142 @@ void AliRDHFCutsD0toKpi::SetStandardCutsPP2010() {
 
 }
 
+//___________________________________________________________________________
+void AliRDHFCutsD0toKpi::SetStandardCutsPP2010vsMult() {
+  //
+  // Cuts for 2010 pp 7 TeV data analysis in Ntracklets bins (vs multiplicity)
+  //
+  SetName("D0toKpiCuts");
+  SetTitle("Cuts for D0 analysis in 2010-data pp 7 TeV vs multiplicity");
+
+  //
+  // Track cuts
+  //
+  AliESDtrackCuts* esdTrackCuts=new AliESDtrackCuts();
+  esdTrackCuts->SetRequireSigmaToVertex(kFALSE);
+  //default
+  esdTrackCuts->SetRequireTPCRefit(kTRUE);
+  esdTrackCuts->SetRequireITSRefit(kTRUE);
+  esdTrackCuts->SetEtaRange(-0.8,0.8);
+  esdTrackCuts->SetClusterRequirementITS(AliESDtrackCuts::kSPD,
+                                        AliESDtrackCuts::kAny); 
+ // default is kBoth, otherwise kAny
+  esdTrackCuts->SetMinDCAToVertexXY(0.);
+  esdTrackCuts->SetPtRange(0.3,1.e10);
+
+  AddTrackCuts(esdTrackCuts);
+  delete esdTrackCuts;
+  esdTrackCuts=NULL;
+
+
+  //
+  // Cut values per pt bin
+  //
+  const Int_t nvars=11;
+  const Int_t nptbins=14;
+  Float_t* ptbins;
+  ptbins=new Float_t[nptbins+1];
+  ptbins[0]=0.;
+  ptbins[1]=0.5;
+  ptbins[2]=1.;
+  ptbins[3]=2.;
+  ptbins[4]=3.;
+  ptbins[5]=4.;
+  ptbins[6]=5.;
+  ptbins[7]=6.;
+  ptbins[8]=7.;
+  ptbins[9]=8.;
+  ptbins[10]=12.;
+  ptbins[11]=16.;
+  ptbins[12]=20.;
+  ptbins[13]=24.;
+  ptbins[14]=9999.;
+
+  SetPtBins(nptbins+1,ptbins);
+  
+  //setting cut values
+  Float_t cutsMatrixD0toKpiStand[nptbins][nvars]={{0.400,350.*1E-4,0.8,0.3,0.3,1000.*1E-4,1000.*1E-4,-5000.*1E-8,0.80,0.,3.},/* pt<0.5*/
+                                                  {0.400,350.*1E-4,0.8,0.3,0.3,1000.*1E-4,1000.*1E-4,-5000.*1E-8,0.80,0.,3.},/* 0.5<pt<1*/
+                                                  {0.400,300.*1E-4,0.8,0.4,0.4,1000.*1E-4,1000.*1E-4,-30000.*1E-8,0.80,0.,4.},/* 1<pt<2 */
+                                                  {0.400,300.*1E-4,0.8,0.7,0.7,1000.*1E-4,1000.*1E-4,-20000.*1E-8,0.85,0.,4.},/* 2<pt<3 */
+                                                  {0.400,300.*1E-4,0.8,0.7,0.7,1000.*1E-4,1000.*1E-4,-15000.*1E-8,0.85,0.,4.},/* 3<pt<4 */
+                                                  {0.400,300.*1E-4,0.75,0.7,0.7,1000.*1E-4,1000.*1E-4,-10000.*1E-8,0.875,0.,0.},/* 4<pt<5 */
+                                                  {0.400,300.*1E-4,0.75,0.7,0.7,1000.*1E-4,1000.*1E-4,-10000.*1E-8,0.875,0.,0.},/* 5<pt<6 */
+                                                  {0.400,400.*1E-4,0.8,0.7,0.7,1000.*1E-4,1000.*1E-4,-8000.*1E-8,0.85,0.,0.},/* 6<pt<7 */
+                                                  {0.400,400.*1E-4,0.8,0.7,0.7,1000.*1E-4,1000.*1E-4,-8000.*1E-8,0.85,0.,0.},/* 7<pt<8 */
+                                                 {0.400,0.06,0.85,0.7,0.7,1000.*1E-4,1000.*1E-4,-0.00001,0.85,0.,0.},/* 8<pt<12 */
+                                                 {0.400,0.09,1.0,0.7,0.7,9999.,9999.,0.,0.,0.,0.},/* 12<pt<16 */
+                                                 {0.400,0.09,1.0,0.7,0.7,9999.,9999.,0.,0.,0.,0.},/* 16<pt<20 */
+                                                 {0.400,0.09,1.0,0.7,0.7,9999.,9999.,0.,0.,0.,0.},/* 20<pt<24 */
+                                                 {0.400,0.09,1.0,0.7,0.7,9999.,9999.,0.,0.,0.,0.}};/* pt>24 */
+  
+  
+  //CREATE TRANSPOSE MATRIX...REVERSE INDICES as required by AliRDHFCuts
+  Float_t **cutsMatrixTransposeStand=new Float_t*[nvars];
+  for(Int_t iv=0;iv<nvars;iv++)cutsMatrixTransposeStand[iv]=new Float_t[nptbins];
+  
+  for (Int_t ibin=0;ibin<nptbins;ibin++){
+    for (Int_t ivar = 0; ivar<nvars; ivar++){
+      cutsMatrixTransposeStand[ivar][ibin]=cutsMatrixD0toKpiStand[ibin][ivar];      
+    }
+  }
+  
+  SetCuts(nvars,nptbins,cutsMatrixTransposeStand);
+  for(Int_t iv=0;iv<nvars;iv++) delete [] cutsMatrixTransposeStand[iv];
+  delete [] cutsMatrixTransposeStand;
+
+  SetUseSpecialCuts(kTRUE);
+  SetMaximumPtSpecialCuts(8.);
+
+  //Do recalculate the vertex
+  SetRemoveDaughtersFromPrim(kTRUE);
+
+  //
+  // Pid settings
+  //
+  Bool_t pidflag=kTRUE;
+  SetUsePID(pidflag);
+  if(pidflag) cout<<"PID is used"<<endl;
+  else cout<<"PID is not used"<<endl;
+  //
+  AliAODPidHF* pidObj=new AliAODPidHF();
+  Int_t mode=1;
+  const Int_t nlims=2;
+  Double_t plims[nlims]={0.6,0.8}; //TPC limits in momentum [GeV/c]
+  Bool_t compat=kTRUE; //effective only for this mode
+  Bool_t asym=kTRUE;
+  Double_t sigmas[5]={2.,1.,0.,3.,0.}; //to be checked and to be modified with new implementation of setters by Rossella
+  pidObj->SetAsym(asym);// if you want to use the asymmetric bands in TPC
+  pidObj->SetMatch(mode);
+  pidObj->SetPLimit(plims,nlims);
+  pidObj->SetSigma(sigmas);
+  pidObj->SetCompat(compat);
+  pidObj->SetTPC(kTRUE);
+  pidObj->SetTOF(kTRUE);
+  pidObj->SetPCompatTOF(1.5);
+  pidObj->SetSigmaForTPCCompat(3.);
+  pidObj->SetSigmaForTOFCompat(3.);
+  pidObj->SetOldPid(kTRUE);
+  //  pidObj->SetOldPid(kFALSE);
+  SetPidHF(pidObj);
+
+  SetUseDefaultPID(kFALSE); //to use the AliAODPidHF
+
+  SetLowPt(kFALSE);
+
+  //activate pileup rejection (for pp)
+  SetOptPileup(AliRDHFCuts::kRejectPileupEvent);
+
+
+  PrintAll();
+
+  delete pidObj;
+  pidObj=NULL;
+  delete [] ptbins;
+  ptbins=NULL;
+
+  return;
+}
 
 //---------------------------------------------------------------------------
 void AliRDHFCutsD0toKpi::SetStandardCutsPP2011_276TeV() {
@@ -1359,6 +1494,8 @@ void AliRDHFCutsD0toKpi::SetStandardCutsPP2011_276TeV() {
   esdTrackCuts->SetMinDCAToVertexXYPtDep("0.0075*TMath::Max(0.,(1-TMath::Floor(TMath::Abs(pt)/2.)))");
 
   AddTrackCuts(esdTrackCuts);
+  delete esdTrackCuts;
+  esdTrackCuts=NULL;
 
 
   const Int_t nvars=11;
@@ -1501,6 +1638,8 @@ void AliRDHFCutsD0toKpi::SetStandardCutsPbPb2010() {
 
 
   AddTrackCuts(esdTrackCuts);
+  delete esdTrackCuts;
+  esdTrackCuts=NULL;
   // SPD k FIRST for Pt<3 GeV/c
   SetSelectCandTrackSPDFirst(kTRUE, 3); 
 
@@ -1637,6 +1776,8 @@ void AliRDHFCutsD0toKpi::SetStandardCutsPbPb2010Peripherals() {
 
 
   AddTrackCuts(esdTrackCuts);
+  delete esdTrackCuts;
+  esdTrackCuts=NULL;
   // SPD k FIRST for Pt<3 GeV/c
   SetSelectCandTrackSPDFirst(kTRUE, 3); 
 
@@ -1765,281 +1906,189 @@ Int_t AliRDHFCutsD0toKpi::IsSelectedCombPID(AliAODRecoDecayHF* d)
 
   if (!fUsePID || !d) return 3;
 
-  switch (fBayesianStrategy) {
-
-  case kBayesWeightNoFilter: {
-    CalculateBayesianWeights(d);   //If weighted, no filtering: Calculate weights, return as unidentified
-    return 3;
-    break;
-  }
-  case kBayesSimple: {
-    return IsSelectedSimpleBayesianPID(d);   // If simple, go to simple method
-    break;
-  }
-  case(kBayesWeight): {
-    break;
-  }
-  case(kBayesMomentum): {
-    break;  //If weighted or momentum method, stay in this function
-  }
-
+  
+  if (fBayesianStrategy == kBayesWeightNoFilter) {
+     //WeightNoFilter: Accept all particles (no PID cut) but fill mass histos with weights in task
+     CalculateBayesianWeights(d);
+     return 3;
   }
 
-  //  Int_t isD0D0barPID[2]={1,2};
-
-
 
 
   Int_t isD0 = 0;
   Int_t isD0bar = 0;
   Int_t returnvalue = 0;
 
-  Int_t isPosKaon = 0, isNegKaon = 0;
-  //   Int_t combinedPID[2][2];// CONVENTION: [daught][isK,IsPi]; [0][0]=(prong 1, isK)=value [0][1]=(prong 1, isPi)=value;
-  //                                                                                                 same for prong 2
-  //                                               values convention  0 = not identified (but compatible) || No PID (->hasPID flag)
-  //                                                                  1 = identified
-  // PID search:   pion (TPC) or not K (TOF), Kaon hypothesis for both
-  // Initial hypothesis: unknown (but compatible)
-
-
-  //   combinedPID[0][0] = 0; // First daughter, kaon
-  //   combinedPID[0][1] = 0; // First daughter, pion
-  //   combinedPID[1][0] = 0; // Second daughter, kaon
-  //   combinedPID[1][1] = 0; // Second daughter, pion
+  Int_t isPosKaon = 0, isNegKaon = 0, isPosPion = 0, isNegPion = 0;
 
+  //Bayesian methods used here check for ID of kaon, and whether it is positive or negative.
+  
+  
   Bool_t checkPIDInfo[2] = {kTRUE, kTRUE};
   AliAODTrack *aodtrack[2] = {(AliAODTrack*)d->GetDaughter(0), (AliAODTrack*)d->GetDaughter(1)};
 
   if ((aodtrack[0]->Charge()*aodtrack[1]->Charge()) != -1) {
     return 0;  //reject if charges not opposite
   }
+   Double_t momentumpositive=0., momentumnegative=0.;  //Used in "prob > prior" method
   for (Int_t daught = 0; daught < 2; daught++) {
     //Loop over prongs
 
-    if (!(fPidHF->CheckStatus(aodtrack[daught], "TPC")) && !(fPidHF->CheckStatus(aodtrack[daught], "TOF"))) {
-      checkPIDInfo[daught] = kFALSE;
-      continue;
-    }
-    CalculateBayesianWeights(d);
-    //      Double_t prob0[AliPID::kSPECIES];
-    //      fPidHF->GetPidCombined()->ComputeProbabilities(aodtrack[daught], fPidHF->GetPidResponse(), prob0);
-    ///Three possible Bayesian methods: Pick one (only first implemented so far. TODO: implement other two methods)
-    ///A: Standard max. probability method   
-    if (TMath::MaxElement(AliPID::kSPECIES, fWeightsPositive) == fWeightsPositive[AliPID::kKaon]) { //If highest probability lies with kaon
-      isPosKaon = 1;  //flag [daught] as a kaon
-    }
-
-    if (TMath::MaxElement(AliPID::kSPECIES, fWeightsNegative) == fWeightsNegative[AliPID::kKaon]) { //If highest probability lies with kaon
-      isNegKaon = 1;  //flag [daught] as a kaon
-    }
-    ///B: Method based on probability greater than prior
-    /* Double_t kaonpriorprob;
-       TH1D *kaonprior = fPidHF->GetPidCombined()->GetPriorDistribution(AliPID::kKaon);
-       kaonpriorprob =   kaonprior->GetBinContent(kaonprior->FindBin(aodtrack[daught]->P()));
-
-       if (prob0[AliPID::kKaon] > kaonpriorprob){combinedPID[daught][0] = 1;}
-
-    */
-    ///C: Probability threshold (TODO: Add necessary threshold variables; currently unused)
-    /*
-      if(prob0[AliPID::kKaon] >= fMinProbKaon) {combinedPID[daught][0] = 1;}
-
-    */
-
-  }
-   
-   
-  //Loop over daughters ends here
-
-  if (!checkPIDInfo[0] && !checkPIDInfo[1]) {
-    return 0;      //Reject if both daughters lack both TPC and TOF info
-  }
-  //Momentum-based selection
-
-
-  if (isNegKaon && isPosKaon) { // If both are kaons, reject
-    isD0 = 0;
-    isD0bar = 0;
-  } else if (isNegKaon) {       //If negative kaon present, D0
-    isD0 = 1;
-  } else if (isPosKaon) {       //If positive kaon present, D0bar
-    isD0bar = 1;
-  } else {                      //If neither ID'd as kaon, subject to extra tests
-    isD0 = 1;
-    isD0bar = 1;
-    if (aodtrack[0]->P() < 1.5 && (fPidHF->CheckStatus(aodtrack[0], "TOF"))) { //If it's low momentum and not ID'd as a kaon, we assume it must be a pion
-      if (aodtrack[0]->Charge() == -1) {
-       isD0 = 0;  //Check charge and reject based on pion hypothesis
-      }
-      if (aodtrack[0]->Charge() == 1) {
-       isD0bar = 0;
-      }
-    }
-    if (aodtrack[1]->P() < 1.5 && (fPidHF->CheckStatus(aodtrack[1], "TOF"))) {
-      if (aodtrack[1]->Charge() == -1) {
-       isD0 = 0;
+      if (aodtrack[daught]->Charge() == -1) {
+         momentumnegative = aodtrack[daught]->P();
       }
-      if (aodtrack[1]->Charge() == 1) {
-       isD0bar = 0;
+      if (aodtrack[daught]->Charge() == +1) {
+         momentumpositive = aodtrack[daught]->P();
       }
-    }
-  }
-
-   
-
-  if (isD0 && isD0bar) {
-    returnvalue = 3;
-  }
-  if (isD0 && !isD0bar) {
-    returnvalue = 1;
-  }
-  if (!isD0 && isD0bar) {
-    returnvalue = 2;
-  }
-  if (!isD0 && !isD0bar) {
-    returnvalue = 0;
-  }
-
-  return returnvalue;
-
-
+         if (!(fPidHF->CheckStatus(aodtrack[daught], "TPC")) && !(fPidHF->CheckStatus(aodtrack[daught], "TOF"))) {
+            checkPIDInfo[daught] = kFALSE;
+            continue;
+         }
 
-  /*
+   }
 
-  //OLD CODE
-  if (combinedPID[0][0] == 1 && combinedPID[1][1] == 1) {          //First track is kaon, second is pion
-  if (aodtrack[0]->Charge() == -1) isD0 = 1;            //Kaon negative: particle is D0
-  else isD0bar = 1;                        //Kaon positive: particle is anti
-  } else if (combinedPID[1][0] == 1 && combinedPID[0][1] == 1) {         //First track is pion, second is kaon
-  if (aodtrack[1]->Charge() == -1) isD0 = 1;              //Kaon negative: particle is D0
-  else isD0bar = 1;                        //Kaon positive: particle is anti
-  }
+   //Loop over daughters ends here
 
-  else if (combinedPID[0][0] == 1 && combinedPID[1][0] == 1) {
-  isD0 = 1;   //If both are kaons, leave to top. cuts
-  isD0bar = 1;
-  }
+   if (!checkPIDInfo[0] && !checkPIDInfo[1]) {
+      return 0;      //Reject if both daughters lack both TPC and TOF info
+   }
 
 
-  else {
-  isD0 = 1;
-  isD0bar = 1;
-  if (combinedPID[0][0] == 0 && aodtrack[0]->P() < 1.5 && (fPidHF->CheckStatus(aodtrack[0], "TOF"))) { //If it's low momentum and not ID'd as a kaon, we assume it must be a pion
-  if (aodtrack[0]->Charge() == -1) {
-  isD0 = 0;  //Check charge and reject depending on pion hypothesis
-  }
-  if (aodtrack[0]->Charge() == 1) {
-  isD0bar = 0;
-  }
-  }
-  if (combinedPID[1][0] == 0 && aodtrack[1]->P() < 1.5 && (fPidHF->CheckStatus(aodtrack[1], "TOF"))) {
-  if (aodtrack[1]->Charge() == -1) {
-  isD0 = 0;
-  }
-  if (aodtrack[1]->Charge() == 1) {
-  isD0bar = 0;
-  }
-  }
-  }*/
-  /*    if(combinedPID[0][1] == 1 && combinedPID[1][1] == 1) {         //If both are ID'd as pions:
-       if ((aodtrack[0]->P()) <1.0 && (aodtrack[1]->P() < 1.0)) {     //Check both momenta. If both low, we require one kaon.
-       isD0 == 0; isD0bar == 0;}              //Both pions, both low momenta: reject
-       else {isD0 == 1; isD0bar == 1;}              //Both pions, both high momenta: Allow as both (filter with topological cuts)
+   CalculateBayesianWeights(d);        //Calculates all Bayesian probabilities for both positive and negative tracks
+    //      Double_t prob0[AliPID::kSPECIES];
+    //      fPidHF->GetPidCombined()->ComputeProbabilities(aodtrack[daught], fPidHF->GetPidResponse(), prob0);
+   ///Three possible Bayesian probability cuts: Picked using SetBayesianCondition(int).
+   switch (fBayesianCondition) {
+      ///A: Standard max. probability method (accept most likely species) 
+      case kMaxProb:
+       if (TMath::MaxElement(AliPID::kSPECIES, fWeightsPositive) == fWeightsPositive[AliPID::kKaon]) { //If highest probability lies with kaon
+         isPosKaon = 1;  //flag [daught] as a kaon
        }
-  */
-
-
-}
-
-
-//---------------------------------------------------------------------------
-Int_t AliRDHFCutsD0toKpi::IsSelectedSimpleBayesianPID(AliAODRecoDecayHF* d)
-
-{
-  //Apply Bayesian PID selection; "simple" method.
-  //Looks for a kaon via max. probability. If neither daughter is a kaon, candidate is rejected.
-
-  Int_t isPosKaon = 0, isNegKaon = 0;
-  Int_t returnvalue;          //0 = rejected, 1 = D0, 2 = D0bar. This version will not output 3 (both).
-  Bool_t checkPIDInfo[2] = {kTRUE, kTRUE};
-
-  AliAODTrack *aodtrack[2] = {(AliAODTrack*)d->GetDaughter(0), (AliAODTrack*)d->GetDaughter(1)};
-  if ((aodtrack[0]->Charge() * aodtrack[1]->Charge()) != -1) {
-    return 0;  //Reject if charges do not oppose
-  }
-  for (Int_t daught = 0; daught < 2; daught++) {
-    //Loop over prongs to check PID information
-
-    if (!(fPidHF->CheckStatus(aodtrack[daught], "TPC")) && !(fPidHF->CheckStatus(aodtrack[daught], "TOF"))) {
-      checkPIDInfo[daught] = kFALSE;
-    }
-
-  }
-  //Loop over daughters ends here.
-
-  if (!checkPIDInfo[0] && !checkPIDInfo[1]) {
-    return 0;
-  }
-
-  CalculateBayesianWeights(d);
-
-   
-  // identify kaon
-  ///Three possible Bayesian methods: Pick one (only first implemented so far. TODO: implement other two methods)
-  ///A: Standard max. probability method
-
+       
+       if (TMath::MaxElement(AliPID::kSPECIES, fWeightsPositive) == fWeightsPositive[AliPID::kPion]) { //If highest probability lies with pion
+         isPosPion = 1;  //flag [daught] as a pion
+       }            
+       
+       if (TMath::MaxElement(AliPID::kSPECIES, fWeightsNegative) == fWeightsNegative[AliPID::kKaon]) { //If highest probability lies with kaon
+         isNegKaon = 1;  //flag [daught] as a kaon
+       }
+       
+       if (TMath::MaxElement(AliPID::kSPECIES, fWeightsNegative) == fWeightsNegative[AliPID::kPion]) { //If highest probability lies with kaon
+         isNegPion = 1;  //flag [daught] as a pion
+       }
+       
+       
+       break;
+       ///B: Accept if probability greater than prior
+   case kAbovePrior:
+     
+     if (fWeightsNegative[AliPID::kKaon] > (fPidHF->GetPidCombined()->GetPriorDistribution(AliPID::kKaon)->
+                                           GetBinContent(fPidHF->GetPidCombined()->GetPriorDistribution(AliPID::kKaon)->FindBin(momentumnegative)))) {  //Retrieves relevant prior, gets value at momentum
+       isNegKaon = 1;
+     }
+     if (fWeightsPositive[AliPID::kKaon] > (fPidHF->GetPidCombined()->GetPriorDistribution(AliPID::kKaon)->
+                                           GetBinContent(fPidHF->GetPidCombined()->GetPriorDistribution(AliPID::kKaon)->FindBin(momentumpositive)))) {  //Retrieves relevant prior, gets value at momentum
+       isPosKaon = 1;
+     }
+     
+     break;
+     
+     ///C: Accept if probability greater than user-defined threshold
+      case kThreshold:
+         if (fWeightsNegative[AliPID::kKaon] > fProbThreshold) {
+            isNegKaon = 1;
+         }
+         if (fWeightsNegative[AliPID::kPion] > fProbThreshold) {
+            isNegPion = 1;
+         }
+         
+         if (fWeightsPositive[AliPID::kKaon] > fProbThreshold) {
+            isPosKaon = 1;
+         }
+        
+         if (fWeightsPositive[AliPID::kPion] > fProbThreshold) {
+            isPosPion = 1;
+         }
+
+         break;
+   }
    
-  if (TMath::MaxElement(AliPID::kSPECIES, fWeightsPositive) == fWeightsPositive[AliPID::kKaon]) { //If highest probability lies with kaon
-    isPosKaon = 1;  //flag [daught] as a kaon
-  }
-
-  if (TMath::MaxElement(AliPID::kSPECIES, fWeightsNegative) == fWeightsNegative[AliPID::kKaon]) { //If highest probability lies with kaon
-    isNegKaon = 1;  //flag [daught] as a kaon
-  }
-
-
-  ///B: Method based on probability greater than prior
-  /* Double_t kaonpriorprob;
-     TH1D *kaonprior = fPidHF->GetPidCombined()->GetPriorDistribution(AliPID::kKaon);
-     kaonpriorprob =   kaonprior->GetBinContent(kaonprior->FindBin(aodtrack[daught]->P()));
-
-     if (prob0[AliPID::kKaon] > kaonpriorprob){combinedPID[daught][0] = 1;}
-
-  */
-  ///C: Probability threshold (TODO: Add necessary threshold variables; currently unused)
-  /*
-    if(prob0[AliPID::kKaon] >= fMinProbKaon) {combinedPID[daught][0] = 1;}
-
-  */
-
+     
+     //Momentum-based selection (also applied to filtered weighted method)
+     
+     if (fBayesianStrategy == kBayesMomentum || fBayesianCondition == kBayesWeight) {
+         if (isNegKaon && isPosKaon) { // If both are kaons, reject
+            isD0 = 1;
+            isD0bar = 1;
+         } else if (isNegKaon && isPosPion) {       //If negative kaon present, D0
+            isD0 = 1;
+         } else if (isPosKaon && isNegPion) {       //If positive kaon present, D0bar
+            isD0bar = 1;
+         } else {                      //If neither ID'd as kaon, subject to extra tests
+            isD0 = 1;
+            isD0bar = 1;
+            if (aodtrack[0]->P() < 1.5 && (fPidHF->CheckStatus(aodtrack[0], "TOF"))) { //If it's low momentum and not ID'd as a kaon, we assume it must be a pion
+               if (aodtrack[0]->Charge() == -1) {
+            isD0 = 0;  //Check charge and reject based on pion hypothesis
+               }
+               if (aodtrack[0]->Charge() == 1) {
+            isD0bar = 0;
+               }
+            }
+            if (aodtrack[1]->P() < 1.5 && (fPidHF->CheckStatus(aodtrack[1], "TOF"))) {
+               if (aodtrack[1]->Charge() == -1) {
+            isD0 = 0;
+               }
+               if (aodtrack[1]->Charge() == 1) {
+            isD0bar = 0;
+               }
+            }
+         }
+
+            
+
+         if (isD0 && isD0bar) {
+            returnvalue = 3;
+         }
+         if (isD0 && !isD0bar) {
+            returnvalue = 1;
+         }
+         if (!isD0 && isD0bar) {
+            returnvalue = 2;
+         }
+         if (!isD0 && !isD0bar) {
+            returnvalue = 0;
+         }
+     }
 
+    //Simple Bayesian
+    if (fBayesianStrategy == kBayesSimple) {
+       
+         if (isPosKaon && isNegKaon)   {  //If both are ID'd as kaons, accept as possible
+               returnvalue = 3;
+            } else if (isNegKaon && isPosPion)   {     //If negative kaon, D0
+               returnvalue = 1;
+            } else if (isPosKaon && isNegPion)   {     //If positive kaon, D0-bar
+               returnvalue = 2;
+            } else if (isPosPion && isNegPion)   {
+               returnvalue = 0;  //If neither kaon, reject
+            } else {returnvalue = 0;}  //default
+            
+    }
+    
+  return returnvalue;
 
 
-  if (isPosKaon && isNegKaon)   {  //If both are ID'd as kaons, reject
-    returnvalue = 0;
-  } else if (isNegKaon)   {     //If negative kaon, D0
-    returnvalue = 1;
-  } else if (isPosKaon)   {     //If positive kaon, D0-bar
-    returnvalue = 2;
-  } else {
-    returnvalue = 0;  //If neither kaon, reject
-  }
 
-  return returnvalue;
 }
 
 
 
-
 //---------------------------------------------------------------------------
 void AliRDHFCutsD0toKpi::CalculateBayesianWeights(AliAODRecoDecayHF* d)
 
 {
   //Function to compute weights for Bayesian method
-  Double_t *prob0;
-  prob0 = new Double_t[AliPID::kSPECIES];
-
 
   AliAODTrack *aodtrack[2] = {(AliAODTrack*)d->GetDaughter(0), (AliAODTrack*)d->GetDaughter(1)};
   if ((aodtrack[0]->Charge() * aodtrack[1]->Charge()) != -1) {
@@ -2053,16 +2102,14 @@ void AliRDHFCutsD0toKpi::CalculateBayesianWeights(AliAODRecoDecayHF* d)
     }
 
 
-    // identify kaon
-    fPidHF->GetPidCombined()->ComputeProbabilities(aodtrack[daught], fPidHF->GetPidResponse(), prob0);
-
+    // identify kaon, define weights
     if (aodtrack[daught]->Charge() == +1) {
-      SetWeightsPositive(prob0);
+    fPidHF->GetPidCombined()->ComputeProbabilities(aodtrack[daught], fPidHF->GetPidResponse(), fWeightsPositive);
     }
+    
     if (aodtrack[daught]->Charge() == -1) {
-      SetWeightsNegative(prob0);
+    fPidHF->GetPidCombined()->ComputeProbabilities(aodtrack[daught], fPidHF->GetPidResponse(), fWeightsNegative);
     }
-
   }
 }