+++ /dev/null
-//------------------------------------------------------------------------------
-// Implementation of AliPerformanceEff class. It keeps information from
-// comparison of reconstructed and MC particle tracks. In addtion,
-// it keeps selection cuts used during comparison. The comparison
-// information is stored in the ROOT histograms. Analysis of these
-// histograms can be done by using Analyse() class function. The result of
-// the analysis (histograms/graphs) are stored in the folder which is
-// a data member of AliPerformanceEff.
-//
-// Author: J.Otwinowski 04/02/2008
-//------------------------------------------------------------------------------
-
-/*
-
- // after running comparison task, read the file, and get component
- gROOT->LoadMacro("$ALICE_ROOT/PWGPP/Macros/LoadMyLibs.C");
- LoadMyLibs();
- TFile f("Output.root");
- AliPerformanceEff * compObj = (AliPerformanceEff*)coutput->FindObject("AliPerformanceEff");
-
- // Analyse comparison data
- compObj->Analyse();
-
- // the output histograms/graphs will be stored in the folder "folderEff"
- compObj->GetAnalysisFolder()->ls("*");
-
- // user can save whole comparison object (or only folder with anlysed histograms)
- // in the seperate output file (e.g.)
- TFile fout("Analysed_Eff.root","recreate");
- compObj->Write(); // compObj->GetAnalysisFolder()->Write();
- fout.Close();
-
-*/
-
-#include <TAxis.h>
-#include <TH1D.h>
-#include "THnSparse.h"
-
-//
-#include "AliESDtrack.h"
-#include "AliRecInfoCuts.h"
-#include "AliMCInfoCuts.h"
-#include "AliLog.h"
-#include "AliESDVertex.h"
-#include "AliExternalTrackParam.h"
-#include "AliTracker.h"
-#include "AliESDEvent.h"
-#include "AliMCEvent.h"
-#include "AliMCParticle.h"
-#include "AliHeader.h"
-#include "AliGenEventHeader.h"
-#include "AliStack.h"
-#include "AliPerformanceEff.h"
-
-using namespace std;
-
-ClassImp(AliPerformanceEff)
-
-//_____________________________________________________________________________
-AliPerformanceEff::AliPerformanceEff(const Char_t* name, const Char_t* title, Int_t analysisMode, Bool_t hptGenerator): AliPerformanceObject(name,title),
-
- // histograms
- fEffHisto(0),
- fEffSecHisto(0),
-
- // Cuts
- fCutsRC(0),
- fCutsMC(0),
-
- // histogram folder
- fAnalysisFolder(0)
-{
- // named constructor
- //
- SetAnalysisMode(analysisMode);
- SetHptGenerator(hptGenerator);
-
- Init();
-}
-
-
-//_____________________________________________________________________________
-AliPerformanceEff::~AliPerformanceEff()
-{
-// destructor
-
- if(fEffHisto) delete fEffHisto; fEffHisto=0;
- if(fEffSecHisto) delete fEffSecHisto; fEffSecHisto=0;
- if(fAnalysisFolder) delete fAnalysisFolder; fAnalysisFolder=0;
-}
-
-//_____________________________________________________________________________
-void AliPerformanceEff::Init()
-{
- // Init histograms
- //
-
- // set pt bins
- Int_t nPtBins = 50;
- Double_t ptMin = 1.e-2, ptMax = 20.;
-
- Double_t *binsPt = 0;
-
- if (IsHptGenerator()) {
- ptMax = 100.;
- }
- binsPt = CreateLogAxis(nPtBins,ptMin,ptMax);
-
- /*
- Int_t nPtBins = 31;
- 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.};
- Double_t ptMin = 0., ptMax = 10.;
-
- if(IsHptGenerator() == kTRUE) {
- nPtBins = 100;
- ptMin = 0.; ptMax = 100.;
- }
- */
-
- //mceta:mcphi:mcpt:pid:recStatus:findable:charge
- Int_t binsEffHisto[9]={30,144,nPtBins,5,2,2,3,fgkMaxClones+1,fgkMaxFakes+1};
- Double_t minEffHisto[9]={-1.5,0.,ptMin,0.,0.,0.,-1.5,0,0};
- Double_t maxEffHisto[9]={ 1.5,2.*TMath::Pi(), ptMax,5.,2.,2.,1.5,fgkMaxClones+1,fgkMaxFakes+1};
-
- fEffHisto = new THnSparseF("fEffHisto","mceta:mcphi:mcpt:pid:recStatus:findable:charge:nclones:nfakes",9,binsEffHisto,minEffHisto,maxEffHisto);
- fEffHisto->SetBinEdges(2,binsPt);
-
- fEffHisto->GetAxis(0)->SetTitle("#eta_{mc}");
- fEffHisto->GetAxis(1)->SetTitle("#phi_{mc} (rad)");
- fEffHisto->GetAxis(2)->SetTitle("p_{Tmc} (GeV/c)");
- fEffHisto->GetAxis(3)->SetTitle("pid");
- fEffHisto->GetAxis(4)->SetTitle("recStatus");
- fEffHisto->GetAxis(5)->SetTitle("findable");
- fEffHisto->GetAxis(6)->SetTitle("charge");
- fEffHisto->GetAxis(7)->SetTitle("nClones");
- fEffHisto->GetAxis(8)->SetTitle("nFakes");
- fEffHisto->Sumw2();
-
- //mceta:mcphi:mcpt:pid:recStatus:findable:mcR:mother_phi:mother_eta:charge
- Int_t binsEffSecHisto[12]={30,60,nPtBins,5,2,2,100,60,30,3,fgkMaxClones+1,fgkMaxFakes+1};
- Double_t minEffSecHisto[12]={-1.5,0.,ptMin,0.,0.,0.,0.,0.,-1.5,-1.5,0,0};
- 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};
-
- fEffSecHisto = new THnSparseF("fEffSecHisto","mceta:mcphi:mcpt:pid:recStatus:findable:mcR:mother_phi:mother_eta:charge:nclones:nfakes",12,binsEffSecHisto,minEffSecHisto,maxEffSecHisto);
- fEffSecHisto->SetBinEdges(2,binsPt);
-
- fEffSecHisto->GetAxis(0)->SetTitle("#eta_{mc}");
- fEffSecHisto->GetAxis(1)->SetTitle("#phi_{mc} (rad)");
- fEffSecHisto->GetAxis(2)->SetTitle("p_{Tmc} (GeV/c)");
- fEffSecHisto->GetAxis(3)->SetTitle("pid");
- fEffSecHisto->GetAxis(4)->SetTitle("recStatus");
- fEffSecHisto->GetAxis(5)->SetTitle("findable");
- fEffSecHisto->GetAxis(6)->SetTitle("mcR (cm)");
- fEffSecHisto->GetAxis(7)->SetTitle("mother_phi (rad)");
- fEffSecHisto->GetAxis(8)->SetTitle("mother_eta");
- fEffSecHisto->GetAxis(9)->SetTitle("charge");
- fEffSecHisto->GetAxis(10)->SetTitle("nClones");
- fEffSecHisto->GetAxis(11)->SetTitle("nFakes");
- fEffSecHisto->Sumw2();
-
- // init cuts
- if(!fCutsMC)
- AliDebug(AliLog::kError, "ERROR: Cannot find AliMCInfoCuts object");
- if(!fCutsRC)
- AliDebug(AliLog::kError, "ERROR: Cannot find AliRecInfoCuts object");
-
- // init folder
- fAnalysisFolder = CreateFolder("folderEff","Analysis Efficiency Folder");
-}
-
-//_____________________________________________________________________________
-void AliPerformanceEff::ProcessTPC(AliMCEvent* const mcEvent, AliESDEvent *const esdEvent)
-{
- // Fill TPC only efficiency comparison information
- if(!esdEvent) return;
- if(!mcEvent) return;
-
- AliStack *stack = mcEvent->Stack();
- if (!stack) {
- AliDebug(AliLog::kError, "Stack not available");
- return;
- }
-
- Int_t *labelsRec = NULL;
- labelsRec = new Int_t[esdEvent->GetNumberOfTracks()];
- if(!labelsRec)
- {
- Printf("Cannot create labelsRec");
- return;
- }
- for(Int_t i=0;i<esdEvent->GetNumberOfTracks();i++) { labelsRec[i] = 0; }
-
- Int_t *labelsAllRec = NULL;
- labelsAllRec = new Int_t[esdEvent->GetNumberOfTracks()];
- if(!labelsAllRec) {
- delete [] labelsRec;
- Printf("Cannot create labelsAllRec");
- return;
- }
- for(Int_t i=0;i<esdEvent->GetNumberOfTracks();i++) { labelsAllRec[i] = 0; }
-
- // loop over rec. tracks
- AliESDtrack *track=0;
- for (Int_t iTrack = 0; iTrack < esdEvent->GetNumberOfTracks(); iTrack++)
- {
- track = esdEvent->GetTrack(iTrack);
- if(!track) continue;
- if(track->Charge()==0) continue;
-
- // if not fUseKinkDaughters don't use tracks with kink index > 0
- if(!fUseKinkDaughters && track->GetKinkIndex(0) > 0) continue;
-
- //Int_t label = TMath::Abs(track->GetLabel());
- Int_t label = track->GetTPCLabel(); //Use TPC-only label for TPC-only efficiency analysis
- labelsAllRec[iTrack]=label;
-
- // TPC only
- if(IsRecTPC(track) != 0) {
- labelsRec[iTrack]=label;
- }
-
- }
-
- //
- // MC histograms for efficiency studies
- //
- Int_t nPart = stack->GetNtrack();
- //Int_t nPart = stack->GetNprimary();
- for (Int_t iMc = 0; iMc < nPart; ++iMc)
- {
- if (iMc == 0) continue; //Cannot distinguish between track or fake track
- TParticle* particle = stack->Particle(iMc);
- if (!particle) continue;
- if (!particle->GetPDG()) continue;
- if (particle->GetPDG()->Charge() == 0.0) continue;
-
- // physical primary
- Bool_t prim = stack->IsPhysicalPrimary(iMc);
- if(!prim) continue;
-
- // --- check for double filling in stack
- // use only particles with no daughters in the list of primaries
- Int_t nDaughters = 0;// particle->GetNDaughters();
-
- for( Int_t iDaught=0; iDaught < particle->GetNDaughters(); iDaught++ ) {
- if( particle->GetDaughter(iDaught) < stack->GetNprimary() )
- nDaughters++;
- }
-
- if( nDaughters > 0 )
- continue;
- // --- check for double filling in stack
-
- /*Bool_t findable = kFALSE;
- for(Int_t iRec=0; iRec<esdEvent->GetNumberOfTracks(); ++iRec)
- {
- // check findable
- if(iMc == labelsAllRec[iRec])
- {
- findable = IsFindable(mcEvent,iMc);
- break;
- }
- }*/
- Bool_t findable = IsFindable(mcEvent,iMc);
-
- Bool_t recStatus = kFALSE;
- Int_t nClones = 0, nFakes = 0;
- for(Int_t iRec=0; iRec<esdEvent->GetNumberOfTracks(); ++iRec)
- {
- // check reconstructed
- if(iMc == labelsRec[iRec])
- {
- if (recStatus && nClones < fgkMaxClones) nClones++;
- recStatus = kTRUE;
- }
- //In order to relate the fake track to track parameters, we assign it to the best matching ESD track.
- if (labelsRec[iRec] < 0 && -labelsRec[iRec] == iMc && nFakes < fgkMaxFakes) nFakes++;
- }
-
- // Only 5 charged particle species (e,mu,pi,K,p)
- if (fCutsMC->IsPdgParticle(TMath::Abs(particle->GetPdgCode())) == kFALSE) continue;
-
- // transform Pdg to Pid
- Int_t pid = TransformToPID(particle);
-
- Float_t mceta = particle->Eta();
- Float_t mcphi = particle->Phi();
- if(mcphi<0) mcphi += 2.*TMath::Pi();
- Float_t mcpt = particle->Pt();
- Float_t charge = 0.;
- if (particle->GetPDG()->Charge() < 0) charge = -1.;
- else if (particle->GetPDG()->Charge() > 0) charge = 1.;
-
- // Fill histograms
- 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)};
- fEffHisto->Fill(vEffHisto);
- }
- if(labelsRec) delete [] labelsRec; labelsRec = 0;
- if(labelsAllRec) delete [] labelsAllRec; labelsAllRec = 0;
-}
-
-//_____________________________________________________________________________
-void AliPerformanceEff::ProcessTPCSec(AliMCEvent* const mcEvent, AliESDEvent *const esdEvent)
-{
- // Fill TPC only efficiency comparison information for secondaries
-
- if(!esdEvent) return;
-
- Int_t *labelsRecSec = NULL;
- labelsRecSec = new Int_t[esdEvent->GetNumberOfTracks()];
- if(!labelsRecSec)
- {
- Printf("Cannot create labelsRecSec");
- return;
- }
- for(Int_t i=0;i<esdEvent->GetNumberOfTracks();i++) { labelsRecSec[i] = 0; }
-
- Int_t *labelsAllRecSec = NULL;
- labelsAllRecSec = new Int_t[esdEvent->GetNumberOfTracks()];
- if(!labelsAllRecSec) {
- delete [] labelsRecSec;
- Printf("Cannot create labelsAllRecSec");
- return;
- }
- for(Int_t i=0;i<esdEvent->GetNumberOfTracks();i++) { labelsAllRecSec[i] = 0; }
-
- // loop over rec. tracks
- AliESDtrack *track=0;
- Int_t multAll=0, multRec=0;
- for (Int_t iTrack = 0; iTrack < esdEvent->GetNumberOfTracks(); iTrack++)
- {
- track = esdEvent->GetTrack(iTrack);
- if(!track) continue;
- if(track->Charge()==0) continue;
-
- // if not fUseKinkDaughters don't use tracks with kink index > 0
- if(!fUseKinkDaughters && track->GetKinkIndex(0) > 0) continue;
-
- //Int_t label = TMath::Abs(track->GetLabel());
- Int_t label = track->GetTPCLabel(); //Use TPC-only label for TPC-only efficiency analysis
- labelsAllRecSec[multAll]=label;
- multAll++;
-
- // TPC only
- if(IsRecTPC(track) != 0) {
- labelsRecSec[multRec]=label;
- multRec++;
- }
- }
-
- //
- // MC histograms for efficiency studies
- //
- if(mcEvent) {
-
- AliStack *stack = mcEvent->Stack();
- if (stack) {
-
- Int_t nPart = stack->GetNtrack();
- //Int_t nPart = stack->GetNprimary();
- for (Int_t iMc = 0; iMc < nPart; ++iMc)
- {
- if (iMc == 0) continue; //Cannot distinguish between track or fake track
- TParticle* particle = stack->Particle(iMc);
- if (!particle) continue;
- if (!particle->GetPDG()) continue;
- if (particle->GetPDG()->Charge() == 0.0) continue;
-
- // physical primary
- Bool_t prim = stack->IsPhysicalPrimary(iMc);
-
- // only secondaries which can be reconstructed at TPC
- if(prim) continue;
-
- //Float_t radius = TMath::Sqrt(particle->Vx()*particle->Vx()+particle->Vy()*particle->Vy()+particle->Vz()*particle->Vz());
- //if(radius > fCutsMC->GetMaxR()) continue;
-
- // only secondary electrons from gamma conversion
- //if( TMath::Abs(particle->GetPdgCode())!=fCutsMC->GetEM() || particle->GetUniqueID() != 5) continue;
-
- /*Bool_t findable = kFALSE;
- for(Int_t iRec=0; iRec<multAll; ++iRec)
- {
- // check findable
- if(iMc == labelsAllRecSec[iRec])
- {
- findable = IsFindable(mcEvent,iMc);
- break;
- }
- }*/
- Bool_t findable = IsFindable(mcEvent,iMc);
-
- Bool_t recStatus = kFALSE;
- Int_t nClones = 0, nFakes = 0;
- for(Int_t iRec=0; iRec<multRec; ++iRec)
- {
- // check reconstructed
- if(iMc == labelsRecSec[iRec])
- {
- if (recStatus && nClones < fgkMaxClones) nClones++;
- recStatus = kTRUE;
- }
- //In order to relate the fake track to track parameters, we assign it to the best matching ESD track.
- if (labelsRecSec[iRec] < 0 && -labelsRecSec[iRec] == iMc && nFakes < fgkMaxFakes) nFakes++;
- }
-
- // Only 5 charged particle species (e,mu,pi,K,p)
- if (fCutsMC->IsPdgParticle(TMath::Abs(particle->GetPdgCode())) == kFALSE) continue;
-
- // transform Pdg to Pid
- Int_t pid = TransformToPID(particle);
-
- Float_t mceta = particle->Eta();
- Float_t mcphi = particle->Phi();
- if(mcphi<0) mcphi += 2.*TMath::Pi();
- Float_t mcpt = particle->Pt();
- Float_t mcR = particle->R();
-
- // get info about mother
- Int_t motherLabel = particle->GetMother(0);
- if(motherLabel < 0) continue;
- TParticle *mother = stack->Particle(motherLabel);
- if(!mother) continue;
-
- Float_t mother_eta = mother->Eta();
- Float_t mother_phi = mother->Phi();
- if(mother_phi<0) mother_phi += 2.*TMath::Pi();
-
- Float_t charge = 0.;
- if (particle->GetPDG()->Charge() < 0) charge = -1.;
- else if (particle->GetPDG()->Charge() > 0) charge = 1.;
-
- // Fill histograms
- 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) };
- fEffSecHisto->Fill(vEffSecHisto);
- }
- }
- }
-
- if(labelsRecSec) delete [] labelsRecSec; labelsRecSec = 0;
- if(labelsAllRecSec) delete [] labelsAllRecSec; labelsAllRecSec = 0;
-}
-
-
-
-
-//_____________________________________________________________________________
-void AliPerformanceEff::ProcessTPCITS(AliMCEvent* const mcEvent, AliESDEvent *const esdEvent)
-{
- // Fill efficiency comparison information
-
- if(!esdEvent) return;
- if(!mcEvent) return;
-
- AliStack *stack = mcEvent->Stack();
- if (!stack) {
- AliDebug(AliLog::kError, "Stack not available");
- return;
- }
-
- Int_t *labelsRecTPCITS = NULL;
- labelsRecTPCITS = new Int_t[esdEvent->GetNumberOfTracks()];
- if(!labelsRecTPCITS)
- {
- Printf("Cannot create labelsRecTPCITS");
- return;
- }
- for(Int_t i=0;i<esdEvent->GetNumberOfTracks();i++) { labelsRecTPCITS[i] = 0; }
-
- Int_t *labelsAllRecTPCITS = NULL;
- labelsAllRecTPCITS = new Int_t[esdEvent->GetNumberOfTracks()];
- if(!labelsAllRecTPCITS) {
- delete [] labelsRecTPCITS;
- Printf("Cannot create labelsAllRecTPCITS");
- return;
- }
- for(Int_t i=0;i<esdEvent->GetNumberOfTracks();i++) { labelsAllRecTPCITS[i] = 0; }
-
- // loop over rec. tracks
- AliESDtrack *track=0;
- for (Int_t iTrack = 0; iTrack < esdEvent->GetNumberOfTracks(); iTrack++)
- {
- track = esdEvent->GetTrack(iTrack);
- if(!track) continue;
- if(track->Charge()==0) continue;
-
- // if not fUseKinkDaughters don't use tracks with kink index > 0
- if(!fUseKinkDaughters && track->GetKinkIndex(0) > 0) continue;
-
- //Int_t label = TMath::Abs(track->GetLabel());
- Int_t label = track->GetLabel(); //Use global label for combined efficiency analysis
- labelsAllRecTPCITS[iTrack]=label;
-
- // iTPC+ITS
- if(IsRecTPCITS(track) != 0)
- labelsRecTPCITS[iTrack]=label;
- }
-
- //
- // MC histograms for efficiency studies
- //
- //Int_t nPart = stack->GetNtrack();
- Int_t nPart = stack->GetNprimary();
- for (Int_t iMc = 0; iMc < nPart; ++iMc)
- {
- if (iMc == 0) continue; //Cannot distinguish between track or fake track
- TParticle* particle = stack->Particle(iMc);
- if (!particle) continue;
- if (!particle->GetPDG()) continue;
- if (particle->GetPDG()->Charge() == 0.0) continue;
-
- // physical primary
- Bool_t prim = stack->IsPhysicalPrimary(iMc);
- if(!prim) continue;
-
- /*Bool_t findable = kFALSE;
- for(Int_t iRec=0; iRec<esdEvent->GetNumberOfTracks(); ++iRec)
- {
- // check findable
- if(iMc == labelsAllRecTPCITS[iRec])
- {
- findable = IsFindable(mcEvent,iMc);
- break;
- }
- }*/
- Bool_t findable = IsFindable(mcEvent,iMc);
-
- Bool_t recStatus = kFALSE;
- Int_t nClones = 0, nFakes = 0;
- for(Int_t iRec=0; iRec<esdEvent->GetNumberOfTracks(); ++iRec)
- {
- // check reconstructed
- if(iMc == labelsRecTPCITS[iRec])
- {
- if (recStatus && nClones < fgkMaxClones) nClones++;
- recStatus = kTRUE;
- }
- //In order to relate the fake track to track parameters, we assign it to the best matching ESD track.
- if (labelsRecTPCITS[iRec] < 0 && -labelsRecTPCITS[iRec] == iMc && nFakes < fgkMaxFakes) nFakes++;
- }
-
- // Only 5 charged particle species (e,mu,pi,K,p)
- if (fCutsMC->IsPdgParticle(TMath::Abs(particle->GetPdgCode())) == kFALSE) continue;
-
- // transform Pdg to Pid
- Int_t pid = TransformToPID(particle);
-
- Float_t mceta = particle->Eta();
- Float_t mcphi = particle->Phi();
- if(mcphi<0) mcphi += 2.*TMath::Pi();
- Float_t mcpt = particle->Pt();
-
- Float_t charge = 0.;
- if (particle->GetPDG()->Charge() < 0) charge = -1.;
- else if (particle->GetPDG()->Charge() > 0) charge = 1.;
-
- // Fill histograms
- 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)};
- fEffHisto->Fill(vEffHisto);
- }
-
- if(labelsRecTPCITS) delete [] labelsRecTPCITS; labelsRecTPCITS = 0;
- if(labelsAllRecTPCITS) delete [] labelsAllRecTPCITS; labelsAllRecTPCITS = 0;
-}
-
-//_____________________________________________________________________________
-void AliPerformanceEff::ProcessConstrained(AliMCEvent* const mcEvent, AliESDEvent *const esdEvent)
-{
- // Process comparison information
- if(!esdEvent) return;
- if(!mcEvent) return;
-
- AliStack *stack = mcEvent->Stack();
- if (!stack) {
- AliDebug(AliLog::kError, "Stack not available");
- return;
- }
-
- Int_t *labelsRecConstrained = NULL;
- labelsRecConstrained = new Int_t[esdEvent->GetNumberOfTracks()];
- if(!labelsRecConstrained)
- {
- Printf("Cannot create labelsRecConstrained");
- return;
- }
- for(Int_t i=0;i<esdEvent->GetNumberOfTracks();i++) { labelsRecConstrained[i] = 0; }
-
- Int_t *labelsAllRecConstrained = NULL;
- labelsAllRecConstrained = new Int_t[esdEvent->GetNumberOfTracks()];
- if(!labelsAllRecConstrained) {
- delete [] labelsRecConstrained;
- Printf("Cannot create labelsAllRecConstrained");
- return;
- }
- for(Int_t i=0;i<esdEvent->GetNumberOfTracks();i++) { labelsAllRecConstrained[i] = 0; }
-
- // loop over rec. tracks
- AliESDtrack *track=0;
- for (Int_t iTrack = 0; iTrack < esdEvent->GetNumberOfTracks(); iTrack++)
- {
- track = esdEvent->GetTrack(iTrack);
- if(!track) continue;
- if(track->Charge()==0) continue;
-
- // if not fUseKinkDaughters don't use tracks with kink index > 0
- if(!fUseKinkDaughters && track->GetKinkIndex(0) > 0) continue;
-
- //Int_t label = TMath::Abs(track->GetLabel());
- Int_t label = track->GetLabel();
- labelsAllRecConstrained[iTrack]=label;
-
- // Constrained
- if(IsRecConstrained(track) != 0)
- labelsRecConstrained[iTrack]=label;
-
- }
-
- //
- // MC histograms for efficiency studies
- //
-
-
- //Int_t nPart = stack->GetNtrack();
- Int_t nPart = stack->GetNprimary();
- for (Int_t iMc = 0; iMc < nPart; ++iMc)
- {
- if (iMc == 0) continue; //Cannot distinguish between track or fake track
- TParticle* particle = stack->Particle(iMc);
- if (!particle) continue;
- if (!particle->GetPDG()) continue;
- if (particle->GetPDG()->Charge() == 0.0) continue;
-
- // physical primary
- Bool_t prim = stack->IsPhysicalPrimary(iMc);
- if(!prim) continue;
-
- /*Bool_t findable = kFALSE;
- for(Int_t iRec=0; iRec<esdEvent->GetNumberOfTracks(); ++iRec)
- {
- // check findable
- if(iMc == labelsAllRecConstrained[iRec])
- {
- findable = IsFindable(mcEvent,iMc);
- break;
- }
- }*/
- Bool_t findable = IsFindable(mcEvent,iMc);
-
- Bool_t recStatus = kFALSE;
- Int_t nClones = 0, nFakes = 0;
- for(Int_t iRec=0; iRec<esdEvent->GetNumberOfTracks(); ++iRec)
- {
- // check reconstructed
- if(iMc == labelsRecConstrained[iRec])
- {
- if (recStatus && nClones < fgkMaxClones) nClones++;
- recStatus = kTRUE;
- }
- //In order to relate the fake track to track parameters, we assign it to the best matching ESD track.
- if (labelsRecConstrained[iRec] < 0 && -labelsRecConstrained[iRec] == iMc && nFakes < fgkMaxFakes) nFakes++;
- }
-
- // Only 5 charged particle species (e,mu,pi,K,p)
- if (fCutsMC->IsPdgParticle(TMath::Abs(particle->GetPdgCode())) == kFALSE) continue;
-
- // transform Pdg to Pid
- Int_t pid = TransformToPID(particle);
-
- Float_t mceta = particle->Eta();
- Float_t mcphi = particle->Phi();
- if(mcphi<0) mcphi += 2.*TMath::Pi();
- Float_t mcpt = particle->Pt();
-
- Float_t charge = 0.;
- if (particle->GetPDG()->Charge() < 0) charge = -1.;
- else if (particle->GetPDG()->Charge() > 0) charge = 1.;
-
- // Fill histograms
- 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) };
- fEffHisto->Fill(vEffHisto);
- }
-
- if(labelsRecConstrained) delete [] labelsRecConstrained; labelsRecConstrained = 0;
- if(labelsAllRecConstrained) delete [] labelsAllRecConstrained; labelsAllRecConstrained = 0;
-}
-
-//_____________________________________________________________________________
-void AliPerformanceEff::Exec(AliMCEvent* const mcEvent, AliESDEvent *const esdEvent, AliESDfriend *const esdFriend, const Bool_t bUseMC, const Bool_t bUseESDfriend)
-{
- // Process comparison information
- //
- if(!esdEvent)
- {
- Error("Exec","esdEvent not available");
- return;
- }
- AliHeader* header = 0;
- AliGenEventHeader* genHeader = 0;
- AliStack* stack = 0;
- TArrayF vtxMC(3);
-
- if(bUseMC)
- {
- if(!mcEvent) {
- Error("Exec","mcEvent not available");
- return;
- }
- // get MC event header
- header = mcEvent->Header();
- if (!header) {
- Error("Exec","Header not available");
- return;
- }
- // MC particle stack
- stack = mcEvent->Stack();
- if (!stack) {
- Error("Exec","Stack not available");
- return;
- }
- // get MC vertex
- genHeader = header->GenEventHeader();
- if (!genHeader) {
- Error("Exec","Could not retrieve genHeader from Header");
- return;
- }
- genHeader->PrimaryVertex(vtxMC);
- }
- else {
- Error("Exec","MC information required!");
- return;
- }
-
- // use ESD friends
- if(bUseESDfriend) {
- if(!esdFriend) {
- Error("Exec","esdFriend not available");
- return;
- }
- }
-
- //
- // Process events
- //
- if(GetAnalysisMode() == 0) ProcessTPC(mcEvent,esdEvent);
- else if(GetAnalysisMode() == 1) ProcessTPCITS(mcEvent,esdEvent);
- else if(GetAnalysisMode() == 2) ProcessConstrained(mcEvent,esdEvent);
- else if(GetAnalysisMode() == 5) ProcessTPCSec(mcEvent,esdEvent);
- else {
- printf("ERROR: AnalysisMode %d \n",fAnalysisMode);
- return;
- }
-}
-
-//_____________________________________________________________________________
-Int_t AliPerformanceEff::TransformToPID(TParticle *particle)
-{
-// transform Pdg to Pid
-// Pdg convension is different for hadrons and leptons
-// (e.g. K+/K- = 321/-321; e+/e- = -11/11 )
-
- Int_t pid = -1;
- if( TMath::Abs(particle->GetPdgCode())==fCutsMC->GetEM() ) pid = 0;
- if( TMath::Abs(particle->GetPdgCode())==fCutsMC->GetMuM() ) pid = 1;
- if( TMath::Abs(particle->GetPdgCode())==fCutsMC->GetPiP() ) pid = 2;
- if( TMath::Abs(particle->GetPdgCode())==fCutsMC->GetKP() ) pid = 3;
- if( TMath::Abs(particle->GetPdgCode())==fCutsMC->GetProt() ) pid = 4;
-
-return pid;
-}
-
-//_____________________________________________________________________________
-Bool_t AliPerformanceEff::IsFindable(const AliMCEvent *mcEvent, Int_t label)
-{
-//
-// Findfindable tracks
-//
-if(!mcEvent) return kFALSE;
-
- AliMCParticle *mcParticle = (AliMCParticle*) mcEvent->GetTrack(label);
- if(!mcParticle) return kFALSE;
-
- Int_t counter;
- Float_t tpcTrackLength = mcParticle->GetTPCTrackLength(AliTracker::GetBz(),0.05,counter,3.0);
- //printf("tpcTrackLength %f \n", tpcTrackLength);
-
-return (tpcTrackLength>fCutsMC->GetMinTrackLength());
-}
-
-//_____________________________________________________________________________
-Bool_t AliPerformanceEff::IsRecTPC(AliESDtrack *esdTrack)
-{
-//
-// Check whether track is reconstructed in TPC
-//
-if(!esdTrack) return kFALSE;
-
- const AliExternalTrackParam *track = esdTrack->GetTPCInnerParam();
- if(!track) return kFALSE;
-
- Float_t dca[2], cov[3]; // dca_xy, dca_z, sigma_xy, sigma_xy_z, sigma_z
- esdTrack->GetImpactParametersTPC(dca,cov);
-
- Bool_t recStatus = kFALSE;
- if(esdTrack->GetTPCNcls()>fCutsRC->GetMinNClustersTPC()) recStatus = kTRUE;
-
- /*
- if( TMath::Abs(dca[0])<fCutsRC->GetMaxDCAToVertexXY() &&
- TMath::Abs(dca[1])<fCutsRC->GetMaxDCAToVertexZ())
- {
- recStatus = kTRUE;
- }
- */
-
-return recStatus;
-}
-
-//_____________________________________________________________________________
-Bool_t AliPerformanceEff::IsRecTPCITS(AliESDtrack *esdTrack)
-{
-//
-// Check whether track is reconstructed in TPCITS
-//
-if(!esdTrack) return kFALSE;
-
- Float_t dca[2], cov[3]; // dca_xy, dca_z, sigma_xy, sigma_xy_z, sigma_z
- esdTrack->GetImpactParameters(dca,cov);
-
- Bool_t recStatus = kFALSE;
-
- if ((esdTrack->GetStatus()&AliESDtrack::kTPCrefit)==0) return kFALSE; // TPC refit
- if ((esdTrack->GetStatus()&AliESDtrack::kITSrefit)==0) return kFALSE; // ITS refit
- if (esdTrack->GetITSclusters(0)<fCutsRC->GetMinNClustersITS()) return kFALSE; // min. nb. ITS clusters
- //if ((esdTrack->GetStatus()&AliESDtrack::kITSrefit)==0) return kFALSE; // ITS refit
- //Int_t clusterITS[200];
- //if(esdTrack->GetITSclusters(clusterITS)<fCutsRC->GetMinNClustersITS()) return kFALSE; // min. nb. ITS clusters
-
- recStatus = kTRUE;
- /*
- if(TMath::Abs(dca[0])<fCutsRC->GetMaxDCAToVertexXY() &&
- TMath::Abs(dca[1])<fCutsRC->GetMaxDCAToVertexZ())
- {
- recStatus = kTRUE;
- }
- */
-
-return recStatus;
-}
-
-//_____________________________________________________________________________
-Bool_t AliPerformanceEff::IsRecConstrained(AliESDtrack *esdTrack)
-{
-//
-// Check whether track is reconstructed in IsRecConstrained
-//
- if(!esdTrack) return kFALSE;
-
- const AliExternalTrackParam * track = esdTrack->GetConstrainedParam();
- if(!track) return kFALSE;
-
- Float_t dca[2], cov[3]; // dca_xy, dca_z, sigma_xy, sigma_xy_z, sigma_z
- esdTrack->GetImpactParameters(dca,cov);
- //Int_t label = TMath::Abs(esdTrack->GetLabel());
-
- Bool_t recStatus = kFALSE;
-
- if ((esdTrack->GetStatus()&AliESDtrack::kTPCrefit)==0) return kFALSE; // TPC refit
- if (esdTrack->GetTPCNcls()<fCutsRC->GetMinNClustersTPC()) return kFALSE; // min. nb. TPC clusters
- Int_t clusterITS[200];
- if(esdTrack->GetITSclusters(clusterITS)<fCutsRC->GetMinNClustersITS()) return kFALSE; // min. nb. ITS clusters
-
- recStatus = kTRUE;
-
- /*
- if(TMath::Abs(dca[0])<fCutsRC->GetMaxDCAToVertexXY() &&
- TMath::Abs(dca[1])<fCutsRC->GetMaxDCAToVertexZ())
- {
- recStatus = kTRUE;
- }
- */
-
-return recStatus;
-}
-
-//_____________________________________________________________________________
-Long64_t AliPerformanceEff::Merge(TCollection* const list)
-{
- // Merge list of objects (needed by PROOF)
-
- if (!list)
- return 0;
-
- if (list->IsEmpty())
- return 1;
-
- TIterator* iter = list->MakeIterator();
- TObject* obj = 0;
-
- // collection of generated histograms
-
- Int_t count=0;
- while((obj = iter->Next()) != 0)
- {
- AliPerformanceEff* entry = dynamic_cast<AliPerformanceEff*>(obj);
- if (entry == 0) continue;
-
- fEffHisto->Add(entry->fEffHisto);
- fEffSecHisto->Add(entry->fEffSecHisto);
- count++;
- }
-
-return count;
-}
-
-//_____________________________________________________________________________
-void AliPerformanceEff::Analyse()
-{
- // Analyse comparison information and store output histograms
- // in the folder "folderEff"
- //
- TH1::AddDirectory(kFALSE);
- TObjArray *aFolderObj = new TObjArray;
- if(!aFolderObj) return;
- // char title[256];
-
- //
- // efficiency vs pt
- //
-
- if(GetAnalysisMode() != 5) {
-
- fEffHisto->GetAxis(0)->SetRangeUser(-0.9,0.89); // eta range
- fEffHisto->GetAxis(2)->SetRangeUser(0.1,19.99); // pt range // FIXME maybe remove since range is defined in THnSparse
-
- // rec efficiency vs pt
- fEffHisto->GetAxis(3)->SetRangeUser(0.,3.99); // reconstructed
-
- fEffHisto->GetAxis(6)->SetRangeUser(-2.,2.); // charge all
- aFolderObj->Add(AddHistoEff(2, "ptRecEff", "rec. efficiency", 0));
- aFolderObj->Add(AddHistoEff(2, "ptClone", "clone rate", 1));
- aFolderObj->Add(AddHistoEff(2, "ptFake", "fake rate", 2));
-
- fEffHisto->GetAxis(6)->SetRangeUser(-2.,0.); // charge negativ
- aFolderObj->Add(AddHistoEff(2, "ptRecEffNeg", "rec. efficiency neg.", 0));
-
- fEffHisto->GetAxis(6)->SetRangeUser(0.,2.); // charge positiv
- aFolderObj->Add(AddHistoEff(2, "ptRecEffPos", "rec. efficiency pos.", 0));
-
- // rec efficiency vs pid vs pt
- fEffHisto->GetAxis(3)->SetRange(3,3); // pions
-
- fEffHisto->GetAxis(6)->SetRangeUser(-2.,2.); // charge all
- aFolderObj->Add(AddHistoEff(2, "ptRecEffPi", "rec. efficiency (pions)", 0));
-
- fEffHisto->GetAxis(6)->SetRangeUser(-2.,0.); // charge negativ
- aFolderObj->Add(AddHistoEff(2, "ptRecEffPiNeg", "rec. efficiency (pions) neg.", 0));
-
- fEffHisto->GetAxis(6)->SetRangeUser(0.,2.); // charge positiv
- aFolderObj->Add(AddHistoEff(2, "ptRecEffPiPos", "rec. efficiency (pions) pos.", 0));
-
-
- fEffHisto->GetAxis(3)->SetRange(4,4); // kaons
-
- fEffHisto->GetAxis(6)->SetRangeUser(-2.,2.); // charge all
- aFolderObj->Add(AddHistoEff(2, "ptRecEffK", "rec. efficiency (kaons)", 0));
-
- fEffHisto->GetAxis(6)->SetRangeUser(-2.,0.); // charge negativ
- aFolderObj->Add(AddHistoEff(2, "ptRecEffKNeg", "rec. efficiency (kaons) neg.", 0));
-
- fEffHisto->GetAxis(6)->SetRangeUser(0.,2.); // charge positiv
- aFolderObj->Add(AddHistoEff(2, "ptRecEffKPos", "rec. efficiency (kaons) pos.", 0));
-
-
- fEffHisto->GetAxis(3)->SetRange(5,5); // protons
-
- fEffHisto->GetAxis(6)->SetRangeUser(-2.,2.); // charge all
- aFolderObj->Add(AddHistoEff(2, "ptRecEffP", "rec. efficiency (protons)", 0));
-
- fEffHisto->GetAxis(6)->SetRangeUser(-2.,0.); // charge negativ
- aFolderObj->Add(AddHistoEff(2, "ptRecEffPNeg", "rec. efficiency (protons) neg.", 0));
-
- fEffHisto->GetAxis(6)->SetRangeUser(0.,2.); // charge positiv
- aFolderObj->Add(AddHistoEff(2, "ptRecEffPPos", "rec. efficiency (protons) pos.", 0));
-
- // findable efficiency vs pt
- fEffHisto->GetAxis(3)->SetRangeUser(0.,4.);
- fEffHisto->GetAxis(5)->SetRange(2,2); // findable
-
- fEffHisto->GetAxis(6)->SetRangeUser(-2.,2.); // charge all
- aFolderObj->Add(AddHistoEff(2, "ptRecEffF", "rec. efficiency (findable)", 0));
- aFolderObj->Add(AddHistoEff(2, "ptCloneF", "clone rate (findable)", 1));
- aFolderObj->Add(AddHistoEff(2, "ptFakeF", "fake rate (findable)", 2));
-
-
- fEffHisto->GetAxis(6)->SetRangeUser(-2.,0.); // charge negativ
- aFolderObj->Add(AddHistoEff(2, "ptRecEffFNeg", "rec. efficiency (findable) neg.", 0));
-
- fEffHisto->GetAxis(6)->SetRangeUser(0.,2.); // charge positiv
- aFolderObj->Add(AddHistoEff(2, "ptRecEffFPos", "rec. efficiency (findable) pos.", 0));
-
- //
- // efficiency vs eta
- //
-
- fEffHisto->GetAxis(0)->SetRangeUser(-1.5,1.49); // eta range
- fEffHisto->GetAxis(2)->SetRangeUser(0.1,19.99); // pt range
- fEffHisto->GetAxis(5)->SetRangeUser(0.,1.0); // all
-
- // rec efficiency vs eta
- fEffHisto->GetAxis(3)->SetRangeUser(0.,4.); // reconstructed
-
- fEffHisto->GetAxis(6)->SetRangeUser(-2.,2.); // charge all
- aFolderObj->Add(AddHistoEff(0, "etaRecEff", "rec. efficiency", 0));
- aFolderObj->Add(AddHistoEff(0, "etaClone", "clone rate", 1));
- aFolderObj->Add(AddHistoEff(0, "etaFake", "fake rate", 2));
-
-
- fEffHisto->GetAxis(6)->SetRangeUser(-2.,0.); // charge negativ
- aFolderObj->Add(AddHistoEff(0, "etaRecEffNeg", "rec. efficiency neg.", 0));
-
- fEffHisto->GetAxis(6)->SetRangeUser(0.,2.); // charge positiv
- aFolderObj->Add(AddHistoEff(0, "etaRecEffPos", "rec. efficiency pos.", 0));
-
- // rec efficiency vs pid vs eta
- fEffHisto->GetAxis(3)->SetRange(3,3); // pions
-
- fEffHisto->GetAxis(6)->SetRangeUser(-2.,2.); // charge all
- aFolderObj->Add(AddHistoEff(0, "etaRecEffPi", "rec. efficiency (pions)", 0));
-
- fEffHisto->GetAxis(6)->SetRangeUser(-2.,0.); // charge negativ
- aFolderObj->Add(AddHistoEff(0, "etaRecEffPiNeg", "rec. efficiency (pions) neg.", 0));
-
- fEffHisto->GetAxis(6)->SetRangeUser(0.,2.); // charge positiv
- aFolderObj->Add(AddHistoEff(0, "etaRecEffPiPos", "rec. efficiency (pions) pos.", 0));
-
-
- fEffHisto->GetAxis(3)->SetRange(4,4); // kaons
-
- fEffHisto->GetAxis(6)->SetRangeUser(-2.,2.); // charge all
- aFolderObj->Add(AddHistoEff(0, "etaRecEffK", "rec. efficiency (kaons)", 0));
-
- fEffHisto->GetAxis(6)->SetRangeUser(-2.,0.); // charge negativ
- aFolderObj->Add(AddHistoEff(0, "etaRecEffKNeg", "rec. efficiency (kaons) neg.", 0));
-
- fEffHisto->GetAxis(6)->SetRangeUser(0.,2.); // charge positiv
- aFolderObj->Add(AddHistoEff(0, "etaRecEffKPos", "rec. efficiency (kaons) pos.", 0));
-
-
- fEffHisto->GetAxis(3)->SetRange(5,5); // protons
-
- fEffHisto->GetAxis(6)->SetRangeUser(-2.,2.); // charge all
- aFolderObj->Add(AddHistoEff(0, "etaRecEffP", "rec. efficiency (protons)", 0));
-
- fEffHisto->GetAxis(6)->SetRangeUser(-2.,0.); // charge negativ
- aFolderObj->Add(AddHistoEff(0, "etaRecEffPNeg", "rec. efficiency (protons) neg.", 0));
-
- fEffHisto->GetAxis(6)->SetRangeUser(0.,2.); // charge positiv
- aFolderObj->Add(AddHistoEff(0, "etaRecEffPPos", "rec. efficiency (protons) pos.", 0));
-
- // findable efficiency vs eta
- fEffHisto->GetAxis(3)->SetRangeUser(0.,4.);
- fEffHisto->GetAxis(5)->SetRange(2,2); // findable
-
- fEffHisto->GetAxis(6)->SetRangeUser(-2.,2.); // charge all
- aFolderObj->Add(AddHistoEff(0, "etaRecEffF", "rec. efficiency (findable)", 0));
- aFolderObj->Add(AddHistoEff(0, "etaCloneF", "clone rate (findable)", 1));
- aFolderObj->Add(AddHistoEff(0, "etaFakeF", "fake rate (findable)", 2));
-
-
- fEffHisto->GetAxis(6)->SetRangeUser(-2.,0.); // charge negativ
- aFolderObj->Add(AddHistoEff(0, "etaRecEffFNeg", "rec. efficiency (findable) neg.", 0));
-
- fEffHisto->GetAxis(6)->SetRangeUser(0.,2.); // charge positiv
- aFolderObj->Add(AddHistoEff(0, "etaRecEffFPos", "rec. efficiency (findable) pos.", 0));
-
- //
- // efficiency vs phi
- //
-
- fEffHisto->GetAxis(0)->SetRangeUser(-0.9,0.89); // eta range
- fEffHisto->GetAxis(2)->SetRangeUser(0.1,19.99); // pt range
- fEffHisto->GetAxis(5)->SetRangeUser(0.,1.); // all
-
- // rec efficiency vs phi
- fEffHisto->GetAxis(3)->SetRangeUser(0.,4.); // reconstructed
-
- fEffHisto->GetAxis(6)->SetRangeUser(-2.,2.); // charge all
- aFolderObj->Add(AddHistoEff(1, "phiRecEff", "rec. efficiency", 0));
- aFolderObj->Add(AddHistoEff(1, "phiClone", "clone rate", 1));
- aFolderObj->Add(AddHistoEff(1, "phiFake", "fake rate", 2));
-
- fEffHisto->GetAxis(6)->SetRangeUser(-2.,0.); // charge negativ
- aFolderObj->Add(AddHistoEff(1, "phiRecEffNeg", "rec. efficiency neg.", 0));
-
- fEffHisto->GetAxis(6)->SetRangeUser(0.,2.); // charge positiv
- aFolderObj->Add(AddHistoEff(1, "phiRecEffPos", "rec. efficiency pos.", 0));
-
- // rec efficiency vs pid vs phi
- fEffHisto->GetAxis(3)->SetRange(3,3); // pions
-
- fEffHisto->GetAxis(6)->SetRangeUser(-2.,2.); // charge all
- aFolderObj->Add(AddHistoEff(1, "phiRecEffPi", "rec. efficiency (pions)", 0));
-
- fEffHisto->GetAxis(6)->SetRangeUser(-2.,0.); // charge negativ
- aFolderObj->Add(AddHistoEff(1, "phiRecEffPiNeg", "rec. efficiency (pions) neg.", 0));
-
- fEffHisto->GetAxis(6)->SetRangeUser(0.,2.); // charge positiv
- aFolderObj->Add(AddHistoEff(1, "phiRecEffPiPos", "rec. efficiency (pions) pos.", 0));
-
-
- fEffHisto->GetAxis(3)->SetRange(4,4); // kaons
-
- fEffHisto->GetAxis(6)->SetRangeUser(-2.,2.); // charge all
- aFolderObj->Add(AddHistoEff(1, "phiRecEffK", "rec. efficiency (kaons)", 0));
-
- fEffHisto->GetAxis(6)->SetRangeUser(-2.,0.); // charge negativ
- aFolderObj->Add(AddHistoEff(1, "phiRecEffKNeg", "rec. efficiency (kaons) neg.", 0));
-
- fEffHisto->GetAxis(6)->SetRangeUser(0.,2.); // charge positiv
- aFolderObj->Add(AddHistoEff(1, "phiRecEffKPos", "rec. efficiency (kaons) pos.", 0));
-
-
- fEffHisto->GetAxis(3)->SetRange(5,5); // protons
-
- fEffHisto->GetAxis(6)->SetRangeUser(-2.,2.); // charge all
- aFolderObj->Add(AddHistoEff(1, "phiRecEffP", "rec. efficiency (protons)", 0));
-
- fEffHisto->GetAxis(6)->SetRangeUser(-2.,0.); // charge negativ
- aFolderObj->Add(AddHistoEff(1, "phiRecEffPNeg", "rec. efficiency (protons) neg.", 0));
-
- fEffHisto->GetAxis(6)->SetRangeUser(0.,2.); // charge positiv
- aFolderObj->Add(AddHistoEff(1, "phiRecEffPPos", "rec. efficiency (protons) pos.", 0));
-
- // findable efficiency vs phi
- fEffHisto->GetAxis(3)->SetRangeUser(0.,4.);
- fEffHisto->GetAxis(5)->SetRange(2,2); // findable
-
- fEffHisto->GetAxis(6)->SetRangeUser(-2.,2.); // charge all
- aFolderObj->Add(AddHistoEff(1, "phiRecEffF", "rec. efficiency (findable)", 0));
- aFolderObj->Add(AddHistoEff(1, "phiCloneF", "clone rate (findable)", 1));
- aFolderObj->Add(AddHistoEff(1, "phiFakeF", "fake rate (findable)", 2));
-
-
- fEffHisto->GetAxis(6)->SetRangeUser(-2.,0.); // charge negativ
- aFolderObj->Add(AddHistoEff(1, "phiRecEffFNeg", "rec. efficiency (findable) neg.", 0));
-
- fEffHisto->GetAxis(6)->SetRangeUser(0.,2.); // charge positiv
- aFolderObj->Add(AddHistoEff(1, "phiRecEffFPos", "rec. efficiency (findable) pos.", 0));
- }
- else {
- //
- Float_t minEta=-1.5, maxEta=1.49;
- Float_t minR=0.0, maxR=150.0;
- Float_t minPt=0.10, maxPt=19.99;
-
- // mother eta range
- fEffSecHisto->GetAxis(8)->SetRangeUser(minEta,maxEta);
-
- // particle creation radius range
- fEffSecHisto->GetAxis(6)->SetRangeUser(minR,maxR);
-
- //
- fEffSecHisto->GetAxis(0)->SetRangeUser(minEta,maxEta);
- fEffSecHisto->GetAxis(2)->SetRangeUser(minPt,maxPt);
-
- // rec efficiency vs pt
-
- aFolderObj->Add(AddHistoEff(2, "ptRecEff", "rec. efficiency", 0, 1));
- aFolderObj->Add(AddHistoEff(2, "ptClone", "clone rate", 1, 1));
- aFolderObj->Add(AddHistoEff(2, "ptFake", "fake rate", 2, 1));
-
- // rec efficiency vs pid vs pt
-
- fEffSecHisto->GetAxis(3)->SetRange(1,1); // electrons
- aFolderObj->Add(AddHistoEff(2, "ptRecEffEle", "rec. efficiency (electrons)", 0, 1));
-
- fEffSecHisto->GetAxis(3)->SetRange(3,3); // pions
- aFolderObj->Add(AddHistoEff(2, "ptRecEffPi", "rec. efficiency (pions)", 0, 1));
-
- fEffSecHisto->GetAxis(3)->SetRange(4,4); // kaons
- aFolderObj->Add(AddHistoEff(2, "ptRecEffK", "rec. efficiency (kaons)", 0, 1));
-
- fEffSecHisto->GetAxis(3)->SetRange(5,5); // protons
- aFolderObj->Add(AddHistoEff(2, "ptRecEffP", "rec. efficiency (protons)", 0, 1));
-
- fEffSecHisto->GetAxis(3)->SetRangeUser(0.,4.);
-
- // findable efficiency vs pt
-
- fEffSecHisto->GetAxis(5)->SetRange(2,2); // findable
- aFolderObj->Add(AddHistoEff(2, "ptRecEffF", "rec. efficiency (findable)", 0, 1));
- aFolderObj->Add(AddHistoEff(2, "ptCloneF", "clone rate (findable)", 1, 1));
- aFolderObj->Add(AddHistoEff(2, "ptFakeF", "fake rate (findable)", 2, 1));
-
- //
- // efficiency vs eta
- //
- fEffSecHisto->GetAxis(2)->SetRangeUser(minPt,maxPt);
- fEffSecHisto->GetAxis(4)->SetRangeUser(0.,1.); // all
- fEffSecHisto->GetAxis(5)->SetRangeUser(0.,1.); // all
-
- aFolderObj->Add(AddHistoEff(0, "etaRecEff", "rec. efficiency", 0, 1));
- aFolderObj->Add(AddHistoEff(0, "etaClone", "clone rate", 1, 1));
- aFolderObj->Add(AddHistoEff(0, "etaFake", "fake rate", 2, 1));
-
- // rec efficiency vs pid vs eta
- fEffSecHisto->GetAxis(3)->SetRange(1,1); // electrons
- aFolderObj->Add(AddHistoEff(0, "etaRecEffEle", "rec. efficiency (electrons)", 0, 1));
-
- fEffSecHisto->GetAxis(3)->SetRange(3,3); // pions
- aFolderObj->Add(AddHistoEff(0, "etaRecEffPi", "rec. efficiency (pions)", 0, 1));
-
- fEffSecHisto->GetAxis(3)->SetRange(4,4); // kaons
- aFolderObj->Add(AddHistoEff(0, "etaRecEffK", "rec. efficiency (kaons)", 0, 1));
-
- fEffSecHisto->GetAxis(3)->SetRange(5,5); // protons
- aFolderObj->Add(AddHistoEff(0, "etaRecEffP", "rec. efficiency (protons)", 0, 1));
-
- fEffSecHisto->GetAxis(3)->SetRangeUser(0.,4.);
-
- // findable efficiency vs eta
-
- fEffSecHisto->GetAxis(5)->SetRange(2,2); // findable
- aFolderObj->Add(AddHistoEff(0, "etaRecEffF", "rec. efficiency (findable)", 0, 1));
- aFolderObj->Add(AddHistoEff(0, "etaCloneF", "clone rate (findable)", 1, 1));
- aFolderObj->Add(AddHistoEff(0, "etaFakeF", "fake rate (findable)", 2, 1));
-
- //
- // efficiency vs phi
- //
-
- fEffSecHisto->GetAxis(0)->SetRangeUser(minEta,maxEta);
- fEffSecHisto->GetAxis(2)->SetRangeUser(minPt,maxPt);
-
- fEffSecHisto->GetAxis(4)->SetRangeUser(0.,1.); // all
- fEffSecHisto->GetAxis(5)->SetRangeUser(0.,1.); // all
-
- aFolderObj->Add(AddHistoEff(1, "phiRecEff", "rec. efficiency", 0, 1));
- aFolderObj->Add(AddHistoEff(1, "phiClone", "clone rate", 1, 1));
- aFolderObj->Add(AddHistoEff(1, "phiFake", "fake rate", 2, 1));
-
- // rec efficiency vs pid vs phi
- fEffSecHisto->GetAxis(3)->SetRange(1,1);
- aFolderObj->Add(AddHistoEff(1, "phiRecEffEle", "rec. efficiency (electrons)", 0, 1));
-
- fEffSecHisto->GetAxis(3)->SetRange(3,3); // pions
- aFolderObj->Add(AddHistoEff(1, "phiRecEffPi", "rec. efficiency (pions)", 0, 1));
-
- fEffSecHisto->GetAxis(3)->SetRange(4,4); // kaons
- aFolderObj->Add(AddHistoEff(1, "phiRecEffK", "rec. efficiency (kaons)", 0, 1));
-
- fEffSecHisto->GetAxis(3)->SetRange(5,5); // protons
- aFolderObj->Add(AddHistoEff(1, "phiRecEffP", "rec. efficiency (protons)", 0, 1));
-
- fEffSecHisto->GetAxis(3)->SetRangeUser(0.,4.);
-
- // findable efficiency vs phi
-
- fEffSecHisto->GetAxis(5)->SetRange(2,2); // findable
- aFolderObj->Add(AddHistoEff(1, "phiRecEffF", "rec. efficiency (findable)", 0, 1));
- aFolderObj->Add(AddHistoEff(1, "phiCloneF", "clone rate (findable)", 1, 1));
- aFolderObj->Add(AddHistoEff(1, "phiFakeF", "fake rate (findable)", 2, 1));
- }
-
- for (Int_t i = 0;i < fEffHisto->GetNdimensions();i++)
- {
- fEffHisto->GetAxis(i)->SetRange(1,0); //Reset Range
- }
-
- // export objects to analysis folder
- fAnalysisFolder = ExportToFolder(aFolderObj);
-
- // delete only TObjArray
- if(aFolderObj) delete aFolderObj;
-}
-
-//_____________________________________________________________________________
-TFolder* AliPerformanceEff::ExportToFolder(TObjArray * array)
-{
- // recreate folder avery time and export objects to new one
- //
- AliPerformanceEff * comp=this;
- TFolder *folder = comp->GetAnalysisFolder();
-
- TString name, title;
- TFolder *newFolder = 0;
- Int_t i = 0;
- Int_t size = array->GetSize();
-
- if(folder) {
- // get name and title from old folder
- name = folder->GetName();
- title = folder->GetTitle();
-
- // delete old one
- delete folder;
-
- // create new one
- newFolder = CreateFolder(name.Data(),title.Data());
- newFolder->SetOwner();
-
- // add objects to folder
- while(i < size) {
- newFolder->Add(array->At(i));
- i++;
- }
- }
-
-return newFolder;
-}
-
-
-//_____________________________________________________________________________
-TFolder* AliPerformanceEff::CreateFolder(TString name,TString title) {
-// create folder for analysed histograms
-//
-TFolder *folder = 0;
- folder = new TFolder(name.Data(),title.Data());
-
- return folder;
-}
-
-TH1D* WeightedProjection(THnSparseF* src, Int_t axis, Int_t nWeights, Int_t* weightCoords)
-{
- THnSparseF* tmp = (THnSparseF*) src->Clone();
- Int_t i;
- for (i = 0;i < tmp->GetNbins();i++)
- {
- Int_t coords[12];
- tmp->GetBinContent(i, coords);
- Int_t weight = 0, j;
- for (j = 0;j < nWeights;j++)
- {
- //The coordinates range from 1 to maxClones / maxFakes + 1, so we have to subtract one
- weight += coords[weightCoords[j]] - 1;
- }
- tmp->SetBinContent(i, weight);
- }
-
- TH1D* ret = tmp->Projection(axis);
- delete tmp;
- return(ret);
-}
-
-//_____________________________________________________________________________
-TH1D* AliPerformanceEff::AddHistoEff(Int_t axis, const Char_t *name, const Char_t* vsTitle,
- const Int_t type, const Int_t secondary) {
- // Create and add rec efficiency vs pt, eta, phi
-
- char title[256];
-
- TH1D *recc = NULL;
-
- THnSparseF* EffHisto = secondary ? fEffSecHisto : fEffHisto;
-
- Int_t axis_clone = secondary ? 10 : 7;
- Int_t axis_fake = secondary ? 11 : 8;
- Int_t axis_all[3] = {4, axis_clone, axis_fake};
-
-
- if (type == 0) // Efficiency
- {
- EffHisto->GetAxis(4)->SetRange(1.,2.); // all
- TH1D *all = EffHisto->Projection(axis);
-
- EffHisto->GetAxis(4)->SetRange(2.,2.); // reconstructed
- TH1D *rec = EffHisto->Projection(axis);
- recc = (TH1D*)rec->Clone();
-
- if(recc)
- {
- recc->Divide(rec,all,1,1,"B");
- recc->GetYaxis()->SetTitle("efficiency");
- }
- }
- else if (type == 1) // Clone Rate
- {
- EffHisto->GetAxis(4)->SetRange(2.,2.); // reconstructed
- TH1D *all = WeightedProjection(EffHisto, axis, 3, axis_all);
-
- EffHisto->GetAxis(4)->SetRange(2.,2.); // reconstructed
- TH1D *clone = WeightedProjection(EffHisto, axis, 1, &axis_clone);
- recc = (TH1D*) clone->Clone();
-
- if(recc)
- {
- recc->Divide(clone,all,1,1,"B");
- recc->GetYaxis()->SetTitle("clone rate");
- }
- }
- else if (type == 2) // Fake Rate
- {
- EffHisto->GetAxis(4)->SetRange(2.,2.); // reconstructed
- TH1D *all = WeightedProjection(EffHisto, axis, 3, axis_all);
-
- EffHisto->GetAxis(4)->SetRange(2.,2.); // reconstructed
- TH1D *fake = WeightedProjection(EffHisto, axis, 1, &axis_fake);
- recc = (TH1D*) fake->Clone();
-
- if(recc)
- {
- recc->Divide(fake,all,1,1,"B");
- recc->GetYaxis()->SetTitle("fake rate");
- }
- }
-
- EffHisto->GetAxis(4)->SetRange(1,0); //Reset Range
-
- if (recc) { // Coverity fix
- recc->SetName(name);
-
- recc->GetXaxis()->SetTitle(fEffHisto->GetAxis(axis)->GetTitle());
-
- snprintf(title,256,"%s vs %s",vsTitle, fEffHisto->GetAxis(axis)->GetTitle());
- recc->SetTitle(title);
-
- if (axis == 2 ) recc->SetBit(TH1::kLogX);
- }
- return recc;
-}
git://git.uio.no / u / mrichter / AliRoot.git / blobdiff