Added methods for finding a given module in the DDL map (F. Prino)
[u/mrichter/AliRoot.git] / STEER / AliESDv0Analysis.C
CommitLineData
e23730c7 1//********************************************************************
2// Example (very naive for the moment) of the data analysis
3// using the ESD classes.
4// It demonstrates the idea of the "combined PID"
5// applied to the Lambda0 reconstruction.
6// Origin: Iouri Belikov, CERN, Jouri.Belikov@cern.ch
7//********************************************************************
8
9#if !defined( __CINT__) || defined(__MAKECINT__)
10 #include <Riostream.h>
3088e2dd 11 #include <TROOT.h>
8b462fd8 12 #include <TTree.h>
9d0d65be 13 #include <TFile.h>
14 #include <TH1F.h>
15 #include <TCanvas.h>
e23730c7 16
86ad5fcb 17 #include "AliESDEvent.h"
b4eb042a 18 #include "AliESDv0.h"
e23730c7 19#endif
20
3088e2dd 21extern TROOT *gROOT;
22
9d0d65be 23Int_t AliESDv0Analysis(const Char_t *dir=".") {
24 TH1F *hm=(TH1F*)gROOT->FindObject("hm");
25 if (!hm) {
26 hm=new TH1F("hm","Lambda+LambdaBar Effective Mass",60,1.065,1.165);
27 hm->SetXTitle("Mass (GeV/c**2)");
28 }
29 Char_t fname[100];
30 sprintf(fname,"%s/AliESDs.root",dir);
31 TFile *ef=TFile::Open(fname);
32 if (!ef||!ef->IsOpen()) {cerr<<"Can't AliESDs.root !\n"; return 1;}
33 cerr<<"\n****** "<<fname<<" ******\n";
e23730c7 34
4514157e 35 AliESDEvent* event = new AliESDEvent();
9d0d65be 36
37
8b462fd8 38 TTree* tree = (TTree*) ef->Get("esdTree");
39 if (!tree) {cerr<<"no ESD tree found\n"; return 1;};
86ad5fcb 40 event->ReadFromTree(tree);
e23730c7 41
e23730c7 42 Int_t rc=0,n=0;
e23730c7 43
44 //****** Tentative particle type "concentrations"
9d0d65be 45 Double_t c[5]={0.0, 0.0, 1, 0, 1};
3088e2dd 46 AliPID pid;
47 pid.SetPriors(c);
e23730c7 48
49 //******* The loop over events
4ca848eb 50 while (tree->GetEvent(n)) {
e23730c7 51
52 cerr<<"Processing event number : "<<n++<<endl;
53
e23730c7 54 Int_t nv0=event->GetNumberOfV0s();
55 cerr<<"Number of ESD v0s : "<<nv0<<endl;
56
57 while (nv0--) {
58 AliESDv0 *v0=event->GetV0(nv0);
9d0d65be 59
60 Int_t protonIdx=v0->GetPindex();
61 Int_t pionIdx =v0->GetNindex();
62
63 v0->ChangeMassHypothesis(3122);
64 Double_t mass=v0->GetEffMass();
4ca848eb 65 if (mass>1.17) { //check also the LambdaBar hypothesis
9d0d65be 66 v0->ChangeMassHypothesis(-3122);
67 mass=v0->GetEffMass();
68 if (mass>1.17) continue;
69 Int_t tmp=protonIdx; protonIdx=pionIdx; pionIdx=tmp;
70 }
71
72 AliESDtrack *protonTrk=event->GetTrack(protonIdx);
73 AliESDtrack *pionTrk =event->GetTrack(pionIdx);
74
75 if (protonTrk->GetP()<0.5) continue;
76
77 // Check if the "proton track" is a proton
78 if ((protonTrk->GetStatus()&AliESDtrack::kESDpid)!=0) {
79 Double_t r[10]; protonTrk->GetESDpid(r);
3088e2dd 80 pid.SetProbabilities(r);
4ca848eb 81 Double_t pp=pid.GetProbability(AliPID::kProton);
82 if (pp < pid.GetProbability(AliPID::kElectron)) continue;
83 if (pp < pid.GetProbability(AliPID::kMuon)) continue;
84 if (pp < pid.GetProbability(AliPID::kPion)) continue;
85 if (pp < pid.GetProbability(AliPID::kKaon)) continue;
e23730c7 86 }
9d0d65be 87
88 //Check if the "pion track" is a pion
89 if ((pionTrk->GetStatus()&AliESDtrack::kESDpid)!=0) {
90 Double_t r[10]; pionTrk->GetESDpid(r);
3088e2dd 91 pid.SetProbabilities(r);
4ca848eb 92 Double_t ppi=pid.GetProbability(AliPID::kPion);
93 if (ppi < pid.GetProbability(AliPID::kElectron)) continue;
94 if (ppi < pid.GetProbability(AliPID::kMuon)) continue;
95 if (ppi < pid.GetProbability(AliPID::kKaon)) continue;
96 if (ppi < pid.GetProbability(AliPID::kProton)) continue;
97 }
9d0d65be 98
e23730c7 99 hm->Fill(mass);
100 }
9d0d65be 101
e23730c7 102 }
103
8b462fd8 104 delete event;
4514157e 105 delete tree;
8b462fd8 106 ef->Close();
107
9d0d65be 108 TCanvas *c1=(TCanvas*)gROOT->FindObject("c1");
109 if (!c1) {
110 c1=new TCanvas();
111 }
e23730c7 112
9d0d65be 113 if (hm->GetEntries()>100) hm->Fit("gaus","","",1.11,1.12);
114 else hm->Draw();
e23730c7 115
9d0d65be 116 c1->Update();
e23730c7 117
118 return rc;
119}