[u/mrichter/AliRoot.git] / PWGGA / EMCALJetTasks / UserTasks / AliAnalysisTaskEmcalJetHMEC.cxx

#include "TChain.h"
#include "TTree.h"
#include "TList.h"
#include "TH1F.h"
#include "TH2F.h"
#include "THnSparse.h"
#include "TCanvas.h"
#include <TClonesArray.h>
#include <TParticle.h>
#include "AliVTrack.h"
#include "TParameter.h"

#include "AliAODEvent.h"
#include "AliAnalysisTask.h"
#include "AliAnalysisManager.h"

#include "AliESDEvent.h"
#include "AliESDInputHandler.h"
#include "AliESDCaloCluster.h"
#include "AliESDVertex.h"
#include "AliCentrality.h"
#include "AliAODJet.h"
#include "AliEmcalJet.h"
#include "AliESDtrackCuts.h"
#include "AliAnalysisTaskEmcalJetHMEC.h"
#include "TVector3.h"

ClassImp(AliAnalysisTaskEmcalJetHMEC)

//________________________________________________________________________
AliAnalysisTaskEmcalJetHMEC::AliAnalysisTaskEmcalJetHMEC() : 
  AliAnalysisTaskSE(),
  fTracksName("tracks"),
  fJetsName("jets"),
  fPhimin(-10), 
  fPhimax(10),
  fEtamin(-0.9), 
  fEtamax(0.9),
  fAreacut(0.0),
  fESD(0), 
  fOutputList(0), 
  fHistTrackPt(0),
  fHistCentrality(0), 
  fHistJetEtaPhi(0), 
  fHistTrackEtaPhi(0), 
  fHistJetHEtaPhi(0) 
{
  // Default Constructor

  for(Int_t icent = 0; icent<6; ++icent){
    fHistJetPt[icent]=0;
    fHistJetPtBias[icent]=0;
    fHistJetPtTT[icent]=0;
    for(Int_t iptjet = 0; iptjet<3; ++iptjet){
      for(Int_t ieta = 0; ieta<3; ++ieta){	
	fHistJetH[icent][iptjet][ieta]=0;
	fHistJetHBias[icent][iptjet][ieta]=0;
	fHistJetHTT[icent][iptjet][ieta]=0;

      }
    }
  }

}
//________________________________________________________________________
AliAnalysisTaskEmcalJetHMEC::AliAnalysisTaskEmcalJetHMEC(const char *name) : 
  AliAnalysisTaskSE(name),
  fTracksName("tracks"),
  fJetsName("jets"),
  fPhimin(-10), 
  fPhimax(10),
  fEtamin(-0.9), 
  fEtamax(0.9),
  fAreacut(0.0),
  fESD(0), 
  fOutputList(0), 
  fHistTrackPt(0),
  fHistCentrality(0), 
  fHistJetEtaPhi(0), 
  fHistTrackEtaPhi(0), 
  fHistJetHEtaPhi(0) 
{
  // Constructor
  for(Int_t icent = 0; icent<6; ++icent){
    fHistJetPt[icent]=0;
    fHistJetPtBias[icent]=0;
    fHistJetPtTT[icent]=0;
    for(Int_t iptjet = 0; iptjet<5; ++iptjet){
      for(Int_t ieta = 0; ieta<3; ++ieta){	
	fHistJetH[icent][iptjet][ieta]=0;
	fHistJetHBias[icent][iptjet][ieta]=0;
	fHistJetHTT[icent][iptjet][ieta]=0;
      }
    }
  }

  DefineInput(0, TChain::Class());
  DefineOutput(1, TList::Class());

}

//________________________________________________________________________
void AliAnalysisTaskEmcalJetHMEC::UserCreateOutputObjects()
{
  // Called once

 
  AliVEventHandler* handler = AliAnalysisManager::GetAnalysisManager()->GetInputEventHandler();
  if (!handler) {
    AliError("Input handler not available!");
    return;
  }

  OpenFile(1);
  fOutputList = new TList();
  fOutputList->SetOwner();

  // Create histograms

  fHistTrackPt = new TH1F("fHistTrackPt", "P_{T} distribution", 1000, 0.0, 100.0);


  fHistCentrality = new TH1F("fHistCentrality","centrality",100,0,100);

  fHistJetEtaPhi = new TH2F("fHistJetEtaPhi","Jet eta-phi",900,-1.8,1.8,640,-3.2,3.2);
  fHistTrackEtaPhi = new TH2F("fHistTrackEtaPhi","Track eta-phi",900,-1.8,1.8,640,-3.2,3.2);
  fHistJetHEtaPhi = new TH2F("fHistJetHEtaPhi","Jet-Hadron deta-dphi",900,-1.8,1.8,640,-1.6,4.8);

  char name[200];

 
  for(Int_t icent = 0; icent<6; ++icent){
    sprintf(name,"fHistJetPt_%i",icent);   
    fHistJetPt[icent] = new TH1F(name,name,200,0,200);
    fOutputList->Add(fHistJetPt[icent]);

    sprintf(name,"fHistJetPtBias_%i",icent);   
    fHistJetPtBias[icent] = new TH1F(name,name,200,0,200);
    fOutputList->Add(fHistJetPtBias[icent]);

    sprintf(name,"fHistJetPtTT_%i",icent);   
    fHistJetPtTT[icent] = new TH1F(name,name,200,0,200);
    fOutputList->Add(fHistJetPtTT[icent]);

    for(Int_t iptjet = 0; iptjet<5; ++iptjet){
      for(Int_t ieta = 0; ieta<3; ++ieta){	
	sprintf(name,"fHistJetH_%i_%i_%i",icent,iptjet,ieta);   
	fHistJetH[icent][iptjet][ieta]=new TH2F(name,name,64,-0.5*TMath::Pi(),1.5*TMath::Pi(),300,0,30);
	fOutputList->Add(fHistJetH[icent][iptjet][ieta]);

	sprintf(name,"fHistJetHBias_%i_%i_%i",icent,iptjet,ieta);   
	fHistJetHBias[icent][iptjet][ieta]=new TH2F(name,name,64,-0.5*TMath::Pi(),1.5*TMath::Pi(),300,0,30);
	fOutputList->Add(fHistJetHBias[icent][iptjet][ieta]);

	sprintf(name,"fHistJetHTT_%i_%i_%i",icent,iptjet,ieta);   
	fHistJetHTT[icent][iptjet][ieta]=new TH2F(name,name,64,-0.5*TMath::Pi(),1.5*TMath::Pi(),300,0,30);
	fOutputList->Add(fHistJetHTT[icent][iptjet][ieta]);

      }
    }
  }


  fOutputList->Add(fHistTrackPt);
  fOutputList->Add(fHistCentrality);
  fOutputList->Add(fHistJetEtaPhi);
  fOutputList->Add(fHistTrackEtaPhi);
  fOutputList->Add(fHistJetHEtaPhi);


  PostData(1, fOutputList);

}

//________________________________________________________________________

Double_t AliAnalysisTaskEmcalJetHMEC:: RelativePhi(Double_t mphi,Double_t vphi) {


  if (vphi < -1*TMath::Pi()) vphi += (2*TMath::Pi());
  else if (vphi > TMath::Pi()) vphi -= (2*TMath::Pi());
  if (mphi < -1*TMath::Pi()) mphi += (2*TMath::Pi());
  else if (mphi > TMath::Pi()) mphi -= (2*TMath::Pi());
  Double_t dphi = vphi-mphi;
  if (dphi < -1*TMath::Pi()) dphi += (2*TMath::Pi());
  else if (dphi > TMath::Pi()) dphi -= (2*TMath::Pi());

  return dphi;//dphi in [-Pi, Pi]                                                                                                    
}


//________________________________________________________________________
Int_t AliAnalysisTaskEmcalJetHMEC::GetCentBin(Double_t cent) const 
{
  // Get centrality bin.

  Int_t centbin = -1;
  if (cent>=0 && cent<10)
    centbin = 0;
  else if (cent>=10 && cent<20)
    centbin = 1;
  else if (cent>=20 && cent<30)
    centbin = 2;
  else if (cent>=30 && cent<40)
    centbin = 3;
  else if (cent>=40 && cent<50)
    centbin = 4;
  else if (cent>=50 && cent<90)
    centbin = 5;
  return centbin;
}


//________________________________________________________________________
Int_t AliAnalysisTaskEmcalJetHMEC::GetEtaBin(Double_t eta) const 
{
  // Get eta bin for histos.

  Int_t etabin = -1;
  if (TMath::Abs(eta)<=0.4)
    etabin = 0;
  else if (TMath::Abs(eta)>0.4 && TMath::Abs(eta)<0.8)
    etabin = 1;
  else if (TMath::Abs(eta)>=0.8)
    etabin = 2;
  return etabin;
}
//________________________________________________________________________
Int_t AliAnalysisTaskEmcalJetHMEC::GetpTjetBin(Double_t pt) const 
{
  // Get jet pt  bin for histos.

  Int_t ptbin = -1;
  if (pt>=10 && pt<15)
    ptbin = 0;
  else if (pt>=15 && pt<20)
    ptbin = 1;
  else if (pt>=20 && pt<25)
    ptbin = 2;
  else if (pt>=25 && pt<30)
    ptbin = 3;
  else if (pt>=30)
    ptbin = 4;


  return ptbin;
}


//________________________________________________________________________
void AliAnalysisTaskEmcalJetHMEC::UserExec(Option_t *) 
{


  // Main loop called for each event
 // esd or aod mode
  Bool_t esdMode = kTRUE;
  if (dynamic_cast<AliAODEvent*>(InputEvent()))
    esdMode = kFALSE;


  if (esdMode) {
    // optimization in case autobranch loading is off
    AliAnalysisManager *am = AliAnalysisManager::GetAnalysisManager();
    if (fTracksName == "Tracks")
      am->LoadBranch("Tracks");
  }
 

  //get centrality
  TList *list = InputEvent()->GetList();
  AliCentrality *centrality = InputEvent()->GetCentrality() ;
  Double_t fcent=-1; 
  if(centrality)
    fcent = centrality->GetCentralityPercentile("V0M");
  else
    fcent=99;//probably pp data

  if (fcent<0) {
    AliError(Form("Centrality negative: %f", fcent));
    return;
  }


  fHistCentrality->Fill(fcent);
  Int_t centbin = GetCentBin(fcent);

  TClonesArray *jets = 0;
  TClonesArray *tracks = 0;

  tracks = dynamic_cast<TClonesArray*>(list->FindObject(fTracksName));
  if (!tracks) {
    AliError(Form("Pointer to tracks %s == 0", fTracksName.Data() ));
    return;
  }
  const Int_t Ntracks=tracks->GetEntries();

  jets= dynamic_cast<TClonesArray*>(list->FindObject(fJetsName));
  const Int_t Njets = jets->GetEntries();
 
  //Leticia's loop to find hardest track

  Int_t iTT=-1;
  Double_t ptmax=-10;

  for (Int_t iTracks = 0; iTracks < Ntracks; iTracks++) 
    {
      AliVTrack* track = static_cast<AliVTrack*>(tracks->At(iTracks));
      if (!track) {
	printf("ERROR: Could not receive track %d\n", iTracks);
	continue;
      }
      
      if(TMath::Abs(track->Eta())>0.9) continue;
      if(track->Pt()<0.15)continue;
      //iCount++;
      if(track->Pt()>ptmax){
	ptmax=track->Pt();
	iTT=iTracks;
      }
    }

  
  Int_t ijethi=-1;

  Double_t highestjetpt=0.0;

  Int_t passedTTcut=0;

  for (Int_t ijet = 0; ijet < Njets; ijet++)
    {
 
      AliEmcalJet *jet = static_cast<AliEmcalJet*>(jets->At(ijet));
      
      if (!jet)
	continue;
      
      //pt,eta,phi,centrality
      float jetphi = jet->Phi();
      if (jetphi>TMath::Pi())
	jetphi = jetphi-2*TMath::Pi();
      
      if ((jet->Phi()<fPhimin)||(jet->Phi()>fPhimax))
	continue;
      if ((jet->Eta()<fEtamin)||(jet->Eta()>fEtamax))
	continue;
      //fHistAreavsRawPt[centbin]->Fill(jet->Pt(),jet->Area());
      if (jet->Area()<fAreacut)
	continue;
      //prevents 0 area jets from sneaking by when area cut == 0
      if (jet->Area()==0)
	continue;

      Double_t jetPt = jet->Pt();

      if(highestjetpt<jetPt){
	ijethi=ijet;
	highestjetpt=jetPt;
      }

    }


  //Only look at the highest pT jet in the event

  if(ijethi>-1){
    AliEmcalJet *jet = static_cast<AliEmcalJet*>(jets->At(ijethi)); 
    
    float jetphi = jet->Phi();
    Double_t jetPt = jet->Pt();

      fHistJetPt[centbin]->Fill(jet->Pt());

      if ((jet->MaxTrackPt()>6) || (jet->MaxClusterPt()>6))
	fHistJetPtBias[centbin]->Fill(jet->Pt());


      fHistJetEtaPhi->Fill(jet->Eta(),jetphi);

      //fHistDeltaPtvsArea->Fill(jetPt,jet->Area());

      if(iTT>0){
	AliVTrack* TT = static_cast<AliVTrack*>(tracks->At(iTT));
	if(TMath::Abs(jetphi-TT->Phi()-TMath::Pi())<0.6) passedTTcut=1;
	else passedTTcut=0;
      }

      if(passedTTcut)
	fHistJetPtTT[centbin]->Fill(jet->Pt());
   

  if (highestjetpt>10) {
   
    for (Int_t iTracks = 0; iTracks < Ntracks; iTracks++) 
      {
	AliVTrack* track = static_cast<AliVTrack*>(tracks->At(iTracks));
	if (!track) {
	  printf("ERROR: Could not receive track %d\n", iTracks);
	  continue;
	}
	
	  if(TMath::Abs(track->Eta())>0.9) continue;

	  fHistTrackPt->Fill(track->Pt());
	  	  
	  if (track->Pt()<0.15)
	    continue;
	  
	  Double_t trackphi = track->Phi();
	  if (trackphi > TMath::Pi())
	    trackphi = trackphi-2*TMath::Pi();

	  Double_t tracketa=track->Eta();
	  Double_t jeteta=jet->Eta();

	  Double_t deta=tracketa-jeteta;
	  Int_t ieta=GetEtaBin(deta);

	  //Jet pt, track pt, dPhi,deta,fcent
	  Double_t dphijh = RelativePhi(jetphi,trackphi);
	  if (dphijh < -0.5*TMath::Pi())
	    dphijh= dphijh+ 2*TMath::Pi();

	  Int_t iptjet=-1;
	  iptjet=GetpTjetBin(jetPt);

	  fHistJetH[centbin][iptjet][ieta]->Fill(dphijh,track->Pt());
	  fHistJetHEtaPhi->Fill(deta,dphijh);
	  fHistTrackEtaPhi->Fill(tracketa,trackphi);
	  if ((jet->MaxTrackPt()>6) || (jet->MaxClusterPt()>6))
	    fHistJetHBias[centbin][iptjet][ieta]->Fill(dphijh,track->Pt());
	  
	  if(passedTTcut)
	  fHistJetHTT[centbin][iptjet][ieta]->Fill(dphijh,track->Pt());


	} //track loop
  }//jet pt cut
  }

  PostData(1, fOutputList);
}      

//________________________________________________________________________
void AliAnalysisTaskEmcalJetHMEC::Terminate(Option_t *) 
{
  //just terminate

}
Commit	Line	Data
1101e468	1	#include "TChain.h"
	2	#include "TTree.h"
	3	#include "TList.h"
	4	#include "TH1F.h"
	5	#include "TH2F.h"
	6	#include "THnSparse.h"
	7	#include "TCanvas.h"
	8	#include <TClonesArray.h>
	9	#include <TParticle.h>
	10	#include "AliVTrack.h"
	11	#include "TParameter.h"
	12
	13	#include "AliAODEvent.h"
	14	#include "AliAnalysisTask.h"
	15	#include "AliAnalysisManager.h"
	16
	17	#include "AliESDEvent.h"
	18	#include "AliESDInputHandler.h"
	19	#include "AliESDCaloCluster.h"
	20	#include "AliESDVertex.h"
	21	#include "AliCentrality.h"
	22	#include "AliAODJet.h"
	23	#include "AliEmcalJet.h"
	24	#include "AliESDtrackCuts.h"
	25	#include "AliAnalysisTaskEmcalJetHMEC.h"
	26	#include "TVector3.h"
	27
	28	ClassImp(AliAnalysisTaskEmcalJetHMEC)
	29
	30	//________________________________________________________________________
	31	AliAnalysisTaskEmcalJetHMEC::AliAnalysisTaskEmcalJetHMEC() :
	32	AliAnalysisTaskSE(),
	33	fTracksName("tracks"),
	34	fJetsName("jets"),
	35	fPhimin(-10),
	36	fPhimax(10),
	37	fEtamin(-0.9),
	38	fEtamax(0.9),
	39	fAreacut(0.0),
	40	fESD(0),
	41	fOutputList(0),
	42	fHistTrackPt(0),
	43	fHistCentrality(0),
	44	fHistJetEtaPhi(0),
	45	fHistTrackEtaPhi(0),
	46	fHistJetHEtaPhi(0)
	47	{
	48	// Default Constructor
	49
	50	for(Int_t icent = 0; icent<6; ++icent){
	51	fHistJetPt[icent]=0;
81941f29	52	fHistJetPtBias[icent]=0;
81941f29	53	fHistJetPtTT[icent]=0;
1101e468	54	for(Int_t iptjet = 0; iptjet<3; ++iptjet){
	55	for(Int_t ieta = 0; ieta<3; ++ieta){
	56	fHistJetH[icent][iptjet][ieta]=0;
	57	fHistJetHBias[icent][iptjet][ieta]=0;
81941f29	58	fHistJetHTT[icent][iptjet][ieta]=0;
81941f29	59
1101e468	60	}
	61	}
	62	}
	63
	64	}
	65	//________________________________________________________________________
	66	AliAnalysisTaskEmcalJetHMEC::AliAnalysisTaskEmcalJetHMEC(const char *name) :
	67	AliAnalysisTaskSE(name),
	68	fTracksName("tracks"),
	69	fJetsName("jets"),
	70	fPhimin(-10),
	71	fPhimax(10),
	72	fEtamin(-0.9),
	73	fEtamax(0.9),
	74	fAreacut(0.0),
	75	fESD(0),
	76	fOutputList(0),
	77	fHistTrackPt(0),
	78	fHistCentrality(0),
	79	fHistJetEtaPhi(0),
	80	fHistTrackEtaPhi(0),
	81	fHistJetHEtaPhi(0)
	82	{
	83	// Constructor
	84	for(Int_t icent = 0; icent<6; ++icent){
	85	fHistJetPt[icent]=0;
	86	fHistJetPtBias[icent]=0;
81941f29	87	fHistJetPtTT[icent]=0;
81941f29	88	for(Int_t iptjet = 0; iptjet<5; ++iptjet){
1101e468	89	for(Int_t ieta = 0; ieta<3; ++ieta){
	90	fHistJetH[icent][iptjet][ieta]=0;
	91	fHistJetHBias[icent][iptjet][ieta]=0;
81941f29	92	fHistJetHTT[icent][iptjet][ieta]=0;
1101e468	93	}
	94	}
	95	}
	96
	97	DefineInput(0, TChain::Class());
	98	DefineOutput(1, TList::Class());
	99
	100	}
	101
	102	//________________________________________________________________________
	103	void AliAnalysisTaskEmcalJetHMEC::UserCreateOutputObjects()
	104	{
	105	// Called once
	106
	107
	108	AliVEventHandler* handler = AliAnalysisManager::GetAnalysisManager()->GetInputEventHandler();
	109	if (!handler) {
	110	AliError("Input handler not available!");
	111	return;
	112	}
	113
	114	OpenFile(1);
	115	fOutputList = new TList();
	116	fOutputList->SetOwner();
	117
	118	// Create histograms
	119
	120	fHistTrackPt = new TH1F("fHistTrackPt", "P_{T} distribution", 1000, 0.0, 100.0);
	121
	122
	123
	124	fHistCentrality = new TH1F("fHistCentrality","centrality",100,0,100);
	125
	126	fHistJetEtaPhi = new TH2F("fHistJetEtaPhi","Jet eta-phi",900,-1.8,1.8,640,-3.2,3.2);
	127	fHistTrackEtaPhi = new TH2F("fHistTrackEtaPhi","Track eta-phi",900,-1.8,1.8,640,-3.2,3.2);
	128	fHistJetHEtaPhi = new TH2F("fHistJetHEtaPhi","Jet-Hadron deta-dphi",900,-1.8,1.8,640,-1.6,4.8);
	129
	130	char name[200];
	131
	132
	133	for(Int_t icent = 0; icent<6; ++icent){
	134	sprintf(name,"fHistJetPt_%i",icent);
	135	fHistJetPt[icent] = new TH1F(name,name,200,0,200);
	136	fOutputList->Add(fHistJetPt[icent]);
	137
	138	sprintf(name,"fHistJetPtBias_%i",icent);
	139	fHistJetPtBias[icent] = new TH1F(name,name,200,0,200);
	140	fOutputList->Add(fHistJetPtBias[icent]);
	141
81941f29	142	sprintf(name,"fHistJetPtTT_%i",icent);
	143	fHistJetPtTT[icent] = new TH1F(name,name,200,0,200);
	144	fOutputList->Add(fHistJetPtTT[icent]);
	145
	146	for(Int_t iptjet = 0; iptjet<5; ++iptjet){
1101e468	147	for(Int_t ieta = 0; ieta<3; ++ieta){
	148	sprintf(name,"fHistJetH_%i_%i_%i",icent,iptjet,ieta);
	149	fHistJetH[icent][iptjet][ieta]=new TH2F(name,name,64,-0.5TMath::Pi(),1.5TMath::Pi(),300,0,30);
	150	fOutputList->Add(fHistJetH[icent][iptjet][ieta]);
	151
	152	sprintf(name,"fHistJetHBias_%i_%i_%i",icent,iptjet,ieta);
	153	fHistJetHBias[icent][iptjet][ieta]=new TH2F(name,name,64,-0.5TMath::Pi(),1.5TMath::Pi(),300,0,30);
	154	fOutputList->Add(fHistJetHBias[icent][iptjet][ieta]);
	155
81941f29	156	sprintf(name,"fHistJetHTT_%i_%i_%i",icent,iptjet,ieta);
	157	fHistJetHTT[icent][iptjet][ieta]=new TH2F(name,name,64,-0.5TMath::Pi(),1.5TMath::Pi(),300,0,30);
	158	fOutputList->Add(fHistJetHTT[icent][iptjet][ieta]);
	159
1101e468	160	}
	161	}
	162	}
	163
	164
	165
	166	fOutputList->Add(fHistTrackPt);
	167	fOutputList->Add(fHistCentrality);
	168	fOutputList->Add(fHistJetEtaPhi);
	169	fOutputList->Add(fHistTrackEtaPhi);
	170	fOutputList->Add(fHistJetHEtaPhi);
	171
	172
	173	PostData(1, fOutputList);
	174
	175	}
	176
	177	//________________________________________________________________________
	178
	179	Double_t AliAnalysisTaskEmcalJetHMEC:: RelativePhi(Double_t mphi,Double_t vphi) {
	180
	181
	182	if (vphi < -1TMath::Pi()) vphi += (2TMath::Pi());
	183	else if (vphi > TMath::Pi()) vphi -= (2*TMath::Pi());
	184	if (mphi < -1TMath::Pi()) mphi += (2TMath::Pi());
	185	else if (mphi > TMath::Pi()) mphi -= (2*TMath::Pi());
	186	Double_t dphi = vphi-mphi;
	187	if (dphi < -1TMath::Pi()) dphi += (2TMath::Pi());
	188	else if (dphi > TMath::Pi()) dphi -= (2*TMath::Pi());
	189
	190	return dphi;//dphi in [-Pi, Pi]
	191	}
	192
	193
	194	//________________________________________________________________________
	195	Int_t AliAnalysisTaskEmcalJetHMEC::GetCentBin(Double_t cent) const
	196	{
	197	// Get centrality bin.
	198
	199	Int_t centbin = -1;
	200	if (cent>=0 && cent<10)
	201	centbin = 0;
	202	else if (cent>=10 && cent<20)
	203	centbin = 1;
	204	else if (cent>=20 && cent<30)
	205	centbin = 2;
	206	else if (cent>=30 && cent<40)
	207	centbin = 3;
	208	else if (cent>=40 && cent<50)
	209	centbin = 4;
	210	else if (cent>=50 && cent<90)
	211	centbin = 5;
	212	return centbin;
	213	}
	214
	215
	216	//________________________________________________________________________
	217	Int_t AliAnalysisTaskEmcalJetHMEC::GetEtaBin(Double_t eta) const
	218	{
	219	// Get eta bin for histos.
	220
	221	Int_t etabin = -1;
	222	if (TMath::Abs(eta)<=0.4)
	223	etabin = 0;
224	else if (TMath::Abs(eta)>0.4 && TMath::Abs(eta)<0.8)
225	etabin = 1;
226	else if (TMath::Abs(eta)>=0.8)
227	etabin = 2;
228	return etabin;
229	}
230	//________________________________________________________________________
81941f29	231	Int_t AliAnalysisTaskEmcalJetHMEC::GetpTjetBin(Double_t pt) const
1101e468	232	{
81941f29	233	// Get jet pt bin for histos.
	234
	235	Int_t ptbin = -1;
	236	if (pt>=10 && pt<15)
	237	ptbin = 0;
	238	else if (pt>=15 && pt<20)
	239	ptbin = 1;
	240	else if (pt>=20 && pt<25)
	241	ptbin = 2;
	242	else if (pt>=25 && pt<30)
	243	ptbin = 3;
	244	else if (pt>=30)
	245	ptbin = 4;
	246
	247
	248	return ptbin;
1101e468	249	}
	250
	251
	252	//________________________________________________________________________
	253	void AliAnalysisTaskEmcalJetHMEC::UserExec(Option_t *)
	254	{
	255
	256
	257	// Main loop called for each event
	258	// esd or aod mode
	259	Bool_t esdMode = kTRUE;
	260	if (dynamic_cast<AliAODEvent*>(InputEvent()))
	261	esdMode = kFALSE;
	262
	263
	264	if (esdMode) {
	265	// optimization in case autobranch loading is off
	266	AliAnalysisManager *am = AliAnalysisManager::GetAnalysisManager();
	267	if (fTracksName == "Tracks")
	268	am->LoadBranch("Tracks");
	269	}
	270
	271
	272	//get centrality
	273	TList *list = InputEvent()->GetList();
	274	AliCentrality *centrality = InputEvent()->GetCentrality() ;
	275	Double_t fcent=-1;
	276	if(centrality)
	277	fcent = centrality->GetCentralityPercentile("V0M");
	278	else
	279	fcent=99;//probably pp data
	280
	281	if (fcent<0) {
	282	AliError(Form("Centrality negative: %f", fcent));
	283	return;
	284	}
	285
	286
	287	fHistCentrality->Fill(fcent);
	288	Int_t centbin = GetCentBin(fcent);
	289
	290	TClonesArray *jets = 0;
	291	TClonesArray *tracks = 0;
	292
	293	tracks = dynamic_cast<TClonesArray*>(list->FindObject(fTracksName));
	294	if (!tracks) {
	295	AliError(Form("Pointer to tracks %s == 0", fTracksName.Data() ));
	296	return;
	297	}
	298	const Int_t Ntracks=tracks->GetEntries();
	299
	300	jets= dynamic_cast<TClonesArray*>(list->FindObject(fJetsName));
	301	const Int_t Njets = jets->GetEntries();
	302
81941f29	303	//Leticia's loop to find hardest track
	304
	305	Int_t iTT=-1;
	306	Double_t ptmax=-10;
	307
	308	for (Int_t iTracks = 0; iTracks < Ntracks; iTracks++)
	309	{
	310	AliVTrack* track = static_cast<AliVTrack*>(tracks->At(iTracks));
	311	if (!track) {
	312	printf("ERROR: Could not receive track %d\n", iTracks);
	313	continue;
	314	}
	315
	316	if(TMath::Abs(track->Eta())>0.9) continue;
	317	if(track->Pt()<0.15)continue;
	318	//iCount++;
	319	if(track->Pt()>ptmax){
	320	ptmax=track->Pt();
	321	iTT=iTracks;
	322	}
	323	}
1101e468	324
1101e468	325
81941f29	326	Int_t ijethi=-1;
	327
	328	Double_t highestjetpt=0.0;
	329
	330	Int_t passedTTcut=0;
1101e468	331
	332	for (Int_t ijet = 0; ijet < Njets; ijet++)
	333	{
81941f29	334
1101e468	335	AliEmcalJet jet = static_cast<AliEmcalJet>(jets->At(ijet));
	336
	337	if (!jet)
	338	continue;
	339
	340	//pt,eta,phi,centrality
	341	float jetphi = jet->Phi();
	342	if (jetphi>TMath::Pi())
	343	jetphi = jetphi-2*TMath::Pi();
	344
	345	if ((jet->Phi()<fPhimin)\|\|(jet->Phi()>fPhimax))
	346	continue;
	347	if ((jet->Eta()<fEtamin)\|\|(jet->Eta()>fEtamax))
	348	continue;
	349	//fHistAreavsRawPt[centbin]->Fill(jet->Pt(),jet->Area());
	350	if (jet->Area()<fAreacut)
	351	continue;
	352	//prevents 0 area jets from sneaking by when area cut == 0
	353	if (jet->Area()==0)
	354	continue;
	355
1101e468	356	Double_t jetPt = jet->Pt();
1101e468	357
59d5ae99	358	if(highestjetpt<jetPt){
	359	ijethi=ijet;
	360	highestjetpt=jetPt;
	361	}
	362
	363	}
	364
	365
	366
	367	//Only look at the highest pT jet in the event
	368
	369	if(ijethi>-1){
	370	AliEmcalJet jet = static_cast<AliEmcalJet>(jets->At(ijethi));
	371
	372	float jetphi = jet->Phi();
	373	Double_t jetPt = jet->Pt();
1101e468	374
	375	fHistJetPt[centbin]->Fill(jet->Pt());
	376
81941f29	377	if ((jet->MaxTrackPt()>6) \|\| (jet->MaxClusterPt()>6))
	378	fHistJetPtBias[centbin]->Fill(jet->Pt());
	379
	380
1101e468	381	fHistJetEtaPhi->Fill(jet->Eta(),jetphi);
	382
	383	//fHistDeltaPtvsArea->Fill(jetPt,jet->Area());
	384
81941f29	385	if(iTT>0){
	386	AliVTrack* TT = static_cast<AliVTrack*>(tracks->At(iTT));
	387	if(TMath::Abs(jetphi-TT->Phi()-TMath::Pi())<0.6) passedTTcut=1;
	388	else passedTTcut=0;
	389	}
1101e468	390
81941f29	391	if(passedTTcut)
81941f29	392	fHistJetPtTT[centbin]->Fill(jet->Pt());
59d5ae99	393
81941f29	394
81941f29	395	if (highestjetpt>10) {
81941f29	396
81941f29	397	for (Int_t iTracks = 0; iTracks < Ntracks; iTracks++)
	398	{
	399	AliVTrack* track = static_cast<AliVTrack*>(tracks->At(iTracks));
	400	if (!track) {
	401	printf("ERROR: Could not receive track %d\n", iTracks);
	402	continue;
	403	}
	404
1101e468	405	if(TMath::Abs(track->Eta())>0.9) continue;
	406
	407	fHistTrackPt->Fill(track->Pt());
	408
81941f29	409	if (track->Pt()<0.15)
1101e468	410	continue;
	411
	412	Double_t trackphi = track->Phi();
	413	if (trackphi > TMath::Pi())
	414	trackphi = trackphi-2*TMath::Pi();
	415
	416	Double_t tracketa=track->Eta();
	417	Double_t jeteta=jet->Eta();
	418
	419	Double_t deta=tracketa-jeteta;
	420	Int_t ieta=GetEtaBin(deta);
	421
	422	//Jet pt, track pt, dPhi,deta,fcent
	423	Double_t dphijh = RelativePhi(jetphi,trackphi);
	424	if (dphijh < -0.5*TMath::Pi())
	425	dphijh= dphijh+ 2*TMath::Pi();
	426
	427	Int_t iptjet=-1;
	428	iptjet=GetpTjetBin(jetPt);
	429
	430	fHistJetH[centbin][iptjet][ieta]->Fill(dphijh,track->Pt());
	431	fHistJetHEtaPhi->Fill(deta,dphijh);
	432	fHistTrackEtaPhi->Fill(tracketa,trackphi);
81941f29	433	if ((jet->MaxTrackPt()>6) \|\| (jet->MaxClusterPt()>6))
1101e468	434	fHistJetHBias[centbin][iptjet][ieta]->Fill(dphijh,track->Pt());
81941f29	435
81941f29	436	if(passedTTcut)
59d5ae99	437	fHistJetHTT[centbin][iptjet][ieta]->Fill(dphijh,track->Pt());
1101e468	438
1101e468	439
1101e468	440	} //track loop
59d5ae99	441	}//jet pt cut
	442	}
	443
1101e468	444	PostData(1, fOutputList);
	445	}
	446
	447	//________________________________________________________________________
	448	void AliAnalysisTaskEmcalJetHMEC::Terminate(Option_t *)
	449	{
	450	//just terminate
	451
	452	}
	453
	454