3 #include "AlidNdEtaAnalysisESDSelector.h"
13 #include <AliESDVertex.h>
16 #include "esdTrackCuts/AliESDtrackCuts.h"
17 #include "dNdEta/dNdEtaAnalysis.h"
18 #include "dNdEta/AlidNdEtaCorrection.h"
19 #include "AliPWG0Helper.h"
21 ClassImp(AlidNdEtaAnalysisESDSelector)
23 AlidNdEtaAnalysisESDSelector::AlidNdEtaAnalysisESDSelector() :
25 fdNdEtaAnalysisMBVtx(0),
32 // Constructor. Initialization of pointers
35 AliLog::SetClassDebugLevel("AlidNdEtaAnalysisESDSelector", AliLog::kDebug);
38 AlidNdEtaAnalysisESDSelector::~AlidNdEtaAnalysisESDSelector()
44 // histograms are in the output list and deleted when the output
45 // list is deleted by the TSelector dtor
48 void AlidNdEtaAnalysisESDSelector::Begin(TTree* tree)
52 ReadUserObjects(tree);
55 void AlidNdEtaAnalysisESDSelector::ReadUserObjects(TTree* tree)
57 // read the user objects, called from slavebegin and begin
59 if (!fEsdTrackCuts && fInput)
60 fEsdTrackCuts = dynamic_cast<AliESDtrackCuts*> (fInput->FindObject("AliESDtrackCuts"));
62 if (!fEsdTrackCuts && tree)
63 fEsdTrackCuts = dynamic_cast<AliESDtrackCuts*> (tree->GetUserInfo()->FindObject("AliESDtrackCuts"));
66 AliDebug(AliLog::kError, "ERROR: Could not read EsdTrackCuts from input list.");
69 if (!fdNdEtaCorrection && fInput)
70 fdNdEtaCorrection = dynamic_cast<AlidNdEtaCorrection*> (fInput->FindObject("dndeta_correction"));
72 if (!fdNdEtaCorrection && tree)
73 fdNdEtaCorrection = dynamic_cast<AlidNdEtaCorrection*> (tree->GetUserInfo()->FindObject("dndeta_correction"));
75 if (!fdNdEtaCorrection)
76 AliDebug(AliLog::kError, "ERROR: Could not read dndeta_correction from input list.");
79 void AlidNdEtaAnalysisESDSelector::SlaveBegin(TTree* tree)
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).
85 AliSelector::SlaveBegin(tree);
87 ReadUserObjects(tree);
89 fdNdEtaAnalysisMBVtx = new dNdEtaAnalysis("dndeta_mbvtx", "dndeta_mbvtx");
90 fdNdEtaAnalysisMB = new dNdEtaAnalysis("dndeta_mb", "dndeta_mb");
91 fdNdEtaAnalysis = new dNdEtaAnalysis("dndeta", "dndeta");
94 void AlidNdEtaAnalysisESDSelector::Init(TTree* tree)
96 // read the user objects
98 AliSelector::Init(tree);
100 // Enable only the needed branches
103 tree->SetBranchStatus("*", 0);
104 tree->SetBranchStatus("fTriggerMask", 1);
105 tree->SetBranchStatus("fSPDVertex*", 1);
107 AliESDtrackCuts::EnableNeededBranches(tree);
111 Bool_t AlidNdEtaAnalysisESDSelector::Process(Long64_t entry)
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.
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.
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.
128 // Assuming that fTree is the pointer to the TChain being processed,
129 // use fTree->GetTree()->GetEntry(entry).
131 if (AliSelector::Process(entry) == kFALSE)
134 // Check prerequisites
137 AliDebug(AliLog::kError, "ESD branch not available");
143 AliDebug(AliLog::kError, "fESDTrackCuts not available");
147 if (!fdNdEtaCorrection)
149 AliDebug(AliLog::kError, "fdNdEtaCorrection not available");
153 if (AliPWG0Helper::IsEventTriggered(fESD) == kFALSE)
155 AliDebug(AliLog::kDebug+1, Form("Skipping event %d because it was not triggered", (Int_t) entry));
159 if (AliPWG0Helper::IsVertexReconstructed(fESD) == kFALSE)
161 AliDebug(AliLog::kDebug+1, Form("Skipping event %d because its vertex was not reconstructed", (Int_t) entry));
165 // ########################################################
166 // get the EDS vertex
167 const AliESDVertex* vtxESD = fESD->GetVertex();
171 // get number of "good" tracks
172 TObjArray* list = fEsdTrackCuts->GetAcceptedTracks(fESD);
173 Int_t nGoodTracks = list->GetEntries();
175 Float_t vertexRecoCorr = fdNdEtaCorrection->GetVertexRecoCorrection(vtx[2], nGoodTracks);
176 if (vertexRecoCorr <= 0)
178 AliDebug(AliLog::kError, Form("INFO: Skipping event because vertexRecoCorr is <= 0 (%f)", vertexRecoCorr));
183 Float_t triggerCorr = fdNdEtaCorrection->GetTriggerBiasCorrection(vtx[2], nGoodTracks);
184 if (triggerCorr <= 0)
186 AliDebug(AliLog::kError, Form("INFO: Skipping event because triggerCorr is <= 0 (%f)", triggerCorr));
191 // loop over esd tracks
192 for (Int_t t=0; t<nGoodTracks; t++)
194 AliESDtrack* esdTrack = dynamic_cast<AliESDtrack*> (list->At(t));
197 AliDebug(AliLog::kError, Form("ERROR: Could not retrieve track %d.", t));
202 esdTrack->GetConstrainedPxPyPz(p); // ### TODO should be okay because we have a vertex, however GetInnerPxPyPy / GetOuterPxPyPy also exist
205 Float_t theta = vector.Theta();
206 Float_t eta = -TMath::Log(TMath::Tan(theta/2.));
207 Float_t pt = vector.Pt();
209 Float_t track2particleCorr = fdNdEtaCorrection->GetTrack2ParticleCorrection(vtx[2], eta, pt);
211 Float_t weight = track2particleCorr * vertexRecoCorr * triggerCorr;
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));
218 fdNdEtaAnalysisMBVtx->FillTrack(vtx[2], eta, pt, track2particleCorr);
219 fdNdEtaAnalysisMB->FillTrack(vtx[2], eta, pt, track2particleCorr * vertexRecoCorr);
220 fdNdEtaAnalysis->FillTrack(vtx[2], eta, pt, weight);
221 } // end of track loop
226 // for event count per vertex
227 fdNdEtaAnalysisMBVtx->FillEvent(vtx[2], 1);
228 fdNdEtaAnalysisMB->FillEvent(vtx[2], vertexRecoCorr);
229 fdNdEtaAnalysis->FillEvent(vtx[2], vertexRecoCorr * triggerCorr);
234 void AlidNdEtaAnalysisESDSelector::SlaveTerminate()
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.
240 AliSelector::SlaveTerminate();
242 // Add the histograms to the output on each slave server
245 AliDebug(AliLog::kError, Form("ERROR: Output list not initialized."));
249 // Add the objects to the output list and set them to 0, so that the destructor does not delete them.
251 fOutput->Add(fdNdEtaAnalysis);
254 fOutput->Add(fdNdEtaAnalysisMB);
255 fdNdEtaAnalysisMB = 0;
257 fOutput->Add(fdNdEtaAnalysisMBVtx);
258 fdNdEtaAnalysisMBVtx = 0;
261 void AlidNdEtaAnalysisESDSelector::Terminate()
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.
267 AliSelector::Terminate();
269 fdNdEtaAnalysis = dynamic_cast<dNdEtaAnalysis*> (fOutput->FindObject("dndeta"));
270 fdNdEtaAnalysisMB = dynamic_cast<dNdEtaAnalysis*> (fOutput->FindObject("dndeta_mb"));
271 fdNdEtaAnalysisMBVtx = dynamic_cast<dNdEtaAnalysis*> (fOutput->FindObject("dndeta_mbvtx"));
273 if (!fdNdEtaAnalysis || !fdNdEtaAnalysisMB || !fdNdEtaAnalysisMBVtx)
275 AliDebug(AliLog::kError, Form("ERROR: Histograms not available %p %p %p", (void*) fdNdEtaAnalysis, (void*) fdNdEtaAnalysisMB, (void*) fdNdEtaAnalysisMBVtx));
280 fdNdEtaAnalysis->Finish(fdNdEtaCorrection, 0.3);
282 if (fdNdEtaAnalysisMB)
283 fdNdEtaAnalysisMB->Finish(fdNdEtaCorrection, 0.3);
285 if (fdNdEtaAnalysisMBVtx)
286 fdNdEtaAnalysisMBVtx->Finish(fdNdEtaCorrection, 0.3);
288 TFile* fout = new TFile("analysis_esd.root","RECREATE");
291 fdNdEtaAnalysis->SaveHistograms();
293 if (fdNdEtaAnalysisMB)
294 fdNdEtaAnalysisMB->SaveHistograms();
296 if (fdNdEtaAnalysisMBVtx)
297 fdNdEtaAnalysisMBVtx->SaveHistograms();
300 fEsdTrackCuts->SaveHistograms("esd_tracks_cuts");
302 if (fdNdEtaCorrection)
303 fdNdEtaCorrection->SaveHistograms();