1 //------------------------------------------------------------------------------
2 // Implementation of AliPerformanceEff class. It keeps information from
3 // comparison of reconstructed and MC particle tracks. In addtion,
4 // it keeps selection cuts used during comparison. The comparison
5 // information is stored in the ROOT histograms. Analysis of these
6 // histograms can be done by using Analyse() class function. The result of
7 // the analysis (histograms/graphs) are stored in the folder which is
8 // a data member of AliPerformanceEff.
10 // Author: J.Otwinowski 04/02/2008
11 //------------------------------------------------------------------------------
15 // after running comparison task, read the file, and get component
16 gROOT->LoadMacro("$ALICE_ROOT/PWGPP/Macros/LoadMyLibs.C");
18 TFile f("Output.root");
19 AliPerformanceEff * compObj = (AliPerformanceEff*)coutput->FindObject("AliPerformanceEff");
21 // Analyse comparison data
24 // the output histograms/graphs will be stored in the folder "folderEff"
25 compObj->GetAnalysisFolder()->ls("*");
27 // user can save whole comparison object (or only folder with anlysed histograms)
28 // in the seperate output file (e.g.)
29 TFile fout("Analysed_Eff.root","recreate");
30 compObj->Write(); // compObj->GetAnalysisFolder()->Write();
37 #include "THnSparse.h"
40 #include "AliESDtrack.h"
41 #include "AliRecInfoCuts.h"
42 #include "AliMCInfoCuts.h"
44 #include "AliESDVertex.h"
45 #include "AliExternalTrackParam.h"
46 #include "AliTracker.h"
47 #include "AliESDEvent.h"
48 #include "AliMCEvent.h"
49 #include "AliMCParticle.h"
50 #include "AliHeader.h"
51 #include "AliGenEventHeader.h"
53 #include "AliPerformanceEff.h"
57 ClassImp(AliPerformanceEff)
59 //_____________________________________________________________________________
60 AliPerformanceEff::AliPerformanceEff(const Char_t* name, const Char_t* title, Int_t analysisMode, Bool_t hptGenerator): AliPerformanceObject(name,title),
75 SetAnalysisMode(analysisMode);
76 SetHptGenerator(hptGenerator);
82 //_____________________________________________________________________________
83 AliPerformanceEff::~AliPerformanceEff()
87 if(fEffHisto) delete fEffHisto; fEffHisto=0;
88 if(fEffSecHisto) delete fEffSecHisto; fEffSecHisto=0;
89 if(fAnalysisFolder) delete fAnalysisFolder; fAnalysisFolder=0;
92 //_____________________________________________________________________________
93 void AliPerformanceEff::Init()
100 Double_t ptMin = 1.e-2, ptMax = 20.;
102 Double_t *binsPt = 0;
104 if (IsHptGenerator()) {
107 binsPt = CreateLogAxis(nPtBins,ptMin,ptMax);
111 Double_t binsPt[32] = {0.,0.05,0.1,0.15,0.2,0.25,0.3,0.35,0.4,0.45,0.5,0.55,0.6,0.7,0.8,0.9,1.0,1.2,1.4,1.6,1.8,2.0,2.25,2.5,2.75,3.,3.5,4.,5.,6.,8.,10.};
112 Double_t ptMin = 0., ptMax = 10.;
114 if(IsHptGenerator() == kTRUE) {
116 ptMin = 0.; ptMax = 100.;
120 //mceta:mcphi:mcpt:pid:recStatus:findable:charge
121 Int_t binsEffHisto[9]={30,144,nPtBins,5,2,2,3,fgkMaxClones+1,fgkMaxFakes+1};
122 Double_t minEffHisto[9]={-1.5,0.,ptMin,0.,0.,0.,-1.5,0,0};
123 Double_t maxEffHisto[9]={ 1.5,2.*TMath::Pi(), ptMax,5.,2.,2.,1.5,fgkMaxClones+1,fgkMaxFakes+1};
125 fEffHisto = new THnSparseF("fEffHisto","mceta:mcphi:mcpt:pid:recStatus:findable:charge:nclones:nfakes",9,binsEffHisto,minEffHisto,maxEffHisto);
126 fEffHisto->SetBinEdges(2,binsPt);
128 fEffHisto->GetAxis(0)->SetTitle("#eta_{mc}");
129 fEffHisto->GetAxis(1)->SetTitle("#phi_{mc} (rad)");
130 fEffHisto->GetAxis(2)->SetTitle("p_{Tmc} (GeV/c)");
131 fEffHisto->GetAxis(3)->SetTitle("pid");
132 fEffHisto->GetAxis(4)->SetTitle("recStatus");
133 fEffHisto->GetAxis(5)->SetTitle("findable");
134 fEffHisto->GetAxis(6)->SetTitle("charge");
135 fEffHisto->GetAxis(7)->SetTitle("nClones");
136 fEffHisto->GetAxis(8)->SetTitle("nFakes");
139 //mceta:mcphi:mcpt:pid:recStatus:findable:mcR:mother_phi:mother_eta:charge
140 Int_t binsEffSecHisto[12]={30,60,nPtBins,5,2,2,100,60,30,3,fgkMaxClones+1,fgkMaxFakes+1};
141 Double_t minEffSecHisto[12]={-1.5,0.,ptMin,0.,0.,0.,0.,0.,-1.5,-1.5,0,0};
142 Double_t maxEffSecHisto[12]={ 1.5,2.*TMath::Pi(), ptMax,5.,2.,2.,200,2.*TMath::Pi(),1.5,1.5,fgkMaxClones+1,fgkMaxFakes+1};
144 fEffSecHisto = new THnSparseF("fEffSecHisto","mceta:mcphi:mcpt:pid:recStatus:findable:mcR:mother_phi:mother_eta:charge:nclones:nfakes",12,binsEffSecHisto,minEffSecHisto,maxEffSecHisto);
145 fEffSecHisto->SetBinEdges(2,binsPt);
147 fEffSecHisto->GetAxis(0)->SetTitle("#eta_{mc}");
148 fEffSecHisto->GetAxis(1)->SetTitle("#phi_{mc} (rad)");
149 fEffSecHisto->GetAxis(2)->SetTitle("p_{Tmc} (GeV/c)");
150 fEffSecHisto->GetAxis(3)->SetTitle("pid");
151 fEffSecHisto->GetAxis(4)->SetTitle("recStatus");
152 fEffSecHisto->GetAxis(5)->SetTitle("findable");
153 fEffSecHisto->GetAxis(6)->SetTitle("mcR (cm)");
154 fEffSecHisto->GetAxis(7)->SetTitle("mother_phi (rad)");
155 fEffSecHisto->GetAxis(8)->SetTitle("mother_eta");
156 fEffSecHisto->GetAxis(9)->SetTitle("charge");
157 fEffSecHisto->GetAxis(10)->SetTitle("nClones");
158 fEffSecHisto->GetAxis(11)->SetTitle("nFakes");
159 fEffSecHisto->Sumw2();
163 AliDebug(AliLog::kError, "ERROR: Cannot find AliMCInfoCuts object");
165 AliDebug(AliLog::kError, "ERROR: Cannot find AliRecInfoCuts object");
168 fAnalysisFolder = CreateFolder("folderEff","Analysis Efficiency Folder");
171 //_____________________________________________________________________________
172 void AliPerformanceEff::ProcessTPC(AliMCEvent* const mcEvent, AliESDEvent *const esdEvent)
174 // Fill TPC only efficiency comparison information
175 if(!esdEvent) return;
178 AliStack *stack = mcEvent->Stack();
180 AliDebug(AliLog::kError, "Stack not available");
184 Int_t *labelsRec = NULL;
185 labelsRec = new Int_t[esdEvent->GetNumberOfTracks()];
188 Printf("Cannot create labelsRec");
191 for(Int_t i=0;i<esdEvent->GetNumberOfTracks();i++) { labelsRec[i] = 0; }
193 Int_t *labelsAllRec = NULL;
194 labelsAllRec = new Int_t[esdEvent->GetNumberOfTracks()];
197 Printf("Cannot create labelsAllRec");
200 for(Int_t i=0;i<esdEvent->GetNumberOfTracks();i++) { labelsAllRec[i] = 0; }
202 // loop over rec. tracks
203 AliESDtrack *track=0;
204 for (Int_t iTrack = 0; iTrack < esdEvent->GetNumberOfTracks(); iTrack++)
206 track = esdEvent->GetTrack(iTrack);
208 if(track->Charge()==0) continue;
210 // if not fUseKinkDaughters don't use tracks with kink index > 0
211 if(!fUseKinkDaughters && track->GetKinkIndex(0) > 0) continue;
213 //Int_t label = TMath::Abs(track->GetLabel());
214 Int_t label = track->GetTPCLabel(); //Use TPC-only label for TPC-only efficiency analysis
215 labelsAllRec[iTrack]=label;
218 if(IsRecTPC(track) != 0) {
219 labelsRec[iTrack]=label;
225 // MC histograms for efficiency studies
227 Int_t nPart = stack->GetNtrack();
228 //Int_t nPart = stack->GetNprimary();
229 for (Int_t iMc = 0; iMc < nPart; ++iMc)
231 if (iMc == 0) continue; //Cannot distinguish between track or fake track
232 TParticle* particle = stack->Particle(iMc);
233 if (!particle) continue;
234 if (!particle->GetPDG()) continue;
235 if (particle->GetPDG()->Charge() == 0.0) continue;
238 Bool_t prim = stack->IsPhysicalPrimary(iMc);
241 // --- check for double filling in stack
242 // use only particles with no daughters in the list of primaries
243 Int_t nDaughters = 0;// particle->GetNDaughters();
245 for( Int_t iDaught=0; iDaught < particle->GetNDaughters(); iDaught++ ) {
246 if( particle->GetDaughter(iDaught) < stack->GetNprimary() )
252 // --- check for double filling in stack
254 /*Bool_t findable = kFALSE;
255 for(Int_t iRec=0; iRec<esdEvent->GetNumberOfTracks(); ++iRec)
258 if(iMc == labelsAllRec[iRec])
260 findable = IsFindable(mcEvent,iMc);
264 Bool_t findable = IsFindable(mcEvent,iMc);
266 Bool_t recStatus = kFALSE;
267 Int_t nClones = 0, nFakes = 0;
268 for(Int_t iRec=0; iRec<esdEvent->GetNumberOfTracks(); ++iRec)
270 // check reconstructed
271 if(iMc == labelsRec[iRec])
273 if (recStatus && nClones < fgkMaxClones) nClones++;
276 //In order to relate the fake track to track parameters, we assign it to the best matching ESD track.
277 if (labelsRec[iRec] < 0 && -labelsRec[iRec] == iMc && nFakes < fgkMaxFakes) nFakes++;
280 // Only 5 charged particle species (e,mu,pi,K,p)
281 if (fCutsMC->IsPdgParticle(TMath::Abs(particle->GetPdgCode())) == kFALSE) continue;
283 // transform Pdg to Pid
284 Int_t pid = TransformToPID(particle);
286 Float_t mceta = particle->Eta();
287 Float_t mcphi = particle->Phi();
288 if(mcphi<0) mcphi += 2.*TMath::Pi();
289 Float_t mcpt = particle->Pt();
291 if (particle->GetPDG()->Charge() < 0) charge = -1.;
292 else if (particle->GetPDG()->Charge() > 0) charge = 1.;
295 Double_t vEffHisto[9] = {mceta, mcphi, mcpt, static_cast<Double_t>(pid), static_cast<Double_t>(recStatus), static_cast<Double_t>(findable), static_cast<Double_t>(charge), static_cast<Double_t>(nClones), static_cast<Double_t>(nFakes)};
296 fEffHisto->Fill(vEffHisto);
298 if(labelsRec) delete [] labelsRec; labelsRec = 0;
299 if(labelsAllRec) delete [] labelsAllRec; labelsAllRec = 0;
302 //_____________________________________________________________________________
303 void AliPerformanceEff::ProcessTPCSec(AliMCEvent* const mcEvent, AliESDEvent *const esdEvent)
305 // Fill TPC only efficiency comparison information for secondaries
307 if(!esdEvent) return;
309 Int_t *labelsRecSec = NULL;
310 labelsRecSec = new Int_t[esdEvent->GetNumberOfTracks()];
313 Printf("Cannot create labelsRecSec");
316 for(Int_t i=0;i<esdEvent->GetNumberOfTracks();i++) { labelsRecSec[i] = 0; }
318 Int_t *labelsAllRecSec = NULL;
319 labelsAllRecSec = new Int_t[esdEvent->GetNumberOfTracks()];
320 if(!labelsAllRecSec) {
321 delete [] labelsRecSec;
322 Printf("Cannot create labelsAllRecSec");
325 for(Int_t i=0;i<esdEvent->GetNumberOfTracks();i++) { labelsAllRecSec[i] = 0; }
327 // loop over rec. tracks
328 AliESDtrack *track=0;
329 Int_t multAll=0, multRec=0;
330 for (Int_t iTrack = 0; iTrack < esdEvent->GetNumberOfTracks(); iTrack++)
332 track = esdEvent->GetTrack(iTrack);
334 if(track->Charge()==0) continue;
336 // if not fUseKinkDaughters don't use tracks with kink index > 0
337 if(!fUseKinkDaughters && track->GetKinkIndex(0) > 0) continue;
339 //Int_t label = TMath::Abs(track->GetLabel());
340 Int_t label = track->GetTPCLabel(); //Use TPC-only label for TPC-only efficiency analysis
341 labelsAllRecSec[multAll]=label;
345 if(IsRecTPC(track) != 0) {
346 labelsRecSec[multRec]=label;
352 // MC histograms for efficiency studies
356 AliStack *stack = mcEvent->Stack();
359 Int_t nPart = stack->GetNtrack();
360 //Int_t nPart = stack->GetNprimary();
361 for (Int_t iMc = 0; iMc < nPart; ++iMc)
363 if (iMc == 0) continue; //Cannot distinguish between track or fake track
364 TParticle* particle = stack->Particle(iMc);
365 if (!particle) continue;
366 if (!particle->GetPDG()) continue;
367 if (particle->GetPDG()->Charge() == 0.0) continue;
370 Bool_t prim = stack->IsPhysicalPrimary(iMc);
372 // only secondaries which can be reconstructed at TPC
375 //Float_t radius = TMath::Sqrt(particle->Vx()*particle->Vx()+particle->Vy()*particle->Vy()+particle->Vz()*particle->Vz());
376 //if(radius > fCutsMC->GetMaxR()) continue;
378 // only secondary electrons from gamma conversion
379 //if( TMath::Abs(particle->GetPdgCode())!=fCutsMC->GetEM() || particle->GetUniqueID() != 5) continue;
381 /*Bool_t findable = kFALSE;
382 for(Int_t iRec=0; iRec<multAll; ++iRec)
385 if(iMc == labelsAllRecSec[iRec])
387 findable = IsFindable(mcEvent,iMc);
391 Bool_t findable = IsFindable(mcEvent,iMc);
393 Bool_t recStatus = kFALSE;
394 Int_t nClones = 0, nFakes = 0;
395 for(Int_t iRec=0; iRec<multRec; ++iRec)
397 // check reconstructed
398 if(iMc == labelsRecSec[iRec])
400 if (recStatus && nClones < fgkMaxClones) nClones++;
403 //In order to relate the fake track to track parameters, we assign it to the best matching ESD track.
404 if (labelsRecSec[iRec] < 0 && -labelsRecSec[iRec] == iMc && nFakes < fgkMaxFakes) nFakes++;
407 // Only 5 charged particle species (e,mu,pi,K,p)
408 if (fCutsMC->IsPdgParticle(TMath::Abs(particle->GetPdgCode())) == kFALSE) continue;
410 // transform Pdg to Pid
411 Int_t pid = TransformToPID(particle);
413 Float_t mceta = particle->Eta();
414 Float_t mcphi = particle->Phi();
415 if(mcphi<0) mcphi += 2.*TMath::Pi();
416 Float_t mcpt = particle->Pt();
417 Float_t mcR = particle->R();
419 // get info about mother
420 Int_t motherLabel = particle->GetMother(0);
421 if(motherLabel < 0) continue;
422 TParticle *mother = stack->Particle(motherLabel);
423 if(!mother) continue;
425 Float_t mother_eta = mother->Eta();
426 Float_t mother_phi = mother->Phi();
427 if(mother_phi<0) mother_phi += 2.*TMath::Pi();
430 if (particle->GetPDG()->Charge() < 0) charge = -1.;
431 else if (particle->GetPDG()->Charge() > 0) charge = 1.;
434 Double_t vEffSecHisto[12] = { mceta, mcphi, mcpt, static_cast<Double_t>(pid), static_cast<Double_t>(recStatus), static_cast<Double_t>(findable), mcR, mother_phi, mother_eta, static_cast<Double_t>(charge), static_cast<Double_t>(nClones), static_cast<Double_t>(nFakes) };
435 fEffSecHisto->Fill(vEffSecHisto);
440 if(labelsRecSec) delete [] labelsRecSec; labelsRecSec = 0;
441 if(labelsAllRecSec) delete [] labelsAllRecSec; labelsAllRecSec = 0;
447 //_____________________________________________________________________________
448 void AliPerformanceEff::ProcessTPCITS(AliMCEvent* const mcEvent, AliESDEvent *const esdEvent)
450 // Fill efficiency comparison information
452 if(!esdEvent) return;
455 AliStack *stack = mcEvent->Stack();
457 AliDebug(AliLog::kError, "Stack not available");
461 Int_t *labelsRecTPCITS = NULL;
462 labelsRecTPCITS = new Int_t[esdEvent->GetNumberOfTracks()];
465 Printf("Cannot create labelsRecTPCITS");
468 for(Int_t i=0;i<esdEvent->GetNumberOfTracks();i++) { labelsRecTPCITS[i] = 0; }
470 Int_t *labelsAllRecTPCITS = NULL;
471 labelsAllRecTPCITS = new Int_t[esdEvent->GetNumberOfTracks()];
472 if(!labelsAllRecTPCITS) {
473 delete [] labelsRecTPCITS;
474 Printf("Cannot create labelsAllRecTPCITS");
477 for(Int_t i=0;i<esdEvent->GetNumberOfTracks();i++) { labelsAllRecTPCITS[i] = 0; }
479 // loop over rec. tracks
480 AliESDtrack *track=0;
481 for (Int_t iTrack = 0; iTrack < esdEvent->GetNumberOfTracks(); iTrack++)
483 track = esdEvent->GetTrack(iTrack);
485 if(track->Charge()==0) continue;
487 // if not fUseKinkDaughters don't use tracks with kink index > 0
488 if(!fUseKinkDaughters && track->GetKinkIndex(0) > 0) continue;
490 //Int_t label = TMath::Abs(track->GetLabel());
491 Int_t label = track->GetLabel(); //Use global label for combined efficiency analysis
492 labelsAllRecTPCITS[iTrack]=label;
495 if(IsRecTPCITS(track) != 0)
496 labelsRecTPCITS[iTrack]=label;
500 // MC histograms for efficiency studies
502 //Int_t nPart = stack->GetNtrack();
503 Int_t nPart = stack->GetNprimary();
504 for (Int_t iMc = 0; iMc < nPart; ++iMc)
506 if (iMc == 0) continue; //Cannot distinguish between track or fake track
507 TParticle* particle = stack->Particle(iMc);
508 if (!particle) continue;
509 if (!particle->GetPDG()) continue;
510 if (particle->GetPDG()->Charge() == 0.0) continue;
513 Bool_t prim = stack->IsPhysicalPrimary(iMc);
516 /*Bool_t findable = kFALSE;
517 for(Int_t iRec=0; iRec<esdEvent->GetNumberOfTracks(); ++iRec)
520 if(iMc == labelsAllRecTPCITS[iRec])
522 findable = IsFindable(mcEvent,iMc);
526 Bool_t findable = IsFindable(mcEvent,iMc);
528 Bool_t recStatus = kFALSE;
529 Int_t nClones = 0, nFakes = 0;
530 for(Int_t iRec=0; iRec<esdEvent->GetNumberOfTracks(); ++iRec)
532 // check reconstructed
533 if(iMc == labelsRecTPCITS[iRec])
535 if (recStatus && nClones < fgkMaxClones) nClones++;
538 //In order to relate the fake track to track parameters, we assign it to the best matching ESD track.
539 if (labelsRecTPCITS[iRec] < 0 && -labelsRecTPCITS[iRec] == iMc && nFakes < fgkMaxFakes) nFakes++;
542 // Only 5 charged particle species (e,mu,pi,K,p)
543 if (fCutsMC->IsPdgParticle(TMath::Abs(particle->GetPdgCode())) == kFALSE) continue;
545 // transform Pdg to Pid
546 Int_t pid = TransformToPID(particle);
548 Float_t mceta = particle->Eta();
549 Float_t mcphi = particle->Phi();
550 if(mcphi<0) mcphi += 2.*TMath::Pi();
551 Float_t mcpt = particle->Pt();
554 if (particle->GetPDG()->Charge() < 0) charge = -1.;
555 else if (particle->GetPDG()->Charge() > 0) charge = 1.;
558 Double_t vEffHisto[9] = { mceta, mcphi, mcpt, static_cast<Double_t>(pid), static_cast<Double_t>(recStatus), static_cast<Double_t>(findable), static_cast<Double_t>(charge), static_cast<Double_t>(nClones), static_cast<Double_t>(nFakes)};
559 fEffHisto->Fill(vEffHisto);
562 if(labelsRecTPCITS) delete [] labelsRecTPCITS; labelsRecTPCITS = 0;
563 if(labelsAllRecTPCITS) delete [] labelsAllRecTPCITS; labelsAllRecTPCITS = 0;
566 //_____________________________________________________________________________
567 void AliPerformanceEff::ProcessConstrained(AliMCEvent* const mcEvent, AliESDEvent *const esdEvent)
569 // Process comparison information
570 if(!esdEvent) return;
573 AliStack *stack = mcEvent->Stack();
575 AliDebug(AliLog::kError, "Stack not available");
579 Int_t *labelsRecConstrained = NULL;
580 labelsRecConstrained = new Int_t[esdEvent->GetNumberOfTracks()];
581 if(!labelsRecConstrained)
583 Printf("Cannot create labelsRecConstrained");
586 for(Int_t i=0;i<esdEvent->GetNumberOfTracks();i++) { labelsRecConstrained[i] = 0; }
588 Int_t *labelsAllRecConstrained = NULL;
589 labelsAllRecConstrained = new Int_t[esdEvent->GetNumberOfTracks()];
590 if(!labelsAllRecConstrained) {
591 delete [] labelsRecConstrained;
592 Printf("Cannot create labelsAllRecConstrained");
595 for(Int_t i=0;i<esdEvent->GetNumberOfTracks();i++) { labelsAllRecConstrained[i] = 0; }
597 // loop over rec. tracks
598 AliESDtrack *track=0;
599 for (Int_t iTrack = 0; iTrack < esdEvent->GetNumberOfTracks(); iTrack++)
601 track = esdEvent->GetTrack(iTrack);
603 if(track->Charge()==0) continue;
605 // if not fUseKinkDaughters don't use tracks with kink index > 0
606 if(!fUseKinkDaughters && track->GetKinkIndex(0) > 0) continue;
608 //Int_t label = TMath::Abs(track->GetLabel());
609 Int_t label = track->GetLabel();
610 labelsAllRecConstrained[iTrack]=label;
613 if(IsRecConstrained(track) != 0)
614 labelsRecConstrained[iTrack]=label;
619 // MC histograms for efficiency studies
623 //Int_t nPart = stack->GetNtrack();
624 Int_t nPart = stack->GetNprimary();
625 for (Int_t iMc = 0; iMc < nPart; ++iMc)
627 if (iMc == 0) continue; //Cannot distinguish between track or fake track
628 TParticle* particle = stack->Particle(iMc);
629 if (!particle) continue;
630 if (!particle->GetPDG()) continue;
631 if (particle->GetPDG()->Charge() == 0.0) continue;
634 Bool_t prim = stack->IsPhysicalPrimary(iMc);
637 /*Bool_t findable = kFALSE;
638 for(Int_t iRec=0; iRec<esdEvent->GetNumberOfTracks(); ++iRec)
641 if(iMc == labelsAllRecConstrained[iRec])
643 findable = IsFindable(mcEvent,iMc);
647 Bool_t findable = IsFindable(mcEvent,iMc);
649 Bool_t recStatus = kFALSE;
650 Int_t nClones = 0, nFakes = 0;
651 for(Int_t iRec=0; iRec<esdEvent->GetNumberOfTracks(); ++iRec)
653 // check reconstructed
654 if(iMc == labelsRecConstrained[iRec])
656 if (recStatus && nClones < fgkMaxClones) nClones++;
659 //In order to relate the fake track to track parameters, we assign it to the best matching ESD track.
660 if (labelsRecConstrained[iRec] < 0 && -labelsRecConstrained[iRec] == iMc && nFakes < fgkMaxFakes) nFakes++;
663 // Only 5 charged particle species (e,mu,pi,K,p)
664 if (fCutsMC->IsPdgParticle(TMath::Abs(particle->GetPdgCode())) == kFALSE) continue;
666 // transform Pdg to Pid
667 Int_t pid = TransformToPID(particle);
669 Float_t mceta = particle->Eta();
670 Float_t mcphi = particle->Phi();
671 if(mcphi<0) mcphi += 2.*TMath::Pi();
672 Float_t mcpt = particle->Pt();
675 if (particle->GetPDG()->Charge() < 0) charge = -1.;
676 else if (particle->GetPDG()->Charge() > 0) charge = 1.;
679 Double_t vEffHisto[9] = { mceta, mcphi, mcpt, static_cast<Double_t>(pid), static_cast<Double_t>(recStatus), static_cast<Double_t>(findable), static_cast<Double_t>(charge), static_cast<Double_t>(nClones), static_cast<Double_t>(nFakes) };
680 fEffHisto->Fill(vEffHisto);
683 if(labelsRecConstrained) delete [] labelsRecConstrained; labelsRecConstrained = 0;
684 if(labelsAllRecConstrained) delete [] labelsAllRecConstrained; labelsAllRecConstrained = 0;
687 //_____________________________________________________________________________
688 void AliPerformanceEff::Exec(AliMCEvent* const mcEvent, AliESDEvent *const esdEvent, AliESDfriend *const esdFriend, const Bool_t bUseMC, const Bool_t bUseESDfriend)
690 // Process comparison information
694 Error("Exec","esdEvent not available");
697 AliHeader* header = 0;
698 AliGenEventHeader* genHeader = 0;
705 Error("Exec","mcEvent not available");
708 // get MC event header
709 header = mcEvent->Header();
711 Error("Exec","Header not available");
715 stack = mcEvent->Stack();
717 Error("Exec","Stack not available");
721 genHeader = header->GenEventHeader();
723 Error("Exec","Could not retrieve genHeader from Header");
726 genHeader->PrimaryVertex(vtxMC);
729 Error("Exec","MC information required!");
736 Error("Exec","esdFriend not available");
744 if(GetAnalysisMode() == 0) ProcessTPC(mcEvent,esdEvent);
745 else if(GetAnalysisMode() == 1) ProcessTPCITS(mcEvent,esdEvent);
746 else if(GetAnalysisMode() == 2) ProcessConstrained(mcEvent,esdEvent);
747 else if(GetAnalysisMode() == 5) ProcessTPCSec(mcEvent,esdEvent);
749 printf("ERROR: AnalysisMode %d \n",fAnalysisMode);
754 //_____________________________________________________________________________
755 Int_t AliPerformanceEff::TransformToPID(TParticle *particle)
757 // transform Pdg to Pid
758 // Pdg convension is different for hadrons and leptons
759 // (e.g. K+/K- = 321/-321; e+/e- = -11/11 )
762 if( TMath::Abs(particle->GetPdgCode())==fCutsMC->GetEM() ) pid = 0;
763 if( TMath::Abs(particle->GetPdgCode())==fCutsMC->GetMuM() ) pid = 1;
764 if( TMath::Abs(particle->GetPdgCode())==fCutsMC->GetPiP() ) pid = 2;
765 if( TMath::Abs(particle->GetPdgCode())==fCutsMC->GetKP() ) pid = 3;
766 if( TMath::Abs(particle->GetPdgCode())==fCutsMC->GetProt() ) pid = 4;
771 //_____________________________________________________________________________
772 Bool_t AliPerformanceEff::IsFindable(const AliMCEvent *mcEvent, Int_t label)
775 // Findfindable tracks
777 if(!mcEvent) return kFALSE;
779 AliMCParticle *mcParticle = (AliMCParticle*) mcEvent->GetTrack(label);
780 if(!mcParticle) return kFALSE;
783 Float_t tpcTrackLength = mcParticle->GetTPCTrackLength(AliTracker::GetBz(),0.05,counter,3.0);
784 //printf("tpcTrackLength %f \n", tpcTrackLength);
786 return (tpcTrackLength>fCutsMC->GetMinTrackLength());
789 //_____________________________________________________________________________
790 Bool_t AliPerformanceEff::IsRecTPC(AliESDtrack *esdTrack)
793 // Check whether track is reconstructed in TPC
795 if(!esdTrack) return kFALSE;
797 const AliExternalTrackParam *track = esdTrack->GetTPCInnerParam();
798 if(!track) return kFALSE;
800 Float_t dca[2], cov[3]; // dca_xy, dca_z, sigma_xy, sigma_xy_z, sigma_z
801 esdTrack->GetImpactParametersTPC(dca,cov);
803 Bool_t recStatus = kFALSE;
804 if(esdTrack->GetTPCNcls()>fCutsRC->GetMinNClustersTPC()) recStatus = kTRUE;
807 if( TMath::Abs(dca[0])<fCutsRC->GetMaxDCAToVertexXY() &&
808 TMath::Abs(dca[1])<fCutsRC->GetMaxDCAToVertexZ())
817 //_____________________________________________________________________________
818 Bool_t AliPerformanceEff::IsRecTPCITS(AliESDtrack *esdTrack)
821 // Check whether track is reconstructed in TPCITS
823 if(!esdTrack) return kFALSE;
825 Float_t dca[2], cov[3]; // dca_xy, dca_z, sigma_xy, sigma_xy_z, sigma_z
826 esdTrack->GetImpactParameters(dca,cov);
828 Bool_t recStatus = kFALSE;
830 if ((esdTrack->GetStatus()&AliESDtrack::kTPCrefit)==0) return kFALSE; // TPC refit
831 if ((esdTrack->GetStatus()&AliESDtrack::kITSrefit)==0) return kFALSE; // ITS refit
832 if (esdTrack->GetITSclusters(0)<fCutsRC->GetMinNClustersITS()) return kFALSE; // min. nb. ITS clusters
833 //if ((esdTrack->GetStatus()&AliESDtrack::kITSrefit)==0) return kFALSE; // ITS refit
834 //Int_t clusterITS[200];
835 //if(esdTrack->GetITSclusters(clusterITS)<fCutsRC->GetMinNClustersITS()) return kFALSE; // min. nb. ITS clusters
839 if(TMath::Abs(dca[0])<fCutsRC->GetMaxDCAToVertexXY() &&
840 TMath::Abs(dca[1])<fCutsRC->GetMaxDCAToVertexZ())
849 //_____________________________________________________________________________
850 Bool_t AliPerformanceEff::IsRecConstrained(AliESDtrack *esdTrack)
853 // Check whether track is reconstructed in IsRecConstrained
855 if(!esdTrack) return kFALSE;
857 const AliExternalTrackParam * track = esdTrack->GetConstrainedParam();
858 if(!track) return kFALSE;
860 Float_t dca[2], cov[3]; // dca_xy, dca_z, sigma_xy, sigma_xy_z, sigma_z
861 esdTrack->GetImpactParameters(dca,cov);
862 //Int_t label = TMath::Abs(esdTrack->GetLabel());
864 Bool_t recStatus = kFALSE;
866 if ((esdTrack->GetStatus()&AliESDtrack::kTPCrefit)==0) return kFALSE; // TPC refit
867 if (esdTrack->GetTPCNcls()<fCutsRC->GetMinNClustersTPC()) return kFALSE; // min. nb. TPC clusters
868 Int_t clusterITS[200];
869 if(esdTrack->GetITSclusters(clusterITS)<fCutsRC->GetMinNClustersITS()) return kFALSE; // min. nb. ITS clusters
874 if(TMath::Abs(dca[0])<fCutsRC->GetMaxDCAToVertexXY() &&
875 TMath::Abs(dca[1])<fCutsRC->GetMaxDCAToVertexZ())
884 //_____________________________________________________________________________
885 Long64_t AliPerformanceEff::Merge(TCollection* const list)
887 // Merge list of objects (needed by PROOF)
895 TIterator* iter = list->MakeIterator();
898 // collection of generated histograms
901 while((obj = iter->Next()) != 0)
903 AliPerformanceEff* entry = dynamic_cast<AliPerformanceEff*>(obj);
904 if (entry == 0) continue;
906 fEffHisto->Add(entry->fEffHisto);
907 fEffSecHisto->Add(entry->fEffSecHisto);
914 //_____________________________________________________________________________
915 void AliPerformanceEff::Analyse()
917 // Analyse comparison information and store output histograms
918 // in the folder "folderEff"
920 TH1::AddDirectory(kFALSE);
921 TObjArray *aFolderObj = new TObjArray;
922 if(!aFolderObj) return;
929 if(GetAnalysisMode() != 5) {
931 fEffHisto->GetAxis(0)->SetRangeUser(-0.9,0.89); // eta range
932 fEffHisto->GetAxis(2)->SetRangeUser(0.1,19.99); // pt range // FIXME maybe remove since range is defined in THnSparse
934 // rec efficiency vs pt
935 fEffHisto->GetAxis(3)->SetRangeUser(0.,3.99); // reconstructed
937 fEffHisto->GetAxis(6)->SetRangeUser(-2.,2.); // charge all
938 aFolderObj->Add(AddHistoEff(2, "ptRecEff", "rec. efficiency", 0));
939 aFolderObj->Add(AddHistoEff(2, "ptClone", "clone rate", 1));
940 aFolderObj->Add(AddHistoEff(2, "ptFake", "fake rate", 2));
942 fEffHisto->GetAxis(6)->SetRangeUser(-2.,0.); // charge negativ
943 aFolderObj->Add(AddHistoEff(2, "ptRecEffNeg", "rec. efficiency neg.", 0));
945 fEffHisto->GetAxis(6)->SetRangeUser(0.,2.); // charge positiv
946 aFolderObj->Add(AddHistoEff(2, "ptRecEffPos", "rec. efficiency pos.", 0));
948 // rec efficiency vs pid vs pt
949 fEffHisto->GetAxis(3)->SetRange(3,3); // pions
951 fEffHisto->GetAxis(6)->SetRangeUser(-2.,2.); // charge all
952 aFolderObj->Add(AddHistoEff(2, "ptRecEffPi", "rec. efficiency (pions)", 0));
954 fEffHisto->GetAxis(6)->SetRangeUser(-2.,0.); // charge negativ
955 aFolderObj->Add(AddHistoEff(2, "ptRecEffPiNeg", "rec. efficiency (pions) neg.", 0));
957 fEffHisto->GetAxis(6)->SetRangeUser(0.,2.); // charge positiv
958 aFolderObj->Add(AddHistoEff(2, "ptRecEffPiPos", "rec. efficiency (pions) pos.", 0));
961 fEffHisto->GetAxis(3)->SetRange(4,4); // kaons
963 fEffHisto->GetAxis(6)->SetRangeUser(-2.,2.); // charge all
964 aFolderObj->Add(AddHistoEff(2, "ptRecEffK", "rec. efficiency (kaons)", 0));
966 fEffHisto->GetAxis(6)->SetRangeUser(-2.,0.); // charge negativ
967 aFolderObj->Add(AddHistoEff(2, "ptRecEffKNeg", "rec. efficiency (kaons) neg.", 0));
969 fEffHisto->GetAxis(6)->SetRangeUser(0.,2.); // charge positiv
970 aFolderObj->Add(AddHistoEff(2, "ptRecEffKPos", "rec. efficiency (kaons) pos.", 0));
973 fEffHisto->GetAxis(3)->SetRange(5,5); // protons
975 fEffHisto->GetAxis(6)->SetRangeUser(-2.,2.); // charge all
976 aFolderObj->Add(AddHistoEff(2, "ptRecEffP", "rec. efficiency (protons)", 0));
978 fEffHisto->GetAxis(6)->SetRangeUser(-2.,0.); // charge negativ
979 aFolderObj->Add(AddHistoEff(2, "ptRecEffPNeg", "rec. efficiency (protons) neg.", 0));
981 fEffHisto->GetAxis(6)->SetRangeUser(0.,2.); // charge positiv
982 aFolderObj->Add(AddHistoEff(2, "ptRecEffPPos", "rec. efficiency (protons) pos.", 0));
984 // findable efficiency vs pt
985 fEffHisto->GetAxis(3)->SetRangeUser(0.,4.);
986 fEffHisto->GetAxis(5)->SetRange(2,2); // findable
988 fEffHisto->GetAxis(6)->SetRangeUser(-2.,2.); // charge all
989 aFolderObj->Add(AddHistoEff(2, "ptRecEffF", "rec. efficiency (findable)", 0));
990 aFolderObj->Add(AddHistoEff(2, "ptCloneF", "clone rate (findable)", 1));
991 aFolderObj->Add(AddHistoEff(2, "ptFakeF", "fake rate (findable)", 2));
994 fEffHisto->GetAxis(6)->SetRangeUser(-2.,0.); // charge negativ
995 aFolderObj->Add(AddHistoEff(2, "ptRecEffFNeg", "rec. efficiency (findable) neg.", 0));
997 fEffHisto->GetAxis(6)->SetRangeUser(0.,2.); // charge positiv
998 aFolderObj->Add(AddHistoEff(2, "ptRecEffFPos", "rec. efficiency (findable) pos.", 0));
1001 // efficiency vs eta
1004 fEffHisto->GetAxis(0)->SetRangeUser(-1.5,1.49); // eta range
1005 fEffHisto->GetAxis(2)->SetRangeUser(0.1,19.99); // pt range
1006 fEffHisto->GetAxis(5)->SetRangeUser(0.,1.0); // all
1008 // rec efficiency vs eta
1009 fEffHisto->GetAxis(3)->SetRangeUser(0.,4.); // reconstructed
1011 fEffHisto->GetAxis(6)->SetRangeUser(-2.,2.); // charge all
1012 aFolderObj->Add(AddHistoEff(0, "etaRecEff", "rec. efficiency", 0));
1013 aFolderObj->Add(AddHistoEff(0, "etaClone", "clone rate", 1));
1014 aFolderObj->Add(AddHistoEff(0, "etaFake", "fake rate", 2));
1017 fEffHisto->GetAxis(6)->SetRangeUser(-2.,0.); // charge negativ
1018 aFolderObj->Add(AddHistoEff(0, "etaRecEffNeg", "rec. efficiency neg.", 0));
1020 fEffHisto->GetAxis(6)->SetRangeUser(0.,2.); // charge positiv
1021 aFolderObj->Add(AddHistoEff(0, "etaRecEffPos", "rec. efficiency pos.", 0));
1023 // rec efficiency vs pid vs eta
1024 fEffHisto->GetAxis(3)->SetRange(3,3); // pions
1026 fEffHisto->GetAxis(6)->SetRangeUser(-2.,2.); // charge all
1027 aFolderObj->Add(AddHistoEff(0, "etaRecEffPi", "rec. efficiency (pions)", 0));
1029 fEffHisto->GetAxis(6)->SetRangeUser(-2.,0.); // charge negativ
1030 aFolderObj->Add(AddHistoEff(0, "etaRecEffPiNeg", "rec. efficiency (pions) neg.", 0));
1032 fEffHisto->GetAxis(6)->SetRangeUser(0.,2.); // charge positiv
1033 aFolderObj->Add(AddHistoEff(0, "etaRecEffPiPos", "rec. efficiency (pions) pos.", 0));
1036 fEffHisto->GetAxis(3)->SetRange(4,4); // kaons
1038 fEffHisto->GetAxis(6)->SetRangeUser(-2.,2.); // charge all
1039 aFolderObj->Add(AddHistoEff(0, "etaRecEffK", "rec. efficiency (kaons)", 0));
1041 fEffHisto->GetAxis(6)->SetRangeUser(-2.,0.); // charge negativ
1042 aFolderObj->Add(AddHistoEff(0, "etaRecEffKNeg", "rec. efficiency (kaons) neg.", 0));
1044 fEffHisto->GetAxis(6)->SetRangeUser(0.,2.); // charge positiv
1045 aFolderObj->Add(AddHistoEff(0, "etaRecEffKPos", "rec. efficiency (kaons) pos.", 0));
1048 fEffHisto->GetAxis(3)->SetRange(5,5); // protons
1050 fEffHisto->GetAxis(6)->SetRangeUser(-2.,2.); // charge all
1051 aFolderObj->Add(AddHistoEff(0, "etaRecEffP", "rec. efficiency (protons)", 0));
1053 fEffHisto->GetAxis(6)->SetRangeUser(-2.,0.); // charge negativ
1054 aFolderObj->Add(AddHistoEff(0, "etaRecEffPNeg", "rec. efficiency (protons) neg.", 0));
1056 fEffHisto->GetAxis(6)->SetRangeUser(0.,2.); // charge positiv
1057 aFolderObj->Add(AddHistoEff(0, "etaRecEffPPos", "rec. efficiency (protons) pos.", 0));
1059 // findable efficiency vs eta
1060 fEffHisto->GetAxis(3)->SetRangeUser(0.,4.);
1061 fEffHisto->GetAxis(5)->SetRange(2,2); // findable
1063 fEffHisto->GetAxis(6)->SetRangeUser(-2.,2.); // charge all
1064 aFolderObj->Add(AddHistoEff(0, "etaRecEffF", "rec. efficiency (findable)", 0));
1065 aFolderObj->Add(AddHistoEff(0, "etaCloneF", "clone rate (findable)", 1));
1066 aFolderObj->Add(AddHistoEff(0, "etaFakeF", "fake rate (findable)", 2));
1069 fEffHisto->GetAxis(6)->SetRangeUser(-2.,0.); // charge negativ
1070 aFolderObj->Add(AddHistoEff(0, "etaRecEffFNeg", "rec. efficiency (findable) neg.", 0));
1072 fEffHisto->GetAxis(6)->SetRangeUser(0.,2.); // charge positiv
1073 aFolderObj->Add(AddHistoEff(0, "etaRecEffFPos", "rec. efficiency (findable) pos.", 0));
1076 // efficiency vs phi
1079 fEffHisto->GetAxis(0)->SetRangeUser(-0.9,0.89); // eta range
1080 fEffHisto->GetAxis(2)->SetRangeUser(0.1,19.99); // pt range
1081 fEffHisto->GetAxis(5)->SetRangeUser(0.,1.); // all
1083 // rec efficiency vs phi
1084 fEffHisto->GetAxis(3)->SetRangeUser(0.,4.); // reconstructed
1086 fEffHisto->GetAxis(6)->SetRangeUser(-2.,2.); // charge all
1087 aFolderObj->Add(AddHistoEff(1, "phiRecEff", "rec. efficiency", 0));
1088 aFolderObj->Add(AddHistoEff(1, "phiClone", "clone rate", 1));
1089 aFolderObj->Add(AddHistoEff(1, "phiFake", "fake rate", 2));
1091 fEffHisto->GetAxis(6)->SetRangeUser(-2.,0.); // charge negativ
1092 aFolderObj->Add(AddHistoEff(1, "phiRecEffNeg", "rec. efficiency neg.", 0));
1094 fEffHisto->GetAxis(6)->SetRangeUser(0.,2.); // charge positiv
1095 aFolderObj->Add(AddHistoEff(1, "phiRecEffPos", "rec. efficiency pos.", 0));
1097 // rec efficiency vs pid vs phi
1098 fEffHisto->GetAxis(3)->SetRange(3,3); // pions
1100 fEffHisto->GetAxis(6)->SetRangeUser(-2.,2.); // charge all
1101 aFolderObj->Add(AddHistoEff(1, "phiRecEffPi", "rec. efficiency (pions)", 0));
1103 fEffHisto->GetAxis(6)->SetRangeUser(-2.,0.); // charge negativ
1104 aFolderObj->Add(AddHistoEff(1, "phiRecEffPiNeg", "rec. efficiency (pions) neg.", 0));
1106 fEffHisto->GetAxis(6)->SetRangeUser(0.,2.); // charge positiv
1107 aFolderObj->Add(AddHistoEff(1, "phiRecEffPiPos", "rec. efficiency (pions) pos.", 0));
1110 fEffHisto->GetAxis(3)->SetRange(4,4); // kaons
1112 fEffHisto->GetAxis(6)->SetRangeUser(-2.,2.); // charge all
1113 aFolderObj->Add(AddHistoEff(1, "phiRecEffK", "rec. efficiency (kaons)", 0));
1115 fEffHisto->GetAxis(6)->SetRangeUser(-2.,0.); // charge negativ
1116 aFolderObj->Add(AddHistoEff(1, "phiRecEffKNeg", "rec. efficiency (kaons) neg.", 0));
1118 fEffHisto->GetAxis(6)->SetRangeUser(0.,2.); // charge positiv
1119 aFolderObj->Add(AddHistoEff(1, "phiRecEffKPos", "rec. efficiency (kaons) pos.", 0));
1122 fEffHisto->GetAxis(3)->SetRange(5,5); // protons
1124 fEffHisto->GetAxis(6)->SetRangeUser(-2.,2.); // charge all
1125 aFolderObj->Add(AddHistoEff(1, "phiRecEffP", "rec. efficiency (protons)", 0));
1127 fEffHisto->GetAxis(6)->SetRangeUser(-2.,0.); // charge negativ
1128 aFolderObj->Add(AddHistoEff(1, "phiRecEffPNeg", "rec. efficiency (protons) neg.", 0));
1130 fEffHisto->GetAxis(6)->SetRangeUser(0.,2.); // charge positiv
1131 aFolderObj->Add(AddHistoEff(1, "phiRecEffPPos", "rec. efficiency (protons) pos.", 0));
1133 // findable efficiency vs phi
1134 fEffHisto->GetAxis(3)->SetRangeUser(0.,4.);
1135 fEffHisto->GetAxis(5)->SetRange(2,2); // findable
1137 fEffHisto->GetAxis(6)->SetRangeUser(-2.,2.); // charge all
1138 aFolderObj->Add(AddHistoEff(1, "phiRecEffF", "rec. efficiency (findable)", 0));
1139 aFolderObj->Add(AddHistoEff(1, "phiCloneF", "clone rate (findable)", 1));
1140 aFolderObj->Add(AddHistoEff(1, "phiFakeF", "fake rate (findable)", 2));
1143 fEffHisto->GetAxis(6)->SetRangeUser(-2.,0.); // charge negativ
1144 aFolderObj->Add(AddHistoEff(1, "phiRecEffFNeg", "rec. efficiency (findable) neg.", 0));
1146 fEffHisto->GetAxis(6)->SetRangeUser(0.,2.); // charge positiv
1147 aFolderObj->Add(AddHistoEff(1, "phiRecEffFPos", "rec. efficiency (findable) pos.", 0));
1151 Float_t minEta=-1.5, maxEta=1.49;
1152 Float_t minR=0.0, maxR=150.0;
1153 Float_t minPt=0.10, maxPt=19.99;
1156 fEffSecHisto->GetAxis(8)->SetRangeUser(minEta,maxEta);
1158 // particle creation radius range
1159 fEffSecHisto->GetAxis(6)->SetRangeUser(minR,maxR);
1162 fEffSecHisto->GetAxis(0)->SetRangeUser(minEta,maxEta);
1163 fEffSecHisto->GetAxis(2)->SetRangeUser(minPt,maxPt);
1165 // rec efficiency vs pt
1167 aFolderObj->Add(AddHistoEff(2, "ptRecEff", "rec. efficiency", 0, 1));
1168 aFolderObj->Add(AddHistoEff(2, "ptClone", "clone rate", 1, 1));
1169 aFolderObj->Add(AddHistoEff(2, "ptFake", "fake rate", 2, 1));
1171 // rec efficiency vs pid vs pt
1173 fEffSecHisto->GetAxis(3)->SetRange(1,1); // electrons
1174 aFolderObj->Add(AddHistoEff(2, "ptRecEffEle", "rec. efficiency (electrons)", 0, 1));
1176 fEffSecHisto->GetAxis(3)->SetRange(3,3); // pions
1177 aFolderObj->Add(AddHistoEff(2, "ptRecEffPi", "rec. efficiency (pions)", 0, 1));
1179 fEffSecHisto->GetAxis(3)->SetRange(4,4); // kaons
1180 aFolderObj->Add(AddHistoEff(2, "ptRecEffK", "rec. efficiency (kaons)", 0, 1));
1182 fEffSecHisto->GetAxis(3)->SetRange(5,5); // protons
1183 aFolderObj->Add(AddHistoEff(2, "ptRecEffP", "rec. efficiency (protons)", 0, 1));
1185 fEffSecHisto->GetAxis(3)->SetRangeUser(0.,4.);
1187 // findable efficiency vs pt
1189 fEffSecHisto->GetAxis(5)->SetRange(2,2); // findable
1190 aFolderObj->Add(AddHistoEff(2, "ptRecEffF", "rec. efficiency (findable)", 0, 1));
1191 aFolderObj->Add(AddHistoEff(2, "ptCloneF", "clone rate (findable)", 1, 1));
1192 aFolderObj->Add(AddHistoEff(2, "ptFakeF", "fake rate (findable)", 2, 1));
1195 // efficiency vs eta
1197 fEffSecHisto->GetAxis(2)->SetRangeUser(minPt,maxPt);
1198 fEffSecHisto->GetAxis(4)->SetRangeUser(0.,1.); // all
1199 fEffSecHisto->GetAxis(5)->SetRangeUser(0.,1.); // all
1201 aFolderObj->Add(AddHistoEff(0, "etaRecEff", "rec. efficiency", 0, 1));
1202 aFolderObj->Add(AddHistoEff(0, "etaClone", "clone rate", 1, 1));
1203 aFolderObj->Add(AddHistoEff(0, "etaFake", "fake rate", 2, 1));
1205 // rec efficiency vs pid vs eta
1206 fEffSecHisto->GetAxis(3)->SetRange(1,1); // electrons
1207 aFolderObj->Add(AddHistoEff(0, "etaRecEffEle", "rec. efficiency (electrons)", 0, 1));
1209 fEffSecHisto->GetAxis(3)->SetRange(3,3); // pions
1210 aFolderObj->Add(AddHistoEff(0, "etaRecEffPi", "rec. efficiency (pions)", 0, 1));
1212 fEffSecHisto->GetAxis(3)->SetRange(4,4); // kaons
1213 aFolderObj->Add(AddHistoEff(0, "etaRecEffK", "rec. efficiency (kaons)", 0, 1));
1215 fEffSecHisto->GetAxis(3)->SetRange(5,5); // protons
1216 aFolderObj->Add(AddHistoEff(0, "etaRecEffP", "rec. efficiency (protons)", 0, 1));
1218 fEffSecHisto->GetAxis(3)->SetRangeUser(0.,4.);
1220 // findable efficiency vs eta
1222 fEffSecHisto->GetAxis(5)->SetRange(2,2); // findable
1223 aFolderObj->Add(AddHistoEff(0, "etaRecEffF", "rec. efficiency (findable)", 0, 1));
1224 aFolderObj->Add(AddHistoEff(0, "etaCloneF", "clone rate (findable)", 1, 1));
1225 aFolderObj->Add(AddHistoEff(0, "etaFakeF", "fake rate (findable)", 2, 1));
1228 // efficiency vs phi
1231 fEffSecHisto->GetAxis(0)->SetRangeUser(minEta,maxEta);
1232 fEffSecHisto->GetAxis(2)->SetRangeUser(minPt,maxPt);
1234 fEffSecHisto->GetAxis(4)->SetRangeUser(0.,1.); // all
1235 fEffSecHisto->GetAxis(5)->SetRangeUser(0.,1.); // all
1237 aFolderObj->Add(AddHistoEff(1, "phiRecEff", "rec. efficiency", 0, 1));
1238 aFolderObj->Add(AddHistoEff(1, "phiClone", "clone rate", 1, 1));
1239 aFolderObj->Add(AddHistoEff(1, "phiFake", "fake rate", 2, 1));
1241 // rec efficiency vs pid vs phi
1242 fEffSecHisto->GetAxis(3)->SetRange(1,1);
1243 aFolderObj->Add(AddHistoEff(1, "phiRecEffEle", "rec. efficiency (electrons)", 0, 1));
1245 fEffSecHisto->GetAxis(3)->SetRange(3,3); // pions
1246 aFolderObj->Add(AddHistoEff(1, "phiRecEffPi", "rec. efficiency (pions)", 0, 1));
1248 fEffSecHisto->GetAxis(3)->SetRange(4,4); // kaons
1249 aFolderObj->Add(AddHistoEff(1, "phiRecEffK", "rec. efficiency (kaons)", 0, 1));
1251 fEffSecHisto->GetAxis(3)->SetRange(5,5); // protons
1252 aFolderObj->Add(AddHistoEff(1, "phiRecEffP", "rec. efficiency (protons)", 0, 1));
1254 fEffSecHisto->GetAxis(3)->SetRangeUser(0.,4.);
1256 // findable efficiency vs phi
1258 fEffSecHisto->GetAxis(5)->SetRange(2,2); // findable
1259 aFolderObj->Add(AddHistoEff(1, "phiRecEffF", "rec. efficiency (findable)", 0, 1));
1260 aFolderObj->Add(AddHistoEff(1, "phiCloneF", "clone rate (findable)", 1, 1));
1261 aFolderObj->Add(AddHistoEff(1, "phiFakeF", "fake rate (findable)", 2, 1));
1264 for (Int_t i = 0;i < fEffHisto->GetNdimensions();i++)
1266 fEffHisto->GetAxis(i)->SetRange(1,0); //Reset Range
1269 // export objects to analysis folder
1270 fAnalysisFolder = ExportToFolder(aFolderObj);
1272 // delete only TObjArray
1273 if(aFolderObj) delete aFolderObj;
1276 //_____________________________________________________________________________
1277 TFolder* AliPerformanceEff::ExportToFolder(TObjArray * array)
1279 // recreate folder avery time and export objects to new one
1281 AliPerformanceEff * comp=this;
1282 TFolder *folder = comp->GetAnalysisFolder();
1284 TString name, title;
1285 TFolder *newFolder = 0;
1287 Int_t size = array->GetSize();
1290 // get name and title from old folder
1291 name = folder->GetName();
1292 title = folder->GetTitle();
1298 newFolder = CreateFolder(name.Data(),title.Data());
1299 newFolder->SetOwner();
1301 // add objects to folder
1303 newFolder->Add(array->At(i));
1312 //_____________________________________________________________________________
1313 TFolder* AliPerformanceEff::CreateFolder(TString name,TString title) {
1314 // create folder for analysed histograms
1316 TFolder *folder = 0;
1317 folder = new TFolder(name.Data(),title.Data());
1322 TH1D* WeightedProjection(THnSparseF* src, Int_t axis, Int_t nWeights, Int_t* weightCoords)
1324 THnSparseF* tmp = (THnSparseF*) src->Clone();
1326 for (i = 0;i < tmp->GetNbins();i++)
1329 tmp->GetBinContent(i, coords);
1330 Int_t weight = 0, j;
1331 for (j = 0;j < nWeights;j++)
1333 //The coordinates range from 1 to maxClones / maxFakes + 1, so we have to subtract one
1334 weight += coords[weightCoords[j]] - 1;
1336 tmp->SetBinContent(i, weight);
1339 TH1D* ret = tmp->Projection(axis);
1344 //_____________________________________________________________________________
1345 TH1D* AliPerformanceEff::AddHistoEff(Int_t axis, const Char_t *name, const Char_t* vsTitle,
1346 const Int_t type, const Int_t secondary) {
1347 // Create and add rec efficiency vs pt, eta, phi
1353 THnSparseF* EffHisto = secondary ? fEffSecHisto : fEffHisto;
1355 Int_t axis_clone = secondary ? 10 : 7;
1356 Int_t axis_fake = secondary ? 11 : 8;
1357 Int_t axis_all[3] = {4, axis_clone, axis_fake};
1360 if (type == 0) // Efficiency
1362 EffHisto->GetAxis(4)->SetRange(1.,2.); // all
1363 TH1D *all = EffHisto->Projection(axis);
1365 EffHisto->GetAxis(4)->SetRange(2.,2.); // reconstructed
1366 TH1D *rec = EffHisto->Projection(axis);
1367 recc = (TH1D*)rec->Clone();
1371 recc->Divide(rec,all,1,1,"B");
1372 recc->GetYaxis()->SetTitle("efficiency");
1375 else if (type == 1) // Clone Rate
1377 EffHisto->GetAxis(4)->SetRange(2.,2.); // reconstructed
1378 TH1D *all = WeightedProjection(EffHisto, axis, 3, axis_all);
1380 EffHisto->GetAxis(4)->SetRange(2.,2.); // reconstructed
1381 TH1D *clone = WeightedProjection(EffHisto, axis, 1, &axis_clone);
1382 recc = (TH1D*) clone->Clone();
1386 recc->Divide(clone,all,1,1,"B");
1387 recc->GetYaxis()->SetTitle("clone rate");
1390 else if (type == 2) // Fake Rate
1392 EffHisto->GetAxis(4)->SetRange(2.,2.); // reconstructed
1393 TH1D *all = WeightedProjection(EffHisto, axis, 3, axis_all);
1395 EffHisto->GetAxis(4)->SetRange(2.,2.); // reconstructed
1396 TH1D *fake = WeightedProjection(EffHisto, axis, 1, &axis_fake);
1397 recc = (TH1D*) fake->Clone();
1401 recc->Divide(fake,all,1,1,"B");
1402 recc->GetYaxis()->SetTitle("fake rate");
1406 EffHisto->GetAxis(4)->SetRange(1,0); //Reset Range
1408 if (recc) { // Coverity fix
1409 recc->SetName(name);
1411 recc->GetXaxis()->SetTitle(fEffHisto->GetAxis(axis)->GetTitle());
1413 snprintf(title,256,"%s vs %s",vsTitle, fEffHisto->GetAxis(axis)->GetTitle());
1414 recc->SetTitle(title);
1416 if (axis == 2 ) recc->SetBit(TH1::kLogX);