f9771dda8eaa135b2b2ef707b025e79794f2d700
[u/mrichter/AliRoot.git] / PWG0 / dNdEta / AlidNdEtaAnalysisESDSelector.cxx
1 /* $Id$ */
2
3 #include "AlidNdEtaAnalysisESDSelector.h"
4
5 #include <TStyle.h>
6 #include <TSystem.h>
7 #include <TCanvas.h>
8 #include <TVector3.h>
9 #include <TChain.h>
10 #include <TFile.h>
11
12 #include <AliLog.h>
13 #include <AliESDVertex.h>
14 #include <AliESD.h>
15
16 #include "esdTrackCuts/AliESDtrackCuts.h"
17 #include "dNdEta/dNdEtaAnalysis.h"
18 #include "dNdEta/AlidNdEtaCorrection.h"
19 #include "AliPWG0Helper.h"
20
21 ClassImp(AlidNdEtaAnalysisESDSelector)
22
23 AlidNdEtaAnalysisESDSelector::AlidNdEtaAnalysisESDSelector() :
24   AliSelector(),
25   fdNdEtaAnalysisMBVtx(0),
26   fdNdEtaAnalysisMB(0),
27   fdNdEtaAnalysis(0),
28   fEsdTrackCuts(0),
29   fdNdEtaCorrection(0)
30 {
31   //
32   // Constructor. Initialization of pointers
33   //
34
35   AliLog::SetClassDebugLevel("AlidNdEtaAnalysisESDSelector", AliLog::kDebug);
36 }
37
38 AlidNdEtaAnalysisESDSelector::~AlidNdEtaAnalysisESDSelector()
39 {
40   //
41   // Destructor
42   //
43
44   // histograms are in the output list and deleted when the output
45   // list is deleted by the TSelector dtor
46 }
47
48 void AlidNdEtaAnalysisESDSelector::Begin(TTree* tree)
49 {
50   // Begin function
51
52   ReadUserObjects(tree);
53 }
54
55 void AlidNdEtaAnalysisESDSelector::ReadUserObjects(TTree* tree)
56 {
57   // read the user objects, called from slavebegin and begin
58
59   if (!fEsdTrackCuts && fInput)
60     fEsdTrackCuts = dynamic_cast<AliESDtrackCuts*> (fInput->FindObject("AliESDtrackCuts"));
61
62   if (!fEsdTrackCuts && tree)
63     fEsdTrackCuts = dynamic_cast<AliESDtrackCuts*> (tree->GetUserInfo()->FindObject("AliESDtrackCuts"));
64
65   if (!fEsdTrackCuts)
66      AliDebug(AliLog::kError, "ERROR: Could not read EsdTrackCuts from input list.");
67
68
69   if (!fdNdEtaCorrection && fInput)
70     fdNdEtaCorrection = dynamic_cast<AlidNdEtaCorrection*> (fInput->FindObject("dndeta_correction"));
71
72   if (!fdNdEtaCorrection && tree)
73     fdNdEtaCorrection = dynamic_cast<AlidNdEtaCorrection*> (tree->GetUserInfo()->FindObject("dndeta_correction"));
74
75   if (!fdNdEtaCorrection)
76      AliDebug(AliLog::kError, "ERROR: Could not read dndeta_correction from input list.");
77 }
78
79 void AlidNdEtaAnalysisESDSelector::SlaveBegin(TTree* tree)
80 {
81   // The SlaveBegin() function is called after the Begin() function.
82   // When running with PROOF SlaveBegin() is called on each slave server.
83   // The tree argument is deprecated (on PROOF 0 is passed).
84
85   AliSelector::SlaveBegin(tree);
86
87   ReadUserObjects(tree);
88
89   fdNdEtaAnalysisMBVtx = new dNdEtaAnalysis("dndeta_mbvtx", "dndeta_mbvtx");
90   fdNdEtaAnalysisMB = new dNdEtaAnalysis("dndeta_mb", "dndeta_mb");
91   fdNdEtaAnalysis = new dNdEtaAnalysis("dndeta", "dndeta");
92 }
93
94 void AlidNdEtaAnalysisESDSelector::Init(TTree* tree)
95 {
96   // read the user objects
97
98   AliSelector::Init(tree);
99 }
100
101 Bool_t AlidNdEtaAnalysisESDSelector::Process(Long64_t entry)
102 {
103   // The Process() function is called for each entry in the tree (or possibly
104   // keyed object in the case of PROOF) to be processed. The entry argument
105   // specifies which entry in the currently loaded tree is to be processed.
106   // It can be passed to either TTree::GetEntry() or TBranch::GetEntry()
107   // to read either all or the required parts of the data. When processing
108   // keyed objects with PROOF, the object is already loaded and is available
109   // via the fObject pointer.
110   //
111   // This function should contain the "body" of the analysis. It can contain
112   // simple or elaborate selection criteria, run algorithms on the data
113   // of the event and typically fill histograms.
114
115   // WARNING when a selector is used with a TChain, you must use
116   //  the pointer to the current TTree to call GetEntry(entry).
117   //  The entry is always the local entry number in the current tree.
118   //  Assuming that fTree is the pointer to the TChain being processed,
119   //  use fTree->GetTree()->GetEntry(entry).
120
121   if (AliSelector::Process(entry) == kFALSE)
122     return kFALSE;
123
124   // Check prerequisites
125   if (!fESD)
126   {
127     AliDebug(AliLog::kError, "ESD branch not available");
128     return kFALSE;
129   }
130
131   if (!fEsdTrackCuts)
132   {
133     AliDebug(AliLog::kError, "fESDTrackCuts not available");
134     return kFALSE;
135   }
136
137   if (!fdNdEtaCorrection)
138   {
139     AliDebug(AliLog::kError, "fdNdEtaCorrection not available");
140     return kFALSE;
141   }
142
143   if (AliPWG0Helper::IsEventTriggered(fESD) == kFALSE)
144     return kTRUE;
145
146   if (AliPWG0Helper::IsVertexReconstructed(fESD) == kFALSE)
147     return kTRUE;
148
149   // ########################################################
150   // get the EDS vertex
151   const AliESDVertex* vtxESD = fESD->GetVertex();
152   Double_t vtx[3];
153   vtxESD->GetXYZ(vtx);
154
155   // get number of "good" tracks
156   TObjArray* list = fEsdTrackCuts->GetAcceptedTracks(fESD);
157   Int_t nGoodTracks = list->GetEntries();
158
159   Float_t vertexRecoCorr = fdNdEtaCorrection->GetVertexRecoCorrection(vtx[2], nGoodTracks);
160   if (vertexRecoCorr <= 0)
161   {
162     AliDebug(AliLog::kError, Form("INFO: Skipping event because vertexRecoCorr is <= 0 (%f)", vertexRecoCorr));
163     delete list;
164     return kTRUE;
165   }
166
167   Float_t triggerCorr = fdNdEtaCorrection->GetTriggerCorrection(vtx[2], nGoodTracks);
168   if (triggerCorr <= 0)
169   {
170     AliDebug(AliLog::kError, Form("INFO: Skipping event because triggerCorr is <= 0 (%f)", triggerCorr));
171     delete list;
172     return kTRUE;
173   }
174
175   // loop over esd tracks
176   for (Int_t t=0; t<nGoodTracks; t++)
177   {
178     AliESDtrack* esdTrack = dynamic_cast<AliESDtrack*> (list->At(t));
179     if (!esdTrack)
180     {
181       AliDebug(AliLog::kError, Form("ERROR: Could not retrieve track %d.", t));
182       continue;
183     }
184
185     Double_t p[3];
186     esdTrack->GetConstrainedPxPyPz(p); // ### TODO or GetInnerPxPyPy / GetOuterPxPyPy
187     TVector3 vector(p);
188
189     Float_t theta = vector.Theta();
190     Float_t eta   = -TMath::Log(TMath::Tan(theta/2.));
191     Float_t pt = vector.Pt();
192
193     Float_t track2particleCorr = fdNdEtaCorrection->GetTrack2ParticleCorrection(vtx[2], eta, pt);
194
195     Float_t weight = track2particleCorr * vertexRecoCorr * triggerCorr;
196     if (weight <= 0)
197     {
198       AliDebug(AliLog::kError, Form("INFO: Skipping track because weight is <= 0 (track %d, weight %f) (vtx %f, eta %f, pt %f)", t, weight, vtx[2], eta, pt));
199       continue;
200     }
201
202     fdNdEtaAnalysisMBVtx->FillTrack(vtx[2], eta, pt, track2particleCorr);
203     fdNdEtaAnalysisMB->FillTrack(vtx[2], eta, pt, track2particleCorr * vertexRecoCorr);
204     fdNdEtaAnalysis->FillTrack(vtx[2], eta, pt, weight);
205   } // end of track loop
206
207   delete list;
208   list = 0;
209
210   // for event count per vertex
211   fdNdEtaAnalysisMBVtx->FillEvent(vtx[2], 1);
212   fdNdEtaAnalysisMB->FillEvent(vtx[2], vertexRecoCorr);
213   fdNdEtaAnalysis->FillEvent(vtx[2], vertexRecoCorr * triggerCorr);
214
215   return kTRUE;
216 }
217
218 void AlidNdEtaAnalysisESDSelector::SlaveTerminate()
219 {
220   // The SlaveTerminate() function is called after all entries or objects
221   // have been processed. When running with PROOF SlaveTerminate() is called
222   // on each slave server.
223
224   AliSelector::SlaveTerminate();
225
226   // Add the histograms to the output on each slave server
227   if (!fOutput)
228   {
229     AliDebug(AliLog::kError, Form("ERROR: Output list not initialized."));
230     return;
231   }
232
233   fOutput->Add(fdNdEtaAnalysis);
234   fOutput->Add(fdNdEtaAnalysisMB);
235   fOutput->Add(fdNdEtaAnalysisMBVtx);
236 }
237
238 void AlidNdEtaAnalysisESDSelector::Terminate()
239 {
240   // The Terminate() function is the last function to be called during
241   // a query. It always runs on the client, it can be used to present
242   // the results graphically or save the results to file.
243
244   AliSelector::Terminate();
245
246   fdNdEtaAnalysis = dynamic_cast<dNdEtaAnalysis*> (fOutput->FindObject("dndeta"));
247   fdNdEtaAnalysisMB = dynamic_cast<dNdEtaAnalysis*> (fOutput->FindObject("dndeta_mb"));
248   fdNdEtaAnalysisMBVtx = dynamic_cast<dNdEtaAnalysis*> (fOutput->FindObject("dndeta_mbvtx"));
249
250   if (!fdNdEtaAnalysis || !fdNdEtaAnalysisMB || !fdNdEtaAnalysisMBVtx)
251   {
252     AliDebug(AliLog::kError, Form("ERROR: Histograms not available %p %p %p", (void*) fdNdEtaAnalysis, (void*) fdNdEtaAnalysisMB, (void*) fdNdEtaAnalysisMBVtx));
253     return;
254   }
255
256   if (fdNdEtaAnalysis)
257     fdNdEtaAnalysis->Finish(fdNdEtaCorrection, 0.3);
258
259   if (fdNdEtaAnalysisMB)
260     fdNdEtaAnalysisMB->Finish(fdNdEtaCorrection, 0.3);
261
262   if (fdNdEtaAnalysisMBVtx)
263     fdNdEtaAnalysisMBVtx->Finish(fdNdEtaCorrection, 0.3);
264
265   TFile* fout = new TFile("analysis_esd.root","RECREATE");
266
267   if (fdNdEtaAnalysis)
268     fdNdEtaAnalysis->SaveHistograms();
269
270   if (fdNdEtaAnalysisMB)
271     fdNdEtaAnalysisMB->SaveHistograms();
272
273   if (fdNdEtaAnalysisMBVtx)
274     fdNdEtaAnalysisMBVtx->SaveHistograms();
275
276   if (fEsdTrackCuts)
277     fEsdTrackCuts->SaveHistograms("esd_tracks_cuts");
278
279   fout->Write();
280   fout->Close();
281 }