STEER/AliESDanalysis.C

   1 //********************************************************************
   2 //     Example (very naive for the moment) of the data analysis
   3 //                    using the ESD classes
   4 //     Origin: Iouri Belikov, CERN, Jouri.Belikov@cern.ch
   5 //********************************************************************
   6
   7 #if !defined( __CINT__) || defined(__MAKECINT__)
   8   #include <Riostream.h>
   9   #include <TTree.h>
  10   #include "TFile.h"
  11   #include "TH1F.h"
  12   #include "TCanvas.h"
  13   #include "TStyle.h"
  14   #include "TStopwatch.h"
  15
  16   #include "AliESD.h"
  17 #endif
  18
  19 extern TStyle *gStyle;
  20
  21 Int_t AliESDanalysis() {
  22    TStopwatch timer;
  23
  24    gStyle->SetOptStat(111110);
  25    gStyle->SetOptFit(1);
  26
  27    Double_t V0mass=0.497672, V0width=0.020, V0window=0.05;
  28    Double_t mmin=V0mass-V0window, mmax=V0mass+V0window;
  29    TH1F *hm =new TH1F("hm","K0s",40, mmin, mmax);
  30    hm->SetXTitle("Mass (GeV/c**2)"); hm->SetLineColor(2);
  31    TH1F *hp =new TH1F("hp","Momentum of the positive daughter",40, 0, 2);
  32    hp->SetXTitle("P (GeV/c)"); hp->SetLineColor(4);
  33
  34 //****** File with the ESD
  35    TFile *ef=TFile::Open("AliESDs.root");
  36    if (!ef || !ef->IsOpen()) {cerr<<"Can't AliESDs.root !\n"; return 1;}
  37    AliESD* event = new AliESD;
  38    TTree* tree = (TTree*) ef->Get("esdTree");
  39    if (!tree) {cerr<<"no ESD tree found\n"; return 1;};
  40    tree->SetBranchAddress("ESD", &event);
  41
  42    Int_t n=0;
  43
  44 //******* The loop over events
  45    while (tree->GetEvent(n)) {
  46      cout<<endl<<"Processing event number : "<<n++<<endl;
  47
  48      Int_t ntrk=event->GetNumberOfTracks();
  49      cout<<"Number of ESD tracks : "<<ntrk<<endl;
  50      Int_t nv0=event->GetNumberOfV0s();
  51      cout<<"Number of ESD V0s : "<<nv0<<endl;
  52      Int_t ncas=event->GetNumberOfCascades();
  53      cout<<"Number of ESD cascades : "<<ncas<<endl;
  54
  55      //****** The loop over tracks
  56      Int_t nk=0;
  57      while (ntrk--) {
  58         AliESDtrack *track=event->GetTrack(ntrk);
  59         UInt_t status=track->GetStatus();
  60
  61         //select only tracks with the "combined PID"
  62         if ((status&AliESDtrack::kESDpid)==0) continue;
  63
  64         Double_t w[10]; track->GetESDpid(w);
  65         //count only "Kaon-like" tracks
  66         if (w[3]>w[4] && w[3]>w[2] && w[3]>w[1] && w[3]>w[0]) nk++;
  67      }
  68      cout<<"Number of \"Kaon-like\" tracks : "<<nk<<endl;
  69
  70      //****** The loop over V0s
  71      while (nv0--) {
  72         AliESDv0 *v0=event->GetV0(nv0);
  73         v0->ChangeMassHypothesis(310); // K0s
  74         Double_t mass=v0->GetEffMass();
  75         hm->Fill(mass);
  76
  77         Int_t pidx=v0->GetPindex();               // now let's get an access
  78         AliESDtrack *track=event->GetTrack(pidx); // to the positive daughter
  79         Double_t p=track->GetP();
  80         hp->Fill(p);
  81      }
  82
  83      //****** The loop over cascades
  84      while (ncas--) {
  85         AliESDcascade *cas=event->GetCascade(ncas);
  86         Double_t q; //"quality" of the associated Lambda
  87         cas->ChangeMassHypothesis(q,3312); // Xi-
  88         // Here you do something with your Xis
  89         //  ...
  90         // You can get the access to the daughters
  91      }
  92
  93    }
  94
  95    delete event;
  96    ef->Close();
  97
  98    timer.Stop(); timer.Print();
  99
 100    TCanvas *c1=new TCanvas("c1","",0,0,600,1200);
 101    c1->Divide(1,2);
 102
 103    c1->cd(1);
 104    hm->Fit("gaus","","",V0mass-V0width,V0mass+V0width);
 105
 106    c1->cd(2);
 107    hp->Fit("expo","","",0.3,2);
 108
 109    return 0;
 110 }