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
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) :
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
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;
}
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);
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);
}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;
esdTrackCuts->SetPtRange(0.3,1.e10);
AddTrackCuts(esdTrackCuts);
-
+ delete esdTrackCuts;
+ esdTrackCuts=NULL;
const Int_t nptbins =14;
const Double_t ptmax = 9999.;
}
+//___________________________________________________________________________
+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() {
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;
AddTrackCuts(esdTrackCuts);
+ delete esdTrackCuts;
+ esdTrackCuts=NULL;
// SPD k FIRST for Pt<3 GeV/c
SetSelectCandTrackSPDFirst(kTRUE, 3);
AddTrackCuts(esdTrackCuts);
+ delete esdTrackCuts;
+ esdTrackCuts=NULL;
// SPD k FIRST for Pt<3 GeV/c
SetSelectCandTrackSPDFirst(kTRUE, 3);
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) {
}
- // 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);
}
-
}
}