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