[u/mrichter/AliRoot.git] / PWGGA / EMCALJetTasks / AliAnalysisTaskEmcalJetSpectra.cxx

#include <TChain.h>
#include <TTree.h>
#include <TList.h>
#include <TH1F.h>
#include <TH2F.h>
#include <TCanvas.h>
#include <TClonesArray.h>
#include <TParticle.h>
#include <TLorentzVector.h>
#include <TVector3.h>
#include <TParameter.h>

#include "AliAnalysisTask.h"
#include "AliAnalysisManager.h"
#include "AliEmcalJet.h"
#include "AliVCluster.h"

#include "AliESDEvent.h"
#include "AliAODEvent.h"
#include "AliESDInputHandler.h"
#include "AliCentrality.h"

#include "AliAnalysisTaskEmcalJetSpectra.h"

ClassImp(AliAnalysisTaskEmcalJetSpectra)

//________________________________________________________________________
AliAnalysisTaskEmcalJetSpectra::AliAnalysisTaskEmcalJetSpectra() 
  : AliAnalysisTaskSE(), fESD(0), fOutputList(0), fHistCentrality(0), fHistJetArea(0), fHistJetMaxPt(0), fHistJetZ(0), fHistJetNEF(0), fHistJetPtvsCent(0),  fHistJetPtM3vsCent(0),  fHistLeadingJetPtvsCent(0), fHistLeadingJetPtM3vsCent(0),  fHistJetAreavsCent(0),  fHistJetMaxPtvsCent(0),  fHistJetZvsCent(0),  fHistJetNEFvsCent(0), fHistNjetvsCent(0),  fHistJetPtvsNtrack(0),  fHistJetAreavsNtrack(0),  fHistJetMaxPtvsNtrack(0),  fHistJetZvsNtrack(0),  fHistJetNEFvsNtrack(0), fHistNjetvsNtrack(0), fHistDeltaRho12vsCent(0), fHistDeltaRho13vsCent(0), fHistDeltaRho23vsCent(0), fHistDeltaJetPt12vsCent(0), fHistDeltaJetPt13vsCent(0), fHistDeltaJetPt23vsCent(0), fHistRho1vsCent(0), fHistRho2vsCent(0), fHistRho3vsCent(0),
    fTracksName("tracks"),
    fJetsName("jets"),
    fClustersName("clusters"),
    fRhos1Name("fArrRhos"),
    fRhos2Name("fArrRhos"),
    fRhos3Name("fArrRhos"),
    phimin(-10), phimax(10),
    etamin(-0.9), etamax(0.9),
    areacut(0.0)
 {
  // Constructor

   for (Int_t i = 0;i<6;++i){
     fHistRawJetPt[i]       = 0;
     fHistAreavsRawPt[i]    = 0;
     for (Int_t j = 0;j<4;++j){
       fHistNEFvsPt[i][j]   = 0;
       fHistZvsPt[i][j]     = 0;
       fHistZchvsPt[i][j]   = 0;
       fHistZemvsPt[i][j]   = 0;
       fHistAreavsPt[i][j]  = 0;
       fHistJetPt[i][j]     = 0;
       fHistNconsvsPt[i][j] = 0;
       fHistJetPt5[i][j]    = 0;
       fHistJetPt6[i][j]    = 0;
       fHistJetPt7[i][j]    = 0;
       fHistJetPt8[i][j]    = 0;
     }
   }
   
   DefineInput(0, TChain::Class());
   DefineOutput(1, TList::Class());
}
//________________________________________________________________________
AliAnalysisTaskEmcalJetSpectra::AliAnalysisTaskEmcalJetSpectra(const char *name) 
  : AliAnalysisTaskSE(name), fESD(0), fOutputList(0), fHistCentrality(0), fHistJetArea(0), fHistJetMaxPt(0), fHistJetZ(0), fHistJetNEF(0), fHistJetPtvsCent(0),  fHistJetPtM3vsCent(0),  fHistLeadingJetPtvsCent(0), fHistLeadingJetPtM3vsCent(0),  fHistJetAreavsCent(0),  fHistJetMaxPtvsCent(0),  fHistJetZvsCent(0),  fHistJetNEFvsCent(0), fHistNjetvsCent(0),  fHistJetPtvsNtrack(0),  fHistJetAreavsNtrack(0),  fHistJetMaxPtvsNtrack(0),  fHistJetZvsNtrack(0),  fHistJetNEFvsNtrack(0), fHistNjetvsNtrack(0), fHistDeltaRho12vsCent(0), fHistDeltaRho13vsCent(0), fHistDeltaRho23vsCent(0), fHistDeltaJetPt12vsCent(0), fHistDeltaJetPt13vsCent(0), fHistDeltaJetPt23vsCent(0), fHistRho1vsCent(0), fHistRho2vsCent(0), fHistRho3vsCent(0),
    fTracksName("tracks"),
    fJetsName("jets"),
    fClustersName("clusters"),
    fRhos1Name("fArrRhos"),
    fRhos2Name("fArrRhos"),
    fRhos3Name("fArrRhos"),
    phimin(-10), phimax(10),
    etamin(-0.9), etamax(0.9),
    areacut(0.0)
{
  // Constructor
   for (Int_t i = 0;i<6;++i){
     fHistRawJetPt[i]       = 0;
     fHistAreavsRawPt[i]    = 0;
     for (Int_t j = 0;j<4;++j){
       fHistNEFvsPt[i][j]   = 0;
       fHistZvsPt[i][j]     = 0;
       fHistZchvsPt[i][j]   = 0;
       fHistZemvsPt[i][j]   = 0;
       fHistAreavsPt[i][j]  = 0;
       fHistJetPt[i][j]     = 0;
       fHistNconsvsPt[i][j] = 0;
       fHistJetPt5[i][j]    = 0;
       fHistJetPt6[i][j]    = 0;
       fHistJetPt7[i][j]    = 0;
       fHistJetPt8[i][j]    = 0;
     }
   }
  
  DefineInput(0, TChain::Class());
  DefineOutput(1, TList::Class());
}

//________________________________________________________________________
void AliAnalysisTaskEmcalJetSpectra::UserCreateOutputObjects()
{

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

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

  fHistCentrality             = new TH1F("Centrality",            "Centrality",            101, -1,  100);
  fHistJetPtvsCent            = new TH2F("JetPtvsCent",           "JetPtvsCent",           101, -1,  100,   200,  0,   500);
  fHistJetPtM3vsCent          = new TH2F("JetPtM3vsCent",         "JetPtM3vsCent",         101, -1,  100,   200, -250, 250);
  fHistLeadingJetPtvsCent     = new TH2F("LeadingJetPtvsCent",    "LeadingJetPtvsCent",    101, -1,  100,   200,  0,   500);
  fHistLeadingJetPtM3vsCent   = new TH2F("LeadingJetPtM3vsCent",  "LeadingJetPtM3vsCent",  101, -1,  100,   200, -250, 250);
  fHistJetAreavsCent          = new TH2F("JetAreavsCent",         "JetAreavsCent",         101, -1,  100,   100,  0,   1.0);
  fHistJetMaxPtvsCent         = new TH2F("JetMaxPtvsCent",        "JetMaxPtvsCent",        101, -1,  100,   200,  0,   100);
  fHistJetZvsCent             = new TH2F("JetZvsCent",            "JetZvsCent",            101, -1,  100,   100,  0,   1.0);
  fHistJetNEFvsCent           = new TH2F("JetNEFvsCent",          "JetNEFvsCent",          101, -1,  100,   100,  0,   1.0);
  fHistNjetvsCent             = new TH2F("NjetvsCent",            "NjetvsCent",            101, -1,  100,   100,  0,   100);

  fHistJetPtvsNtrack          = new TH2F("JetPtvsNtrack",         "JetPtvsNtrack",         500,  0,  2500,  200,  0,   500);
  fHistJetAreavsNtrack        = new TH2F("JetAreavsNtrack",       "JetAreavsNtrack",       500,  0,  2500,  100,  0,   1.0);
  fHistJetMaxPtvsNtrack       = new TH2F("JetMaxPtvsNtrack",      "JetMaxPtvsNtrack",      500,  0,  2500,  200,  0,   100);
  fHistJetZvsNtrack           = new TH2F("JetZvsNtrack",          "JetZvsNtrack",          500,  0,  2500,  100,  0,   1.0);
  fHistJetNEFvsNtrack         = new TH2F("JetNEFvsNtrack",        "JetNEFvsNtrack",        500,  0,  2500,  100,  0,   1.0);
  fHistNjetvsNtrack           = new TH2F("NjetvsNtrack",          "rNjetvsNtrack",         500,  0,  2500,  100,  0,   100);
 
  fHistJetArea               = new TH1F("JetArea",                "Jet Area",              100, 0.0,    1.0);
  fHistJetMaxPt              = new TH1F("JetMaxPt",               "Jet pt max track",      100, 0.0,  100.0);
  fHistJetZ                  = new TH1F("JetZ",                   "Jet Z",                 100, 0.0,    1.0);
  fHistJetNEF                = new TH1F("JetNEF",                 "Jet NEF",               100, 0.0,    1.0);
  fHistDeltaRho12vsCent      = new TH2F("DeltaRho12vsCent",       "DeltaRho12vsCent",      100, 0.0,  100.0, 500,-250,250);
  fHistDeltaRho13vsCent      = new TH2F("DeltaRho13vsCent",       "DeltaRho13vsCent",      100, 0.0,  100.0, 500,-250,250);
  fHistDeltaRho23vsCent      = new TH2F("DeltaRho23vsCent",       "DeltaRho23vsCent",      100, 0.0,  100.0, 500,-250,250);
  fHistDeltaJetPt12vsCent    = new TH2F("DeltaJetPt12vsCent",     "DeltaJetPt12vsCent",    100, 0.0,  100.0, 500,-250,250);
  fHistDeltaJetPt13vsCent    = new TH2F("DeltaJetPt13vsCent",     "DeltaJetPt13vsCent",    100, 0.0,  100.0, 500,-250,250);
  fHistDeltaJetPt23vsCent    = new TH2F("DeltaJetPt23vsCent",     "DeltaJetPt23vsCent",    100, 0.0,  100.0, 500,-250,250);
  fHistRho1vsCent            = new TH2F("Rho1vsCent",             "Rho1vsCent",            100, 0.0, 100.0, 500, 0, 500);
  fHistRho2vsCent            = new TH2F("Rho2vsCent",             "Rho2vsCent",            100, 0.0, 100.0, 500, 0, 500);
  fHistRho3vsCent            = new TH2F("Rho3vsCent",             "Rho3vsCent",            100, 0.0, 100.0, 500, 0, 500);


  for (Int_t i = 0;i<6;++i){
      TString name00(Form("fHistRawJetPt_%i",i));
      fHistRawJetPt[i] = new TH1F(name00,name00,250,0,500);
      fOutputList->Add(fHistRawJetPt[i]);
      TString name01(Form("fHistAreavsRawPt_%i",i));
      fHistAreavsRawPt[i] = new TH2F(name01,name01,250,0,500,100,0,1);
      fOutputList->Add(fHistAreavsRawPt[i]);
      for (Int_t j = 0;j<4;j++){
	TString name0(Form("fHistNEFvsPt_%i_%i",i,j));
	fHistNEFvsPt[i][j] = new TH2F(name0,name0,250,-250,250,200,0,2);
	fOutputList->Add(fHistNEFvsPt[i][j]);
	TString name1(Form("fHistZvsPt_%i_%i",i,j));
	fHistZvsPt[i][j] = new TH2F(name1,name1,250,-250,250,200,0,2);
	fOutputList->Add(fHistZvsPt[i][j]);
	TString name2(Form("fHistZchvsPt_%i_%i",i,j));
	fHistZchvsPt[i][j] = new TH2F(name2,name2,250,-250,250,200,0,2);
	fOutputList->Add(fHistZchvsPt[i][j]);
	TString name3(Form("fHistZemvsPt_%i_%i",i,j));
	fHistZemvsPt[i][j] = new TH2F(name3,name3,250,-250,250,200,0,2);
	fOutputList->Add(fHistZemvsPt[i][j]);
	TString name4(Form("fHistAreavsPt_%i_%i",i,j));
	fHistAreavsPt[i][j] = new TH2F(name4,name4,250,-250,250,100,0,1);
	fOutputList->Add(fHistAreavsPt[i][j]);
	TString name5(Form("fHistJetPt_%i_%i",i,j));
	fHistJetPt[i][j] = new TH1F(name5,name5,250,-250,250);
	fOutputList->Add(fHistJetPt[i][j]);
	TString name6(Form("fHistNconsvsPt_%i_%i",i,j));
	fHistNconsvsPt[i][j] = new TH2F(name6,name6,250,-250,250,500,0,500);
	fOutputList->Add(fHistNconsvsPt[i][j]);
	TString name7(Form("fHistJetPt5_%i_%i",i,j));
	fHistJetPt5[i][j] = new TH1F(name7,name7,250,-250,250);
	fOutputList->Add(fHistJetPt5[i][j]);
	TString name8(Form("fHistJetPt6_%i_%i",i,j));
	fHistJetPt6[i][j] = new TH1F(name8,name8,250,-250,250);
	fOutputList->Add(fHistJetPt6[i][j]);
	TString name9(Form("fHistJetPt7_%i_%i",i,j));
	fHistJetPt7[i][j] = new TH1F(name9,name9,250,-250,250);
	fOutputList->Add(fHistJetPt7[i][j]);
	TString name10(Form("fHistJetPt8_%i_%i",i,j));
	fHistJetPt8[i][j] = new TH1F(name10,name10,250,-250,250);
	fOutputList->Add(fHistJetPt8[i][j]);
      }
  }
  fOutputList->Add(fHistCentrality);
  fOutputList->Add(fHistJetArea);
  fOutputList->Add(fHistJetMaxPt);
  fOutputList->Add(fHistJetZ);
  fOutputList->Add(fHistJetNEF);
  fOutputList->Add(fHistJetPtvsCent);
  fOutputList->Add(fHistJetPtM3vsCent);
  fOutputList->Add(fHistLeadingJetPtvsCent);
  fOutputList->Add(fHistLeadingJetPtM3vsCent);
  fOutputList->Add(fHistJetAreavsCent);
  fOutputList->Add(fHistJetMaxPtvsCent);
  fOutputList->Add(fHistJetZvsCent); 
  fOutputList->Add(fHistJetNEFvsCent);
  fOutputList->Add(fHistNjetvsCent);
  
  fOutputList->Add(fHistJetPtvsNtrack);
  fOutputList->Add(fHistJetAreavsNtrack);
  fOutputList->Add(fHistJetMaxPtvsNtrack);
  fOutputList->Add(fHistJetZvsNtrack);
  fOutputList->Add(fHistJetNEFvsNtrack);
  fOutputList->Add(fHistNjetvsNtrack);
  fOutputList->Add(fHistDeltaRho12vsCent);
  fOutputList->Add(fHistDeltaRho13vsCent);
  fOutputList->Add(fHistDeltaRho23vsCent);
  fOutputList->Add(fHistDeltaJetPt12vsCent);
  fOutputList->Add(fHistDeltaJetPt13vsCent);
  fOutputList->Add(fHistDeltaJetPt23vsCent);
  fOutputList->Add(fHistRho1vsCent);
  fOutputList->Add(fHistRho2vsCent);
  fOutputList->Add(fHistRho3vsCent);
  
   PostData(1, fOutputList);
  
}
//________________________________________________________________________

Int_t AliAnalysisTaskEmcalJetSpectra::GetCentBin(Double_t cent) const {
  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;
}
//________________________________________________________________________
void AliAnalysisTaskEmcalJetSpectra::UserExec(Option_t *) 
{
  // Main loop, called for each event.

  // esd or aod mode
  Bool_t esdMode = kTRUE;
  if (dynamic_cast<AliAODEvent*>(InputEvent()))
    esdMode = kFALSE;

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

  // get centrality 
  Double_t fCent = -1; 
  TList *l = InputEvent()->GetList();
  AliCentrality *centrality = InputEvent()->GetCentrality() ;
  if (centrality)
    fCent = centrality->GetCentralityPercentile("V0M");
  else
    fCent=99; // probably pp data
  if (fCent<0) {
    AliError(Form("Centrality negative: %f", fCent));
    return;
  }

  Int_t centbin = GetCentBin(fCent);

  //don't want to analyze events 90-100
  if (centbin < 0)
    return;

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

  tracks = dynamic_cast<TClonesArray*>(l->FindObject(fTracksName));
  if (!tracks) {
    AliError(Form("Pointer to tracks %s == 0", fTracksName.Data() ));
    return;
  }
    
  jets = dynamic_cast<TClonesArray*>(l->FindObject(fJetsName));
  if (!jets) {
    AliError(Form("Pointer to tracks %s == 0", fTracksName.Data() ));
    return;
  }

  Double_t rho1 = -1;
  TParameter<Double_t> *Rho1Param = dynamic_cast<TParameter<Double_t>*>(InputEvent()->FindListObject(fRhos1Name));
  if (Rho1Param)
    rho1 = Rho1Param->GetVal();

  Double_t rho2 = -1;
  TParameter<Double_t> *Rho2Param = dynamic_cast<TParameter<Double_t>*>(InputEvent()->FindListObject(fRhos2Name));
  if (Rho2Param)
    rho2 = Rho2Param->GetVal();

  Double_t rho3 = -1;
  TParameter<Double_t> *Rho3Param = dynamic_cast<TParameter<Double_t>*>(InputEvent()->FindListObject(fRhos3Name));
  if (Rho3Param)
    rho3 = Rho3Param->GetVal();

  if (rho1>0)
    fHistRho1vsCent->Fill(fCent,rho1);
  if (rho2>0)
    fHistRho2vsCent->Fill(fCent,rho2);
  if (rho3>0)
    fHistRho3vsCent->Fill(fCent,rho3);

  if (( rho1>0 ) && ( rho2>0 )) 
    fHistDeltaRho12vsCent->Fill(fCent,rho1-rho2);
  if (( rho1>0 ) && ( rho3>0 )) 
    fHistDeltaRho13vsCent->Fill(fCent,rho1-rho3);
  if (( rho2>0 ) && ( rho3>0 )) 
    fHistDeltaRho23vsCent->Fill(fCent,rho2-rho3);
  fHistCentrality->Fill(fCent);

  const Int_t Ntracks = tracks->GetEntries();
  Double_t lJetPt = -500;
  Double_t lJetPtun = -500;
    
  const Int_t Njets = jets->GetEntries();
  Int_t NjetAcc = 0;

 
  for (Int_t iJets = 0; iJets < Njets; ++iJets) {
    
    AliEmcalJet *jet = static_cast<AliEmcalJet*>(jets->At(iJets));
    if (!jet)
      continue; 
     if ((jet->Phi()<phimin)||(jet->Phi()>phimax))
      continue;
    if ((jet->Eta()<etamin)||(jet->Eta()>etamax))
      continue;
    fHistAreavsRawPt[centbin]->Fill(jet->Pt(),jet->Area());
    if (jet->Area()<areacut)
      continue;
    //prevents 0 area jets from sneaking by when area cut == 0
    if (jet->Area()==0)
      continue;
    fHistRawJetPt[centbin]->Fill(jet->Pt());
    
    NjetAcc++;
    fHistJetArea->Fill(jet->Area());
    fHistJetMaxPt->Fill(jet->MaxTrackPt());
    float Z; 
    float jetPt1 = -500;
    float jetPt2 = -500;
    float jetPt3 = -500;
    if (rho1 > 0)
      jetPt1 = jet->Pt()-jet->Area()*rho1;
    if (rho2 > 0)
      jetPt2 = jet->Pt()-jet->Area()*rho2;
    if (rho3 > 0)
      jetPt3 = jet->Pt()-jet->Area()*rho3;
  if (( rho1>0 ) && ( rho2>0 )) 
    fHistDeltaJetPt12vsCent->Fill(fCent,jetPt1-jetPt2);
  if (( rho1>0 ) && ( rho3>0 )) 
    fHistDeltaJetPt13vsCent->Fill(fCent,jetPt1-jetPt3);
  if (( rho2>0 ) && ( rho3>0 )) 
    fHistDeltaJetPt23vsCent->Fill(fCent,jetPt2-jetPt3);

    if (lJetPt < jetPt2){
      lJetPt= jetPt2;
      lJetPtun = jet->Pt();
    }

    if (jet->MaxTrackPt()>jet->MaxClusterPt()){
      Z = jet->MaxTrackPt();
      fHistJetMaxPtvsCent->Fill(fCent,jet->MaxTrackPt());
      fHistJetMaxPtvsNtrack->Fill(Ntracks,jet->MaxTrackPt());
    }
    else{
      Z = jet->MaxClusterPt();
      fHistJetMaxPtvsCent->Fill(fCent,jet->MaxClusterPt());
      fHistJetMaxPtvsNtrack->Fill(Ntracks,jet->MaxClusterPt());
    }

    fHistJetZ->Fill(jet->MaxTrackPt()/jetPt2);
    
    fHistJetNEF->Fill(jet->NEF());
    fHistJetPtvsCent->Fill(fCent,jet->Pt());
    fHistJetPtM3vsCent->Fill(fCent,jetPt2);
    fHistJetAreavsCent->Fill(fCent,jet->Area());
    fHistJetZvsCent->Fill(fCent,Z);
    fHistJetZvsNtrack->Fill(Ntracks,Z);
  
    fHistNEFvsPt[centbin][0]->Fill(jet->Pt(),jet->NEF());
    fHistZvsPt[centbin][0]->Fill(jet->Pt(),Z/jet->Pt());
    fHistZchvsPt[centbin][0]->Fill(jet->Pt(),jet->MaxTrackPt()/jet->Pt());
    fHistZemvsPt[centbin][0]->Fill(jet->Pt(),jet->MaxClusterPt()/jet->Pt());
    fHistAreavsPt[centbin][0]->Fill(jet->Pt(),jet->Area());
    fHistJetPt[centbin][0]->Fill(jet->Pt());
    fHistNconsvsPt[centbin][0]->Fill(jet->Pt(),jet->N());
    if ((jet->MaxTrackPt()>5) || (jet->MaxClusterPt()>5))
      fHistJetPt5[centbin][0]->Fill(jet->Pt());
    if ((jet->MaxTrackPt()>6) || (jet->MaxClusterPt()>6))
      fHistJetPt6[centbin][0]->Fill(jet->Pt());
    if ((jet->MaxTrackPt()>7) || (jet->MaxClusterPt()>7))
      fHistJetPt7[centbin][0]->Fill(jet->Pt());
    if ((jet->MaxTrackPt()>8) || (jet->MaxClusterPt()>8))
      fHistJetPt8[centbin][0]->Fill(jet->Pt());

    fHistNEFvsPt[centbin][1]->Fill(jetPt1,jet->NEF());
    fHistZvsPt[centbin][1]->Fill(jetPt1,Z/jetPt1);
    fHistZchvsPt[centbin][1]->Fill(jetPt1,jet->MaxTrackPt()/jetPt1);
    fHistZemvsPt[centbin][1]->Fill(jetPt1,jet->MaxClusterPt()/jetPt1);
    fHistAreavsPt[centbin][1]->Fill(jetPt1,jet->Area());
    fHistJetPt[centbin][1]->Fill(jetPt1);
    fHistNconsvsPt[centbin][1]->Fill(jetPt1,jet->N());
    if ((jet->MaxTrackPt()>5) || (jet->MaxClusterPt()>5))
      fHistJetPt5[centbin][1]->Fill(jetPt1);
    if ((jet->MaxTrackPt()>6) || (jet->MaxClusterPt()>6))
      fHistJetPt6[centbin][1]->Fill(jetPt1);
    if ((jet->MaxTrackPt()>7) || (jet->MaxClusterPt()>7))
      fHistJetPt7[centbin][1]->Fill(jetPt1);
    if ((jet->MaxTrackPt()>8) || (jet->MaxClusterPt()>8))
      fHistJetPt8[centbin][1]->Fill(jetPt1);

    fHistNEFvsPt[centbin][2]->Fill(jetPt2,jet->NEF());
    fHistZvsPt[centbin][2]->Fill(jetPt2,Z/jetPt2);
    fHistZchvsPt[centbin][2]->Fill(jetPt2,jet->MaxTrackPt()/jetPt2);
    fHistZemvsPt[centbin][2]->Fill(jetPt2,jet->MaxClusterPt()/jetPt2);
    fHistAreavsPt[centbin][2]->Fill(jetPt2,jet->Area());
    fHistJetPt[centbin][2]->Fill(jetPt2);
    fHistNconsvsPt[centbin][2]->Fill(jetPt2,jet->N());
    if ((jet->MaxTrackPt()>5) || (jet->MaxClusterPt()>5))
      fHistJetPt5[centbin][2]->Fill(jetPt2);
    if ((jet->MaxTrackPt()>6) || (jet->MaxClusterPt()>6))
      fHistJetPt6[centbin][2]->Fill(jetPt2);
    if ((jet->MaxTrackPt()>7) || (jet->MaxClusterPt()>7))
      fHistJetPt7[centbin][2]->Fill(jetPt2);
    if ((jet->MaxTrackPt()>8) || (jet->MaxClusterPt()>8))
      fHistJetPt8[centbin][2]->Fill(jetPt2);

    fHistNEFvsPt[centbin][3]->Fill(jetPt3,jet->NEF());
    fHistZvsPt[centbin][3]->Fill(jetPt3,Z/jetPt3);
    fHistZchvsPt[centbin][3]->Fill(jetPt3,jet->MaxTrackPt()/jetPt3);
    fHistZemvsPt[centbin][3]->Fill(jetPt3,jet->MaxClusterPt()/jetPt3);
    fHistAreavsPt[centbin][3]->Fill(jetPt3,jet->Area());
    fHistJetPt[centbin][3]->Fill(jetPt3);
    fHistNconsvsPt[centbin][3]->Fill(jetPt3,jet->N());
    if ((jet->MaxTrackPt()>5) || (jet->MaxClusterPt()>5))
      fHistJetPt5[centbin][3]->Fill(jetPt3);
    if ((jet->MaxTrackPt()>6) || (jet->MaxClusterPt()>6))
      fHistJetPt6[centbin][3]->Fill(jetPt3);
    if ((jet->MaxTrackPt()>7) || (jet->MaxClusterPt()>7))
      fHistJetPt7[centbin][3]->Fill(jetPt3);
    if ((jet->MaxTrackPt()>8) || (jet->MaxClusterPt()>8))
      fHistJetPt8[centbin][3]->Fill(jetPt3);
    
    fHistJetNEFvsCent->Fill(fCent,jet->NEF());
    fHistJetPtvsNtrack->Fill(Ntracks,jetPt2);
    fHistJetAreavsNtrack->Fill(Ntracks,jet->Area());
    fHistJetNEFvsNtrack->Fill(Ntracks,jet->NEF());
  }

  fHistLeadingJetPtvsCent->Fill(fCent,lJetPtun);
  fHistLeadingJetPtM3vsCent->Fill(fCent,lJetPt);
  fHistNjetvsCent->Fill(fCent,NjetAcc);
  fHistNjetvsNtrack->Fill(Ntracks,NjetAcc);
  
  PostData(1, fOutputList);
}      

//________________________________________________________________________
void AliAnalysisTaskEmcalJetSpectra::Terminate(Option_t *) 
{

}
Commit	Line	Data
74ec9ac9	1	#include <TChain.h>
	2	#include <TTree.h>
	3	#include <TList.h>
	4	#include <TH1F.h>
	5	#include <TH2F.h>
	6	#include <TCanvas.h>
020052e4	7	#include <TClonesArray.h>
	8	#include <TParticle.h>
	9	#include <TLorentzVector.h>
6172360f	10	#include <TVector3.h>
74ec9ac9	11	#include <TParameter.h>
020052e4	12
	13	#include "AliAnalysisTask.h"
	14	#include "AliAnalysisManager.h"
	15	#include "AliEmcalJet.h"
	16	#include "AliVCluster.h"
	17
	18	#include "AliESDEvent.h"
	19	#include "AliAODEvent.h"
	20	#include "AliESDInputHandler.h"
	21	#include "AliCentrality.h"
	22
	23	#include "AliAnalysisTaskEmcalJetSpectra.h"
	24
	25	ClassImp(AliAnalysisTaskEmcalJetSpectra)
	26
	27	//________________________________________________________________________
	28	AliAnalysisTaskEmcalJetSpectra::AliAnalysisTaskEmcalJetSpectra()
6172360f	29	: AliAnalysisTaskSE(), fESD(0), fOutputList(0), fHistCentrality(0), fHistJetArea(0), fHistJetMaxPt(0), fHistJetZ(0), fHistJetNEF(0), fHistJetPtvsCent(0), fHistJetPtM3vsCent(0), fHistLeadingJetPtvsCent(0), fHistLeadingJetPtM3vsCent(0), fHistJetAreavsCent(0), fHistJetMaxPtvsCent(0), fHistJetZvsCent(0), fHistJetNEFvsCent(0), fHistNjetvsCent(0), fHistJetPtvsNtrack(0), fHistJetAreavsNtrack(0), fHistJetMaxPtvsNtrack(0), fHistJetZvsNtrack(0), fHistJetNEFvsNtrack(0), fHistNjetvsNtrack(0), fHistDeltaRho12vsCent(0), fHistDeltaRho13vsCent(0), fHistDeltaRho23vsCent(0), fHistDeltaJetPt12vsCent(0), fHistDeltaJetPt13vsCent(0), fHistDeltaJetPt23vsCent(0), fHistRho1vsCent(0), fHistRho2vsCent(0), fHistRho3vsCent(0),
020052e4	30	fTracksName("tracks"),
	31	fJetsName("jets"),
	32	fClustersName("clusters"),
6172360f	33	fRhos1Name("fArrRhos"),
	34	fRhos2Name("fArrRhos"),
	35	fRhos3Name("fArrRhos"),
74ec9ac9	36	phimin(-10), phimax(10),
	37	etamin(-0.9), etamax(0.9),
	38	areacut(0.0)
020052e4	39	{
	40	// Constructor
	41
	42	for (Int_t i = 0;i<6;++i){
	43	fHistRawJetPt[i] = 0;
	44	fHistAreavsRawPt[i] = 0;
74ec9ac9	45	for (Int_t j = 0;j<4;++j){
020052e4	46	fHistNEFvsPt[i][j] = 0;
	47	fHistZvsPt[i][j] = 0;
	48	fHistZchvsPt[i][j] = 0;
	49	fHistZemvsPt[i][j] = 0;
	50	fHistAreavsPt[i][j] = 0;
	51	fHistJetPt[i][j] = 0;
	52	fHistNconsvsPt[i][j] = 0;
6172360f	53	fHistJetPt5[i][j] = 0;
	54	fHistJetPt6[i][j] = 0;
	55	fHistJetPt7[i][j] = 0;
	56	fHistJetPt8[i][j] = 0;
020052e4	57	}
	58	}
	59
	60	DefineInput(0, TChain::Class());
	61	DefineOutput(1, TList::Class());
	62	}
	63	//________________________________________________________________________
	64	AliAnalysisTaskEmcalJetSpectra::AliAnalysisTaskEmcalJetSpectra(const char *name)
6172360f	65	: AliAnalysisTaskSE(name), fESD(0), fOutputList(0), fHistCentrality(0), fHistJetArea(0), fHistJetMaxPt(0), fHistJetZ(0), fHistJetNEF(0), fHistJetPtvsCent(0), fHistJetPtM3vsCent(0), fHistLeadingJetPtvsCent(0), fHistLeadingJetPtM3vsCent(0), fHistJetAreavsCent(0), fHistJetMaxPtvsCent(0), fHistJetZvsCent(0), fHistJetNEFvsCent(0), fHistNjetvsCent(0), fHistJetPtvsNtrack(0), fHistJetAreavsNtrack(0), fHistJetMaxPtvsNtrack(0), fHistJetZvsNtrack(0), fHistJetNEFvsNtrack(0), fHistNjetvsNtrack(0), fHistDeltaRho12vsCent(0), fHistDeltaRho13vsCent(0), fHistDeltaRho23vsCent(0), fHistDeltaJetPt12vsCent(0), fHistDeltaJetPt13vsCent(0), fHistDeltaJetPt23vsCent(0), fHistRho1vsCent(0), fHistRho2vsCent(0), fHistRho3vsCent(0),
020052e4	66	fTracksName("tracks"),
	67	fJetsName("jets"),
	68	fClustersName("clusters"),
6172360f	69	fRhos1Name("fArrRhos"),
	70	fRhos2Name("fArrRhos"),
	71	fRhos3Name("fArrRhos"),
74ec9ac9	72	phimin(-10), phimax(10),
	73	etamin(-0.9), etamax(0.9),
	74	areacut(0.0)
020052e4	75	{
020052e4	76	// Constructor
6172360f	77	for (Int_t i = 0;i<6;++i){
	78	fHistRawJetPt[i] = 0;
	79	fHistAreavsRawPt[i] = 0;
	80	for (Int_t j = 0;j<4;++j){
	81	fHistNEFvsPt[i][j] = 0;
020052e4	82	fHistZvsPt[i][j] = 0;
	83	fHistZchvsPt[i][j] = 0;
	84	fHistZemvsPt[i][j] = 0;
	85	fHistAreavsPt[i][j] = 0;
	86	fHistJetPt[i][j] = 0;
	87	fHistNconsvsPt[i][j] = 0;
6172360f	88	fHistJetPt5[i][j] = 0;
	89	fHistJetPt6[i][j] = 0;
	90	fHistJetPt7[i][j] = 0;
	91	fHistJetPt8[i][j] = 0;
020052e4	92	}
	93	}
	94
	95	DefineInput(0, TChain::Class());
	96	DefineOutput(1, TList::Class());
	97	}
	98
	99	//________________________________________________________________________
	100	void AliAnalysisTaskEmcalJetSpectra::UserCreateOutputObjects()
	101	{
	102
	103	AliVEventHandler* handler = AliAnalysisManager::GetAnalysisManager()->GetInputEventHandler();
	104	if (!handler) {
	105	AliError("Input handler not available!");
	106	return;
	107	}
	108
020052e4	109	OpenFile(1);
	110	fOutputList = new TList();
	111	fOutputList->SetOwner();
	112
	113	fHistCentrality = new TH1F("Centrality", "Centrality", 101, -1, 100);
	114	fHistJetPtvsCent = new TH2F("JetPtvsCent", "JetPtvsCent", 101, -1, 100, 200, 0, 500);
	115	fHistJetPtM3vsCent = new TH2F("JetPtM3vsCent", "JetPtM3vsCent", 101, -1, 100, 200, -250, 250);
	116	fHistLeadingJetPtvsCent = new TH2F("LeadingJetPtvsCent", "LeadingJetPtvsCent", 101, -1, 100, 200, 0, 500);
	117	fHistLeadingJetPtM3vsCent = new TH2F("LeadingJetPtM3vsCent", "LeadingJetPtM3vsCent", 101, -1, 100, 200, -250, 250);
	118	fHistJetAreavsCent = new TH2F("JetAreavsCent", "JetAreavsCent", 101, -1, 100, 100, 0, 1.0);
	119	fHistJetMaxPtvsCent = new TH2F("JetMaxPtvsCent", "JetMaxPtvsCent", 101, -1, 100, 200, 0, 100);
	120	fHistJetZvsCent = new TH2F("JetZvsCent", "JetZvsCent", 101, -1, 100, 100, 0, 1.0);
	121	fHistJetNEFvsCent = new TH2F("JetNEFvsCent", "JetNEFvsCent", 101, -1, 100, 100, 0, 1.0);
	122	fHistNjetvsCent = new TH2F("NjetvsCent", "NjetvsCent", 101, -1, 100, 100, 0, 100);
	123
	124	fHistJetPtvsNtrack = new TH2F("JetPtvsNtrack", "JetPtvsNtrack", 500, 0, 2500, 200, 0, 500);
	125	fHistJetAreavsNtrack = new TH2F("JetAreavsNtrack", "JetAreavsNtrack", 500, 0, 2500, 100, 0, 1.0);
	126	fHistJetMaxPtvsNtrack = new TH2F("JetMaxPtvsNtrack", "JetMaxPtvsNtrack", 500, 0, 2500, 200, 0, 100);
	127	fHistJetZvsNtrack = new TH2F("JetZvsNtrack", "JetZvsNtrack", 500, 0, 2500, 100, 0, 1.0);
	128	fHistJetNEFvsNtrack = new TH2F("JetNEFvsNtrack", "JetNEFvsNtrack", 500, 0, 2500, 100, 0, 1.0);
	129	fHistNjetvsNtrack = new TH2F("NjetvsNtrack", "rNjetvsNtrack", 500, 0, 2500, 100, 0, 100);
	130
	131	fHistJetArea = new TH1F("JetArea", "Jet Area", 100, 0.0, 1.0);
	132	fHistJetMaxPt = new TH1F("JetMaxPt", "Jet pt max track", 100, 0.0, 100.0);
	133	fHistJetZ = new TH1F("JetZ", "Jet Z", 100, 0.0, 1.0);
	134	fHistJetNEF = new TH1F("JetNEF", "Jet NEF", 100, 0.0, 1.0);
	135	fHistDeltaRho12vsCent = new TH2F("DeltaRho12vsCent", "DeltaRho12vsCent", 100, 0.0, 100.0, 500,-250,250);
	136	fHistDeltaRho13vsCent = new TH2F("DeltaRho13vsCent", "DeltaRho13vsCent", 100, 0.0, 100.0, 500,-250,250);
	137	fHistDeltaRho23vsCent = new TH2F("DeltaRho23vsCent", "DeltaRho23vsCent", 100, 0.0, 100.0, 500,-250,250);
6172360f	138	fHistDeltaJetPt12vsCent = new TH2F("DeltaJetPt12vsCent", "DeltaJetPt12vsCent", 100, 0.0, 100.0, 500,-250,250);
	139	fHistDeltaJetPt13vsCent = new TH2F("DeltaJetPt13vsCent", "DeltaJetPt13vsCent", 100, 0.0, 100.0, 500,-250,250);
	140	fHistDeltaJetPt23vsCent = new TH2F("DeltaJetPt23vsCent", "DeltaJetPt23vsCent", 100, 0.0, 100.0, 500,-250,250);
020052e4	141	fHistRho1vsCent = new TH2F("Rho1vsCent", "Rho1vsCent", 100, 0.0, 100.0, 500, 0, 500);
	142	fHistRho2vsCent = new TH2F("Rho2vsCent", "Rho2vsCent", 100, 0.0, 100.0, 500, 0, 500);
	143	fHistRho3vsCent = new TH2F("Rho3vsCent", "Rho3vsCent", 100, 0.0, 100.0, 500, 0, 500);
	144
	145
	146	for (Int_t i = 0;i<6;++i){
74ec9ac9	147	TString name00(Form("fHistRawJetPt_%i",i));
74ec9ac9	148	fHistRawJetPt[i] = new TH1F(name00,name00,250,0,500);
020052e4	149	fOutputList->Add(fHistRawJetPt[i]);
74ec9ac9	150	TString name01(Form("fHistAreavsRawPt_%i",i));
74ec9ac9	151	fHistAreavsRawPt[i] = new TH2F(name01,name01,250,0,500,100,0,1);
020052e4	152	fOutputList->Add(fHistAreavsRawPt[i]);
74ec9ac9	153	for (Int_t j = 0;j<4;j++){
020052e4	154	TString name0(Form("fHistNEFvsPt_%i_%i",i,j));
74ec9ac9	155	fHistNEFvsPt[i][j] = new TH2F(name0,name0,250,-250,250,200,0,2);
020052e4	156	fOutputList->Add(fHistNEFvsPt[i][j]);
020052e4	157	TString name1(Form("fHistZvsPt_%i_%i",i,j));
74ec9ac9	158	fHistZvsPt[i][j] = new TH2F(name1,name1,250,-250,250,200,0,2);
020052e4	159	fOutputList->Add(fHistZvsPt[i][j]);
020052e4	160	TString name2(Form("fHistZchvsPt_%i_%i",i,j));
74ec9ac9	161	fHistZchvsPt[i][j] = new TH2F(name2,name2,250,-250,250,200,0,2);
020052e4	162	fOutputList->Add(fHistZchvsPt[i][j]);
020052e4	163	TString name3(Form("fHistZemvsPt_%i_%i",i,j));
74ec9ac9	164	fHistZemvsPt[i][j] = new TH2F(name3,name3,250,-250,250,200,0,2);
020052e4	165	fOutputList->Add(fHistZemvsPt[i][j]);
	166	TString name4(Form("fHistAreavsPt_%i_%i",i,j));
	167	fHistAreavsPt[i][j] = new TH2F(name4,name4,250,-250,250,100,0,1);
	168	fOutputList->Add(fHistAreavsPt[i][j]);
	169	TString name5(Form("fHistJetPt_%i_%i",i,j));
	170	fHistJetPt[i][j] = new TH1F(name5,name5,250,-250,250);
	171	fOutputList->Add(fHistJetPt[i][j]);
	172	TString name6(Form("fHistNconsvsPt_%i_%i",i,j));
74ec9ac9	173	fHistNconsvsPt[i][j] = new TH2F(name6,name6,250,-250,250,500,0,500);
020052e4	174	fOutputList->Add(fHistNconsvsPt[i][j]);
6172360f	175	TString name7(Form("fHistJetPt5_%i_%i",i,j));
	176	fHistJetPt5[i][j] = new TH1F(name7,name7,250,-250,250);
	177	fOutputList->Add(fHistJetPt5[i][j]);
	178	TString name8(Form("fHistJetPt6_%i_%i",i,j));
	179	fHistJetPt6[i][j] = new TH1F(name8,name8,250,-250,250);
	180	fOutputList->Add(fHistJetPt6[i][j]);
	181	TString name9(Form("fHistJetPt7_%i_%i",i,j));
	182	fHistJetPt7[i][j] = new TH1F(name9,name9,250,-250,250);
	183	fOutputList->Add(fHistJetPt7[i][j]);
	184	TString name10(Form("fHistJetPt8_%i_%i",i,j));
	185	fHistJetPt8[i][j] = new TH1F(name10,name10,250,-250,250);
	186	fOutputList->Add(fHistJetPt8[i][j]);
020052e4	187	}
	188	}
	189	fOutputList->Add(fHistCentrality);
	190	fOutputList->Add(fHistJetArea);
	191	fOutputList->Add(fHistJetMaxPt);
	192	fOutputList->Add(fHistJetZ);
	193	fOutputList->Add(fHistJetNEF);
	194	fOutputList->Add(fHistJetPtvsCent);
	195	fOutputList->Add(fHistJetPtM3vsCent);
	196	fOutputList->Add(fHistLeadingJetPtvsCent);
	197	fOutputList->Add(fHistLeadingJetPtM3vsCent);
	198	fOutputList->Add(fHistJetAreavsCent);
	199	fOutputList->Add(fHistJetMaxPtvsCent);
	200	fOutputList->Add(fHistJetZvsCent);
	201	fOutputList->Add(fHistJetNEFvsCent);
	202	fOutputList->Add(fHistNjetvsCent);
	203
	204	fOutputList->Add(fHistJetPtvsNtrack);
	205	fOutputList->Add(fHistJetAreavsNtrack);
	206	fOutputList->Add(fHistJetMaxPtvsNtrack);
	207	fOutputList->Add(fHistJetZvsNtrack);
	208	fOutputList->Add(fHistJetNEFvsNtrack);
	209	fOutputList->Add(fHistNjetvsNtrack);
	210	fOutputList->Add(fHistDeltaRho12vsCent);
	211	fOutputList->Add(fHistDeltaRho13vsCent);
	212	fOutputList->Add(fHistDeltaRho23vsCent);
6172360f	213	fOutputList->Add(fHistDeltaJetPt12vsCent);
	214	fOutputList->Add(fHistDeltaJetPt13vsCent);
	215	fOutputList->Add(fHistDeltaJetPt23vsCent);
020052e4	216	fOutputList->Add(fHistRho1vsCent);
	217	fOutputList->Add(fHistRho2vsCent);
	218	fOutputList->Add(fHistRho3vsCent);
	219
	220	PostData(1, fOutputList);
	221
	222	}
	223	//________________________________________________________________________
	224
	225	Int_t AliAnalysisTaskEmcalJetSpectra::GetCentBin(Double_t cent) const {
	226	Int_t centbin = -1;
	227	if (cent>=0 && cent<10)
	228	centbin = 0;
	229	else if (cent>=10 && cent<20)
	230	centbin = 1;
	231	else if (cent>=20 && cent<30)
	232	centbin = 2;
	233	else if (cent>=30 && cent<40)
	234	centbin = 3;
	235	else if (cent>=40 && cent<50)
	236	centbin = 4;
	237	else if (cent>=50 && cent<90)
	238	centbin = 5;
	239	return centbin;
	240	}
	241	//________________________________________________________________________
	242	void AliAnalysisTaskEmcalJetSpectra::UserExec(Option_t *)
	243	{
	244	// Main loop, called for each event.
	245
	246	// esd or aod mode
	247	Bool_t esdMode = kTRUE;
	248	if (dynamic_cast<AliAODEvent*>(InputEvent()))
	249	esdMode = kFALSE;
	250
	251	// optimization in case autobranch loading is off
	252	AliAnalysisManager *am = AliAnalysisManager::GetAnalysisManager();
	253	if (fTracksName == "Tracks")
	254	am->LoadBranch("Tracks");
	255
	256	// get centrality
	257	Double_t fCent = -1;
	258	TList *l = InputEvent()->GetList();
	259	AliCentrality *centrality = InputEvent()->GetCentrality() ;
	260	if (centrality)
	261	fCent = centrality->GetCentralityPercentile("V0M");
	262	else
	263	fCent=99; // probably pp data
	264	if (fCent<0) {
	265	AliError(Form("Centrality negative: %f", fCent));
	266	return;
	267	}
	268
	269	Int_t centbin = GetCentBin(fCent);
	270
	271	//don't want to analyze events 90-100
	272	if (centbin < 0)
	273	return;
	274
	275	TClonesArray *jets = 0;
	276	TClonesArray *tracks = 0;
020052e4	277
	278	tracks = dynamic_cast<TClonesArray*>(l->FindObject(fTracksName));
	279	if (!tracks) {
	280	AliError(Form("Pointer to tracks %s == 0", fTracksName.Data() ));
	281	return;
	282	}
	283
	284	jets = dynamic_cast<TClonesArray*>(l->FindObject(fJetsName));
	285	if (!jets) {
	286	AliError(Form("Pointer to tracks %s == 0", fTracksName.Data() ));
	287	return;
	288	}
	289
6172360f	290	Double_t rho1 = -1;
	291	TParameter<Double_t> Rho1Param = dynamic_cast<TParameter<Double_t>>(InputEvent()->FindListObject(fRhos1Name));
	292	if (Rho1Param)
	293	rho1 = Rho1Param->GetVal();
74ec9ac9	294
6172360f	295	Double_t rho2 = -1;
	296	TParameter<Double_t> Rho2Param = dynamic_cast<TParameter<Double_t>>(InputEvent()->FindListObject(fRhos2Name));
	297	if (Rho2Param)
	298	rho2 = Rho2Param->GetVal();
020052e4	299
6172360f	300	Double_t rho3 = -1;
	301	TParameter<Double_t> Rho3Param = dynamic_cast<TParameter<Double_t>>(InputEvent()->FindListObject(fRhos3Name));
	302	if (Rho3Param)
	303	rho3 = Rho3Param->GetVal();
020052e4	304
6172360f	305	if (rho1>0)
020052e4	306	fHistRho1vsCent->Fill(fCent,rho1);
	307	if (rho2>0)
	308	fHistRho2vsCent->Fill(fCent,rho2);
	309	if (rho3>0)
	310	fHistRho3vsCent->Fill(fCent,rho3);
	311
	312	if (( rho1>0 ) && ( rho2>0 ))
	313	fHistDeltaRho12vsCent->Fill(fCent,rho1-rho2);
	314	if (( rho1>0 ) && ( rho3>0 ))
	315	fHistDeltaRho13vsCent->Fill(fCent,rho1-rho3);
	316	if (( rho2>0 ) && ( rho3>0 ))
	317	fHistDeltaRho23vsCent->Fill(fCent,rho2-rho3);
	318	fHistCentrality->Fill(fCent);
	319
	320	const Int_t Ntracks = tracks->GetEntries();
	321	Double_t lJetPt = -500;
	322	Double_t lJetPtun = -500;
	323
	324	const Int_t Njets = jets->GetEntries();
	325	Int_t NjetAcc = 0;
	326
	327
	328	for (Int_t iJets = 0; iJets < Njets; ++iJets) {
	329
	330	AliEmcalJet jet = static_cast<AliEmcalJet>(jets->At(iJets));
	331	if (!jet)
	332	continue;
	333	if ((jet->Phi()<phimin)\|\|(jet->Phi()>phimax))
	334	continue;
	335	if ((jet->Eta()<etamin)\|\|(jet->Eta()>etamax))
	336	continue;
	337	fHistAreavsRawPt[centbin]->Fill(jet->Pt(),jet->Area());
	338	if (jet->Area()<areacut)
	339	continue;
	340	//prevents 0 area jets from sneaking by when area cut == 0
	341	if (jet->Area()==0)
	342	continue;
	343	fHistRawJetPt[centbin]->Fill(jet->Pt());
	344
	345	NjetAcc++;
	346	fHistJetArea->Fill(jet->Area());
	347	fHistJetMaxPt->Fill(jet->MaxTrackPt());
	348	float Z;
	349	float jetPt1 = -500;
	350	float jetPt2 = -500;
	351	float jetPt3 = -500;
	352	if (rho1 > 0)
	353	jetPt1 = jet->Pt()-jet->Area()*rho1;
	354	if (rho2 > 0)
	355	jetPt2 = jet->Pt()-jet->Area()*rho2;
	356	if (rho3 > 0)
	357	jetPt3 = jet->Pt()-jet->Area()*rho3;
6172360f	358	if (( rho1>0 ) && ( rho2>0 ))
	359	fHistDeltaJetPt12vsCent->Fill(fCent,jetPt1-jetPt2);
	360	if (( rho1>0 ) && ( rho3>0 ))
	361	fHistDeltaJetPt13vsCent->Fill(fCent,jetPt1-jetPt3);
	362	if (( rho2>0 ) && ( rho3>0 ))
	363	fHistDeltaJetPt23vsCent->Fill(fCent,jetPt2-jetPt3);
020052e4	364
	365	if (lJetPt < jetPt2){
	366	lJetPt= jetPt2;
	367	lJetPtun = jet->Pt();
	368	}
	369
	370	if (jet->MaxTrackPt()>jet->MaxClusterPt()){
6172360f	371	Z = jet->MaxTrackPt();
020052e4	372	fHistJetMaxPtvsCent->Fill(fCent,jet->MaxTrackPt());
	373	fHistJetMaxPtvsNtrack->Fill(Ntracks,jet->MaxTrackPt());
	374	}
	375	else{
6172360f	376	Z = jet->MaxClusterPt();
020052e4	377	fHistJetMaxPtvsCent->Fill(fCent,jet->MaxClusterPt());
	378	fHistJetMaxPtvsNtrack->Fill(Ntracks,jet->MaxClusterPt());
	379	}
	380
	381	fHistJetZ->Fill(jet->MaxTrackPt()/jetPt2);
	382
	383	fHistJetNEF->Fill(jet->NEF());
	384	fHistJetPtvsCent->Fill(fCent,jet->Pt());
	385	fHistJetPtM3vsCent->Fill(fCent,jetPt2);
	386	fHistJetAreavsCent->Fill(fCent,jet->Area());
	387	fHistJetZvsCent->Fill(fCent,Z);
	388	fHistJetZvsNtrack->Fill(Ntracks,Z);
	389
74ec9ac9	390	fHistNEFvsPt[centbin][0]->Fill(jet->Pt(),jet->NEF());
6172360f	391	fHistZvsPt[centbin][0]->Fill(jet->Pt(),Z/jet->Pt());
	392	fHistZchvsPt[centbin][0]->Fill(jet->Pt(),jet->MaxTrackPt()/jet->Pt());
	393	fHistZemvsPt[centbin][0]->Fill(jet->Pt(),jet->MaxClusterPt()/jet->Pt());
74ec9ac9	394	fHistAreavsPt[centbin][0]->Fill(jet->Pt(),jet->Area());
	395	fHistJetPt[centbin][0]->Fill(jet->Pt());
	396	fHistNconsvsPt[centbin][0]->Fill(jet->Pt(),jet->N());
6172360f	397	if ((jet->MaxTrackPt()>5) \|\| (jet->MaxClusterPt()>5))
	398	fHistJetPt5[centbin][0]->Fill(jet->Pt());
	399	if ((jet->MaxTrackPt()>6) \|\| (jet->MaxClusterPt()>6))
	400	fHistJetPt6[centbin][0]->Fill(jet->Pt());
	401	if ((jet->MaxTrackPt()>7) \|\| (jet->MaxClusterPt()>7))
	402	fHistJetPt7[centbin][0]->Fill(jet->Pt());
	403	if ((jet->MaxTrackPt()>8) \|\| (jet->MaxClusterPt()>8))
	404	fHistJetPt8[centbin][0]->Fill(jet->Pt());
74ec9ac9	405
74ec9ac9	406	fHistNEFvsPt[centbin][1]->Fill(jetPt1,jet->NEF());
6172360f	407	fHistZvsPt[centbin][1]->Fill(jetPt1,Z/jetPt1);
	408	fHistZchvsPt[centbin][1]->Fill(jetPt1,jet->MaxTrackPt()/jetPt1);
	409	fHistZemvsPt[centbin][1]->Fill(jetPt1,jet->MaxClusterPt()/jetPt1);
74ec9ac9	410	fHistAreavsPt[centbin][1]->Fill(jetPt1,jet->Area());
	411	fHistJetPt[centbin][1]->Fill(jetPt1);
	412	fHistNconsvsPt[centbin][1]->Fill(jetPt1,jet->N());
6172360f	413	if ((jet->MaxTrackPt()>5) \|\| (jet->MaxClusterPt()>5))
	414	fHistJetPt5[centbin][1]->Fill(jetPt1);
	415	if ((jet->MaxTrackPt()>6) \|\| (jet->MaxClusterPt()>6))
	416	fHistJetPt6[centbin][1]->Fill(jetPt1);
	417	if ((jet->MaxTrackPt()>7) \|\| (jet->MaxClusterPt()>7))
	418	fHistJetPt7[centbin][1]->Fill(jetPt1);
	419	if ((jet->MaxTrackPt()>8) \|\| (jet->MaxClusterPt()>8))
	420	fHistJetPt8[centbin][1]->Fill(jetPt1);
74ec9ac9	421
74ec9ac9	422	fHistNEFvsPt[centbin][2]->Fill(jetPt2,jet->NEF());
6172360f	423	fHistZvsPt[centbin][2]->Fill(jetPt2,Z/jetPt2);
	424	fHistZchvsPt[centbin][2]->Fill(jetPt2,jet->MaxTrackPt()/jetPt2);
	425	fHistZemvsPt[centbin][2]->Fill(jetPt2,jet->MaxClusterPt()/jetPt2);
74ec9ac9	426	fHistAreavsPt[centbin][2]->Fill(jetPt2,jet->Area());
	427	fHistJetPt[centbin][2]->Fill(jetPt2);
	428	fHistNconsvsPt[centbin][2]->Fill(jetPt2,jet->N());
6172360f	429	if ((jet->MaxTrackPt()>5) \|\| (jet->MaxClusterPt()>5))
	430	fHistJetPt5[centbin][2]->Fill(jetPt2);
	431	if ((jet->MaxTrackPt()>6) \|\| (jet->MaxClusterPt()>6))
	432	fHistJetPt6[centbin][2]->Fill(jetPt2);
	433	if ((jet->MaxTrackPt()>7) \|\| (jet->MaxClusterPt()>7))
	434	fHistJetPt7[centbin][2]->Fill(jetPt2);
	435	if ((jet->MaxTrackPt()>8) \|\| (jet->MaxClusterPt()>8))
	436	fHistJetPt8[centbin][2]->Fill(jetPt2);
74ec9ac9	437
74ec9ac9	438	fHistNEFvsPt[centbin][3]->Fill(jetPt3,jet->NEF());
6172360f	439	fHistZvsPt[centbin][3]->Fill(jetPt3,Z/jetPt3);
	440	fHistZchvsPt[centbin][3]->Fill(jetPt3,jet->MaxTrackPt()/jetPt3);
	441	fHistZemvsPt[centbin][3]->Fill(jetPt3,jet->MaxClusterPt()/jetPt3);
74ec9ac9	442	fHistAreavsPt[centbin][3]->Fill(jetPt3,jet->Area());
	443	fHistJetPt[centbin][3]->Fill(jetPt3);
	444	fHistNconsvsPt[centbin][3]->Fill(jetPt3,jet->N());
6172360f	445	if ((jet->MaxTrackPt()>5) \|\| (jet->MaxClusterPt()>5))
	446	fHistJetPt5[centbin][3]->Fill(jetPt3);
	447	if ((jet->MaxTrackPt()>6) \|\| (jet->MaxClusterPt()>6))
	448	fHistJetPt6[centbin][3]->Fill(jetPt3);
	449	if ((jet->MaxTrackPt()>7) \|\| (jet->MaxClusterPt()>7))
	450	fHistJetPt7[centbin][3]->Fill(jetPt3);
	451	if ((jet->MaxTrackPt()>8) \|\| (jet->MaxClusterPt()>8))
	452	fHistJetPt8[centbin][3]->Fill(jetPt3);
020052e4	453
	454	fHistJetNEFvsCent->Fill(fCent,jet->NEF());
	455	fHistJetPtvsNtrack->Fill(Ntracks,jetPt2);
	456	fHistJetAreavsNtrack->Fill(Ntracks,jet->Area());
	457	fHistJetNEFvsNtrack->Fill(Ntracks,jet->NEF());
	458	}
	459
	460	fHistLeadingJetPtvsCent->Fill(fCent,lJetPtun);
	461	fHistLeadingJetPtM3vsCent->Fill(fCent,lJetPt);
	462	fHistNjetvsCent->Fill(fCent,NjetAcc);
	463	fHistNjetvsNtrack->Fill(Ntracks,NjetAcc);
	464
	465	PostData(1, fOutputList);
	466	}
	467
	468	//________________________________________________________________________
	469	void AliAnalysisTaskEmcalJetSpectra::Terminate(Option_t *)
	470	{
	471
	472	}
	473
	474