[u/mrichter/AliRoot.git] / PWGGA / EMCALTasks / AliAnalysisTaskSAJF.cxx

// $Id$

#include <TChain.h>
#include <TClonesArray.h>
#include <TH1F.h>
#include <TH2F.h>
#include <TList.h>
#include <TLorentzVector.h>
#include <TParticle.h>

#include "AliAnalysisManager.h"
#include "AliCentrality.h"
#include "AliVCluster.h"
#include "AliESDtrack.h"
#include "AliEmcalJet.h"
#include "AliAODTrack.h"
#include "AliEmcalJet.h"
#include "AliVEventHandler.h"
#include "AliPicoTrack.h"

#include "AliAnalysisTaskSAJF.h"

ClassImp(AliAnalysisTaskSAJF)

//________________________________________________________________________
AliAnalysisTaskSAJF::AliAnalysisTaskSAJF() : 
  AliAnalysisTaskSE("AliAnalysisTaskSAJF"),
  fOutput(0),
  fTracksName("Tracks"),
  fCaloName("CaloClusters"),
  fJetsName("Jets"),
  fTrgClusName("ClustersL1GAMMAFEE"),
  fTracks(0),
  fCaloClusters(0),
  fJets(0),
  fTrgClusters(0),
  fCent(0),
  fHistCentrality(0),
  Ptbins(400),
  Ptlow(0),
  Ptup(200),
  Ebins(400),
  Elow(0),
  Eup(200)
{
  // Default constructor.
  fBranchNames="ESD:AliESDRun.,AliESDHeader.,PrimaryVertex.";

  for (Int_t i = 0; i < 4; i++) {
    fHistJetsE[i] = 0;
    fHistJetsNE[i] = 0;
    fHistJetsNEF[i] = 0;
    fHistJetsZ[i] = 0;
    fHistLeadingJetE[i] = 0;
    fHistTracksPtLJ[i] = 0;
    fHistClusELJ[i] = 0;
    fHistTracksPtBkg[i] = 0;
    fHistClusEBkg[i] = 0;
  }
}

//________________________________________________________________________
AliAnalysisTaskSAJF::AliAnalysisTaskSAJF(const char *name) : 
  AliAnalysisTaskSE(name),
  fOutput(0),
  fTracksName("Tracks"),
  fCaloName("CaloClusters"),
  fJetsName("Jets"),
  fTrgClusName("ClustersL1GAMMAFEE"),
  fTracks(0),
  fCaloClusters(0),
  fJets(0),
  fTrgClusters(0),
  fCent(0),
  fHistCentrality(0),
  Ptbins(400),
  Ptlow(0),
  Ptup(200),
  Ebins(400),
  Elow(0),
  Eup(200)
{
  // Standard constructor.

  fBranchNames="ESD:AliESDRun.,AliESDHeader.,PrimaryVertex.";

  for (Int_t i = 0; i < 4; i++) {
    fHistJetsE[i] = 0;
    fHistJetsNE[i] = 0;
    fHistJetsNEF[i] = 0;
    fHistJetsZ[i] = 0;
    fHistLeadingJetE[i] = 0;
    fHistTracksPtLJ[i] = 0;
    fHistClusELJ[i] = 0;
    fHistTracksPtBkg[i] = 0;
    fHistClusEBkg[i] = 0;
  }
}

//________________________________________________________________________
AliAnalysisTaskSAJF::~AliAnalysisTaskSAJF()
{
  // Destructor
}

//________________________________________________________________________
void AliAnalysisTaskSAJF::UserCreateOutputObjects()
{
  // Create histograms
  
  fOutput = new TList();
  fOutput->SetOwner();  // IMPORTANT!
  
  fHistCentrality = new TH1F("fHistCentrality","Event centrality distribution", Ebins, 0, 100);
  fHistCentrality->GetXaxis()->SetTitle("Centrality (%)");
  fHistCentrality->GetYaxis()->SetTitle("counts");
  fOutput->Add(fHistCentrality);

  TString histname;

  for (Int_t i = 0; i < 4; i++) {
    histname = "fHistJetsE_";
    histname += i;
    fHistJetsE[i] = new TH1F(histname.Data(), histname.Data(), Ebins, Elow, Eup);
    fHistJetsE[i]->GetXaxis()->SetTitle("E [GeV]");
    fHistJetsE[i]->GetYaxis()->SetTitle("counts");
    fOutput->Add(fHistJetsE[i]);
    
    histname = "fHistJetsNE_";
    histname += i;
    fHistJetsNE[i] = new TH1F(histname.Data(), histname.Data(), Ebins, Elow, Eup);
    fHistJetsNE[i]->GetXaxis()->SetTitle("E [GeV]");
    fHistJetsNE[i]->GetYaxis()->SetTitle("counts");
    fOutput->Add(fHistJetsNE[i]);
    
    histname = "fHistJetsNEF_";
    histname += i;
    fHistJetsNEF[i] = new TH1F(histname.Data(), histname.Data(), Ebins, 0, 1.2);
    fHistJetsNEF[i]->GetXaxis()->SetTitle("NEF");
    fHistJetsNEF[i]->GetYaxis()->SetTitle("counts");
    fOutput->Add(fHistJetsNEF[i]);

    histname = "fHistJetsZ_";
    histname += i;
    fHistJetsZ[i] = new TH1F(histname.Data(), histname.Data(), Ebins, 0, 1.2);
    fHistJetsZ[i]->GetXaxis()->SetTitle("Z");
    fHistJetsZ[i]->GetYaxis()->SetTitle("counts");
    fOutput->Add(fHistJetsZ[i]);

    histname = "fHistLeadingJetE_";
    histname += i;
    fHistLeadingJetE[i] = new TH1F(histname.Data(), histname.Data(), Ebins, Elow, Eup);
    fHistLeadingJetE[i]->GetXaxis()->SetTitle("E [GeV]");
    fHistLeadingJetE[i]->GetYaxis()->SetTitle("counts");
    fOutput->Add(fHistLeadingJetE[i]);
    
    histname = "fHistTracksPtLJ_";
    histname += i;
    fHistTracksPtLJ[i] = new TH1F(histname.Data(), histname.Data(), Ptbins, Ptlow, Ptup);
    fHistTracksPtLJ[i]->GetXaxis()->SetTitle("P_{T} [GeV/c]");
    fHistTracksPtLJ[i]->GetYaxis()->SetTitle("counts");
    fOutput->Add(fHistTracksPtLJ[i]);
    
    histname = "fHistClusELJ_";
    histname += i;
    fHistClusELJ[i] = new TH1F(histname.Data(), histname.Data(), Ebins, Elow, Eup);
    fHistClusELJ[i]->GetXaxis()->SetTitle("E [GeV]");
    fHistClusELJ[i]->GetYaxis()->SetTitle("counts");
    fOutput->Add(fHistClusELJ[i]);

    histname = "fHistTracksPtBkg_";
    histname += i;
    fHistTracksPtBkg[i] = new TH1F(histname.Data(), histname.Data(), Ptbins, Ptlow, Ptup);
    fHistTracksPtBkg[i]->GetXaxis()->SetTitle("P_{T} [GeV/c]");
    fHistTracksPtBkg[i]->GetYaxis()->SetTitle("counts");
    fOutput->Add(fHistTracksPtBkg[i]);
    
    histname = "fHistClusEBkg_";
    histname += i;
    fHistClusEBkg[i] = new TH1F(histname.Data(), histname.Data(), Ebins, Elow, Eup);
    fHistClusEBkg[i]->GetXaxis()->SetTitle("E [GeV]");
    fHistClusEBkg[i]->GetYaxis()->SetTitle("counts");
    fOutput->Add(fHistClusEBkg[i]);
  }

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

void AliAnalysisTaskSAJF::RetrieveEventObjects()
{
  fCaloClusters =  dynamic_cast<TClonesArray*>(InputEvent()->FindListObject(fCaloName));
  if (!fCaloClusters) {
    AliWarning(Form("Could not retrieve clusters!")); 
  }

  fTracks = dynamic_cast<TClonesArray*>(InputEvent()->FindListObject(fTracksName));
  if (!fTracks) {
    AliWarning(Form("Could not retrieve tracks!")); 
  }

  fJets = dynamic_cast<TClonesArray*>(InputEvent()->FindListObject(fJetsName));
  if (!fJets) {
    AliWarning(Form("Could not retrieve jets!")); 
  }

  if (strcmp(fTrgClusName,"")) {
    fTrgClusters =  dynamic_cast<TClonesArray*>(InputEvent()->FindListObject(fTrgClusName));
    if (!fTrgClusters) {
      AliWarning(Form("Could not retrieve trigger clusters!")); 
    }
  }

  fCent = InputEvent()->GetCentrality();
}

AliVTrack* AliAnalysisTaskSAJF::GetTrack(const Int_t i) const
{
  if (fTracks)
    return dynamic_cast<AliVTrack*>(fTracks->At(i));
  else
    return 0;
}

Int_t AliAnalysisTaskSAJF::GetNumberOfTracks() const
{
  if (fTracks)
    return fTracks->GetEntriesFast();
  else
    return 0;
}

AliVCluster* AliAnalysisTaskSAJF::GetCaloCluster(const Int_t i) const
{ 
  if (fCaloClusters)
    return dynamic_cast<AliVCluster*>(fCaloClusters->At(i));
  else
    return 0;
}

Int_t AliAnalysisTaskSAJF::GetNumberOfCaloClusters() const
{ 
  if (fCaloClusters)
    return fCaloClusters->GetEntriesFast();
  else
    return 0;
}

AliEmcalJet* AliAnalysisTaskSAJF::GetJet(const Int_t i) const
{
  if (fJets)
    return dynamic_cast<AliEmcalJet*>(fJets->At(i));
  else
    return 0;
}

Int_t AliAnalysisTaskSAJF::GetNumberOfJets() const
{
  if (fJets)
    return fJets->GetEntriesFast();
  else
    return 0;
}

AliVCluster* AliAnalysisTaskSAJF::GetTrgCluster(const Int_t i) const
{
  if (fTrgClusters)
    return dynamic_cast<AliVCluster*>(fTrgClusters->At(i));
  else
    return 0;
}

Int_t AliAnalysisTaskSAJF::GetNumberOfTrgClusters() const
{
  if (fTrgClusters)
    return fTrgClusters->GetEntriesFast();
  else
    return 0;
}

void AliAnalysisTaskSAJF::FillHistograms()
{
  Float_t cent = fCent->GetCentralityPercentile("V0M");

  fHistCentrality->Fill(cent);

  Int_t centbin=-1;
  if(cent>=0 && cent<10) centbin=0;
  else if(cent>=10 && cent<30) centbin=1;
  else if(cent>=30 && cent<50) centbin=2;
  else if(cent>=50 && cent<=100) centbin=3;

  // Jet loop
  Int_t njets =  GetNumberOfJets();
  //cout << njets << " jets" << endl;
  Float_t maxJetEnergy = 0;
  Int_t maxJetIndex = -1;
  for (Int_t ij = 0; ij < njets; ij++) {

    AliEmcalJet* jet = GetJet(ij);

    if (!jet) {
      printf("ERROR: Could not receive jet %d\n", ij);
      continue;
    }  
    
    if (jet->E() <= 0)
      continue;

    fHistJetsE[centbin]->Fill(jet->E());
    fHistJetsNEF[centbin]->Fill(jet->NEF());
    fHistJetsNE[centbin]->Fill(jet->E() * jet->NEF());

    //cout << "****** jet id = " << ij << " ******" << endl;

    for (Int_t it = 0; it < jet->GetNumberOfTracks(); it++) {
      Int_t trackid = jet->TrackAt(it);
      AliVTrack *track = GetTrack(trackid);
      //cout << "track id = " << trackid << endl;
      if (track)
	fHistJetsZ[centbin]->Fill(track->Pt() / jet->E());
    }
    for (Int_t ic = 0; ic < jet->GetNumberOfClusters(); ic++) {
      Int_t clusterid = jet->ClusterAt(ic);
      //cout << "cluster id = " << clusterid << endl;
      AliVCluster *cluster = GetCaloCluster(clusterid);
      if (cluster)
	fHistJetsZ[centbin]->Fill(cluster->E() / jet->E());
    }

    if (maxJetEnergy < jet->E()) {
      maxJetEnergy = jet->E();
      maxJetIndex = ij;
    }
  } //jet loop 

  
  if (!(maxJetEnergy > 0) || maxJetIndex < 0) 
    return;

  fHistLeadingJetE[centbin]->Fill(maxJetEnergy);
  
  AliEmcalJet* jet = GetJet(maxJetIndex);
  if (!jet)
    return;
  
  jet->SortConstituents();
  
  // Cluster loop
  Int_t clusJetId = 0;
  Int_t nclusters =  GetNumberOfCaloClusters();
  for (Int_t iClusters = 0; iClusters < nclusters; iClusters++) {
    AliVCluster* cluster = GetCaloCluster(iClusters);
    if (!cluster) {
      printf("ERROR: Could not receive cluster %d\n", iClusters);
      continue;
    }  
    
    if (!(cluster->IsEMCAL())) continue;
    
    if (jet->GetNumberOfClusters() > 0) {
      if (jet->ClusterAt(clusJetId) < iClusters) {
	if (clusJetId < jet->GetNumberOfClusters() - 1)
	  clusJetId++;
      }
      
      if (jet->ClusterAt(clusJetId) == iClusters) {
	fHistClusELJ[centbin]->Fill(cluster->E());
      }
      else {
	fHistClusEBkg[centbin]->Fill(cluster->E());
      }
    }
    else {
      fHistClusEBkg[centbin]->Fill(cluster->E());
    }
  } //cluster loop 
  
  
  // Track loop 
  Int_t trackJetId = 0;
  Int_t ntracks = GetNumberOfTracks();
  for(Int_t iTracks = 0; iTracks < ntracks; iTracks++) {
    AliVTrack* track = GetTrack(iTracks);         
    if(!track) {
      AliError(Form("ERROR: Could not retrieve track %d",iTracks)); 
      continue; 
    }
    
    if (!AcceptTrack(track)) continue;
    
    if (jet->GetNumberOfTracks() > 0) {
      if (jet->TrackAt(trackJetId) < iTracks) {
	if (trackJetId < jet->GetNumberOfTracks() - 1)
	  trackJetId++;
      }
      
      if (jet->TrackAt(trackJetId) == iTracks) {
	fHistTracksPtLJ[centbin]->Fill(track->Pt());
      }
      else {
	fHistTracksPtBkg[centbin]->Fill(track->Pt());
      }
    }
    else {
      fHistTracksPtBkg[centbin]->Fill(track->Pt());
    }
  }
}


//________________________________________________________________________
Bool_t AliAnalysisTaskSAJF::AcceptTrack(AliVTrack* /*track*/)
{
  return kTRUE;
}

//________________________________________________________________________
void AliAnalysisTaskSAJF::UserExec(Option_t *) 
{
  // Main loop, called for each event.
  // Add jets to event if not yet there

  RetrieveEventObjects();

  FillHistograms();
    
  // information for this iteration of the UserExec in the container
  PostData(1, fOutput);
}

//________________________________________________________________________
void AliAnalysisTaskSAJF::Terminate(Option_t *) 
{
  // Called once at the end of the analysis.
}
Commit	Line	Data
00514d01	1	// $Id$
25283b37	2
	3	#include <TChain.h>
	4	#include <TClonesArray.h>
	5	#include <TH1F.h>
	6	#include <TH2F.h>
	7	#include <TList.h>
	8	#include <TLorentzVector.h>
	9	#include <TParticle.h>
	10
	11	#include "AliAnalysisManager.h"
	12	#include "AliCentrality.h"
f0a0fd33	13	#include "AliVCluster.h"
25283b37	14	#include "AliESDtrack.h"
25283b37	15	#include "AliEmcalJet.h"
25283b37	16	#include "AliAODTrack.h"
25283b37	17	#include "AliEmcalJet.h"
25283b37	18	#include "AliVEventHandler.h"
35789a2d	19	#include "AliPicoTrack.h"
25283b37	20
00514d01	21	#include "AliAnalysisTaskSAJF.h"
25283b37	22
00514d01	23	ClassImp(AliAnalysisTaskSAJF)
25283b37	24
25283b37	25	//________________________________________________________________________
00514d01	26	AliAnalysisTaskSAJF::AliAnalysisTaskSAJF() :
00514d01	27	AliAnalysisTaskSE("AliAnalysisTaskSAJF"),
25283b37	28	fOutput(0),
	29	fTracksName("Tracks"),
	30	fCaloName("CaloClusters"),
	31	fJetsName("Jets"),
e82e282c	32	fTrgClusName("ClustersL1GAMMAFEE"),
25283b37	33	fTracks(0),
	34	fCaloClusters(0),
	35	fJets(0),
e82e282c	36	fTrgClusters(0),
f0a0fd33	37	fCent(0),
	38	fHistCentrality(0),
	39	Ptbins(400),
25283b37	40	Ptlow(0),
f0a0fd33	41	Ptup(200),
f0a0fd33	42	Ebins(400),
25283b37	43	Elow(0),
f0a0fd33	44	Eup(200)
25283b37	45	{
	46	// Default constructor.
	47	fBranchNames="ESD:AliESDRun.,AliESDHeader.,PrimaryVertex.";
f0a0fd33	48
	49	for (Int_t i = 0; i < 4; i++) {
	50	fHistJetsE[i] = 0;
	51	fHistJetsNE[i] = 0;
	52	fHistJetsNEF[i] = 0;
	53	fHistJetsZ[i] = 0;
	54	fHistLeadingJetE[i] = 0;
	55	fHistTracksPtLJ[i] = 0;
	56	fHistClusELJ[i] = 0;
	57	fHistTracksPtBkg[i] = 0;
	58	fHistClusEBkg[i] = 0;
	59	}
25283b37	60	}
	61
	62	//________________________________________________________________________
00514d01	63	AliAnalysisTaskSAJF::AliAnalysisTaskSAJF(const char *name) :
00514d01	64	AliAnalysisTaskSE(name),
25283b37	65	fOutput(0),
	66	fTracksName("Tracks"),
	67	fCaloName("CaloClusters"),
	68	fJetsName("Jets"),
e82e282c	69	fTrgClusName("ClustersL1GAMMAFEE"),
25283b37	70	fTracks(0),
	71	fCaloClusters(0),
	72	fJets(0),
e82e282c	73	fTrgClusters(0),
f0a0fd33	74	fCent(0),
	75	fHistCentrality(0),
	76	Ptbins(400),
25283b37	77	Ptlow(0),
f0a0fd33	78	Ptup(200),
f0a0fd33	79	Ebins(400),
25283b37	80	Elow(0),
f0a0fd33	81	Eup(200)
25283b37	82	{
	83	// Standard constructor.
	84
25283b37	85	fBranchNames="ESD:AliESDRun.,AliESDHeader.,PrimaryVertex.";
25283b37	86
f0a0fd33	87	for (Int_t i = 0; i < 4; i++) {
	88	fHistJetsE[i] = 0;
	89	fHistJetsNE[i] = 0;
	90	fHistJetsNEF[i] = 0;
	91	fHistJetsZ[i] = 0;
	92	fHistLeadingJetE[i] = 0;
	93	fHistTracksPtLJ[i] = 0;
	94	fHistClusELJ[i] = 0;
	95	fHistTracksPtBkg[i] = 0;
	96	fHistClusEBkg[i] = 0;
	97	}
25283b37	98	}
	99
	100	//________________________________________________________________________
00514d01	101	AliAnalysisTaskSAJF::~AliAnalysisTaskSAJF()
25283b37	102	{
	103	// Destructor
	104	}
	105
	106	//________________________________________________________________________
00514d01	107	void AliAnalysisTaskSAJF::UserCreateOutputObjects()
25283b37	108	{
f0a0fd33	109	// Create histograms
25283b37	110
	111	fOutput = new TList();
	112	fOutput->SetOwner(); // IMPORTANT!
	113
f0a0fd33	114	fHistCentrality = new TH1F("fHistCentrality","Event centrality distribution", Ebins, 0, 100);
	115	fHistCentrality->GetXaxis()->SetTitle("Centrality (%)");
	116	fHistCentrality->GetYaxis()->SetTitle("counts");
	117	fOutput->Add(fHistCentrality);
25283b37	118
f0a0fd33	119	TString histname;
	120
	121	for (Int_t i = 0; i < 4; i++) {
	122	histname = "fHistJetsE_";
	123	histname += i;
	124	fHistJetsE[i] = new TH1F(histname.Data(), histname.Data(), Ebins, Elow, Eup);
	125	fHistJetsE[i]->GetXaxis()->SetTitle("E [GeV]");
	126	fHistJetsE[i]->GetYaxis()->SetTitle("counts");
	127	fOutput->Add(fHistJetsE[i]);
25283b37	128
f0a0fd33	129	histname = "fHistJetsNE_";
	130	histname += i;
	131	fHistJetsNE[i] = new TH1F(histname.Data(), histname.Data(), Ebins, Elow, Eup);
	132	fHistJetsNE[i]->GetXaxis()->SetTitle("E [GeV]");
	133	fHistJetsNE[i]->GetYaxis()->SetTitle("counts");
	134	fOutput->Add(fHistJetsNE[i]);
25283b37	135
f0a0fd33	136	histname = "fHistJetsNEF_";
	137	histname += i;
	138	fHistJetsNEF[i] = new TH1F(histname.Data(), histname.Data(), Ebins, 0, 1.2);
	139	fHistJetsNEF[i]->GetXaxis()->SetTitle("NEF");
	140	fHistJetsNEF[i]->GetYaxis()->SetTitle("counts");
	141	fOutput->Add(fHistJetsNEF[i]);
	142
	143	histname = "fHistJetsZ_";
	144	histname += i;
	145	fHistJetsZ[i] = new TH1F(histname.Data(), histname.Data(), Ebins, 0, 1.2);
	146	fHistJetsZ[i]->GetXaxis()->SetTitle("Z");
	147	fHistJetsZ[i]->GetYaxis()->SetTitle("counts");
	148	fOutput->Add(fHistJetsZ[i]);
	149
	150	histname = "fHistLeadingJetE_";
	151	histname += i;
	152	fHistLeadingJetE[i] = new TH1F(histname.Data(), histname.Data(), Ebins, Elow, Eup);
	153	fHistLeadingJetE[i]->GetXaxis()->SetTitle("E [GeV]");
	154	fHistLeadingJetE[i]->GetYaxis()->SetTitle("counts");
	155	fOutput->Add(fHistLeadingJetE[i]);
25283b37	156
f0a0fd33	157	histname = "fHistTracksPtLJ_";
	158	histname += i;
	159	fHistTracksPtLJ[i] = new TH1F(histname.Data(), histname.Data(), Ptbins, Ptlow, Ptup);
	160	fHistTracksPtLJ[i]->GetXaxis()->SetTitle("P_{T} [GeV/c]");
	161	fHistTracksPtLJ[i]->GetYaxis()->SetTitle("counts");
	162	fOutput->Add(fHistTracksPtLJ[i]);
	163
	164	histname = "fHistClusELJ_";
	165	histname += i;
	166	fHistClusELJ[i] = new TH1F(histname.Data(), histname.Data(), Ebins, Elow, Eup);
	167	fHistClusELJ[i]->GetXaxis()->SetTitle("E [GeV]");
	168	fHistClusELJ[i]->GetYaxis()->SetTitle("counts");
	169	fOutput->Add(fHistClusELJ[i]);
	170
	171	histname = "fHistTracksPtBkg_";
	172	histname += i;
	173	fHistTracksPtBkg[i] = new TH1F(histname.Data(), histname.Data(), Ptbins, Ptlow, Ptup);
	174	fHistTracksPtBkg[i]->GetXaxis()->SetTitle("P_{T} [GeV/c]");
	175	fHistTracksPtBkg[i]->GetYaxis()->SetTitle("counts");
	176	fOutput->Add(fHistTracksPtBkg[i]);
	177
	178	histname = "fHistClusEBkg_";
	179	histname += i;
	180	fHistClusEBkg[i] = new TH1F(histname.Data(), histname.Data(), Ebins, Elow, Eup);
	181	fHistClusEBkg[i]->GetXaxis()->SetTitle("E [GeV]");
	182	fHistClusEBkg[i]->GetYaxis()->SetTitle("counts");
	183	fOutput->Add(fHistClusEBkg[i]);
e82e282c	184	}
e82e282c	185
25283b37	186	PostData(1, fOutput); // Post data for ALL output slots >0 here, to get at least an empty histogram
	187	}
	188
00514d01	189	void AliAnalysisTaskSAJF::RetrieveEventObjects()
25283b37	190	{
25283b37	191	fCaloClusters = dynamic_cast<TClonesArray*>(InputEvent()->FindListObject(fCaloName));
f0a0fd33	192	if (!fCaloClusters) {
	193	AliWarning(Form("Could not retrieve clusters!"));
	194	}
25283b37	195
f0a0fd33	196	fTracks = dynamic_cast<TClonesArray*>(InputEvent()->FindListObject(fTracksName));
25283b37	197	if (!fTracks) {
e82e282c	198	AliWarning(Form("Could not retrieve tracks!"));
25283b37	199	}
f0a0fd33	200
f0a0fd33	201	fJets = dynamic_cast<TClonesArray*>(InputEvent()->FindListObject(fJetsName));
25283b37	202	if (!fJets) {
e82e282c	203	AliWarning(Form("Could not retrieve jets!"));
e82e282c	204	}
f0a0fd33	205
	206	if (strcmp(fTrgClusName,"")) {
	207	fTrgClusters = dynamic_cast<TClonesArray*>(InputEvent()->FindListObject(fTrgClusName));
	208	if (!fTrgClusters) {
	209	AliWarning(Form("Could not retrieve trigger clusters!"));
	210	}
25283b37	211	}
f0a0fd33	212
f0a0fd33	213	fCent = InputEvent()->GetCentrality();
25283b37	214	}
25283b37	215
00514d01	216	AliVTrack* AliAnalysisTaskSAJF::GetTrack(const Int_t i) const
25283b37	217	{
e82e282c	218	if (fTracks)
	219	return dynamic_cast<AliVTrack*>(fTracks->At(i));
	220	else
	221	return 0;
25283b37	222	}
25283b37	223
00514d01	224	Int_t AliAnalysisTaskSAJF::GetNumberOfTracks() const
25283b37	225	{
e82e282c	226	if (fTracks)
	227	return fTracks->GetEntriesFast();
	228	else
	229	return 0;
25283b37	230	}
25283b37	231
00514d01	232	AliVCluster* AliAnalysisTaskSAJF::GetCaloCluster(const Int_t i) const
e82e282c	233	{
	234	if (fCaloClusters)
	235	return dynamic_cast<AliVCluster*>(fCaloClusters->At(i));
	236	else
	237	return 0;
25283b37	238	}
25283b37	239
00514d01	240	Int_t AliAnalysisTaskSAJF::GetNumberOfCaloClusters() const
e82e282c	241	{
	242	if (fCaloClusters)
	243	return fCaloClusters->GetEntriesFast();
	244	else
	245	return 0;
25283b37	246	}
25283b37	247
00514d01	248	AliEmcalJet* AliAnalysisTaskSAJF::GetJet(const Int_t i) const
25283b37	249	{
e82e282c	250	if (fJets)
	251	return dynamic_cast<AliEmcalJet*>(fJets->At(i));
	252	else
	253	return 0;
25283b37	254	}
25283b37	255
00514d01	256	Int_t AliAnalysisTaskSAJF::GetNumberOfJets() const
25283b37	257	{
e82e282c	258	if (fJets)
	259	return fJets->GetEntriesFast();
	260	else
	261	return 0;
	262	}
	263
00514d01	264	AliVCluster* AliAnalysisTaskSAJF::GetTrgCluster(const Int_t i) const
e82e282c	265	{
	266	if (fTrgClusters)
	267	return dynamic_cast<AliVCluster*>(fTrgClusters->At(i));
	268	else
	269	return 0;
	270	}
	271
00514d01	272	Int_t AliAnalysisTaskSAJF::GetNumberOfTrgClusters() const
e82e282c	273	{
	274	if (fTrgClusters)
	275	return fTrgClusters->GetEntriesFast();
	276	else
	277	return 0;
25283b37	278	}
25283b37	279
00514d01	280	void AliAnalysisTaskSAJF::FillHistograms()
25283b37	281	{
f0a0fd33	282	Float_t cent = fCent->GetCentralityPercentile("V0M");
25283b37	283
f0a0fd33	284	fHistCentrality->Fill(cent);
e82e282c	285
f0a0fd33	286	Int_t centbin=-1;
	287	if(cent>=0 && cent<10) centbin=0;
	288	else if(cent>=10 && cent<30) centbin=1;
	289	else if(cent>=30 && cent<50) centbin=2;
	290	else if(cent>=50 && cent<=100) centbin=3;
25283b37	291
	292	// Jet loop
	293	Int_t njets = GetNumberOfJets();
	294	//cout << njets << " jets" << endl;
f0a0fd33	295	Float_t maxJetEnergy = 0;
f0a0fd33	296	Int_t maxJetIndex = -1;
25283b37	297	for (Int_t ij = 0; ij < njets; ij++) {
f0a0fd33	298
25283b37	299	AliEmcalJet* jet = GetJet(ij);
f0a0fd33	300
25283b37	301	if (!jet) {
	302	printf("ERROR: Could not receive jet %d\n", ij);
	303	continue;
	304	}
	305
f0a0fd33	306	if (jet->E() <= 0)
	307	continue;
	308
	309	fHistJetsE[centbin]->Fill(jet->E());
	310	fHistJetsNEF[centbin]->Fill(jet->NEF());
	311	fHistJetsNE[centbin]->Fill(jet->E() * jet->NEF());
	312
	313	//cout << "**** jet id = " << ij << " ****" << endl;
	314
35789a2d	315	for (Int_t it = 0; it < jet->GetNumberOfTracks(); it++) {
	316	Int_t trackid = jet->TrackAt(it);
	317	AliVTrack *track = GetTrack(trackid);
f0a0fd33	318	//cout << "track id = " << trackid << endl;
35789a2d	319	if (track)
f0a0fd33	320	fHistJetsZ[centbin]->Fill(track->Pt() / jet->E());
35789a2d	321	}
	322	for (Int_t ic = 0; ic < jet->GetNumberOfClusters(); ic++) {
	323	Int_t clusterid = jet->ClusterAt(ic);
f0a0fd33	324	//cout << "cluster id = " << clusterid << endl;
35789a2d	325	AliVCluster *cluster = GetCaloCluster(clusterid);
35789a2d	326	if (cluster)
f0a0fd33	327	fHistJetsZ[centbin]->Fill(cluster->E() / jet->E());
	328	}
	329
	330	if (maxJetEnergy < jet->E()) {
	331	maxJetEnergy = jet->E();
	332	maxJetIndex = ij;
35789a2d	333	}
25283b37	334	} //jet loop
e82e282c	335
f0a0fd33	336
	337	if (!(maxJetEnergy > 0) \|\| maxJetIndex < 0)
	338	return;
	339
	340	fHistLeadingJetE[centbin]->Fill(maxJetEnergy);
	341
	342	AliEmcalJet* jet = GetJet(maxJetIndex);
	343	if (!jet)
	344	return;
	345
	346	jet->SortConstituents();
	347
	348	// Cluster loop
	349	Int_t clusJetId = 0;
	350	Int_t nclusters = GetNumberOfCaloClusters();
	351	for (Int_t iClusters = 0; iClusters < nclusters; iClusters++) {
	352	AliVCluster* cluster = GetCaloCluster(iClusters);
e82e282c	353	if (!cluster) {
	354	printf("ERROR: Could not receive cluster %d\n", iClusters);
	355	continue;
	356	}
	357
	358	if (!(cluster->IsEMCAL())) continue;
f0a0fd33	359
	360	if (jet->GetNumberOfClusters() > 0) {
	361	if (jet->ClusterAt(clusJetId) < iClusters) {
	362	if (clusJetId < jet->GetNumberOfClusters() - 1)
	363	clusJetId++;
	364	}
	365
	366	if (jet->ClusterAt(clusJetId) == iClusters) {
	367	fHistClusELJ[centbin]->Fill(cluster->E());
	368	}
	369	else {
	370	fHistClusEBkg[centbin]->Fill(cluster->E());
	371	}
	372	}
	373	else {
	374	fHistClusEBkg[centbin]->Fill(cluster->E());
	375	}
e82e282c	376	} //cluster loop
e82e282c	377
f0a0fd33	378
	379	// Track loop
	380	Int_t trackJetId = 0;
	381	Int_t ntracks = GetNumberOfTracks();
	382	for(Int_t iTracks = 0; iTracks < ntracks; iTracks++) {
	383	AliVTrack* track = GetTrack(iTracks);
	384	if(!track) {
	385	AliError(Form("ERROR: Could not retrieve track %d",iTracks));
	386	continue;
25283b37	387	}
f0a0fd33	388
	389	if (!AcceptTrack(track)) continue;
	390
	391	if (jet->GetNumberOfTracks() > 0) {
	392	if (jet->TrackAt(trackJetId) < iTracks) {
	393	if (trackJetId < jet->GetNumberOfTracks() - 1)
	394	trackJetId++;
	395	}
25283b37	396
f0a0fd33	397	if (jet->TrackAt(trackJetId) == iTracks) {
	398	fHistTracksPtLJ[centbin]->Fill(track->Pt());
	399	}
	400	else {
	401	fHistTracksPtBkg[centbin]->Fill(track->Pt());
	402	}
	403	}
	404	else {
	405	fHistTracksPtBkg[centbin]->Fill(track->Pt());
25283b37	406	}
25283b37	407	}
f0a0fd33	408	}
	409
	410
	411	//________________________________________________________________________
00514d01	412	Bool_t AliAnalysisTaskSAJF::AcceptTrack(AliVTrack* /track/)
f0a0fd33	413	{
f0a0fd33	414	return kTRUE;
25283b37	415	}
	416
	417	//________________________________________________________________________
00514d01	418	void AliAnalysisTaskSAJF::UserExec(Option_t *)
25283b37	419	{
	420	// Main loop, called for each event.
	421	// Add jets to event if not yet there
	422
	423	RetrieveEventObjects();
	424
	425	FillHistograms();
	426
	427	// information for this iteration of the UserExec in the container
	428	PostData(1, fOutput);
	429	}
	430
	431	//________________________________________________________________________
00514d01	432	void AliAnalysisTaskSAJF::Terminate(Option_t *)
25283b37	433	{
	434	// Called once at the end of the analysis.
	435	}