From b52f13298cdb3fa5c756990d5e4241134339bc0f Mon Sep 17 00:00:00 2001 From: mverweij Date: Tue, 7 Oct 2014 01:56:45 -0700 Subject: [PATCH] add Dcal dijet performance task (Rosi) --- PWGJE/CMakelibPWGJEEMCALJetTasks.pkg | 1 + .../AliAnalysisTaskDcalDijetPerf.cxx | 397 ++++++++++++++++++ .../UserTasks/AliAnalysisTaskDcalDijetPerf.h | 59 +++ .../macros/AddTaskDcalDijetPerf.C | 127 ++++++ PWGJE/PWGJEEMCALJetTasksLinkDef.h | 1 + 5 files changed, 585 insertions(+) create mode 100644 PWGJE/EMCALJetTasks/UserTasks/AliAnalysisTaskDcalDijetPerf.cxx create mode 100644 PWGJE/EMCALJetTasks/UserTasks/AliAnalysisTaskDcalDijetPerf.h create mode 100644 PWGJE/EMCALJetTasks/macros/AddTaskDcalDijetPerf.C diff --git a/PWGJE/CMakelibPWGJEEMCALJetTasks.pkg b/PWGJE/CMakelibPWGJEEMCALJetTasks.pkg index 1467606a7b5..dd6591d32b5 100644 --- a/PWGJE/CMakelibPWGJEEMCALJetTasks.pkg +++ b/PWGJE/CMakelibPWGJEEMCALJetTasks.pkg @@ -53,6 +53,7 @@ set ( SRCS EMCALJetTasks/AliJetTriggerSelectionTask.cxx EMCALJetTasks/UserTasks/AliAnalysisTaskCLQA.cxx EMCALJetTasks/UserTasks/AliAnalysisTaskChargedJetsPA.cxx + EMCALJetTasks/UserTasks/AliAnalysisTaskDcalDijetPerf.cxx EMCALJetTasks/UserTasks/AliAnalysisTaskDeltaPtJEmb.cxx EMCALJetTasks/UserTasks/AliAnalysisTaskEmcalBadCells.cxx EMCALJetTasks/UserTasks/AliAnalysisTaskEmcalDiJetBase.cxx diff --git a/PWGJE/EMCALJetTasks/UserTasks/AliAnalysisTaskDcalDijetPerf.cxx b/PWGJE/EMCALJetTasks/UserTasks/AliAnalysisTaskDcalDijetPerf.cxx new file mode 100644 index 00000000000..73b4a72c503 --- /dev/null +++ b/PWGJE/EMCALJetTasks/UserTasks/AliAnalysisTaskDcalDijetPerf.cxx @@ -0,0 +1,397 @@ +// $Id$ +// +// Dcal dijet performance task +// +// Author: R. Reed + +#include +#include +#include +#include +#include +#include +#include + +#include "AliVCluster.h" +#include "AliAODCaloCluster.h" +#include "AliESDCaloCluster.h" +#include "AliVTrack.h" +#include "AliEmcalJet.h" +#include "AliRhoParameter.h" +#include "AliLog.h" +#include "AliJetContainer.h" +#include "AliParticleContainer.h" +#include "AliClusterContainer.h" +#include "AliPicoTrack.h" + +#include "AliAnalysisTaskDcalDijetPerf.h" + +ClassImp(AliAnalysisTaskDcalDijetPerf) + +//________________________________________________________________________ +AliAnalysisTaskDcalDijetPerf::AliAnalysisTaskDcalDijetPerf() : + AliAnalysisTaskEmcalJet("AliAnalysisTaskDcalDijetPerf", kTRUE), + fHistTracksPt(0), + fHistTracksEtaPhi(0), + fHistClustersPt(0), + fHistLeadingJetPt(0), + fHistJetsPhiEta(0), + fHistJetsPtArea(0), + fHistJetsPtLeadHad(0), + fHistJetsCorrPtArea(0), + fHistJet1(0), + fHistJet1m(0), + fHistJet1nm(0), + fHistJet2(0), + fHistJet1to2(0), + fJetsCont(0), + fJetsCont2(0), + fTracksCont(0), + fCaloClustersCont(0) + +{ + // Default constructor. + + fHistTracksPt = new TH1*[fNcentBins]; + fHistTracksEtaPhi = new TH2*[fNcentBins]; + fHistClustersPt = new TH1*[fNcentBins]; + fHistLeadingJetPt = new TH1*[fNcentBins]; + fHistJetsPhiEta = new TH2*[fNcentBins]; + fHistJetsPtArea = new TH2*[fNcentBins]; + fHistJetsPtLeadHad = new TH2*[fNcentBins]; + fHistJetsCorrPtArea = new TH2*[fNcentBins]; + + for (Int_t i = 0; i < fNcentBins; i++) { + fHistTracksPt[i] = 0; + fHistTracksEtaPhi[i] = 0; + fHistClustersPt[i] = 0; + fHistLeadingJetPt[i] = 0; + fHistJetsPhiEta[i] = 0; + fHistJetsPtArea[i] = 0; + fHistJetsPtLeadHad[i] = 0; + fHistJetsCorrPtArea[i] = 0; + } + fHistJet1 = 0; + fHistJet1m = 0; + fHistJet1nm = 0; + fHistJet2 = 0; + fHistJet1to2 = 0; + SetMakeGeneralHistograms(kTRUE); +} + +//________________________________________________________________________ +AliAnalysisTaskDcalDijetPerf::AliAnalysisTaskDcalDijetPerf(const char *name) : + AliAnalysisTaskEmcalJet(name, kTRUE), + fHistTracksPt(0), + fHistTracksEtaPhi(0), + fHistClustersPt(0), + fHistLeadingJetPt(0), + fHistJetsPhiEta(0), + fHistJetsPtArea(0), + fHistJetsPtLeadHad(0), + fHistJetsCorrPtArea(0), + fHistJet1(0), + fHistJet1m(0), + fHistJet1nm(0), + fHistJet2(0), + fHistJet1to2(0), + fJetsCont(0), + fJetsCont2(0), + fTracksCont(0), + fCaloClustersCont(0) +{ + // Standard constructor. + + fHistTracksPt = new TH1*[fNcentBins]; + fHistTracksEtaPhi = new TH2*[fNcentBins]; + fHistClustersPt = new TH1*[fNcentBins]; + fHistLeadingJetPt = new TH1*[fNcentBins]; + fHistJetsPhiEta = new TH2*[fNcentBins]; + fHistJetsPtArea = new TH2*[fNcentBins]; + fHistJetsPtLeadHad = new TH2*[fNcentBins]; + fHistJetsCorrPtArea = new TH2*[fNcentBins]; + + for (Int_t i = 0; i < fNcentBins; i++) { + fHistTracksPt[i] = 0; + fHistTracksEtaPhi[i] = 0; + fHistClustersPt[i] = 0; + fHistLeadingJetPt[i] = 0; + fHistJetsPhiEta[i] = 0; + fHistJetsPtArea[i] = 0; + fHistJetsPtLeadHad[i] = 0; + fHistJetsCorrPtArea[i] = 0; + } + fHistJet1 = 0; + fHistJet1m = 0; + fHistJet1nm = 0; + fHistJet2 = 0; + fHistJet1to2 = 0; + + SetMakeGeneralHistograms(kTRUE); +} + +//________________________________________________________________________ +AliAnalysisTaskDcalDijetPerf::~AliAnalysisTaskDcalDijetPerf() +{ + // Destructor. +} + +//________________________________________________________________________ +void AliAnalysisTaskDcalDijetPerf::UserCreateOutputObjects() +{ + // Create user output. + + AliAnalysisTaskEmcalJet::UserCreateOutputObjects(); + + fJetsCont = GetJetContainer(0); + fJetsCont2 = GetJetContainer(1); + if(fJetsCont) { //get particles and clusters connected to jets + fTracksCont = fJetsCont->GetParticleContainer(); + fCaloClustersCont = fJetsCont->GetClusterContainer(); + } else { //no jets, just analysis tracks and clusters + fTracksCont = GetParticleContainer(0); + fCaloClustersCont = GetClusterContainer(0); + } + fTracksCont->SetClassName("AliVTrack"); + //fCaloClustersCont->SetClassName("AliVCluster"); + TString histname; + + for (Int_t i = 0; i < fNcentBins; i++) { + if (fParticleCollArray.GetEntriesFast()>0) { + histname = "fHistTracksPt_"; + histname += i; + fHistTracksPt[i] = new TH1F(histname.Data(), histname.Data(), fNbins / 2, fMinBinPt, fMaxBinPt / 2); + fHistTracksPt[i]->GetXaxis()->SetTitle("p_{T,track} (GeV/c)"); + fHistTracksPt[i]->GetYaxis()->SetTitle("counts"); + fOutput->Add(fHistTracksPt[i]); + histname = "fHistTracksEtaPhi_"; + histname += i; + fHistTracksEtaPhi[i] = new TH2F(histname.Data(), histname.Data(), fNbins, -0.7, 0.7, fNbins, 0, TMath::Pi()*2); + fHistTracksEtaPhi[i]->GetXaxis()->SetTitle("eta"); + fHistTracksEtaPhi[i]->GetYaxis()->SetTitle("phi"); + fOutput->Add(fHistTracksEtaPhi[i]); + } + + if (fClusterCollArray.GetEntriesFast()>0) { + histname = "fHistClustersPt_"; + histname += i; + fHistClustersPt[i] = new TH1F(histname.Data(), histname.Data(), fNbins / 2, fMinBinPt, fMaxBinPt / 2); + fHistClustersPt[i]->GetXaxis()->SetTitle("p_{T,clus} (GeV/c)"); + fHistClustersPt[i]->GetYaxis()->SetTitle("counts"); + fOutput->Add(fHistClustersPt[i]); + } + + if (fJetCollArray.GetEntriesFast()>0) { + histname = "fHistLeadingJetPt_"; + histname += i; + fHistLeadingJetPt[i] = new TH1F(histname.Data(), histname.Data(), fNbins, fMinBinPt, fMaxBinPt); + fHistLeadingJetPt[i]->GetXaxis()->SetTitle("p_{T}^{raw} (GeV/c)"); + fHistLeadingJetPt[i]->GetYaxis()->SetTitle("counts"); + fOutput->Add(fHistLeadingJetPt[i]); + + histname = "fHistJetsPhiEta_"; + histname += i; + fHistJetsPhiEta[i] = new TH2F(histname.Data(), histname.Data(), 50, -1, 1, 101, 0, TMath::Pi()*2 + TMath::Pi()/200); + fHistJetsPhiEta[i]->GetXaxis()->SetTitle("#eta"); + fHistJetsPhiEta[i]->GetYaxis()->SetTitle("#phi"); + fOutput->Add(fHistJetsPhiEta[i]); + + histname = "fHistJetsPtArea_"; + histname += i; + fHistJetsPtArea[i] = new TH2F(histname.Data(), histname.Data(), fNbins, fMinBinPt, fMaxBinPt, 30, 0, 3); + fHistJetsPtArea[i]->GetXaxis()->SetTitle("p_{T}^{raw} (GeV/c)"); + fHistJetsPtArea[i]->GetYaxis()->SetTitle("area"); + fOutput->Add(fHistJetsPtArea[i]); + + histname = "fHistJetsPtLeadHad_"; + histname += i; + fHistJetsPtLeadHad[i] = new TH2F(histname.Data(), histname.Data(), fNbins, fMinBinPt, fMaxBinPt, fNbins / 2, fMinBinPt, fMaxBinPt / 2); + fHistJetsPtLeadHad[i]->GetXaxis()->SetTitle("p_{T}^{raw} (GeV/c)"); + fHistJetsPtLeadHad[i]->GetYaxis()->SetTitle("p_{T,lead} (GeV/c)"); + fHistJetsPtLeadHad[i]->GetZaxis()->SetTitle("counts"); + fOutput->Add(fHistJetsPtLeadHad[i]); + + if (!(GetJetContainer()->GetRhoName().IsNull())) { + histname = "fHistJetsCorrPtArea_"; + histname += i; + fHistJetsCorrPtArea[i] = new TH2F(histname.Data(), histname.Data(), fNbins*2, -fMaxBinPt, fMaxBinPt, 30, 0, 3); + fHistJetsCorrPtArea[i]->GetXaxis()->SetTitle("p_{T}^{corr} [GeV/c]"); + fHistJetsCorrPtArea[i]->GetYaxis()->SetTitle("area"); + fOutput->Add(fHistJetsCorrPtArea[i]); + } + } + } + + Int_t nbins[] = {150,100,100,100,150}; + Double_t xmin[] = { 0,-0.7, 0, 0, 0}; + Double_t xmax[] = {150, 0.7,TMath::TwoPi(), 1, 150}; + fHistJet1 = new THnSparseF("Jets1Collection","Jets1Collection",5,nbins,xmin,xmax); + fOutput->Add(fHistJet1); + fHistJet1m = new THnSparseF("Jets1CollectionMatched","Jets1Collection",5,nbins,xmin,xmax); + fOutput->Add(fHistJet1m); + fHistJet1nm = new THnSparseF("Jets1CollectionNotMatched","Jets1Collection",5,nbins,xmin,xmax); + fOutput->Add(fHistJet1nm); + fHistJet2 = new THnSparseF("Jets2Collection","Jets2Collection",5,nbins,xmin,xmax); + fOutput->Add(fHistJet2); + + Int_t nbins2[] = {150,100,100,100,150,100,100,100,100}; + Double_t xmin2[] = {0,-0.7,0,0,0,-0.7,0,0,0}; + Double_t xmax2[] = {150,0.7,6.28,1,150,0.7,6.28,1,0.2}; + fHistJet1to2 = new THnSparseF("Jets1to2Collection","Jets1to2Collection",9,nbins2,xmin2,xmax2); + fOutput->Add(fHistJet1to2); + + PostData(1, fOutput); // Post data for ALL output slots > 0 here. +} + +//________________________________________________________________________ +Bool_t AliAnalysisTaskDcalDijetPerf::FillHistograms() +{ + // Fill histograms. + + if (fTracksCont) { + AliVTrack *track = static_cast(fTracksCont->GetNextAcceptParticle(0)); + while(track) { + fHistTracksPt[fCentBin]->Fill(track->Pt()); + fHistTracksEtaPhi[fCentBin]->Fill(track->Eta(),track->Phi()); + track = static_cast(fTracksCont->GetNextAcceptParticle()); + } + } + + if (fCaloClustersCont) { + AliVCluster *cluster = fCaloClustersCont->GetNextAcceptCluster(0); + while(cluster) { + TLorentzVector nPart; + cluster->GetMomentum(nPart, fVertex); + fHistClustersPt[fCentBin]->Fill(nPart.Pt()); + + cluster = fCaloClustersCont->GetNextAcceptCluster(); + } + } + + int N1 = 0; + if (fJetsCont) { + AliEmcalJet *jet = fJetsCont->GetNextAcceptJet(0); + while(jet) { + Float_t NEFpT = jet->Pt()*jet->NEF(); + N1++; + //cout<<"loop 1 jet 1 has label, pt,eta,phi,NEF = "<GetLabel()<<" "<Pt()<<" "<Eta()<<" "<Phi()<<" "<NEF()<<" and NEFpT "<Fill(jet->Pt(), jet->Area()); + fHistJetsPhiEta[fCentBin]->Fill(jet->Eta(), jet->Phi()); + Float_t ptLeading = fJetsCont->GetLeadingHadronPt(jet); + fHistJetsPtLeadHad[fCentBin]->Fill(jet->Pt(), ptLeading); + + //110 degrees in the azimuthal angle |eta|<0.7 + //EMCal Tower 0.0143 x 0.0143 + //16x16 = 0.2288 x 0.2288 -> R = 0.1144 + + Double_t jetarray[] = {jet->Pt(),jet->Eta(),jet->Phi(),jet->NEF(),NEFpT}; + //cout<<"Filling "<Pt()<<" "<Eta()<<" "<Phi()<<" "<NEF()<Fill(jetarray); + if (fHistJetsCorrPtArea[fCentBin]) { + Float_t corrPt = jet->Pt() - fJetsCont->GetRhoVal() * jet->Area(); + fHistJetsCorrPtArea[fCentBin]->Fill(corrPt, jet->Area()); + } + jet = fJetsCont->GetNextAcceptJet(); + } + + jet = fJetsCont->GetLeadingJet(); + if(jet) fHistLeadingJetPt[fCentBin]->Fill(jet->Pt()); + } + // cout<<"DONE LOOP 1"<GetNextAcceptJet(0); + while(jet){ + Float_t NEFpT = jet->Pt()*jet->NEF(); + //cout<<"loop 2 jet 2 has label, pt,eta,phi,NEF = "<GetLabel()<<" "<Pt()<<" "<Eta()<<" "<Phi()<<" "<NEF()<Pt(),jet->Eta(),jet->Phi(),jet->NEF(),NEFpT}; + //cout<<"Filling "<Pt()<<" "<Eta()<<" "<Phi()<<" "<NEF()<Fill(jetarray); + // cout<<" we have a jet! wiht pt = "<Pt()<GetNextAcceptJet(); + } + } + //<<"DONE LOOP 2"<GetNJets()<<" cont1 jets and "<GetNJets()<<" cont2 jets"<GetNextAcceptJet(0); + AliEmcalJet *jet2 = fJetsCont2->GetNextAcceptJet(0); + while(jet1){ + bool ismatched = false; + Float_t NEFpT1 = jet1->Pt()*jet1->NEF(); + Double_t jetarray1[] = {jet1->Pt(),jet1->Eta(),jet1->Phi(),jet1->NEF(),NEFpT1}; + while(jet2){ + N1N2++; + //cout<<"jet 1 has label, pt,eta,phi,NEF = "<GetLabel()<<" "<Pt()<<" "<Eta()<<" "<Phi()<<" "<NEF()<Pt()<<" "<Eta()<<" "<Phi()<<" "<NEF()<Eta()-jet2->Eta(); + Double_t dphi = RelativePhi(jet1->Phi(),jet2->Phi()); + Double_t deta2 = deta*deta; + Double_t dphi2 = dphi*dphi; + Double_t dR = pow(deta2+dphi2,0.5); + //cout<<"dR is "<Pt(),jet1->Eta(),jet1->Phi(),jet1->NEF(),jet2->Pt(),jet2->Eta(),jet2->Phi(),jet2->NEF(),dR}; + if (dR<0.2){ + N1N2m++; + fHistJet1to2->Fill(jetarray); + ismatched = true; + } + jet2 = fJetsCont2->GetNextAcceptJet(); + } + if (ismatched) + fHistJet1m->Fill(jetarray1); + else + fHistJet1nm->Fill(jetarray1); + jet1 = fJetsCont->GetNextAcceptJet(); + jet2 = fJetsCont2->GetNextAcceptJet(0); + } + } + + return kTRUE; +} + +//________________________________________________________________________ +Float_t AliAnalysisTaskDcalDijetPerf:: RelativePhi(Double_t mphi,Double_t vphi) const +{ + if (vphi < -1*TMath::Pi()) vphi += (2*TMath::Pi()); + else if (vphi > TMath::Pi()) vphi -= (2*TMath::Pi()); + if (mphi < -1*TMath::Pi()) mphi += (2*TMath::Pi()); + else if (mphi > TMath::Pi()) mphi -= (2*TMath::Pi()); + double dphi = mphi-vphi; + if (dphi < -1*TMath::Pi()) dphi += (2*TMath::Pi()); + else if (dphi > TMath::Pi()) dphi -= (2*TMath::Pi()); + + return dphi;//dphi in [-Pi, Pi] +} + + +//________________________________________________________________________ +void AliAnalysisTaskDcalDijetPerf::ExecOnce() { + + AliAnalysisTaskEmcalJet::ExecOnce(); + + if (fJetsCont && fJetsCont->GetArray() == 0) fJetsCont = 0; + if (fTracksCont && fTracksCont->GetArray() == 0) fTracksCont = 0; + if (fCaloClustersCont && fCaloClustersCont->GetArray() == 0) fCaloClustersCont = 0; + +} + +//________________________________________________________________________ +Bool_t AliAnalysisTaskDcalDijetPerf::Run() +{ + // Run analysis code here, if needed. It will be executed before FillHistograms(). + + return kTRUE; // If return kFALSE FillHistogram() will NOT be executed. +} + +//________________________________________________________________________ +void AliAnalysisTaskDcalDijetPerf::Terminate(Option_t *) +{ + // Called once at the end of the analysis. +} diff --git a/PWGJE/EMCALJetTasks/UserTasks/AliAnalysisTaskDcalDijetPerf.h b/PWGJE/EMCALJetTasks/UserTasks/AliAnalysisTaskDcalDijetPerf.h new file mode 100644 index 00000000000..457cc46ffc1 --- /dev/null +++ b/PWGJE/EMCALJetTasks/UserTasks/AliAnalysisTaskDcalDijetPerf.h @@ -0,0 +1,59 @@ +#ifndef ALIANALYSISTASKDCALDIJETPERF_H +#define ALIANALYSISTASKDCALDIJETPERF_H + +// $Id$ + +class TH1; +class TH2; +class TH3; +class THnSparse; +class AliJetContainer; +class AliParticleContainer; +class AliClusterContainer; + +#include "AliAnalysisTaskEmcalJet.h" + +class AliAnalysisTaskDcalDijetPerf : public AliAnalysisTaskEmcalJet { + public: + + AliAnalysisTaskDcalDijetPerf(); + AliAnalysisTaskDcalDijetPerf(const char *name); + virtual ~AliAnalysisTaskDcalDijetPerf(); + + void UserCreateOutputObjects(); + void Terminate(Option_t *option); + + protected: + Float_t RelativePhi(Double_t mphi,Double_t vphi) const; + void ExecOnce(); + Bool_t FillHistograms() ; + Bool_t Run() ; + + // General histograms + TH1 **fHistTracksPt; //!Track pt spectrum + TH2 **fHistTracksEtaPhi; //!Track eta phi + TH1 **fHistClustersPt; //!Cluster pt spectrum + TH1 **fHistLeadingJetPt; //!Leading jet pt spectrum + TH2 **fHistJetsPhiEta; //!Phi-Eta distribution of jets + TH2 **fHistJetsPtArea; //!Jet pt vs. area + TH2 **fHistJetsPtLeadHad; //!Jet pt vs. leading hadron + TH2 **fHistJetsCorrPtArea; //!Jet pt - bkg vs. area + + THnSparse *fHistJet1; //!jet collection 1 + THnSparse *fHistJet1m; //!jet collection 1 matched + THnSparse *fHistJet1nm; //!jet collection 1 unmatched + THnSparse *fHistJet2; //!jet collection 2 + THnSparse *fHistJet1to2; //!jet collection 1 and 2 + + AliJetContainer *fJetsCont; //!Jets + AliJetContainer *fJetsCont2; //!Jets + AliParticleContainer *fTracksCont; //!Tracks + AliClusterContainer *fCaloClustersCont; //!Clusters + + private: + AliAnalysisTaskDcalDijetPerf(const AliAnalysisTaskDcalDijetPerf&); // not implemented + AliAnalysisTaskDcalDijetPerf &operator=(const AliAnalysisTaskDcalDijetPerf&); // not implemented + + ClassDef(AliAnalysisTaskDcalDijetPerf, 1) +}; +#endif diff --git a/PWGJE/EMCALJetTasks/macros/AddTaskDcalDijetPerf.C b/PWGJE/EMCALJetTasks/macros/AddTaskDcalDijetPerf.C new file mode 100644 index 00000000000..9c0ec509fee --- /dev/null +++ b/PWGJE/EMCALJetTasks/macros/AddTaskDcalDijetPerf.C @@ -0,0 +1,127 @@ +// $Id$ + +AliAnalysisTaskDcalDijetPerf* AddTaskDcalDijetPerf( + const char *ntracks = "Tracks", + const char *nclusters = "CaloClusters", + const char *njets = "Jets", + const char *njets2 = "Jets2", + const char *nrho = "Rho", + Int_t nCentBins = 1, + Double_t jetradius = 0.2, + Double_t jetradius2 = 0.2, + Double_t jetptcut = 1, + Double_t jetareacut = 0.6, + const char *type = "TPC", + Int_t leadhadtype = 0, + const char *taskname = "AliAnalysisTaskDcalDijetPerf" +) +{ + // Get the pointer to the existing analysis manager via the static access method. + //============================================================================== + AliAnalysisManager *mgr = AliAnalysisManager::GetAnalysisManager(); + if (!mgr) + { + ::Error("AddTaskDcalDijetPerf", "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("AddTaskDcalDijetPerf", "This task requires an input event handler"); + return NULL; + } + + //------------------------------------------------------- + // Init the task and do settings + //------------------------------------------------------- + + TString name(taskname); + if (strcmp(njets,"")) { + name += "_"; + name += njets; + } + if (strcmp(njets,"")) { + name += "_"; + name += njets2; + } + if (strcmp(nrho,"")) { + name += "_"; + name += nrho; + } + if (strcmp(type,"")) { + name += "_"; + name += type; + } + + Printf("name: %s",name.Data()); + + AliAnalysisTaskDcalDijetPerf* jetTask = new AliAnalysisTaskDcalDijetPerf(name); + jetTask->SetCentRange(0.,100.); + jetTask->SetNCentBins(nCentBins); + + AliParticleContainer *trackCont = jetTask->AddParticleContainer(ntracks); + trackCont->SetClassName("AliVTrack"); + AliClusterContainer *clusterCont = jetTask->AddClusterContainer(nclusters); + + TString strType(type); + AliJetContainer *jetCont = jetTask->AddJetContainer(njets,strType,jetradius); + AliJetContainer *jetCont2 = jetTask->AddJetContainer(njets2,strType,jetradius2); + if(jetCont) { + jetCont->SetRhoName(nrho); + jetCont->ConnectParticleContainer(trackCont); + jetCont->ConnectClusterContainer(clusterCont); + //jetCont->SetZLeadingCut(0.98,0.98); + //jetCont->SetPercAreaCut(0.6); + jetCont->SetJetPtCut(jetptcut); + jetCont->SetLeadingHadronType(leadhadtype); + } + if(jetCont2) { + jetCont2->SetRhoName(nrho); + jetCont2->ConnectParticleContainer(trackCont); + jetCont2->ConnectClusterContainer(clusterCont); + //jetCont->SetZLeadingCut(0.98,0.98); + //jetCont->SetPercAreaCut(0.6); + jetCont2->SetJetPtCut(jetptcut); + jetCont2->SetLeadingHadronType(leadhadtype); + } + + //------------------------------------------------------- + // Final settings, pass to manager and set the containers + //------------------------------------------------------- + + mgr->AddTask(jetTask); + + // 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 (jetTask, 0, cinput1 ); + mgr->ConnectOutput (jetTask, 1, coutput1 ); + + return jetTask; +} + +AliAnalysisTaskDcalDijetPerf* AddTaskDcalDijetPerf( AliEmcalJetTask* jetFinderTask, + Int_t nCentBins = 1, + Double_t jetareacut = 0.6, + const char *type = "EMCAL", + Int_t leadhadtype = 0, + const char *taskname = "AliAnalysisTaskDcalDijetPerf" + ) +{ + const char* ntracks = jetFinderTask->GetTracksName(); + const char* nclusters = jetFinderTask->GetClusName(); + const char* njets = jetFinderTask->GetJetsName(); + const char* nrho = jetFinderTask->GetRhoName(); + Double_t jetradius = jetFinderTask->GetRadius(); + Double_t jetptcut = jetFinderTask->GetMinJetPt(); + + AliAnalysisTaskDcalDijetPerf* jetTask = AddTaskDcalDijetPerf(ntracks , nclusters, njets, nrho, nCentBins, jetradius, jetptcut, jetareacut, type, leadhadtype, taskname); + + return jetTask; +} diff --git a/PWGJE/PWGJEEMCALJetTasksLinkDef.h b/PWGJE/PWGJEEMCALJetTasksLinkDef.h index b96ee9cb4d6..e43cc4fdb8c 100644 --- a/PWGJE/PWGJEEMCALJetTasksLinkDef.h +++ b/PWGJE/PWGJEEMCALJetTasksLinkDef.h @@ -34,6 +34,7 @@ // user tasks #pragma link C++ class AliAnalysisTaskCLQA+; #pragma link C++ class AliAnalysisTaskChargedJetsPA; +#pragma link C++ class AliAnalysisTaskDcalDijetPerf; #pragma link C++ class AliAnalysisTaskDeltaPtJEmb+; #pragma link C++ class AliAnalysisTaskEmcalBadCells+; #pragma link C++ class AliAnalysisTaskEmcalDiJetBase+; -- 2.43.0