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