Adding a setter for the probabiliy densities (A.Mastroserio)
authorbelikov <belikov@f7af4fe6-9843-0410-8265-dc069ae4e863>
Wed, 9 May 2007 08:41:05 +0000 (08:41 +0000)
committerbelikov <belikov@f7af4fe6-9843-0410-8265-dc069ae4e863>
Wed, 9 May 2007 08:41:05 +0000 (08:41 +0000)
STEER/AliESDv0Analysis.C
STEER/AliPID.cxx
STEER/AliPID.h

index ee7406e937c3daa0610f40b3ec69087380008b05..c5967693a6c43d51c9880581be996529a4d2c935 100644 (file)
@@ -8,6 +8,7 @@
 
 #if !defined( __CINT__) || defined(__MAKECINT__)
   #include <Riostream.h>
 
 #if !defined( __CINT__) || defined(__MAKECINT__)
   #include <Riostream.h>
+  #include <TROOT.h>
   #include <TTree.h>
   #include <TFile.h>
   #include <TH1F.h>
   #include <TTree.h>
   #include <TFile.h>
   #include <TH1F.h>
@@ -17,6 +18,8 @@
 
 #endif
 
 
 #endif
 
+extern TROOT *gROOT;
+
 Int_t AliESDv0Analysis(const Char_t *dir=".") { 
    TH1F *hm=(TH1F*)gROOT->FindObject("hm");
    if (!hm) {
 Int_t AliESDv0Analysis(const Char_t *dir=".") { 
    TH1F *hm=(TH1F*)gROOT->FindObject("hm");
    if (!hm) {
@@ -40,7 +43,8 @@ Int_t AliESDv0Analysis(const Char_t *dir=".") {
 
    //****** Tentative particle type "concentrations"
    Double_t c[5]={0.0, 0.0, 1, 0, 1};
 
    //****** Tentative particle type "concentrations"
    Double_t c[5]={0.0, 0.0, 1, 0, 1};
-   AliPID::SetPriors(c);
+   AliPID pid;
+   pid.SetPriors(c);
 
    //******* The loop over events
     while (tree->GetEvent(n)) {
 
    //******* The loop over events
     while (tree->GetEvent(n)) {
@@ -73,7 +77,7 @@ Int_t AliESDv0Analysis(const Char_t *dir=".") {
        // Check if the "proton track" is a proton
        if ((protonTrk->GetStatus()&AliESDtrack::kESDpid)!=0) {
         Double_t r[10]; protonTrk->GetESDpid(r);
        // Check if the "proton track" is a proton
        if ((protonTrk->GetStatus()&AliESDtrack::kESDpid)!=0) {
         Double_t r[10]; protonTrk->GetESDpid(r);
-         AliPID pid(r);
+         pid.SetProbabilities(r);
          Double_t pp=pid.GetProbability(AliPID::kProton);
          if (pp < pid.GetProbability(AliPID::kElectron)) continue;
          if (pp < pid.GetProbability(AliPID::kMuon)) continue;
          Double_t pp=pid.GetProbability(AliPID::kProton);
          if (pp < pid.GetProbability(AliPID::kElectron)) continue;
          if (pp < pid.GetProbability(AliPID::kMuon)) continue;
@@ -84,7 +88,7 @@ Int_t AliESDv0Analysis(const Char_t *dir=".") {
        //Check if the "pion track" is a pion
        if ((pionTrk->GetStatus()&AliESDtrack::kESDpid)!=0) {
         Double_t r[10]; pionTrk->GetESDpid(r);
        //Check if the "pion track" is a pion
        if ((pionTrk->GetStatus()&AliESDtrack::kESDpid)!=0) {
         Double_t r[10]; pionTrk->GetESDpid(r);
-         AliPID pid(r);
+         pid.SetProbabilities(r);
          Double_t ppi=pid.GetProbability(AliPID::kPion);
          if (ppi < pid.GetProbability(AliPID::kElectron)) continue;
          if (ppi < pid.GetProbability(AliPID::kMuon)) continue;
          Double_t ppi=pid.GetProbability(AliPID::kPion);
          if (ppi < pid.GetProbability(AliPID::kElectron)) continue;
          if (ppi < pid.GetProbability(AliPID::kMuon)) continue;
index a22c09a446b3948834dd29d40a6986a42b62bbf5..4238ac1080ec918c90b856f523a506cc48749e05 100644 (file)
@@ -158,6 +158,19 @@ AliPID::AliPID(const AliPID& pid) :
     fProbDensity[i] = pid.fProbDensity[i];
 }
 
     fProbDensity[i] = pid.fProbDensity[i];
 }
 
+//_______________________________________________________________________
+void AliPID::SetProbabilities(const Double_t* probDensity, Bool_t charged) 
+{
+  //
+  // Set the probability densities
+  //
+  for (Int_t i = 0; i < kSPECIES; i++) 
+    fProbDensity[i] = probDensity[i];
+
+  for (Int_t i = kSPECIES; i < kSPECIESN; i++) 
+    fProbDensity[i] = ((charged) ? 0 : probDensity[i]);
+}
+
 //_______________________________________________________________________
 AliPID& AliPID::operator = (const AliPID& pid)
 {
 //_______________________________________________________________________
 AliPID& AliPID::operator = (const AliPID& pid)
 {
index 9c9494e4dfe5aef86ee17c3b7eb13fb2161b4ecd..755cd5956f8cee408a3324f5436fcbe33e40ca49 100644 (file)
@@ -57,6 +57,9 @@ class AliPID : public TObject {
   EParticleType        GetMostProbable(const Double_t* prior) const;
   EParticleType        GetMostProbable() const;
   
   EParticleType        GetMostProbable(const Double_t* prior) const;
   EParticleType        GetMostProbable() const;
   
+  void                 SetProbabilities(const Double_t* probabilities,
+                                        Bool_t charged = kTRUE);
+
   static void          SetPriors(const Double_t* prior,
                                 Bool_t charged = kTRUE);
   static void          SetPrior(EParticleType iType, Double_t prior);
   static void          SetPriors(const Double_t* prior,
                                 Bool_t charged = kTRUE);
   static void          SetPrior(EParticleType iType, Double_t prior);