#include <TClonesArray.h>
#include <TH1F.h>
#include <TH2F.h>
+#include <TH3F.h>
#include <TList.h>
#include <TLorentzVector.h>
#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 "AliAnalysisTaskEmcalJetSample.h"
//________________________________________________________________________
AliAnalysisTaskEmcalJetSample::AliAnalysisTaskEmcalJetSample() :
- AliAnalysisTaskEmcalJet("AliAnalysisTaskEmcalJetSample", kTRUE)
+ AliAnalysisTaskEmcalJet("AliAnalysisTaskEmcalJetSample", kTRUE),
+ fHistTracksPt(0),
+ fHistClustersPt(0),
+ fHistLeadingJetPt(0),
+ fHistJetsPhiEta(0),
+ fHistJetsPtArea(0),
+ fHistJetsPtLeadHad(0),
+ fHistJetsCorrPtArea(0),
+ fHistPtDEtaDPhiTrackClus(0),
+ fHistPtDEtaDPhiClusTrack(0),
+ fHistClustDx(0),
+ fHistClustDz(0),
+ fJetsCont(0),
+ fTracksCont(0),
+ fCaloClustersCont(0)
{
// Default constructor.
- for (Int_t i = 0; i < 4; i++) {
+ fHistTracksPt = new TH1*[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;
fHistClustersPt[i] = 0;
fHistLeadingJetPt[i] = 0;
//________________________________________________________________________
AliAnalysisTaskEmcalJetSample::AliAnalysisTaskEmcalJetSample(const char *name) :
- AliAnalysisTaskEmcalJet(name, kTRUE)
+ AliAnalysisTaskEmcalJet(name, kTRUE),
+ fHistTracksPt(0),
+ fHistClustersPt(0),
+ fHistLeadingJetPt(0),
+ fHistJetsPhiEta(0),
+ fHistJetsPtArea(0),
+ fHistJetsPtLeadHad(0),
+ fHistJetsCorrPtArea(0),
+ fHistPtDEtaDPhiTrackClus(0),
+ fHistPtDEtaDPhiClusTrack(0),
+ fHistClustDx(0),
+ fHistClustDz(0),
+ fJetsCont(0),
+ fTracksCont(0),
+ fCaloClustersCont(0)
{
// Standard constructor.
- for (Int_t i = 0; i < 4; i++) {
+ fHistTracksPt = new TH1*[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;
fHistClustersPt[i] = 0;
fHistLeadingJetPt[i] = 0;
AliAnalysisTaskEmcalJet::UserCreateOutputObjects();
+ fJetsCont = GetJetContainer(0);
+ 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);
+ }
+ if(fTracksCont) fTracksCont->SetClassName("AliVTrack");
+ if(fCaloClustersCont) fCaloClustersCont->SetClassName("AliVCluster");
+
TString histname;
- for (Int_t i = 0; i < 4; i++) {
- if (!fTracksName.IsNull()) {
+ 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);
fOutput->Add(fHistTracksPt[i]);
}
- if (!fCaloName.IsNull()) {
+ if (fClusterCollArray.GetEntriesFast()>0) {
histname = "fHistClustersPt_";
histname += i;
fHistClustersPt[i] = new TH1F(histname.Data(), histname.Data(), fNbins / 2, fMinBinPt, fMaxBinPt / 2);
fOutput->Add(fHistClustersPt[i]);
}
- if (!fJetsName.IsNull()) {
+ if (fJetCollArray.GetEntriesFast()>0) {
histname = "fHistLeadingJetPt_";
histname += i;
fHistLeadingJetPt[i] = new TH1F(histname.Data(), histname.Data(), fNbins, fMinBinPt, fMaxBinPt);
fHistJetsPtLeadHad[i]->GetZaxis()->SetTitle("counts");
fOutput->Add(fHistJetsPtLeadHad[i]);
- if (!fRhoName.IsNull()) {
+ if (!(GetJetContainer()->GetRhoName().IsNull())) {
histname = "fHistJetsCorrPtArea_";
histname += i;
fHistJetsCorrPtArea[i] = new TH2F(histname.Data(), histname.Data(), fNbins*2, -fMaxBinPt, fMaxBinPt, 30, 0, 3);
}
}
}
+
+ histname = "fHistPtDEtaDPhiTrackClus";
+ fHistPtDEtaDPhiTrackClus = new TH3F(histname.Data(),Form("%s;#it{p}_{T}^{track};#Delta#eta;#Delta#varphi",histname.Data()),100,0.,100.,100,-0.1,0.1,100,-0.1,0.1);
+ fOutput->Add(fHistPtDEtaDPhiTrackClus);
+
+ histname = "fHistPtDEtaDPhiClusTrack";
+ fHistPtDEtaDPhiClusTrack = new TH3F(histname.Data(),Form("%s;#it{p}_{T}^{clus};#Delta#eta;#Delta#varphi",histname.Data()),100,0.,100.,100,-0.1,0.1,100,-0.1,0.1);
+ fOutput->Add(fHistPtDEtaDPhiClusTrack);
+
+ fHistClustDx = new TH1F("fHistClustDx","fHistClustDx;Dx",1000,0.,1.);
+ fOutput->Add(fHistClustDx);
+
+ fHistClustDz = new TH1F("fHistClustDz","fHistClustDz;Dz",1000,0.,1.);
+ fOutput->Add(fHistClustDz);
+
PostData(1, fOutput); // Post data for ALL output slots > 0 here.
}
{
// Fill histograms.
- if (fTracks) {
- const Int_t ntracks = fTracks->GetEntriesFast();
-
- for (Int_t it = 0; it < ntracks; it++) {
- AliVTrack *track = static_cast<AliVTrack*>(fTracks->At(it));
-
- if (!track) {
- AliError(Form("Could not receive track %d", it));
- continue;
- }
-
- if (!AcceptTrack(track))
- continue;
-
+ if (fTracksCont) {
+ AliVTrack *track = static_cast<AliVTrack*>(fTracksCont->GetNextAcceptParticle(0));
+ while(track) {
fHistTracksPt[fCentBin]->Fill(track->Pt());
+ track = static_cast<AliVTrack*>(fTracksCont->GetNextAcceptParticle());
}
}
-
- if (fCaloClusters) {
- const Int_t nclusters = fCaloClusters->GetEntriesFast();
-
- for (Int_t ic = 0; ic < nclusters; ic++) {
- AliVCluster *cluster = static_cast<AliVCluster*>(fCaloClusters->At(ic));
-
- if (!cluster) {
- AliError(Form("Could not receive cluster %d", ic));
- continue;
- }
+ if (fCaloClustersCont) {
+ AliVCluster *cluster = fCaloClustersCont->GetNextAcceptCluster(0);
+ while(cluster) {
TLorentzVector nPart;
cluster->GetMomentum(nPart, fVertex);
fHistClustersPt[fCentBin]->Fill(nPart.Pt());
+ Double_t dx = cluster->GetTrackDx();
+ Double_t dz = cluster->GetTrackDz();
+ fHistClustDx->Fill(dx);
+ fHistClustDz->Fill(dz);
+ cluster = fCaloClustersCont->GetNextAcceptCluster();
}
}
- if (fJets) {
- static Int_t sortedJets[9999] = {-1};
- GetSortedArray(sortedJets, fJets);
-
- if (sortedJets[0]>=0) {
- AliEmcalJet* leadJet = static_cast<AliEmcalJet*>(fJets->At(sortedJets[0]));
- if (leadJet)
- fHistLeadingJetPt[fCentBin]->Fill(leadJet->Pt());
- else
- AliError("Could not retrieve leading jet!");
- }
-
- const Int_t njets = fJets->GetEntriesFast();
- for (Int_t ij = 0; ij < njets; ij++) {
-
- AliEmcalJet* jet = static_cast<AliEmcalJet*>(fJets->At(ij));
- if (!jet) {
- AliError(Form("Could not receive jet %d", ij));
- continue;
- }
-
- if (!AcceptJet(jet))
- continue;
+ if (fJetsCont) {
+ AliEmcalJet *jet = fJetsCont->GetNextAcceptJet(0);
+ while(jet) {
fHistJetsPtArea[fCentBin]->Fill(jet->Pt(), jet->Area());
fHistJetsPhiEta[fCentBin]->Fill(jet->Eta(), jet->Phi());
- Float_t ptLeading = GetLeadingHadronPt(jet);
+ Float_t ptLeading = fJetsCont->GetLeadingHadronPt(jet);
fHistJetsPtLeadHad[fCentBin]->Fill(jet->Pt(), ptLeading);
- if (fRho) {
- Float_t corrPt = jet->Pt() - fRhoVal * jet->Area();
+ 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());
}
+ CheckClusTrackMatching();
+
return kTRUE;
}
+//________________________________________________________________________
+void AliAnalysisTaskEmcalJetSample::CheckClusTrackMatching()
+{
+
+ if(!fTracksCont || !fCaloClustersCont)
+ return;
+
+ Double_t deta = 999;
+ Double_t dphi = 999;
+
+ //Get closest cluster to track
+ AliVTrack *track = static_cast<AliVTrack*>(fTracksCont->GetNextAcceptParticle(0));
+ while(track) {
+ //Get matched cluster
+ Int_t emc1 = track->GetEMCALcluster();
+ if(fCaloClustersCont && emc1>=0) {
+ AliVCluster *clusMatch = fCaloClustersCont->GetCluster(emc1);
+ if(clusMatch) {
+ AliPicoTrack::GetEtaPhiDiff(track, clusMatch, dphi, deta);
+ fHistPtDEtaDPhiTrackClus->Fill(track->Pt(),deta,dphi);
+ }
+ }
+ track = static_cast<AliVTrack*>(fTracksCont->GetNextAcceptParticle());
+ }
+
+ //Get closest track to cluster
+ AliVCluster *cluster = fCaloClustersCont->GetNextAcceptCluster(0);
+ while(cluster) {
+ TLorentzVector nPart;
+ cluster->GetMomentum(nPart, fVertex);
+ fHistClustersPt[fCentBin]->Fill(nPart.Pt());
+
+ //Get matched track
+ AliVTrack *mt = NULL;
+ AliAODCaloCluster *acl = dynamic_cast<AliAODCaloCluster*>(cluster);
+ if(acl) {
+ if(acl->GetNTracksMatched()>1)
+ mt = static_cast<AliVTrack*>(acl->GetTrackMatched(0));
+ }
+ else {
+ AliESDCaloCluster *ecl = dynamic_cast<AliESDCaloCluster*>(cluster);
+ Int_t im = ecl->GetTrackMatchedIndex();
+ if(fTracksCont && im>=0) {
+ mt = static_cast<AliVTrack*>(fTracksCont->GetParticle(im));
+ }
+ }
+ if(mt) {
+ AliPicoTrack::GetEtaPhiDiff(mt, cluster, dphi, deta);
+ fHistPtDEtaDPhiClusTrack->Fill(nPart.Pt(),deta,dphi);
+
+ /* //debugging
+ if(mt->IsEMCAL()) {
+ Int_t emc1 = mt->GetEMCALcluster();
+ Printf("current id: %d emc1: %d",fCaloClustersCont->GetCurrentID(),emc1);
+ AliVCluster *clm = fCaloClustersCont->GetCluster(emc1);
+ AliPicoTrack::GetEtaPhiDiff(mt, clm, dphi, deta);
+ Printf("deta: %f dphi: %f",deta,dphi);
+ }
+ */
+ }
+ cluster = fCaloClustersCont->GetNextAcceptCluster();
+ }
+}
+
+//________________________________________________________________________
+void AliAnalysisTaskEmcalJetSample::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 AliAnalysisTaskEmcalJetSample::Run()
{