[u/mrichter/AliRoot.git] / PWGJE / EMCALJetTasks / UserTasks / AliAnalysisTaskSAJF.cxx

// $Id$
//
// Jet analysis task.
//
// Author: S.Aiola

#include <TClonesArray.h>
#include <TH2F.h>
#include <TH3F.h>
#include <THnSparse.h>
#include <TList.h>
#include <TLorentzVector.h>

#include "AliVCluster.h"
#include "AliVParticle.h"
#include "AliEmcalJet.h"
#include "AliRhoParameter.h"
#include "AliLog.h"

#include "AliAnalysisTaskSAJF.h"

ClassImp(AliAnalysisTaskSAJF)

//________________________________________________________________________
AliAnalysisTaskSAJF::AliAnalysisTaskSAJF() : 
  AliAnalysisTaskEmcalJetDev("AliAnalysisTaskSAJF", kTRUE),
  fHistoType(1),
  fHistJetObservables(0)

{
  // Default constructor.

  for (Int_t i = 0; i < 4; i++) {
    fHistTracksJetPt[i] = 0;
    fHistClustersJetPt[i] = 0;
    fHistTracksPtDist[i] = 0;
    fHistClustersPtDist[i] = 0;

    fHistJetPtEtaPhi[i] = 0;
    fHistJetPtArea[i] = 0;
    fHistJetPtEP[i] = 0;
    fHistJetPtNEF[i] = 0;
    fHistJetPtZ[i] = 0;
    fHistJetPtLeadingPartPt[i] = 0;
    fHistJetCorrPtEtaPhi[i] = 0;
    fHistJetCorrPtArea[i] = 0;
    fHistJetCorrPtEP[i] = 0;
    fHistJetCorrPtNEF[i] = 0;
    fHistJetCorrPtZ[i] = 0;
    fHistJetCorrPtLeadingPartPt[i] = 0;
    fHistJetPtCorrPt[i] = 0;
    fHistJetPtMCPt[i] = 0;
    fHistJetMCPtCorrPt[i] = 0;
  }

  SetMakeGeneralHistograms(kTRUE);
}

//________________________________________________________________________
AliAnalysisTaskSAJF::AliAnalysisTaskSAJF(const char *name) : 
  AliAnalysisTaskEmcalJetDev(name, kTRUE),
  fHistoType(1),
  fHistJetObservables(0)
{
  // Standard constructor.

  for (Int_t i = 0; i < 4; i++) {
    fHistTracksJetPt[i] = 0;
    fHistClustersJetPt[i] = 0;
    fHistTracksPtDist[i] = 0;
    fHistClustersPtDist[i] = 0;

    fHistJetPtEtaPhi[i] = 0;
    fHistJetPtArea[i] = 0;
    fHistJetPtEP[i] = 0;
    fHistJetPtNEF[i] = 0;
    fHistJetPtZ[i] = 0;
    fHistJetPtLeadingPartPt[i] = 0;
    fHistJetCorrPtEtaPhi[i] = 0;
    fHistJetCorrPtArea[i] = 0;
    fHistJetCorrPtEP[i] = 0;
    fHistJetCorrPtNEF[i] = 0;
    fHistJetCorrPtZ[i] = 0;
    fHistJetCorrPtLeadingPartPt[i] = 0;
    fHistJetPtCorrPt[i] = 0;
    fHistJetPtMCPt[i] = 0;
    fHistJetMCPtCorrPt[i] = 0;
  }

  SetMakeGeneralHistograms(kTRUE);
}

//________________________________________________________________________
void AliAnalysisTaskSAJF::AllocateTHnSparse()
{
    TString title[20]= {""};
  Int_t nbins[20]  = {0};
  Double_t min[20] = {0.};
  Double_t max[20] = {0.};
  Int_t dim = 0;

  if (fForceBeamType != kpp) {
    title[dim] = "Centrality (%)";
    nbins[dim] = 22;
    min[dim] = -5;
    max[dim] = 105;
    dim++;
    
    title[dim] = "#phi_{jet} - #psi_{RP}";
    nbins[dim] = 100;
    min[dim] = 0;
    max[dim] = TMath::Pi();
    dim++;
  }

  title[dim] = "#eta";
  nbins[dim] = 100;
  min[dim] = -1;
  max[dim] = 1;
  dim++;

  title[dim] = "#phi_{jet} (rad)";
  nbins[dim] = 201;
  min[dim] = 0;
  max[dim] = 2*TMath::Pi()*nbins[dim]/(nbins[dim]-1);
  dim++;

  title[dim] = "p_{T} (GeV/c)";
  nbins[dim] = fNbins;
  min[dim] = fMinBinPt;
  max[dim] = fMaxBinPt;
  dim++;

  if (fIsEmbedded) {
    title[dim] = "p_{T}^{MC} (GeV/c)";
    nbins[dim] = fNbins;
    min[dim] = fMinBinPt;
    max[dim] = fMaxBinPt;
    dim++;
  }

  if (!GetRhoName().IsNull()) {
    title[dim] = "p_{T}^{corr} (GeV/c)";
    nbins[dim] = fNbins*2;
    min[dim] = -fMaxBinPt;
    max[dim] = fMaxBinPt;
    dim++;
  }

  title[dim] = "A_{jet}";
  nbins[dim] = 150;
  min[dim] = 0;
  max[dim] = 1.5;
  dim++;

  title[dim] = "NEF";
  nbins[dim] = 102;
  min[dim] = 0;
  max[dim] = 1.02;
  dim++;

  title[dim] = "Z";
  nbins[dim] = 102;
  min[dim] = 0;
  max[dim] = 1.02;
  dim++;

  title[dim] = "No. of constituents";
  nbins[dim] = 250;
  min[dim] = -0.5;
  max[dim] = 249.5;
  dim++;

  title[dim] = "p_{T,particle}^{leading} (GeV/c)";
  nbins[dim] = 120;
  min[dim] = 0;
  max[dim] = 120;
  dim++;

  fHistJetObservables = new THnSparseD("fHistJetObservables","fHistJetObservables",dim,nbins,min,max);
  fOutput->Add(fHistJetObservables);
  for (Int_t i = 0; i < dim; i++)
    fHistJetObservables->GetAxis(i)->SetTitle(title[i]);
}

//________________________________________________________________________
void AliAnalysisTaskSAJF::AllocateTHX()
{
  for (Int_t i = 0; i < 4; i++) {
    TString histname;

    histname = "fHistJetPtEtaPhi_";
    histname += i;
    fHistJetPtEtaPhi[i] = new TH3F(histname.Data(), histname.Data(), fNbins, fMinBinPt, fMaxBinPt, 20, -1, 1, 41, 0, 2*TMath::Pi()*41/40);
    fHistJetPtEtaPhi[i]->GetXaxis()->SetTitle("p_{T} (GeV/c)");
    fHistJetPtEtaPhi[i]->GetYaxis()->SetTitle("#eta");
    fHistJetPtEtaPhi[i]->GetZaxis()->SetTitle("#phi_{jet} (rad)");
    fOutput->Add(fHistJetPtEtaPhi[i]);
      
    histname = "fHistJetPtArea_";
    histname += i;
    fHistJetPtArea[i] = new TH2F(histname.Data(), histname.Data(), fNbins, fMinBinPt, fMaxBinPt, 150, 0, 1.5);
    fHistJetPtArea[i]->GetXaxis()->SetTitle("p_{T} (GeV/c)");
    fHistJetPtArea[i]->GetYaxis()->SetTitle("A_{jet}");
    fHistJetPtArea[i]->GetZaxis()->SetTitle("counts");
    fOutput->Add(fHistJetPtArea[i]);

    histname = "fHistJetPtEP_";
    histname += i;
    fHistJetPtEP[i] = new TH2F(histname.Data(), histname.Data(), fNbins, fMinBinPt, fMaxBinPt, 100, 0, TMath::Pi());
    fHistJetPtEP[i]->GetXaxis()->SetTitle("p_{T} (GeV/c)");
    fHistJetPtEP[i]->GetYaxis()->SetTitle("#phi_{jet} - #psi_{RP}");
    fHistJetPtEP[i]->GetZaxis()->SetTitle("counts");
    fOutput->Add(fHistJetPtEP[i]);

    histname = "fHistJetPtNEF_";
    histname += i;
    fHistJetPtNEF[i] = new TH2F(histname.Data(), histname.Data(), fNbins, fMinBinPt, fMaxBinPt, 102, 0, 1.02);
    fHistJetPtNEF[i]->GetXaxis()->SetTitle("p_{T} (GeV/c)");
    fHistJetPtNEF[i]->GetYaxis()->SetTitle("NEF");
    fHistJetPtNEF[i]->GetZaxis()->SetTitle("counts");
    fOutput->Add(fHistJetPtNEF[i]);

    histname = "fHistJetPtZ_";
    histname += i;
    fHistJetPtZ[i] = new TH2F(histname.Data(), histname.Data(), fNbins, fMinBinPt, fMaxBinPt, 102, 0, 1.02);
    fHistJetPtZ[i]->GetXaxis()->SetTitle("p_{T} (GeV/c)");
    fHistJetPtZ[i]->GetYaxis()->SetTitle("z");
    fHistJetPtZ[i]->GetZaxis()->SetTitle("counts");
    fOutput->Add(fHistJetPtZ[i]);

    histname = "fHistJetPtLeadingPartPt_";
    histname += i;
    fHistJetPtLeadingPartPt[i] = new TH2F(histname.Data(), histname.Data(), fNbins, fMinBinPt, fMaxBinPt, 120, 0, 120);
    fHistJetPtLeadingPartPt[i]->GetXaxis()->SetTitle("p_{T} (GeV/c)");
    fHistJetPtLeadingPartPt[i]->GetYaxis()->SetTitle("p_{T,particle}^{leading} (GeV/c)");
    fHistJetPtLeadingPartPt[i]->GetZaxis()->SetTitle("counts");
    fOutput->Add(fHistJetPtLeadingPartPt[i]);

    if (!GetRhoName().IsNull()) {
      histname = "fHistJetCorrPtEtaPhi_";
      histname += i;
      fHistJetCorrPtEtaPhi[i] = new TH3F(histname.Data(), histname.Data(), fNbins*2, -fMaxBinPt, fMaxBinPt, 20, -1, 1, 41, 0, 2*TMath::Pi()*201/200);
      fHistJetCorrPtEtaPhi[i]->GetXaxis()->SetTitle("p_{T} (GeV/c)");
      fHistJetCorrPtEtaPhi[i]->GetYaxis()->SetTitle("#eta");
      fHistJetCorrPtEtaPhi[i]->GetZaxis()->SetTitle("#phi_{jet} (rad)");
      fOutput->Add(fHistJetCorrPtEtaPhi[i]);
      
      histname = "fHistJetCorrPtArea_";
      histname += i;
      fHistJetCorrPtArea[i] = new TH2F(histname.Data(), histname.Data(), fNbins*2, -fMaxBinPt, fMaxBinPt, 150, 0, 1.5);
      fHistJetCorrPtArea[i]->GetXaxis()->SetTitle("p_{T}^{corr} (GeV/c)");
      fHistJetCorrPtArea[i]->GetYaxis()->SetTitle("A_{jet}");
      fHistJetCorrPtArea[i]->GetZaxis()->SetTitle("counts");
      fOutput->Add(fHistJetCorrPtArea[i]);

      histname = "fHistJetCorrPtEP_";
      histname += i;
      fHistJetCorrPtEP[i] = new TH2F(histname.Data(), histname.Data(), fNbins*2, -fMaxBinPt, fMaxBinPt, 100, 0, TMath::Pi());
      fHistJetCorrPtEP[i]->GetXaxis()->SetTitle("p_{T}^{corr} (GeV/c)");
      fHistJetCorrPtEP[i]->GetYaxis()->SetTitle("#phi_{jet} - #psi_{RP}");
      fHistJetCorrPtEP[i]->GetZaxis()->SetTitle("counts");
      fOutput->Add(fHistJetCorrPtEP[i]);

      histname = "fHistJetCorrPtNEF_";
      histname += i;
      fHistJetCorrPtNEF[i] = new TH2F(histname.Data(), histname.Data(), fNbins*2, -fMaxBinPt, fMaxBinPt, 102, 0, 1.02);
      fHistJetCorrPtNEF[i]->GetXaxis()->SetTitle("p_{T}^{corr} (GeV/c)");
      fHistJetCorrPtNEF[i]->GetYaxis()->SetTitle("NEF");
      fHistJetCorrPtNEF[i]->GetZaxis()->SetTitle("counts");
      fOutput->Add(fHistJetCorrPtNEF[i]);

      histname = "fHistJetCorrPtZ_";
      histname += i;
      fHistJetCorrPtZ[i] = new TH2F(histname.Data(), histname.Data(), fNbins*2, -fMaxBinPt, fMaxBinPt, 102, 0, 1.02);
      fHistJetCorrPtZ[i]->GetXaxis()->SetTitle("p_{T}^{corr} (GeV/c)");
      fHistJetCorrPtZ[i]->GetYaxis()->SetTitle("z");
      fHistJetCorrPtZ[i]->GetZaxis()->SetTitle("counts");
      fOutput->Add(fHistJetCorrPtZ[i]);

      histname = "fHistJetCorrPtLeadingPartPt_";
      histname += i;
      fHistJetCorrPtLeadingPartPt[i] = new TH2F(histname.Data(), histname.Data(), fNbins*2, -fMaxBinPt, fMaxBinPt, 120, 0, 120);
      fHistJetCorrPtLeadingPartPt[i]->GetXaxis()->SetTitle("p_{T}^{corr} (GeV/c)");
      fHistJetCorrPtLeadingPartPt[i]->GetYaxis()->SetTitle("p_{T,particle}^{leading} (GeV/c)");
      fHistJetCorrPtLeadingPartPt[i]->GetZaxis()->SetTitle("counts");
      fOutput->Add(fHistJetCorrPtLeadingPartPt[i]);

      histname = "fHistJetPtCorrPt_";
      histname += i;
      fHistJetPtCorrPt[i] = new TH2F(histname.Data(), histname.Data(), fNbins, fMinBinPt, fMaxBinPt, fNbins*2, -fMaxBinPt, fMaxBinPt);
      fHistJetPtCorrPt[i]->GetXaxis()->SetTitle("p_{T} (GeV/c)");
      fHistJetPtCorrPt[i]->GetYaxis()->SetTitle("p_{T}^{corr} (GeV/c)");
      fHistJetPtCorrPt[i]->GetZaxis()->SetTitle("counts");
      fOutput->Add(fHistJetPtCorrPt[i]);

      if (fIsEmbedded) {
	histname = "fHistJetMCPtCorrPt_";
	histname += i;
	fHistJetMCPtCorrPt[i] = new TH2F(histname.Data(), histname.Data(), fNbins, fMinBinPt, fMaxBinPt, fNbins*2, -fMaxBinPt, fMaxBinPt);
	fHistJetMCPtCorrPt[i]->GetXaxis()->SetTitle("p_{T}^{MC} (GeV/c)");
	fHistJetMCPtCorrPt[i]->GetYaxis()->SetTitle("p_{T}^{corr} (GeV/c)");
	fHistJetMCPtCorrPt[i]->GetZaxis()->SetTitle("counts");
	fOutput->Add(fHistJetMCPtCorrPt[i]);
      }
    }

    if (fIsEmbedded) {
      histname = "fHistJetPtMCPt_";
      histname += i;
      fHistJetPtMCPt[i] = new TH2F(histname.Data(), histname.Data(), fNbins, fMinBinPt, fMaxBinPt, fNbins, fMinBinPt, fMaxBinPt);
      fHistJetPtMCPt[i]->GetXaxis()->SetTitle("p_{T} (GeV/c)");
      fHistJetPtMCPt[i]->GetYaxis()->SetTitle("p_{T}^{MC} (GeV/c)");
      fHistJetPtMCPt[i]->GetZaxis()->SetTitle("counts");
      fOutput->Add(fHistJetPtMCPt[i]);
    }
  }
}

//________________________________________________________________________
void AliAnalysisTaskSAJF::UserCreateOutputObjects()
{
  // Create user output.

  AliAnalysisTaskEmcalJetDev::UserCreateOutputObjects();

  if (fHistoType == 0)
    AllocateTHX();
  else
    AllocateTHnSparse();

  for (Int_t i = 0; i < 4; i++) {
    TString histname;

    if (fParticleCollArray.GetEntriesFast()>0) {
      histname = "fHistTracksJetPt_";
      histname += i;
      fHistTracksJetPt[i] = new TH2F(histname.Data(), histname.Data(), fNbins / 2, fMinBinPt, fMaxBinPt / 2, fNbins, fMinBinPt, fMaxBinPt);
      fHistTracksJetPt[i]->GetXaxis()->SetTitle("p_{T,track} (GeV/c)");
      fHistTracksJetPt[i]->GetYaxis()->SetTitle("p_{T,jet} (GeV/c)");
      fHistTracksJetPt[i]->GetZaxis()->SetTitle("counts");
      fOutput->Add(fHistTracksJetPt[i]);
      
      histname = "fHistTracksPtDist_";
      histname += i;
      fHistTracksPtDist[i] = new TH2F(histname.Data(), histname.Data(), fNbins / 2, fMinBinPt, fMaxBinPt / 2, 100, 0, 5);
      fHistTracksPtDist[i]->GetXaxis()->SetTitle("p_{T,track} (GeV/c)");
      fHistTracksPtDist[i]->GetYaxis()->SetTitle("d");
      fHistTracksPtDist[i]->GetZaxis()->SetTitle("counts");
      fOutput->Add(fHistTracksPtDist[i]);
    }

    if (fClusterCollArray.GetEntriesFast()>0) {
      histname = "fHistClustersJetPt_";
      histname += i;
      fHistClustersJetPt[i] = new TH2F(histname.Data(), histname.Data(), fNbins / 2, fMinBinPt, fMaxBinPt / 2, fNbins, fMinBinPt, fMaxBinPt);
      fHistClustersJetPt[i]->GetXaxis()->SetTitle("p_{T,clus} (GeV/c)");
      fHistClustersJetPt[i]->GetYaxis()->SetTitle("p_{T,jet} (GeV/c)");
      fHistClustersJetPt[i]->GetZaxis()->SetTitle("counts");
      fOutput->Add(fHistClustersJetPt[i]);

      histname = "fHistClustersPtDist_";
      histname += i;
      fHistClustersPtDist[i] = new TH2F(histname.Data(), histname.Data(), fNbins / 2, fMinBinPt, fMaxBinPt / 2, 100, 0, 5);
      fHistClustersPtDist[i]->GetXaxis()->SetTitle("p_{T,clus} (GeV/c)");
      fHistClustersPtDist[i]->GetYaxis()->SetTitle("d");
      fHistClustersPtDist[i]->GetZaxis()->SetTitle("counts");
      fOutput->Add(fHistClustersPtDist[i]);
    }
  }

  PostData(1, fOutput); // Post data for ALL output slots >0 here, to get at least an empty histogram

}

//________________________________________________________________________
Bool_t AliAnalysisTaskSAJF::FillHistograms()
{
  // Fill histograms.

  if (!fJets) {
    AliError(Form("%s - Jet array not provided, returning...", GetName()));
    return kFALSE;
  }

  for (Int_t ij = 0; ij < fJets->GetEntriesFast(); ij++) {

    AliEmcalJet* jet = static_cast<AliEmcalJet*>(fJets->At(ij));

    if (!jet) {
      AliError(Form("Could not receive jet %d", ij));
      continue;
    }

    if (!AcceptJet(jet))
      continue;

    Float_t ptLeading = GetLeadingHadronPt(jet);
    Float_t corrPt = jet->Pt() - fRhoVal * jet->Area();

    // Fill THnSparse
    Double_t ep = jet->Phi() - fEPV0;
    while (ep < 0) ep += TMath::Pi();
    while (ep >= TMath::Pi()) ep -= TMath::Pi();

    FillJetHisto(fCent, ep, jet->Eta(), jet->Phi(), jet->Pt(), jet->MCPt(), corrPt, jet->Area(), 
		 jet->NEF(), ptLeading/jet->Pt(), jet->GetNumberOfConstituents(), ptLeading);

    if (fTracks) {
      for (Int_t it = 0; it < jet->GetNumberOfTracks(); it++) {
	AliVParticle *track = jet->TrackAt(it, fTracks);
	if (track) {
	  fHistTracksJetPt[fCentBin]->Fill(track->Pt(), jet->Pt());
	  Double_t dist = TMath::Sqrt((track->Eta() - jet->Eta()) * (track->Eta() - jet->Eta()) + (track->Phi() - jet->Phi()) * (track->Phi() - jet->Phi()));
	  fHistTracksPtDist[fCentBin]->Fill(track->Pt(), dist);
	}
      }
    }

    if (fCaloClusters) {
      for (Int_t ic = 0; ic < jet->GetNumberOfClusters(); ic++) {
	AliVCluster *cluster = jet->ClusterAt(ic, fCaloClusters);
	
	if (cluster) {
	  TLorentzVector nPart;
	  cluster->GetMomentum(nPart, fVertex);

	  fHistClustersJetPt[fCentBin]->Fill(nPart.Pt(), jet->Pt());
	  Double_t dist = TMath::Sqrt((nPart.Eta() - jet->Eta()) * (nPart.Eta() - jet->Eta()) + (nPart.Phi() - jet->Phi()) * (nPart.Phi() - jet->Phi()));
	  fHistClustersPtDist[fCentBin]->Fill(nPart.Pt(), dist);
	}
      }
    }
  } //jet loop 

  return kTRUE;
}

//________________________________________________________________________
void AliAnalysisTaskSAJF::FillJetHisto(Double_t cent, Double_t ep, Double_t eta, Double_t phi, Double_t pt, Double_t MCpt, Double_t corrpt, Double_t area, 
				       Double_t NEF, Double_t z, Int_t n, Double_t leadingpt)
{
  if (fHistoType == 0) {
    fHistJetPtEtaPhi[fCentBin]->Fill(pt,eta,phi);
    fHistJetPtArea[fCentBin]->Fill(pt,area);
    fHistJetPtEP[fCentBin]->Fill(pt,ep);
    fHistJetPtNEF[fCentBin]->Fill(pt,NEF);
    fHistJetPtZ[fCentBin]->Fill(pt,z);
    fHistJetPtLeadingPartPt[fCentBin]->Fill(pt,leadingpt);
    if (fHistJetCorrPtEtaPhi[fCentBin]) {
      fHistJetCorrPtEtaPhi[fCentBin]->Fill(corrpt,eta,phi);
      fHistJetCorrPtArea[fCentBin]->Fill(corrpt,area);
      fHistJetCorrPtEP[fCentBin]->Fill(corrpt,ep);
      fHistJetCorrPtNEF[fCentBin]->Fill(corrpt,NEF);
      fHistJetCorrPtZ[fCentBin]->Fill(corrpt,z);
      fHistJetCorrPtLeadingPartPt[fCentBin]->Fill(corrpt,leadingpt);
      fHistJetPtCorrPt[fCentBin]->Fill(pt,corrpt);
      if (fIsEmbedded)
	fHistJetMCPtCorrPt[fCentBin]->Fill(MCpt,corrpt);
    }
    if (fIsEmbedded)
      fHistJetPtMCPt[fCentBin]->Fill(pt,MCpt);
  }
  else {
 
    Double_t contents[20]={0};

    for (Int_t i = 0; i < fHistJetObservables->GetNdimensions(); i++) {
      TString title(fHistJetObservables->GetAxis(i)->GetTitle());
      if (title=="Centrality (%)")
	contents[i] = cent;
      else if (title=="#phi_{jet} - #psi_{RP}")
	contents[i] = ep;
      else if (title=="#eta")
	contents[i] = eta;
      else if (title=="#phi_{jet} (rad)")
	contents[i] = phi;
      else if (title=="p_{T} (GeV/c)")
	contents[i] = pt;
      else if (title=="p_{T}^{MC} (GeV/c)")
	contents[i] = MCpt;
      else if (title=="p_{T}^{corr} (GeV/c)")
	contents[i] = corrpt;
      else if (title=="A_{jet}")
	contents[i] = area;
      else if (title=="NEF")
	contents[i] = NEF;
      else if (title=="Z")
	contents[i] = z;
      else if (title=="No. of constituents")
	contents[i] = n;
      else if (title=="p_{T,particle}^{leading} (GeV/c)")
	contents[i] = leadingpt;
      else 
	AliWarning(Form("Unable to fill dimension %s!",title.Data()));
    }

    fHistJetObservables->Fill(contents);
  }
}
Commit	Line	Data
00514d01	1	// $Id$
91f4b7c5	2	//
297edd60	3	// Jet analysis task.
91f4b7c5	4	//
297edd60	5	// Author: S.Aiola
25283b37	6
25283b37	7	#include <TClonesArray.h>
25283b37	8	#include <TH2F.h>
85d48773	9	#include <TH3F.h>
43032ce2	10	#include <THnSparse.h>
25283b37	11	#include <TList.h>
25283b37	12	#include <TLorentzVector.h>
25283b37	13
f0a0fd33	14	#include "AliVCluster.h"
df43b607	15	#include "AliVParticle.h"
25283b37	16	#include "AliEmcalJet.h"
e44e8726	17	#include "AliRhoParameter.h"
55264f20	18	#include "AliLog.h"
25283b37	19
00514d01	20	#include "AliAnalysisTaskSAJF.h"
25283b37	21
00514d01	22	ClassImp(AliAnalysisTaskSAJF)
25283b37	23
25283b37	24	//________________________________________________________________________
00514d01	25	AliAnalysisTaskSAJF::AliAnalysisTaskSAJF() :
7cd832c7	26	AliAnalysisTaskEmcalJetDev("AliAnalysisTaskSAJF", kTRUE),
85d48773	27	fHistoType(1),
43032ce2	28	fHistJetObservables(0)
226f511d	29
25283b37	30	{
25283b37	31	// Default constructor.
f0a0fd33	32
f0a0fd33	33	for (Int_t i = 0; i < 4; i++) {
a487deae	34	fHistTracksJetPt[i] = 0;
a487deae	35	fHistClustersJetPt[i] = 0;
7030f36f	36	fHistTracksPtDist[i] = 0;
7030f36f	37	fHistClustersPtDist[i] = 0;
85d48773	38
	39	fHistJetPtEtaPhi[i] = 0;
	40	fHistJetPtArea[i] = 0;
	41	fHistJetPtEP[i] = 0;
	42	fHistJetPtNEF[i] = 0;
	43	fHistJetPtZ[i] = 0;
	44	fHistJetPtLeadingPartPt[i] = 0;
	45	fHistJetCorrPtEtaPhi[i] = 0;
	46	fHistJetCorrPtArea[i] = 0;
	47	fHistJetCorrPtEP[i] = 0;
	48	fHistJetCorrPtNEF[i] = 0;
	49	fHistJetCorrPtZ[i] = 0;
	50	fHistJetCorrPtLeadingPartPt[i] = 0;
	51	fHistJetPtCorrPt[i] = 0;
	52	fHistJetPtMCPt[i] = 0;
	53	fHistJetMCPtCorrPt[i] = 0;
f0a0fd33	54	}
a487deae	55
a487deae	56	SetMakeGeneralHistograms(kTRUE);
25283b37	57	}
	58
	59	//________________________________________________________________________
00514d01	60	AliAnalysisTaskSAJF::AliAnalysisTaskSAJF(const char *name) :
7cd832c7	61	AliAnalysisTaskEmcalJetDev(name, kTRUE),
85d48773	62	fHistoType(1),
43032ce2	63	fHistJetObservables(0)
25283b37	64	{
	65	// Standard constructor.
	66
f0a0fd33	67	for (Int_t i = 0; i < 4; i++) {
a487deae	68	fHistTracksJetPt[i] = 0;
a487deae	69	fHistClustersJetPt[i] = 0;
7030f36f	70	fHistTracksPtDist[i] = 0;
7030f36f	71	fHistClustersPtDist[i] = 0;
85d48773	72
	73	fHistJetPtEtaPhi[i] = 0;
	74	fHistJetPtArea[i] = 0;
	75	fHistJetPtEP[i] = 0;
	76	fHistJetPtNEF[i] = 0;
	77	fHistJetPtZ[i] = 0;
	78	fHistJetPtLeadingPartPt[i] = 0;
	79	fHistJetCorrPtEtaPhi[i] = 0;
	80	fHistJetCorrPtArea[i] = 0;
	81	fHistJetCorrPtEP[i] = 0;
	82	fHistJetCorrPtNEF[i] = 0;
	83	fHistJetCorrPtZ[i] = 0;
	84	fHistJetCorrPtLeadingPartPt[i] = 0;
	85	fHistJetPtCorrPt[i] = 0;
	86	fHistJetPtMCPt[i] = 0;
	87	fHistJetMCPtCorrPt[i] = 0;
f0a0fd33	88	}
a487deae	89
a487deae	90	SetMakeGeneralHistograms(kTRUE);
25283b37	91	}
25283b37	92
25283b37	93	//________________________________________________________________________
85d48773	94	void AliAnalysisTaskSAJF::AllocateTHnSparse()
25283b37	95	{
85d48773	96	TString title[20]= {""};
43032ce2	97	Int_t nbins[20] = {0};
	98	Double_t min[20] = {0.};
	99	Double_t max[20] = {0.};
	100	Int_t dim = 0;
	101
	102	if (fForceBeamType != kpp) {
	103	title[dim] = "Centrality (%)";
85d48773	104	nbins[dim] = 22;
	105	min[dim] = -5;
	106	max[dim] = 105;
43032ce2	107	dim++;
43032ce2	108
3fe08cdb	109	title[dim] = "#phi_{jet} - #psi_{RP}";
85d48773	110	nbins[dim] = 100;
	111	min[dim] = 0;
	112	max[dim] = TMath::Pi();
43032ce2	113	dim++;
43032ce2	114	}
624bef5b	115
3fe08cdb	116	title[dim] = "#eta";
85d48773	117	nbins[dim] = 100;
3fe08cdb	118	min[dim] = -1;
43032ce2	119	max[dim] = 1;
	120	dim++;
	121
43032ce2	122	title[dim] = "#phi_{jet} (rad)";
85d48773	123	nbins[dim] = 201;
43032ce2	124	min[dim] = 0;
3fe08cdb	125	max[dim] = 2TMath::Pi()nbins[dim]/(nbins[dim]-1);
43032ce2	126	dim++;
	127
	128	title[dim] = "p_{T} (GeV/c)";
	129	nbins[dim] = fNbins;
85d48773	130	min[dim] = fMinBinPt;
85d48773	131	max[dim] = fMaxBinPt;
43032ce2	132	dim++;
43032ce2	133
43032ce2	134	if (fIsEmbedded) {
	135	title[dim] = "p_{T}^{MC} (GeV/c)";
	136	nbins[dim] = fNbins;
85d48773	137	min[dim] = fMinBinPt;
85d48773	138	max[dim] = fMaxBinPt;
43032ce2	139	dim++;
43032ce2	140	}
59f16b27	141
7cd832c7	142	if (!GetRhoName().IsNull()) {
43032ce2	143	title[dim] = "p_{T}^{corr} (GeV/c)";
43032ce2	144	nbins[dim] = fNbins*2;
85d48773	145	min[dim] = -fMaxBinPt;
85d48773	146	max[dim] = fMaxBinPt;
43032ce2	147	dim++;
43032ce2	148	}
05077f28	149
3fe08cdb	150	title[dim] = "A_{jet}";
85d48773	151	nbins[dim] = 150;
3fe08cdb	152	min[dim] = 0;
85d48773	153	max[dim] = 1.5;
43032ce2	154	dim++;
	155
	156	title[dim] = "NEF";
85d48773	157	nbins[dim] = 102;
43032ce2	158	min[dim] = 0;
85d48773	159	max[dim] = 1.02;
43032ce2	160	dim++;
	161
	162	title[dim] = "Z";
85d48773	163	nbins[dim] = 102;
43032ce2	164	min[dim] = 0;
85d48773	165	max[dim] = 1.02;
3fe08cdb	166	dim++;
	167
	168	title[dim] = "No. of constituents";
	169	nbins[dim] = 250;
	170	min[dim] = -0.5;
	171	max[dim] = 249.5;
43032ce2	172	dim++;
	173
	174	title[dim] = "p_{T,particle}^{leading} (GeV/c)";
85d48773	175	nbins[dim] = 120;
43032ce2	176	min[dim] = 0;
85d48773	177	max[dim] = 120;
43032ce2	178	dim++;
	179
	180	fHistJetObservables = new THnSparseD("fHistJetObservables","fHistJetObservables",dim,nbins,min,max);
3fe08cdb	181	fOutput->Add(fHistJetObservables);
43032ce2	182	for (Int_t i = 0; i < dim; i++)
43032ce2	183	fHistJetObservables->GetAxis(i)->SetTitle(title[i]);
85d48773	184	}
	185
	186	//________________________________________________________________________
	187	void AliAnalysisTaskSAJF::AllocateTHX()
	188	{
	189	for (Int_t i = 0; i < 4; i++) {
	190	TString histname;
	191
	192	histname = "fHistJetPtEtaPhi_";
	193	histname += i;
	194	fHistJetPtEtaPhi[i] = new TH3F(histname.Data(), histname.Data(), fNbins, fMinBinPt, fMaxBinPt, 20, -1, 1, 41, 0, 2TMath::Pi()41/40);
	195	fHistJetPtEtaPhi[i]->GetXaxis()->SetTitle("p_{T} (GeV/c)");
	196	fHistJetPtEtaPhi[i]->GetYaxis()->SetTitle("#eta");
	197	fHistJetPtEtaPhi[i]->GetZaxis()->SetTitle("#phi_{jet} (rad)");
	198	fOutput->Add(fHistJetPtEtaPhi[i]);
	199
	200	histname = "fHistJetPtArea_";
	201	histname += i;
	202	fHistJetPtArea[i] = new TH2F(histname.Data(), histname.Data(), fNbins, fMinBinPt, fMaxBinPt, 150, 0, 1.5);
	203	fHistJetPtArea[i]->GetXaxis()->SetTitle("p_{T} (GeV/c)");
	204	fHistJetPtArea[i]->GetYaxis()->SetTitle("A_{jet}");
	205	fHistJetPtArea[i]->GetZaxis()->SetTitle("counts");
	206	fOutput->Add(fHistJetPtArea[i]);
	207
	208	histname = "fHistJetPtEP_";
	209	histname += i;
	210	fHistJetPtEP[i] = new TH2F(histname.Data(), histname.Data(), fNbins, fMinBinPt, fMaxBinPt, 100, 0, TMath::Pi());
	211	fHistJetPtEP[i]->GetXaxis()->SetTitle("p_{T} (GeV/c)");
	212	fHistJetPtEP[i]->GetYaxis()->SetTitle("#phi_{jet} - #psi_{RP}");
	213	fHistJetPtEP[i]->GetZaxis()->SetTitle("counts");
	214	fOutput->Add(fHistJetPtEP[i]);
	215
	216	histname = "fHistJetPtNEF_";
	217	histname += i;
	218	fHistJetPtNEF[i] = new TH2F(histname.Data(), histname.Data(), fNbins, fMinBinPt, fMaxBinPt, 102, 0, 1.02);
	219	fHistJetPtNEF[i]->GetXaxis()->SetTitle("p_{T} (GeV/c)");
	220	fHistJetPtNEF[i]->GetYaxis()->SetTitle("NEF");
	221	fHistJetPtNEF[i]->GetZaxis()->SetTitle("counts");
	222	fOutput->Add(fHistJetPtNEF[i]);
	223
	224	histname = "fHistJetPtZ_";
	225	histname += i;
	226	fHistJetPtZ[i] = new TH2F(histname.Data(), histname.Data(), fNbins, fMinBinPt, fMaxBinPt, 102, 0, 1.02);
	227	fHistJetPtZ[i]->GetXaxis()->SetTitle("p_{T} (GeV/c)");
	228	fHistJetPtZ[i]->GetYaxis()->SetTitle("z");
	229	fHistJetPtZ[i]->GetZaxis()->SetTitle("counts");
	230	fOutput->Add(fHistJetPtZ[i]);
	231
	232	histname = "fHistJetPtLeadingPartPt_";
	233	histname += i;
	234	fHistJetPtLeadingPartPt[i] = new TH2F(histname.Data(), histname.Data(), fNbins, fMinBinPt, fMaxBinPt, 120, 0, 120);
	235	fHistJetPtLeadingPartPt[i]->GetXaxis()->SetTitle("p_{T} (GeV/c)");
	236	fHistJetPtLeadingPartPt[i]->GetYaxis()->SetTitle("p_{T,particle}^{leading} (GeV/c)");
	237	fHistJetPtLeadingPartPt[i]->GetZaxis()->SetTitle("counts");
	238	fOutput->Add(fHistJetPtLeadingPartPt[i]);
	239
7cd832c7	240	if (!GetRhoName().IsNull()) {
85d48773	241	histname = "fHistJetCorrPtEtaPhi_";
	242	histname += i;
	243	fHistJetCorrPtEtaPhi[i] = new TH3F(histname.Data(), histname.Data(), fNbins2, -fMaxBinPt, fMaxBinPt, 20, -1, 1, 41, 0, 2TMath::Pi()*201/200);
	244	fHistJetCorrPtEtaPhi[i]->GetXaxis()->SetTitle("p_{T} (GeV/c)");
	245	fHistJetCorrPtEtaPhi[i]->GetYaxis()->SetTitle("#eta");
	246	fHistJetCorrPtEtaPhi[i]->GetZaxis()->SetTitle("#phi_{jet} (rad)");
	247	fOutput->Add(fHistJetCorrPtEtaPhi[i]);
	248
	249	histname = "fHistJetCorrPtArea_";
	250	histname += i;
	251	fHistJetCorrPtArea[i] = new TH2F(histname.Data(), histname.Data(), fNbins*2, -fMaxBinPt, fMaxBinPt, 150, 0, 1.5);
	252	fHistJetCorrPtArea[i]->GetXaxis()->SetTitle("p_{T}^{corr} (GeV/c)");
	253	fHistJetCorrPtArea[i]->GetYaxis()->SetTitle("A_{jet}");
	254	fHistJetCorrPtArea[i]->GetZaxis()->SetTitle("counts");
	255	fOutput->Add(fHistJetCorrPtArea[i]);
	256
	257	histname = "fHistJetCorrPtEP_";
	258	histname += i;
	259	fHistJetCorrPtEP[i] = new TH2F(histname.Data(), histname.Data(), fNbins*2, -fMaxBinPt, fMaxBinPt, 100, 0, TMath::Pi());
	260	fHistJetCorrPtEP[i]->GetXaxis()->SetTitle("p_{T}^{corr} (GeV/c)");
	261	fHistJetCorrPtEP[i]->GetYaxis()->SetTitle("#phi_{jet} - #psi_{RP}");
	262	fHistJetCorrPtEP[i]->GetZaxis()->SetTitle("counts");
	263	fOutput->Add(fHistJetCorrPtEP[i]);
	264
	265	histname = "fHistJetCorrPtNEF_";
	266	histname += i;
	267	fHistJetCorrPtNEF[i] = new TH2F(histname.Data(), histname.Data(), fNbins*2, -fMaxBinPt, fMaxBinPt, 102, 0, 1.02);
	268	fHistJetCorrPtNEF[i]->GetXaxis()->SetTitle("p_{T}^{corr} (GeV/c)");
	269	fHistJetCorrPtNEF[i]->GetYaxis()->SetTitle("NEF");
	270	fHistJetCorrPtNEF[i]->GetZaxis()->SetTitle("counts");
	271	fOutput->Add(fHistJetCorrPtNEF[i]);
	272
	273	histname = "fHistJetCorrPtZ_";
	274	histname += i;
	275	fHistJetCorrPtZ[i] = new TH2F(histname.Data(), histname.Data(), fNbins*2, -fMaxBinPt, fMaxBinPt, 102, 0, 1.02);
	276	fHistJetCorrPtZ[i]->GetXaxis()->SetTitle("p_{T}^{corr} (GeV/c)");
	277	fHistJetCorrPtZ[i]->GetYaxis()->SetTitle("z");
	278	fHistJetCorrPtZ[i]->GetZaxis()->SetTitle("counts");
	279	fOutput->Add(fHistJetCorrPtZ[i]);
	280
	281	histname = "fHistJetCorrPtLeadingPartPt_";
	282	histname += i;
	283	fHistJetCorrPtLeadingPartPt[i] = new TH2F(histname.Data(), histname.Data(), fNbins*2, -fMaxBinPt, fMaxBinPt, 120, 0, 120);
	284	fHistJetCorrPtLeadingPartPt[i]->GetXaxis()->SetTitle("p_{T}^{corr} (GeV/c)");
	285	fHistJetCorrPtLeadingPartPt[i]->GetYaxis()->SetTitle("p_{T,particle}^{leading} (GeV/c)");
	286	fHistJetCorrPtLeadingPartPt[i]->GetZaxis()->SetTitle("counts");
	287	fOutput->Add(fHistJetCorrPtLeadingPartPt[i]);
	288
	289	histname = "fHistJetPtCorrPt_";
	290	histname += i;
	291	fHistJetPtCorrPt[i] = new TH2F(histname.Data(), histname.Data(), fNbins, fMinBinPt, fMaxBinPt, fNbins*2, -fMaxBinPt, fMaxBinPt);
	292	fHistJetPtCorrPt[i]->GetXaxis()->SetTitle("p_{T} (GeV/c)");
	293	fHistJetPtCorrPt[i]->GetYaxis()->SetTitle("p_{T}^{corr} (GeV/c)");
	294	fHistJetPtCorrPt[i]->GetZaxis()->SetTitle("counts");
	295	fOutput->Add(fHistJetPtCorrPt[i]);
	296
	297	if (fIsEmbedded) {
	298	histname = "fHistJetMCPtCorrPt_";
	299	histname += i;
	300	fHistJetMCPtCorrPt[i] = new TH2F(histname.Data(), histname.Data(), fNbins, fMinBinPt, fMaxBinPt, fNbins*2, -fMaxBinPt, fMaxBinPt);
	301	fHistJetMCPtCorrPt[i]->GetXaxis()->SetTitle("p_{T}^{MC} (GeV/c)");
	302	fHistJetMCPtCorrPt[i]->GetYaxis()->SetTitle("p_{T}^{corr} (GeV/c)");
	303	fHistJetMCPtCorrPt[i]->GetZaxis()->SetTitle("counts");
	304	fOutput->Add(fHistJetMCPtCorrPt[i]);
305	}
306	}
307
308	if (fIsEmbedded) {
309	histname = "fHistJetPtMCPt_";
310	histname += i;
311	fHistJetPtMCPt[i] = new TH2F(histname.Data(), histname.Data(), fNbins, fMinBinPt, fMaxBinPt, fNbins, fMinBinPt, fMaxBinPt);
312	fHistJetPtMCPt[i]->GetXaxis()->SetTitle("p_{T} (GeV/c)");
313	fHistJetPtMCPt[i]->GetYaxis()->SetTitle("p_{T}^{MC} (GeV/c)");
314	fHistJetPtMCPt[i]->GetZaxis()->SetTitle("counts");
315	fOutput->Add(fHistJetPtMCPt[i]);
316	}
317	}
318	}
319
320	//________________________________________________________________________
321	void AliAnalysisTaskSAJF::UserCreateOutputObjects()
322	{
323	// Create user output.
324
7cd832c7	325	AliAnalysisTaskEmcalJetDev::UserCreateOutputObjects();
85d48773	326
	327	if (fHistoType == 0)
	328	AllocateTHX();
	329	else
	330	AllocateTHnSparse();
43032ce2	331
	332	for (Int_t i = 0; i < 4; i++) {
	333	TString histname;
	334
7cd832c7	335	if (fParticleCollArray.GetEntriesFast()>0) {
43032ce2	336	histname = "fHistTracksJetPt_";
05077f28	337	histname += i;
43032ce2	338	fHistTracksJetPt[i] = new TH2F(histname.Data(), histname.Data(), fNbins / 2, fMinBinPt, fMaxBinPt / 2, fNbins, fMinBinPt, fMaxBinPt);
	339	fHistTracksJetPt[i]->GetXaxis()->SetTitle("p_{T,track} (GeV/c)");
	340	fHistTracksJetPt[i]->GetYaxis()->SetTitle("p_{T,jet} (GeV/c)");
	341	fHistTracksJetPt[i]->GetZaxis()->SetTitle("counts");
	342	fOutput->Add(fHistTracksJetPt[i]);
	343
	344	histname = "fHistTracksPtDist_";
05077f28	345	histname += i;
43032ce2	346	fHistTracksPtDist[i] = new TH2F(histname.Data(), histname.Data(), fNbins / 2, fMinBinPt, fMaxBinPt / 2, 100, 0, 5);
	347	fHistTracksPtDist[i]->GetXaxis()->SetTitle("p_{T,track} (GeV/c)");
	348	fHistTracksPtDist[i]->GetYaxis()->SetTitle("d");
	349	fHistTracksPtDist[i]->GetZaxis()->SetTitle("counts");
	350	fOutput->Add(fHistTracksPtDist[i]);
05077f28	351	}
05077f28	352
7cd832c7	353	if (fClusterCollArray.GetEntriesFast()>0) {
a487deae	354	histname = "fHistClustersJetPt_";
e44e8726	355	histname += i;
a487deae	356	fHistClustersJetPt[i] = new TH2F(histname.Data(), histname.Data(), fNbins / 2, fMinBinPt, fMaxBinPt / 2, fNbins, fMinBinPt, fMaxBinPt);
	357	fHistClustersJetPt[i]->GetXaxis()->SetTitle("p_{T,clus} (GeV/c)");
	358	fHistClustersJetPt[i]->GetYaxis()->SetTitle("p_{T,jet} (GeV/c)");
	359	fHistClustersJetPt[i]->GetZaxis()->SetTitle("counts");
	360	fOutput->Add(fHistClustersJetPt[i]);
7030f36f	361
	362	histname = "fHistClustersPtDist_";
	363	histname += i;
	364	fHistClustersPtDist[i] = new TH2F(histname.Data(), histname.Data(), fNbins / 2, fMinBinPt, fMaxBinPt / 2, 100, 0, 5);
	365	fHistClustersPtDist[i]->GetXaxis()->SetTitle("p_{T,clus} (GeV/c)");
	366	fHistClustersPtDist[i]->GetYaxis()->SetTitle("d");
	367	fHistClustersPtDist[i]->GetZaxis()->SetTitle("counts");
	368	fOutput->Add(fHistClustersPtDist[i]);
e44e8726	369	}
e82e282c	370	}
e82e282c	371
25283b37	372	PostData(1, fOutput); // Post data for ALL output slots >0 here, to get at least an empty histogram
25283b37	373
2bee31e9	374	}
2bee31e9	375
55264f20	376	//________________________________________________________________________
6fd5039f	377	Bool_t AliAnalysisTaskSAJF::FillHistograms()
25283b37	378	{
16d143bd	379	// Fill histograms.
16d143bd	380
a487deae	381	if (!fJets) {
	382	AliError(Form("%s - Jet array not provided, returning...", GetName()));
	383	return kFALSE;
	384	}
	385
3fe08cdb	386	for (Int_t ij = 0; ij < fJets->GetEntriesFast(); ij++) {
6fd5039f	387
3fe08cdb	388	AliEmcalJet* jet = static_cast<AliEmcalJet*>(fJets->At(ij));
6fd5039f	389
	390	if (!jet) {
	391	AliError(Form("Could not receive jet %d", ij));
	392	continue;
43032ce2	393	}
6fd5039f	394
3fe08cdb	395	if (!AcceptJet(jet))
3fe08cdb	396	continue;
226f511d	397
83888eef	398	Float_t ptLeading = GetLeadingHadronPt(jet);
43032ce2	399	Float_t corrPt = jet->Pt() - fRhoVal * jet->Area();
e44e8726	400
43032ce2	401	// Fill THnSparse
85d48773	402	Double_t ep = jet->Phi() - fEPV0;
	403	while (ep < 0) ep += TMath::Pi();
	404	while (ep >= TMath::Pi()) ep -= TMath::Pi();
	405
	406	FillJetHisto(fCent, ep, jet->Eta(), jet->Phi(), jet->Pt(), jet->MCPt(), corrPt, jet->Area(),
3fe08cdb	407	jet->NEF(), ptLeading/jet->Pt(), jet->GetNumberOfConstituents(), ptLeading);
a825589f	408
6e8d91c9	409	if (fTracks) {
	410	for (Int_t it = 0; it < jet->GetNumberOfTracks(); it++) {
	411	AliVParticle *track = jet->TrackAt(it, fTracks);
a825589f	412	if (track) {
a487deae	413	fHistTracksJetPt[fCentBin]->Fill(track->Pt(), jet->Pt());
7030f36f	414	Double_t dist = TMath::Sqrt((track->Eta() - jet->Eta()) * (track->Eta() - jet->Eta()) + (track->Phi() - jet->Phi()) * (track->Phi() - jet->Phi()));
7030f36f	415	fHistTracksPtDist[fCentBin]->Fill(track->Pt(), dist);
a825589f	416	}
6e8d91c9	417	}
35789a2d	418	}
a55e4f1d	419
6e8d91c9	420	if (fCaloClusters) {
	421	for (Int_t ic = 0; ic < jet->GetNumberOfClusters(); ic++) {
	422	AliVCluster *cluster = jet->ClusterAt(ic, fCaloClusters);
	423
	424	if (cluster) {
	425	TLorentzVector nPart;
	426	cluster->GetMomentum(nPart, fVertex);
43032ce2	427
a487deae	428	fHistClustersJetPt[fCentBin]->Fill(nPart.Pt(), jet->Pt());
7030f36f	429	Double_t dist = TMath::Sqrt((nPart.Eta() - jet->Eta()) * (nPart.Eta() - jet->Eta()) + (nPart.Phi() - jet->Phi()) * (nPart.Phi() - jet->Phi()));
7030f36f	430	fHistClustersPtDist[fCentBin]->Fill(nPart.Pt(), dist);
6e8d91c9	431	}
55264f20	432	}
f0a0fd33	433	}
25283b37	434	} //jet loop
3fe08cdb	435
	436	return kTRUE;
	437	}
	438
	439	//________________________________________________________________________
	440	void AliAnalysisTaskSAJF::FillJetHisto(Double_t cent, Double_t ep, Double_t eta, Double_t phi, Double_t pt, Double_t MCpt, Double_t corrpt, Double_t area,
	441	Double_t NEF, Double_t z, Int_t n, Double_t leadingpt)
	442	{
85d48773	443	if (fHistoType == 0) {
	444	fHistJetPtEtaPhi[fCentBin]->Fill(pt,eta,phi);
	445	fHistJetPtArea[fCentBin]->Fill(pt,area);
	446	fHistJetPtEP[fCentBin]->Fill(pt,ep);
	447	fHistJetPtNEF[fCentBin]->Fill(pt,NEF);
	448	fHistJetPtZ[fCentBin]->Fill(pt,z);
	449	fHistJetPtLeadingPartPt[fCentBin]->Fill(pt,leadingpt);
7cd832c7	450	if (fHistJetCorrPtEtaPhi[fCentBin]) {
85d48773	451	fHistJetCorrPtEtaPhi[fCentBin]->Fill(corrpt,eta,phi);
	452	fHistJetCorrPtArea[fCentBin]->Fill(corrpt,area);
	453	fHistJetCorrPtEP[fCentBin]->Fill(corrpt,ep);
	454	fHistJetCorrPtNEF[fCentBin]->Fill(corrpt,NEF);
	455	fHistJetCorrPtZ[fCentBin]->Fill(corrpt,z);
	456	fHistJetCorrPtLeadingPartPt[fCentBin]->Fill(corrpt,leadingpt);
	457	fHistJetPtCorrPt[fCentBin]->Fill(pt,corrpt);
	458	if (fIsEmbedded)
	459	fHistJetMCPtCorrPt[fCentBin]->Fill(MCpt,corrpt);
	460	}
	461	if (fIsEmbedded)
	462	fHistJetPtMCPt[fCentBin]->Fill(pt,MCpt);
3fe08cdb	463	}
85d48773	464	else {
	465
	466	Double_t contents[20]={0};
	467
	468	for (Int_t i = 0; i < fHistJetObservables->GetNdimensions(); i++) {
	469	TString title(fHistJetObservables->GetAxis(i)->GetTitle());
	470	if (title=="Centrality (%)")
	471	contents[i] = cent;
	472	else if (title=="#phi_{jet} - #psi_{RP}")
	473	contents[i] = ep;
	474	else if (title=="#eta")
	475	contents[i] = eta;
	476	else if (title=="#phi_{jet} (rad)")
	477	contents[i] = phi;
	478	else if (title=="p_{T} (GeV/c)")
	479	contents[i] = pt;
	480	else if (title=="p_{T}^{MC} (GeV/c)")
	481	contents[i] = MCpt;
	482	else if (title=="p_{T}^{corr} (GeV/c)")
	483	contents[i] = corrpt;
	484	else if (title=="A_{jet}")
	485	contents[i] = area;
	486	else if (title=="NEF")
	487	contents[i] = NEF;
	488	else if (title=="Z")
	489	contents[i] = z;
	490	else if (title=="No. of constituents")
	491	contents[i] = n;
	492	else if (title=="p_{T,particle}^{leading} (GeV/c)")
	493	contents[i] = leadingpt;
	494	else
	495	AliWarning(Form("Unable to fill dimension %s!",title.Data()));
	496	}
3fe08cdb	497
85d48773	498	fHistJetObservables->Fill(contents);
85d48773	499	}
55264f20	500	}