+++ /dev/null
-// $Id$
-//
-// Jet deltaPt task.
-//
-// Author: S.Aiola
-
-#include <TClonesArray.h>
-#include <TH1F.h>
-#include <TH2F.h>
-#include <TH3F.h>
-#include <TList.h>
-#include <TLorentzVector.h>
-#include <TRandom3.h>
-
-#include "AliVCluster.h"
-#include "AliVParticle.h"
-#include "AliVTrack.h"
-#include "AliEmcalJet.h"
-#include "AliRhoParameter.h"
-#include "AliLog.h"
-#include "AliJetContainer.h"
-#include "AliParticleContainer.h"
-#include "AliClusterContainer.h"
-
-#include "AliAnalysisTaskDeltaPt.h"
-
-ClassImp(AliAnalysisTaskDeltaPt)
-
-//________________________________________________________________________
-AliAnalysisTaskDeltaPt::AliAnalysisTaskDeltaPt() :
- AliAnalysisTaskEmcalJet("AliAnalysisTaskDeltaPt", kTRUE),
- fMCJetPtThreshold(1),
- fMinRC2LJ(-1),
- fRCperEvent(-1),
- fConeRadius(0.2),
- fConeMinEta(-0.9),
- fConeMaxEta(0.9),
- fConeMinPhi(0),
- fConeMaxPhi(TMath::Pi()*2),
- fJetsCont(0),
- fTracksCont(0),
- fCaloClustersCont(0),
- fEmbJetsCont(0),
- fEmbTracksCont(0),
- fEmbCaloClustersCont(0),
- fRandTracksCont(0),
- fRandCaloClustersCont(0),
- fHistRhovsCent(0),
- fHistRCPhiEta(0),
- fHistRCPt(0),
- fHistRCPtExLJ(0),
- fHistRCPtExPartialLJ(0),
- fHistRCPtRand(0),
- fHistRhoVSRCPt(0),
- fHistDeltaPtRCvsEP(0),
- fHistDeltaPtRCExLJ(0),
- fHistDeltaPtRCExPartialLJ(0),
- fHistDeltaPtRCRand(0),
- fHistEmbJetsPtArea(0),
- fHistEmbJetsCorrPtArea(0),
- fHistEmbPartPtvsJetPt(0),
- fHistEmbPartPtvsJetCorrPt(0),
- fHistJetPtvsJetCorrPt(0),
- fHistDistLeadPart2JetAxis(0),
- fHistEmbBkgArea(0),
- fHistRhoVSEmbBkg(0),
- fHistDeltaPtEmbArea(0),
- fHistDeltaPtEmbvsEP(0),
- fHistRCPtExLJVSDPhiLJ(0),
- fHistRCPtExPartialLJVSDPhiLJ(0),
- fHistEmbJetsPhiEta(0),
- fHistLeadPartPhiEta(0)
-{
- // Default constructor.
-
- fHistRCPt = 0;
- fHistRCPtExLJ = 0;
- fHistRCPtExPartialLJ = 0;
- fHistRCPtRand = 0;
- fHistRhoVSRCPt = 0;
- fHistDeltaPtRCvsEP = 0;
- fHistDeltaPtRCExLJ = 0;
- fHistDeltaPtRCExPartialLJ = 0;
- fHistDeltaPtRCRand = 0;
- fHistEmbJetsPtArea = 0;
- fHistEmbJetsCorrPtArea = 0;
- fHistEmbPartPtvsJetPt = 0;
- fHistEmbPartPtvsJetCorrPt = 0;
- fHistJetPtvsJetCorrPt = 0;
- fHistDistLeadPart2JetAxis = 0;
- fHistEmbBkgArea = 0;
- fHistRhoVSEmbBkg = 0;
- fHistDeltaPtEmbArea = 0;
- fHistDeltaPtEmbvsEP = 0;
-
- SetMakeGeneralHistograms(kTRUE);
-}
-
-//________________________________________________________________________
-AliAnalysisTaskDeltaPt::AliAnalysisTaskDeltaPt(const char *name) :
- AliAnalysisTaskEmcalJet(name, kTRUE),
- fMCJetPtThreshold(1),
- fMinRC2LJ(-1),
- fRCperEvent(-1),
- fConeRadius(0.2),
- fConeMinEta(-0.9),
- fConeMaxEta(0.9),
- fConeMinPhi(0),
- fConeMaxPhi(TMath::Pi()*2),
- fJetsCont(0),
- fTracksCont(0),
- fCaloClustersCont(0),
- fEmbJetsCont(0),
- fEmbTracksCont(0),
- fEmbCaloClustersCont(0),
- fRandTracksCont(0),
- fRandCaloClustersCont(0),
- fHistRhovsCent(0),
- fHistRCPhiEta(0),
- fHistRCPt(0),
- fHistRCPtExLJ(0),
- fHistRCPtExPartialLJ(0),
- fHistRCPtRand(0),
- fHistRhoVSRCPt(0),
- fHistDeltaPtRCvsEP(0),
- fHistDeltaPtRCExLJ(0),
- fHistDeltaPtRCExPartialLJ(0),
- fHistDeltaPtRCRand(0),
- fHistEmbJetsPtArea(0),
- fHistEmbJetsCorrPtArea(0),
- fHistEmbPartPtvsJetPt(0),
- fHistEmbPartPtvsJetCorrPt(0),
- fHistJetPtvsJetCorrPt(0),
- fHistDistLeadPart2JetAxis(0),
- fHistEmbBkgArea(0),
- fHistRhoVSEmbBkg(0),
- fHistDeltaPtEmbArea(0),
- fHistDeltaPtEmbvsEP(0),
- fHistRCPtExLJVSDPhiLJ(0),
- fHistRCPtExPartialLJVSDPhiLJ(0),
- fHistEmbJetsPhiEta(0),
- fHistLeadPartPhiEta(0)
-{
- // Standard constructor.
-
- fHistRCPt = 0;
- fHistRCPtExLJ = 0;
- fHistRCPtExPartialLJ = 0;
- fHistRCPtRand = 0;
- fHistRhoVSRCPt = 0;
- fHistDeltaPtRCvsEP = 0;
- fHistDeltaPtRCExLJ = 0;
- fHistDeltaPtRCExPartialLJ = 0;
- fHistDeltaPtRCRand = 0;
- fHistEmbJetsPtArea = 0;
- fHistEmbJetsCorrPtArea = 0;
- fHistEmbPartPtvsJetPt = 0;
- fHistEmbPartPtvsJetCorrPt = 0;
- fHistJetPtvsJetCorrPt = 0;
- fHistDistLeadPart2JetAxis = 0;
- fHistEmbBkgArea = 0;
- fHistRhoVSEmbBkg = 0;
- fHistDeltaPtEmbArea = 0;
- fHistDeltaPtEmbvsEP = 0;
-
- SetMakeGeneralHistograms(kTRUE);
-}
-
-//________________________________________________________________________
-void AliAnalysisTaskDeltaPt::AllocateHistogramArrays()
-{
- fHistRCPt = new TH1*[fNcentBins];
- fHistRCPtExLJ = new TH1*[fNcentBins];
- fHistRCPtExPartialLJ = new TH1*[fNcentBins];
- fHistRCPtRand = new TH1*[fNcentBins];
- fHistRhoVSRCPt = new TH2*[fNcentBins];
- fHistDeltaPtRCvsEP = new TH2*[fNcentBins];
- fHistDeltaPtRCExLJ = new TH1*[fNcentBins];
- fHistDeltaPtRCExPartialLJ = new TH1*[fNcentBins];
- fHistDeltaPtRCRand = new TH1*[fNcentBins];
- fHistEmbJetsPtArea = new TH3*[fNcentBins];
- fHistEmbJetsCorrPtArea = new TH3*[fNcentBins];
- fHistEmbPartPtvsJetPt = new TH2*[fNcentBins];
- fHistEmbPartPtvsJetCorrPt = new TH2*[fNcentBins];
- fHistJetPtvsJetCorrPt = new TH2*[fNcentBins];
- fHistDistLeadPart2JetAxis = new TH1*[fNcentBins];
- fHistEmbBkgArea = new TH2*[fNcentBins];
- fHistRhoVSEmbBkg = new TH2*[fNcentBins];
- fHistDeltaPtEmbArea = new TH2*[fNcentBins];
- fHistDeltaPtEmbvsEP = new TH2*[fNcentBins];
-
- for (Int_t i = 0; i < fNcentBins; i++) {
- fHistRCPt[i] = 0;
- fHistRCPtExLJ[i] = 0;
- fHistRCPtExPartialLJ[i] = 0;
- fHistRCPtRand[i] = 0;
- fHistRhoVSRCPt[i] = 0;
- fHistDeltaPtRCvsEP[i] = 0;
- fHistDeltaPtRCExLJ[i] = 0;
- fHistDeltaPtRCExPartialLJ[i] = 0;
- fHistDeltaPtRCRand[i] = 0;
- fHistEmbJetsPtArea[i] = 0;
- fHistEmbJetsCorrPtArea[i] = 0;
- fHistEmbPartPtvsJetPt[i] = 0;
- fHistEmbPartPtvsJetCorrPt[i] = 0;
- fHistJetPtvsJetCorrPt[i] = 0;
- fHistDistLeadPart2JetAxis[i] = 0;
- fHistEmbBkgArea[i] = 0;
- fHistRhoVSEmbBkg[i] = 0;
- fHistDeltaPtEmbArea[i] = 0;
- fHistDeltaPtEmbvsEP[i] = 0;
- }
-}
-
-//________________________________________________________________________
-void AliAnalysisTaskDeltaPt::UserCreateOutputObjects()
-{
- // Create user output.
-
- AliAnalysisTaskEmcalJet::UserCreateOutputObjects();
-
- AllocateHistogramArrays();
-
- fHistRhovsCent = new TH2F("fHistRhovsCent", "fHistRhovsCent", 101, -1, 100, fNbins, 0, fMaxBinPt*2);
- fHistRhovsCent->GetXaxis()->SetTitle("Centrality (%)");
- fHistRhovsCent->GetYaxis()->SetTitle("#rho (GeV/c * rad^{-1})");
- fOutput->Add(fHistRhovsCent);
-
- fJetsCont = GetJetContainer("Jets");
- fTracksCont = GetParticleContainer("Tracks");
- fCaloClustersCont = GetClusterContainer("CaloClusters");
- fEmbJetsCont = GetJetContainer("EmbJets");
- fEmbTracksCont = GetParticleContainer("EmbTracks");
- fEmbCaloClustersCont = GetClusterContainer("EmbCaloClusters");
- fRandTracksCont = GetParticleContainer("RandTracks");
- fRandCaloClustersCont = GetClusterContainer("RandCaloClusters");
-
- if (fTracksCont || fCaloClustersCont) {
- fHistRCPhiEta = new TH2F("fHistRCPhiEta","fHistRCPhiEta", 100, -1, 1, 201, 0, TMath::Pi() * 2.01);
- fHistRCPhiEta->GetXaxis()->SetTitle("#eta");
- fHistRCPhiEta->GetYaxis()->SetTitle("#phi");
- fOutput->Add(fHistRCPhiEta);
-
- if (fJetsCont) {
- fHistRCPtExLJVSDPhiLJ = new TH2F("fHistRCPtExLJVSDPhiLJ","fHistRCPtExLJVSDPhiLJ", fNbins, fMinBinPt, fMaxBinPt, 128, -1.6, 4.8);
- fHistRCPtExLJVSDPhiLJ->GetXaxis()->SetTitle("#it{p}_{T} (GeV/#it{c})");
- fHistRCPtExLJVSDPhiLJ->GetYaxis()->SetTitle("#Delta#phi");
- fOutput->Add(fHistRCPtExLJVSDPhiLJ);
-
- fHistRCPtExPartialLJVSDPhiLJ = new TH2F("fHistRCPtExPartialLJVSDPhiLJ","fHistRCPtExPartialLJVSDPhiLJ", fNbins, fMinBinPt, fMaxBinPt, 128, -1.6, 4.8);
- fHistRCPtExPartialLJVSDPhiLJ->GetXaxis()->SetTitle("#it{p}_{T} (GeV/#it{c})");
- fHistRCPtExPartialLJVSDPhiLJ->GetYaxis()->SetTitle("#Delta#phi");
- fOutput->Add(fHistRCPtExPartialLJVSDPhiLJ);
- }
- }
-
- if (fEmbJetsCont) {
- fHistEmbJetsPhiEta = new TH2F("fHistEmbJetsPhiEta","fHistEmbJetsPhiEta", 100, -1, 1, 201, 0, TMath::Pi() * 2.01);
- fHistEmbJetsPhiEta->GetXaxis()->SetTitle("#eta");
- fHistEmbJetsPhiEta->GetYaxis()->SetTitle("#phi");
- fOutput->Add(fHistEmbJetsPhiEta);
-
- fHistLeadPartPhiEta = new TH2F("fHistLeadPartPhiEta","fHistLeadPartPhiEta", 100, -1, 1, 201, 0, TMath::Pi() * 2.01);
- fHistLeadPartPhiEta->GetXaxis()->SetTitle("#eta");
- fHistLeadPartPhiEta->GetYaxis()->SetTitle("#phi");
- fOutput->Add(fHistLeadPartPhiEta);
- }
-
- TString histname;
-
- const Int_t nbinsZ = 12;
- Double_t binsZ[nbinsZ+1] = {0,1,2,3,4,5,6,7,8,9,10,20,1000};
-
- Double_t *binsPt = GenerateFixedBinArray(fNbins, fMinBinPt, fMaxBinPt);
- Double_t *binsCorrPt = GenerateFixedBinArray(fNbins*2, -fMaxBinPt, fMaxBinPt);
- Double_t *binsArea = GenerateFixedBinArray(50, 0, 2);
-
- for (Int_t i = 0; i < fNcentBins; i++) {
- if (fTracksCont || fCaloClustersCont) {
- histname = "fHistRCPt_";
- histname += i;
- fHistRCPt[i] = new TH1F(histname.Data(), histname.Data(), fNbins, fMinBinPt, fMaxBinPt * 2);
- fHistRCPt[i]->GetXaxis()->SetTitle("#it{p}_{T} (GeV/#it{c})");
- fHistRCPt[i]->GetYaxis()->SetTitle("counts");
- fOutput->Add(fHistRCPt[i]);
-
- histname = "fHistRhoVSRCPt_";
- histname += i;
- fHistRhoVSRCPt[i] = new TH2F(histname.Data(), histname.Data(), fNbins, fMinBinPt, fMaxBinPt, fNbins, fMinBinPt, fMaxBinPt);
- fHistRhoVSRCPt[i]->GetXaxis()->SetTitle("A#rho (GeV/#it{c})");
- fHistRhoVSRCPt[i]->GetYaxis()->SetTitle("#it{p}_{T} (GeV/#it{c})");
- fOutput->Add(fHistRhoVSRCPt[i]);
-
- histname = "fHistDeltaPtRCvsEP_";
- histname += i;
- fHistDeltaPtRCvsEP[i] = new TH2F(histname.Data(), histname.Data(), 101, 0, TMath::Pi()*1.01, fNbins * 2, -fMaxBinPt, fMaxBinPt);
- fHistDeltaPtRCvsEP[i]->GetXaxis()->SetTitle("#phi_{RC} - #psi_{RP}");
- fHistDeltaPtRCvsEP[i]->GetYaxis()->SetTitle("#delta#it{p}_{T}^{RC} (GeV/#it{c})");
- fHistDeltaPtRCvsEP[i]->GetZaxis()->SetTitle("counts");
- fOutput->Add(fHistDeltaPtRCvsEP[i]);
-
- if (fJetsCont) {
- histname = "fHistRCPtExLJ_";
- histname += i;
- fHistRCPtExLJ[i] = new TH1F(histname.Data(), histname.Data(), fNbins, fMinBinPt, fMaxBinPt * 2);
- fHistRCPtExLJ[i]->GetXaxis()->SetTitle("#it{p}_{T}^{RC} (GeV/#it{c})");
- fHistRCPtExLJ[i]->GetYaxis()->SetTitle("counts");
- fOutput->Add(fHistRCPtExLJ[i]);
-
- histname = "fHistDeltaPtRCExLJ_";
- histname += i;
- fHistDeltaPtRCExLJ[i] = new TH1F(histname.Data(), histname.Data(), fNbins * 2, -fMaxBinPt, fMaxBinPt);
- fHistDeltaPtRCExLJ[i]->GetXaxis()->SetTitle("#delta#it{p}_{T}^{RC} (GeV/#it{c})");
- fHistDeltaPtRCExLJ[i]->GetYaxis()->SetTitle("counts");
- fOutput->Add(fHistDeltaPtRCExLJ[i]);
-
- histname = "fHistRCPtExPartialLJ_";
- histname += i;
- fHistRCPtExPartialLJ[i] = new TH1F(histname.Data(), histname.Data(), fNbins, fMinBinPt, fMaxBinPt * 2);
- fHistRCPtExPartialLJ[i]->GetXaxis()->SetTitle("#it{p}_{T}^{RC} (GeV/#it{c})");
- fHistRCPtExPartialLJ[i]->GetYaxis()->SetTitle("counts");
- fOutput->Add(fHistRCPtExPartialLJ[i]);
-
- histname = "fHistDeltaPtRCExPartialLJ_";
- histname += i;
- fHistDeltaPtRCExPartialLJ[i] = new TH1F(histname.Data(), histname.Data(), fNbins * 2, -fMaxBinPt, fMaxBinPt);
- fHistDeltaPtRCExPartialLJ[i]->GetXaxis()->SetTitle("#delta#it{p}_{T}^{RC} (GeV/#it{c})");
- fHistDeltaPtRCExPartialLJ[i]->GetYaxis()->SetTitle("counts");
- fOutput->Add(fHistDeltaPtRCExPartialLJ[i]);
- }
- }
-
- if (fRandTracksCont || fRandCaloClustersCont) {
- histname = "fHistRCPtRand_";
- histname += i;
- fHistRCPtRand[i] = new TH1F(histname.Data(), histname.Data(), fNbins, fMinBinPt, fMaxBinPt * 2);
- fHistRCPtRand[i]->GetXaxis()->SetTitle("#it{p}_{T}^{RC} (GeV/#it{c})");
- fHistRCPtRand[i]->GetYaxis()->SetTitle("counts");
- fOutput->Add(fHistRCPtRand[i]);
-
- histname = "fHistDeltaPtRCRand_";
- histname += i;
- fHistDeltaPtRCRand[i] = new TH1F(histname.Data(), histname.Data(), fNbins * 2, -fMaxBinPt, fMaxBinPt);
- fHistDeltaPtRCRand[i]->GetXaxis()->SetTitle("#delta#it{p}_{T}^{RC} (GeV/#it{c})");
- fHistDeltaPtRCRand[i]->GetYaxis()->SetTitle("counts");
- fOutput->Add(fHistDeltaPtRCRand[i]);
- }
-
- if (fEmbJetsCont) {
- histname = "fHistEmbJetsPtArea_";
- histname += i;
- fHistEmbJetsPtArea[i] = new TH3F(histname.Data(), histname.Data(), 50, binsArea, fNbins, binsPt, nbinsZ, binsZ);
- fHistEmbJetsPtArea[i]->GetXaxis()->SetTitle("area");
- fHistEmbJetsPtArea[i]->GetYaxis()->SetTitle("#it{p}_{T,jet}^{emb,raw} (GeV/#it{c})");
- fOutput->Add(fHistEmbJetsPtArea[i]);
-
- histname = "fHistEmbJetsCorrPtArea_";
- histname += i;
- fHistEmbJetsCorrPtArea[i] = new TH3F(histname.Data(), histname.Data(), 50, binsArea, fNbins * 2, binsCorrPt, nbinsZ, binsZ);
- fHistEmbJetsCorrPtArea[i]->GetXaxis()->SetTitle("area");
- fHistEmbJetsCorrPtArea[i]->GetYaxis()->SetTitle("#it{p}_{T,jet}^{emb,corr} (GeV/#it{c})");
- fOutput->Add(fHistEmbJetsCorrPtArea[i]);
-
- histname = "fHistEmbPartPtvsJetPt_";
- histname += i;
- fHistEmbPartPtvsJetPt[i] = new TH2F(histname.Data(), histname.Data(), fNbins, fMinBinPt, fMaxBinPt, fNbins, fMinBinPt, fMaxBinPt);
- fHistEmbPartPtvsJetPt[i]->GetXaxis()->SetTitle("#sum#it{p}_{T,const}^{emb} (GeV/#it{c})");
- fHistEmbPartPtvsJetPt[i]->GetYaxis()->SetTitle("#it{p}_{T,jet}^{emb} (GeV/#it{c})");
- fHistEmbPartPtvsJetPt[i]->GetZaxis()->SetTitle("counts");
- fOutput->Add(fHistEmbPartPtvsJetPt[i]);
-
- histname = "fHistEmbPartPtvsJetCorrPt_";
- histname += i;
- fHistEmbPartPtvsJetCorrPt[i] = new TH2F(histname.Data(), histname.Data(),
- fNbins, fMinBinPt, fMaxBinPt, fNbins*2, -fMaxBinPt, fMaxBinPt);
- fHistEmbPartPtvsJetCorrPt[i]->GetXaxis()->SetTitle("#sum#it{p}_{T,const}^{emb} (GeV/#it{c})");
- fHistEmbPartPtvsJetCorrPt[i]->GetYaxis()->SetTitle("#it{p}_{T,jet}^{emb} - A#rho (GeV/#it{c})");
- fHistEmbPartPtvsJetCorrPt[i]->GetZaxis()->SetTitle("counts");
- fOutput->Add(fHistEmbPartPtvsJetCorrPt[i]);
-
- histname = "fHistJetPtvsJetCorrPt_";
- histname += i;
- fHistJetPtvsJetCorrPt[i] = new TH2F(histname.Data(), histname.Data(),
- fNbins, fMinBinPt, fMaxBinPt, fNbins*2, -fMaxBinPt, fMaxBinPt);
- fHistJetPtvsJetCorrPt[i]->GetXaxis()->SetTitle("#it{p}_{T,jet}^{emb} (GeV/#it{c})");
- fHistJetPtvsJetCorrPt[i]->GetYaxis()->SetTitle("#it{p}_{T,jet}^{emb} - A#rho (GeV/#it{c})");
- fHistJetPtvsJetCorrPt[i]->GetZaxis()->SetTitle("counts");
- fOutput->Add(fHistJetPtvsJetCorrPt[i]);
-
- histname = "fHistDistLeadPart2JetAxis_";
- histname += i;
- fHistDistLeadPart2JetAxis[i] = new TH1F(histname.Data(), histname.Data(), 50, 0, 0.5);
- fHistDistLeadPart2JetAxis[i]->GetXaxis()->SetTitle("distance");
- fHistDistLeadPart2JetAxis[i]->GetYaxis()->SetTitle("counts");
- fOutput->Add(fHistDistLeadPart2JetAxis[i]);
-
- histname = "fHistEmbBkgArea_";
- histname += i;
- fHistEmbBkgArea[i] = new TH2F(histname.Data(), histname.Data(), 50, 0, 2, fNbins, fMinBinPt, fMaxBinPt);
- fHistEmbBkgArea[i]->GetXaxis()->SetTitle("area");
- fHistEmbBkgArea[i]->GetYaxis()->SetTitle("#it{p}_{T,jet}^{emb} - #sum#it{p}_{T,const}^{emb} (GeV/#it{c})");
- fOutput->Add(fHistEmbBkgArea[i]);
-
- histname = "fHistRhoVSEmbBkg_";
- histname += i;
- fHistRhoVSEmbBkg[i] = new TH2F(histname.Data(), histname.Data(), fNbins, fMinBinPt, fMaxBinPt, fNbins, fMinBinPt, fMaxBinPt);
- fHistRhoVSEmbBkg[i]->GetXaxis()->SetTitle("A#rho (GeV/#it{c})");
- fHistRhoVSEmbBkg[i]->GetYaxis()->SetTitle("#it{p}_{T,jet}^{emb} - #sum#it{p}_{T,const}^{emb} (GeV/#it{c})");
- fOutput->Add(fHistRhoVSEmbBkg[i]);
-
- histname = "fHistDeltaPtEmbArea_";
- histname += i;
- fHistDeltaPtEmbArea[i] = new TH2F(histname.Data(), histname.Data(),
- 50, 0, 2, fNbins * 2, -fMaxBinPt, fMaxBinPt);
- fHistDeltaPtEmbArea[i]->GetXaxis()->SetTitle("area");
- fHistDeltaPtEmbArea[i]->GetYaxis()->SetTitle("#delta#it{p}_{T}^{emb} (GeV/#it{c})");
- fHistDeltaPtEmbArea[i]->GetZaxis()->SetTitle("counts");
- fOutput->Add(fHistDeltaPtEmbArea[i]);
-
- histname = "fHistDeltaPtEmbvsEP_";
- histname += i;
- fHistDeltaPtEmbvsEP[i] = new TH2F(histname.Data(), histname.Data(), 101, 0, TMath::Pi()*1.01, fNbins * 2, -fMaxBinPt, fMaxBinPt);
- fHistDeltaPtEmbvsEP[i]->GetXaxis()->SetTitle("#phi_{jet} - #Psi_{EP}");
- fHistDeltaPtEmbvsEP[i]->GetYaxis()->SetTitle("#delta#it{p}_{T}^{emb} (GeV/#it{c})");
- fHistDeltaPtEmbvsEP[i]->GetZaxis()->SetTitle("counts");
- fOutput->Add(fHistDeltaPtEmbvsEP[i]);
- }
- }
-
- delete[] binsPt;
- delete[] binsCorrPt;
- delete[] binsArea;
-
- PostData(1, fOutput); // Post data for ALL output slots >0 here, to get at least an empty histogram
-}
-
-//________________________________________________________________________
-Bool_t AliAnalysisTaskDeltaPt::FillHistograms()
-{
- // Fill histograms.
-
- fHistRhovsCent->Fill(fCent, fRhoVal);
-
- // ************
- // Random cones
- // _________________________________
-
- const Float_t rcArea = fConeRadius * fConeRadius * TMath::Pi();
- Float_t RCpt = 0;
- Float_t RCeta = 0;
- Float_t RCphi = 0;
-
- if (fTracksCont || fCaloClustersCont) {
-
- for (Int_t i = 0; i < fRCperEvent; i++) {
- // Simple random cones
- RCpt = 0;
- RCeta = 0;
- RCphi = 0;
- GetRandomCone(RCpt, RCeta, RCphi, fTracksCont, fCaloClustersCont, 0);
- if (RCpt > 0) {
- fHistRCPhiEta->Fill(RCeta, RCphi);
- fHistRhoVSRCPt[fCentBin]->Fill(fRhoVal * rcArea, RCpt);
-
- fHistRCPt[fCentBin]->Fill(RCpt);
-
- Double_t ep = RCphi - fEPV0;
- while (ep < 0) ep += TMath::Pi();
- while (ep >= TMath::Pi()) ep -= TMath::Pi();
-
- fHistDeltaPtRCvsEP[fCentBin]->Fill(ep, RCpt - rcArea * fRhoVal);
- }
-
- if (fJetsCont) {
-
- // Random cones far from leading jet
- AliEmcalJet* jet = fJetsCont->GetLeadingJet("rho");
-
- RCpt = 0;
- RCeta = 0;
- RCphi = 0;
- GetRandomCone(RCpt, RCeta, RCphi, fTracksCont, fCaloClustersCont, jet);
- if (RCpt > 0) {
- if (jet) {
- Float_t dphi = RCphi - jet->Phi();
- if (dphi > 4.8) dphi -= TMath::Pi() * 2;
- if (dphi < -1.6) dphi += TMath::Pi() * 2;
- fHistRCPtExLJVSDPhiLJ->Fill(RCpt, dphi);
- }
- fHistRCPtExLJ[fCentBin]->Fill(RCpt);
- fHistDeltaPtRCExLJ[fCentBin]->Fill(RCpt - rcArea * fRhoVal);
- }
-
- //partial exclusion
- if(fBeamType == kpA) {
-
- RCpt = 0;
- RCeta = 0;
- RCphi = 0;
- GetRandomCone(RCpt, RCeta, RCphi, fTracksCont, fCaloClustersCont, jet, kTRUE);
-
- if (RCpt > 0) {
- if (jet) {
- Float_t dphi = RCphi - jet->Phi();
- if (dphi > 4.8) dphi -= TMath::Pi() * 2;
- if (dphi < -1.6) dphi += TMath::Pi() * 2;
- fHistRCPtExPartialLJVSDPhiLJ->Fill(RCpt, dphi);
- }
- fHistRCPtExPartialLJ[fCentBin]->Fill(RCpt);
- fHistDeltaPtRCExPartialLJ[fCentBin]->Fill(RCpt - rcArea * fRhoVal);
- }
- }
- }
- }
- }
-
- // Random cones with randomized particles
- if (fRandTracksCont || fRandCaloClustersCont) {
- RCpt = 0;
- RCeta = 0;
- RCphi = 0;
- GetRandomCone(RCpt, RCeta, RCphi, fRandTracksCont, fRandCaloClustersCont, 0);
- if (RCpt > 0) {
- fHistRCPtRand[fCentBin]->Fill(RCpt);
- fHistDeltaPtRCRand[fCentBin]->Fill(RCpt - rcArea * fRhoVal);
- }
- }
-
- // ************
- // Embedding
- // _________________________________
-
- if (fEmbJetsCont) {
-
- AliEmcalJet *embJet = NextEmbeddedJet(kTRUE);
-
- while (embJet != 0) {
- TLorentzVector mom;
- fEmbJetsCont->GetLeadingHadronMomentum(mom,embJet);
-
- Double_t distLeading2Jet = TMath::Sqrt((embJet->Eta() - mom.Eta()) * (embJet->Eta() - mom.Eta()) + (embJet->Phi() - mom.Phi()) * (embJet->Phi() - mom.Phi()));
-
- fHistEmbPartPtvsJetPt[fCentBin]->Fill(embJet->MCPt(), embJet->Pt());
- fHistEmbPartPtvsJetCorrPt[fCentBin]->Fill(embJet->MCPt(), embJet->Pt() - embJet->Area() * fRhoVal);
- fHistLeadPartPhiEta->Fill(mom.Eta(), mom.Phi());
- fHistDistLeadPart2JetAxis[fCentBin]->Fill(distLeading2Jet);
-
- fHistEmbJetsPtArea[fCentBin]->Fill(embJet->Area(), embJet->Pt(), mom.Pt());
- fHistEmbJetsCorrPtArea[fCentBin]->Fill(embJet->Area(), embJet->Pt() - fRhoVal * embJet->Area(), mom.Pt());
- fHistEmbJetsPhiEta->Fill(embJet->Eta(), embJet->Phi());
- fHistJetPtvsJetCorrPt[fCentBin]->Fill(embJet->Pt(), embJet->Pt() - fRhoVal * embJet->Area());
-
- fHistEmbBkgArea[fCentBin]->Fill(embJet->Area(), embJet->Pt() - embJet->MCPt());
- fHistRhoVSEmbBkg[fCentBin]->Fill(fRhoVal * embJet->Area(), embJet->Pt() - embJet->MCPt());
- fHistDeltaPtEmbArea[fCentBin]->Fill(embJet->Area(), embJet->Pt() - embJet->Area() * fRhoVal - embJet->MCPt());
-
- Double_t ep = embJet->Phi() - fEPV0;
- while (ep < 0) ep += TMath::Pi();
- while (ep >= TMath::Pi()) ep -= TMath::Pi();
-
- fHistDeltaPtEmbvsEP[fCentBin]->Fill(ep, embJet->Pt() - embJet->Area() * fRhoVal - embJet->MCPt());
-
- embJet = NextEmbeddedJet();
- }
- }
-
- return kTRUE;
-}
-
-//________________________________________________________________________
-AliEmcalJet* AliAnalysisTaskDeltaPt::NextEmbeddedJet(Bool_t reset)
-{
- // Get the next accepted embedded jet.
-
- Int_t i = reset ? 0 : -1;
-
- AliEmcalJet* jet = fEmbJetsCont->GetNextAcceptJet(i);
- while (jet && jet->MCPt() < fMCJetPtThreshold) jet = fEmbJetsCont->GetNextAcceptJet();
-
- return jet;
-}
-
-//________________________________________________________________________
-void AliAnalysisTaskDeltaPt::GetRandomCone(Float_t &pt, Float_t &eta, Float_t &phi,
- AliParticleContainer* tracks, AliClusterContainer* clusters,
- AliEmcalJet *jet, Bool_t bPartialExclusion) const
-{
- // Get rigid cone.
-
- eta = -999;
- phi = -999;
- pt = 0;
-
- if (!tracks && !clusters)
- return;
-
- Float_t LJeta = 999;
- Float_t LJphi = 999;
-
- if (jet) {
- LJeta = jet->Eta();
- LJphi = jet->Phi();
- }
-
- Float_t maxEta = fConeMaxEta;
- Float_t minEta = fConeMinEta;
- Float_t maxPhi = fConeMaxPhi;
- Float_t minPhi = fConeMinPhi;
-
- if (maxPhi > TMath::Pi() * 2) maxPhi = TMath::Pi() * 2;
- if (minPhi < 0) minPhi = 0;
-
- Float_t dLJ = 0;
- Int_t repeats = 0;
- Bool_t reject = kTRUE;
- do {
- eta = gRandom->Rndm() * (maxEta - minEta) + minEta;
- phi = gRandom->Rndm() * (maxPhi - minPhi) + minPhi;
- dLJ = TMath::Sqrt((LJeta - eta) * (LJeta - eta) + (LJphi - phi) * (LJphi - phi));
-
- if(bPartialExclusion) {
- reject = kFALSE;
-
- TRandom3 rnd;
- rnd.SetSeed(0);
-
- Double_t ncoll = GetNColl();
-
- Double_t prob = 0.;
- if(ncoll>0)
- prob = 1./ncoll;
-
- if(rnd.Rndm()<=prob) reject = kTRUE; //reject cone
- }
-
- repeats++;
- } while (dLJ < fMinRC2LJ && repeats < 999 && reject);
-
- if (repeats == 999) {
- AliWarning(Form("%s: Could not get random cone!", GetName()));
- return;
- }
-
- if (clusters) {
- AliVCluster* cluster = clusters->GetNextAcceptCluster(0);
- while (cluster) {
- TLorentzVector nPart;
- cluster->GetMomentum(nPart, const_cast<Double_t*>(fVertex));
-
- Float_t cluseta = nPart.Eta();
- Float_t clusphi = nPart.Phi();
-
- if (TMath::Abs(clusphi - phi) > TMath::Abs(clusphi - phi + 2 * TMath::Pi()))
- clusphi += 2 * TMath::Pi();
- if (TMath::Abs(clusphi - phi) > TMath::Abs(clusphi - phi - 2 * TMath::Pi()))
- clusphi -= 2 * TMath::Pi();
-
- Float_t d = TMath::Sqrt((cluseta - eta) * (cluseta - eta) + (clusphi - phi) * (clusphi - phi));
- if (d <= fConeRadius)
- pt += nPart.Pt();
-
- cluster = clusters->GetNextAcceptCluster();
- }
- }
-
- if (tracks) {
- AliVParticle* track = tracks->GetNextAcceptParticle(0);
- while(track) {
- Float_t tracketa = track->Eta();
- Float_t trackphi = track->Phi();
-
- if (TMath::Abs(trackphi - phi) > TMath::Abs(trackphi - phi + 2 * TMath::Pi()))
- trackphi += 2 * TMath::Pi();
- if (TMath::Abs(trackphi - phi) > TMath::Abs(trackphi - phi - 2 * TMath::Pi()))
- trackphi -= 2 * TMath::Pi();
-
- Float_t d = TMath::Sqrt((tracketa - eta) * (tracketa - eta) + (trackphi - phi) * (trackphi - phi));
- if (d <= fConeRadius)
- pt += track->Pt();
-
- track = tracks->GetNextAcceptParticle();
- }
- }
-}
-
-//________________________________________________________________________
-void AliAnalysisTaskDeltaPt::SetConeEtaPhiEMCAL()
-{
- // Set default cuts for full cones
-
- SetConeEtaLimits(-0.7+fConeRadius,0.7-fConeRadius);
- SetConePhiLimits(1.4+fConeRadius,TMath::Pi()-fConeRadius);
-}
-
-//________________________________________________________________________
-void AliAnalysisTaskDeltaPt::SetConeEtaPhiTPC()
-{
- // Set default cuts for charged cones
-
- SetConeEtaLimits(-0.9+fConeRadius, 0.9-fConeRadius);
- SetConePhiLimits(-10, 10);
-}
-
-//________________________________________________________________________
-void AliAnalysisTaskDeltaPt::ExecOnce()
-{
- // Initialize the analysis.
-
- AliAnalysisTaskEmcalJet::ExecOnce();
-
- if (fTracksCont && fTracksCont->GetArray() == 0) fTracksCont = 0;
- if (fCaloClustersCont && fCaloClustersCont->GetArray() == 0) fCaloClustersCont = 0;
- if (fEmbTracksCont && fEmbTracksCont->GetArray() == 0) fEmbTracksCont = 0;
- if (fEmbCaloClustersCont && fEmbCaloClustersCont->GetArray() == 0) fEmbCaloClustersCont = 0;
- if (fRandTracksCont && fRandTracksCont->GetArray() == 0) fRandTracksCont = 0;
- if (fRandCaloClustersCont && fRandCaloClustersCont->GetArray() == 0) fRandCaloClustersCont = 0;
- if (fJetsCont && fJetsCont->GetArray() == 0) fJetsCont = 0;
- if (fEmbJetsCont && fEmbJetsCont->GetArray() == 0) fEmbJetsCont = 0;
-
- if (fRCperEvent < 0) {
- Double_t area = (fConeMaxEta - fConeMinEta) * (fConeMaxPhi - fConeMinPhi);
- Double_t rcArea = TMath::Pi() * fConeRadius * fConeRadius;
- fRCperEvent = TMath::FloorNint(area / rcArea - 0.5);
- if (fRCperEvent == 0)
- fRCperEvent = 1;
- }
-
- if (fMinRC2LJ < 0)
- fMinRC2LJ = fConeRadius * 1.5;
-
- const Float_t maxDist = TMath::Max(fConeMaxPhi - fConeMinPhi, fConeMaxEta - fConeMinEta) / 2;
- if (fMinRC2LJ > maxDist) {
- AliWarning(Form("The parameter fMinRC2LJ = %f is too large for the considered acceptance. "
- "Will use fMinRC2LJ = %f", fMinRC2LJ, maxDist));
- fMinRC2LJ = maxDist;
- }
-}
-
-//________________________________________________________________________
-Double_t AliAnalysisTaskDeltaPt::GetNColl() const {
- // Get NColl - returns value of corresponding bin
- // only works for pA
- // values taken from V0A slicing https://twiki.cern.ch/twiki/bin/viewauth/ALICE/PACentStudies#Tables_with_centrality_bins_from
-
- if(fBeamType == kpA) {
-
- const Int_t nNCollBins = 7;
-
- Double_t centMin[nNCollBins] = {0.,5.,10.,20.,40.,60.,80.};
- Double_t centMax[nNCollBins] = {5.,10.,20.,40.,60.,80.,100.};
-
- Double_t nColl[nNCollBins] = {14.7,13.,11.7,9.38,6.49,3.96,1.52};
-
- for(Int_t i = 0; i<nNCollBins; i++) {
- if(fCent>=centMin[i] && fCent<centMax[i])
- return nColl[i];
- }
-
- return -1.;
- }
- else {
- AliWarning(Form("%s: Only works for pA analysis. Returning -1",GetName()));
- return -1.;
- }
-}
git://git.uio.no / u / mrichter / AliRoot.git / blobdiff