[u/mrichter/AliRoot.git] / PWGJE / EMCALJetTasks / AliAnalysisTaskEmcalJetSample.cxx

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

#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"

ClassImp(AliAnalysisTaskEmcalJetSample)

//________________________________________________________________________
AliAnalysisTaskEmcalJetSample::AliAnalysisTaskEmcalJetSample() : 
  AliAnalysisTaskEmcalJet("AliAnalysisTaskEmcalJetSample", kTRUE),
  fHistTracksPt(0),
  fHistClustersPt(0),
  fHistLeadingJetPt(0),
  fHistJetsPhiEta(0),
  fHistJetsPtArea(0),
  fHistJetsPtLeadHad(0),
  fHistJetsCorrPtArea(0),
  fHistPtDEtaDPhiTrackClus(0),
  fHistPtDEtaDPhiClusTrack(0),
  fJetsCont(0),
  fTracksCont(0),
  fCaloClustersCont(0)

{
  // Default constructor.

  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;
    fHistJetsPhiEta[i] = 0;
    fHistJetsPtArea[i] = 0;
    fHistJetsPtLeadHad[i] = 0;
    fHistJetsCorrPtArea[i] = 0;
  }

  SetMakeGeneralHistograms(kTRUE);
}

//________________________________________________________________________
AliAnalysisTaskEmcalJetSample::AliAnalysisTaskEmcalJetSample(const char *name) : 
  AliAnalysisTaskEmcalJet(name, kTRUE),
  fHistTracksPt(0),
  fHistClustersPt(0),
  fHistLeadingJetPt(0),
  fHistJetsPhiEta(0),
  fHistJetsPtArea(0),
  fHistJetsPtLeadHad(0),
  fHistJetsCorrPtArea(0),
  fHistPtDEtaDPhiTrackClus(0),
  fHistPtDEtaDPhiClusTrack(0),
  fJetsCont(0),
  fTracksCont(0),
  fCaloClustersCont(0)
{
  // Standard constructor.

  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;
    fHistJetsPhiEta[i] = 0;
    fHistJetsPtArea[i] = 0;
    fHistJetsPtLeadHad[i] = 0;
    fHistJetsCorrPtArea[i] = 0;
  }

  SetMakeGeneralHistograms(kTRUE);
}

//________________________________________________________________________
AliAnalysisTaskEmcalJetSample::~AliAnalysisTaskEmcalJetSample()
{
  // Destructor.
}

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

  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);
  }
  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]);
    }

    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]);
      }
    }
  }

  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);

  PostData(1, fOutput); // Post data for ALL output slots > 0 here.
}

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

  if (fTracksCont) {
    AliVTrack *track = static_cast<AliVTrack*>(fTracksCont->GetNextAcceptParticle(0)); 
    while(track) {
      fHistTracksPt[fCentBin]->Fill(track->Pt()); 
      track = static_cast<AliVTrack*>(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();
    }
  }

  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 = fJetsCont->GetLeadingHadronPt(jet);
      fHistJetsPtLeadHad[fCentBin]->Fill(jet->Pt(), ptLeading);

      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()
{
  // Run analysis code here, if needed. It will be executed before FillHistograms().

  return kTRUE;  // If return kFALSE FillHistogram() will NOT be executed.
}

//________________________________________________________________________
void AliAnalysisTaskEmcalJetSample::Terminate(Option_t *) 
{
  // Called once at the end of the analysis.
}
Commit	Line	Data
8628b70c	1	// $Id$
9f52c61f	2	//
	3	// Jet sample analysis task.
	4	//
	5	// Author: S.Aiola
	6
	7	#include <TClonesArray.h>
	8	#include <TH1F.h>
	9	#include <TH2F.h>
0925db21	10	#include <TH3F.h>
9f52c61f	11	#include <TList.h>
	12	#include <TLorentzVector.h>
	13
	14	#include "AliVCluster.h"
4bb8c6df	15	#include "AliAODCaloCluster.h"
0925db21	16	#include "AliESDCaloCluster.h"
9f52c61f	17	#include "AliVTrack.h"
	18	#include "AliEmcalJet.h"
	19	#include "AliRhoParameter.h"
	20	#include "AliLog.h"
9239b066	21	#include "AliJetContainer.h"
8612dfc8	22	#include "AliParticleContainer.h"
8612dfc8	23	#include "AliClusterContainer.h"
0925db21	24	#include "AliPicoTrack.h"
9f52c61f	25
	26	#include "AliAnalysisTaskEmcalJetSample.h"
	27
	28	ClassImp(AliAnalysisTaskEmcalJetSample)
	29
	30	//________________________________________________________________________
	31	AliAnalysisTaskEmcalJetSample::AliAnalysisTaskEmcalJetSample() :
8612dfc8	32	AliAnalysisTaskEmcalJet("AliAnalysisTaskEmcalJetSample", kTRUE),
4bb8c6df	33	fHistTracksPt(0),
	34	fHistClustersPt(0),
	35	fHistLeadingJetPt(0),
	36	fHistJetsPhiEta(0),
	37	fHistJetsPtArea(0),
	38	fHistJetsPtLeadHad(0),
	39	fHistJetsCorrPtArea(0),
0925db21	40	fHistPtDEtaDPhiTrackClus(0),
0925db21	41	fHistPtDEtaDPhiClusTrack(0),
8612dfc8	42	fJetsCont(0),
	43	fTracksCont(0),
	44	fCaloClustersCont(0)
9f52c61f	45
	46	{
	47	// Default constructor.
	48
4bb8c6df	49	fHistTracksPt = new TH1*[fNcentBins];
	50	fHistClustersPt = new TH1*[fNcentBins];
	51	fHistLeadingJetPt = new TH1*[fNcentBins];
	52	fHistJetsPhiEta = new TH2*[fNcentBins];
	53	fHistJetsPtArea = new TH2*[fNcentBins];
	54	fHistJetsPtLeadHad = new TH2*[fNcentBins];
	55	fHistJetsCorrPtArea = new TH2*[fNcentBins];
	56
	57	for (Int_t i = 0; i < fNcentBins; i++) {
9f52c61f	58	fHistTracksPt[i] = 0;
	59	fHistClustersPt[i] = 0;
	60	fHistLeadingJetPt[i] = 0;
	61	fHistJetsPhiEta[i] = 0;
	62	fHistJetsPtArea[i] = 0;
	63	fHistJetsPtLeadHad[i] = 0;
	64	fHistJetsCorrPtArea[i] = 0;
	65	}
	66
	67	SetMakeGeneralHistograms(kTRUE);
	68	}
	69
	70	//________________________________________________________________________
	71	AliAnalysisTaskEmcalJetSample::AliAnalysisTaskEmcalJetSample(const char *name) :
8612dfc8	72	AliAnalysisTaskEmcalJet(name, kTRUE),
4bb8c6df	73	fHistTracksPt(0),
	74	fHistClustersPt(0),
	75	fHistLeadingJetPt(0),
	76	fHistJetsPhiEta(0),
	77	fHistJetsPtArea(0),
	78	fHistJetsPtLeadHad(0),
	79	fHistJetsCorrPtArea(0),
0925db21	80	fHistPtDEtaDPhiTrackClus(0),
0925db21	81	fHistPtDEtaDPhiClusTrack(0),
8612dfc8	82	fJetsCont(0),
	83	fTracksCont(0),
	84	fCaloClustersCont(0)
9f52c61f	85	{
	86	// Standard constructor.
	87
4bb8c6df	88	fHistTracksPt = new TH1*[fNcentBins];
	89	fHistClustersPt = new TH1*[fNcentBins];
	90	fHistLeadingJetPt = new TH1*[fNcentBins];
	91	fHistJetsPhiEta = new TH2*[fNcentBins];
	92	fHistJetsPtArea = new TH2*[fNcentBins];
	93	fHistJetsPtLeadHad = new TH2*[fNcentBins];
	94	fHistJetsCorrPtArea = new TH2*[fNcentBins];
9239b066	95
4bb8c6df	96	for (Int_t i = 0; i < fNcentBins; i++) {
9f52c61f	97	fHistTracksPt[i] = 0;
	98	fHistClustersPt[i] = 0;
	99	fHistLeadingJetPt[i] = 0;
	100	fHistJetsPhiEta[i] = 0;
	101	fHistJetsPtArea[i] = 0;
	102	fHistJetsPtLeadHad[i] = 0;
	103	fHistJetsCorrPtArea[i] = 0;
	104	}
	105
	106	SetMakeGeneralHistograms(kTRUE);
	107	}
	108
	109	//________________________________________________________________________
	110	AliAnalysisTaskEmcalJetSample::~AliAnalysisTaskEmcalJetSample()
	111	{
	112	// Destructor.
	113	}
	114
	115	//________________________________________________________________________
	116	void AliAnalysisTaskEmcalJetSample::UserCreateOutputObjects()
	117	{
	118	// Create user output.
	119
	120	AliAnalysisTaskEmcalJet::UserCreateOutputObjects();
	121
4bb8c6df	122	fJetsCont = GetJetContainer(0);
	123	if(fJetsCont) { //get particles and clusters connected to jets
	124	fTracksCont = fJetsCont->GetParticleContainer();
	125	fCaloClustersCont = fJetsCont->GetClusterContainer();
	126	} else { //no jets, just analysis tracks and clusters
	127	fTracksCont = GetParticleContainer(0);
	128	fCaloClustersCont = GetClusterContainer(0);
	129	}
	130	fTracksCont->SetClassName("AliVTrack");
05f22af9	131	fCaloClustersCont->SetClassName("AliVCluster");
8612dfc8	132
9f52c61f	133	TString histname;
9f52c61f	134
4bb8c6df	135	for (Int_t i = 0; i < fNcentBins; i++) {
9239b066	136	if (fParticleCollArray.GetEntriesFast()>0) {
9f52c61f	137	histname = "fHistTracksPt_";
	138	histname += i;
	139	fHistTracksPt[i] = new TH1F(histname.Data(), histname.Data(), fNbins / 2, fMinBinPt, fMaxBinPt / 2);
	140	fHistTracksPt[i]->GetXaxis()->SetTitle("p_{T,track} (GeV/c)");
	141	fHistTracksPt[i]->GetYaxis()->SetTitle("counts");
	142	fOutput->Add(fHistTracksPt[i]);
	143	}
	144
9239b066	145	if (fClusterCollArray.GetEntriesFast()>0) {
9f52c61f	146	histname = "fHistClustersPt_";
	147	histname += i;
	148	fHistClustersPt[i] = new TH1F(histname.Data(), histname.Data(), fNbins / 2, fMinBinPt, fMaxBinPt / 2);
	149	fHistClustersPt[i]->GetXaxis()->SetTitle("p_{T,clus} (GeV/c)");
	150	fHistClustersPt[i]->GetYaxis()->SetTitle("counts");
	151	fOutput->Add(fHistClustersPt[i]);
	152	}
	153
9239b066	154	if (fJetCollArray.GetEntriesFast()>0) {
9f52c61f	155	histname = "fHistLeadingJetPt_";
	156	histname += i;
	157	fHistLeadingJetPt[i] = new TH1F(histname.Data(), histname.Data(), fNbins, fMinBinPt, fMaxBinPt);
	158	fHistLeadingJetPt[i]->GetXaxis()->SetTitle("p_{T}^{raw} (GeV/c)");
	159	fHistLeadingJetPt[i]->GetYaxis()->SetTitle("counts");
	160	fOutput->Add(fHistLeadingJetPt[i]);
	161
	162	histname = "fHistJetsPhiEta_";
	163	histname += i;
	164	fHistJetsPhiEta[i] = new TH2F(histname.Data(), histname.Data(), 50, -1, 1, 101, 0, TMath::Pi()*2 + TMath::Pi()/200);
	165	fHistJetsPhiEta[i]->GetXaxis()->SetTitle("#eta");
	166	fHistJetsPhiEta[i]->GetYaxis()->SetTitle("#phi");
	167	fOutput->Add(fHistJetsPhiEta[i]);
	168
	169	histname = "fHistJetsPtArea_";
	170	histname += i;
	171	fHistJetsPtArea[i] = new TH2F(histname.Data(), histname.Data(), fNbins, fMinBinPt, fMaxBinPt, 30, 0, 3);
	172	fHistJetsPtArea[i]->GetXaxis()->SetTitle("p_{T}^{raw} (GeV/c)");
	173	fHistJetsPtArea[i]->GetYaxis()->SetTitle("area");
	174	fOutput->Add(fHistJetsPtArea[i]);
	175
	176	histname = "fHistJetsPtLeadHad_";
	177	histname += i;
	178	fHistJetsPtLeadHad[i] = new TH2F(histname.Data(), histname.Data(), fNbins, fMinBinPt, fMaxBinPt, fNbins / 2, fMinBinPt, fMaxBinPt / 2);
	179	fHistJetsPtLeadHad[i]->GetXaxis()->SetTitle("p_{T}^{raw} (GeV/c)");
	180	fHistJetsPtLeadHad[i]->GetYaxis()->SetTitle("p_{T,lead} (GeV/c)");
	181	fHistJetsPtLeadHad[i]->GetZaxis()->SetTitle("counts");
	182	fOutput->Add(fHistJetsPtLeadHad[i]);
	183
9239b066	184	if (!(GetJetContainer()->GetRhoName().IsNull())) {
9f52c61f	185	histname = "fHistJetsCorrPtArea_";
	186	histname += i;
	187	fHistJetsCorrPtArea[i] = new TH2F(histname.Data(), histname.Data(), fNbins*2, -fMaxBinPt, fMaxBinPt, 30, 0, 3);
	188	fHistJetsCorrPtArea[i]->GetXaxis()->SetTitle("p_{T}^{corr} [GeV/c]");
	189	fHistJetsCorrPtArea[i]->GetYaxis()->SetTitle("area");
	190	fOutput->Add(fHistJetsCorrPtArea[i]);
	191	}
	192	}
	193	}
0925db21	194
	195	histname = "fHistPtDEtaDPhiTrackClus";
	196	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);
	197	fOutput->Add(fHistPtDEtaDPhiTrackClus);
	198
	199	histname = "fHistPtDEtaDPhiClusTrack";
	200	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);
	201	fOutput->Add(fHistPtDEtaDPhiClusTrack);
	202
9f52c61f	203	PostData(1, fOutput); // Post data for ALL output slots > 0 here.
	204	}
	205
	206	//________________________________________________________________________
	207	Bool_t AliAnalysisTaskEmcalJetSample::FillHistograms()
	208	{
	209	// Fill histograms.
	210
8612dfc8	211	if (fTracksCont) {
4bb8c6df	212	AliVTrack track = static_cast<AliVTrack>(fTracksCont->GetNextAcceptParticle(0));
8612dfc8	213	while(track) {
9f52c61f	214	fHistTracksPt[fCentBin]->Fill(track->Pt());
4bb8c6df	215	track = static_cast<AliVTrack*>(fTracksCont->GetNextAcceptParticle());
9f52c61f	216	}
	217	}
	218
8612dfc8	219	if (fCaloClustersCont) {
0925db21	220	AliVCluster *cluster = fCaloClustersCont->GetNextAcceptCluster(0);
8612dfc8	221	while(cluster) {
9f52c61f	222	TLorentzVector nPart;
	223	cluster->GetMomentum(nPart, fVertex);
	224	fHistClustersPt[fCentBin]->Fill(nPart.Pt());
0925db21	225
0925db21	226	cluster = fCaloClustersCont->GetNextAcceptCluster();
9f52c61f	227	}
	228	}
	229
8612dfc8	230	if (fJetsCont) {
	231	AliEmcalJet *jet = fJetsCont->GetNextAcceptJet(0);
	232	while(jet) {
9239b066	233
9f52c61f	234	fHistJetsPtArea[fCentBin]->Fill(jet->Pt(), jet->Area());
	235	fHistJetsPhiEta[fCentBin]->Fill(jet->Eta(), jet->Phi());
	236
22339f55	237	Float_t ptLeading = fJetsCont->GetLeadingHadronPt(jet);
9f52c61f	238	fHistJetsPtLeadHad[fCentBin]->Fill(jet->Pt(), ptLeading);
9f52c61f	239
388f9cd5	240	if (fHistJetsCorrPtArea[fCentBin]) {
	241	Float_t corrPt = jet->Pt() - fJetsCont->GetRhoVal() * jet->Area();
	242	fHistJetsCorrPtArea[fCentBin]->Fill(corrPt, jet->Area());
	243	}
8612dfc8	244	jet = fJetsCont->GetNextAcceptJet();
9f52c61f	245	}
8612dfc8	246
8612dfc8	247	jet = fJetsCont->GetLeadingJet();
6b840003	248	if(jet) fHistLeadingJetPt[fCentBin]->Fill(jet->Pt());
9f52c61f	249	}
9f52c61f	250
0925db21	251	CheckClusTrackMatching();
0925db21	252
9f52c61f	253	return kTRUE;
	254	}
	255
0925db21	256	//________________________________________________________________________
	257	void AliAnalysisTaskEmcalJetSample::CheckClusTrackMatching()
	258	{
	259
	260	if(!fTracksCont \|\| !fCaloClustersCont)
	261	return;
	262
	263	Double_t deta = 999;
	264	Double_t dphi = 999;
	265
	266	//Get closest cluster to track
	267	AliVTrack track = static_cast<AliVTrack>(fTracksCont->GetNextAcceptParticle(0));
	268	while(track) {
	269	//Get matched cluster
	270	Int_t emc1 = track->GetEMCALcluster();
	271	if(fCaloClustersCont && emc1>=0) {
	272	AliVCluster *clusMatch = fCaloClustersCont->GetCluster(emc1);
	273	if(clusMatch) {
	274	AliPicoTrack::GetEtaPhiDiff(track, clusMatch, dphi, deta);
	275	fHistPtDEtaDPhiTrackClus->Fill(track->Pt(),deta,dphi);
	276	}
	277	}
	278	track = static_cast<AliVTrack*>(fTracksCont->GetNextAcceptParticle());
	279	}
	280
	281	//Get closest track to cluster
	282	AliVCluster *cluster = fCaloClustersCont->GetNextAcceptCluster(0);
	283	while(cluster) {
	284	TLorentzVector nPart;
	285	cluster->GetMomentum(nPart, fVertex);
	286	fHistClustersPt[fCentBin]->Fill(nPart.Pt());
	287
	288	//Get matched track
	289	AliVTrack *mt = NULL;
	290	AliAODCaloCluster acl = dynamic_cast<AliAODCaloCluster>(cluster);
	291	if(acl) {
	292	if(acl->GetNTracksMatched()>1)
	293	mt = static_cast<AliVTrack*>(acl->GetTrackMatched(0));
	294	}
	295	else {
	296	AliESDCaloCluster ecl = dynamic_cast<AliESDCaloCluster>(cluster);
	297	Int_t im = ecl->GetTrackMatchedIndex();
	298	if(fTracksCont && im>=0) {
	299	mt = static_cast<AliVTrack*>(fTracksCont->GetParticle(im));
	300	}
	301	}
	302	if(mt) {
	303	AliPicoTrack::GetEtaPhiDiff(mt, cluster, dphi, deta);
	304	fHistPtDEtaDPhiClusTrack->Fill(nPart.Pt(),deta,dphi);
	305
	306	/* //debugging
	307	if(mt->IsEMCAL()) {
	308	Int_t emc1 = mt->GetEMCALcluster();
	309	Printf("current id: %d emc1: %d",fCaloClustersCont->GetCurrentID(),emc1);
	310	AliVCluster *clm = fCaloClustersCont->GetCluster(emc1);
	311	AliPicoTrack::GetEtaPhiDiff(mt, clm, dphi, deta);
	312	Printf("deta: %f dphi: %f",deta,dphi);
	313	}
	314	*/
	315	}
	316	cluster = fCaloClustersCont->GetNextAcceptCluster();
	317	}
	318	}
	319
8612dfc8	320	//________________________________________________________________________
	321	void AliAnalysisTaskEmcalJetSample::ExecOnce() {
	322
	323	AliAnalysisTaskEmcalJet::ExecOnce();
	324
	325	if (fJetsCont && fJetsCont->GetArray() == 0) fJetsCont = 0;
	326	if (fTracksCont && fTracksCont->GetArray() == 0) fTracksCont = 0;
	327	if (fCaloClustersCont && fCaloClustersCont->GetArray() == 0) fCaloClustersCont = 0;
	328
	329	}
	330
9f52c61f	331	//________________________________________________________________________
	332	Bool_t AliAnalysisTaskEmcalJetSample::Run()
	333	{
	334	// Run analysis code here, if needed. It will be executed before FillHistograms().
	335
	336	return kTRUE; // If return kFALSE FillHistogram() will NOT be executed.
	337	}
	338
	339	//________________________________________________________________________
	340	void AliAnalysisTaskEmcalJetSample::Terminate(Option_t *)
	341	{
	342	// Called once at the end of the analysis.
	343	}