[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("AliAODCaloCluster");

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

      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>
bced96f7	10	#include <TH3F.h>
9f52c61f	11	#include <TList.h>
	12	#include <TLorentzVector.h>
	13
	14	#include "AliVCluster.h"
a7626402	15	#include "AliAODCaloCluster.h"
bced96f7	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"
bced96f7	24	#include "AliPicoTrack.h"
9f52c61f	25
	26	#include "AliAnalysisTaskEmcalJetSample.h"
	27
	28	ClassImp(AliAnalysisTaskEmcalJetSample)
	29
	30	//________________________________________________________________________
	31	AliAnalysisTaskEmcalJetSample::AliAnalysisTaskEmcalJetSample() :
8612dfc8	32	AliAnalysisTaskEmcalJet("AliAnalysisTaskEmcalJetSample", kTRUE),
a7626402	33	fHistTracksPt(0),
	34	fHistClustersPt(0),
	35	fHistLeadingJetPt(0),
	36	fHistJetsPhiEta(0),
	37	fHistJetsPtArea(0),
	38	fHistJetsPtLeadHad(0),
	39	fHistJetsCorrPtArea(0),
bced96f7	40	fHistPtDEtaDPhiTrackClus(0),
bced96f7	41	fHistPtDEtaDPhiClusTrack(0),
8612dfc8	42	fJetsCont(0),
	43	fTracksCont(0),
	44	fCaloClustersCont(0)
9f52c61f	45
	46	{
	47	// Default constructor.
	48
a7626402	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),
a7626402	73	fHistTracksPt(0),
	74	fHistClustersPt(0),
	75	fHistLeadingJetPt(0),
	76	fHistJetsPhiEta(0),
	77	fHistJetsPtArea(0),
	78	fHistJetsPtLeadHad(0),
	79	fHistJetsCorrPtArea(0),
bced96f7	80	fHistPtDEtaDPhiTrackClus(0),
bced96f7	81	fHistPtDEtaDPhiClusTrack(0),
8612dfc8	82	fJetsCont(0),
	83	fTracksCont(0),
	84	fCaloClustersCont(0)
9f52c61f	85	{
	86	// Standard constructor.
	87
a7626402	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
a7626402	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
a7626402	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");
	131	fCaloClustersCont->SetClassName("AliAODCaloCluster");
8612dfc8	132
9f52c61f	133	TString histname;
9f52c61f	134
a7626402	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	}
bced96f7	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) {
a7626402	212	AliVTrack track = static_cast<AliVTrack>(fTracksCont->GetNextAcceptParticle(0));
8612dfc8	213	while(track) {
9f52c61f	214	fHistTracksPt[fCentBin]->Fill(track->Pt());
a7626402	215	track = static_cast<AliVTrack*>(fTracksCont->GetNextAcceptParticle());
9f52c61f	216	}
	217	}
	218
8612dfc8	219	if (fCaloClustersCont) {
bced96f7	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());
bced96f7	225
bced96f7	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
8612dfc8	240	Float_t corrPt = jet->Pt() - fJetsCont->GetRhoVal() * jet->Area();
	241	fHistJetsCorrPtArea[fCentBin]->Fill(corrPt, jet->Area());
	242
	243	jet = fJetsCont->GetNextAcceptJet();
9f52c61f	244	}
8612dfc8	245
8612dfc8	246	jet = fJetsCont->GetLeadingJet();
6b840003	247	if(jet) fHistLeadingJetPt[fCentBin]->Fill(jet->Pt());
9f52c61f	248	}
9f52c61f	249
bced96f7	250	CheckClusTrackMatching();
bced96f7	251
9f52c61f	252	return kTRUE;
	253	}
	254
bced96f7	255	//________________________________________________________________________
	256	void AliAnalysisTaskEmcalJetSample::CheckClusTrackMatching()
	257	{
	258
	259	if(!fTracksCont \|\| !fCaloClustersCont)
	260	return;
	261
	262	Double_t deta = 999;
	263	Double_t dphi = 999;
	264
	265	//Get closest cluster to track
	266	AliVTrack track = static_cast<AliVTrack>(fTracksCont->GetNextAcceptParticle(0));
	267	while(track) {
	268	//Get matched cluster
	269	Int_t emc1 = track->GetEMCALcluster();
	270	if(fCaloClustersCont && emc1>=0) {
	271	AliVCluster *clusMatch = fCaloClustersCont->GetCluster(emc1);
	272	if(clusMatch) {
	273	AliPicoTrack::GetEtaPhiDiff(track, clusMatch, dphi, deta);
	274	fHistPtDEtaDPhiTrackClus->Fill(track->Pt(),deta,dphi);
	275	}
	276	}
	277	track = static_cast<AliVTrack*>(fTracksCont->GetNextAcceptParticle());
	278	}
	279
	280	//Get closest track to cluster
	281	AliVCluster *cluster = fCaloClustersCont->GetNextAcceptCluster(0);
	282	while(cluster) {
	283	TLorentzVector nPart;
	284	cluster->GetMomentum(nPart, fVertex);
	285	fHistClustersPt[fCentBin]->Fill(nPart.Pt());
	286
	287	//Get matched track
	288	AliVTrack *mt = NULL;
	289	AliAODCaloCluster acl = dynamic_cast<AliAODCaloCluster>(cluster);
	290	if(acl) {
	291	if(acl->GetNTracksMatched()>1)
	292	mt = static_cast<AliVTrack*>(acl->GetTrackMatched(0));
	293	}
	294	else {
	295	AliESDCaloCluster ecl = dynamic_cast<AliESDCaloCluster>(cluster);
	296	Int_t im = ecl->GetTrackMatchedIndex();
	297	if(fTracksCont && im>=0) {
	298	mt = static_cast<AliVTrack*>(fTracksCont->GetParticle(im));
	299	}
	300	}
	301	if(mt) {
	302	AliPicoTrack::GetEtaPhiDiff(mt, cluster, dphi, deta);
	303	fHistPtDEtaDPhiClusTrack->Fill(nPart.Pt(),deta,dphi);
	304
	305	/* //debugging
	306	if(mt->IsEMCAL()) {
	307	Int_t emc1 = mt->GetEMCALcluster();
	308	Printf("current id: %d emc1: %d",fCaloClustersCont->GetCurrentID(),emc1);
	309	AliVCluster *clm = fCaloClustersCont->GetCluster(emc1);
	310	AliPicoTrack::GetEtaPhiDiff(mt, clm, dphi, deta);
	311	Printf("deta: %f dphi: %f",deta,dphi);
	312	}
	313	*/
	314	}
	315	cluster = fCaloClustersCont->GetNextAcceptCluster();
	316	}
	317	}
	318
8612dfc8	319	//________________________________________________________________________
	320	void AliAnalysisTaskEmcalJetSample::ExecOnce() {
	321
	322	AliAnalysisTaskEmcalJet::ExecOnce();
	323
	324	if (fJetsCont && fJetsCont->GetArray() == 0) fJetsCont = 0;
	325	if (fTracksCont && fTracksCont->GetArray() == 0) fTracksCont = 0;
	326	if (fCaloClustersCont && fCaloClustersCont->GetArray() == 0) fCaloClustersCont = 0;
	327
	328	}
	329
9f52c61f	330	//________________________________________________________________________
	331	Bool_t AliAnalysisTaskEmcalJetSample::Run()
	332	{
	333	// Run analysis code here, if needed. It will be executed before FillHistograms().
	334
	335	return kTRUE; // If return kFALSE FillHistogram() will NOT be executed.
	336	}
	337
	338	//________________________________________________________________________
	339	void AliAnalysisTaskEmcalJetSample::Terminate(Option_t *)
	340	{
	341	// Called once at the end of the analysis.
	342	}