//_________________________________________________________________________________________________
AliRsnCutPID::AliRsnCutPID() :
- AliRsnCut(AliRsnCut::kDaughter),
- fPerfect(kFALSE),
- fUseDefault(kTRUE)
+ AliRsnCut(),
+ fPerfect(kFALSE),
+ fUseDefault(kTRUE)
{
//
// Default constructor.
// to which we want to compare this object.
//
- Int_t i;
-
- for (i = 0; i < kDetectors; i++)
- {
- fUseDetector[i] = kFALSE;
- fPtThreshold[i] = 0.0;
- fGoAboveThreshold[i] = kTRUE;
- }
-
- for (i = 0; i < AliPID::kSPECIES; i++)
- {
- fWeight[i] = 0.0;
- fPrior[i] = 1.0;
- }
+ Int_t i;
+
+ for (i = 0; i < kDetectors; i++) {
+ fUseDetector[i] = kFALSE;
+ fPtThreshold[i] = 0.0;
+ fGoAboveThreshold[i] = kTRUE;
+ }
+
+ for (i = 0; i < AliPID::kSPECIES; i++) {
+ fWeight[i] = 0.0;
+ fPrior[i] = 1.0;
+ }
+
+ SetTargetType(AliRsnTarget::kDaughter);
}
//_________________________________________________________________________________________________
AliRsnCutPID::AliRsnCutPID(const char *name, AliPID::EParticleType pid, Double_t probMin, Bool_t perfectPID) :
- AliRsnCut(name, AliRsnCut::kDaughter, (Int_t)pid),
- fPerfect(perfectPID),
- fUseDefault(kTRUE)
+ AliRsnCut(name, AliRsnCut::kDaughter, (Int_t)pid),
+ fPerfect(perfectPID),
+ fUseDefault(kTRUE)
{
//
// Default constructor.
// to which we want to compare this object.
//
- Int_t i;
-
- for (i = 0; i < kDetectors; i++)
- {
- fUseDetector[i] = kFALSE;
- fPtThreshold[i] = 0.0;
- fGoAboveThreshold[i] = kTRUE;
- }
-
- for (i = 0; i < AliPID::kSPECIES; i++)
- {
- fWeight[i] = 0.0;
- fPrior[i] = 1.0;
- }
-
- fMinD = probMin;
- fMaxD = 1.000001;
+ Int_t i;
+
+ for (i = 0; i < kDetectors; i++) {
+ fUseDetector[i] = kFALSE;
+ fPtThreshold[i] = 0.0;
+ fGoAboveThreshold[i] = kTRUE;
+ }
+
+ for (i = 0; i < AliPID::kSPECIES; i++) {
+ fWeight[i] = 0.0;
+ fPrior[i] = 1.0;
+ }
+
+ fMinD = probMin;
+ fMaxD = 1.000001;
+
+ SetTargetType(AliRsnTarget::kDaughter);
}
//_____________________________________________________________________________
-Bool_t AliRsnCutPID::CheckThreshold(EDetector det, Double_t value)
+Bool_t AliRsnCutPID::CheckThreshold(EDetector det, Double_t value) const
{
//
-// Checks if the passed value (track pT) stays in the
+// Checks if the passed value (track pT) stays in the
// interval where the detector should be accepted
//
- if (!CheckBounds(det)) return kFALSE;
-
- if (fGoAboveThreshold[det]) return (value >= fPtThreshold[det]);
- else return (value <= fPtThreshold[det]);
+ if (!CheckBounds(det)) return kFALSE;
+
+ if (fGoAboveThreshold[det]) return (value >= fPtThreshold[det]);
+ else return (value <= fPtThreshold[det]);
}
//_________________________________________________________________________________________________
// The return value tells if the operation was successful.
//
- Int_t i, j;
- Bool_t useDefault = fUseDefault;
- Bool_t perfectPID = fPerfect;
- if (perfectPID && !daughter->GetRefMC()) perfectPID = kFALSE;
- if (!daughter->GetRefESDtrack()) useDefault = kTRUE;
- if (!daughter->GetRefESDtrack() && !daughter->GetRefAODtrack()) return kFALSE;
-
- // if perfect PID ise required, this overcomes all
- // in this case the weight of the correct species is set to 1
- // and the others to 0
- // of course this happens only if there is a reference MC
- if (perfectPID)
- {
- j = TMath::Abs(daughter->GetRefMC()->Particle()->GetPdgCode());
- for (i = 0; i < AliPID::kSPECIES; i++)
- {
- if (AliPID::ParticleCode((AliPID::EParticleType)i) == j) fWeight[i] = 1.0;
- else fWeight[i] = 0.0;
- }
- return kTRUE;
- }
-
- // if default weights are (or need to be) used,
- // they are taken directly and function exits
- if (useDefault)
- {
- if (daughter->GetRefESDtrack())
- daughter->GetRefESDtrack()->GetESDpid(fWeight);
- else
- {
- for (i = 0; i < AliPID::kSPECIES; i++)
- fWeight[i] = daughter->GetRefAODtrack()->PID()[i];
- }
- return kTRUE;
- }
-
- // if we arrive here, this means that we have an ESD track
- // and we want to customize the PID
- AliESDtrack *track = daughter->GetRefESDtrack();
- Double_t w[kDetectors][AliPID::kSPECIES];
- track->GetITSpid(w[kITS]);
- track->GetTPCpid(w[kTPC]);
- track->GetTRDpid(w[kTRD]);
- track->GetTOFpid(w[kTOF]);
- track->GetHMPIDpid(w[kHMPID]);
-
- // if a detector is excluded or the track has not the right pT
- // all related weights are set to 1 in order not to contribute
- // to final multiplication
- for (i = 0; i < kDetectors; i++)
- {
- if (!fUseDetector[i] || !CheckThreshold((EDetector)i, track->Pt()))
- {
- for (j = 0; j < AliPID::kSPECIES; j++) {
- w[i][j] = 1.0;
+ Int_t i, j;
+ Bool_t useDefault = fUseDefault;
+ Bool_t perfectPID = fPerfect;
+ if (perfectPID && !daughter->GetRefMC()) return kFALSE;
+ if (!daughter->Ref2ESDtrack()) useDefault = kTRUE;
+ if (!daughter->Ref2ESDtrack() && !daughter->Ref2AODtrack()) return kFALSE;
+
+ // if perfect PID ise required,
+ // compare the PDG code of the type stored in 'fMinI' of the cut
+ // and that of the particle which is checked, looking at its MC
+ if (perfectPID) {
+ i = TMath::Abs(AliPID::ParticleCode(fMinI));
+ j = daughter->GetPDG();
+ return (i == j);
+ }
+
+ // if default weights are (or need to be) used,
+ // they are taken directly and function exits
+ if (useDefault) {
+ if (daughter->Ref2ESDtrack())
+ daughter->Ref2ESDtrack()->GetESDpid(fWeight);
+ else {
+ for (i = 0; i < AliPID::kSPECIES; i++)
+ fWeight[i] = daughter->Ref2AODtrack()->PID()[i];
+ }
+ return kTRUE;
+ }
+
+ // if we arrive here, this means that we have an ESD track
+ // and we want to customize the PID
+ AliESDtrack *track = daughter->Ref2ESDtrack();
+ Double_t w[kDetectors][AliPID::kSPECIES];
+ track->GetITSpid(w[kITS]);
+ track->GetTPCpid(w[kTPC]);
+ track->GetTRDpid(w[kTRD]);
+ track->GetTOFpid(w[kTOF]);
+ track->GetHMPIDpid(w[kHMPID]);
+
+ // if a detector is excluded or the track has not the right pT
+ // all related weights are set to 1 in order not to contribute
+ // to final multiplication
+ for (i = 0; i < kDetectors; i++) {
+ if (!fUseDetector[i] || !CheckThreshold((EDetector)i, track->Pt())) {
+ for (j = 0; j < AliPID::kSPECIES; j++) {
+ w[i][j] = 1.0;
+ }
+ }
+ }
+
+ // multiplicate all weights to compute final one
+ for (i = 0; i < AliPID::kSPECIES; i++) {
+ fWeight[i] = w[kITS][i] * w[kTPC][i] * w[kTRD][i] * w[kTOF][i] * w[kHMPID][i];
+ }
+
+ return kTRUE;
+}
+
+//_________________________________________________________________________________________________
+Int_t AliRsnCutPID::RealisticPID(AliRsnDaughter * const daughter, Double_t &prob)
+{
+//
+// Combines the weights collected from chosen detectors with the priors
+// and gives the corresponding particle with the largest probability,
+// in terms of the AliPID particle type enumeration.
+// The argument, passed by reference, gives the corresponding probability,
+// since the cut could require that it is larger than a threshold.
+//
+
+ // try to compute the weights
+ if (!ComputeWeights(daughter)) {
+ prob = -1.0;
+ return AliPID::kUnknown;
+ }
+
+ // combine with priors and normalize
+ Int_t i;
+ Double_t sum = 0.0, w[AliPID::kSPECIES];
+ for (i = 0; i < AliPID::kSPECIES; i++) {
+ w[i] = fWeight[i] * fPrior[i];
+ sum += w[i];
+ }
+ if (sum <= 0.0) {
+ AliError("Sum = 0");
+ prob = -1.0;
+ return AliPID::kUnknown;
+ }
+ for (i = 0; i < AliPID::kSPECIES; i++) w[i] /= sum;
+
+ // find the largest probability and related PID
+ Int_t ibest = 0;
+ prob = w[0];
+ for (i = 1; i < AliPID::kSPECIES; i++) {
+ if (w[i] > prob) {
+ prob = w[i];
+ ibest = i;
}
- }
- }
-
- // multiplicate all weights to compute final one
- for (i = 0; i < AliPID::kSPECIES; i++)
- {
- fWeight[i] = w[kITS][i] * w[kTPC][i] * w[kTRD][i] * w[kTOF][i] * w[kHMPID][i];
- }
-
- return kTRUE;
+ }
+
+ // return the value, while the probability
+ // will be consequentially set
+ return ibest;
+}
+
+//_________________________________________________________________________________________________
+Int_t AliRsnCutPID::PerfectPID(AliRsnDaughter * const daughter)
+{
+//
+// If MC is present, retrieve the particle corresponding to the used track
+// (using the fRefMC data member) and then return the true particle type
+// taken from the PDG code of the reference particle itself, converted
+// into the enumeration from AliPID object.
+//
+
+ // works only if the MC is present
+ if (!daughter->GetRefMC()) return AliPID::kUnknown;
+
+ // get the PDG code of the particle
+ Int_t pdg = daughter->GetPDG();
+
+ // loop over all species listed in AliPID to find the match
+ Int_t i;
+ for (i = 0; i < AliPID::kSPECIES; i++) {
+ if (AliPID::ParticleCode(i) == pdg) return i;
+ }
+
+ return AliPID::kUnknown;
}
//_________________________________________________________________________________________________
-Bool_t AliRsnCutPID::IsSelected(TObject *obj1, TObject* /*obj2*/)
+Bool_t AliRsnCutPID::IsSelected(TObject *object)
{
//
// Cut checker.
//
- // coherence check
- if (!AliRsnCut::TargetOK(obj1))
- {
- AliError(Form("Wrong target. Skipping cut", GetName()));
- return kTRUE;
- }
-
- // convert the object into the unique correct type
- AliRsnDaughter *daughter = dynamic_cast<AliRsnDaughter*>(obj1);
-
- // try to compute the weights
- if (!ComputeWeights(daughter)) return kFALSE;
-
- // combine with priors and get the majority
- Int_t i;
- Double_t sum = 0.0, w[AliPID::kSPECIES];
- for (i = 0; i < AliPID::kSPECIES; i++)
- {
- w[i] = fWeight[i] * fPrior[i];
- sum += w[i];
- }
- if (sum <= 0.0)
- {
- AliError("Sum = 0");
- return kFALSE;
- }
- for (i = 0; i < AliPID::kSPECIES; i++) w[i] /= sum;
-
- // find the largest probability and related PID
- // and assign them to the mother class members which
- // are checked for the cut
- fCutValueI = 0;
- fCutValueD = w[0];
- for (i = 1; i < AliPID::kSPECIES; i++)
- {
- if (w[i] > fCutValueD)
- {
- fCutValueD = w[i];
- fCutValueI = i;
- }
- }
-
- // if the best probability is too small, the cut is failed anyway
- if (!OkRangeD()) return kFALSE;
-
- // if the best probability is OK, the cut is passed
- // if it correspond to the right particle
- return OkValue();
+ // convert the object into the unique correct type
+
+ if (!TargetOK(object)) {
+ AliError(Form("[%s]: this cut works only with AliRsnDaughter objects", GetName()));
+ return kTRUE;
+ }
+
+ // if target is OK, do a dynamic cast
+ AliRsnDaughter *daughter = fDaughter;
+
+ // depending on the PID type, recalls the appropriate method:
+ // in case of perfect PID, checks only if the PID type is
+ // corresponding to the request in the cut (fMinI)
+ // while in case of realistic PID checks also the probability
+ // to be within the required interval
+ if (fPerfect && daughter) {
+ fCutValueI = PerfectPID(daughter);
+ return OkValueI();
+ } else if (daughter) {
+ fCutValueI = RealisticPID(daughter, fCutValueD);
+ return OkValueI() && OkRangeD();
+ } else
+ return kFALSE;
}
+//__________________________________________________________________________________________________
void AliRsnCutPID::IncludeDetector(EDetector det, Double_t threshold, Bool_t goAbove)
{
//
// By default the threshold is zero and detector is always used.
//
- if (!CheckBounds(det)) return;
-
- fUseDetector[det] = kTRUE;
- fPtThreshold[det] = threshold;
- fGoAboveThreshold[det] = goAbove;
+ if (!CheckBounds(det)) return;
+
+ fUseDetector[det] = kTRUE;
+ fPtThreshold[det] = threshold;
+ fGoAboveThreshold[det] = goAbove;
}