AliPhysics  8d00e07 (8d00e07)
AliAnalysisTaskPWGJEQA.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 <cstring>
17 
18 #include <TClonesArray.h>
19 #include <TH1F.h>
20 #include <TH2F.h>
21 #include <THnSparse.h>
22 #include <TMath.h>
23 #include <TString.h>
24 #include <Riostream.h>
25 #include <TChain.h>
26 #include <TFile.h>
27 
28 #include <AliVEventHandler.h>
29 #include <AliAnalysisManager.h>
30 #include "AliParticleContainer.h"
31 #include "AliClusterContainer.h"
32 #include "AliTrackContainer.h"
33 #include "AliCentrality.h"
34 #include "AliVCluster.h"
35 #include "AliVParticle.h"
36 #include "AliVTrack.h"
37 #include "AliLog.h"
38 #include "AliEMCALGeometry.h"
39 #include "AliEMCALGeoParams.h"
40 #include "AliESDtrack.h"
41 #include "AliAODTrack.h"
42 #include "AliAODMCParticle.h"
43 #include "AliMCParticleContainer.h"
44 #include "AliEmcalJet.h"
45 #include "AliJetContainer.h"
46 #include "AliMCEvent.h"
47 #include "AliGenPythiaEventHeader.h"
48 #include "AliAODMCHeader.h"
49 #include "AliPHOSGeometry.h"
50 #include "AliOADBContainer.h"
51 
52 #include "AliAnalysisTaskPWGJEQA.h"
53 
55 ClassImp(AliAnalysisTaskPWGJEQA);
57 
58 //________________________________________________________________________
61  AliAnalysisTaskEmcalJet("AliAnalysisTaskPWGJEQA", kTRUE),
62  fCellEnergyCut(0.05),
63  fMaxPt(250),
64  fNTotTracks(0),
65  fLeadingTrack(),
66  fDoTrackQA(kTRUE),
67  fDoCaloQA(kTRUE),
68  fDoJetQA(kTRUE),
69  fDoEventQA(kTRUE),
70  fGeneratorLevelName(),
71  fDetectorLevelName(),
72  fRejectOutlierEvents(kFALSE),
73  fIsPtHard(kFALSE),
74  fGeneratorLevel(0),
75  fDetectorLevel(0),
76  fNPtHistBins(0),
77  fPtHistBins(0),
78  fNEtaHistBins(0),
79  fEtaHistBins(0),
80  fNPhiHistBins(0),
81  fPhiHistBins(0),
82  fNCentHistBins(0),
83  fCentHistBins(0),
84  fNPtRelDiffHistBins(0),
85  fPtRelDiffHistBins(0),
86  fNPtResHistBins(0),
87  fPtResHistBins(0),
88  fNIntegerHistBins(0),
89  fIntegerHistBins(0),
90  fPHOSGeo(nullptr),
91  fHistManager("AliAnalysisTaskPWGJEQA")
92 {
93  // Default constructor.
94  memset(fNTotClusters, 0, sizeof(Int_t)*3);
97 }
98 
99 //________________________________________________________________________
104  AliAnalysisTaskEmcalJet(name, kTRUE),
105  fCellEnergyCut(0.05),
106  fMaxPt(250),
107  fNTotTracks(0),
108  fLeadingTrack(),
109  fDoTrackQA(kTRUE),
110  fDoCaloQA(kTRUE),
111  fDoJetQA(kTRUE),
112  fDoEventQA(kTRUE),
115  fRejectOutlierEvents(kFALSE),
116  fIsPtHard(kFALSE),
117  fGeneratorLevel(0),
118  fDetectorLevel(0),
119  fNPtHistBins(0),
120  fPtHistBins(0),
121  fNEtaHistBins(0),
122  fEtaHistBins(0),
123  fNPhiHistBins(0),
124  fPhiHistBins(0),
125  fNCentHistBins(0),
126  fCentHistBins(0),
129  fNPtResHistBins(0),
130  fPtResHistBins(0),
132  fIntegerHistBins(0),
133  fPHOSGeo(nullptr),
134  fHistManager(name)
135 {
136  // Standard
137  memset(fNTotClusters, 0, sizeof(Int_t)*3);
140 }
141 
142 //________________________________________________________________________
144 {
145  // Destructor
146 }
147 
148 //________________________________________________________________________
150 {
151  // Create histograms
153 
154  // Set track container pointers
157 
158  // Allocate histograms for tracks, cells, clusters, jets
164 
165  TIter next(fHistManager.GetListOfHistograms());
166  TObject* obj = 0;
167  while ((obj = next())) {
168  fOutput->Add(obj);
169  }
170 
171  PostData(1, fOutput); // Post data for ALL output slots >0 here, to get at least an empty histogram
172 }
173 
174 //________________________________________________________________________
176 
178  if (fGeneratorLevel) {
181  }
182 }
183 
184 //________________________________________________________________________
186 
187  TString histname;
188  TString title;
189 
190  if (!fCaloCellsName.IsNull()) {
191 
192  histname = TString::Format("%s/fHistCellEnergy", fCaloCellsName.Data());
193  title = histname + ";#it{E}_{cell} (GeV);counts";
194  fHistManager.CreateTH1(histname.Data(), title.Data(), 500,0, 100);
195 
196  histname = TString::Format("%s/fProfCellAbsIdEnergy", fCaloCellsName.Data());
197  title = histname + ";cell absId;<#it{E}_{cell}> (GeV)";
198  fHistManager.CreateTProfile(histname.Data(), title.Data(), 18000,0,18000);
199 
200  histname = TString::Format("%s/fHistCellTime", fCaloCellsName.Data());
201  title = histname + ";#it{t}_{cell} (s);counts";
202  fHistManager.CreateTH1(histname.Data(), title.Data(), 500,-5e-6, 5e-6);
203 
204  histname = TString::Format("%s/fProfCellAbsIdTime", fCaloCellsName.Data());
205  title = histname + ";cell absId;<#it{t}_{cell}> (s)";
206  fHistManager.CreateTProfile(histname.Data(), title.Data(), 18000,0,18000);
207 
208  histname = TString::Format("%s/fHistCellEvsTime", fCaloCellsName.Data());
209  title = histname + ";<#it{t}_{cell}> (s);<#it{E}_{cell}> (GeV)";
210  fHistManager.CreateTH2(histname.Data(), title.Data(), 1000, -1e-6, 1e-6, 200, 0, 100);
211 
212  }
213 }
214 
215 //________________________________________________________________________
217 
218  AliEmcalContainer* cont = 0;
219  TString histname;
220  TString title;
221  Int_t nPtBins = 2*fMaxPt;
222 
223  TIter nextClusColl(&fClusterCollArray);
224  while ((cont = static_cast<AliEmcalContainer*>(nextClusColl()))) {
225 
226  histname = TString::Format("%s/fHistClusterRejectionReason", cont->GetArrayName().Data());
227  title = histname + ";Rejection reason;#it{E}_{clus} (GeV/);counts";
228  TH2* hist = fHistManager.CreateTH2(histname.Data(), title.Data(), 32, 0, 32, 50, 0, 250);
229  SetRejectionReasonLabels(hist->GetXaxis());
230 
231  histname = TString::Format("%s/hPhosNCellsVsEnergy", cont->GetArrayName().Data());
232  title = histname + ";N cells;#it{E}_{clus} (GeV/);counts";
233  fHistManager.CreateTH2(histname.Data(), title.Data(), 40, 0, 40, nPtBins, 0, fMaxPt);
234 
235  histname = TString::Format("%s/hPhosM02VsEnergy", cont->GetArrayName().Data());
236  title = histname + ";M02;#it{E}_{clus} (GeV/);counts";
237  fHistManager.CreateTH2(histname.Data(), title.Data(), 100, 0, 20, nPtBins, 0, fMaxPt);
238 
239  histname = TString::Format("%s/hPhosCellIdVsEnergy", cont->GetArrayName().Data());
240  title = histname + ";Cell ID;#it{E}_{clus} (GeV/);counts";
241  fHistManager.CreateTH2(histname.Data(), title.Data(), 20000, 0, 20000, nPtBins, 0, fMaxPt);
242 
243  const Int_t nEmcalSM = 20;
244  for (Int_t sm = 0; sm < nEmcalSM; sm++) {
245  histname = TString::Format("%s/BySM/hEmcalClusEnergy_SM%d", cont->GetArrayName().Data(), sm);
246  title = histname + ";#it{E}_{cluster} (GeV);counts";
247  fHistManager.CreateTH1(histname.Data(), title.Data(), nPtBins, 0, fMaxPt);
248  }
249 
250  for (Int_t sm = 1; sm < 5; sm++) {
251  histname = TString::Format("%s/BySM/hPhosClusEnergy_SM%d", cont->GetArrayName().Data(), sm);
252  title = histname + ";#it{E}_{cluster} (GeV);counts";
253  fHistManager.CreateTH1(histname.Data(), title.Data(), nPtBins, 0, fMaxPt);
254  }
255 
256  // Cluster THnSparse
257  Int_t dim = 0;
258  TString title[20];
259  Int_t nbins[20] = {0};
260  Double_t min[30] = {0.};
261  Double_t max[30] = {0.};
262  Double_t *binEdges[20] = {0};
263 
265  title[dim] = "Centrality %";
266  nbins[dim] = fNCentHistBins;
267  binEdges[dim] = fCentHistBins;
268  min[dim] = fCentHistBins[0];
269  max[dim] = fCentHistBins[fNCentHistBins];
270  dim++;
271  }
272 
273  title[dim] = "#it{E}_{clus} (GeV)";
274  nbins[dim] = fNPtHistBins;
275  binEdges[dim] = fPtHistBins;
276  min[dim] = fPtHistBins[0];
277  max[dim] = fPtHistBins[fNPtHistBins];
278  dim++;
279 
280  title[dim] = "#eta";
281  nbins[dim] = fNEtaHistBins;
282  binEdges[dim] = fEtaHistBins;
283  min[dim] = fEtaHistBins[0];
284  max[dim] = fEtaHistBins[fNEtaHistBins];
285  dim++;
286 
287  title[dim] = "#phi";
288  nbins[dim] = fNPhiHistBins;
289  binEdges[dim] = fPhiHistBins;
290  min[dim] = fPhiHistBins[0];
291  max[dim] = fPhiHistBins[fNPhiHistBins];
292  dim++;
293 
294  title[dim] = "cluster type";
295  nbins[dim] = 3;
296  binEdges[dim] = fIntegerHistBins;
297  min[dim] = fIntegerHistBins[0];
298  max[dim] = fIntegerHistBins[3];
299  dim++;
300 
301  TString thnname = TString::Format("%s/clusterObservables", cont->GetArrayName().Data());
302  THnSparse* hn = fHistManager.CreateTHnSparse(thnname.Data(), thnname.Data(), dim, nbins, min, max);
303  for (Int_t i = 0; i < dim; i++) {
304  hn->GetAxis(i)->SetTitle(title[i]);
305  hn->SetBinEdges(i, binEdges[i]);
306  }
307  }
308 }
309 
310 //________________________________________________________________________
312 
313  AliJetContainer* jets = 0;
314  TIter nextJetColl(&fJetCollArray);
315  while ((jets = static_cast<AliJetContainer*>(nextJetColl()))) {
316 
317  // Allocate THnSparse
318 
319  TString axisTitle[30]= {""};
320  Int_t nbins[30] = {0};
321  Double_t min[30] = {0.};
322  Double_t max[30] = {0.};
323  Double_t *binEdges[20] = {0};
324  Int_t dim = 0;
325 
326  if (fForceBeamType != kpp) {
327  axisTitle[dim] = "Centrality (%)";
328  nbins[dim] = fNCentHistBins;
329  binEdges[dim] = fCentHistBins;
330  min[dim] = fCentHistBins[0];
331  max[dim] = fCentHistBins[fNCentHistBins];
332  dim++;
333  }
334 
335  axisTitle[dim] = "#eta_{jet}";
336  nbins[dim] = fNEtaHistBins;
337  binEdges[dim] = fEtaHistBins;
338  min[dim] = fEtaHistBins[0];
339  max[dim] = fEtaHistBins[fNEtaHistBins];
340  dim++;
341 
342  axisTitle[dim] = "#phi_{jet} (rad)";
343  nbins[dim] = fNPhiHistBins;
344  binEdges[dim] = fPhiHistBins;
345  min[dim] = fPhiHistBins[0];
346  max[dim] = fPhiHistBins[fNPhiHistBins];
347  dim++;
348 
349  axisTitle[dim] = "#it{p}_{T} (GeV/#it{c})";
350  nbins[dim] = TMath::CeilNint(fMaxPt/3);
351  min[dim] = 0;
352  max[dim] = fMaxPt;
353  binEdges[dim] = GenerateFixedBinArray(nbins[dim], min[dim], max[dim]);
354  dim++;
355 
356  if (fForceBeamType != kpp) {
357  axisTitle[dim] = "#it{p}_{T}^{corr} (GeV/#it{c})";
358  nbins[dim] = TMath::CeilNint(fMaxPt/3);
359  min[dim] = -fMaxPt/2 + 25;
360  max[dim] = fMaxPt/2 + 25;
361  binEdges[dim] = GenerateFixedBinArray(nbins[dim], min[dim], max[dim]);
362  dim++;
363  }
364 
365  axisTitle[dim] = "#it{p}_{T,particle}^{leading} (GeV/#it{c})";
366  nbins[dim] = fMaxPt/5;
367  min[dim] = 0;
368  max[dim] = 150;
369  binEdges[dim] = GenerateFixedBinArray(nbins[dim], min[dim], max[dim]);
370  dim++;
371 
372  TString thnname = TString::Format("%s/fHistJetObservables", jets->GetArrayName().Data());
373  THnSparse* hn = fHistManager.CreateTHnSparse(thnname.Data(), thnname.Data(), dim, nbins, min, max);
374  for (Int_t i = 0; i < dim; i++) {
375  hn->GetAxis(i)->SetTitle(axisTitle[i]);
376  hn->SetBinEdges(i, binEdges[i]);
377  }
378 
379  // Allocate other jet histograms
380  TString histname;
381  TString title;
382  histname = TString::Format("%s/fHistJetRejectionReason", jets->GetArrayName().Data());
383  title = histname + ";Rejection reason;#it{p}_{T,jet} (GeV/#it{c});counts";
384  TH2* hist = fHistManager.CreateTH2(histname.Data(), title.Data(), 32, 0, 32, 50, 0, 250);
385  SetRejectionReasonLabels(hist->GetXaxis());
386 
387  if (!jets->GetRhoName().IsNull()) {
388  histname = TString::Format("%s/fHistRhoVsCent", jets->GetArrayName().Data());
389  title = histname + ";Centrality (%);#rho (GeV/#it{c});counts";
390  fHistManager.CreateTH2(histname.Data(), title.Data(), 101, 0, 101, 100, 0, 500);
391  }
392  }
393 }
394 
395 //________________________________________________________________________
397 
398  Int_t nPtBins = 2*fMaxPt;
399 
400  Int_t dim = 0;
401  TString axistitle[40];
402  Int_t nbins[40] = {0};
403  Double_t min[40] = {0};
404  Double_t max[40] = {0};
405  Double_t *binEdges[20] = {0};
406 
408  axistitle[dim] = "Centrality %";
409  nbins[dim] = fNCentHistBins;
410  binEdges[dim] = fCentHistBins;
411  min[dim] = fCentHistBins[0];
412  max[dim] = fCentHistBins[fNCentHistBins];
413  dim++;
414  }
415 
416  if (fParticleCollArray.GetEntriesFast()>0) {
417  axistitle[dim] = "No. of tracks";
419  nbins[dim] = 50;
420  min[dim] = -0.5;
421  max[dim] = 10000-0.5;
422  binEdges[dim] = GenerateFixedBinArray(nbins[dim], min[dim], max[dim]);
423  }
424  else {
425  nbins[dim] = 50;
426  min[dim] = 0;
427  max[dim] = 200;
428  binEdges[dim] = GenerateFixedBinArray(nbins[dim], min[dim], max[dim]);
429  }
430  dim++;
431 
432  axistitle[dim] = "#it{p}_{T,track}^{leading} (GeV/c)";
433  nbins[dim] = fMaxPt/5;
434  min[dim] = 0;
435  max[dim] = 150;
436  binEdges[dim] = GenerateFixedBinArray(nbins[dim], min[dim], max[dim]);
437  dim++;
438  }
439 
440  if (fClusterCollArray.GetEntriesFast()>0) {
441  axistitle[dim] = "No. of clusters";
442 
444  nbins[dim] = 50;
445  min[dim] = -0.5;
446  max[dim] = 4000-0.5;
447  binEdges[dim] = GenerateFixedBinArray(nbins[dim], min[dim], max[dim]);
448  }
449  else {
450  nbins[dim] = 50;
451  min[dim] = 0;
452  max[dim] = 200;
453  binEdges[dim] = GenerateFixedBinArray(nbins[dim], min[dim], max[dim]);
454  }
455  dim++;
456 
457  axistitle[dim] = "#it{E}_{cluster}^{leading} (GeV)";
458  nbins[dim] = fMaxPt/5;
459  min[dim] = 0;
460  max[dim] = 150;
461  binEdges[dim] = GenerateFixedBinArray(nbins[dim], min[dim], max[dim]);
462  dim++;
463  }
464 
465  THnSparse* hn = fHistManager.CreateTHnSparse("eventQA","eventQA",dim,nbins,min,max);
466  for (Int_t i = 0; i < dim; i++) {
467  hn->GetAxis(i)->SetTitle(axistitle[i]);
468  hn->SetBinEdges(i, binEdges[i]);
469  }
470 
471  if (fIsPtHard) {
472  TString histname = "hPtHard";
473  TString title = histname + ";#it{p}_{T,hard} (GeV/c);counts";
474  fHistManager.CreateTH1(histname.Data(), title.Data(), nPtBins, 0, fMaxPt);
475 
476  TH1* hTrials = fHistManager.CreateTH1("hNtrials", "hNtrials", 1, 0, 1);
477  hTrials->GetXaxis()->SetBinLabel(1,"#sum{ntrials}");
478 
479  TH1* hxsection = fHistManager.CreateTH1("hXsec", "hXsec", 1, 0, 1);
480  hxsection->GetXaxis()->SetBinLabel(1,"<#sigma>");
481  }
482 
483  fHistEventRejection->GetXaxis()->SetBinLabel(15,"PtHardBinJetOutlier");
484 
485 }
486 
487 //________________________________________________________________________
489 {
490  Int_t dim = 0;
491  TString title[20];
492  Int_t nbins[20] = {0};
493  Double_t min[30] = {0.};
494  Double_t max[30] = {0.};
495  Double_t *binEdges[20] = {0};
496 
498  title[dim] = "Centrality %";
499  nbins[dim] = fNCentHistBins;
500  binEdges[dim] = fCentHistBins;
501  min[dim] = fCentHistBins[0];
502  max[dim] = fCentHistBins[fNCentHistBins];
503  dim++;
504  }
505 
506  title[dim] = "#it{p}_{T} (GeV/#it{c})";
507  nbins[dim] = fNPtHistBins;
508  binEdges[dim] = fPtHistBins;
509  min[dim] = fPtHistBins[0];
510  max[dim] = fPtHistBins[fNPtHistBins];
511  dim++;
512 
513  title[dim] = "#eta";
514  nbins[dim] = fNEtaHistBins;
515  binEdges[dim] = fEtaHistBins;
516  min[dim] = fEtaHistBins[0];
517  max[dim] = fEtaHistBins[fNEtaHistBins];
518  dim++;
519 
520  title[dim] = "#phi";
521  nbins[dim] = fNPhiHistBins;
522  binEdges[dim] = fPhiHistBins;
523  min[dim] = fPhiHistBins[0];
524  max[dim] = fPhiHistBins[fNPhiHistBins];
525  dim++;
526 
527  title[dim] = "track type";
528  nbins[dim] = 3;
529  binEdges[dim] = fIntegerHistBins;
530  min[dim] = fIntegerHistBins[0];
531  max[dim] = fIntegerHistBins[3];
532  dim++;
533 
534  title[dim] = "#sigma(#it{p}_{T}) / #it{p}_{T}";
535  nbins[dim] = fNPtResHistBins;
536  binEdges[dim] = fPtResHistBins;
537  min[dim] = fPtResHistBins[0];
538  max[dim] = fPtResHistBins[fNPtResHistBins];
539  dim++;
540 
541  THnSparse* hn = fHistManager.CreateTHnSparse("tracks","tracks", dim, nbins, min, max);
542  for (Int_t i = 0; i < dim; i++) {
543  hn->GetAxis(i)->SetTitle(title[i]);
544  hn->SetBinEdges(i, binEdges[i]);
545  }
546 
547 }
548 
549 //________________________________________________________________________
551 {
552  Int_t dim = 0;
553  TString title[20];
554  Int_t nbins[20] = {0};
555  Double_t min[30] = {0.};
556  Double_t max[30] = {0.};
557  Double_t *binEdges[20] = {0};
558 
560  title[dim] = "Centrality %";
561  nbins[dim] = fNCentHistBins;
562  binEdges[dim] = fCentHistBins;
563  min[dim] = fCentHistBins[0];
564  max[dim] = fCentHistBins[fNCentHistBins];
565  dim++;
566  }
567 
568  title[dim] = "#it{p}_{T} (GeV/#it{c})";
569  nbins[dim] = fNPtHistBins;
570  binEdges[dim] = fPtHistBins;
571  min[dim] = fPtHistBins[0];
572  max[dim] = fPtHistBins[fNPtHistBins];
573  dim++;
574 
575  title[dim] = "#eta";
576  nbins[dim] = fNEtaHistBins;
577  binEdges[dim] = fEtaHistBins;
578  min[dim] = fEtaHistBins[0];
579  max[dim] = fEtaHistBins[fNEtaHistBins];
580  dim++;
581 
582  title[dim] = "#phi";
583  nbins[dim] = fNPhiHistBins;
584  binEdges[dim] = fPhiHistBins;
585  min[dim] = fPhiHistBins[0];
586  max[dim] = fPhiHistBins[fNPhiHistBins];
587  dim++;
588 
589  title[dim] = "Findable";
590  nbins[dim] = 2;
591  binEdges[dim] = fIntegerHistBins;
592  min[dim] = fIntegerHistBins[0];
593  max[dim] = fIntegerHistBins[2];
594  dim++;
595 
596  THnSparse* hn = fHistManager.CreateTHnSparse("tracks_PhysPrim","tracks_PhysPrim", dim, nbins, min, max);
597  for (Int_t i = 0; i < dim; i++) {
598  hn->GetAxis(i)->SetTitle(title[i]);
599  hn->SetBinEdges(i, binEdges[i]);
600  }
601 }
602 
603 //________________________________________________________________________
605 {
606  Int_t dim = 0;
607  TString title[20];
608  Int_t nbins[20] = {0};
609  Double_t min[30] = {0.};
610  Double_t max[30] = {0.};
611  Double_t *binEdges[20] = {0};
612 
614  title[dim] = "Centrality %";
615  nbins[dim] = fNCentHistBins;
616  binEdges[dim] = fCentHistBins;
617  min[dim] = fCentHistBins[0];
618  max[dim] = fCentHistBins[fNCentHistBins];
619  dim++;
620  }
621 
622  title[dim] = "#it{p}_{T}^{gen} (GeV/#it{c})";
623  nbins[dim] = fNPtHistBins;
624  binEdges[dim] = fPtHistBins;
625  min[dim] = fPtHistBins[0];
626  max[dim] = fPtHistBins[fNPtHistBins];
627  dim++;
628 
629  title[dim] = "#eta^{gen}";
630  nbins[dim] = fNEtaHistBins;
631  binEdges[dim] = fEtaHistBins;
632  min[dim] = fEtaHistBins[0];
633  max[dim] = fEtaHistBins[fNEtaHistBins];
634  dim++;
635 
636  title[dim] = "#phi^{gen}";
637  nbins[dim] = fNPhiHistBins;
638  binEdges[dim] = fPhiHistBins;
639  min[dim] = fPhiHistBins[0];
640  max[dim] = fPhiHistBins[fNPhiHistBins];
641  dim++;
642 
643  title[dim] = "#it{p}_{T}^{det} (GeV/#it{c})";
644  nbins[dim] = fNPtHistBins;
645  binEdges[dim] = fPtHistBins;
646  min[dim] = fPtHistBins[0];
647  max[dim] = fPtHistBins[fNPtHistBins];
648  dim++;
649 
650  title[dim] = "#eta^{det}";
651  nbins[dim] = fNEtaHistBins;
652  binEdges[dim] = fEtaHistBins;
653  min[dim] = fEtaHistBins[0];
654  max[dim] = fEtaHistBins[fNEtaHistBins];
655  dim++;
656 
657  title[dim] = "#phi^{det}";
658  nbins[dim] = fNPhiHistBins;
659  binEdges[dim] = fPhiHistBins;
660  min[dim] = fPhiHistBins[0];
661  max[dim] = fPhiHistBins[fNPhiHistBins];
662  dim++;
663 
664  title[dim] = "(#it{p}_{T}^{gen} - #it{p}_{T}^{det}) / #it{p}_{T}^{gen}";
665  nbins[dim] = fNPtRelDiffHistBins;
666  binEdges[dim] = fPtRelDiffHistBins;
667  min[dim] = fPtRelDiffHistBins[0];
669  dim++;
670 
671  title[dim] = "track type";
672  nbins[dim] = 3;
673  binEdges[dim] = fIntegerHistBins;
674  min[dim] = fIntegerHistBins[0];
675  max[dim] = fIntegerHistBins[3];
676  dim++;
677 
678  THnSparse* hn = fHistManager.CreateTHnSparse("tracks_Matched","tracks_Matched", dim, nbins, min, max);
679  for (Int_t i = 0; i < dim; i++) {
680  hn->GetAxis(i)->SetTitle(title[i]);
681  hn->SetBinEdges(i, binEdges[i]);
682  }
683 }
684 
685 
686 //________________________________________________________________________
689 {
690  fNPtHistBins = 82;
693  GenerateFixedBinArray(7, 0.3, 1, fPtHistBins+6);
694  GenerateFixedBinArray(10, 1, 3, fPtHistBins+13);
695  GenerateFixedBinArray(14, 3, 10, fPtHistBins+23);
696  GenerateFixedBinArray(10, 10, 20, fPtHistBins+37);
697  GenerateFixedBinArray(15, 20, 50, fPtHistBins+47);
698  GenerateFixedBinArray(20, 50, 150, fPtHistBins+62);
699 
700  fNEtaHistBins = 70;
703 
704  fNPhiHistBins = 101;
706  GenerateFixedBinArray(fNPhiHistBins, 0, TMath::Pi() * 2.02, fPhiHistBins);
707 
708  fNCentHistBins = 4;
710  fCentHistBins[0] = 0;
711  fCentHistBins[1] = 10;
712  fCentHistBins[2] = 30;
713  fCentHistBins[3] = 50;
714  fCentHistBins[4] = 90;
715 
716  fNPtResHistBins = 100;
719  GenerateFixedBinArray(25, 0.05, 0.10, fPtResHistBins+50);
720  GenerateFixedBinArray(15, 0.10, 0.20, fPtResHistBins+75);
721  GenerateFixedBinArray(10, 0.20, 0.50, fPtResHistBins+90);
722 
723  fNPtRelDiffHistBins = 50;
726 
727  fNIntegerHistBins = 3;
730 }
731 
732 //________________________________________________________________________
734 {
735  if (fClusterCollArray.GetEntriesFast() == 0 && fCaloCellsName.IsNull()) {
736  fNeedEmcalGeom = kFALSE;
737  }
738  else {
739  fNeedEmcalGeom = kTRUE;
740  }
741 
743 
744  // Load the PHOS geometry
745  fPHOSGeo = AliPHOSGeometry::GetInstance();
746  if (fPHOSGeo) {
747  AliInfo("Found instance of PHOS geometry!");
748  }
749  else {
750  AliInfo("Creating PHOS geometry!");
751  Int_t runNum = InputEvent()->GetRunNumber();
752  if(runNum<209122) //Run1
753  fPHOSGeo = AliPHOSGeometry::GetInstance("IHEP");
754  else
755  fPHOSGeo = AliPHOSGeometry::GetInstance("Run2");
756 
757  if (fPHOSGeo) {
758  AliOADBContainer geomContainer("phosGeo");
759  geomContainer.InitFromFile("$ALICE_PHYSICS/OADB/PHOS/PHOSMCGeometry.root","PHOSMCRotationMatrixes");
760  TObjArray* matrixes = (TObjArray*)geomContainer.GetObject(runNum,"PHOSRotationMatrixes");
761  for(Int_t mod=0; mod<6; mod++) {
762  if(!matrixes->At(mod)) continue;
763  fPHOSGeo->SetMisalMatrix(((TGeoHMatrix*)matrixes->At(mod)),mod);
764  printf(".........Adding Matrix(%d), geo=%p\n",mod,fPHOSGeo);
765  ((TGeoHMatrix*)matrixes->At(mod))->Print();
766  }
767  }
768  }
769 
770 }
771 
772 //________________________________________________________________________
774 {
775  if (!AliAnalysisTaskEmcalJet::RetrieveEventObjects()) return kFALSE;
776 
777  // If Pt-hard production, get the Pt-hard of the event, and have possibility to reject the event for jet outliers
778  if (fIsPtHard) {
779  AliGenPythiaEventHeader* pygen = nullptr;
780  if (MCEvent()) {
781  pygen = dynamic_cast<AliGenPythiaEventHeader*>(MCEvent()->GenEventHeader());
782  if (!pygen) {
783  // Check if AOD
784  AliAODMCHeader* aodMCH = dynamic_cast<AliAODMCHeader*>(InputEvent()->FindListObject(AliAODMCHeader::StdBranchName()));
785  if (aodMCH) {
786  for (UInt_t i = 0;i<aodMCH->GetNCocktailHeaders();i++) {
787  pygen = dynamic_cast<AliGenPythiaEventHeader*>(aodMCH->GetCocktailHeader(i));
788  if (pygen) break;
789  }
790  }
791  }
792  }
793  if (pygen) {
794  fPtHard = pygen->GetPtHard();
795  //fNTrials = pygen->Trials();
796  //fXsection = pygen->GetXsection();
797 
798  // reject outlier events, where the jet Pt is much larger than Pt-hard
799  if (fRejectOutlierEvents) {
800  AliTLorentzVector jet;
801  Int_t nTriggerJets = pygen->NTriggerJets();
802  Float_t tmpjet[]={0,0,0,0};
803  for (Int_t ijet = 0; ijet< nTriggerJets; ijet++) {
804  pygen->TriggerJet(ijet, tmpjet);
805  jet.SetPxPyPzE(tmpjet[0],tmpjet[1],tmpjet[2],tmpjet[3]);
806 
807  if (jet.Pt() > 4. * fPtHard) {
808  //AliInfo(Form("Reject jet event with: pT Hard %2.2f, pycell jet pT %2.2f, rejection factor %1.1f\n", fPtHard, jet.Pt(), 4.));
809  if (fGeneralHistograms) fHistEventRejection->Fill("PtHardBinJetOutlier",1);
810  return kFALSE;
811  }
812  }
813  }
814  }
815  }
816 
817  return kTRUE;
818 }
819 
827 //________________________________________________________________________
829 {
830  if (fIsPtHard) {
831  TTree *tree = AliAnalysisManager::GetAnalysisManager()->GetTree();
832  if (!tree) {
833  AliError(Form("%s - UserNotify: No current tree!",GetName()));
834  return kFALSE;
835  }
836 
837  Float_t xsection = 0;
838  Float_t trials = 0;
839  Int_t pthardbin = 0;
840 
841  TFile *curfile = tree->GetCurrentFile();
842  if (!curfile) {
843  AliError(Form("%s - UserNotify: No current file!",GetName()));
844  return kFALSE;
845  }
846 
847  TChain *chain = dynamic_cast<TChain*>(tree);
848  if (chain) tree = chain->GetTree();
849 
850  PythiaInfoFromFile(curfile->GetName(), xsection, trials, pthardbin);
851 
852  fHistManager.FillTH1("hNtrials", "#sum{ntrials}", trials);
853  fHistManager.FillTH1("hXsec", "<#sigma>", xsection);
854  }
855  return kTRUE;
856 }
857 
858 //________________________________________________________________________
860 {
866 
867  return kTRUE;
868 }
869 
870 //________________________________________________________________________
872 
873  fNTotTracks = 0;
874  fLeadingTrack.SetPxPyPzE(0,0,0,0);
875 
876  AliVTrack* track;
877  for (auto trackIterator : fDetectorLevel->accepted_momentum() ) {
878 
879  fNTotTracks++;
880  if (fLeadingTrack.Pt() < trackIterator.first.Pt()) fLeadingTrack = trackIterator.first;
881 
882  track = trackIterator.second;
883  Byte_t type = fDetectorLevel->GetTrackType(track);
884  if (type <= 2) {
885  Double_t sigma = 0;
886 
887  if (fIsEsd) {
888 
889  AliESDtrack *esdTrack = dynamic_cast<AliESDtrack*>(track);
890  if (esdTrack) sigma = TMath::Sqrt(esdTrack->GetSigma1Pt2());
891 
892  } else { // AOD
893 
894  AliAODTrack *aodtrack = dynamic_cast<AliAODTrack*>(track);
895  if(!aodtrack) AliFatal("Not a standard AOD");
896 
897  AliExternalTrackParam exParam;
898 
899  //get covariance matrix
900  Double_t cov[21] = {0,};
901  aodtrack->GetCovMatrix(cov);
902  Double_t pxpypz[3] = {0,};
903  aodtrack->PxPyPz(pxpypz);
904  Double_t xyz[3] = {0,};
905  aodtrack->GetXYZ(xyz);
906  Short_t sign = aodtrack->Charge();
907  exParam.Set(xyz,pxpypz,cov,sign);
908 
909  sigma = TMath::Sqrt(exParam.GetSigma1Pt2());
910 
911  }
912 
913  Int_t label = TMath::Abs(track->GetLabel());
914 
915  FillDetectorLevelTHnSparse(fCent, track->Eta(), track->Phi(), track->Pt(), sigma, type);
916 
917  if (fGeneratorLevel && label > 0) {
918  AliAODMCParticle *part = fGeneratorLevel->GetAcceptMCParticleWithLabel(label);
919  if (part) {
920  if (part->GetGeneratorIndex() == 0) {
921  Int_t pdg = TMath::Abs(part->PdgCode());
922  // select charged pions, protons, kaons, electrons, muons
923  if (pdg == 211 || pdg == 2212 || pdg == 321 || pdg == 11 || pdg == 13) {
924  FillMatchedParticlesTHnSparse(fCent, part->Eta(), part->Phi(), part->Pt(), track->Eta(), track->Phi(), track->Pt(), type);
925  }
926  }
927  }
928  }
929  }
930  else {
931  AliError(Form("Track %d has type %d not recognized!", fDetectorLevel->GetCurrentID(), type));
932  }
933  }
934 
935  if (fGeneratorLevel) {
936  AliAODMCParticle* part;
937  for (auto partIterator : fGeneratorLevel->accepted_momentum() ) {
938  part = partIterator.second;
939 
940  Byte_t findable = 0;
941 
942  Int_t pdg = TMath::Abs(part->PdgCode());
943  // select charged pions, protons, kaons, electrons, muons
944  if (pdg == 211 || pdg == 2212 || pdg == 321 || pdg == 11 || pdg == 13) findable = 1;
945 
946  FillGeneratorLevelTHnSparse(fCent, part->Eta(), part->Phi(), part->Pt(), findable);
947  }
948  }
949 }
950 
951 //________________________________________________________________________
953 
954  if (!fCaloCells) return;
955 
956  TString histname_energy = TString::Format("%s/fHistCellEnergy", fCaloCellsName.Data());
957  TString histname_time = TString::Format("%s/fHistCellTime", fCaloCellsName.Data());
958  TString histname_energyProf = TString::Format("%s/fProfCellAbsIdEnergy", fCaloCellsName.Data());
959  TString histname_timeProf = TString::Format("%s/fProfCellAbsIdTime", fCaloCellsName.Data());
960  TString histname_EnergyvsTime = TString::Format("%s/fHistCellEvsTime", fCaloCellsName.Data());
961 
962  const Int_t ncells = fCaloCells->GetNumberOfCells();
963  for (Int_t pos = 0; pos < ncells; pos++) {
964  Float_t amp = fCaloCells->GetAmplitude(pos);
965  Float_t time = fCaloCells->GetTime(pos);
966  Int_t absId = fCaloCells->GetCellNumber(pos);
967 
968  if (amp < fCellEnergyCut) continue;
969 
970  fHistManager.FillTH1(histname_energy, amp);
971  fHistManager.FillTH1(histname_time, time);
972 
973  fHistManager.FillProfile(histname_energyProf, absId, amp);
974  fHistManager.FillProfile(histname_timeProf, absId, time);
975 
976  fHistManager.FillTH2(histname_EnergyvsTime, time, amp);
977  }
978 }
979 
980 //________________________________________________________________________
982 
983  memset(fNTotClusters, 0, sizeof(Int_t)*3);
984  for (Int_t i = 0; i < 3; i++) fLeadingCluster[i].SetPxPyPzE(0,0,0,0);
985 
986  TString histname;
987  AliClusterContainer* clusters = 0;
988  TIter nextClusColl(&fClusterCollArray);
989  while ((clusters = static_cast<AliClusterContainer*>(nextClusColl()))) {
990  // Cluster loop
992  for (AliClusterIterableMomentumContainer::iterator it = itcont.begin(); it != itcont.end(); it++) {
993 
994  UInt_t rejectionReason = 0;
995  if (!clusters->AcceptCluster(it.current_index(), rejectionReason)) {
996  histname = TString::Format("%s/fHistClusterRejectionReason", clusters->GetArrayName().Data());
997  fHistManager.FillTH2(histname, clusters->GetRejectionReasonBitPosition(rejectionReason), it->first.E());
998  continue;
999  }
1000 
1001  // Determine cluster type (EMCal/DCal/Phos) and record relevant eventQA info
1002  ClusterType clusType = kNA;
1003  if (it->second->IsEMCAL()) {
1004  Double_t phi = it->first.Phi_0_2pi();
1005  Int_t isDcal = Int_t(phi > fgkEMCalDCalPhiDivide);
1006  if (isDcal == 0) {
1007  clusType = kEMCal;
1008  } else if (isDcal == 1) {
1009  clusType = kDCal;
1010  }
1011 
1012  if (fLeadingCluster[isDcal].E() < it->first.E()) fLeadingCluster[isDcal] = it->first;
1013  fNTotClusters[isDcal]++;
1014 
1015  Int_t sm = fGeom->GetSuperModuleNumber(it->second->GetCellAbsId(0));
1016  if (sm >=0 && sm < 20) {
1017  histname = TString::Format("%s/BySM/hEmcalClusEnergy_SM%d", clusters->GetArrayName().Data(), sm);
1018  fHistManager.FillTH1(histname, it->second->E());
1019  }
1020  else {
1021  AliError(Form("Supermodule %d does not exist!", sm));
1022  }
1023 
1024  } else if (it->second->GetType() == AliVCluster::kPHOSNeutral){
1025 
1026  Int_t nCells = it->second->GetNCells();
1027  Double_t M02 = it->second->GetM02();
1028  Double_t energy = it->second->E();
1029 
1030  histname = TString::Format("%s/hPhosNCellsVsEnergy", clusters->GetArrayName().Data());
1031  fHistManager.FillTH2(histname, nCells, energy);
1032 
1033  histname = TString::Format("%s/hPhosM02VsEnergy", clusters->GetArrayName().Data());
1034  fHistManager.FillTH2(histname, M02, energy);
1035 
1036  clusType = kPHOS;
1037  fNTotClusters[2]++;
1038  if (fLeadingCluster[2].E() < it->first.E()) fLeadingCluster[2] = it->first;
1039 
1040  // Fill Phos spectra by module
1041  Int_t relid[4];
1042  if (fPHOSGeo) {
1043  fPHOSGeo->AbsToRelNumbering(it->second->GetCellAbsId(0), relid);
1044  Int_t sm = relid[0];
1045  if (sm >=1 && sm < 5) {
1046  histname = TString::Format("%s/BySM/hPhosClusEnergy_SM%d", clusters->GetArrayName().Data(), sm);
1047  fHistManager.FillTH1(histname, energy);
1048  }
1049  else {
1050  AliError(Form("Supermodule %d does not exist!", sm));
1051  }
1052  }
1053 
1054  // Loop through cells in the cluster, and histogram their energy
1055  histname = TString::Format("%s/hPhosCellIdVsEnergy", clusters->GetArrayName().Data());
1056  for (Int_t i=0; i < nCells; i++) {
1057  fHistManager.FillTH2(histname, it->second->GetCellAbsId(i), energy);
1058  }
1059 
1060  }
1061 
1062  Double_t contents[30]={0};
1063  histname = TString::Format("%s/clusterObservables", clusters->GetArrayName().Data());
1064  THnSparse* histClusterObservables = static_cast<THnSparse*>(fHistManager.FindObject(histname));
1065  if (!histClusterObservables) return;
1066  for (Int_t i = 0; i < histClusterObservables->GetNdimensions(); i++) {
1067  TString title(histClusterObservables->GetAxis(i)->GetTitle());
1068  if (title=="Centrality %")
1069  contents[i] = fCent;
1070  else if (title=="#it{E}_{clus} (GeV)")
1071  contents[i] = it->first.E();
1072  else if (title=="#eta")
1073  contents[i] = it->first.Eta();
1074  else if (title=="#phi")
1075  contents[i] = it->first.Phi_0_2pi();
1076  else if (title=="cluster type")
1077  contents[i] = clusType;
1078  else
1079  AliWarning(Form("Unable to fill dimension %s!",title.Data()));
1080  }
1081  histClusterObservables->Fill(contents);
1082 
1083  }
1084  }
1085 }
1086 
1087 //________________________________________________________________________
1089 
1090  TString histname;
1091  AliJetContainer* jets = 0;
1092  TIter nextJetColl(&fJetCollArray);
1093  while ((jets = static_cast<AliJetContainer*>(nextJetColl()))) {
1094  Double_t rhoVal = 0;
1095  if (jets->GetRhoParameter()) {
1096  rhoVal = jets->GetRhoVal();
1097  histname = TString::Format("%s/fHistRhoVsCent", jets->GetArrayName().Data());
1098  fHistManager.FillTH2(histname.Data(), fCent, rhoVal);
1099  }
1100 
1101  for (auto jet : jets->all()) {
1102 
1103  UInt_t rejectionReason = 0;
1104  if (!jets->AcceptJet(jet, rejectionReason)) {
1105  histname = TString::Format("%s/fHistJetRejectionReason", jets->GetArrayName().Data());
1106  fHistManager.FillTH2(histname.Data(), jets->GetRejectionReasonBitPosition(rejectionReason), jet->Pt());
1107  continue;
1108  }
1109 
1110  Float_t ptLeading = jets->GetLeadingHadronPt(jet);
1111  Float_t corrPt = jet->Pt() - rhoVal * jet->Area();
1112 
1113  TLorentzVector leadPart;
1114  jets->GetLeadingHadronMomentum(leadPart, jet);
1115 
1116  Double_t contents[30]={0};
1117  histname = TString::Format("%s/fHistJetObservables", jets->GetArrayName().Data());
1118  THnSparse* histJetObservables = static_cast<THnSparse*>(fHistManager.FindObject(histname));
1119  if (!histJetObservables) return;
1120  for (Int_t i = 0; i < histJetObservables->GetNdimensions(); i++) {
1121  TString title(histJetObservables->GetAxis(i)->GetTitle());
1122  if (title=="Centrality (%)")
1123  contents[i] = fCent;
1124  else if (title=="#eta_{jet}")
1125  contents[i] = jet->Eta();
1126  else if (title=="#phi_{jet} (rad)")
1127  contents[i] = jet->Phi_0_2pi();
1128  else if (title=="#it{p}_{T} (GeV/#it{c})")
1129  contents[i] = jet->Pt();
1130  else if (title=="#it{p}_{T}^{MC} (GeV/#it{c})")
1131  contents[i] = jet->MCPt();
1132  else if (title=="#it{p}_{T}^{corr} (GeV/#it{c})")
1133  contents[i] = corrPt;
1134  else if (title=="#it{p}_{T,particle}^{leading} (GeV/#it{c})")
1135  contents[i] = ptLeading;
1136  else
1137  AliWarning(Form("Unable to fill dimension %s!",title.Data()));
1138  }
1139  histJetObservables->Fill(contents);
1140 
1141  } //jet loop
1142  }
1143 }
1144 
1145 //________________________________________________________________________
1147 
1148  EventQA_t eventQA;
1149  for (Int_t i = 0; i < 3; i++) {
1150  eventQA.fMaxCluster[i] = fLeadingCluster[i];
1151  }
1152  eventQA.fCent = fCent;
1153  eventQA.fNTracks = fNTotTracks;
1154  eventQA.fNClusters[0] = fNTotClusters[0];
1155  eventQA.fNClusters[1] = fNTotClusters[1];
1156  eventQA.fNClusters[2] = fNTotClusters[2];
1157  eventQA.fMaxTrack = fLeadingTrack;
1158 
1159  Int_t globalNclusters = eventQA.fNClusters[0] + eventQA.fNClusters[1] + eventQA.fNClusters[2];
1160  AliTLorentzVector globalMaxCluster;
1161  for (Int_t i = 0; i < 3; i++) {
1162  if (globalMaxCluster.E() < eventQA.fMaxCluster[i].E()) globalMaxCluster = eventQA.fMaxCluster[i];
1163  }
1164 
1165  THnSparse* histEventQA = static_cast<THnSparse*>(fHistManager.FindObject("eventQA"));
1166  Double_t contents[40]={0};
1167  for (Int_t i = 0; i < histEventQA->GetNdimensions(); i++) {
1168  TString title(histEventQA->GetAxis(i)->GetTitle());
1169  if (title=="Centrality %")
1170  contents[i] = eventQA.fCent;
1171  else if (title=="No. of tracks")
1172  contents[i] = eventQA.fNTracks;
1173  else if (title=="No. of clusters")
1174  contents[i] = globalNclusters;
1175  else if (title=="#it{p}_{T,track}^{leading} (GeV/c)")
1176  contents[i] = eventQA.fMaxTrack.Pt();
1177  else if (title=="#it{E}_{cluster}^{leading} (GeV)")
1178  contents[i] = globalMaxCluster.E();
1179  else
1180  AliWarning(Form("Unable to fill dimension %s!",title.Data()));
1181  }
1182 
1183  histEventQA->Fill(contents);
1184 
1185  // Fill pythia pt hard histograms, if applicable
1186  if (fPtHard > 1e-6) {
1187  fHistManager.FillTH1("hPtHard", fPtHard);
1188  }
1189 }
1190 
1191 //________________________________________________________________________
1193  Double_t sigma1OverPt, Byte_t trackType)
1194 {
1195  Double_t contents[20]={0};
1196 
1197  THnSparse* thnTracks = static_cast<THnSparse*>(fHistManager.FindObject("tracks"));
1198  if (!thnTracks) return;
1199  for (Int_t i = 0; i < thnTracks->GetNdimensions(); i++) {
1200  TString title(thnTracks->GetAxis(i)->GetTitle());
1201  if (title=="Centrality %")
1202  contents[i] = cent;
1203  else if (title=="#it{p}_{T} (GeV/#it{c})")
1204  contents[i] = trackPt;
1205  else if (title=="#eta")
1206  contents[i] = trackEta;
1207  else if (title=="#phi")
1208  contents[i] = trackPhi;
1209  else if (title=="#sigma(1/#it{p}_{T}) (GeV/#it{c})^{-1}")
1210  contents[i] = sigma1OverPt;
1211  else if (title=="#sigma(#it{p}_{T}) / #it{p}_{T}")
1212  contents[i] = sigma1OverPt*trackPt;
1213  else if (title=="track type")
1214  contents[i] = trackType;
1215  else
1216  AliWarning(Form("Unable to fill dimension %s of histogram %s!", title.Data(), thnTracks->GetName()));
1217  }
1218 
1219  thnTracks->Fill(contents);
1220 }
1221 
1222 //________________________________________________________________________
1223 void AliAnalysisTaskPWGJEQA::FillGeneratorLevelTHnSparse(Double_t cent, Double_t partEta, Double_t partPhi, Double_t partPt, Byte_t findable)
1224 {
1225  Double_t contents[20]={0};
1226 
1227  THnSparse* thnTracks_PhysPrim = static_cast<THnSparse*>(fHistManager.FindObject("tracks_PhysPrim"));
1228  if (!thnTracks_PhysPrim) return;
1229  for (Int_t i = 0; i < thnTracks_PhysPrim->GetNdimensions(); i++) {
1230  TString title(thnTracks_PhysPrim->GetAxis(i)->GetTitle());
1231  if (title=="Centrality %")
1232  contents[i] = cent;
1233  else if (title=="#it{p}_{T} (GeV/#it{c})")
1234  contents[i] = partPt;
1235  else if (title=="#eta")
1236  contents[i] = partEta;
1237  else if (title=="#phi")
1238  contents[i] = partPhi;
1239  else if (title=="Findable")
1240  contents[i] = findable;
1241  else
1242  AliWarning(Form("Unable to fill dimension %s of histogram %s!", title.Data(), thnTracks_PhysPrim->GetName()));
1243  }
1244 
1245  thnTracks_PhysPrim->Fill(contents);
1246 }
1247 
1248 //________________________________________________________________________
1250  Double_t trackEta, Double_t trackPhi, Double_t trackPt, Byte_t trackType)
1251 {
1252  Double_t contents[20]={0};
1253 
1254  THnSparse* thnTracks_Matched = static_cast<THnSparse*>(fHistManager.FindObject("tracks_Matched"));
1255  if (!thnTracks_Matched) return;
1256  for (Int_t i = 0; i < thnTracks_Matched->GetNdimensions(); i++) {
1257  TString title(thnTracks_Matched->GetAxis(i)->GetTitle());
1258  if (title=="Centrality %")
1259  contents[i] = cent;
1260  else if (title=="#it{p}_{T}^{gen} (GeV/#it{c})")
1261  contents[i] = partPt;
1262  else if (title=="#eta^{gen}")
1263  contents[i] = partEta;
1264  else if (title=="#phi^{gen}")
1265  contents[i] = partPhi;
1266  else if (title=="#it{p}_{T}^{det} (GeV/#it{c})")
1267  contents[i] = trackPt;
1268  else if (title=="#eta^{det}")
1269  contents[i] = trackEta;
1270  else if (title=="#phi^{det}")
1271  contents[i] = trackPhi;
1272  else if (title=="(#it{p}_{T}^{gen} - #it{p}_{T}^{det}) / #it{p}_{T}^{gen}")
1273  contents[i] = (partPt - trackPt) / partPt;
1274  else if (title=="(#it{p}_{T}^{gen} - #it{p}_{T}^{det}) / #it{p}_{T}^{det}")
1275  contents[i] = (partPt - trackPt) / trackPt;
1276  else if (title=="track type")
1277  contents[i] = (Double_t)trackType;
1278  else
1279  AliWarning(Form("Unable to fill dimension %s of histogram %s!", title.Data(), thnTracks_Matched->GetName()));
1280  }
1281 
1282  thnTracks_Matched->Fill(contents);
1283 }
Int_t pdg
Double_t * fEtaHistBins
! eta bins
Double_t * fPtResHistBins
! pt res bins
void Print(std::ostream &o, const char *name, Double_t dT, Double_t dVM, Double_t alldT, Double_t alldVM)
Definition: PlotSysInfo.C:121
TObjArray fClusterCollArray
cluster collection array
Bool_t fDoEventQA
Set whether to enable event QA.
Double_t GetRhoVal() const
const TString & GetRhoName() const
Int_t fNEtaHistBins
! number of eta bins
double Double_t
Definition: External.C:58
const char * title
Definition: MakeQAPdf.C:27
Int_t fNPtResHistBins
! number of pt res bins
AliTrackContainer * fDetectorLevel
! detector level container
Int_t fNPhiHistBins
! number of phi bins
Bool_t RetrieveEventObjects()
Retrieve common objects from event.
const AliMCParticleIterableMomentumContainer accepted_momentum() const
bidirectional stl iterator over the EMCAL iterable container
energy
Definition: HFPtSpectrum.C:44
void FillTH2(const char *hname, double x, double y, double weight=1., Option_t *opt="")
Fill a 2D histogram within the container.
AliAnalysisTaskPWGJEQA()
Default constructor for ROOT I/O purposes.
Bool_t fGeneralHistograms
whether or not it should fill some general histograms
THistManager fHistManager
Histogram manager.
Bool_t fRejectOutlierEvents
flag to reject pythia pt-hard jet outlier events
virtual Bool_t AcceptJet(Int_t i, UInt_t &rejectionReason) const
Int_t fNTotClusters[3]
!Total number of accepted clusters in current event (DCal/EMCal)
TString fDetectorLevelName
detector level container name
Bool_t FillHistograms()
Function filling histograms.
Bool_t fDoTrackQA
Set whether to enable track QA.
TString part
use mixed event to constrain combinatorial background
Definition: InvMassFit.C:52
TObjArray fParticleCollArray
particle/track collection array
Double_t * sigma
const Int_t nPtBins
Double_t * fPtHistBins
! pt bins
void GetLeadingHadronMomentum(TLorentzVector &mom, const AliEmcalJet *jet) const
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
void CreateTProfile(const char *name, const char *title, int nbinsX, double xmin, double xmax, Option_t *opt="")
Create a new TProfile within the container.
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
AliEMCALGeometry * fGeom
!emcal geometry
AliPHOSGeometry * fPHOSGeo
! phos geometry
void FillDetectorLevelTHnSparse(Double_t cent, Double_t trackEta, Double_t trackPhi, Double_t trackPt, Double_t sigma1OverPt, Byte_t trackType)
virtual AliAODMCParticle * GetAcceptMCParticleWithLabel(Int_t lab)
Int_t fNTotTracks
!Total number of accepted tracks in current event
TH1 * CreateTH1(const char *name, const char *title, int nbins, double xmin, double xmax, Option_t *opt="")
Create a new TH1 within the container.
Float_t fCellEnergyCut
Energy cell cut.
BeamType fForceBeamType
forced beam type
virtual Bool_t AcceptCluster(Int_t i, UInt_t &rejectionReason) const
Double_t * fPtRelDiffHistBins
! pt relative difference bins
void FillProfile(const char *name, double x, double y, double weight=1.)
Declaration of class AliAnalysisTaskPWGJEQA.
Double_t fCent
!event centrality
Int_t fNPtHistBins
! number of pt bins
TString fCaloCellsName
name of calo cell collection
Int_t fNIntegerHistBins
! number of integer bins
void FillTH1(const char *hname, double x, double weight=1., Option_t *opt="")
Fill a 1D histogram within the container.
Bool_t PythiaInfoFromFile(const char *currFile, Float_t &fXsec, Float_t &fTrials, Int_t &pthard)
Loading PYTHIA information from external cross section file into the task.
TObjArray fJetCollArray
jet collection array
short Short_t
Definition: External.C:23
AliVCaloCells * fCaloCells
!cells
AliRhoParameter * GetRhoParameter()
Int_t fNCentHistBins
! number of cent bins
void ExecOnce()
Perform steps needed to initialize the analysis.
Double_t * fIntegerHistBins
! integer bins
const AliClusterIterableMomentumContainer all_momentum() const
void FillMatchedParticlesTHnSparse(Double_t cent, Double_t partEta, Double_t partPhi, Double_t partPt, Double_t trackEta, Double_t trackPhi, Double_t trackPt, Byte_t trackType)
Float_t fPtHard
!event -hard
void FillGeneratorLevelTHnSparse(Double_t cent, Double_t partEta, Double_t partPhi, Double_t partPt, Byte_t findable)
static Double_t * GenerateFixedBinArray(Int_t n, Double_t min, Double_t max)
AliTLorentzVector fLeadingTrack
!Leading track in current event
AliEmcalList * fOutput
!output list
AliMCParticleContainer * GetMCParticleContainer(Int_t i=0) const
Char_t GetTrackType(const AliVTrack *track) const
Definition: External.C:220
Bool_t fIsEsd
!whether it&#39;s an ESD analysis
void GenerateHistoBins()
Generate histogram binning arrays.
AliTrackContainer * GetTrackContainer(Int_t i=0) const
Bool_t fDoCaloQA
Set whether to enable cell/cluster QA.
TH1 * fHistEventRejection
!book keep reasons for rejecting event
void SetMakeGeneralHistograms(Bool_t g)
This is a task used to do basic PWGJE QA on tracks, clusters, and jets.
Base task in the EMCAL jet framework.
Bool_t fIsPtHard
flag to enable pt-hard histos and make available outlier cut
const AliTrackIterableMomentumContainer accepted_momentum() const
Bool_t fDoJetQA
Set whether to enable jet QA.
TString fGeneratorLevelName
generator level container name
void SetRejectionReasonLabels(TAxis *axis)
void UserCreateOutputObjects()
Main initialization function on the worker.
Float_t fMaxPt
Histogram pt limit.
AliMCParticleContainer * fGeneratorLevel
! generator level container
const Int_t nbins
bool Bool_t
Definition: External.C:53
Double_t * fCentHistBins
! cent bins
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.
Bool_t fNeedEmcalGeom
whether or not the task needs the emcal geometry
Int_t fNPtRelDiffHistBins
! number of pt relative difference bins
Container for jet within the EMCAL jet framework.
Definition: External.C:196
Double_t * fPhiHistBins
! phi bins
AliTLorentzVector fLeadingCluster[3]
!Leading cluster in current event (EMCal/DCal)
const AliJetIterableContainer all() const
static Double_t fgkEMCalDCalPhiDivide
phi value used to distinguish between DCal and EMCal