c068884b28a6977ed0dd28b97188ccd96b50f0a6
[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   if (fdNdEtaAnalysis)
48   {
49     delete fdNdEtaAnalysis;
50     fdNdEtaAnalysis = 0;
51   }
52
53   if (fdNdEtaAnalysisMB)
54   {
55     delete fdNdEtaAnalysisMB;
56     fdNdEtaAnalysisMB = 0;
57   }
58
59   if (fdNdEtaAnalysisMBVtx)
60   {
61     delete fdNdEtaAnalysisMBVtx;
62     fdNdEtaAnalysisMBVtx = 0;
63   }
64
65   /*if (fEsdTrackCuts)
66   {
67     delete fEsdTrackCuts;
68     fEsdTrackCuts = 0;
69   }
70
71   if (fdNdEtaCorrection)
72   {
73     delete fdNdEtaCorrection;
74     fdNdEtaCorrection = 0;
75   }*/
76 }
77
78 void AlidNdEtaAnalysisESDSelector::Begin(TTree* tree)
79 {
80   // Begin function
81
82   ReadUserObjects(tree);
83 }
84
85 void AlidNdEtaAnalysisESDSelector::ReadUserObjects(TTree* tree)
86 {
87   // read the user objects, called from slavebegin and begin
88
89   if (!fEsdTrackCuts && fInput)
90     fEsdTrackCuts = dynamic_cast<AliESDtrackCuts*> (fInput->FindObject("AliESDtrackCuts"));
91
92   if (!fEsdTrackCuts && tree)
93     fEsdTrackCuts = dynamic_cast<AliESDtrackCuts*> (tree->GetUserInfo()->FindObject("AliESDtrackCuts"));
94
95   if (!fEsdTrackCuts)
96      AliDebug(AliLog::kError, "ERROR: Could not read EsdTrackCuts from input list.");
97
98
99   if (!fdNdEtaCorrection && fInput)
100     fdNdEtaCorrection = dynamic_cast<AlidNdEtaCorrection*> (fInput->FindObject("dndeta_correction"));
101
102   if (!fdNdEtaCorrection && tree)
103     fdNdEtaCorrection = dynamic_cast<AlidNdEtaCorrection*> (tree->GetUserInfo()->FindObject("dndeta_correction"));
104
105   if (!fdNdEtaCorrection)
106      AliDebug(AliLog::kError, "ERROR: Could not read dndeta_correction from input list.");
107 }
108
109 void AlidNdEtaAnalysisESDSelector::SlaveBegin(TTree* tree)
110 {
111   // The SlaveBegin() function is called after the Begin() function.
112   // When running with PROOF SlaveBegin() is called on each slave server.
113   // The tree argument is deprecated (on PROOF 0 is passed).
114
115   AliSelector::SlaveBegin(tree);
116
117   ReadUserObjects(tree);
118
119   fdNdEtaAnalysisMBVtx = new dNdEtaAnalysis("dndeta_mbvtx", "dndeta_mbvtx");
120   fdNdEtaAnalysisMB = new dNdEtaAnalysis("dndeta_mb", "dndeta_mb");
121   fdNdEtaAnalysis = new dNdEtaAnalysis("dndeta", "dndeta");
122 }
123
124 void AlidNdEtaAnalysisESDSelector::Init(TTree* tree)
125 {
126   // read the user objects
127
128   AliSelector::Init(tree);
129 }
130
131 Bool_t AlidNdEtaAnalysisESDSelector::Process(Long64_t entry)
132 {
133   // The Process() function is called for each entry in the tree (or possibly
134   // keyed object in the case of PROOF) to be processed. The entry argument
135   // specifies which entry in the currently loaded tree is to be processed.
136   // It can be passed to either TTree::GetEntry() or TBranch::GetEntry()
137   // to read either all or the required parts of the data. When processing
138   // keyed objects with PROOF, the object is already loaded and is available
139   // via the fObject pointer.
140   //
141   // This function should contain the "body" of the analysis. It can contain
142   // simple or elaborate selection criteria, run algorithms on the data
143   // of the event and typically fill histograms.
144
145   // WARNING when a selector is used with a TChain, you must use
146   //  the pointer to the current TTree to call GetEntry(entry).
147   //  The entry is always the local entry number in the current tree.
148   //  Assuming that fTree is the pointer to the TChain being processed,
149   //  use fTree->GetTree()->GetEntry(entry).
150
151   if (AliSelector::Process(entry) == kFALSE)
152     return kFALSE;
153
154   // Check prerequisites
155   if (!fESD)
156   {
157     AliDebug(AliLog::kError, "ESD branch not available");
158     return kFALSE;
159   }
160
161   if (!fEsdTrackCuts)
162   {
163     AliDebug(AliLog::kError, "fESDTrackCuts not available");
164     return kFALSE;
165   }
166
167   if (!fdNdEtaCorrection)
168   {
169     AliDebug(AliLog::kError, "fdNdEtaCorrection not available");
170     return kFALSE;
171   }
172
173   if (AliPWG0Helper::IsEventTriggered(fESD) == kFALSE)
174     return kTRUE;
175
176   if (AliPWG0Helper::IsVertexReconstructed(fESD) == kFALSE)
177     return kTRUE;
178
179   // ########################################################
180   // get the EDS vertex
181   const AliESDVertex* vtxESD = fESD->GetVertex();
182   Double_t vtx[3];
183   vtxESD->GetXYZ(vtx);
184
185   // get number of "good" tracks
186   TObjArray* list = fEsdTrackCuts->GetAcceptedTracks(fESD);
187   Int_t nGoodTracks = list->GetEntries();
188
189   Float_t vertexRecoCorr = fdNdEtaCorrection->GetVertexRecoCorrection(vtx[2], nGoodTracks);
190   if (vertexRecoCorr <= 0)
191   {
192     AliDebug(AliLog::kError, Form("INFO: Skipping event because vertexRecoCorr is <= 0 (%f)", vertexRecoCorr));
193     delete list;
194     return kTRUE;
195   }
196
197   Float_t triggerCorr = fdNdEtaCorrection->GetTriggerCorrection(vtx[2], nGoodTracks);
198   if (triggerCorr <= 0)
199   {
200     AliDebug(AliLog::kError, Form("INFO: Skipping event because triggerCorr is <= 0 (%f)", triggerCorr));
201     delete list;
202     return kTRUE;
203   }
204
205   // loop over esd tracks
206   for (Int_t t=0; t<nGoodTracks; t++)
207   {
208     AliESDtrack* esdTrack = dynamic_cast<AliESDtrack*> (list->At(t));
209     if (!esdTrack)
210     {
211       AliDebug(AliLog::kError, Form("ERROR: Could not retrieve track %d.", t));
212       continue;
213     }
214
215     Double_t p[3];
216     esdTrack->GetConstrainedPxPyPz(p); // ### TODO should be okay because we have a vertex, however GetInnerPxPyPy / GetOuterPxPyPy also exist
217     TVector3 vector(p);
218
219     Float_t theta = vector.Theta();
220     Float_t eta   = -TMath::Log(TMath::Tan(theta/2.));
221     Float_t pt = vector.Pt();
222
223     Float_t track2particleCorr = fdNdEtaCorrection->GetTrack2ParticleCorrection(vtx[2], eta, pt);
224
225     Float_t weight = track2particleCorr * vertexRecoCorr * triggerCorr;
226     if (weight <= 0)
227     {
228       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));
229       continue;
230     }
231
232     fdNdEtaAnalysisMBVtx->FillTrack(vtx[2], eta, pt, track2particleCorr);
233     fdNdEtaAnalysisMB->FillTrack(vtx[2], eta, pt, track2particleCorr * vertexRecoCorr);
234     fdNdEtaAnalysis->FillTrack(vtx[2], eta, pt, weight);
235   } // end of track loop
236
237   delete list;
238   list = 0;
239
240   // for event count per vertex
241   fdNdEtaAnalysisMBVtx->FillEvent(vtx[2], 1);
242   fdNdEtaAnalysisMB->FillEvent(vtx[2], vertexRecoCorr);
243   fdNdEtaAnalysis->FillEvent(vtx[2], vertexRecoCorr * triggerCorr);
244
245   return kTRUE;
246 }
247
248 void AlidNdEtaAnalysisESDSelector::SlaveTerminate()
249 {
250   // The SlaveTerminate() function is called after all entries or objects
251   // have been processed. When running with PROOF SlaveTerminate() is called
252   // on each slave server.
253
254   AliSelector::SlaveTerminate();
255
256   // Add the histograms to the output on each slave server
257   if (!fOutput)
258   {
259     AliDebug(AliLog::kError, Form("ERROR: Output list not initialized."));
260     return;
261   }
262
263   // Add the objects to the output list and set them to 0, so that the destructor does not delete them.
264
265   fOutput->Add(fdNdEtaAnalysis);
266   fdNdEtaAnalysis = 0;
267
268   fOutput->Add(fdNdEtaAnalysisMB);
269   fdNdEtaAnalysisMB = 0;
270
271   fOutput->Add(fdNdEtaAnalysisMBVtx);
272   fdNdEtaAnalysisMBVtx = 0;
273 }
274
275 void AlidNdEtaAnalysisESDSelector::Terminate()
276 {
277   // The Terminate() function is the last function to be called during
278   // a query. It always runs on the client, it can be used to present
279   // the results graphically or save the results to file.
280
281   AliSelector::Terminate();
282
283   fdNdEtaAnalysis = dynamic_cast<dNdEtaAnalysis*> (fOutput->FindObject("dndeta"));
284   fdNdEtaAnalysisMB = dynamic_cast<dNdEtaAnalysis*> (fOutput->FindObject("dndeta_mb"));
285   fdNdEtaAnalysisMBVtx = dynamic_cast<dNdEtaAnalysis*> (fOutput->FindObject("dndeta_mbvtx"));
286
287   if (!fdNdEtaAnalysis || !fdNdEtaAnalysisMB || !fdNdEtaAnalysisMBVtx)
288   {
289     AliDebug(AliLog::kError, Form("ERROR: Histograms not available %p %p %p", (void*) fdNdEtaAnalysis, (void*) fdNdEtaAnalysisMB, (void*) fdNdEtaAnalysisMBVtx));
290     return;
291   }
292
293   if (fdNdEtaAnalysis)
294     fdNdEtaAnalysis->Finish(fdNdEtaCorrection, 0.3);
295
296   if (fdNdEtaAnalysisMB)
297     fdNdEtaAnalysisMB->Finish(fdNdEtaCorrection, 0.3);
298
299   if (fdNdEtaAnalysisMBVtx)
300     fdNdEtaAnalysisMBVtx->Finish(fdNdEtaCorrection, 0.3);
301
302   TFile* fout = new TFile("analysis_esd.root","RECREATE");
303
304   if (fdNdEtaAnalysis)
305     fdNdEtaAnalysis->SaveHistograms();
306
307   if (fdNdEtaAnalysisMB)
308     fdNdEtaAnalysisMB->SaveHistograms();
309
310   if (fdNdEtaAnalysisMBVtx)
311     fdNdEtaAnalysisMBVtx->SaveHistograms();
312
313   if (fEsdTrackCuts)
314     fEsdTrackCuts->SaveHistograms("esd_tracks_cuts");
315
316   if (fdNdEtaCorrection)
317     fdNdEtaCorrection->SaveHistograms();
318
319   fout->Write();
320   fout->Close();
321 }