PWGJE/EMCALJetTasks/UserTasks/AliAnalysisTaskEmcalJetSpectraMECpA.cxx

   1 // $Id$
   2 //
   3 // Jet spectrum task.
   4 //
   5 // Author: R.Reed, M.Connors
   6
   7 #include "AliAnalysisTaskEmcalJetSpectraMECpA.h"
   8
   9 #include <TCanvas.h>
  10 #include <TChain.h>
  11 #include <TClonesArray.h>
  12 #include <TH1F.h>
  13 #include <TH2F.h>
  14 #include <TH3F.h>
  15 #include <THnSparse.h>
  16 #include <TList.h>
  17 #include <TLorentzVector.h>
  18 #include <TParameter.h>
  19 #include <TParticle.h>
  20 #include <TTree.h>
  21 #include <TVector3.h>
  22
  23 #include "AliAODEvent.h"
  24 #include "AliAnalysisManager.h"
  25 #include "AliAnalysisTask.h"
  26 #include "AliCentrality.h"
  27 #include "AliESDEvent.h"
  28 #include "AliESDInputHandler.h"
  29 #include "AliEmcalJet.h"
  30 #include "AliVCluster.h"
  31 #include "AliRhoParameter.h"
  32 #include "AliEmcalParticle.h"
  33
  34 ClassImp(AliAnalysisTaskEmcalJetSpectraMECpA)
  35
  36 //________________________________________________________________________
  37 AliAnalysisTaskEmcalJetSpectraMECpA::AliAnalysisTaskEmcalJetSpectraMECpA() :
  38   AliAnalysisTaskEmcalJet("spectra",kFALSE),
  39   fHistRhovsCent(0),
  40   fHistRhoScvsCent(0),
  41   fHistNjetvsCent(0)
  42 {
  43   // Default constructor.
  44   for (Int_t i = 0;i<7;++i){
  45     fHistJetPtvsTrackPt[i]      = 0;
  46     fHistJetPtScvsTrackPt[i]      = 0;
  47     fHistRawJetPtvsTrackPt[i]   = 0;
  48     fHistTrackPt[i]             = 0;
  49     fHistEP0[i]                 = 0;
  50     fHistEP0A[i]                = 0;
  51     fHistEP0C[i]                = 0;
  52     fHistEPAvsC[i]              = 0;
  53     fHistJetPtvsdEP[i]          = 0;
  54     fHistJetPtvsdEPBias[i]      = 0;
  55     fHistJetPtvsEP[i]          = 0;
  56     fHistJetPtvsEPBias[i]      = 0;
  57     fHistRhovsEP[i]             = 0;
  58     fHistJetPtEtaPhi[i]         = 0;
  59
  60   }
  61   SetMakeGeneralHistograms(kTRUE);
  62 }
  63
  64 //________________________________________________________________________
  65 AliAnalysisTaskEmcalJetSpectraMECpA::AliAnalysisTaskEmcalJetSpectraMECpA(const char *name) :
  66   AliAnalysisTaskEmcalJet(name,kTRUE),
  67   fHistRhovsCent(0),
  68   fHistRhoScvsCent(0),
  69   fHistNjetvsCent(0)
  70  {
  71    for (Int_t i = 0;i<7;++i){
  72     fHistJetPtvsTrackPt[i]      = 0;
  73     fHistJetPtScvsTrackPt[i]      = 0;
  74     fHistRawJetPtvsTrackPt[i]   = 0;
  75     fHistTrackPt[i]             = 0;
  76     fHistEP0[i]                 = 0;
  77     fHistEP0A[i]                = 0;
  78     fHistEP0C[i]                = 0;
  79     fHistEPAvsC[i]              = 0;
  80     fHistJetPtvsdEP[i]          = 0;
  81     fHistJetPtvsdEPBias[i]      = 0;
  82     fHistJetPtvsEP[i]          = 0;
  83     fHistJetPtvsEPBias[i]      = 0;
  84     fHistRhovsEP[i]             = 0;
  85     fHistJetPtEtaPhi[i]         = 0;
  86    }
  87    SetMakeGeneralHistograms(kTRUE);
  88  }
  89
  90 //________________________________________________________________________
  91 void AliAnalysisTaskEmcalJetSpectraMECpA::UserCreateOutputObjects()
  92 {
  93   if (! fCreateHisto)
  94     return;
  95   AliAnalysisTaskEmcalJet::UserCreateOutputObjects();
  96
  97   fHistRhovsCent             = new TH2F("RhovsCent",              "RhovsCent",             100, 0.0, 100.0, 500, 0, 500);
  98   fHistRhoScvsCent             = new TH2F("RhoScvsCent",              "RhoScvsCent",       100, 0.0, 100.0, 500, 0, 500);
  99   fHistNjetvsCent            = new TH2F("NjetvsCent",             "NjetvsCent",            100, 0.0, 100.0, 100, 0, 100);
 100
 101   TString name;
 102   TString title;
 103   for (Int_t i = 0;i<7;++i){
 104     name = TString(Form("JetPtvsTrackPt_%i",i));
 105     title = TString(Form("Jet pT vs Leading Track pT cent bin %i",i));
 106     fHistJetPtvsTrackPt[i] = new TH2F(name,title,1000,-500,500,100,0,100);
 107     fOutput->Add(fHistJetPtvsTrackPt[i]);
 108
 109     name = TString(Form("JetPtScvsTrackPt_%i",i));
 110     title = TString(Form("Jet pT Rho_scaled vs Leading Track pT cent bin %i",i));
 111     fHistJetPtScvsTrackPt[i] = new TH2F(name,title,1000,-500,500,100,0,100);
 112     fOutput->Add(fHistJetPtScvsTrackPt[i]);
 113
 114     name = TString(Form("RawJetPtvsTrackPt_%i",i));
 115     title = TString(Form("Raw Jet pT vs Leading Track pT cent bin %i",i));
 116     fHistRawJetPtvsTrackPt[i] = new TH2F(name,title,1000,-500,500,100,0,100);
 117     fOutput->Add(fHistRawJetPtvsTrackPt[i]);
 118     name = TString(Form("TrackPt_%i",i));
 119     title = TString(Form("Track pT cent bin %i",i));
 120     fHistTrackPt[i] = new TH1F(name,title,1000,0,200);
 121     fOutput->Add(fHistTrackPt[i]);
 122
 123     name = TString(Form("EP0_%i",i));
 124     title = TString(Form("EP VZero cent bin %i",i));
 125     fHistEP0[i] = new TH1F(name,title,100,-TMath::Pi(),TMath::Pi());
 126     fOutput->Add(fHistEP0[i]);
 127     name = TString(Form("EP0A_%i",i));
 128     title = TString(Form("EP VZero cent bin %i",i));
 129     fHistEP0A[i] = new TH1F(name,title,100,-TMath::Pi(),TMath::Pi());
 130     fOutput->Add(fHistEP0A[i]);
 131     name = TString(Form("EP0C_%i",i));
 132     title = TString(Form("EP VZero cent bin %i",i));
 133     fHistEP0C[i] = new TH1F(name,title,100,-TMath::Pi(),TMath::Pi());
 134     fOutput->Add(fHistEP0C[i]);
 135     name = TString(Form("EPAvsC_%i",i));
 136     title = TString(Form("EP VZero cent bin %i",i));
 137     fHistEPAvsC[i] = new TH2F(name,title,100,-TMath::Pi(),TMath::Pi(),100,-TMath::Pi(),TMath::Pi());
 138     fOutput->Add(fHistEPAvsC[i]);
 139     name = TString(Form("JetPtvsdEP_%i",i));
 140     title = TString(Form("Jet pt vs dEP cent bin %i",i));
 141     fHistJetPtvsdEP[i] = new TH2F(name,title,1000,-500,500,400,-2*TMath::Pi(),2*TMath::Pi());
 142     fOutput->Add(fHistJetPtvsdEP[i]);
 143     name = TString(Form("JetPtvsdEPBias_%i",i));
 144     title = TString(Form("Bias Jet pt vs dEP cent bin %i",i));
 145     fHistJetPtvsdEPBias[i] = new TH2F(name,title,1000,-500,500,400,-2*TMath::Pi(),2*TMath::Pi());
 146     fOutput->Add(fHistJetPtvsdEPBias[i]);
 147     name = TString(Form("JetPtvsEP_%i",i));
 148     title = TString(Form("Jet pt vs EP cent bin %i",i));
 149     fHistJetPtvsEP[i] = new TH2F(name,title,1000,-500,500,400,-2*TMath::Pi(),2*TMath::Pi());
 150     fOutput->Add(fHistJetPtvsEP[i]);
 151     name = TString(Form("JetPtvsEPBias_%i",i));
 152     title = TString(Form("Bias Jet pt vs EP cent bin %i",i));
 153     fHistJetPtvsEPBias[i] = new TH2F(name,title,1000,-500,500,400,-2*TMath::Pi(),2*TMath::Pi());
 154     fOutput->Add(fHistJetPtvsEPBias[i]);
 155     name = TString(Form("RhovsEP_%i",i));
 156     title = TString(Form("Rho vs EP cent bin %i",i));
 157     fHistRhovsEP[i] = new TH2F(name,title,500,0,500,400,-2*TMath::Pi(),2*TMath::Pi());
 158     fOutput->Add(fHistRhovsEP[i]);
 159
 160     name = TString(Form("JetPtEtaPhi_%i",i));
 161     title = TString(Form("JetPtEtaPhi_%i",i));
 162     fHistJetPtEtaPhi[i] = new TH3F(name,title,300,-100,200,26,-0.8,0.8,400,-2*TMath::Pi(),2*TMath::Pi());
 163     fOutput->Add(fHistJetPtEtaPhi[i]);
 164
 165
 166   }
 167
 168
 169   fOutput->Add(fHistRhovsCent);
 170   fOutput->Add(fHistRhoScvsCent);
 171   fOutput->Add(fHistNjetvsCent);
 172
 173    PostData(1, fOutput);
 174 }
 175
 176 //________________________________________________________________________
 177
 178 Int_t AliAnalysisTaskEmcalJetSpectraMECpA::GetCentBin(Double_t cent) const
 179 {
 180   // Get centrality bin.
 181
 182   Int_t centbin = -1;
 183   if (cent>=0 && cent<10)
 184     centbin = 0;
 185   else if (cent>=10 && cent<20)
 186     centbin = 1;
 187   else if (cent>=20 && cent<30)
 188     centbin = 2;
 189   else if (cent>=30 && cent<40)
 190     centbin = 3;
 191   else if (cent>=40 && cent<50)
 192     centbin = 4;
 193   else if (cent>=50 && cent<101)
 194     centbin = 5;
 195   return centbin;
 196 }
 197
 198 //________________________________________________________________________
 199
 200 Float_t AliAnalysisTaskEmcalJetSpectraMECpA:: RelativePhi(Double_t mphi,Double_t vphi) const
 201 {
 202   if (vphi < -1*TMath::Pi()) vphi += (2*TMath::Pi());
 203   else if (vphi > TMath::Pi()) vphi -= (2*TMath::Pi());
 204   if (mphi < -1*TMath::Pi()) mphi += (2*TMath::Pi());
 205   else if (mphi > TMath::Pi()) mphi -= (2*TMath::Pi());
 206   double dphi = mphi-vphi;
 207   if (dphi < -1*TMath::Pi()) dphi += (2*TMath::Pi());
 208   else if (dphi > TMath::Pi()) dphi -= (2*TMath::Pi());
 209
 210   return dphi;//dphi in [-Pi, Pi]
 211 }
 212
 213
 214 //________________________________________________________________________
 215 Bool_t AliAnalysisTaskEmcalJetSpectraMECpA::Run()
 216 {
 217
 218
 219   Int_t centbin = GetCentBin(fCent);
 220   //for pp analyses we will just use the first centrality bin
 221   if (centbin == -1)
 222     centbin = 6;
 223
 224   //if (!fTracks){
 225   //cout << "Tracks not found: " << fCent << endl;
 226   //return kTRUE;
 227   //}
 228
 229   if (fTracks){
 230     const Int_t nTrack = fTracks->GetEntriesFast();
 231     for (int i = 0;i<nTrack;i++){
 232       AliVParticle *track = static_cast<AliVParticle*>(fTracks->At(i));
 233       if (!track)
 234         continue;
 235       fHistTrackPt[centbin]->Fill(track->Pt());
 236     }
 237   }
 238
 239
 240   fHistEP0[centbin]->Fill(fEPV0);
 241   fHistEP0A[centbin]->Fill(fEPV0A);
 242   fHistEP0C[centbin]->Fill(fEPV0C);
 243   fHistEPAvsC[centbin]->Fill(fEPV0A,fEPV0C);
 244
 245   fRho = GetRhoFromEvent(fRhoName);
 246   fRhoVal = fRho->GetVal();
 247   fHistRhovsCent->Fill(fCent,fRhoVal);
 248   fHistRhovsEP[centbin]->Fill(fRhoVal,fEPV0);
 249
 250
 251   TString fRhoScaledName = fRhoName;
 252   fRhoScaledName.Append("_Scaled");
 253   Double_t fRhoScVal = 0;
 254   AliRhoParameter *fRhoScaled=GetRhoFromEvent(fRhoScaledName);
 255
 256   if(fRhoScaled){
 257     fRhoScVal=fRhoScaled->GetVal();
 258     fHistRhoScvsCent->Fill(fCent,fRhoScVal);
 259   }
 260   const Int_t Njets = fJets->GetEntriesFast();
 261
 262   Int_t NjetAcc = 0;
 263   for (Int_t iJets = 0; iJets < Njets; ++iJets) {
 264      AliEmcalJet *jet = static_cast<AliEmcalJet*>(fJets->At(iJets));
 265      if (!jet)
 266        continue;
 267      if (jet->Area()==0)
 268        continue;
 269      if (jet->Pt()<0.1)
 270        continue;
 271      if (jet->MaxTrackPt()>100)
 272        continue;
 273      if (! AcceptJet(jet))
 274        continue;
 275      //  jets.push_back(jet);
 276      NjetAcc++;
 277      Double_t jetPt = -500;
 278      jetPt = jet->Pt()-jet->Area()*fRhoVal;
 279      fHistJetPtvsTrackPt[centbin]->Fill(jetPt,jet->MaxTrackPt());
 280      if(centbin<6) fHistJetPtvsTrackPt[6]->Fill(jetPt,jet->MaxTrackPt());
 281
 282      if(fRhoScVal) {
 283        Double_t jetPtSc = -500;
 284        jetPtSc = jet->Pt()-jet->Area()*fRhoScVal;
 285        fHistJetPtScvsTrackPt[centbin]->Fill(jetPtSc,jet->MaxTrackPt());
 286        if(centbin<6) fHistJetPtScvsTrackPt[6]->Fill(jetPtSc,jet->MaxTrackPt());
 287      }
 288
 289      fHistRawJetPtvsTrackPt[centbin]->Fill(jet->Pt(),jet->MaxTrackPt());
 290      if(centbin<6) fHistRawJetPtvsTrackPt[6]->Fill(jet->Pt(),jet->MaxTrackPt());
 291
 292      fHistJetPtEtaPhi[centbin]->Fill(jet->Pt(),jet->Eta(),jet->Phi());
 293
 294      fHistJetPtvsdEP[centbin]->Fill(jetPt,RelativePhi((fEPV0+TMath::Pi()),jet->Phi()));
 295      fHistJetPtvsEP[centbin]->Fill(jetPt,fEPV0);
 296      if (jet->MaxTrackPt()>5.0){
 297        fHistJetPtvsdEPBias[centbin]->Fill(jetPt,RelativePhi((fEPV0+TMath::Pi()),jet->Phi()));
 298        fHistJetPtvsEPBias[centbin]->Fill(jetPt,fEPV0);
 299      }
 300   }
 301
 302   fHistNjetvsCent->Fill(fCent,NjetAcc);
 303   return kTRUE;
 304 }
 305
 306
 307
 308
 309