First update for AOD compatibility (distributions and efficiency), still missing...
[u/mrichter/AliRoot.git] / PWGCF / EBYE / NetParticle / AliAnalysisNetParticleEffCont.cxx
CommitLineData
cb68eb1d 1//-*- Mode: C++ -*-
2
3#include "TMath.h"
4#include "TAxis.h"
9be43c8e 5#include "TDatabasePDG.h"
cb68eb1d 6
7#include "AliESDEvent.h"
8#include "AliESDInputHandler.h"
9#include "AliStack.h"
10#include "AliMCEvent.h"
cb68eb1d 11#include "AliESDtrackCuts.h"
9be43c8e 12#include "AliAODEvent.h"
13#include "AliAODInputHandler.h"
14#include "AliAODMCParticle.h"
cb68eb1d 15
16#include "AliAnalysisNetParticleEffCont.h"
17
18using namespace std;
19
20// Task for NetParticle checks
21// Author: Jochen Thaeder <jochen@thaeder.de>
22
23ClassImp(AliAnalysisNetParticleEffCont)
24
25/*
26 * ---------------------------------------------------------------------------------
27 * Constructor / Destructor
28 * ---------------------------------------------------------------------------------
29 */
30
31//________________________________________________________________________
32AliAnalysisNetParticleEffCont::AliAnalysisNetParticleEffCont() :
33 fHelper(NULL),
34 fPdgCode(2212),
35
36 fESD(NULL),
37 fESDTrackCuts(NULL),
38
9be43c8e 39 fAOD(NULL),
40 fArrayMC(NULL),
41
cb68eb1d 42 fCentralityBin(-1.),
43 fNTracks(0),
44
9be43c8e 45 fAODtrackCutBit(1024),
46
cb68eb1d 47 fStack(NULL),
48 fMCEvent(NULL),
49
50 fLabelsRec(NULL),
51
52 fHnEff(NULL),
53 fHnCont(NULL) {
54 // Constructor
55
56 AliLog::SetClassDebugLevel("AliAnalysisNetParticleEffCont",10);
57}
58
59//________________________________________________________________________
60AliAnalysisNetParticleEffCont::~AliAnalysisNetParticleEffCont() {
61 // Destructor
62
63 if (fLabelsRec[0])
64 delete[] (fLabelsRec[0]);
65 if (fLabelsRec[1])
66 delete[] (fLabelsRec[1]);
67 if (fLabelsRec)
68 delete[] fLabelsRec;
69}
70
71/*
72 * ---------------------------------------------------------------------------------
73 * Public Methods
74 * ---------------------------------------------------------------------------------
75 */
76
77//________________________________________________________________________
9be43c8e 78void AliAnalysisNetParticleEffCont::Initialize(AliESDtrackCuts *cuts , AliAnalysisNetParticleHelper* helper, Int_t trackCutBit) {
cb68eb1d 79
80 // -- ESD track cuts
81 // -------------------
82 fESDTrackCuts = cuts;
83
84 // -- Get Helper
85 // ---------------
86 fHelper = helper;
87
88 // -- Get particle species / pdgCode
89 // -------------------------
90 fPdgCode = AliPID::ParticleCode(fHelper->GetParticleSpecies());
91
92 // -- MC Labels for efficiency
93 // -----------------------------
94 fLabelsRec = new Int_t*[2];
95 for (Int_t ii = 0; ii < 2 ; ++ii)
96 fLabelsRec[ii] = NULL;
97
9be43c8e 98 // -- AOD track filter bit
99 // -------------------------
100 fAODtrackCutBit = trackCutBit;
101
102
cb68eb1d 103 // -- Create THnSparse Histograms
104 // --------------------------------
105 CreateHistograms();
106
107 return;
108}
109
110/*
111 * ---------------------------------------------------------------------------------
112 * Setup/Reset Methods - private
113 * ---------------------------------------------------------------------------------
114 */
115
116//________________________________________________________________________
117Int_t AliAnalysisNetParticleEffCont::SetupEvent(AliESDInputHandler *esdHandler, AliMCEvent *mcEvent) {
118 // -- Setup event
119
120 // -- Reset of event
121 // -------------------
122 ResetEvent();
123
124 // -- Setup of event
125 // -------------------
126 fESD = esdHandler->GetEvent();
127 fNTracks = fESD->GetNumberOfTracks();
128
129 fMCEvent = mcEvent;
130 fStack = fMCEvent->Stack();
131
132 fCentralityBin = fHelper->GetCentralityBin();
133
134 // -- Create label arrays
135 // ------------------------
136 fLabelsRec[0] = new Int_t[fNTracks];
137 if(!fLabelsRec[0]) {
138 AliError("Cannot create fLabelsRec[0]");
139 return -1;
140 }
141
142 fLabelsRec[1] = new Int_t[fNTracks];
143 if(!fLabelsRec[1]) {
144 AliError("Cannot create fLabelsRec[1] for TPC");
145 return -1;
146 }
147
148 for(Int_t ii = 0; ii < fNTracks; ++ii) {
149 fLabelsRec[0][ii] = 0;
150 fLabelsRec[1][ii] = 0;
151 }
152
153 return 0;
154}
155
9be43c8e 156//________________________________________________________________________
157Int_t AliAnalysisNetParticleEffCont::SetupEvent(AliAODInputHandler *aodHandler) {
158 // -- Setup event
159
160 // -- Reset of event
161 // -------------------
162 ResetEvent();
163
164 // -- Setup of event
165 // -------------------
166 fAOD = aodHandler->GetEvent();
167 fNTracks = fAOD->GetNumberOfTracks();
168
169 fArrayMC = dynamic_cast<TClonesArray*>(fAOD->FindListObject(AliAODMCParticle::StdBranchName()));
170 if (!fArrayMC)
171 AliFatal("No array of MC particles found !!!");
172
173 fCentralityBin = fHelper->GetCentralityBin();
174
175 // -- Create label arrays
176 // ------------------------
177 fLabelsRec[0] = new Int_t[fNTracks];
178 if(!fLabelsRec[0]) {
179 AliError("Cannot create fLabelsRec[0]");
180 return -1;
181 }
182
183 fLabelsRec[1] = new Int_t[fNTracks];
184 if(!fLabelsRec[1]) {
185 AliError("Cannot create fLabelsRec[1] for TPC");
186 return -1;
187 }
188
189 for(Int_t ii = 0; ii < fNTracks; ++ii) {
190 fLabelsRec[0][ii] = 0;
191 fLabelsRec[1][ii] = 0;
192 }
193
194 return 0;
195}
196
cb68eb1d 197//________________________________________________________________________
198void AliAnalysisNetParticleEffCont::ResetEvent() {
199 // -- Reset event
200
201 for (Int_t ii = 0; ii < 2 ; ++ii) {
202 if (fLabelsRec[ii])
203 delete[] fLabelsRec[ii];
204 fLabelsRec[ii] = NULL;
205 }
206
207 return;
208}
209
210//________________________________________________________________________
211void AliAnalysisNetParticleEffCont::Process() {
212 // -- Process event
213
214 // -- Setup (clean, create and fill) MC labels
215 // ---------------------------------------------
9be43c8e 216 if(fESD)
217 FillMCLabels();
218 else if(fAOD)
219 FillMCLabelsAOD();
220
cb68eb1d 221 // -- Fill MC histograms for efficiency studies
222 // -----------------------------------------------
9be43c8e 223 if(fESD)
224 FillMCEffHist();
225 else if(fAOD)
226 FillMCEffHistAOD();
cb68eb1d 227
228 return;
229}
230
231/*
232 * ---------------------------------------------------------------------------------
233 * Private Methods
234 * ---------------------------------------------------------------------------------
235 */
236
237//________________________________________________________________________
238void AliAnalysisNetParticleEffCont::CreateHistograms() {
239 // -- Create histograms
240
241 Double_t dCent[2] = {-0.5, 8.5};
242 Int_t iCent = 9;
243
244 Double_t dEta[2] = {-0.9, 0.9};
245 Int_t iEta = Int_t((dEta[1]-dEta[0]) / 0.075) +1 ;
246
247 Double_t dRap[2] = {-0.5, 0.5};
248 Int_t iRap = Int_t((dRap[1]-dRap[0]) / 0.075) +1 ;
249
250 Double_t dPhi[2] = {0.0, TMath::TwoPi()};
251 Int_t iPhi = 42;
252
253 Double_t dPt[2] = {0.1, 3.0};
254 Int_t iPt = Int_t((dPt[1]-dPt[0]) / 0.075);
255
256 Double_t dSign[2] = {-1.5, 1.5};
257 Int_t iSign = 3;
258
259 // ------------------------------------------------------------------
260 // -- Create THnSparseF - Eff
261 // ------------------------------------------------------------------
262
263 // cent: etaMC: yMC: phiMC: ptMC: sign:findable:recStatus:pidStatus: etaRec: phiRec: ptRec:deltaEta: deltaPhi: deltaPt
264 Int_t binHnEff[15] = { iCent, iEta, iRap, iPhi, iPt, iSign, 2, 2 , 2 , iEta, iPhi, iPt, iEta, 2*iPhi+1, 2*iPt+1};
265 Double_t minHnEff[15] = {dCent[0], dEta[0], dRap[0], dPhi[0], dPt[0], dSign[0], -0.5, -0.5, -0.5, dEta[0], dPhi[0], dPt[0], dEta[0], -dPhi[1], -dPt[1]};
266 Double_t maxHnEff[15] = {dCent[1], dEta[1], dRap[0], dPhi[1], dPt[1], dSign[1], 1.5, 1.5, 1.5, dEta[1], dPhi[1], dPt[1], dEta[1], dPhi[1], dPt[1]};
267
268 fHnEff = new THnSparseF("fHnEff", "cent:etaMC:yMC:phiMC:ptMC:sign:findable:recStatus:pidStatus:etaRec:phiRec:ptRec:deltaEta:deltaPhi:deltaPt",
269 15, binHnEff, minHnEff, maxHnEff);
270
271 fHnEff->Sumw2();
272
273 fHnEff->GetAxis(0)->SetTitle("centrality"); // 0-5|5-10|10-20|20-30|30-40|40-50|50-60|60-70|70-80|80-90 --> 10 bins
274 fHnEff->GetAxis(1)->SetTitle("#eta_{MC}"); // eta [-0.9, 0.9]
275 fHnEff->GetAxis(2)->SetTitle("#it{y}_{MC}"); // rapidity [-0.5, 0.5]
276 fHnEff->GetAxis(3)->SetTitle("#varphi_{MC} (rad)"); // phi [ 0. , 2Pi]
277 fHnEff->GetAxis(4)->SetTitle("#it{p}_{T,MC} (GeV/#it{c})"); // pt [ 0.1, 3.0]
278
279 fHnEff->GetAxis(5)->SetTitle("sign"); // -1 | 0 | +1
280 fHnEff->GetAxis(6)->SetTitle("findable"); // 0 not findable | 1 findable
281 fHnEff->GetAxis(7)->SetTitle("recStatus"); // 0 not reconstructed | 1 reconstructed
282 fHnEff->GetAxis(8)->SetTitle("recPid"); // 0 not accepted | 1 accepted
283
284 fHnEff->GetAxis(9)->SetTitle("#eta_{Rec}"); // eta [-0.9, 0.9]
285 fHnEff->GetAxis(10)->SetTitle("#varphi_{Rec} (rad)"); // phi [ 0. , 2Pi]
286 fHnEff->GetAxis(11)->SetTitle("#it{p}_{T,Rec} (GeV/#it{c})"); // pt [ 0.1, 3.0]
287
288 fHnEff->GetAxis(12)->SetTitle("#eta_{MC}-#eta_{Rec}"); // eta [-0.9,0.9]
289 fHnEff->GetAxis(13)->SetTitle("#varphi_{MC}-#varphi_{Rec} (rad)"); // phi [ 0. ,2Pi]
290 fHnEff->GetAxis(14)->SetTitle("#it{p}_{T,MC}-#it{p}_{T,Rec} (GeV/#it{c})"); // pt [ 3.0,3.0]
291
292 fHelper->BinLogAxis(fHnEff, 4);
293 fHelper->BinLogAxis(fHnEff, 11);
294
295 /*
296 >>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>
297
298 creation -> findable -> reconstructed -> pid_TPC+TOF
299 (1) (2) (3) (4)
300 || findable | recStatus | recPid
301 1) all primary probeParticles_MC || - - -
302 2) all findable primary probeParticles_MC || x - -
303 3) all reconstructed primary probeParticles_MC || x x -
304 4) all reconstructed primary probeParticles_MC & recPid_TPC+TOF || x x x
305
306 <<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<
307 */
308
309 // ------------------------------------------------------------------
310 // -- Create THnSparseF - Cont
311 // ------------------------------------------------------------------
312
313 // cent: etaMC: yMC: phiMC: ptMC: sign:contPart: contSign: etaRec: phiRec: ptRec:deltaEta: deltaPhi: deltaPt
314 Int_t binHnCont[14] = { iCent, iEta, iRap, iPhi, iPt, iSign, 8, iSign, iEta, iPhi, iPt, iEta, 2*iPhi+1, 2*iPt+1};
315 Double_t minHnCont[14] = {dCent[0], dEta[0], dRap[0], dPhi[0], dPt[0], dSign[0], 0.5, dSign[0], dEta[0], dPhi[0], dPt[0], dEta[0], -dPhi[1], -dPt[1]};
316 Double_t maxHnCont[14] = {dCent[1], dEta[1], dRap[1], dPhi[1], dPt[1], dSign[1], 8.5, dSign[1], dEta[1], dPhi[1], dPt[1], dEta[1], dPhi[1], dPt[1]};
317 fHnCont = new THnSparseF("fHnCont", "cent:etaMC:yMC:phiMC:ptMC:sign:contPart:contSign:etaRec:phiRec:ptRec:deltaEta:deltaPhi:deltaPt",
318 14, binHnCont, minHnCont, maxHnCont);
319
320 fHnCont->Sumw2();
321
322 fHnCont->GetAxis(0)->SetTitle("centrality"); // 0-5|5-10|10-20|20-30|30-40|40-50|50-60|60-70|70-80|80-90 --> 10 bins
323 fHnCont->GetAxis(1)->SetTitle("#eta_{MC}"); // eta [-0.9,0.9]
324 fHnCont->GetAxis(2)->SetTitle("#it{y}_{MC}"); // rapidity [-0.5, 0.5]
325 fHnCont->GetAxis(3)->SetTitle("#varphi_{MC} (rad)"); // phi [ 0. ,2Pi]
326 fHnCont->GetAxis(4)->SetTitle("#it{p}_{T,MC} (GeV/#it{c})"); // pT [ 0.1,1.3]
327 fHnCont->GetAxis(5)->SetTitle("sign"); // -1 | 0 | +1
328
329 fHnCont->GetAxis(6)->SetTitle("contPart"); // 1 pi | 2 K | 3 p | 4 e | 5 mu | 6 other | 7 p from WeakDecay | 8 p from Material
330 fHnCont->GetAxis(7)->SetTitle("contSign"); // -1 | 0 | +1
331
332 fHnCont->GetAxis(8)->SetTitle("#eta_{Rec}"); // eta [-0.9, 0.9]
333 fHnCont->GetAxis(9)->SetTitle("#varphi_{Rec} (rad)"); // phi [ 0. , 2Pi]
334 fHnCont->GetAxis(10)->SetTitle("#it{p}_{T,Rec} (GeV/#it{c})"); // pt [ 0.1, 3.0]
335
336 fHnCont->GetAxis(11)->SetTitle("#eta_{MC}-#eta_{Rec}"); // eta [-0.9,0.9]
337 fHnCont->GetAxis(12)->SetTitle("#varphi_{MC}-#varphi_{Rec} (rad)"); // phi [ 0. ,2Pi]
338 fHnCont->GetAxis(13)->SetTitle("#it{p}_{T,MC}-#it{p}_{T,Rec} (GeV/#it{c})"); // pt [ 3.0,3.0]
339
340 fHelper->BinLogAxis(fHnCont, 4);
341 fHelper->BinLogAxis(fHnCont, 10);
342
343 // ------------------------------------------------------------------
344
345 return;
346}
347
348//________________________________________________________________________
349void AliAnalysisNetParticleEffCont::FillMCLabels() {
350 // Fill MC labels
351 // Loop over ESD tracks and fill arrays with MC lables
352 // fLabelsRec[0] : all Tracks
353 // fLabelsRec[1] : all Tracks accepted by PID of TPC
354 // Check every accepted track if correctly identified
355 // otherwise check for contamination
356
357 for (Int_t idxTrack = 0; idxTrack < fNTracks; ++idxTrack) {
358 AliESDtrack *track = fESD->GetTrack(idxTrack);
359
360 // -- Check if track is accepted for basic parameters
361 if (!fHelper->IsTrackAcceptedBasicCharged(track))
362 continue;
363
364 // -- Check if accepted
365 if (!fESDTrackCuts->AcceptTrack(track))
366 continue;
367
368 // -- Check if accepted in rapidity window
369 Double_t yP;
370 if (!fHelper->IsTrackAcceptedRapidity(track, yP))
371 continue;
372
373 // -- Check if accepted with thighter DCA cuts
374 if (!fHelper->IsTrackAcceptedDCA(track))
375 continue;
376
377 Int_t label = TMath::Abs(track->GetLabel());
378
379 // -- Fill Label of all reconstructed
380 fLabelsRec[0][idxTrack] = label;
381
382 // -- Check if accepted by PID from TPC or TPC+TOF
383 Double_t pid[2];
384 if (!fHelper->IsTrackAcceptedPID(track, pid))
385 continue;
386
387 // -- Fill Label of all reconstructed && recPid_TPC+TOF
388 fLabelsRec[1][idxTrack] = label;
389
390 // -- Check for contamination and fill contamination THnSparse
391 CheckContTrack(label, track->GetSign(), idxTrack);
392
393 } // for (Int_t idxTrack = 0; idxTrack < fNTracks; ++idxTrack) {
394
395 return;
396}
397
9be43c8e 398//________________________________________________________________________
399void AliAnalysisNetParticleEffCont::FillMCLabelsAOD() {
400 // Fill MC labels for AOD analysis
401 // Loop over ESD tracks and fill arrays with MC lables
402 // fLabelsRec[0] : all Tracks
403 // fLabelsRec[1] : all Tracks accepted by PID of TPC
404 // Check every accepted track if correctly identified
405 // otherwise check for contamination
406
407 for (Int_t idxTrack = 0; idxTrack < fNTracks; ++idxTrack) {
408 AliAODTrack *track = fAOD->GetTrack(idxTrack);
409
410 // -- Check if track is accepted for basic parameters
411 if (!fHelper->IsTrackAcceptedBasicCharged(track))
412 continue;
413
414 // -- Check if accepted
415 if (!track->TestFilterBit(fAODtrackCutBit))
416 continue;
417
418 // -- Check if accepted in rapidity window
419 Double_t yP;
420 if (!fHelper->IsTrackAcceptedRapidity(track, yP))
421 continue;
422
423 // -- Check if accepted with thighter DCA cuts (not yet done for AODs XXX)
424 //if (!fHelper->IsTrackAcceptedDCA(track))
425 // continue;
426
427 Int_t label = TMath::Abs(track->GetLabel());
428
429 // -- Fill Label of all reconstructed
430 fLabelsRec[0][idxTrack] = label;
431
432 // -- Check if accepted by PID from TPC or TPC+TOF
433 Double_t pid[2];
434 if (!fHelper->IsTrackAcceptedPID(track, pid))
435 continue;
436
437 // -- Fill Label of all reconstructed && recPid_TPC+TOF
438 fLabelsRec[1][idxTrack] = label;
439
440 // -- Check for contamination and fill contamination THnSparse
441 CheckContTrackAOD(label, track->Charge(), idxTrack);
442
443 } // for (Int_t idxTrack = 0; idxTrack < fNTracks; ++idxTrack) {
444
445 return;
446}
447
cb68eb1d 448//________________________________________________________________________
449void AliAnalysisNetParticleEffCont::CheckContTrack(Int_t label, Float_t sign, Int_t idxTrack) {
9be43c8e 450 // Check if particle is contamination or correctly identified for ESDs
cb68eb1d 451 // Fill contamination THnSparse
9be43c8e 452
cb68eb1d 453 TParticle* particle = fStack->Particle(label);
454 if (!particle)
455 return;
456
457 Bool_t isSecondaryFromWeakDecay = kFALSE;
458 Bool_t isSecondaryFromMaterial = kFALSE;
459
460 // -- Check if correctly identified
461 if (particle->GetPdgCode() == (sign*fPdgCode)) {
462
463 // Check if is physical primary -> all ok
464 if (fStack->IsPhysicalPrimary(label))
465 return;
466
467 // -- Check if secondaries from material or weak decay
468 isSecondaryFromWeakDecay = fStack->IsSecondaryFromWeakDecay(label);
469 isSecondaryFromMaterial = fStack->IsSecondaryFromMaterial(label);
470 }
471
472 // -- Get MC pdg
473 Float_t contSign = 0.;
474 if (particle->GetPDG()->Charge() == 0.) contSign = 0.;
475 else if (particle->GetPDG()->Charge() < 0.) contSign = -1.;
476 else if (particle->GetPDG()->Charge() > 0.) contSign = 1.;
477
478 // -- Get contaminating particle
479 Float_t contPart = 0;
480 if (isSecondaryFromWeakDecay) contPart = 7; // probeParticle from WeakDecay
481 else if (isSecondaryFromMaterial) contPart = 8; // probeParticle from Material
482 else {
483 if (TMath::Abs(particle->GetPdgCode()) == 211) contPart = 1; // pion
484 else if (TMath::Abs(particle->GetPdgCode()) == 321) contPart = 2; // kaon
485 else if (TMath::Abs(particle->GetPdgCode()) == 2212) contPart = 3; // proton
486 else if (TMath::Abs(particle->GetPdgCode()) == 11) contPart = 4; // electron
487 else if (TMath::Abs(particle->GetPdgCode()) == 13) contPart = 5; // muon
488 else contPart = 6; // other
489 }
490
491 // -- Get Reconstructed values
492 Float_t etaRec = 0.;
493 Float_t phiRec = 0.;
494 Float_t ptRec = 0.;
495
496 AliESDtrack *track = fESD->GetTrack(idxTrack);
497 if (track) {
498 // if no track present (which should not happen)
499 // -> pt = 0. , which is not in the looked at range
500
501 // -- Get Reconstructed eta/phi/pt
502 etaRec = track->Eta();
503 phiRec = track->Phi();
504 ptRec = track->Pt();
505 }
506
507 // -- Fill THnSparse
508 Double_t hnCont[14] = {fCentralityBin, particle->Eta(), particle->Y(), particle->Phi(), particle->Pt(), sign,
509 contPart, contSign, etaRec, phiRec, ptRec,
510 particle->Eta()-etaRec, particle->Phi()-phiRec, particle->Pt()-ptRec};
511 fHnCont->Fill(hnCont);
512}
513
9be43c8e 514//________________________________________________________________________
515void AliAnalysisNetParticleEffCont::CheckContTrackAOD(Int_t label, Float_t sign, Int_t idxTrack) {
516 // Check if particle is contamination or correctly identified for AODs
517 // Fill contamination THnSparse
518
519 AliAODMCParticle* particle = (AliAODMCParticle*)fArrayMC->At(label);
520
521 if (!particle)
522 return;
523
524 Bool_t isSecondaryFromWeakDecay = kFALSE;
525 Bool_t isSecondaryFromMaterial = kFALSE;
526
527 // -- Check if correctly identified
528 if (particle->GetPdgCode() == (sign*fPdgCode)) {
529
530 // Check if is physical primary -> all ok
531 if(particle->IsPhysicalPrimary())
532 return;
533
534 // -- Check if secondaries from material or weak decay
535 isSecondaryFromWeakDecay = particle->IsSecondaryFromWeakDecay();
536 isSecondaryFromMaterial = particle->IsSecondaryFromMaterial();
537
538 }
539
540 // -- Get MC pdg
541 Float_t contSign = 0.;
542 if ((TDatabasePDG::Instance()->GetParticle(particle->PdgCode()))->Charge() == 0.) contSign = 0.;
543 else if ((TDatabasePDG::Instance()->GetParticle(particle->PdgCode()))->Charge() < 0.) contSign = -1.;
544 else if ((TDatabasePDG::Instance()->GetParticle(particle->PdgCode()))->Charge() > 0.) contSign = 1.;
545
546 // -- Get contaminating particle
547 Float_t contPart = 0;
548 if (isSecondaryFromWeakDecay) contPart = 7; // probeParticle from WeakDecay
549 else if (isSecondaryFromMaterial) contPart = 8; // probeParticle from Material
550 else {
551 if (TMath::Abs(particle->GetPdgCode()) == 211) contPart = 1; // pion
552 else if (TMath::Abs(particle->GetPdgCode()) == 321) contPart = 2; // kaon
553 else if (TMath::Abs(particle->GetPdgCode()) == 2212) contPart = 3; // proton
554 else if (TMath::Abs(particle->GetPdgCode()) == 11) contPart = 4; // electron
555 else if (TMath::Abs(particle->GetPdgCode()) == 13) contPart = 5; // muon
556 else contPart = 6; // other
557 }
558
559 // -- Get Reconstructed values
560 Float_t etaRec = 0.;
561 Float_t phiRec = 0.;
562 Float_t ptRec = 0.;
563
564 AliAODTrack *track = fAOD->GetTrack(idxTrack);
565
566 if (track) {
567 // if no track present (which should not happen)
568 // -> pt = 0. , which is not in the looked at range
569
570 // -- Get Reconstructed eta/phi/pt
571 etaRec = track->Eta();
572 phiRec = track->Phi();
573 ptRec = track->Pt();
574 }
575
576 // -- Fill THnSparse
577 Double_t hnCont[14] = {fCentralityBin, particle->Eta(), particle->Y(), particle->Phi(), particle->Pt(), sign,
578 contPart, contSign, etaRec, phiRec, ptRec,
579 particle->Eta()-etaRec, particle->Phi()-phiRec, particle->Pt()-ptRec};
580 fHnCont->Fill(hnCont);
581}
582
cb68eb1d 583//________________________________________________________________________
584void AliAnalysisNetParticleEffCont::FillMCEffHist() {
9be43c8e 585 // Fill efficiency THnSparse for ESDs
cb68eb1d 586
587 Int_t nPart = fStack->GetNprimary();
588
589 for (Int_t idxMC = 0; idxMC < nPart; ++idxMC) {
590 TParticle* particle = fStack->Particle(idxMC);
591
592 // -- Check basic MC properties -> charged physical primary
593 if (!fHelper->IsParticleAcceptedBasicCharged(particle, idxMC))
594 continue;
595
596 // -- Check rapidity window
597 Double_t yMC;
598 if (!fHelper->IsParticleAcceptedRapidity(particle, yMC))
599 continue;
600
601 // -- Check if probeParticle / anti-probeParticle
602 if (TMath::Abs(particle->GetPdgCode()) != fPdgCode)
603 continue;
604
605 // -- Get sign of particle
606 Float_t sign = (particle->GetPdgCode() < 0) ? -1. : 1.;
607
608 // -- Get if particle is findable
609 Float_t findable = Float_t(fHelper->IsParticleFindable(idxMC));
610
611 // -- Get recStatus and pidStatus
612 Float_t recStatus = 0.;
613 Float_t recPid = 0.;
614
615 // -- Get Reconstructed values
616 Float_t etaRec = 0.;
617 Float_t phiRec = 0.;
618 Float_t ptRec = 0.;
619
620 // -- Loop over all labels
621 for (Int_t idxRec=0; idxRec < fNTracks; ++idxRec) {
622 if (idxMC == fLabelsRec[0][idxRec]) {
623 recStatus = 1.;
624
625 if (idxMC == fLabelsRec[1][idxRec])
626 recPid = 1.;
627
628 AliESDtrack *track = fESD->GetTrack(idxRec);
629 if (track) {
630 // if no track present (which should not happen)
631 // -> pt = 0. , which is not in the looked at range
632
633 // -- Get Reconstructed eta/phi/pt
634 etaRec = track->Eta();
635 phiRec = track->Phi();
636 ptRec = track->Pt();
637 }
638 break;
639 }
640 } // for (Int_t idxRec=0; idxRec < fNTracks; ++idxRec) {
641
642 // -- Fill THnSparse
643 Double_t hnEff[15] = {fCentralityBin, particle->Eta(), particle->Y(), particle->Phi(), particle->Pt(), sign,
644 findable, recStatus, recPid, etaRec, phiRec, ptRec,
645 particle->Eta()-etaRec, particle->Phi()-phiRec, particle->Pt()-ptRec};
646 fHnEff->Fill(hnEff);
647
648 } // for (Int_t idxMC = 0; idxMC < nPart; ++idxMC) {
649
650 return;
651}
9be43c8e 652
653//________________________________________________________________________
654void AliAnalysisNetParticleEffCont::FillMCEffHistAOD() {
655 // Fill efficiency THnSparse for AODs
656
657 Int_t nPart = fArrayMC->GetEntriesFast();
658
659 for (Int_t idxMC = 0; idxMC < nPart; ++idxMC) {
660
661 AliAODMCParticle* particle = (AliAODMCParticle*)fArrayMC->At(idxMC);
662
663 // -- Check basic MC properties -> charged physical primary
664 if (!fHelper->IsParticleAcceptedBasicCharged(particle))
665 continue;
666
667 // -- Check rapidity window
668 Double_t yMC;
669 if (!fHelper->IsParticleAcceptedRapidity((TParticle*)particle, yMC))
670 continue;
671
672 // -- Check if probeParticle / anti-probeParticle
673 if (TMath::Abs(particle->GetPdgCode()) != fPdgCode)
674 continue;
675
676 // -- Get sign of particle
677 Float_t sign = (particle->GetPdgCode() < 0) ? -1. : 1.;
678
679 // -- Get if particle is findable
680 Float_t findable = 1.;// Float_t(fHelper->IsParticleFindable(idxMC)); XXX
681
682 // -- Get recStatus and pidStatus
683 Float_t recStatus = 0.;
684 Float_t recPid = 0.;
685
686 // -- Get Reconstructed values
687 Float_t etaRec = 0.;
688 Float_t phiRec = 0.;
689 Float_t ptRec = 0.;
690
691 // -- Loop over all labels
692 for (Int_t idxRec=0; idxRec < fNTracks; ++idxRec) {
693 if (idxMC == fLabelsRec[0][idxRec]) {
694 recStatus = 1.;
695
696 if (idxMC == fLabelsRec[1][idxRec])
697 recPid = 1.;
698
699 AliVTrack *track = NULL;
700 if(fESD)
701 track = fESD->GetTrack(idxRec);
702 else if(fAOD)
703 track = fAOD->GetTrack(idxRec);
704
705 if (track) {
706 // if no track present (which should not happen)
707 // -> pt = 0. , which is not in the looked at range
708
709 // -- Get Reconstructed eta/phi/pt
710 etaRec = track->Eta();
711 phiRec = track->Phi();
712 ptRec = track->Pt();
713 }
714 break;
715 }
716 } // for (Int_t idxRec=0; idxRec < fNTracks; ++idxRec) {
717
718 // -- Fill THnSparse
719 Double_t hnEff[15] = {fCentralityBin, particle->Eta(), particle->Y(), particle->Phi(), particle->Pt(), sign,
720 findable, recStatus, recPid, etaRec, phiRec, ptRec,
721 particle->Eta()-etaRec, particle->Phi()-phiRec, particle->Pt()-ptRec};
722 fHnEff->Fill(hnEff);
723
724 } // for (Int_t idxMC = 0; idxMC < nPart; ++idxMC) {
725
726 return;
727}