AliAODEvent::GetHeader() returns AliVHeader
[u/mrichter/AliRoot.git] / PWGPP / EvTrkSelection / AliCFSingleTrackEfficiencyTask.cxx
CommitLineData
c0757458 1#include "TCanvas.h"
2#include "TParticle.h"
3#include "TH1I.h"
4#include <TDatabasePDG.h>
5
6#include "AliAnalysisDataSlot.h"
7#include "AliAnalysisDataContainer.h"
8#include "AliAnalysisManager.h"
9#include "AliAODEvent.h"
10#include "AliAODMCParticle.h"
11#include "AliAODTrack.h"
12#include "AliAODTracklets.h"
13#include "AliMultiplicity.h"
14#include "AliCFManager.h"
15#include "AliCFCutBase.h"
16#include "AliCFContainer.h"
17#include "AliESDEvent.h"
18#include "AliESDtrack.h"
19#include "AliESDtrackCuts.h"
20#include "AliESDVertex.h"
21#include "AliInputEventHandler.h"
22#include "AliMCEvent.h"
23#include "AliVEvent.h"
24#include "AliAODEvent.h"
25#include "AliVVertex.h"
26#include "AliLog.h"
27#include "AliStack.h"
28#include "AliGenEventHeader.h"
29#include "AliAODMCHeader.h"
2d9f660a 30#include "AliCentrality.h"
c0757458 31
32#include "AliSingleTrackEffCuts.h"
33#include "AliCFSingleTrackEfficiencyTask.h"
34
35
36ClassImp(AliCFSingleTrackEfficiencyTask)
37
38//__________________________________________________________________________
39AliCFSingleTrackEfficiencyTask::AliCFSingleTrackEfficiencyTask() :
40 fReadAODData(0),
41 fCFManager(0x0),
42 fQAHistList(0x0),
43 fTrackCuts(0x0),
44 fMCCuts(0x0),
45 fTriggerMask(AliVEvent::kAny),
46 fSetFilterBit(kFALSE),
47 fbit(0),
2d9f660a 48 fEvalCentrality(kFALSE),
49 fCentralityEstimator("V0M"),
c0757458 50 fHistEventsProcessed(0x0)
51{
52 //
53 //Default constructor
54 //
55}
56
57//___________________________________________________________________________
58AliCFSingleTrackEfficiencyTask::AliCFSingleTrackEfficiencyTask(const Char_t* name,AliESDtrackCuts *trackcuts, AliSingleTrackEffCuts * mccuts) :
59 AliAnalysisTaskSE(name),
60 fReadAODData(0),
61 fCFManager(0x0),
62 fQAHistList(0x0),
63 fTrackCuts(trackcuts),
64 fMCCuts(mccuts),
65 fTriggerMask(AliVEvent::kAny),
66 fSetFilterBit(kFALSE),
67 fbit(0),
2d9f660a 68 fEvalCentrality(kFALSE),
69 fCentralityEstimator("V0M"),
c0757458 70 fHistEventsProcessed(0x0)
71{
72 //
73 // Constructor. Initialization of Inputs and Outputs
74 //
75 Info("AliCFSingleTrackEfficiencyTask","Calling Constructor");
76
77 // Output slot #1 writes into a TList container (nevents histogran)
78 DefineOutput(1,TH1I::Class());
79 // Output slot #2 writes into a TList container (distributions)
80 DefineOutput(2,AliCFContainer::Class());
81 // Output slot #3 writes the QA list
82 DefineOutput(3,TList::Class());
83 // Output slot #3 writes the ESD track cuts used
84 DefineOutput(4,AliESDtrackCuts::Class());
85 // Output slot #4 writes the particle and event selection object
86 DefineOutput(5,AliSingleTrackEffCuts::Class());
87}
88
89//_________________________________________________________________________________________________________________
90AliCFSingleTrackEfficiencyTask& AliCFSingleTrackEfficiencyTask::operator=(const AliCFSingleTrackEfficiencyTask& c)
91{
92 //
93 // Assignment operator
94 //
95 if (this!=&c) {
96 AliAnalysisTaskSE::operator=(c) ;
97
2d9f660a 98 fReadAODData = c.fReadAODData;
c0757458 99 fCFManager = c.fCFManager;
2d9f660a 100 fQAHistList = c.fQAHistList;
c0757458 101
102 if(c.fTrackCuts) { delete fTrackCuts; fTrackCuts = new AliESDtrackCuts(*(c.fTrackCuts)); }
103 if(c.fMCCuts) { delete fMCCuts; fMCCuts = new AliSingleTrackEffCuts(*(c.fMCCuts)); }
104 fTriggerMask = c.fTriggerMask;
105
106 fSetFilterBit = c.fSetFilterBit;
107 fbit = c.fbit;
2d9f660a 108
109 fEvalCentrality = c.fEvalCentrality;
110 fCentralityEstimator = c.fCentralityEstimator;
111
c0757458 112 fHistEventsProcessed = c.fHistEventsProcessed;
113 }
114 return *this;
115}
116
117//________________________________________________________________________________________________________
118AliCFSingleTrackEfficiencyTask::AliCFSingleTrackEfficiencyTask(const AliCFSingleTrackEfficiencyTask& c) :
119 AliAnalysisTaskSE(c),
120 fReadAODData(c.fReadAODData),
121 fCFManager(c.fCFManager),
122 fQAHistList(c.fQAHistList),
123 fTrackCuts(c.fTrackCuts),
124 fMCCuts(c.fMCCuts),
125 fTriggerMask(c.fTriggerMask),
126 fSetFilterBit(c.fSetFilterBit),
127 fbit(c.fbit),
2d9f660a 128 fEvalCentrality(c.fEvalCentrality),
129 fCentralityEstimator(c.fCentralityEstimator),
c0757458 130 fHistEventsProcessed(c.fHistEventsProcessed)
131{
132 //
133 // Copy Constructor
134 //
135}
136
137//___________________________________________________________________________
138AliCFSingleTrackEfficiencyTask::~AliCFSingleTrackEfficiencyTask()
139{
140 //
141 // Destructor
142 //
143 Info("~AliCFSingleTrackEfficiencyTask","Calling Destructor");
144
145 if (fCFManager) delete fCFManager;
146 if (fHistEventsProcessed) delete fHistEventsProcessed;
147 if (fQAHistList) { fQAHistList->Clear(); delete fQAHistList; }
148 if (fTrackCuts) delete fTrackCuts;
149 if (fMCCuts) delete fMCCuts;
150}
151
152//_________________________________________________
153void AliCFSingleTrackEfficiencyTask::Init()
154{
155 //
156 // Initialization, checks + copy cuts
157 //
158 if(!fMCCuts) {
159 AliFatal(" MC Cuts not defined");
160 return;
161 }
162 if(!fTrackCuts) {
163 AliFatal(" Track Cuts not defined");
164 return;
165 }
166
167 AliESDtrackCuts* copyfTrackCuts = new AliESDtrackCuts(*fTrackCuts);
168 const char* nameoutput = GetOutputSlot(4)->GetContainer()->GetName();
169 copyfTrackCuts->SetName(nameoutput);
170 // Post the data
171 PostData(4,copyfTrackCuts);
172
173 AliSingleTrackEffCuts* copyfMCCuts = new AliSingleTrackEffCuts(*fMCCuts);
174 nameoutput = GetOutputSlot(5)->GetContainer()->GetName();
175 copyfMCCuts->SetName(nameoutput);
176 // Post the data
177 PostData(5,copyfMCCuts);
178
2d9f660a 179
180 if(fEvalCentrality) {
181 Bool_t isCentEstimatorOk = kFALSE;
182 TString validEstimators[9] = { "V0M", "V0A", "V0C", "TRK", "TKL", "CL1", "ZNA", "ZNC" "ZPA" };
183 for(Int_t i=0; i<9; i++ ) {
184 if(fCentralityEstimator==validEstimators[i]) isCentEstimatorOk = kTRUE;
185 }
186 if(!isCentEstimatorOk) {
187 AliFatal(Form("Chosen centrality estimator %s is not valid\n",fCentralityEstimator.Data()));
188 return;
189 }
190 }
191
c0757458 192 return;
193}
194
195//_________________________________________________________________
196void AliCFSingleTrackEfficiencyTask::UserExec(Option_t *)
197{
198 //
199 // User Exec
200 //
201
5182a19b 202 Info("UserExec","Start of method") ;
c0757458 203
204 AliVEvent* event = fInputEvent;
205
206 if(!fInputEvent) {
207 AliFatal("NO EVENT FOUND!");
208 return;
209 }
210 if(!fMCEvent) {
211 AliFatal("NO MC INFO FOUND");
212 return;
213 }
214
215 fHistEventsProcessed->Fill(0.5); // # of Event proceed
216 Bool_t IsEventMCSelected = kFALSE;
217 Bool_t isAOD = fInputEvent->IsA()->InheritsFrom("AliAODEvent");
218
219 //
220 // Step 0: MC Gen Event Selection
221 //
222 if(isAOD) {
223 if(!event && AODEvent() && IsStandardAOD()) {
224 // In case there is an AOD handler writing a standard AOD, use the AOD
225 // event in memory rather than the input (ESD) event.
226 event = dynamic_cast<AliAODEvent*> (AODEvent());
227 }
228 IsEventMCSelected = fMCCuts->IsMCEventSelected(event);//AODs
229 } else {
230 IsEventMCSelected = fMCCuts->IsMCEventSelected(fMCEvent);//ESDs
231 }
232
233 // pass to the manager the event info to the cuts that need it
234 if(!isAOD) fCFManager->SetMCEventInfo(fMCEvent);
235 else fCFManager->SetMCEventInfo(event);
236 fCFManager->SetRecEventInfo(event);
237
238 if(!IsEventMCSelected) {
239 AliDebug(3,"MC event not passing the quality criteria \n");
240 PostData(1,fHistEventsProcessed);
241 PostData(2,fCFManager->GetParticleContainer());
242 PostData(3,fQAHistList);
243 return;
244 }
245 fHistEventsProcessed->Fill(1.5); // # of Event after passing MC cuts
246
247
248 //
249 // Step 1-3: Check the MC generated particles
250 //
251 if(!isAOD) CheckESDMCParticles();
252 else CheckAODMCParticles();
253
254 //
255 // Step 4-7: Reconstructed event and track selection
256 //
257 Bool_t isRecoEventOk = fMCCuts->IsRecoEventSelected(event);
258
259 if(isRecoEventOk) {
260 fHistEventsProcessed->Fill(2.5); // # of Event after passing all cuts
261 CheckReconstructedParticles();
262 }
263 else AliDebug(3,"Event not passing quality criteria\n");
264
265
266 PostData(1,fHistEventsProcessed);
267 PostData(2,fCFManager->GetParticleContainer());
268 PostData(3,fQAHistList);
269
270 return;
271}
272
273
274//___________________________________________________________________________
275void AliCFSingleTrackEfficiencyTask::Terminate(Option_t*)
276{
277 //
278 // Terminate
279 //
280
5182a19b 281 Info("Terminate","Start and end of Method");
c0757458 282 AliAnalysisTaskSE::Terminate();
283
284 /*
285 //draw some example plots....
286 AliCFContainer *cont= dynamic_cast<AliCFContainer*> (GetOutputData(2));
287
288 TH1D* h00 = cont->ShowProjection(5,0) ;
289 TH1D* h01 = cont->ShowProjection(5,1) ;
290 TH1D* h02 = cont->ShowProjection(5,2) ;
291 TH1D* h03 = cont->ShowProjection(5,3) ;
292 TH1D* h04 = cont->ShowProjection(5,4) ;
293 TH1D* h05 = cont->ShowProjection(5,5) ;
294 TH1D* h06 = cont->ShowProjection(5,6) ;
2d9f660a 295 TH1D* h07 = cont->ShowProjection(5,7) ;
c0757458 296
297 h00->SetMarkerStyle(23) ;
298 h01->SetMarkerStyle(24) ;
299 h02->SetMarkerStyle(25) ;
300 h03->SetMarkerStyle(26) ;
301 h04->SetMarkerStyle(27) ;
302 h05->SetMarkerStyle(28) ;
303 h06->SetMarkerStyle(28) ;
304
305
306
307 TCanvas * c =new TCanvas("c","",1400,800);
308 c->Divide(4,2);
309
310 c->cd(1);
311 h00->Draw("p");
312 c->cd(2);
313 h01->Draw("p");
314 c->cd(3);
315 h02->Draw("p");
316 c->cd(4);
317 h03->Draw("p");
318 c->cd(5);
319 h04->Draw("p");
320 c->cd(6);
321 h05->Draw("p");
322 c->cd(7);
323 h06->Draw("p");
2d9f660a 324 c->cd(8);
325 h07->Draw("p");
c0757458 326
327 c->SaveAs("plots.eps");
328 */
c0757458 329}
330
331
332//___________________________________________________________________________
333void AliCFSingleTrackEfficiencyTask::UserCreateOutputObjects()
334{
335 //
336 // UserCreateOutputObjects
337 //
338
339 Info("CreateOutputObjects","CreateOutputObjects of task %s", GetName());
340
341 //slot #1
342 OpenFile(1);
343
344 const char* nameoutput=GetOutputSlot(1)->GetContainer()->GetName();
345 // fHistEventsProcessed = new TH1I(nameoutput"fHistEventsProcessed","fHistEventsProcessed",3,0,3) ;
346 fHistEventsProcessed = new TH1I(nameoutput,"fHistEventsProcessed",3,0,3) ;
347 fHistEventsProcessed->GetXaxis()->SetBinLabel(1,"All events");
348 fHistEventsProcessed->GetXaxis()->SetBinLabel(2,"Good MC events");
349 fHistEventsProcessed->GetXaxis()->SetBinLabel(3,"Good Reconstructed events");
350
351 PostData(1,fHistEventsProcessed) ;
352 PostData(2,fCFManager->GetParticleContainer()) ;
353 PostData(3,fQAHistList) ;
354
355 return;
356}
357
358
359//_________________________________________________________________________
360void AliCFSingleTrackEfficiencyTask::CheckESDMCParticles()
361{
362 //
363 // Check ESD generated particles
364 //
365 if (!fMCEvent) {
366 AliFatal("NO MC INFO FOUND");
367 return;
368 }
369
2d9f660a 370 Double_t containerInput[7] = { 0., 0., 0., 0., 0., 0., 0. };
c0757458 371
372 TArrayF vtxPos(3);
373 AliGenEventHeader *genHeader = NULL;
374 genHeader = fMCEvent->GenEventHeader();
375 genHeader->PrimaryVertex(vtxPos);
376 containerInput[4] = vtxPos[2]; // z-vtx position
377
378 Double_t multiplicity = (Double_t)GetNumberOfTrackletsInEtaRange(-1.0,1.0);
379 containerInput[5] = multiplicity; //reconstructed number of tracklets
380
2d9f660a 381 Double_t centrality = -1.;
382 if(fEvalCentrality) centrality = GetCentrality();
383 containerInput[6] = centrality;
384
c0757458 385 // loop on the MC Generated Particles
386 for (Int_t ipart=0; ipart<fMCEvent->GetNumberOfTracks(); ipart++) {
387
388 AliMCParticle *mcPart = (AliMCParticle*)fMCEvent->GetTrack(ipart);
389 containerInput[0] = (Float_t)mcPart->Pt();
390 containerInput[1] = mcPart->Eta();
391 containerInput[2] = mcPart->Phi();
392 containerInput[3] = mcPart->Theta();
393
394 // Step 1. Particle passing through Generation criteria and filling
395 if( !fMCCuts->IsMCParticleGenerated(mcPart) ) {
396 AliDebug(3,"MC Particle not passing through genetations criteria\n");
397 continue;
398 }
399 fCFManager->GetParticleContainer()->Fill(containerInput,kStepMCGenCut);
400
401 // Step 2. Particle passing through Kinematic criteria and filling
402 if( !fMCCuts->IsMCParticleInKineAcceptance(mcPart) ) {
403 AliDebug(3,"MC Particle not in the kine acceptance\n");
404 continue;
405 }
406 fCFManager->GetParticleContainer()->Fill(containerInput,kStepMCKineCut);
407
408 // Step 3. Particle passing through Track ref criteria and filling
409 // did leave signal (enough clusters) on the detector
410 if( !fMCCuts->IsMCParticleInReconstructable(mcPart) ) {
411 AliDebug(3,"MC Particle not in the reconstructible\n");
412 continue;
413 }
414 fCFManager->GetParticleContainer()->Fill(containerInput,kStepMCAccpCut);
415
416 }// end of particle loop
417
418 return;
419}
420
421
422//________________________________________________________________________
423void AliCFSingleTrackEfficiencyTask::CheckAODMCParticles()
424{
425 //
426 // Check AOD generated particles
427 //
428 if (!fInputEvent) {
429 AliFatal("NO EVENT FOUND!");
430 return;
431 }
432
433 AliAODEvent* event = dynamic_cast<AliAODEvent*>(fInputEvent);
434 if(!event && AODEvent() && IsStandardAOD()) {
435 // In case there is an AOD handler writing a standard AOD, use the AOD
436 // event in memory rather than the input (ESD) event.
437 event = dynamic_cast<AliAODEvent*> (AODEvent());
438 }
439 TClonesArray* mcArray = dynamic_cast<TClonesArray*>(event->FindListObject(AliAODMCParticle::StdBranchName()));
440 if (!mcArray) {
441 AliError("Could not find Monte-Carlo in AOD");
442 return;
443 }
444 AliAODMCHeader *mcHeader=NULL;
445 mcHeader = dynamic_cast<AliAODMCHeader*>(event->GetList()->FindObject(AliAODMCHeader::StdBranchName()));
446 if (!mcHeader) {
447 AliError("Could not find MC Header in AOD");
448 return;
449 }
450
2d9f660a 451 Double_t containerInput[7] = { 0., 0., 0., 0., 0., 0., 0. };
c0757458 452 // Set the z-vertex position
453 containerInput[4] = mcHeader->GetVtxZ();
454 // Multiplicity of the event defined as Ntracklets |eta|<1.0
455 Double_t multiplicity = (Double_t)GetNumberOfTrackletsInEtaRange(-1.0,1.0);
456 containerInput[5] = multiplicity;
2d9f660a 457 // Determine the event centrality
458 Double_t centrality = 0.;
459 if(fEvalCentrality) centrality = GetCentrality();
460 containerInput[6] = centrality;
c0757458 461
462 for (Int_t ipart=0; ipart<mcArray->GetEntriesFast(); ipart++) {
463
464 AliAODMCParticle* mcPart = dynamic_cast<AliAODMCParticle*>(mcArray->At(ipart));
465 if (!mcPart){
466 AliError("Failed casting particle from MC array!, Skipping particle");
467 continue;
468 }
469 containerInput[0] = (Float_t)mcPart->Pt();
470 containerInput[1] = mcPart->Eta();
471 containerInput[2] = mcPart->Phi();
472 containerInput[3] = mcPart->Theta();
473
474 // Step 1. Particle passing through Generation criteria and filling
475 if( !fMCCuts->IsMCParticleGenerated(mcPart) ) {
476 AliDebug(3,"MC Particle not passing quality criteria\n");
477 continue;
478 }
479 fCFManager->GetParticleContainer()->Fill(containerInput,kStepMCGenCut);
480
481 // Step 2. Particle passing through Kinematic criteria and filling
482 if( !fMCCuts->IsMCParticleInKineAcceptance(mcPart) ) {
483 AliDebug(3,"MC Particle not in the acceptance\n");
484 continue;
485 }
486 fCFManager->GetParticleContainer()->Fill(containerInput,kStepMCKineCut);
487
488 // Step 3. Particle passing through Track Ref criteria and filling
489 // but no info available for Track ref in AOD fillng same as above
490 fCFManager->GetParticleContainer()->Fill(containerInput,kStepMCAccpCut);
491 }
492
493 return;
494}
495
496
497//_______________________________________________________________________________
498AliESDtrack * AliCFSingleTrackEfficiencyTask::ConvertTrack(AliAODTrack *track)
499{
500 //
501 // Convert an AOD track to an ESD track to apply ESDtrackCuts
502 //
503
504 AliAODEvent *aodEvent = (AliAODEvent*)fInputEvent;
505 const AliAODVertex *primary = aodEvent->GetPrimaryVertex();
506 Double_t pos[3],cov[6];
507 primary->GetXYZ(pos);
508 primary->GetCovarianceMatrix(cov);
509 const AliESDVertex vESD(pos,cov,100.,100);
510
511 AliESDtrack *esdTrack = new AliESDtrack(track);
512 // set the TPC cluster info
513 esdTrack->SetTPCClusterMap(track->GetTPCClusterMap());
514 esdTrack->SetTPCSharedMap(track->GetTPCSharedMap());
515 esdTrack->SetTPCPointsF(track->GetTPCNclsF());
516 // needed to calculate the impact parameters
517 esdTrack->RelateToVertex(&vESD,0.,3.);
518 // std::cout << " primary vtx "<< primary << std::endl;
519 // std::cout << " esdvtx "<< vESD.GetName() << std::endl;
520 // std::cout<< " esdtrack pt "<< esdTrack.Pt() << " and status " << esdTrack.GetStatus() <<endl;
521 // std::cout << " aod track "<< track<< " and status " << track->GetStatus() << std::endl;
522 // std::cout << " esd track "<< esdTrack<< " and status " << esdTrack->GetStatus() << std::endl;
523 return esdTrack;
524}
525
526
527//___________________________________________________________________________
528void AliCFSingleTrackEfficiencyTask::CheckReconstructedParticles()
529{
530 //
531 // Check reconstructed particles
532 //
533
534 AliVEvent* event = fInputEvent;
535 Bool_t isAOD = fInputEvent->IsA()->InheritsFrom("AliAODEvent");
536 if(!event && AODEvent() && IsStandardAOD()) {
537 event = dynamic_cast<AliAODEvent*> (AODEvent());
538 }
c37908f3 539 if (!event) return;
c0757458 540
2d9f660a 541 Double_t containerInput[7] = { 0, 0, 0, 0, 0, 0, 0}; // contains reconstructed quantities
542 Double_t containerInputMC[7] = { 0, 0, 0, 0, 0, 0, 0}; // contains generated quantities
c0757458 543
544 const AliVVertex *vertex = event->GetPrimaryVertex();
545 // set the z-vertex position
546 containerInput[4] = vertex->GetZ();
547 containerInputMC[4] = containerInput[4];
548 // set the event multiplicity as Ntracklets in |eta|<1.0
549 Double_t multiplicity = (Double_t)GetNumberOfTrackletsInEtaRange(-1.0,1.0);
550 containerInput[5] = multiplicity;
551 containerInputMC[5] = multiplicity;
2d9f660a 552 // Determine the event centrality
553 Double_t centrality = 0.;
554 if(fEvalCentrality) centrality = GetCentrality();
555 containerInput[6] = centrality;
556 containerInputMC[6] = centrality;
c0757458 557
558 // Reco tracks track loop
559 AliVParticle* track = NULL;
560 for (Int_t iTrack = 0; iTrack<event->GetNumberOfTracks(); iTrack++) {
561
562 track = event->GetTrack(iTrack);
563 if(!track) {
564 AliDebug(3,Form("Track %d not found",iTrack));
565 continue;
566 }
567
568 Double_t mom[3];
569 track->PxPyPz(mom);
570 Double_t pt = TMath::Sqrt(mom[0]*mom[0]+mom[1]*mom[1]);
571 containerInput[0] = pt;
572 containerInput[1] = track->Eta();
573 containerInput[2] = track->Phi();
574 containerInput[3] = track->Theta();
575
576 //
577 // Step 4. Track that are recostructed, filling
578 //
579 fCFManager->GetParticleContainer()->Fill(containerInput,kStepReconstructed);
580
581 //
582 // Step 5. Track that are recostructed and pass acceptance cuts, filling
583 //
584 if (!fMCCuts->IsRecoParticleKineAcceptance(track)) continue;
585 fCFManager->GetParticleContainer()->Fill(containerInput,kStepRecoKineCuts);
586
587 // is track associated to particle ? if yes + implimenting the physical primary..
588 Int_t label = TMath::Abs(track->GetLabel());
589 if (label<=0) {
590 AliDebug(3,"Particle not matching MC label \n");
591 continue;
592 }
593
594 //
595 // Step 6-7. Track that are recostructed + Kine + Quality criteria, filling
596 //
597
598 // check particle selections at MC level
599 AliVParticle *mcPart = (AliVParticle*)fMCEvent->GetTrack(label);
600 if(!mcPart) continue;
601 containerInputMC[0] = (Float_t)mcPart->Pt();
602 containerInputMC[1] = mcPart->Eta();
603 containerInputMC[2] = mcPart->Phi();
604 containerInputMC[3] = mcPart->Theta();
605
606 if (!fMCCuts->IsMCParticleGenerated(mcPart)) continue;
0ba71228 607 // cout<< "MC matching did work"<<endl;
c0757458 608
609
610 // for filter bit selection
611 AliAODTrack *aodTrack = dynamic_cast<AliAODTrack*>(track);
92fabd03 612 if(isAOD && !aodTrack) continue;
c0757458 613 if(isAOD && fSetFilterBit) if (!aodTrack->TestFilterMask(BIT(fbit))) continue;
0ba71228 614 // cout<<" Filter bit check passed"<<endl;
c0757458 615
616 Bool_t isESDtrack = track->IsA()->InheritsFrom("AliESDtrack");
617 AliESDtrack *tmptrack = NULL;
618 if(isESDtrack) {
619 tmptrack = dynamic_cast<AliESDtrack*>(track); // ESDs
620 } else {
fd6ce779 621 if (aodTrack) tmptrack = ConvertTrack(aodTrack); // AODs
c0757458 622 }
c37908f3 623 if (!tmptrack) continue;
c0757458 624
625 // exclude global constrained and TPC only tracks (filter bits 128 and 512)
626 Int_t id = tmptrack->GetID();
627 if(isAOD && id<0) {
7a88f6d0 628 AliDebug(3,"Track removed bc corresponds to either filter bit 128 or 512 (TPC only tracks)\n");
629 delete tmptrack; tmptrack=NULL;
c0757458 630 continue;
631 }
632
633 // Apply ESD track cuts
634 if( !fTrackCuts->IsSelected(tmptrack) ){
635 AliDebug(3,"Reconstructed track not passing quality criteria\n");
7a88f6d0 636 if(isAOD) { delete tmptrack; tmptrack=NULL; }
c0757458 637 continue;
638 }
0ba71228 639 // cout<<" analysis cuts passed"<<endl;
c0757458 640 fCFManager->GetParticleContainer()->Fill(containerInputMC,kStepReconstructedMC);
641 fCFManager->GetParticleContainer()->Fill(containerInput,kStepRecoQualityCuts);
2d9f660a 642 // cout << " checking particle pid"<<endl;
c0757458 643
644 //
645 // Step8, PID check
646 //
647 if( !fMCCuts->IsRecoParticlePID(track) ){
648 AliDebug(3,"Reconstructed track not passing PID criteria\n");
7a88f6d0 649 if(isAOD) { delete tmptrack; tmptrack=NULL; }
c0757458 650 continue;
651 }
652 fCFManager->GetParticleContainer()->Fill(containerInput,kStepRecoPID);
2d9f660a 653 // cout << " all steps filled"<<endl;
c0757458 654
7a88f6d0 655 if(isAOD) { delete tmptrack; tmptrack=NULL; }
c0757458 656 }
657 return;
658}
659
660//______________________________________________________________________
661Int_t AliCFSingleTrackEfficiencyTask::GetNumberOfTrackletsInEtaRange(Double_t mineta, Double_t maxeta)
662{
663 //
664 // counts tracklets in given eta range
665 //
666
667 AliAODEvent* event = dynamic_cast<AliAODEvent*>(fInputEvent);
668 Bool_t isAOD = fInputEvent->IsA()->InheritsFrom("AliAODEvent");
669 if(!event && AODEvent() && IsStandardAOD()) {
670 event = dynamic_cast<AliAODEvent*> (AODEvent());
671 }
c37908f3 672 if (!event) return -1;
c0757458 673 Int_t count=0;
674
675 if(isAOD) {
676 AliAODTracklets* tracklets = event->GetTracklets();
677 Int_t nTr=tracklets->GetNumberOfTracklets();
678 for(Int_t iTr=0; iTr<nTr; iTr++){
679 Double_t theta=tracklets->GetTheta(iTr);
680 Double_t eta=-TMath::Log(TMath::Tan(theta/2.));
681 if(eta>mineta && eta<maxeta) count++;
682 }
683 } else {
684 AliESDEvent *esdEvent = dynamic_cast<AliESDEvent*>(fInputEvent);
c37908f3 685 if (!esdEvent) return -1;
c0757458 686 const AliMultiplicity *mult = esdEvent->GetMultiplicity();
687 if (mult) {
688 if (mult->GetNumberOfTracklets()>0) {
689 for (Int_t n=0; n<mult->GetNumberOfTracklets(); n++) {
690 Double_t eta = -TMath::Log( TMath::Tan( mult->GetTheta(n) / 2.0 ) );
691 if(eta>mineta && eta<maxeta) count++;
692 }
693 }
694 }
695 }
696
697 return count;
698}
2d9f660a 699
700//______________________________________________________________________
701Double_t AliCFSingleTrackEfficiencyTask::GetCentrality()
702{
703 //
704 // Get centrality
705 //
706
707 Bool_t isAOD = fInputEvent->IsA()->InheritsFrom("AliAODEvent");
708 Double_t cent = -1;
709
710 if(isAOD) {
711 AliAODEvent* aodEvent = dynamic_cast<AliAODEvent*>(fInputEvent);
0d16cd51 712 if(!aodEvent) return cent;
0a918d8d 713 AliAODHeader* header = dynamic_cast<AliAODHeader*>(aodEvent->GetHeader());
714 if(!header) AliFatal("Not a standard AOD");
0d16cd51 715 if(!header) return cent;
2d9f660a 716 AliCentrality *centrality = header->GetCentralityP();
0d16cd51 717 if(!centrality) return cent;
2d9f660a 718 // cout<<" about to get cent perc AOD"<<endl;
719 cent = centrality->GetCentralityPercentile(fCentralityEstimator.Data());
720 } else {
0d16cd51 721 AliESDEvent* esdEvent = dynamic_cast<AliESDEvent*>(fInputEvent);
722 if(!esdEvent) return cent;
723 AliCentrality *centrality = esdEvent->GetCentrality();
724 if(!centrality) return cent;
2d9f660a 725 cent = centrality->GetCentralityPercentile(fCentralityEstimator.Data());
726 }
727
728 return cent;
729}