AliPhysics  aaf9c62 (aaf9c62)
AliAnalysisTaskEmcalJetSpectraQA.cxx
Go to the documentation of this file.
1 /**************************************************************************
2  * Copyright(c) 1998-2016, ALICE Experiment at CERN, All rights reserved. *
3  * *
4  * Author: The ALICE Off-line Project. *
5  * Contributors are mentioned in the code where appropriate. *
6  * *
7  * Permission to use, copy, modify and distribute this software and its *
8  * documentation strictly for non-commercial purposes is hereby granted *
9  * without fee, provided that the above copyright notice appears in all *
10  * copies and that both the copyright notice and this permission notice *
11  * appear in the supporting documentation. The authors make no claims *
12  * about the suitability of this software for any purpose. It is *
13  * provided "as is" without express or implied warranty. *
14  **************************************************************************/
15 
16 #include <TClonesArray.h>
17 #include <TH2F.h>
18 #include <TH3F.h>
19 #include <THnSparse.h>
20 #include <THashList.h>
21 
22 #include <AliAnalysisManager.h>
23 #include <AliVEventHandler.h>
24 
25 #include "AliEmcalJet.h"
26 #include "AliVCluster.h"
27 #include "AliVParticle.h"
28 #include "AliLog.h"
29 #include "AliJetContainer.h"
30 #include "AliClusterContainer.h"
31 #include "AliParticleContainer.h"
32 
34 
35 // Definitions of class AliAnalysisTaskEmcalJetSpectraQA::AliEmcalJetInfo
36 
40  fArea(0),
41  fMCPt(0),
42  fNConstituents(0),
43  fNEF(0),
44  fCent(0),
45  fEP(0),
46  fCorrPt(0),
47  fZ(0),
48  fLeadingPt(0)
49 {
50 }
51 
57  fArea(jet.Area()),
58  fMCPt(jet.MCPt()),
59  fNConstituents(jet.GetNumberOfConstituents()),
60  fNEF(jet.NEF()),
61  fCent(0),
62  fEP(0),
63  fCorrPt(0),
64  fZ(0),
65  fLeadingPt(0)
66 {
67  jet.GetMomentum(*this);
68 }
69 
70 // Definitions of class AliAnalysisTaskEmcalJetSpectraQA
71 
75 
78  AliAnalysisTaskEmcalJetLight("AliAnalysisTaskEmcalJetSpectraQA", kTRUE),
80  fJetEPaxis(kFALSE),
81  fAreaAxis(kTRUE),
82  fPtBinWidth(0.5),
83  fMaxPt(250),
84  fIsEmbedded(kFALSE),
85  fHistManager()
86 
87 {
89 }
90 
95  AliAnalysisTaskEmcalJetLight(name, kTRUE),
97  fJetEPaxis(kFALSE),
98  fAreaAxis(kTRUE),
99  fPtBinWidth(0.5),
100  fMaxPt(250),
101  fIsEmbedded(kFALSE),
102  fHistManager(name)
103 {
105 }
106 
111 {
112  AliError("Tree output not implemented. Falling back to THnSparse output");
113  AllocateTHnSparse(jets);
114 }
115 
120 {
121  Double_t jetRadius = jets->GetJetRadius();
122 
123  TString title[30]= {""};
124  Int_t nbins[30] = {0};
125  Double_t min[30] = {0.};
126  Double_t max[30] = {0.};
127  Int_t dim = 0;
128 
129  Int_t nPtBins = TMath::CeilNint(fMaxPt / fPtBinWidth);
130 
131  if (fForceBeamType != kpp) {
132  title[dim] = "Centrality (%)";
133  nbins[dim] = 20;
134  min[dim] = 0;
135  max[dim] = 100;
136  dim++;
137 
138  if (fJetEPaxis) {
139  title[dim] = "#phi_{jet} - #psi_{EP}";
140  nbins[dim] = nPtBins/5;
141  min[dim] = 0;
142  max[dim] = TMath::Pi();
143  dim++;
144  }
145  }
146 
147  title[dim] = "#eta_{jet}";
148  nbins[dim] = nPtBins/10;
149  min[dim] = -1;
150  max[dim] = 1;
151  dim++;
152 
153  title[dim] = "#phi_{jet} (rad)";
154  nbins[dim] = nPtBins/10*3;
155  min[dim] = 0;
156  max[dim] = 2*TMath::Pi();
157  dim++;
158 
159  title[dim] = "#it{p}_{T} (GeV/#it{c})";
160  nbins[dim] = nPtBins;
161  min[dim] = 0;
162  max[dim] = fMaxPt;
163  dim++;
164 
165  if (fIsEmbedded) {
166  title[dim] = "#it{p}_{T}^{MC} (GeV/#it{c})";
167  nbins[dim] = nPtBins;
168  min[dim] = 0;
169  max[dim] = fMaxPt;
170  dim++;
171  }
172 
173  if (fForceBeamType != kpp) {
174  title[dim] = "#it{p}_{T}^{corr} (GeV/#it{c})";
175  nbins[dim] = nPtBins;
176  min[dim] = -fMaxPt/2;
177  max[dim] = fMaxPt/2;
178  dim++;
179  }
180 
181  if (fAreaAxis) {
182  // area resolution is about 0.01 (w/ ghost area 0.005)
183  // for nPtBins = 250 use bin width 0.01
184  title[dim] = "#it{A}_{jet}";
185  nbins[dim] = TMath::CeilNint(2.0*jetRadius*jetRadius*TMath::Pi() / 0.01 * nPtBins / 250);
186  min[dim] = 0;
187  max[dim] = 2.0*jetRadius*jetRadius*TMath::Pi();
188  dim++;
189  }
190 
191  if (fClusterCollArray.size() > 0 && fParticleCollArray.size() > 0) {
192  title[dim] = "NEF";
193  nbins[dim] = nPtBins/5;
194  min[dim] = 0;
195  max[dim] = 1.0;
196  dim++;
197  }
198 
199  title[dim] = "#it{z}_{leading}";
200  nbins[dim] = nPtBins/5;
201  min[dim] = 0;
202  max[dim] = 1.0;
203  dim++;
204 
205  if (fForceBeamType != kpp) {
206  title[dim] = "No. of constituents";
207  nbins[dim] = 125;
208  min[dim] = 0;
209  max[dim] = 250;
210  dim++;
211  }
212  else {
213  title[dim] = "No. of constituents";
214  nbins[dim] = 50;
215  min[dim] = -0.5;
216  max[dim] = 49.5;
217  dim++;
218  }
219 
220  title[dim] = "#it{p}_{T,particle}^{leading} (GeV/#it{c})";
221  nbins[dim] = nPtBins/10*3;
222  min[dim] = 0;
223  max[dim] = 150;
224  dim++;
225 
226  TString histname = TString::Format("%s/fHistJetObservables", jets->GetArrayName().Data());
227  THnSparse* hn = fHistManager.CreateTHnSparse(histname.Data(), histname.Data(), dim, nbins, min, max);
228  for (Int_t i = 0; i < dim; i++) {
229  hn->GetAxis(i)->SetTitle(title[i]);
230  }
231 }
232 
237 {
238  TString histname;
239  TString title;
240 
241  Int_t nPtBins = TMath::CeilNint(fMaxPt / fPtBinWidth);
242 
243  for (Int_t i = 0; i < GetNCentBins(); i++) {
244  histname = TString::Format("%s/fHistJetPtEtaPhi_%d", jets->GetArrayName().Data(), i);
245  title = histname + ";#it{p}_{T} (GeV/#it{c});#eta;#phi (rad)";
246  fHistManager.CreateTH3(histname.Data(), title.Data(), 20, -1, 1, 41, 0, 2*TMath::Pi()*41/40, nPtBins, 0, fMaxPt);
247 
248  histname = TString::Format("%s/fHistJetPtArea_%d", jets->GetArrayName().Data(), i);
249  title = histname + ";#it{p}_{T} (GeV/#it{c});#it{A}_{jet};counts";
250  fHistManager.CreateTH2(histname.Data(), title.Data(), nPtBins, 0, fMaxPt, 150, 0, 1.5);
251 
252  histname = TString::Format("%s/fHistJetPtEP_%d", jets->GetArrayName().Data(), i);
253  title = histname + ";#it{p}_{T} (GeV/#it{c});#phi_{jet} - #psi_{EP};counts";
254  fHistManager.CreateTH2(histname.Data(), title.Data(), nPtBins, 0, fMaxPt, 100, 0, TMath::Pi());;
255 
256  histname = TString::Format("%s/fHistJetPtNEF_%d", jets->GetArrayName().Data(), i);
257  title = histname + ";#it{p}_{T} (GeV/#it{c});NEF;counts";
258  fHistManager.CreateTH2(histname.Data(), title.Data(), nPtBins, 0, fMaxPt, 102, 0, 1.02);
259 
260  histname = TString::Format("%s/fHistJetPtZ_%d", jets->GetArrayName().Data(), i);
261  title = histname + ";#it{p}_{T} (GeV/#it{c});#it{z}_{leading};counts";
262  fHistManager.CreateTH2(histname.Data(), title.Data(), nPtBins, 0, fMaxPt, 102, 0, 1.02);
263 
264  histname = TString::Format("%s/fHistJetPtLeadingPartPt_%d", jets->GetArrayName().Data(), i);
265  title = histname + ";#it{p}_{T} (GeV/#it{c});#it{p}_{T,particle}^{leading} (GeV/#it{c});counts";
266  fHistManager.CreateTH2(histname.Data(), title.Data(), nPtBins, 0, fMaxPt, 120, 0, 120);
267 
268  if (!jets->GetRhoName().IsNull()) {
269  histname = TString::Format("%s/fHistJetCorrPtEtaPhi_%d", jets->GetArrayName().Data(), i);
270  title = histname + ";#it{p}_{T,corr} (GeV/#it{c});#eta;#phi (rad)";
271  fHistManager.CreateTH3(histname.Data(), title.Data(), 20, -1, 1, 41, 0, 2*TMath::Pi()*201/200, nPtBins*2, -fMaxPt, fMaxPt);
272 
273  histname = TString::Format("%s/fHistJetCorrPtArea_%d", jets->GetArrayName().Data(), i);
274  title = histname + ";#it{p}_{T,corr} (GeV/#it{c});#it{A}_{jet};counts";
275  fHistManager.CreateTH2(histname.Data(), title.Data(), nPtBins, 0, fMaxPt, 150, 0, 1.5);
276 
277  histname = TString::Format("%s/fHistJetCorrPtEP_%d", jets->GetArrayName().Data(), i);
278  title = histname + ";#it{p}_{T,corr} (GeV/#it{c});#phi_{jet} - #psi_{EP};counts";
279  fHistManager.CreateTH2(histname.Data(), title.Data(), nPtBins, 0, fMaxPt, 100, 0, TMath::Pi());;
280 
281  histname = TString::Format("%s/fHistJetCorrPtNEF_%d", jets->GetArrayName().Data(), i);
282  title = histname + ";#it{p}_{T,corr} (GeV/#it{c});NEF;counts";
283  fHistManager.CreateTH2(histname.Data(), title.Data(), nPtBins, 0, fMaxPt, 102, 0, 1.02);
284 
285  histname = TString::Format("%s/fHistJetCorrPtZ_%d", jets->GetArrayName().Data(), i);
286  title = histname + ";#it{p}_{T,corr} (GeV/#it{c});#it{z}_{leading};counts";
287  fHistManager.CreateTH2(histname.Data(), title.Data(), nPtBins, 0, fMaxPt, 102, 0, 1.02);
288 
289  histname = TString::Format("%s/fHistJetCorrPtLeadingPartPt_%d", jets->GetArrayName().Data(), i);
290  title = histname + ";#it{p}_{T,corr} (GeV/#it{c});#it{p}_{T,particle}^{leading} (GeV/#it{c});counts";
291  fHistManager.CreateTH2(histname.Data(), title.Data(), nPtBins, 0, fMaxPt, 120, 0, 120);
292 
293  histname = TString::Format("%s/fHistJetPtCorrPt_%d", jets->GetArrayName().Data(), i);
294  title = histname + ";#it{p}_{T} (GeV/#it{c});#it{p}_{T,corr} (GeV/#it{c});counts";
295  fHistManager.CreateTH2(histname.Data(), title.Data(), nPtBins, 0, fMaxPt, nPtBins*2, -fMaxPt, fMaxPt);
296 
297  if (fIsEmbedded) {
298  histname = TString::Format("%s/fHistJetMCPtCorrPt_%d", jets->GetArrayName().Data(), i);
299  title = histname + ";#it{p}_{T,MC} (GeV/#it{c});#it{p}_{T,corr} (GeV/#it{c});counts";
300  fHistManager.CreateTH2(histname.Data(), title.Data(), nPtBins, 0, fMaxPt, nPtBins*2, -fMaxPt, fMaxPt);
301  }
302  }
303 
304  if (fIsEmbedded) {
305  histname = TString::Format("%s/fHistJetPtMCPt_%d", jets->GetArrayName().Data(), i);
306  title = histname + ";#it{p}_{T} (GeV/#it{c});#it{p}_{T,MC} (GeV/#it{c});counts";
307  fHistManager.CreateTH2(histname.Data(), title.Data(), nPtBins, 0, fMaxPt, nPtBins, 0, fMaxPt);
308  }
309  }
310 }
311 
314 {
316 
317  Int_t maxTracks = 6000;
318  Int_t constituentsNbins = 250;
319  Double_t constituentsMax = 249.5;
320  Double_t maxRho = 500;
321 
322  if (fForceBeamType == kpp) {
323  constituentsNbins = 50;
324  constituentsMax = 49.5;
325  maxRho = 50;
326  maxTracks = 200;
327  }
328  else if (fForceBeamType == kpA) {
329  constituentsNbins = 100;
330  constituentsMax = 99.5;
331  maxRho = 200;
332  maxTracks = 500;
333  }
334 
335  Int_t nPtBins = TMath::CeilNint(fMaxPt / fPtBinWidth);
336 
337  TString histname;
338  TString title;
339 
340  for (auto cont_it : fJetCollArray) {
341  AliJetContainer* jets = cont_it.second;
342  fHistManager.CreateHistoGroup(jets->GetArrayName());
343 
344  switch (fHistoType) {
345  case kTH2:
346  AllocateTHX(jets);
347  break;
348  case kTHnSparse:
349  AllocateTHnSparse(jets);
350  break;
351  case kTTree:
352  AllocateTTree(jets);
353  break;
354  }
355 
356  TString histname;
357 
358  for (Int_t i = 0; i < GetNCentBins(); i++) {
359  if (jets->GetParticleContainer()) {
360  histname = TString::Format("%s/fHistTracksJetPt_%d", jets->GetArrayName().Data(), i);
361  title = histname + ";#it{p}_{T,track} (GeV/#it{c});#it{p}_{T,jet} (GeV/#it{c});counts";
362  fHistManager.CreateTH2(histname.Data(), title.Data(), nPtBins / 2, 0, fMaxPt / 2, nPtBins, 0, fMaxPt);
363 
364  histname = TString::Format("%s/fHistTracksPtDist_%d", jets->GetArrayName().Data(), i);
365  title = histname + ";#it{p}_{T,track} (GeV/#it{c});#it{d};counts";
366  fHistManager.CreateTH2(histname.Data(), title.Data(), nPtBins / 2, 0, fMaxPt / 2, 100, 0, 5);
367 
368  histname = TString::Format("%s/fHistTracksZJetPtJetConst_%d", jets->GetArrayName().Data(), i);
369  title = histname + ";#it{z}_{track} (GeV/#it{c});#it{d};No. of constituents";
370  fHistManager.CreateTH3(histname.Data(), title.Data(), 120, 0.0, 1.2, nPtBins, 0, fMaxPt, constituentsNbins, -0.5, constituentsMax);
371  }
372 
373  if (jets->GetClusterContainer()) {
374  histname = TString::Format("%s/fHistClustersJetPt_%d", jets->GetArrayName().Data(), i);
375  title = histname + ";#it{p}_{T,cluster} (GeV/#it{c});#it{p}_{T,jet} (GeV/#it{c});counts";
376  fHistManager.CreateTH2(histname.Data(), title.Data(), nPtBins / 2, 0, fMaxPt / 2, nPtBins, 0, fMaxPt);
377 
378  histname = TString::Format("%s/fHistClustersPtDist_%d", jets->GetArrayName().Data(), i);
379  title = histname + ";#it{p}_{T,cluster} (GeV/#it{c});#it{d};counts";
380  fHistManager.CreateTH2(histname.Data(), title.Data(), nPtBins / 2, 0, fMaxPt / 2, 100, 0, 5);
381 
382  histname = TString::Format("%s/fHistClustersZJetPtJetConst_%d", jets->GetArrayName().Data(), i);
383  title = histname + ";#it{z}_{cluster} (GeV/#it{c});#it{p}_{T,jet} (GeV/#it{c});No. of constituents";
384  fHistManager.CreateTH3(histname.Data(), title.Data(), 120, 0.0, 1.2, nPtBins, 0, fMaxPt, constituentsNbins, -0.5, constituentsMax);
385  }
386 
387  histname = TString::Format("%s/fHistRejectionReason_%d", jets->GetArrayName().Data(), i);
388  title = histname + ";Rejection reason;#it{p}_{T,jet} (GeV/#it{c});counts";
389  TH2* hist = fHistManager.CreateTH2(histname.Data(), title.Data(), 32, 0, 32, 100, 0, 250);
390  SetRejectionReasonLabels(hist->GetXaxis());
391  }
392 
393  histname = TString::Format("%s/fHistLeadJetPtVsCent", jets->GetArrayName().Data());
394  title = histname + ";Centrality (%);#it{p}_{T,jet} (GeV/#it{c});counts";
395  fHistManager.CreateTH2(histname.Data(), title.Data(), 100, 0, 100, nPtBins, 0, fMaxPt);
396 
397  histname = TString::Format("%s/fHistLeadJetPtVsNTracks", jets->GetArrayName().Data());
398  title = histname + ";no. of tracks;#it{p}_{T,jet} (GeV/#it{c});counts";
399  fHistManager.CreateTH2(histname.Data(), title.Data(), 200, 0, maxTracks, nPtBins, 0, fMaxPt);
400 
401  if (!jets->GetRhoName().IsNull()) {
402  histname = TString::Format("%s/fHistRhoVsCent", jets->GetArrayName().Data());
403  title = histname + ";Centrality (%);#rho (GeV/#it{c} rad^{-1});counts";
404  fHistManager.CreateTH2(histname.Data(), title.Data(), 100, 0, 100, 1000, 0, maxRho);
405 
406  histname = TString::Format("%s/fHistRhoVsNTracks", jets->GetArrayName().Data());
407  title = histname + ";no. of tracks;#rho (GeV/#it{c} rad^{-1});counts";
408  fHistManager.CreateTH2(histname.Data(), title.Data(), 200, 0, maxTracks, 1000, 0, maxRho);
409 
410  histname = TString::Format("%s/fHistRhoVsLeadJetPt", jets->GetArrayName().Data());
411  title = histname + ";#it{p}_{T,jet} (GeV/#it{c});#rho (GeV/#it{c} rad^{-1});counts";
412  fHistManager.CreateTH2(histname.Data(), title.Data(), nPtBins, 0, fMaxPt, 1000, 0, maxRho);
413  }
414  }
415 
416  TIter nextElement(fHistManager.GetListOfHistograms());
417  TObject* obj = 0;
418  while ((obj = nextElement())) fOutput->Add(obj);
419  PostData(1, fOutput);
420 }
421 
426 {
427  TString histname;
428 
429  if (fCentBin < 0) {
430  AliError(Form("fCentBin is %d! fCent = %.3f. Fix the centrality bins to include all possible values of centrality.", fCentBin, fCent));
431  return kFALSE;
432  }
433 
434  for (auto cont_it : fJetCollArray) {
435  AliJetContainer* jets = cont_it.second;
436  Double_t rhoVal = 0;
437  if (jets->GetRhoParameter()) rhoVal = jets->GetRhoVal();
438 
439  Double_t leadJetPt = 0;
440  for (auto jet : jets->accepted()) {
441 
442  UInt_t rejectionReason = 0;
443  if (!jets->AcceptJet(jet, rejectionReason)) {
444  histname = TString::Format("%s/fHistRejectionReason_%d", jets->GetArrayName().Data(), fCentBin);
445  fHistManager.FillTH2(histname.Data(), jets->GetRejectionReasonBitPosition(rejectionReason), jet->Pt());
446  continue;
447  }
448 
449  Float_t ptLeading = jets->GetLeadingHadronPt(jet);
450  Float_t corrPt = jet->Pt() - rhoVal * jet->Area();
451 
452  TLorentzVector leadPart;
453 
454  jets->GetLeadingHadronMomentum(leadPart, jet);
455 
456  // Fill THnSparse
457  Double_t ep = jet->Phi() - fEPV0;
458  while (ep < 0) ep += TMath::Pi();
459  while (ep >= TMath::Pi()) ep -= TMath::Pi();
460 
461  Double_t z = GetParallelFraction(leadPart.Vect(), jet);
462  if (z == 1 || (z > 1 && z - 1 < 1e-3)) z = 0.999; // so that it will contribute to the bin 0.9-1 rather than 1-1.1
463 
464  AliEmcalJetInfo jetInfo(*jet);
465  jetInfo.fCent = fCent;
466  jetInfo.fEP = ep;
467  jetInfo.fCorrPt = corrPt;
468  jetInfo.fZ = z;
469  jetInfo.fLeadingPt = ptLeading;
470 
471  if (jet->Pt() > leadJetPt) leadJetPt = jet->Pt();
472 
473  FillJetHisto(jetInfo, jets);
474 
475  AliParticleContainer* tracks = jets->GetParticleContainer();
476  if (tracks) {
477  for (Int_t it = 0; it < jet->GetNumberOfTracks(); it++) {
478  AliVParticle *track = jet->TrackAt(it, tracks->GetArray());
479  if (track) {
480  histname = TString::Format("%s/fHistTracksJetPt_%d", jets->GetArrayName().Data(), fCentBin);
481  fHistManager.FillTH2(histname.Data(), track->Pt(), jet->Pt());
482 
483  Double_t dphi = TVector2::Phi_0_2pi(track->Phi() - jet->Phi());
484  Double_t deta = track->Eta() - jet->Eta();
485  Double_t dist = TMath::Sqrt(deta * deta + dphi * dphi);
486 
487  histname = TString::Format("%s/fHistTracksPtDist_%d", jets->GetArrayName().Data(), fCentBin);
488  fHistManager.FillTH2(histname.Data(), track->Pt(), dist);
489 
490  histname = TString::Format("%s/fHistTracksZJetPtJetConst_%d", jets->GetArrayName().Data(), fCentBin);
491  fHistManager.FillTH3(histname.Data(), GetParallelFraction(track, jet), jet->Pt(), jet->GetNumberOfConstituents());
492  }
493  }
494  }
495 
496  AliClusterContainer* clusters = jets->GetClusterContainer();
497  if (clusters) {
498  for (Int_t ic = 0; ic < jet->GetNumberOfClusters(); ic++) {
499  AliVCluster *cluster = jet->ClusterAt(ic, clusters->GetArray());
500 
501  if (cluster) {
502  TLorentzVector nPart;
503  if (clusters->GetDefaultClusterEnergy() >=0 && clusters->GetDefaultClusterEnergy() < AliVCluster::kLastUserDefEnergy) {
504  cluster->GetMomentum(nPart, fVertex, (AliVCluster::VCluUserDefEnergy_t)clusters->GetDefaultClusterEnergy());
505  }
506  else {
507  cluster->GetMomentum(nPart, fVertex);
508  }
509 
510  histname = TString::Format("%s/fHistClustersJetPt_%d", jets->GetArrayName().Data(), fCentBin);
511  fHistManager.FillTH2(histname.Data(), nPart.Pt(), jet->Pt());
512 
513  Double_t dphi = TVector2::Phi_0_2pi(nPart.Phi() - jet->Phi());
514  Double_t deta = nPart.Eta() - jet->Eta();
515  Double_t dist = TMath::Sqrt(deta * deta + dphi * dphi);
516 
517  histname = TString::Format("%s/fHistClustersPtDist_%d", jets->GetArrayName().Data(), fCentBin);
518  fHistManager.FillTH2(histname.Data(), nPart.Pt(), dist);
519 
520  histname = TString::Format("%s/fHistClustersZJetPtJetConst_%d", jets->GetArrayName().Data(), fCentBin);
521  fHistManager.FillTH3(histname.Data(), GetParallelFraction(nPart.Vect(), jet), jet->Pt(), jet->GetNumberOfConstituents());
522  }
523  }
524  }
525  } //jet loop
526 
527  Int_t ntracks = 0;
528  for (auto cont : this->fParticleCollArray) ntracks += cont.second->GetNAcceptEntries();
529 
530  histname = TString::Format("%s/fHistLeadJetPtVsCent", jets->GetArrayName().Data());
531  fHistManager.FillTH2(histname.Data(), fCent, leadJetPt);
532 
533  histname = TString::Format("%s/fHistLeadJetPtVsNTracks", jets->GetArrayName().Data());
534  fHistManager.FillTH2(histname.Data(), ntracks, leadJetPt);
535 
536  if (jets->GetRhoParameter()) {
537  histname = TString::Format("%s/fHistRhoVsCent", jets->GetArrayName().Data());
538  fHistManager.FillTH2(histname.Data(), fCent, rhoVal);
539 
540  histname = TString::Format("%s/fHistRhoVsNTracks", jets->GetArrayName().Data());
541  fHistManager.FillTH2(histname.Data(), ntracks, rhoVal);
542 
543  histname = TString::Format("%s/fHistRhoVsLeadJetPt", jets->GetArrayName().Data());
544  fHistManager.FillTH2(histname.Data(), leadJetPt, rhoVal);
545  }
546  }
547  return kTRUE;
548 }
549 
555 {
556  TString histname;
557 
558  histname = TString::Format("%s/fHistJetPtEtaPhi_%d", jets->GetArrayName().Data(), fCentBin);
559  fHistManager.FillTH3(histname.Data(), jet.Eta(), jet.Phi_0_2pi(), jet.Pt());
560 
561  histname = TString::Format("%s/fHistJetPtArea_%d", jets->GetArrayName().Data(), fCentBin);
562  fHistManager.FillTH2(histname.Data(), jet.Pt(), jet.fArea);
563 
564  histname = TString::Format("%s/fHistJetPtEP_%d", jets->GetArrayName().Data(), fCentBin);
565  fHistManager.FillTH2(histname.Data(), jet.Pt(), jet.fEP);
566 
567  histname = TString::Format("%s/fHistJetPtNEF_%d", jets->GetArrayName().Data(), fCentBin);
568  fHistManager.FillTH2(histname.Data(), jet.Pt(), jet.fNEF);
569 
570  histname = TString::Format("%s/fHistJetPtZ_%d", jets->GetArrayName().Data(), fCentBin);
571  fHistManager.FillTH2(histname.Data(), jet.Pt(), jet.fZ);
572 
573  histname = TString::Format("%s/fHistJetPtLeadingPartPt_%d", jets->GetArrayName().Data(), fCentBin);
574  fHistManager.FillTH2(histname.Data(), jet.Pt(), jet.fLeadingPt);
575 
576  if (fIsEmbedded) {
577  histname = TString::Format("%s/fHistJetPtMCPt_%d", jets->GetArrayName().Data(), fCentBin);
578  fHistManager.FillTH2(histname.Data(), jet.Pt(), jet.fMCPt);
579  }
580 
581  if (!jets->GetRhoName().IsNull()) {
582  histname = TString::Format("%s/fHistJetCorrPtEtaPhi_%d", jets->GetArrayName().Data(), fCentBin);
583  fHistManager.FillTH3(histname.Data(), jet.Eta(), jet.Phi_0_2pi(), jet.fCorrPt);
584 
585  histname = TString::Format("%s/fHistJetCorrPtArea_%d", jets->GetArrayName().Data(), fCentBin);
586  fHistManager.FillTH2(histname.Data(), jet.fCorrPt, jet.fArea);
587 
588  histname = TString::Format("%s/fHistJetCorrPtEP_%d", jets->GetArrayName().Data(), fCentBin);
589  fHistManager.FillTH2(histname.Data(), jet.fCorrPt, jet.fEP);
590 
591  histname = TString::Format("%s/fHistJetCorrPtNEF_%d", jets->GetArrayName().Data(), fCentBin);
592  fHistManager.FillTH2(histname.Data(), jet.fCorrPt, jet.fNEF);
593 
594  histname = TString::Format("%s/fHistJetCorrPtZ_%d", jets->GetArrayName().Data(), fCentBin);
595  fHistManager.FillTH2(histname.Data(), jet.fCorrPt, jet.fZ);
596 
597  histname = TString::Format("%s/fHistJetCorrPtLeadingPartPt_%d", jets->GetArrayName().Data(), fCentBin);
598  fHistManager.FillTH2(histname.Data(), jet.fCorrPt, jet.fLeadingPt);
599 
600  histname = TString::Format("%s/fHistJetPtCorrPt_%d", jets->GetArrayName().Data(), fCentBin);
601  fHistManager.FillTH2(histname.Data(), jet.Pt(), jet.fCorrPt);
602 
603  if (fIsEmbedded) {
604  histname = TString::Format("%s/fHistJetMCPtCorrPt_%d", jets->GetArrayName().Data(), fCentBin);
605  fHistManager.FillTH2(histname.Data(), jet.fMCPt, jet.fCorrPt);
606  }
607  }
608 }
609 
615 {
616  AliError("Tree output not implemented. Falling back to THnSparse output");
617  FillTHnSparse(jet, jets);
618 }
619 
625 {
626  TString histname;
627  Double_t contents[30]={0};
628 
629  histname = TString::Format("%s/fHistJetObservables", jets->GetArrayName().Data());
630  THnSparse* histJetObservables = static_cast<THnSparse*>(fHistManager.FindObject(histname));
631 
632  if (!histJetObservables) return;
633 
634  for (Int_t i = 0; i < histJetObservables->GetNdimensions(); i++) {
635  TString title(histJetObservables->GetAxis(i)->GetTitle());
636  if (title=="Centrality (%)")
637  contents[i] = jet.fCent;
638  else if (title=="#phi_{jet} - #psi_{EP}")
639  contents[i] = jet.fEP;
640  else if (title=="#eta_{jet}")
641  contents[i] = jet.Eta();
642  else if (title=="#phi_{jet} (rad)")
643  contents[i] = jet.Phi_0_2pi();
644  else if (title=="#it{p}_{T} (GeV/#it{c})")
645  contents[i] = jet.Pt();
646  else if (title=="#it{p}_{T}^{MC} (GeV/#it{c})")
647  contents[i] = jet.fMCPt;
648  else if (title=="#it{p}_{T}^{corr} (GeV/#it{c})")
649  contents[i] = jet.fCorrPt;
650  else if (title=="#it{A}_{jet}")
651  contents[i] = jet.fArea;
652  else if (title=="NEF")
653  contents[i] = jet.fNEF;
654  else if (title=="#it{z}_{leading}")
655  contents[i] = jet.fZ;
656  else if (title=="No. of constituents")
657  contents[i] = jet.fNConstituents;
658  else if (title=="#it{p}_{T,particle}^{leading} (GeV/#it{c})")
659  contents[i] = jet.fLeadingPt;
660  else
661  AliWarning(Form("Unable to fill dimension %s!",title.Data()));
662  }
663 
664  histJetObservables->Fill(contents);
665 }
666 
672 {
673  switch (fHistoType) {
674  case kTH2:
675  FillTHX(jet, jets);
676  break;
677 
678  case kTHnSparse:
679  FillTHnSparse(jet, jets);
680  break;
681 
682  case kTTree:
683  FillTTree(jet, jets);
684  break;
685  }
686 }
687 
697  Double_t trackPtCut, Double_t clusECut, TString suffix)
698 {
699  // Get the pointer to the existing analysis manager via the static access method.
700  AliAnalysisManager *mgr = AliAnalysisManager::GetAnalysisManager();
701  if (!mgr) {
702  ::Error("AddTaskEmcalJetSpectraQA", "No analysis manager to connect to.");
703  return nullptr;
704  }
705 
706  // Check the analysis type using the event handlers connected to the analysis manager.
707  AliVEventHandler* handler = mgr->GetInputEventHandler();
708  if (!handler) {
709  ::Error("AddTaskEmcalJetSpectraQA", "This task requires an input event handler");
710  return nullptr;
711  }
712 
714 
715  if (handler->InheritsFrom("AliESDInputHandler")) {
716  dataType = kESD;
717  }
718  else if (handler->InheritsFrom("AliAODInputHandler")) {
719  dataType = kAOD;
720  }
721 
722  // Init the task and do settings
723 
724  if (trackName == "usedefault") {
725  if (dataType == kESD) {
726  trackName = "Tracks";
727  }
728  else if (dataType == kAOD) {
729  trackName = "tracks";
730  }
731  else {
732  trackName = "";
733  }
734  }
735 
736  if (clusName == "usedefault") {
737  if (dataType == kESD) {
738  clusName = "CaloClusters";
739  }
740  else if (dataType == kAOD) {
741  clusName = "caloClusters";
742  }
743  else {
744  clusName = "";
745  }
746  }
747 
748  TString name("AliAnalysisTaskEmcalJetSpectraQA");
749  if (strcmp(suffix,"")) {
750  name += "_";
751  name += suffix;
752  }
753 
755  jetTask->SetVzRange(-10,10);
756  jetTask->SetNeedEmcalGeom(kFALSE);
757  AliParticleContainer *partCont = jetTask->AddParticleContainer(trackName.Data());
758  if (partCont) partCont->SetParticlePtCut(trackPtCut);
759 
760  AliClusterContainer *clusterCont = jetTask->AddClusterContainer(clusName.Data());
761  if (clusterCont) {
762  clusterCont->SetClusECut(0.);
763  clusterCont->SetClusPtCut(0.);
764  clusterCont->SetClusHadCorrEnergyCut(clusECut);
765  clusterCont->SetDefaultClusterEnergy(AliVCluster::kHadCorr);
766  }
767 
768  // Final settings, pass to manager and set the containers
769  mgr->AddTask(jetTask);
770 
771  // Create containers for input/output
772  AliAnalysisDataContainer *cinput1 = mgr->GetCommonInputContainer() ;
773  TString contname(name);
774  contname += "_histos";
775  AliAnalysisDataContainer *coutput1 = mgr->CreateContainer(contname.Data(),
776  TList::Class(),AliAnalysisManager::kOutputContainer,
777  Form("%s", AliAnalysisManager::GetCommonFileName()));
778  mgr->ConnectInput (jetTask, 0, cinput1 );
779  mgr->ConnectOutput (jetTask, 1, coutput1 );
780 
781  return jetTask;
782 }
Declaration of class AliAnalysisTaskEmcalJetSpectraQA.
AliClusterContainer * AddClusterContainer(std::string branchName, std::string contName="")
THashList * CreateHistoGroup(const char *groupname)
Create a new group of histograms within a parent group.
void SetParticlePtCut(Double_t cut)
Double_t GetRhoVal() const
const TString & GetRhoName() const
double Double_t
Definition: External.C:58
EDataType_t
Switch for the data type.
const char * title
Definition: MakeQAPdf.C:27
static AliAnalysisTaskEmcalJetSpectraQA * AddTaskEmcalJetSpectraQA(TString ntracks="usedefault", TString nclusters="usedefault", Double_t trackPtCut=0.15, Double_t clusECut=0.30, TString suffix="")
Double_t fEPV0
!event plane V0
Float_t fPtBinWidth
Histogram pt bin width.
void FillTH2(const char *hname, double x, double y, double weight=1., Option_t *opt="")
Fill a 2D histogram within the container.
virtual void FillTHX(const AliEmcalJetInfo &jetInfo, const AliJetContainer *jets)
virtual void FillTHnSparse(const AliEmcalJetInfo &jetInfo, const AliJetContainer *jets)
AliClusterContainer * GetClusterContainer() const
void FillTH3(const char *hname, double x, double y, double z, double weight=1., Option_t *opt="")
Fill a 3D histogram within the container.
virtual Bool_t AcceptJet(Int_t i, UInt_t &rejectionReason) const
void UserCreateOutputObjects()
Overloads base class method. Creates output objects.
Container for particles within the EMCAL framework.
Int_t GetDefaultClusterEnergy() const
virtual void AllocateTTree(const AliJetContainer *jets)
void GetMomentum(TLorentzVector &vec) const
EBeamType_t fForceBeamType
forced beam type
const Int_t nPtBins
AliParticleContainer * GetParticleContainer() const
void GetLeadingHadronMomentum(TLorentzVector &mom, const AliEmcalJet *jet) const
void FillJetHisto(const AliEmcalJetInfo &jetInfo, const AliJetContainer *jets)
TH2 * CreateTH2(const char *name, const char *title, int nbinsx, double xmin, double xmax, int nbinsy, double ymin, double ymax, Option_t *opt="")
Create a new TH2 within the container.
TObject * FindObject(const char *name) const
Find an object inside the container.
int Int_t
Definition: External.C:63
THistManager fHistManager
Histogram manager.
unsigned int UInt_t
Definition: External.C:33
THashList * GetListOfHistograms() const
Get the list of histograms.
Definition: THistManager.h:671
Double_t GetLeadingHadronPt(const AliEmcalJet *jet) const
float Float_t
Definition: External.C:68
Double_t Phi_0_2pi() const
Double_t fVertex[3]
!event vertex
virtual void FillTTree(const AliEmcalJetInfo &jetInfo, const AliJetContainer *jets)
Base task in the EMCAL jet framework (lighter version of AliAnalysisTaskEmcalJet) ...
Int_t fCentBin
!event centrality bin
Bool_t fJetEPaxis
whether a EP-jet axis should be included in the THnSparse
Bool_t fAreaAxis
whether the area axis should be included
AliRhoParameter * GetRhoParameter()
AliParticleContainer * AddParticleContainer(std::string branchName, std::string contName="")
Float_t GetJetRadius() const
Definition: External.C:220
std::map< std::string, AliParticleContainer * > fParticleCollArray
particle/track collection array
void SetVzRange(Double_t min, Double_t max)
Bool_t GetMomentum(TLorentzVector &mom, const AliVCluster *vc, Double_t mass) const
Represent a jet reconstructed using the EMCal jet framework.
Definition: AliEmcalJet.h:51
void SetClusPtCut(Double_t cut)
AliAnalysisTaskEmcalJetSpectraQA()
Default constructor for ROOT I/O purposes.
std::map< std::string, AliJetContainer * > fJetCollArray
jet collection array
const Int_t nbins
const AliJetIterableContainer accepted() const
bool Bool_t
Definition: External.C:53
void SetClusECut(Double_t cut)
void SetDefaultClusterEnergy(Int_t d)
std::map< std::string, AliClusterContainer * > fClusterCollArray
cluster collection array
THnSparse * CreateTHnSparse(const char *name, const char *title, int ndim, const int *nbins, const double *min, const double *max, Option_t *opt="")
Create a new THnSparse within the container.
Container structure for EMCAL clusters.
virtual void AllocateTHnSparse(const AliJetContainer *jets)
virtual void AllocateTHX(const AliJetContainer *jets)
Double_t fCent
!event centrality
Container for jet within the EMCAL jet framework.
Implementation of a task to perform QA on jet spectra.
TH3 * CreateTH3(const char *name, const char *title, int nbinsx, double xmin, double xmax, int nbinsy, double ymin, double ymax, int nbinsz, double zmin, double zmax, Option_t *opt="")
Create a new TH2 within the container.
void SetClusHadCorrEnergyCut(Double_t cut)
static Double_t GetParallelFraction(AliVParticle *part1, AliVParticle *part2)