]>
Commit | Line | Data |
---|---|---|
10ebe68d | 1 | /* $Id$ */ |
2 | ||
3 | #include "AlidNdEtaSystematicsSelector.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 | #include <TH3F.h> | |
12 | #include <TH1F.h> | |
13 | #include <TParticle.h> | |
14 | ||
15 | #include <AliLog.h> | |
16 | #include <AliESD.h> | |
17 | #include <AliGenEventHeader.h> | |
18 | #include <AliStack.h> | |
19 | #include <AliHeader.h> | |
20 | ||
21 | #include "esdTrackCuts/AliESDtrackCuts.h" | |
22 | #include "AliPWG0Helper.h" | |
23 | #include "AlidNdEtaCorrection.h" | |
24 | ||
25 | ClassImp(AlidNdEtaSystematicsSelector) | |
26 | ||
27 | AlidNdEtaSystematicsSelector::AlidNdEtaSystematicsSelector() : | |
28 | AliSelectorRL(), | |
29 | fSecondaries(0), | |
30 | fEsdTrackCuts(0) | |
31 | { | |
32 | // | |
33 | // Constructor. Initialization of pointers | |
34 | // | |
35 | ||
36 | for (Int_t i=0; i<4; ++i) | |
37 | fdNdEtaCorrection[i] = 0; | |
38 | } | |
39 | ||
40 | AlidNdEtaSystematicsSelector::~AlidNdEtaSystematicsSelector() | |
41 | { | |
42 | // | |
43 | // Destructor | |
44 | // | |
45 | } | |
46 | ||
47 | void AlidNdEtaSystematicsSelector::Begin(TTree* tree) | |
48 | { | |
49 | // Begin function | |
50 | ||
51 | AliSelectorRL::Begin(tree); | |
52 | ||
53 | ReadUserObjects(tree); | |
54 | } | |
55 | ||
56 | void AlidNdEtaSystematicsSelector::ReadUserObjects(TTree* tree) | |
57 | { | |
58 | // read the user objects, called from slavebegin and begin | |
59 | ||
60 | if (!fEsdTrackCuts && fInput) | |
61 | fEsdTrackCuts = dynamic_cast<AliESDtrackCuts*> (fInput->FindObject("AliESDtrackCuts")); | |
62 | ||
63 | if (!fEsdTrackCuts && tree) | |
64 | fEsdTrackCuts = dynamic_cast<AliESDtrackCuts*> (tree->GetUserInfo()->FindObject("AliESDtrackCuts")); | |
65 | ||
66 | if (!fEsdTrackCuts) | |
67 | AliDebug(AliLog::kError, "ERROR: Could not read EsdTrackCuts from input list."); | |
68 | } | |
69 | ||
70 | void AlidNdEtaSystematicsSelector::SlaveBegin(TTree* tree) | |
71 | { | |
72 | // The SlaveBegin() function is called after the Begin() function. | |
73 | // When running with PROOF SlaveBegin() is called on each slave server. | |
74 | // The tree argument is deprecated (on PROOF 0 is passed). | |
75 | ||
76 | AliSelector::SlaveBegin(tree); | |
77 | ||
78 | ReadUserObjects(tree); | |
79 | ||
80 | TString option(GetOption()); | |
81 | ||
82 | printf("Running AlidNdEtaSystematicsSelector with options %s\n", option.Data()); | |
83 | ||
84 | if (option.Contains("secondaries")) | |
85 | { | |
86 | fSecondaries = new TH3F("fSecondaries", "fSecondaries;NacceptedTracks;CratioSecondaries;p_{T}", 201, -0.5, 205.5, 16, 0.45, 2.05, 10, 0, 10); | |
87 | } | |
88 | ||
89 | if (option.Contains("particle-composition")) | |
90 | { | |
91 | for (Int_t i=0; i<4; ++i) | |
92 | { | |
93 | TString name; | |
94 | name.Form("correction_%d", i); | |
95 | fdNdEtaCorrection[i] = new AlidNdEtaCorrection(name, name); | |
96 | } | |
97 | } | |
98 | } | |
99 | ||
100 | Bool_t AlidNdEtaSystematicsSelector::Process(Long64_t entry) | |
101 | { | |
102 | // The Process() function is called for each entry in the tree (or possibly | |
103 | // keyed object in the case of PROOF) to be processed. The entry argument | |
104 | // specifies which entry in the currently loaded tree is to be processed. | |
105 | // It can be passed to either TTree::GetEntry() or TBranch::GetEntry() | |
106 | // to read either all or the required parts of the data. When processing | |
107 | // keyed objects with PROOF, the object is already loaded and is available | |
108 | // via the fObject pointer. | |
109 | // | |
110 | // This function should contain the "body" of the analysis. It can contain | |
111 | // simple or elaborate selection criteria, run algorithms on the data | |
112 | // of the event and typically fill histograms. | |
113 | ||
114 | // WARNING when a selector is used with a TChain, you must use | |
115 | // the pointer to the current TTree to call GetEntry(entry). | |
116 | // The entry is always the local entry number in the current tree. | |
117 | // Assuming that fTree is the pointer to the TChain being processed, | |
118 | // use fTree->GetTree()->GetEntry(entry). | |
119 | ||
120 | if (AliSelectorRL::Process(entry) == kFALSE) | |
121 | return kFALSE; | |
122 | ||
123 | // Check prerequisites | |
124 | if (!fESD) | |
125 | { | |
126 | AliDebug(AliLog::kError, "ESD branch not available"); | |
127 | return kFALSE; | |
128 | } | |
129 | ||
130 | if (!fEsdTrackCuts) | |
131 | { | |
132 | AliDebug(AliLog::kError, "fESDTrackCuts not available"); | |
133 | return kFALSE; | |
134 | } | |
135 | ||
136 | AliStack* stack = GetStack(); | |
137 | if (!stack) | |
138 | { | |
139 | AliDebug(AliLog::kError, "Stack not available"); | |
140 | return kFALSE; | |
141 | } | |
142 | ||
143 | TObjArray* list = fEsdTrackCuts->GetAcceptedTracks(fESD); | |
144 | ||
145 | if (fdNdEtaCorrection[0]) | |
146 | FillCorrectionMaps(list); | |
147 | ||
148 | if (fSecondaries) | |
149 | FillSecondaries(list); | |
150 | ||
151 | delete list; | |
152 | list = 0; | |
153 | ||
154 | return kTRUE; | |
155 | } | |
156 | ||
157 | void AlidNdEtaSystematicsSelector::FillCorrectionMaps(TObjArray* listOfTracks) | |
158 | { | |
159 | // fills the correction maps for different particle species | |
160 | ||
161 | AliStack* stack = GetStack(); | |
162 | AliHeader* header = GetHeader(); | |
163 | ||
164 | Bool_t vertexReconstructed = AliPWG0Helper::IsVertexReconstructed(fESD); | |
165 | Bool_t eventTriggered = AliPWG0Helper::IsEventTriggered(fESD); | |
166 | ||
167 | for (Int_t i=0; i<4; ++i) | |
168 | { | |
169 | fdNdEtaCorrection[i]->IncreaseEventCount(); | |
170 | if (eventTriggered) | |
171 | fdNdEtaCorrection[i]->IncreaseTriggeredEventCount(); | |
172 | } | |
173 | ||
174 | // get the MC vertex | |
175 | AliGenEventHeader* genHeader = header->GenEventHeader(); | |
176 | ||
177 | TArrayF vtxMC(3); | |
178 | genHeader->PrimaryVertex(vtxMC); | |
179 | ||
180 | // loop over mc particles | |
181 | Int_t nPrim = stack->GetNprimary(); | |
182 | ||
183 | for (Int_t iMc = 0; iMc < nPrim; ++iMc) | |
184 | { | |
185 | TParticle* particle = stack->Particle(iMc); | |
186 | ||
187 | if (!particle) | |
188 | { | |
189 | AliDebug(AliLog::kError, Form("UNEXPECTED: particle with label %d not found in stack (mc loop).", iMc)); | |
190 | continue; | |
191 | } | |
192 | ||
193 | if (AliPWG0Helper::IsPrimaryCharged(particle, nPrim) == kFALSE) | |
194 | continue; | |
195 | ||
196 | Float_t eta = particle->Eta(); | |
197 | Float_t pt = particle->Pt(); | |
198 | ||
199 | Int_t id = -1; | |
200 | switch (TMath::Abs(particle->GetPdgCode())) | |
201 | { | |
202 | case 211: id = 0; break; | |
203 | case 321: id = 1; break; | |
204 | case 2212: id = 2; break; | |
205 | default: id = 3; break; | |
206 | } | |
207 | ||
208 | if (vertexReconstructed) | |
209 | fdNdEtaCorrection[id]->FillParticle(vtxMC[2], eta, pt); | |
210 | ||
211 | fdNdEtaCorrection[id]->FillParticleAllEvents(eta, pt); | |
212 | if (eventTriggered) | |
213 | fdNdEtaCorrection[id]->FillParticleWhenEventTriggered(eta, pt); | |
214 | }// end of mc particle | |
215 | ||
216 | // loop over esd tracks | |
217 | TIterator* iter = listOfTracks->MakeIterator(); | |
218 | TObject* obj = 0; | |
219 | while ((obj = iter->Next()) != 0) | |
220 | { | |
221 | AliESDtrack* esdTrack = dynamic_cast<AliESDtrack*> (obj); | |
222 | if (!esdTrack) | |
223 | continue; | |
224 | ||
225 | // using the properties of the mc particle | |
226 | Int_t label = TMath::Abs(esdTrack->GetLabel()); | |
227 | if (label == 0) | |
228 | { | |
229 | AliDebug(AliLog::kWarning, Form("WARNING: cannot find corresponding mc part for track %d.", label)); | |
230 | continue; | |
231 | } | |
232 | ||
233 | TParticle* particle = stack->Particle(label); | |
234 | if (!particle) | |
235 | { | |
236 | AliDebug(AliLog::kError, Form("UNEXPECTED: particle with label %d not found in stack (track loop).", label)); | |
237 | continue; | |
238 | } | |
239 | ||
240 | Int_t id = -1; | |
241 | switch (TMath::Abs(particle->GetPdgCode())) | |
242 | { | |
243 | case 211: id = 0; break; | |
244 | case 321: id = 1; break; | |
245 | case 2212: id = 2; break; | |
246 | default: id = 3; break; | |
247 | } | |
248 | ||
249 | if (vertexReconstructed) | |
250 | fdNdEtaCorrection[id]->FillParticleWhenMeasuredTrack(vtxMC[2], particle->Eta(), particle->Pt()); | |
251 | } // end of track loop | |
252 | ||
253 | delete iter; | |
254 | iter = 0; | |
255 | } | |
256 | ||
257 | void AlidNdEtaSystematicsSelector::FillSecondaries(TObjArray* listOfTracks) | |
258 | { | |
259 | // fills the secondary histograms | |
260 | ||
261 | AliStack* stack = GetStack(); | |
262 | ||
263 | TH1* nPrimaries = new TH1F("secondaries_primaries", "secondaries_primaries", fSecondaries->GetZaxis()->GetNbins(), fSecondaries->GetZaxis()->GetXmin(), fSecondaries->GetZaxis()->GetXmax()); | |
264 | TH1* nSecondaries = new TH1F("secondaries_secondaries", "secondaries_secondaries", fSecondaries->GetZaxis()->GetNbins(), fSecondaries->GetZaxis()->GetXmin(), fSecondaries->GetZaxis()->GetXmax()); | |
265 | ||
266 | TIterator* iter = listOfTracks->MakeIterator(); | |
267 | TObject* obj = 0; | |
268 | while ((obj = iter->Next()) != 0) | |
269 | { | |
270 | AliESDtrack* esdTrack = dynamic_cast<AliESDtrack*> (obj); | |
271 | if (!esdTrack) | |
272 | continue; | |
273 | ||
274 | Int_t label = TMath::Abs(esdTrack->GetLabel()); | |
275 | if (label == 0) | |
276 | { | |
277 | AliDebug(AliLog::kWarning, Form("WARNING: cannot find corresponding mc part for track %d.", label)); | |
278 | continue; | |
279 | } | |
280 | ||
281 | TParticle* particle = stack->Particle(label); | |
282 | if (!particle) | |
283 | { | |
284 | AliDebug(AliLog::kError, Form("UNEXPECTED: particle with label %d not found in stack (track loop).", label)); | |
285 | continue; | |
286 | } | |
287 | ||
288 | Int_t nPrim = stack->GetNprimary(); | |
289 | if (label < nPrim) | |
290 | nPrimaries->Fill(particle->Pt()); | |
291 | else | |
292 | nSecondaries->Fill(particle->Pt()); | |
293 | } | |
294 | ||
295 | for (Int_t i=1; i<=nPrimaries->GetNbinsX(); ++i) | |
296 | { | |
297 | Int_t primaries = (Int_t) nPrimaries->GetBinContent(i); | |
298 | Int_t secondaries = (Int_t) nSecondaries->GetBinContent(i); | |
299 | ||
300 | if (primaries + secondaries > 0) | |
301 | { | |
302 | AliDebug(AliLog::kDebug, Form("The ratio between primaries and secondaries is %d:%d = %f", primaries, secondaries, ((secondaries > 0) ? (Double_t) primaries / secondaries : -1))); | |
303 | ||
304 | for (Double_t factor = 0.5; factor < 2.01; factor += 0.1) | |
305 | fSecondaries->Fill((Double_t) primaries + (Double_t) secondaries * factor, factor, nPrimaries->GetBinCenter(i)); | |
306 | } | |
307 | } | |
308 | ||
309 | delete nPrimaries; | |
310 | nPrimaries = 0; | |
311 | ||
312 | delete nSecondaries; | |
313 | nSecondaries = 0; | |
314 | ||
315 | delete iter; | |
316 | iter = 0; | |
317 | } | |
318 | ||
319 | void AlidNdEtaSystematicsSelector::SlaveTerminate() | |
320 | { | |
321 | // The SlaveTerminate() function is called after all entries or objects | |
322 | // have been processed. When running with PROOF SlaveTerminate() is called | |
323 | // on each slave server. | |
324 | ||
325 | AliSelector::SlaveTerminate(); | |
326 | ||
327 | // Add the histograms to the output on each slave server | |
328 | if (!fOutput) | |
329 | { | |
330 | AliDebug(AliLog::kError, Form("ERROR: Output list not initialized.")); | |
331 | return; | |
332 | } | |
333 | ||
334 | if (fSecondaries) | |
335 | fOutput->Add(fSecondaries); | |
336 | ||
337 | for (Int_t i=0; i<4; ++i) | |
338 | if (fdNdEtaCorrection[i]) | |
339 | fOutput->Add(fdNdEtaCorrection[i]); | |
340 | } | |
341 | ||
342 | void AlidNdEtaSystematicsSelector::Terminate() | |
343 | { | |
344 | // The Terminate() function is the last function to be called during | |
345 | // a query. It always runs on the client, it can be used to present | |
346 | // the results graphically or save the results to file. | |
347 | ||
348 | AliSelector::Terminate(); | |
349 | ||
350 | fSecondaries = dynamic_cast<TH3F*> (fOutput->FindObject("fSecondaries")); | |
351 | for (Int_t i=0; i<4; ++i) | |
352 | fdNdEtaCorrection[i] = dynamic_cast<AlidNdEtaCorrection*> (fOutput->FindObject(Form("correction_%d", i))); | |
353 | ||
354 | TFile* fout = TFile::Open("systematics.root", "RECREATE"); | |
355 | ||
356 | if (fEsdTrackCuts) | |
357 | fEsdTrackCuts->SaveHistograms("esd_track_cuts"); | |
358 | ||
359 | if (fSecondaries) | |
360 | fSecondaries->Write(); | |
361 | ||
362 | for (Int_t i=0; i<4; ++i) | |
363 | if (fdNdEtaCorrection[i]) | |
364 | fdNdEtaCorrection[i]->SaveHistograms(); | |
365 | ||
366 | fout->Write(); | |
367 | fout->Close(); | |
368 | ||
369 | if (fSecondaries) | |
370 | { | |
371 | new TCanvas; | |
372 | fSecondaries->Draw(); | |
373 | } | |
374 | } |