From eec7bbb07091054a96c092651e2d2e7e5c0fff15 Mon Sep 17 00:00:00 2001 From: saiola Date: Tue, 8 Jul 2014 20:19:10 -0400 Subject: [PATCH] New tracking efficiency task that uses the EMCal jet framework --- PWG/CMakelibPWGEMCAL.pkg | 1 + PWG/EMCAL/AliEmcalTrackingQATask.cxx | 399 +++++++++++++++++++++++++++ PWG/EMCAL/AliEmcalTrackingQATask.h | 50 ++++ PWG/EMCAL/macros/AddTaskTrackingQA.C | 48 ++++ PWG/PWGEMCALLinkDef.h | 1 + 5 files changed, 499 insertions(+) create mode 100644 PWG/EMCAL/AliEmcalTrackingQATask.cxx create mode 100644 PWG/EMCAL/AliEmcalTrackingQATask.h create mode 100644 PWG/EMCAL/macros/AddTaskTrackingQA.C diff --git a/PWG/CMakelibPWGEMCAL.pkg b/PWG/CMakelibPWGEMCAL.pkg index 5cadd245258..fe98a23dc47 100644 --- a/PWG/CMakelibPWGEMCAL.pkg +++ b/PWG/CMakelibPWGEMCAL.pkg @@ -47,6 +47,7 @@ set ( SRCS EMCAL/AliEmcalPicoTrackMaker.cxx EMCAL/AliEmcalSetupTask.cxx EMCAL/AliEmcalTenderTask.cxx + EMCAL/AliEmcalTrackingQATask.cxx EMCAL/AliEmcalTrackPropagatorTask.cxx EMCAL/AliEmcalTriggerMaker.cxx EMCAL/AliEmcalTriggerPatchInfo.cxx diff --git a/PWG/EMCAL/AliEmcalTrackingQATask.cxx b/PWG/EMCAL/AliEmcalTrackingQATask.cxx new file mode 100644 index 00000000000..f788b38dffb --- /dev/null +++ b/PWG/EMCAL/AliEmcalTrackingQATask.cxx @@ -0,0 +1,399 @@ +// Track QA task (efficiency and pt resolution) +// +// Author: S.Aiola + +#include +#include +#include +#include +#include + +#include "AliPicoTrack.h" +#include "AliAODMCParticle.h" +#include "AliParticleContainer.h" +#include "AliLog.h" + +#include "AliEmcalTrackingQATask.h" + +ClassImp(AliEmcalTrackingQATask) + +//________________________________________________________________________ +AliEmcalTrackingQATask::AliEmcalTrackingQATask() : + AliAnalysisTaskEmcal("AliEmcalTrackingQA", kTRUE), + fGeneratorLevel(0), + fDetectorLevel(0), + fTracksAll(0), + fTracksSelected(0), + fParticlesAllPhysPrim(0), + fParticlesSelected(0), + fParticlesFindable(0), + fParticlesMatched(0) +{ + // Default constructor. + + SetMakeGeneralHistograms(kTRUE); +} + +//________________________________________________________________________ +AliEmcalTrackingQATask::AliEmcalTrackingQATask(const char *name) : + AliAnalysisTaskEmcal(name, kTRUE), + fGeneratorLevel(0), + fDetectorLevel(0), + fTracksAll(0), + fTracksSelected(0), + fParticlesAllPhysPrim(0), + fParticlesSelected(0), + fParticlesFindable(0), + fParticlesMatched(0) +{ + // Standard constructor. + + SetMakeGeneralHistograms(kTRUE); +} + +//________________________________________________________________________ +AliEmcalTrackingQATask::~AliEmcalTrackingQATask() +{ + // Destructor. +} + +//________________________________________________________________________ +void AliEmcalTrackingQATask::UserCreateOutputObjects() +{ + // Create my user objects. + + AliAnalysisTaskEmcal::UserCreateOutputObjects(); + + if (!fCreateHisto) return; + + fTracksAll = new TH3**[fNcentBins]; + fTracksSelected = new TH3**[fNcentBins]; + fParticlesAllPhysPrim = new TH3*[fNcentBins]; + fParticlesSelected = new TH3*[fNcentBins]; + + TString histname; + + for (Int_t i = 0; i < fNcentBins; i++) { + + fTracksAll[i] = new TH3*[3]; + fTracksSelected[i] = new TH3*[3]; + for (Int_t j = 0; j < 3; j++) { + histname = Form("fTracksAll_%d_%d",i,j); + fTracksAll[i][j] = new TH3F(histname,histname, 100, -1, 1, 101, 0, TMath::Pi() * 2.02, fNbins, fMinBinPt, fMaxBinPt); + fTracksAll[i][j]->GetXaxis()->SetTitle("#eta"); + fTracksAll[i][j]->GetYaxis()->SetTitle("#phi"); + fTracksAll[i][j]->GetZaxis()->SetTitle("#it{p}_{T} (GeV/#it{c})"); + fOutput->Add(fTracksAll[i][j]); + + histname = Form("fTracksSelected_%d_%d",i,j); + fTracksSelected[i][j] = new TH3F(histname,histname, 100, -1, 1, 101, 0, TMath::Pi() * 2.02, fNbins, fMinBinPt, fMaxBinPt); + fTracksSelected[i][j]->GetXaxis()->SetTitle("#eta"); + fTracksSelected[i][j]->GetYaxis()->SetTitle("#phi"); + fTracksSelected[i][j]->GetZaxis()->SetTitle("#it{p}_{T} (GeV/#it{c})"); + fOutput->Add(fTracksSelected[i][j]); + } + + histname = Form("fParticlesAllPhysPrim_%d",i); + fParticlesAllPhysPrim[i] = new TH3F(histname,histname, 100, -1, 1, 101, 0, TMath::Pi() * 2.02, fNbins, fMinBinPt, fMaxBinPt); + fParticlesAllPhysPrim[i]->GetXaxis()->SetTitle("#eta"); + fParticlesAllPhysPrim[i]->GetYaxis()->SetTitle("#phi"); + fParticlesAllPhysPrim[i]->GetZaxis()->SetTitle("#it{p}_{T} (GeV/#it{c})"); + fOutput->Add(fParticlesAllPhysPrim[i]); + + histname = Form("fParticlesSelected_%d",i); + fParticlesSelected[i] = new TH3F(histname,histname, 100, -1, 1, 101, 0, TMath::Pi() * 2.02, fNbins, fMinBinPt, fMaxBinPt); + fParticlesSelected[i]->GetXaxis()->SetTitle("#eta"); + fParticlesSelected[i]->GetYaxis()->SetTitle("#phi"); + fParticlesSelected[i]->GetZaxis()->SetTitle("#it{p}_{T} (GeV/#it{c})"); + fOutput->Add(fParticlesSelected[i]); + } + + AllocateFindableParticlesTHnSparse(); + AllocateMatchedParticlesTHnSparse(); + + PostData(1, fOutput); +} + +//________________________________________________________________________ +void AliEmcalTrackingQATask::AllocateFindableParticlesTHnSparse() +{ + Int_t dim = 0; + TString title[20]; + Int_t nbins[20] = {0}; + Double_t min[20] = {0}; + Double_t max[20] = {0}; + + if (fForceBeamType != AliAnalysisTaskEmcal::kpp) { + title[dim] = "Centrality %"; + nbins[dim] = 101; + min[dim] = 0; + max[dim] = 101; + dim++; + } + + title[dim] = "#it{p}_{T} (GeV/#it{c})"; + nbins[dim] = fNbins; + min[dim] = fMinBinPt; + max[dim] = fMaxBinPt; + dim++; + + title[dim] = "#eta"; + nbins[dim] = 100; + min[dim] = -1; + max[dim] = 1; + dim++; + + title[dim] = "#phi"; + nbins[dim] = 101; + min[dim] = 0; + max[dim] = TMath::Pi() * 2.02; + dim++; + + fParticlesFindable = new THnSparseF("fParticlesFindable","fParticlesFindable",dim,nbins,min,max); + for (Int_t i = 0; i < dim; i++) + fParticlesFindable->GetAxis(i)->SetTitle(title[i]); + fOutput->Add(fParticlesFindable); +} + +//________________________________________________________________________ +void AliEmcalTrackingQATask::AllocateMatchedParticlesTHnSparse() +{ + Int_t dim = 0; + TString title[20]; + Int_t nbins[20] = {0}; + Double_t min[20] = {0}; + Double_t max[20] = {0}; + + if (fForceBeamType != AliAnalysisTaskEmcal::kpp) { + title[dim] = "Centrality %"; + nbins[dim] = 101; + min[dim] = 0; + max[dim] = 101; + dim++; + } + + title[dim] = "#it{p}_{T}^{gen} (GeV/#it{c})"; + nbins[dim] = fNbins; + min[dim] = fMinBinPt; + max[dim] = fMaxBinPt; + dim++; + + title[dim] = "#eta^{gen}"; + nbins[dim] = 100; + min[dim] = -1; + max[dim] = 1; + dim++; + + title[dim] = "#phi^{gen}"; + nbins[dim] = 101; + min[dim] = 0; + max[dim] = TMath::Pi() * 2.02; + dim++; + + title[dim] = "#it{p}_{T}^{det} (GeV/#it{c})"; + nbins[dim] = fNbins; + min[dim] = fMinBinPt; + max[dim] = fMaxBinPt; + dim++; + + title[dim] = "#eta^{det}"; + nbins[dim] = 100; + min[dim] = -1; + max[dim] = 1; + dim++; + + title[dim] = "#phi^{det}"; + nbins[dim] = 101; + min[dim] = 0; + max[dim] = TMath::Pi() * 2.02; + dim++; + + title[dim] = "(#it{p}_{T}^{gen} - #it{p}_{T}^{det}) / #it{p}_{T}^{gen}"; + nbins[dim] = fNbins; + min[dim] = -1; + max[dim] = 1; + dim++; + + title[dim] = "track type"; + nbins[dim] = 3; + min[dim] = -0.5; + max[dim] = 2.5; + dim++; + + fParticlesMatched = new THnSparseF("fParticlesMatched","fParticlesMatched",dim,nbins,min,max); + for (Int_t i = 0; i < dim; i++) + fParticlesMatched->GetAxis(i)->SetTitle(title[i]); + fOutput->Add(fParticlesMatched); +} + +//________________________________________________________________________ +void AliEmcalTrackingQATask::SetGeneratorLevelName(const char* name) +{ + if (!fGeneratorLevel) { // first check if the generator level array is set + fGeneratorLevel = static_cast(fParticleCollArray.At(0)); + if (fGeneratorLevel) { // now check if the first collection array has been added already + fGeneratorLevel->SetArrayName(name); + } + else { + fGeneratorLevel = AddParticleContainer(name); + } + fGeneratorLevel->SetClassName("AliAODMCParticle"); + fGeneratorLevel->SelectPhysicalPrimaries(kTRUE); + } + fGeneratorLevel->SetArrayName(name); +} + +//________________________________________________________________________ +void AliEmcalTrackingQATask::SetDetectorLevelName(const char* name) +{ + if (!fGeneratorLevel) { + AliError("Please, first set the generatol level array!"); + return; + } + if (!fDetectorLevel) { // first check if the detector level array is set + fDetectorLevel = static_cast(fParticleCollArray.At(1)); + if (fDetectorLevel) { // now check if the second collection array has been added already + fDetectorLevel->SetArrayName(name); + } + else { + fDetectorLevel = AddParticleContainer(name); + } + fDetectorLevel->SetClassName("AliPicoTrack"); + fDetectorLevel->SelectPhysicalPrimaries(kTRUE); + } + fDetectorLevel->SetArrayName(name); +} + +//________________________________________________________________________ +void AliEmcalTrackingQATask::ExecOnce() +{ + // Init the analysis. + + if (fParticleCollArray.GetEntriesFast() < 2) { + AliFatal("This task needs at least two particle containers!"); + } + + if (!fGeneratorLevel) { + fGeneratorLevel = static_cast(fParticleCollArray.At(0)); + fGeneratorLevel->SetClassName("AliAODMCParticle"); + } + + if (!fDetectorLevel) { + fDetectorLevel = static_cast(fParticleCollArray.At(1)); + fDetectorLevel->SetClassName("AliPicoTrack"); + } + + AliAnalysisTaskEmcal::ExecOnce(); +} + +//________________________________________________________________________ +void AliEmcalTrackingQATask::FillFindableParticlesTHnSparse(Double_t cent, Double_t partEta, Double_t partPhi, Double_t partPt) +{ + Double_t contents[20]={0}; + + for (Int_t i = 0; i < fParticlesFindable->GetNdimensions(); i++) { + TString title(fParticlesFindable->GetAxis(i)->GetTitle()); + if (title=="Centrality %") + contents[i] = cent; + else if (title=="#it{p}_{T} (GeV/#it{c})") + contents[i] = partPt; + else if (title=="#eta") + contents[i] = partEta; + else if (title=="#phi") + contents[i] = partPhi; + else + AliWarning(Form("Unable to fill dimension %s of histogram %s!", title.Data(), fParticlesFindable->GetName())); + } + + fParticlesFindable->Fill(contents); +} + +//________________________________________________________________________ +void AliEmcalTrackingQATask::FillMatchedParticlesTHnSparse(Double_t cent, Double_t partEta, Double_t partPhi, Double_t partPt, + Double_t trackEta, Double_t trackPhi, Double_t trackPt, Byte_t trackType) +{ + Double_t contents[20]={0}; + + for (Int_t i = 0; i < fParticlesMatched->GetNdimensions(); i++) { + TString title(fParticlesMatched->GetAxis(i)->GetTitle()); + if (title=="Centrality %") + contents[i] = cent; + else if (title=="#it{p}_{T}^{gen} (GeV/#it{c})") + contents[i] = partPt; + else if (title=="#eta^{gen}") + contents[i] = partEta; + else if (title=="#phi^{gen}") + contents[i] = partPhi; + else if (title=="#it{p}_{T}^{det} (GeV/#it{c})") + contents[i] = trackPt; + else if (title=="#eta^{det}") + contents[i] = trackEta; + else if (title=="#phi^{det}") + contents[i] = trackPhi; + else if (title=="(#it{p}_{T}^{gen} - #it{p}_{T}^{det}) / #it{p}_{T}^{gen}") + contents[i] = (partPt - trackPt) / partPt; + else if (title=="track type") + contents[i] = trackType; + else + AliWarning(Form("Unable to fill dimension %s of histogram %s!", title.Data(), fParticlesMatched->GetName())); + } + + fParticlesMatched->Fill(contents); +} + +//________________________________________________________________________ +Bool_t AliEmcalTrackingQATask::FillHistograms() +{ + // Fill the histograms. + + AliPicoTrack *track = static_cast(fDetectorLevel->GetNextAcceptParticle(0)); + while (track != 0) { + Byte_t type = track->GetTrackType(); + if (type<= 2 && type >= 0) { + fTracksAll[fCentBin][type]->Fill(track->Eta(), track->Phi(), track->Pt()); + + Int_t label = TMath::Abs(track->GetLabel()); + + if (label==0 || track->GetGeneratorIndex() == 0) { // reject particles generated from other generator in the cocktail but keep fake tracks (label == 0) + fTracksSelected[fCentBin][type]->Fill(track->Eta(), track->Phi(), track->Pt()); + } + + if (label > 0) { + AliAODMCParticle *part = static_cast(fGeneratorLevel->GetAcceptParticleWithLabel(label)); + if (part) { + if (part->GetGeneratorIndex() == 0) { + Int_t pdg = TMath::Abs(part->PdgCode()); + // select charged pions, protons, kaons , electrons, muons + if (pdg == 211 || pdg == 2212 || pdg == 321 || pdg == 11 || pdg == 13) { + FillMatchedParticlesTHnSparse(fCent, part->Eta(), part->Phi(), part->Pt(), track->Eta(), track->Phi(), track->Pt(), type); + } + } + } + } + } + else { + AliError(Form("Track %d has type %d not recognized!", fDetectorLevel->GetCurrentID(), type)); + } + + track = static_cast(fDetectorLevel->GetNextAcceptParticle()); + } + + AliAODMCParticle *part = static_cast(fGeneratorLevel->GetNextAcceptParticle(0)); + while (part != 0) { + fParticlesAllPhysPrim[fCentBin]->Fill(part->Eta(), part->Phi(), part->Pt()); + + if (part->GetGeneratorIndex() == 0) { + fParticlesSelected[fCentBin]->Fill(part->Eta(), part->Phi(), part->Pt()); + + Int_t pdg = TMath::Abs(part->PdgCode()); + // select charged pions, protons, kaons , electrons, muons + if (pdg == 211 || pdg == 2212 || pdg == 321 || pdg == 11 || pdg == 13) { + FillFindableParticlesTHnSparse(fCent, part->Eta(), part->Phi(), part->Pt()); + } + } + + part = static_cast(fGeneratorLevel->GetNextAcceptParticle()); + } + + return kTRUE; +} diff --git a/PWG/EMCAL/AliEmcalTrackingQATask.h b/PWG/EMCAL/AliEmcalTrackingQATask.h new file mode 100644 index 00000000000..50e108a7e41 --- /dev/null +++ b/PWG/EMCAL/AliEmcalTrackingQATask.h @@ -0,0 +1,50 @@ +#ifndef ALIEMCALTRACKINGQATASK_H +#define ALIEMCALTRACKINGQATASK_H + +#include "AliAnalysisTaskEmcal.h" + +class AliParticleContainer; +class THnSparse; +class TH3; + +class AliEmcalTrackingQATask : public AliAnalysisTaskEmcal { + + public: + AliEmcalTrackingQATask(); + AliEmcalTrackingQATask(const char *name); + virtual ~AliEmcalTrackingQATask(); + + void UserCreateOutputObjects(); + void SetGeneratorLevelName(const char* name); + void SetDetectorLevelName(const char* name); + + protected: + Bool_t FillHistograms() ; + void ExecOnce() ; + void AllocateFindableParticlesTHnSparse() ; + void AllocateMatchedParticlesTHnSparse() ; + void FillFindableParticlesTHnSparse(Double_t cent, Double_t partEta, Double_t partPhi, Double_t partPt); + void FillMatchedParticlesTHnSparse(Double_t cent, Double_t partEta, Double_t partPhi, Double_t partPt, + Double_t trackEta, Double_t trackPhi, Double_t trackPt, Byte_t trackType); + + // Task configuration + + // Service fields (non-streamed) + AliParticleContainer* fGeneratorLevel ; //! generator level container + AliParticleContainer* fDetectorLevel ; //! detector level container + + // Histograms + TH3*** fTracksAll ; //! all tracks + TH3*** fTracksSelected ; //! selected tracks (e.g. remove injected signal in HIJING productions) + TH3** fParticlesAllPhysPrim ; //! all physical primary particles + TH3** fParticlesSelected ; //! selected physical primary particles (e.g. remove injected signal in HIJING productions) + THnSparse* fParticlesFindable ; //! findable physical primary particles (use PDG and charge selection) + THnSparse* fParticlesMatched ; //! primary particles matched to detector level tracks + + private: + AliEmcalTrackingQATask(const AliEmcalTrackingQATask&); // not implemented + AliEmcalTrackingQATask &operator=(const AliEmcalTrackingQATask&); // not implemented + + ClassDef(AliEmcalTrackingQATask, 1) // Track QA task (efficiency and pt resolution) +}; +#endif diff --git a/PWG/EMCAL/macros/AddTaskTrackingQA.C b/PWG/EMCAL/macros/AddTaskTrackingQA.C new file mode 100644 index 00000000000..6e86176f65b --- /dev/null +++ b/PWG/EMCAL/macros/AddTaskTrackingQA.C @@ -0,0 +1,48 @@ +AliEmcalTrackingQATask* AddTaskTrackingQA(const char *nGenLev = "mcparticles", + const char *nDetLev = "PicoTracks") +{ + // Get the pointer to the existing analysis manager via the static access method. + //============================================================================== + AliAnalysisManager *mgr = AliAnalysisManager::GetAnalysisManager(); + if (!mgr) + { + ::Error("AddTaskTrackingQA", "No analysis manager to connect to."); + return NULL; + } + + // Check the analysis type using the event handlers connected to the analysis manager. + //============================================================================== + if (!mgr->GetInputEventHandler()) + { + ::Error("AddTaskTrackingQA", "This task requires an input event handler"); + return NULL; + } + + //------------------------------------------------------- + // Init the task and do settings + //------------------------------------------------------- + + TString name(Form("AliEmcalTrackingQATask_%s_%s", nGenLev, nDetLev)); + AliEmcalTrackingQATask *qaTask = new AliEmcalTrackingQATask(name); + qaTask->SetGeneratorLevelName(nGenLev); + qaTask->SetDetectorLevelName(nDetLev); + + //------------------------------------------------------- + // Final settings, pass to manager and set the containers + //------------------------------------------------------- + + mgr->AddTask(qaTask); + + // Create containers for input/output + AliAnalysisDataContainer *cinput1 = mgr->GetCommonInputContainer() ; + + TString contName(name); + contName += "_histos"; + AliAnalysisDataContainer *coutput1 = mgr->CreateContainer(contName.Data(), + TList::Class(),AliAnalysisManager::kOutputContainer, + Form("%s", AliAnalysisManager::GetCommonFileName())); + mgr->ConnectInput (qaTask, 0, cinput1 ); + mgr->ConnectOutput (qaTask, 1, coutput1 ); + + return qaTask; +} diff --git a/PWG/PWGEMCALLinkDef.h b/PWG/PWGEMCALLinkDef.h index fe58f6d58be..9b4c490b3a1 100644 --- a/PWG/PWGEMCALLinkDef.h +++ b/PWG/PWGEMCALLinkDef.h @@ -24,6 +24,7 @@ #pragma link C++ class AliEmcalPicoTrackMaker+; #pragma link C++ class AliEmcalSetupTask+; #pragma link C++ class AliEmcalTenderTask+; +#pragma link C++ class AliEmcalTrackingQATask+; #pragma link C++ class AliEmcalTrackPropagatorTask+; #pragma link C++ class AliEmcalTriggerMaker+; #pragma link C++ class AliEmcalTriggerPatchInfo+; -- 2.43.5