PWGJE/EMCALJetTasks/UserTasks/AliAnalysisTaskSAJF.cxx

   1 //
   2 // Jet analysis task.
   3 //
   4 // Author: S.Aiola
   5
   6 #include <TClonesArray.h>
   7 #include <TH2F.h>
   8 #include <TH3F.h>
   9 #include <THnSparse.h>
  10 #include <TList.h>
  11 #include <TLorentzVector.h>
  12
  13 #include "AliVCluster.h"
  14 #include "AliVParticle.h"
  15 #include "AliEmcalJet.h"
  16 #include "AliRhoParameter.h"
  17 #include "AliLog.h"
  18 #include "AliJetContainer.h"
  19
  20 #include "AliAnalysisTaskSAJF.h"
  21
  22 ClassImp(AliAnalysisTaskSAJF)
  23
  24 //________________________________________________________________________
  25 AliAnalysisTaskSAJF::AliAnalysisTaskSAJF() :
  26   AliAnalysisTaskEmcalJet("AliAnalysisTaskSAJF", kTRUE),
  27   fHistoType(1),
  28   fHistRejectionReason(0),
  29   fHistTracksJetPt(0),
  30   fHistClustersJetPt(0),
  31   fHistTracksPtDist(0),
  32   fHistClustersPtDist(0),
  33   fHistJetObservables(0),
  34   fHistJetPtEtaPhi(0),
  35   fHistJetPtArea(0),
  36   fHistJetPtEP(0),
  37   fHistJetPtNEF(0),
  38   fHistJetPtZ(0),
  39   fHistJetPtLeadingPartPt(0),
  40   fHistJetCorrPtEtaPhi(0),
  41   fHistJetCorrPtArea(0),
  42   fHistJetCorrPtEP(0),
  43   fHistJetCorrPtNEF(0),
  44   fHistJetCorrPtZ(0),
  45   fHistJetCorrPtLeadingPartPt(0),
  46   fHistJetPtCorrPt(0),
  47   fHistJetPtMCPt(0),
  48   fHistJetMCPtCorrPt(0)
  49 {
  50   // Default constructor.
  51
  52   SetMakeGeneralHistograms(kTRUE);
  53 }
  54
  55 //________________________________________________________________________
  56 AliAnalysisTaskSAJF::AliAnalysisTaskSAJF(const char *name) :
  57   AliAnalysisTaskEmcalJet(name, kTRUE),
  58   fHistoType(1),
  59   fHistRejectionReason(0),
  60   fHistTracksJetPt(0),
  61   fHistClustersJetPt(0),
  62   fHistTracksPtDist(0),
  63   fHistClustersPtDist(0),
  64   fHistJetObservables(0),
  65   fHistJetPtEtaPhi(0),
  66   fHistJetPtArea(0),
  67   fHistJetPtEP(0),
  68   fHistJetPtNEF(0),
  69   fHistJetPtZ(0),
  70   fHistJetPtLeadingPartPt(0),
  71   fHistJetCorrPtEtaPhi(0),
  72   fHistJetCorrPtArea(0),
  73   fHistJetCorrPtEP(0),
  74   fHistJetCorrPtNEF(0),
  75   fHistJetCorrPtZ(0),
  76   fHistJetCorrPtLeadingPartPt(0),
  77   fHistJetPtCorrPt(0),
  78   fHistJetPtMCPt(0),
  79   fHistJetMCPtCorrPt(0)
  80 {
  81   // Standard constructor.
  82
  83   SetMakeGeneralHistograms(kTRUE);
  84 }
  85
  86 //________________________________________________________________________
  87 void AliAnalysisTaskSAJF::AllocateTHnSparse()
  88 {
  89   TString title[20]= {""};
  90   Int_t nbins[20]  = {0};
  91   Double_t min[20] = {0.};
  92   Double_t max[20] = {0.};
  93   Int_t dim = 0;
  94
  95   if (fForceBeamType != kpp) {
  96     title[dim] = "Centrality (%)";
  97     nbins[dim] = 22;
  98     min[dim] = -5;
  99     max[dim] = 105;
 100     dim++;
 101
 102     title[dim] = "#phi_{jet} - #psi_{RP}";
 103     nbins[dim] = 100;
 104     min[dim] = 0;
 105     max[dim] = TMath::Pi();
 106     dim++;
 107   }
 108
 109   title[dim] = "#eta_{jet}";
 110   nbins[dim] = 100;
 111   min[dim] = -1;
 112   max[dim] = 1;
 113   dim++;
 114
 115   title[dim] = "#phi_{jet} (rad)";
 116   nbins[dim] = 201;
 117   min[dim] = 0;
 118   max[dim] = 2*TMath::Pi()*nbins[dim]/(nbins[dim]-1);
 119   dim++;
 120
 121   title[dim] = "p_{T} (GeV/c)";
 122   nbins[dim] = fNbins;
 123   min[dim] = fMinBinPt;
 124   max[dim] = fMaxBinPt;
 125   dim++;
 126
 127   if (fIsEmbedded) {
 128     title[dim] = "p_{T}^{MC} (GeV/c)";
 129     nbins[dim] = fNbins;
 130     min[dim] = fMinBinPt;
 131     max[dim] = fMaxBinPt;
 132     dim++;
 133   }
 134
 135   if (!GetRhoName().IsNull()) {
 136     title[dim] = "p_{T}^{corr} (GeV/c)";
 137     nbins[dim] = fNbins*2;
 138     min[dim] = -fMaxBinPt;
 139     max[dim] = fMaxBinPt;
 140     dim++;
 141   }
 142
 143   title[dim] = "A_{jet}";
 144   nbins[dim] = 150;
 145   min[dim] = 0;
 146   max[dim] = 1.5;
 147   dim++;
 148
 149   title[dim] = "NEF";
 150   nbins[dim] = 102;
 151   min[dim] = 0;
 152   max[dim] = 1.02;
 153   dim++;
 154
 155   title[dim] = "Z";
 156   nbins[dim] = 102;
 157   min[dim] = 0;
 158   max[dim] = 1.02;
 159   dim++;
 160
 161   title[dim] = "No. of constituents";
 162   nbins[dim] = 250;
 163   min[dim] = -0.5;
 164   max[dim] = 249.5;
 165   dim++;
 166
 167   title[dim] = "p_{T,particle}^{leading} (GeV/c)";
 168   nbins[dim] = 120;
 169   min[dim] = 0;
 170   max[dim] = 120;
 171   dim++;
 172
 173   fHistJetObservables = new THnSparseD("fHistJetObservables","fHistJetObservables",dim,nbins,min,max);
 174   fOutput->Add(fHistJetObservables);
 175   for (Int_t i = 0; i < dim; i++)
 176     fHistJetObservables->GetAxis(i)->SetTitle(title[i]);
 177 }
 178
 179 //________________________________________________________________________
 180 void AliAnalysisTaskSAJF::AllocateTHX()
 181 {
 182   fHistJetPtEtaPhi = new TH3*[fNcentBins];
 183   fHistJetPtArea = new TH2*[fNcentBins];
 184   fHistJetPtEP = new TH2*[fNcentBins];
 185   fHistJetPtNEF = new TH2*[fNcentBins];
 186   fHistJetPtZ = new TH2*[fNcentBins];
 187   fHistJetPtLeadingPartPt = new TH2*[fNcentBins];
 188   fHistJetCorrPtEtaPhi = new TH3*[fNcentBins];
 189   fHistJetCorrPtArea = new TH2*[fNcentBins];
 190   fHistJetCorrPtEP = new TH2*[fNcentBins];
 191   fHistJetCorrPtNEF = new TH2*[fNcentBins];
 192   fHistJetCorrPtZ = new TH2*[fNcentBins];
 193   fHistJetCorrPtLeadingPartPt = new TH2*[fNcentBins];
 194   fHistJetPtCorrPt = new TH2*[fNcentBins];
 195   fHistJetPtMCPt = new TH2*[fNcentBins];
 196   fHistJetMCPtCorrPt = new TH2*[fNcentBins];
 197
 198   for (Int_t i = 0; i < fNcentBins; i++) {
 199     TString histname;
 200
 201     histname = "fHistJetPtEtaPhi_";
 202     histname += i;
 203     fHistJetPtEtaPhi[i] = new TH3F(histname.Data(), histname.Data(), fNbins, fMinBinPt, fMaxBinPt, 20, -1, 1, 41, 0, 2*TMath::Pi()*41/40);
 204     fHistJetPtEtaPhi[i]->GetXaxis()->SetTitle("p_{T} (GeV/c)");
 205     fHistJetPtEtaPhi[i]->GetYaxis()->SetTitle("#eta");
 206     fHistJetPtEtaPhi[i]->GetZaxis()->SetTitle("#phi_{jet} (rad)");
 207     fOutput->Add(fHistJetPtEtaPhi[i]);
 208
 209     histname = "fHistJetPtArea_";
 210     histname += i;
 211     fHistJetPtArea[i] = new TH2F(histname.Data(), histname.Data(), fNbins, fMinBinPt, fMaxBinPt, 150, 0, 1.5);
 212     fHistJetPtArea[i]->GetXaxis()->SetTitle("p_{T} (GeV/c)");
 213     fHistJetPtArea[i]->GetYaxis()->SetTitle("A_{jet}");
 214     fHistJetPtArea[i]->GetZaxis()->SetTitle("counts");
 215     fOutput->Add(fHistJetPtArea[i]);
 216
 217     histname = "fHistJetPtEP_";
 218     histname += i;
 219     fHistJetPtEP[i] = new TH2F(histname.Data(), histname.Data(), fNbins, fMinBinPt, fMaxBinPt, 100, 0, TMath::Pi());
 220     fHistJetPtEP[i]->GetXaxis()->SetTitle("p_{T} (GeV/c)");
 221     fHistJetPtEP[i]->GetYaxis()->SetTitle("#phi_{jet} - #psi_{RP}");
 222     fHistJetPtEP[i]->GetZaxis()->SetTitle("counts");
 223     fOutput->Add(fHistJetPtEP[i]);
 224
 225     histname = "fHistJetPtNEF_";
 226     histname += i;
 227     fHistJetPtNEF[i] = new TH2F(histname.Data(), histname.Data(), fNbins, fMinBinPt, fMaxBinPt, 102, 0, 1.02);
 228     fHistJetPtNEF[i]->GetXaxis()->SetTitle("p_{T} (GeV/c)");
 229     fHistJetPtNEF[i]->GetYaxis()->SetTitle("NEF");
 230     fHistJetPtNEF[i]->GetZaxis()->SetTitle("counts");
 231     fOutput->Add(fHistJetPtNEF[i]);
 232
 233     histname = "fHistJetPtZ_";
 234     histname += i;
 235     fHistJetPtZ[i] = new TH2F(histname.Data(), histname.Data(), fNbins, fMinBinPt, fMaxBinPt, 102, 0, 1.02);
 236     fHistJetPtZ[i]->GetXaxis()->SetTitle("p_{T} (GeV/c)");
 237     fHistJetPtZ[i]->GetYaxis()->SetTitle("z");
 238     fHistJetPtZ[i]->GetZaxis()->SetTitle("counts");
 239     fOutput->Add(fHistJetPtZ[i]);
 240
 241     histname = "fHistJetPtLeadingPartPt_";
 242     histname += i;
 243     fHistJetPtLeadingPartPt[i] = new TH2F(histname.Data(), histname.Data(), fNbins, fMinBinPt, fMaxBinPt, 120, 0, 120);
 244     fHistJetPtLeadingPartPt[i]->GetXaxis()->SetTitle("p_{T} (GeV/c)");
 245     fHistJetPtLeadingPartPt[i]->GetYaxis()->SetTitle("p_{T,particle}^{leading} (GeV/c)");
 246     fHistJetPtLeadingPartPt[i]->GetZaxis()->SetTitle("counts");
 247     fOutput->Add(fHistJetPtLeadingPartPt[i]);
 248
 249     if (!GetRhoName().IsNull()) {
 250       histname = "fHistJetCorrPtEtaPhi_";
 251       histname += i;
 252       fHistJetCorrPtEtaPhi[i] = new TH3F(histname.Data(), histname.Data(), fNbins*2, -fMaxBinPt, fMaxBinPt, 20, -1, 1, 41, 0, 2*TMath::Pi()*201/200);
 253       fHistJetCorrPtEtaPhi[i]->GetXaxis()->SetTitle("p_{T} (GeV/c)");
 254       fHistJetCorrPtEtaPhi[i]->GetYaxis()->SetTitle("#eta");
 255       fHistJetCorrPtEtaPhi[i]->GetZaxis()->SetTitle("#phi_{jet} (rad)");
 256       fOutput->Add(fHistJetCorrPtEtaPhi[i]);
 257
 258       histname = "fHistJetCorrPtArea_";
 259       histname += i;
 260       fHistJetCorrPtArea[i] = new TH2F(histname.Data(), histname.Data(), fNbins*2, -fMaxBinPt, fMaxBinPt, 150, 0, 1.5);
 261       fHistJetCorrPtArea[i]->GetXaxis()->SetTitle("p_{T}^{corr} (GeV/c)");
 262       fHistJetCorrPtArea[i]->GetYaxis()->SetTitle("A_{jet}");
 263       fHistJetCorrPtArea[i]->GetZaxis()->SetTitle("counts");
 264       fOutput->Add(fHistJetCorrPtArea[i]);
 265
 266       histname = "fHistJetCorrPtEP_";
 267       histname += i;
 268       fHistJetCorrPtEP[i] = new TH2F(histname.Data(), histname.Data(), fNbins*2, -fMaxBinPt, fMaxBinPt, 100, 0, TMath::Pi());
 269       fHistJetCorrPtEP[i]->GetXaxis()->SetTitle("p_{T}^{corr} (GeV/c)");
 270       fHistJetCorrPtEP[i]->GetYaxis()->SetTitle("#phi_{jet} - #psi_{RP}");
 271       fHistJetCorrPtEP[i]->GetZaxis()->SetTitle("counts");
 272       fOutput->Add(fHistJetCorrPtEP[i]);
 273
 274       histname = "fHistJetCorrPtNEF_";
 275       histname += i;
 276       fHistJetCorrPtNEF[i] = new TH2F(histname.Data(), histname.Data(), fNbins*2, -fMaxBinPt, fMaxBinPt, 102, 0, 1.02);
 277       fHistJetCorrPtNEF[i]->GetXaxis()->SetTitle("p_{T}^{corr} (GeV/c)");
 278       fHistJetCorrPtNEF[i]->GetYaxis()->SetTitle("NEF");
 279       fHistJetCorrPtNEF[i]->GetZaxis()->SetTitle("counts");
 280       fOutput->Add(fHistJetCorrPtNEF[i]);
 281
 282       histname = "fHistJetCorrPtZ_";
 283       histname += i;
 284       fHistJetCorrPtZ[i] = new TH2F(histname.Data(), histname.Data(), fNbins*2, -fMaxBinPt, fMaxBinPt, 102, 0, 1.02);
 285       fHistJetCorrPtZ[i]->GetXaxis()->SetTitle("p_{T}^{corr} (GeV/c)");
 286       fHistJetCorrPtZ[i]->GetYaxis()->SetTitle("z");
 287       fHistJetCorrPtZ[i]->GetZaxis()->SetTitle("counts");
 288       fOutput->Add(fHistJetCorrPtZ[i]);
 289
 290       histname = "fHistJetCorrPtLeadingPartPt_";
 291       histname += i;
 292       fHistJetCorrPtLeadingPartPt[i] = new TH2F(histname.Data(), histname.Data(), fNbins*2, -fMaxBinPt, fMaxBinPt, 120, 0, 120);
 293       fHistJetCorrPtLeadingPartPt[i]->GetXaxis()->SetTitle("p_{T}^{corr} (GeV/c)");
 294       fHistJetCorrPtLeadingPartPt[i]->GetYaxis()->SetTitle("p_{T,particle}^{leading} (GeV/c)");
 295       fHistJetCorrPtLeadingPartPt[i]->GetZaxis()->SetTitle("counts");
 296       fOutput->Add(fHistJetCorrPtLeadingPartPt[i]);
 297
 298       histname = "fHistJetPtCorrPt_";
 299       histname += i;
 300       fHistJetPtCorrPt[i] = new TH2F(histname.Data(), histname.Data(), fNbins, fMinBinPt, fMaxBinPt, fNbins*2, -fMaxBinPt, fMaxBinPt);
 301       fHistJetPtCorrPt[i]->GetXaxis()->SetTitle("p_{T} (GeV/c)");
 302       fHistJetPtCorrPt[i]->GetYaxis()->SetTitle("p_{T}^{corr} (GeV/c)");
 303       fHistJetPtCorrPt[i]->GetZaxis()->SetTitle("counts");
 304       fOutput->Add(fHistJetPtCorrPt[i]);
 305
 306       if (fIsEmbedded) {
 307         histname = "fHistJetMCPtCorrPt_";
 308         histname += i;
 309         fHistJetMCPtCorrPt[i] = new TH2F(histname.Data(), histname.Data(), fNbins, fMinBinPt, fMaxBinPt, fNbins*2, -fMaxBinPt, fMaxBinPt);
 310         fHistJetMCPtCorrPt[i]->GetXaxis()->SetTitle("p_{T}^{MC} (GeV/c)");
 311         fHistJetMCPtCorrPt[i]->GetYaxis()->SetTitle("p_{T}^{corr} (GeV/c)");
 312         fHistJetMCPtCorrPt[i]->GetZaxis()->SetTitle("counts");
 313         fOutput->Add(fHistJetMCPtCorrPt[i]);
 314       }
 315     }
 316
 317     if (fIsEmbedded) {
 318       histname = "fHistJetPtMCPt_";
 319       histname += i;
 320       fHistJetPtMCPt[i] = new TH2F(histname.Data(), histname.Data(), fNbins, fMinBinPt, fMaxBinPt, fNbins, fMinBinPt, fMaxBinPt);
 321       fHistJetPtMCPt[i]->GetXaxis()->SetTitle("p_{T} (GeV/c)");
 322       fHistJetPtMCPt[i]->GetYaxis()->SetTitle("p_{T}^{MC} (GeV/c)");
 323       fHistJetPtMCPt[i]->GetZaxis()->SetTitle("counts");
 324       fOutput->Add(fHistJetPtMCPt[i]);
 325     }
 326   }
 327 }
 328
 329 //________________________________________________________________________
 330 void AliAnalysisTaskSAJF::UserCreateOutputObjects()
 331 {
 332   // Create user output.
 333
 334   AliAnalysisTaskEmcalJet::UserCreateOutputObjects();
 335
 336   if (fHistoType == 0)
 337     AllocateTHX();
 338   else
 339     AllocateTHnSparse();
 340
 341   fHistTracksJetPt = new TH2*[fNcentBins];
 342   fHistClustersJetPt = new TH2*[fNcentBins];
 343   fHistTracksPtDist = new TH2*[fNcentBins];
 344   fHistClustersPtDist = new TH2*[fNcentBins];
 345   fHistRejectionReason = new TH2*[fNcentBins];
 346
 347   for (Int_t i = 0; i < fNcentBins; i++) {
 348     TString histname;
 349
 350     if (fParticleCollArray.GetEntriesFast() > 0) {
 351       histname = "fHistTracksJetPt_";
 352       histname += i;
 353       fHistTracksJetPt[i] = new TH2F(histname.Data(), histname.Data(), fNbins / 2, fMinBinPt, fMaxBinPt / 2, fNbins, fMinBinPt, fMaxBinPt);
 354       fHistTracksJetPt[i]->GetXaxis()->SetTitle("p_{T,track} (GeV/c)");
 355       fHistTracksJetPt[i]->GetYaxis()->SetTitle("p_{T,jet} (GeV/c)");
 356       fHistTracksJetPt[i]->GetZaxis()->SetTitle("counts");
 357       fOutput->Add(fHistTracksJetPt[i]);
 358
 359       histname = "fHistTracksPtDist_";
 360       histname += i;
 361       fHistTracksPtDist[i] = new TH2F(histname.Data(), histname.Data(), fNbins / 2, fMinBinPt, fMaxBinPt / 2, 100, 0, 5);
 362       fHistTracksPtDist[i]->GetXaxis()->SetTitle("p_{T,track} (GeV/c)");
 363       fHistTracksPtDist[i]->GetYaxis()->SetTitle("d");
 364       fHistTracksPtDist[i]->GetZaxis()->SetTitle("counts");
 365       fOutput->Add(fHistTracksPtDist[i]);
 366     }
 367
 368     if (fClusterCollArray.GetEntriesFast() > 0) {
 369       histname = "fHistClustersJetPt_";
 370       histname += i;
 371       fHistClustersJetPt[i] = new TH2F(histname.Data(), histname.Data(), fNbins / 2, fMinBinPt, fMaxBinPt / 2, fNbins, fMinBinPt, fMaxBinPt);
 372       fHistClustersJetPt[i]->GetXaxis()->SetTitle("p_{T,clus} (GeV/c)");
 373       fHistClustersJetPt[i]->GetYaxis()->SetTitle("p_{T,jet} (GeV/c)");
 374       fHistClustersJetPt[i]->GetZaxis()->SetTitle("counts");
 375       fOutput->Add(fHistClustersJetPt[i]);
 376
 377       histname = "fHistClustersPtDist_";
 378       histname += i;
 379       fHistClustersPtDist[i] = new TH2F(histname.Data(), histname.Data(), fNbins / 2, fMinBinPt, fMaxBinPt / 2, 100, 0, 5);
 380       fHistClustersPtDist[i]->GetXaxis()->SetTitle("p_{T,clus} (GeV/c)");
 381       fHistClustersPtDist[i]->GetYaxis()->SetTitle("d");
 382       fHistClustersPtDist[i]->GetZaxis()->SetTitle("counts");
 383       fOutput->Add(fHistClustersPtDist[i]);
 384     }
 385
 386     histname = "fHistRejectionReason_";
 387     histname += i;
 388     fHistRejectionReason[i] = new TH2F(histname, histname, 32, 0, 32, 100, 0, 250);
 389     fHistRejectionReason[i]->GetXaxis()->SetTitle("Rejection reason");
 390     fHistRejectionReason[i]->GetYaxis()->SetTitle("p_{T,jet} (GeV/c)");
 391     fHistRejectionReason[i]->GetZaxis()->SetTitle("counts");
 392     SetRejectionReasonLabels(fHistRejectionReason[i]->GetXaxis());
 393     fOutput->Add(fHistRejectionReason[i]);
 394   }
 395
 396   PostData(1, fOutput); // Post data for ALL output slots >0 here, to get at least an empty histogram
 397
 398 }
 399
 400 //________________________________________________________________________
 401 Bool_t AliAnalysisTaskSAJF::FillHistograms()
 402 {
 403   // Fill histograms.
 404
 405   AliJetContainer *jets = static_cast<AliJetContainer*>(fJetCollArray.At(0));
 406
 407   if (!jets) return kFALSE;
 408
 409   AliEmcalJet* jet = 0;
 410
 411   jets->ResetCurrentID();
 412   while ((jet = jets->GetNextJet())) {
 413     if (!jet) {
 414       AliError("Could not receive jet!");
 415       continue;
 416     }
 417
 418     if (!jets->AcceptJet(jet)) {
 419       fHistRejectionReason[fCentBin]->Fill(jets->GetRejectionReasonBitPosition(), jet->Pt());
 420       continue;
 421     }
 422
 423     Float_t ptLeading = GetLeadingHadronPt(jet);
 424     Float_t corrPt = jet->Pt() - fRhoVal * jet->Area();
 425
 426     // Fill THnSparse
 427     Double_t ep = jet->Phi() - fEPV0;
 428     while (ep < 0) ep += TMath::Pi();
 429     while (ep >= TMath::Pi()) ep -= TMath::Pi();
 430
 431     FillJetHisto(fCent, ep, jet->Eta(), jet->Phi(), jet->Pt(), jet->MCPt(), corrPt, jet->Area(),
 432                  jet->NEF(), ptLeading/jet->Pt(), jet->GetNumberOfConstituents(), ptLeading);
 433
 434     if (fTracks) {
 435       for (Int_t it = 0; it < jet->GetNumberOfTracks(); it++) {
 436         AliVParticle *track = jet->TrackAt(it, fTracks);
 437         if (track) {
 438           fHistTracksJetPt[fCentBin]->Fill(track->Pt(), jet->Pt());
 439           Double_t dist = TMath::Sqrt((track->Eta() - jet->Eta()) * (track->Eta() - jet->Eta()) + (track->Phi() - jet->Phi()) * (track->Phi() - jet->Phi()));
 440           fHistTracksPtDist[fCentBin]->Fill(track->Pt(), dist);
 441         }
 442       }
 443     }
 444
 445     if (fCaloClusters) {
 446       for (Int_t ic = 0; ic < jet->GetNumberOfClusters(); ic++) {
 447         AliVCluster *cluster = jet->ClusterAt(ic, fCaloClusters);
 448
 449         if (cluster) {
 450           TLorentzVector nPart;
 451           cluster->GetMomentum(nPart, fVertex);
 452
 453           fHistClustersJetPt[fCentBin]->Fill(nPart.Pt(), jet->Pt());
 454           Double_t dist = TMath::Sqrt((nPart.Eta() - jet->Eta()) * (nPart.Eta() - jet->Eta()) + (nPart.Phi() - jet->Phi()) * (nPart.Phi() - jet->Phi()));
 455           fHistClustersPtDist[fCentBin]->Fill(nPart.Pt(), dist);
 456         }
 457       }
 458     }
 459   } //jet loop
 460
 461   return kTRUE;
 462 }
 463
 464 //________________________________________________________________________
 465 void AliAnalysisTaskSAJF::FillJetHisto(Double_t cent, Double_t ep, Double_t eta, Double_t phi, Double_t pt, Double_t MCpt, Double_t corrpt, Double_t area,
 466                                        Double_t NEF, Double_t z, Int_t n, Double_t leadingpt)
 467 {
 468   if (fHistoType == 0) {
 469     fHistJetPtEtaPhi[fCentBin]->Fill(pt,eta,phi);
 470     fHistJetPtArea[fCentBin]->Fill(pt,area);
 471     fHistJetPtEP[fCentBin]->Fill(pt,ep);
 472     fHistJetPtNEF[fCentBin]->Fill(pt,NEF);
 473     fHistJetPtZ[fCentBin]->Fill(pt,z);
 474     fHistJetPtLeadingPartPt[fCentBin]->Fill(pt,leadingpt);
 475     if (fHistJetCorrPtEtaPhi[fCentBin]) {
 476       fHistJetCorrPtEtaPhi[fCentBin]->Fill(corrpt,eta,phi);
 477       fHistJetCorrPtArea[fCentBin]->Fill(corrpt,area);
 478       fHistJetCorrPtEP[fCentBin]->Fill(corrpt,ep);
 479       fHistJetCorrPtNEF[fCentBin]->Fill(corrpt,NEF);
 480       fHistJetCorrPtZ[fCentBin]->Fill(corrpt,z);
 481       fHistJetCorrPtLeadingPartPt[fCentBin]->Fill(corrpt,leadingpt);
 482       fHistJetPtCorrPt[fCentBin]->Fill(pt,corrpt);
 483       if (fIsEmbedded)
 484         fHistJetMCPtCorrPt[fCentBin]->Fill(MCpt,corrpt);
 485     }
 486     if (fIsEmbedded)
 487       fHistJetPtMCPt[fCentBin]->Fill(pt,MCpt);
 488   }
 489   else {
 490
 491     Double_t contents[20]={0};
 492
 493     for (Int_t i = 0; i < fHistJetObservables->GetNdimensions(); i++) {
 494       TString title(fHistJetObservables->GetAxis(i)->GetTitle());
 495       if (title=="Centrality (%)")
 496         contents[i] = cent;
 497       else if (title=="#phi_{jet} - #psi_{RP}")
 498         contents[i] = ep;
 499       else if (title=="#eta_{jet}")
 500         contents[i] = eta;
 501       else if (title=="#phi_{jet} (rad)")
 502         contents[i] = phi;
 503       else if (title=="p_{T} (GeV/c)")
 504         contents[i] = pt;
 505       else if (title=="p_{T}^{MC} (GeV/c)")
 506         contents[i] = MCpt;
 507       else if (title=="p_{T}^{corr} (GeV/c)")
 508         contents[i] = corrpt;
 509       else if (title=="A_{jet}")
 510         contents[i] = area;
 511       else if (title=="NEF")
 512         contents[i] = NEF;
 513       else if (title=="Z")
 514         contents[i] = z;
 515       else if (title=="No. of constituents")
 516         contents[i] = n;
 517       else if (title=="p_{T,particle}^{leading} (GeV/c)")
 518         contents[i] = leadingpt;
 519       else
 520         AliWarning(Form("Unable to fill dimension %s!",title.Data()));
 521     }
 522
 523     fHistJetObservables->Fill(contents);
 524   }
 525 }