Update User Task

[u/mrichter/AliRoot.git] / STEER / AliESDv0Analysis.C

diff --git a/STEER/AliESDv0Analysis.C b/STEER/AliESDv0Analysis.C
index 102b0586b407ba47950616de3a8b39f715ff1d72..b5d94b49a4a34256ad7f9a6ac4d806366d33a136 100644 (file)
--- a/STEER/AliESDv0Analysis.C
+++ b/STEER/AliESDv0Analysis.C
@@ -8,15 +8,18 @@
 
 #if !defined( __CINT__) || defined(__MAKECINT__)
   #include <Riostream.h>
+  #include <TROOT.h>
   #include <TTree.h>
   #include <TFile.h>
   #include <TH1F.h>
   #include <TCanvas.h>
 
-  #include "AliESD.h"
-
+  #include "AliESDEvent.h"
+  #include "AliESDv0.h"
 #endif
 
+extern TROOT *gROOT;
+
 Int_t AliESDv0Analysis(const Char_t *dir=".") { 
    TH1F *hm=(TH1F*)gROOT->FindObject("hm");
    if (!hm) {
@@ -29,21 +32,22 @@ Int_t AliESDv0Analysis(const Char_t *dir=".") {
    if (!ef||!ef->IsOpen()) {cerr<<"Can't AliESDs.root !\n"; return 1;}
    cerr<<"\n****** "<<fname<<" ******\n";
 
-   AliESD* event = new AliESD;
+   AliESDEvent* event = new AliESDEvent();
 
 
    TTree* tree = (TTree*) ef->Get("esdTree");
    if (!tree) {cerr<<"no ESD tree found\n"; return 1;};
-   tree->SetBranchAddress("ESD", &event);
+   event->ReadFromTree(tree);
 
    Int_t rc=0,n=0;
 
    //****** Tentative particle type "concentrations"
    Double_t c[5]={0.0, 0.0, 1, 0, 1};
+   AliPID pid;
+   pid.SetPriors(c);
 
    //******* The loop over events
-    while (tree->GetEvent(n))
-    {
+    while (tree->GetEvent(n)) {
 
      cerr<<"Processing event number : "<<n++<<endl;
 
@@ -52,13 +56,14 @@ Int_t AliESDv0Analysis(const Char_t *dir=".") {
 
      while (nv0--) {
        AliESDv0 *v0=event->GetV0(nv0);
+       if (v0->GetOnFlyStatus()) continue;
 
        Int_t protonIdx=v0->GetPindex();
        Int_t pionIdx  =v0->GetNindex();
       
        v0->ChangeMassHypothesis(3122);
        Double_t mass=v0->GetEffMass();
-       if (mass>1.17) {  //check also a LambdaBar hypothesis
+       if (mass>1.17) {  //check also the LambdaBar hypothesis
           v0->ChangeMassHypothesis(-3122);
           mass=v0->GetEffMass();
           if (mass>1.17) continue;
@@ -73,36 +78,24 @@ 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);
-         Double_t rcc=0.;
-         Int_t i;
-         for (i=0; i<AliESDtrack::kSPECIES; i++) rcc+=(c[i]*r[i]);
-         if (rcc==0.) continue;
-         //Here we apply Bayes' formula
-         Double_t w[10];
-         for (i=0; i<AliESDtrack::kSPECIES; i++) w[i]=c[i]*r[i]/rcc;
-
-         if (w[4]<w[3]) continue;
-         if (w[4]<w[2]) continue;
-         if (w[4]<w[1]) continue;
-        if (w[4]<w[0]) continue;
+         pid.SetProbabilities(r);
+         Double_t pp=pid.GetProbability(AliPID::kProton);
+         if (pp < pid.GetProbability(AliPID::kElectron)) continue;
+         if (pp < pid.GetProbability(AliPID::kMuon)) continue;
+         if (pp < pid.GetProbability(AliPID::kPion)) continue;
+         if (pp < pid.GetProbability(AliPID::kKaon)) continue;
        }
  
        //Check if the "pion track" is a pion
        if ((pionTrk->GetStatus()&AliESDtrack::kESDpid)!=0) {
         Double_t r[10]; pionTrk->GetESDpid(r);
-         Double_t rcc=0.;
-         Int_t i;
-         for (i=0; i<AliESDtrack::kSPECIES; i++) rcc+=(c[i]*r[i]);
-         if (rcc==0.) continue;
-         //Here we apply Bayes' formula
-         Double_t w[10];
-         for (i=0; i<AliESDtrack::kSPECIES; i++) w[i]=c[i]*r[i]/rcc;
-
-         if (w[2]<w[4]) continue;
-         if (w[2]<w[3]) continue;
-         if (w[2]<w[1]) continue;
-        if (w[2]<w[0]) continue;
-       }
+         pid.SetProbabilities(r);
+         Double_t ppi=pid.GetProbability(AliPID::kPion);
+         if (ppi < pid.GetProbability(AliPID::kElectron)) continue;
+         if (ppi < pid.GetProbability(AliPID::kMuon)) continue;
+         if (ppi < pid.GetProbability(AliPID::kKaon)) continue;
+         if (ppi < pid.GetProbability(AliPID::kProton)) continue;
+        }
 
        hm->Fill(mass);
      } 
@@ -110,6 +103,7 @@ Int_t AliESDv0Analysis(const Char_t *dir=".") {
    }
 
    delete event;
+   delete tree;
    ef->Close();
 
    TCanvas *c1=(TCanvas*)gROOT->FindObject("c1");