AliPhysics  aaf9c62 (aaf9c62)
AliAnalysisTaskEmcalJetPerformance.cxx
Go to the documentation of this file.
1 /**********************************************************************************
2 * Copyright (C) 2016, Copyright Holders of the ALICE Collaboration *
3 * All rights reserved. *
4 * *
5 * Redistribution and use in source and binary forms, with or without *
6 * modification, are permitted provided that the following conditions are met: *
7 * * Redistributions of source code must retain the above copyright *
8 * notice, this list of conditions and the following disclaimer. *
9 * * Redistributions in binary form must reproduce the above copyright *
10 * notice, this list of conditions and the following disclaimer in the *
11 * documentation and/or other materials provided with the distribution. *
12 * * Neither the name of the <organization> nor the *
13 * names of its contributors may be used to endorse or promote products *
14 * derived from this software without specific prior written permission. *
15 * *
16 * THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS "AS IS" AND *
17 * ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE IMPLIED *
18 * WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE ARE *
19 * DISCLAIMED. IN NO EVENT SHALL ALICE COLLABORATION BE LIABLE FOR ANY *
20 * DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES *
21 * (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; *
22 * LOSS OF USE, DATA, OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND *
23 * ON ANY THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT *
24 * (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS *
25 * SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE. *
26 * *********************************************************************************/
27 
28 #include <vector>
29 
30 #include <TClonesArray.h>
31 #include <TH1F.h>
32 #include <TH2F.h>
33 #include <TH3.h>
34 #include <TList.h>
35 #include <THnSparse.h>
36 #include <TRandom3.h>
37 #include <TGrid.h>
38 #include <TFile.h>
39 
40 #include <AliVCluster.h>
41 #include <AliVParticle.h>
42 #include <AliLog.h>
43 #include "AliAnalysisManager.h"
44 #include <AliVEventHandler.h>
45 #include "AliTLorentzVector.h"
46 #include "AliEmcalJet.h"
47 #include "AliRhoParameter.h"
48 #include "AliJetContainer.h"
49 #include "AliParticleContainer.h"
50 #include "AliClusterContainer.h"
51 #include "AliEMCALGeometry.h"
52 #include "AliOADBContainer.h"
54 #include "AliEMCALTriggerPatchInfo.h"
56 #include "AliMCEvent.h"
57 
59 
63 
69  fPlotJetHistograms(kFALSE),
70  fPlotClusterHistograms(kFALSE),
71  fPlotParticleCompositionHistograms(kFALSE),
72  fComputeBackground(kFALSE),
73  fDoTriggerSimulation(kFALSE),
74  fPlotMatchedJetHistograms(kFALSE),
75  fComputeMBDownscaling(kFALSE),
76  fPlotDCal(kFALSE),
77  fDoClosureTest(kFALSE),
78  fMinPt(-100),
79  fMaxPt(250),
80  fNEtaBins(40),
81  fNPhiBins(200),
82  fNCentHistBins(0),
83  fCentHistBins(0),
84  fNPtHistBins(0),
85  fPtHistBins(0),
86  fNM02HistBins(0),
87  fM02HistBins(0),
88  fNEoverPBins(0),
89  fEoverPBins(0),
90  fTrackMatchingDeltaEtaMax(0.015),
91  fTrackMatchingDeltaPhiMax(0.030),
92  fMBUpscaleFactor(1.),
93  fMedianEMCal(0.),
94  fMedianDCal(0.),
95  fkEMCEJE(kFALSE),
96  fDoTriggerResponse(kFALSE),
97  fDoJetMatchingGeometrical(kFALSE),
98  fDoJetMatchingMCFraction(kFALSE),
99  fEmbeddingQA(),
100  fMCJetContainer(nullptr),
101  fDetJetContainer(nullptr),
102  fDetJetContainerPPIntermediate(nullptr),
103  fJetMatchingR(0.),
104  fMinSharedMomentumFraction(0.5),
105  fMCJetMinMatchingPt(0.15),
106  fPlotJetMatchCandThresh(1.),
107  fUseAliEventCuts(kTRUE),
108  fEventCuts(0),
109  fEventCutList(0),
110  fUseManualEventCuts(kFALSE),
111  fGeneratorLevel(0),
112  fHistManager()
113 {
115 }
116 
123  AliAnalysisTaskEmcalJet(name, kTRUE),
124  fPlotJetHistograms(kFALSE),
125  fPlotClusterHistograms(kFALSE),
127  fComputeBackground(kFALSE),
128  fDoTriggerSimulation(kFALSE),
130  fComputeMBDownscaling(kFALSE),
131  fPlotDCal(kFALSE),
132  fDoClosureTest(kFALSE),
133  fMinPt(-100),
134  fMaxPt(250),
135  fNEtaBins(40),
136  fNPhiBins(200),
137  fNCentHistBins(0),
138  fCentHistBins(0),
139  fNPtHistBins(0),
140  fPtHistBins(0),
141  fNM02HistBins(0),
142  fM02HistBins(0),
143  fNEoverPBins(0),
144  fEoverPBins(0),
147  fMBUpscaleFactor(1.),
148  fMedianEMCal(0.),
149  fMedianDCal(0.),
150  fkEMCEJE(kFALSE),
151  fDoTriggerResponse(kFALSE),
153  fDoJetMatchingMCFraction(kFALSE),
154  fEmbeddingQA(),
158  fJetMatchingR(0.),
160  fMCJetMinMatchingPt(0.15),
162  fUseAliEventCuts(kTRUE),
163  fEventCuts(0),
164  fEventCutList(0),
165  fUseManualEventCuts(kFALSE),
166  fGeneratorLevel(0),
167  fHistManager(name)
168 {
170 }
171 
176 {
177 }
178 
183 {
184  fNCentHistBins = 4;
186  fCentHistBins[0] = 0;
187  fCentHistBins[1] = 10;
188  fCentHistBins[2] = 30;
189  fCentHistBins[3] = 50;
190  fCentHistBins[4] = 90;
191 
192  fNPtHistBins = 82;
195  GenerateFixedBinArray(7, 0.3, 1, fPtHistBins+6);
196  GenerateFixedBinArray(10, 1, 3, fPtHistBins+13);
197  GenerateFixedBinArray(14, 3, 10, fPtHistBins+23);
198  GenerateFixedBinArray(10, 10, 20, fPtHistBins+37);
199  GenerateFixedBinArray(15, 20, 50, fPtHistBins+47);
200  GenerateFixedBinArray(20, 50, 150, fPtHistBins+62);
201 
202  fNM02HistBins = 81;
204  GenerateFixedBinArray(35, 0, 0.7, fM02HistBins);
205  GenerateFixedBinArray(6, 0.7, 1., fM02HistBins+35);
206  GenerateFixedBinArray(20, 1., 3., fM02HistBins+41);
207  GenerateFixedBinArray(10, 3., 5., fM02HistBins+61);
208  GenerateFixedBinArray(10, 5., 10., fM02HistBins+71);
209 
210  fNEoverPBins = 47;
212  GenerateFixedBinArray(30, 0, 1.5, fEoverPBins);
213  GenerateFixedBinArray(10, 1.5, 3.5, fEoverPBins+30);
214  GenerateFixedBinArray(7, 3.5, 10.5, fEoverPBins+40);
215 }
216 
222 {
224 
225  // Intialize AliEventCuts
226  if (fUseAliEventCuts) {
227  fEventCutList = new TList();
228  fEventCutList ->SetOwner();
229  fEventCutList ->SetName("EventCutOutput");
230 
231  fEventCuts.OverrideAutomaticTriggerSelection(fOffTrigger);
232  if(fUseManualEventCuts==1)
233  {
234  fEventCuts.SetManualMode();
235  fEventCuts.fMC = false;
236  fEventCuts.SetupLHC15o();
237  fEventCuts.fUseVariablesCorrelationCuts = true;
238  }
239  fEventCuts.AddQAplotsToList(fEventCutList);
240  fOutput->Add(fEventCutList);
241  }
242 
243  // Get the MC particle branch, in case it exists
244  fGeneratorLevel = GetMCParticleContainer("mcparticles");
245 
246  // Get jet containers, in order to check the jet acceptance criteria
247  // For geometrical matching, it is expected that there is one det-level and one truth-level container
248  // For MC fraction matching, it is expected that there is one Pb-Pb det-level, one pp det-level, and one pp truth-level container
250 
252  Printf("Geometrical jet matching enabled.");
253 
254  for (Int_t i=0; i<2; i++) {
255  auto jetCont = GetJetContainer(i);
256  TString jetContName = jetCont->GetName();
257  if (jetContName.Contains("mcparticles")) {
258  fMCJetContainer = jetCont;
259  }
260  else {
261  fDetJetContainer = jetCont;
262  }
263  }
264 
265  }
266  else if (fDoJetMatchingMCFraction) {
267  Printf("MC-fraction jet matching enabled.");
268 
269  for (Int_t i=0; i<3; i++) {
270  auto jetCont = GetJetContainer(i);
271  TString jetContName = jetCont->GetName();
272  if (jetContName.Contains("mcparticles")) {
273  fMCJetContainer = jetCont;
274  }
275  else if (jetContName.Contains("Combined")) {
276  fDetJetContainer = jetCont;
277  }
278  else {
280  }
281  }
282 
283  }
284 
285  if (!fMCJetContainer) {
286  Printf("No MC jet container found!");
287  }
288  Printf("mcJetContainer: %s", fMCJetContainer->GetName());
289 
290  if (!fDetJetContainer) {
291  Printf("No det-level jet container found!");
292  }
293  Printf("det-level JetContainer: %s", fDetJetContainer->GetName());
294 
297  Printf("No intermediate pp det-level jet container found, despite MC-fraction matching enabled!");
298  }
299  Printf("Intermediate pp det-level JetContainer: %s", fDetJetContainerPPIntermediate->GetName());
300  }
301 
302  }
303 
304  // Allocate histograms
305  if (fPlotJetHistograms) {
307  }
310  }
313  }
314  if (fComputeBackground) {
316  }
317  if (fDoTriggerSimulation) {
319  }
322  }
323 
324  // Initialize embedding QA
326  if (embeddingHelper) {
327  bool res = fEmbeddingQA.Initialize();
328  if (res) {
330  }
331  }
332 
333  TIter next(fHistManager.GetListOfHistograms());
334  TObject* obj = 0;
335  while ((obj = next())) {
336  fOutput->Add(obj);
337  }
338 
339  PostData(1, fOutput); // Post data for ALL output slots > 0 here.
340 }
341 
342 /*
343  * This function allocates the histograms for single jets.
344  * A set of histograms is allocated per each jet container.
345  */
347 {
348  TString histname;
349  TString title;
350 
351  Int_t nPtBins1 = TMath::CeilNint(fMaxPt-fMinPt);
352  Int_t nPtBins2 = TMath::CeilNint((fMaxPt-fMinPt)/2);
353  Int_t nPtBins5 = TMath::CeilNint((fMaxPt-fMinPt)/5);
354 
355  AliJetContainer* jets = 0;
356  TIter nextJetColl(&fJetCollArray);
357  while ((jets = static_cast<AliJetContainer*>(nextJetColl()))) {
358 
360  if (!fDoClosureTest) {
361  continue; // don't plot det-level histograms if embedding
362  }
363  }
364 
365  // Jet rejection reason
366  histname = TString::Format("%s/JetHistograms/hJetRejectionReason", jets->GetArrayName().Data());
367  if (fForceBeamType != kpp) {
368  title = histname + ";Centrality (%);Rejection reason;#it{p}_{T,corr} (GeV/#it{c});counts";
369  TH3* hist = fHistManager.CreateTH3(histname.Data(), title.Data(), 10, 0, 100, 32, 0, 32, nPtBins5, fMinPt, fMaxPt);
370  SetRejectionReasonLabels(hist->GetYaxis());
371  } else {
372  title = histname + ";Rejection reason;#it{p}_{T,corr} (GeV/#it{c});counts";
373  TH2* hist = fHistManager.CreateTH2(histname.Data(), title.Data(), 32, 0, 32, nPtBins5, fMinPt, fMaxPt);
374  SetRejectionReasonLabels(hist->GetXaxis());
375  }
376 
377  // Rho vs. Centrality
378  if (!jets->GetRhoName().IsNull()) {
379  histname = TString::Format("%s/JetHistograms/hRhoVsCent", jets->GetArrayName().Data());
380  title = histname + ";Centrality (%);#rho (GeV/#it{c});counts";
381  fHistManager.CreateTH2(histname.Data(), title.Data(), 50, 0, 100, 100, 0, 500);
382  }
383 
384  // (Centrality, pT, NEF)
385  Int_t nbinsx = 20; Int_t minx = 0; Int_t maxx = 100;
386  Int_t nbinsy = nPtBins1; Int_t miny = fMinPt; Int_t maxy = fMaxPt;
387  Int_t nbinsz = 50; Int_t minz = 0; Int_t maxz = 1.;
388 
389  histname = TString::Format("%s/JetHistograms/hNEFVsPtEMCal", jets->GetArrayName().Data());
390  if (fForceBeamType != kpp) {
391  title = histname + ";Centrality (%);#it{p}_{T}^{corr} (GeV/#it{c});NEF";
392  fHistManager.CreateTH3(histname.Data(), title.Data(), nbinsx, minx, maxx, nbinsy, miny, maxy, nbinsz, minz, maxz);
393  }
394  else {
395  title = histname + ";#it{p}_{T}^{corr} (GeV/#it{c});NEF";
396  fHistManager.CreateTH2(histname.Data(), title.Data(), nbinsy, miny, maxy, nbinsz, minz, maxz);
397  }
398 
399  if (fPlotDCal) {
400  histname = TString::Format("%s/JetHistograms/hNEFVsPtDCal", jets->GetArrayName().Data());
401  if (fForceBeamType != kpp) {
402  title = histname + ";Centrality (%);#it{p}_{T}^{corr} (GeV/#it{c});NEF";
403  fHistManager.CreateTH3(histname.Data(), title.Data(), nbinsx, minx, maxx, nbinsy, miny, maxy, nbinsz, minz, maxz);
404  }
405  else {
406  title = histname + ";#it{p}_{T}^{corr} (GeV/#it{c});NEF";
407  fHistManager.CreateTH2(histname.Data(), title.Data(), nbinsy, miny, maxy, nbinsz, minz, maxz);
408  }
409  }
410 
411  // (Centrality, pT upscaled, calo type)
412  if (fComputeMBDownscaling) {
413  histname = TString::Format("%s/JetHistograms/hPtUpscaledMB", jets->GetArrayName().Data());
414  title = histname + ";Centrality (%);#it{p}_{T}^{corr} (GeV/#it{c});type";
415  fHistManager.CreateTH3(histname.Data(), title.Data(), 20, 0, 100, nPtBins2, fMinPt, fMaxPt, 2, -0.5, 1.5, "s");
416  }
417 
418  // pT-leading vs. pT
419  histname = TString::Format("%s/JetHistograms/hPtLeadingVsPt", jets->GetArrayName().Data());
420  if (fForceBeamType != kpp) {
421  title = histname + ";Centrality (%);#it{p}_{T}^{corr} (GeV/#it{c});#it{p}_{T,particle}^{leading} (GeV/#it{c})";
422  fHistManager.CreateTH3(histname.Data(), title.Data(), 10, 0, 100, nPtBins2, fMinPt, fMaxPt, fMaxPt, 0, fMaxPt);
423  }
424  else {
425  title = histname + ";#it{p}_{T}^{corr} (GeV/#it{c});#it{p}_{T,particle}^{leading} (GeV/#it{c})";
426  fHistManager.CreateTH2(histname.Data(), title.Data(), nPtBins2, fMinPt, fMaxPt, nPtBins1, fMinPt, fMaxPt);
427  }
428 
429  // A vs. pT
430  histname = TString::Format("%s/JetHistograms/hAreaVsPt", jets->GetArrayName().Data());
431  if (fForceBeamType != kpp) {
432  title = histname + ";Centrality (%);#it{p}_{T}^{corr} (GeV/#it{c});#it{A}_{jet}";
433  fHistManager.CreateTH3(histname.Data(), title.Data(), 10, 0, 100, nPtBins2, fMinPt, fMaxPt, 50, 0, 0.5);
434  }
435  else {
436  title = histname + ";#it{p}_{T}^{corr} (GeV/#it{c});#it{A}_{jet}";
437  fHistManager.CreateTH2(histname.Data(), title.Data(), nPtBins2, fMinPt, fMaxPt, 50, 0, 0.5);
438  }
439 
440  // (Centrality, pT, z-leading (charged))
441  nbinsx = 20; minx = 0; maxx = 100;
442  nbinsy = nPtBins2; miny = fMinPt; maxy = fMaxPt;
443  nbinsz = 50; minz = 0; maxz = 1.;
444 
445  histname = TString::Format("%s/JetHistograms/hZLeadingVsPtEMCal", jets->GetArrayName().Data());
446  if (fForceBeamType != kpp) {
447  title = histname + ";Centrality (%);#it{p}_{T}^{corr} (GeV/#it{c});#it{z}_{leading}";
448  fHistManager.CreateTH3(histname.Data(), title.Data(), nbinsx, minx, maxx, nbinsy, miny, maxy, nbinsz, minz, maxz);
449  }
450  else {
451  title = histname + ";#it{p}_{T}^{corr} (GeV/#it{c});#it{z}_{leading}";
452  fHistManager.CreateTH2(histname.Data(), title.Data(), nbinsy, miny, maxy, nbinsz, minz, maxz);
453  }
454 
455  if (fPlotDCal) {
456  histname = TString::Format("%s/JetHistograms/hZLeadingVsPtDCal", jets->GetArrayName().Data());
457  if (fForceBeamType != kpp) {
458  title = histname + ";Centrality (%);#it{p}_{T}^{corr} (GeV/#it{c});#it{z}_{leading}";
459  fHistManager.CreateTH3(histname.Data(), title.Data(), nbinsx, minx, maxx, nbinsy, miny, maxy, nbinsz, minz, maxz);
460  }
461  else {
462  title = histname + ";#it{p}_{T}^{corr} (GeV/#it{c});#it{z}_{leading}";
463  fHistManager.CreateTH2(histname.Data(), title.Data(), nbinsy, miny, maxy, nbinsz, minz, maxz);
464  }
465  }
466 
467  // (Centrality, pT, z (charged))
468  nbinsx = 20; minx = 0; maxx = 100;
469  nbinsy = nPtBins2; miny = fMinPt; maxy = fMaxPt;
470  nbinsz = 50; minz = 0; maxz = 1.;
471 
472  histname = TString::Format("%s/JetHistograms/hZVsPtEMCal", jets->GetArrayName().Data());
473  if (fForceBeamType != kpp) {
474  title = histname + ";Centrality (%);#it{p}_{T}^{corr} (GeV/#it{c});#it{z}";
475  fHistManager.CreateTH3(histname.Data(), title.Data(), nbinsx, minx, maxx, nbinsy, miny, maxy, nbinsz, minz, maxz);
476  }
477  else {
478  title = histname + ";#it{p}_{T}^{corr} (GeV/#it{c});#it{z}";
479  fHistManager.CreateTH2(histname.Data(), title.Data(), nbinsy, miny, maxy, nbinsz, minz, maxz);
480  }
481 
482  if (fPlotDCal) {
483  histname = TString::Format("%s/JetHistograms/hZVsPtDCal", jets->GetArrayName().Data());
484  if (fForceBeamType != kpp) {
485  title = histname + ";Centrality (%);#it{p}_{T}^{corr} (GeV/#it{c});#it{z}";
486  fHistManager.CreateTH3(histname.Data(), title.Data(), nbinsx, minx, maxx, nbinsy, miny, maxy, nbinsz, minz, maxz);
487  }
488  else {
489  title = histname + ";#it{p}_{T}^{corr} (GeV/#it{c});#it{z}";
490  fHistManager.CreateTH2(histname.Data(), title.Data(), nbinsy, miny, maxy, nbinsz, minz, maxz);
491  }
492  }
493 
494  // (Centrality, pT, Nconst)
495  nbinsx = 20; minx = 0; maxx = 100;
496  nbinsy = nPtBins2; miny = fMinPt; maxy = fMaxPt;
497  nbinsz = 50; minz = 0; maxz = fMaxPt;
498 
499  histname = TString::Format("%s/JetHistograms/hNConstVsPtEMCal", jets->GetArrayName().Data());
500  if (fForceBeamType != kpp) {
501  title = histname + ";Centrality (%);#it{p}_{T}^{corr} (GeV/#it{c});No. of constituents";
502  fHistManager.CreateTH3(histname.Data(), title.Data(), nbinsx, minx, maxx, nbinsy, miny, maxy, nbinsz, minz, maxz);
503  }
504  else {
505  title = histname + ";#it{p}_{T}^{corr} (GeV/#it{c});No. of constituents";
506  fHistManager.CreateTH2(histname.Data(), title.Data(), nbinsy, miny, maxy, nbinsz, minz, maxz);
507  }
508 
509  if (fPlotDCal) {
510  histname = TString::Format("%s/JetHistograms/hNConstVsPtDCal", jets->GetArrayName().Data());
511  if (fForceBeamType != kpp) {
512  title = histname + ";Centrality (%);#it{p}_{T}^{corr} (GeV/#it{c});No. of constituents";
513  fHistManager.CreateTH3(histname.Data(), title.Data(), nbinsx, minx, maxx, nbinsy, miny, maxy, nbinsz, minz, maxz);
514  }
515  else {
516  title = histname + ";#it{p}_{T}^{corr} (GeV/#it{c});No. of constituents";
517  fHistManager.CreateTH2(histname.Data(), title.Data(), nbinsy, miny, maxy, nbinsz, minz, maxz);
518  }
519  }
520 
521  // (Centrality, pT) for eta<0 and eta>0
522  if (fForceBeamType != kpp) {
523  nbinsx = 20; minx = 0; maxx = 100;
524  nbinsy = nPtBins1; miny = fMinPt; maxy = fMaxPt;
525 
526  histname = TString::Format("%s/JetHistograms/hEtaPosVsPtEMCal", jets->GetArrayName().Data());
527  title = histname + ";Centrality (%);#it{p}_{T}^{corr} (GeV/#it{c})";
528  fHistManager.CreateTH2(histname.Data(), title.Data(), nbinsx, minx, maxx, nbinsy, miny, maxy);
529 
530  histname = TString::Format("%s/JetHistograms/hEtaNegVsPtEMCal", jets->GetArrayName().Data());
531  title = histname + ";Centrality (%);#it{p}_{T}^{corr} (GeV/#it{c})";
532  fHistManager.CreateTH2(histname.Data(), title.Data(), nbinsx, minx, maxx, nbinsy, miny, maxy);
533  }
534 
535  // (Centrality, jet pT, Enonlincorr - Ehadcorr)
536  if (fForceBeamType != kpp) {
537  nbinsx = 20; minx = 0; maxx = 100;
538  nbinsy = nPtBins2; miny = fMinPt; maxy = fMaxPt;
539  nbinsz = nPtBins2; minz = fMinPt; maxz = fMaxPt;
540 
541  histname = TString::Format("%s/JetHistograms/hDeltaEHadCorr", jets->GetArrayName().Data());
542  title = histname + ";Centrality (%);#it{p}_{T}^{corr} (GeV/#it{c});#sum#it{E}_{nonlincorr} - #it{E}_{hadcorr}";
543  fHistManager.CreateTH3(histname.Data(), title.Data(), nbinsx, minx, maxx, nbinsy, miny, maxy, nbinsz, minz, maxz);
544  }
545 
546  // (Median patch energy, calo type, jet pT, centrality)
547  if (fDoTriggerSimulation) {
548  histname = TString::Format("%s/JetHistograms/hMedPatchJet", jets->GetArrayName().Data());
549  title = histname + ";#it{E}_{patch,med};type;#it{p}_{T}^{corr} (GeV/#it{c});Centrality (%)";
550  Int_t nbins5[4] = {100, 2, nPtBins2, 50};
551  Double_t min5[4] = {0,-0.5, fMinPt, 0};
552  Double_t max5[4] = {50,1.5, fMaxPt, 100};
553  fHistManager.CreateTHnSparse(histname.Data(), title.Data(), 4, nbins5, min5, max5);
554  }
555 
556  }
557 
558  // MB downscale factor histogram
559  if (fComputeMBDownscaling) {
560  histname = "Trigger/hMBDownscaleFactor";
561  title = histname + ";Downscale factor;counts";
562  fHistManager.CreateTH1(histname.Data(), title.Data(), 200, 0, 200);
563  }
564 
565 }
566 
567 /*
568  * This function allocates the histograms for the calorimeter performance study.
569  */
571 {
572  TString histname;
573  TString htitle;
574 
575  const Int_t nRcorrBins = 50;
576  Double_t *RcorrBins = GenerateFixedBinArray(nRcorrBins, 0., 1.);
577  const Int_t nCellBins = 30;
578  Double_t *cellBins = GenerateFixedBinArray(nCellBins, -0.5, 29.5);
579  const Int_t nMatchedTrackBins = 5;
580  Double_t *matchedTrackBins = GenerateFixedBinArray(nMatchedTrackBins, -0.5, 4.5);
581  const Int_t nDeltaEtaBins = 60;
582  Double_t *deltaEtaBins = GenerateFixedBinArray(nDeltaEtaBins, -0.015, 0.015);
583 
586 
587  // Plot M02 distribution (centrality, Eclus nonlincorr, M02)
588  histname = "ClusterHistograms/hM02";
589  htitle = histname + ";Centrality (%);#it{E}_{clus} (GeV); M02";
591 
592  // Plot Ncell distribution for M02 > 0.4 and 0.1 < M02 < 0.4 (centrality, Eclus nonlincorr, Ncells)
593  histname = "ClusterHistograms/hNcellsM02G04";
594  htitle = histname + ";Centrality (%);#it{E}_{clus} (GeV); Ncells";
595  fHistManager.CreateTH3(histname.Data(), htitle.Data(), fNCentHistBins, fCentHistBins, fNPtHistBins, fPtHistBins, nCellBins, cellBins);
596 
597  histname = "ClusterHistograms/hNcellsM02L04";
598  htitle = histname + ";Centrality (%);#it{E}_{clus} (GeV); Ncells";
599  fHistManager.CreateTH3(histname.Data(), htitle.Data(), fNCentHistBins, fCentHistBins, fNPtHistBins, fPtHistBins, nCellBins, cellBins);
600 
603 
604  // Plot matched track pT for all clusters, M02 > 0.4 clusters, and 0.1 < M02 < 0.4 clusters (centrality, Eclus nonlincorr, trackPsum)
605  histname = "ClusterHistograms/hMatchedTrackPt";
606  htitle = histname + ";Centrality (%);#it{E}_{clus} (GeV); #Sigma#it{p}_{track} (GeV/c)";
608 
609  histname = "ClusterHistograms/hMatchedTrackPtM02G04";
610  htitle = histname + ";Centrality (%);#it{E}_{clus} (GeV); #Sigma#it{p}_{track} (GeV/c)";
612 
613  histname = "ClusterHistograms/hMatchedTrackPtM02L04";
614  htitle = histname + ";Centrality (%);#it{E}_{clus} (GeV); #Sigma#it{p}_{track} (GeV/c)";
616 
617  // Plot number of matched tracks for all clusters, M02 > 0.4 clusters, and 0.1 < M02 < 0.4 clusters (centrality, Eclus nonlincorr, N matches)
618  histname = "ClusterHistograms/hMatchedTrackN";
619  htitle = histname + ";Centrality (%);#it{E}_{clus} (GeV); N_{tracks}";
620  fHistManager.CreateTH3(histname.Data(), htitle.Data(), fNCentHistBins, fCentHistBins, fNPtHistBins, fPtHistBins, nMatchedTrackBins, matchedTrackBins);
621 
622  histname = "ClusterHistograms/hMatchedTrackNM02G04";
623  htitle = histname + ";Centrality (%);#it{E}_{clus} (GeV); N_{tracks}";
624  fHistManager.CreateTH3(histname.Data(), htitle.Data(), fNCentHistBins, fCentHistBins, fNPtHistBins, fPtHistBins, nMatchedTrackBins, matchedTrackBins);
625 
626  histname = "ClusterHistograms/hMatchedTrackNM02L04";
627  htitle = histname + ";Centrality (%);#it{E}_{clus} (GeV); N_{tracks}";
628  fHistManager.CreateTH3(histname.Data(), htitle.Data(), fNCentHistBins, fCentHistBins, fNPtHistBins, fPtHistBins, nMatchedTrackBins, matchedTrackBins);
629 
630  // Plot M02 distribution for clusters with matched tracks (centrality, Eclus nonlincorr, M02)
631  histname = "ClusterHistograms/hM02Matched";
632  htitle = histname + ";Centrality (%);#it{E}_{clus} (GeV); M02";
634 
635  // Plot M02 distribution for clusters without matched tracks (centrality, Eclus nonlincorr, M02)
636  histname = "ClusterHistograms/hM02Unmatched";
637  htitle = histname + ";Centrality (%);#it{E}_{clus} (GeV); M02";
639 
640  // Plot clus-track deltaEta of matched tracks (deltaEta, Eclus, M02)
641  histname = "ClusterHistograms/hDeltaEtaCentral";
642  htitle = histname + ";#eta_{track} - #eta_{clus};#it{E}_{clus} (GeV); M02";
643  fHistManager.CreateTH3(histname.Data(), htitle.Data(), nDeltaEtaBins, deltaEtaBins, fNPtHistBins, fPtHistBins, fNM02HistBins, fM02HistBins);
644 
645  histname = "ClusterHistograms/hDeltaEtaPeripheral";
646  htitle = histname + ";#eta_{track} - #eta_{clus};#it{E}_{clus} (GeV); M02";
647  fHistManager.CreateTH3(histname.Data(), htitle.Data(), nDeltaEtaBins, deltaEtaBins, fNPtHistBins, fPtHistBins, fNM02HistBins, fM02HistBins);
648 
651 
652  // Plot E/p vs. M02 for 0-10% and 50-90% (Eclus nonlincorr, Eclus nonlincorr / trackPsum, M02)
653  histname = "ClusterHistograms/hEoverPM02Central";
654  htitle = histname + ";#it{E}_{clus} (GeV); #it{E}_{clus} / #Sigma#it{p}_{track} (GeV); M02";
656 
657  histname = "ClusterHistograms/hEoverPM02Peripheral";
658  htitle = histname + ";#it{E}_{clus} (GeV); #it{E}_{clus} / #Sigma#it{p}_{track} (GeV); M02";
660 
663 
664  // Plot Rcorr distribution (centrality, trackPSum, Rcorr = (Enonlincorr - Ehadcorr) / trackPSum)
665  histname = "ClusterHistograms/hRcorrVsCent";
666  htitle = histname + ";Centrality (%);#Sigma#it{p}_{track} (GeV); R_{corr} = #frac{#DeltaE_{clus}}{#Sigmap_{track}}";
667  fHistManager.CreateTH3(histname.Data(), htitle.Data(), fNCentHistBins, fCentHistBins, fNPtHistBins, fPtHistBins, nRcorrBins, RcorrBins);
668 
669  // Plot Rcorr distribution for 0-10% centrality (Eclus nonlincorr, trackPSum, Rcorr)
670  histname = "ClusterHistograms/hRcorr0-10";
671  htitle = histname + ";#it{E}_{clus} (GeV);#Sigma#it{p}_{track} (GeV); R_{corr} = #frac{#DeltaE_{clus}}{#Sigmap_{track}}";
672  fHistManager.CreateTH3(histname.Data(), htitle.Data(), fNPtHistBins, fPtHistBins, fNPtHistBins, fPtHistBins, nRcorrBins, RcorrBins);
673 
674  // Plot Rcorr distribution for 50-90% centrality (Eclus nonlincorr, trackPSum, Rcorr)
675  histname = "ClusterHistograms/hRcorr50-90";
676  htitle = histname + ";#it{E}_{clus} (GeV);#Sigma#it{p}_{track} (GeV); R_{corr} = #frac{#DeltaE_{clus}}{#Sigmap_{track}}";
677  fHistManager.CreateTH3(histname.Data(), htitle.Data(), fNPtHistBins, fPtHistBins, fNPtHistBins, fPtHistBins, nRcorrBins, RcorrBins);
678 
679  // Plot also Rcorr-clus (centrality, trackPSum, Rcorr-clus = (Enonlincorr - Ehadcorr) / Enonlincorr )
680  histname = "ClusterHistograms/hRcorrClusVsCent";
681  htitle = histname + ";Centrality (%);#Sigma#it{p}_{track} (GeV); #frac{#DeltaE_{clus}}{E_{clus}}";
682  fHistManager.CreateTH3(histname.Data(), htitle.Data(), fNCentHistBins, fCentHistBins, fNPtHistBins, fPtHistBins, nRcorrBins, RcorrBins);
683 
684  // Rcorr-clus for 0-10% centrality (Eclus nonlincorr, trackPSum, Rcorr-clus = (Enonlincorr - Ehadcorr) / Enonlincorr )
685  histname = "ClusterHistograms/hRcorrClus0-10";
686  htitle = histname + ";#it{E}_{clus} (GeV);#Sigma#it{p}_{track} (GeV); #frac{#DeltaE_{clus}}{E_{clus}}";
687  fHistManager.CreateTH3(histname.Data(), htitle.Data(), fNPtHistBins, fPtHistBins, fNPtHistBins, fPtHistBins, nRcorrBins, RcorrBins);
688 
689  // Rcorr-clus for 50-90% centrality (Eclus nonlincorr, trackPSum, Rcorr-clus = (Enonlincorr - Ehadcorr) / Enonlincorr )
690  histname = "ClusterHistograms/hRcorrClus50-90";
691  htitle = histname + ";#it{E}_{clus} (GeV);#Sigma#it{p}_{track} (GeV); #frac{#DeltaE_{clus}}{E_{clus}}";
692  fHistManager.CreateTH3(histname.Data(), htitle.Data(), fNPtHistBins, fPtHistBins, fNPtHistBins, fPtHistBins, nRcorrBins, RcorrBins);
693 
694  // Plot total track multiplicity
695  histname = "ClusterHistograms/hTrackMultiplicity";
696  htitle = histname + ";N_{tracks};Centrality (%)";
697  fHistManager.CreateTH2(histname.Data(), htitle.Data(), 1000, 0, 10000, 20, 0, 100);
698 }
699 
700 /*
701  * This function allocates the histograms for the jet composition study.
702  */
704 {
705  TString histname;
706  TString htitle;
707  Int_t nPtBins = TMath::CeilNint(fMaxPt/2);
708 
709  const Int_t nRejBins = 32;
710  Double_t* rejReasonBins = new Double_t[nRejBins+1];
711  GenerateFixedBinArray(nRejBins, 0, nRejBins, rejReasonBins);
712  const Int_t nContributorTypes = 11;
713  Double_t *contributorTypeBins = GenerateFixedBinArray(nContributorTypes, -0.5, 10.5);
714  const Int_t nParticleTypes = 17;
715  Double_t *particleTypeBins = GenerateFixedBinArray(nParticleTypes, -0.5, 16.5);
716 
717  AliEmcalContainer* cont = 0;
718  TIter nextClusColl(&fClusterCollArray);
719  while ((cont = static_cast<AliEmcalContainer*>(nextClusColl()))) {
720 
721  histname = "ClusterHistogramsMC/hClusterRejectionReasonMC";
722  htitle = histname + ";Rejection reason;#it{E}_{clus} (GeV/)";
723  TH2* histMC2 = fHistManager.CreateTH2(histname.Data(), htitle.Data(), nRejBins, rejReasonBins, fNPtHistBins, fPtHistBins);
724  SetRejectionReasonLabels(histMC2->GetXaxis());
725  }
726 
727  // M02 vs. Energy vs. Particle type
728  histname = "ClusterHistogramsMC/hM02VsParticleTypeCentral";
729  htitle = histname + ";M02;#it{E}_{clus} (GeV); Particle type";
730  TH3* hM02VsParticleTypeCentral = fHistManager.CreateTH3(histname.Data(), htitle.Data(), fNM02HistBins, fM02HistBins, fNPtHistBins, fPtHistBins, nParticleTypes, particleTypeBins);
731  SetParticleTypeLabels(hM02VsParticleTypeCentral->GetZaxis());
732 
733  histname = "ClusterHistogramsMC/hM02VsParticleTypePeripheral";
734  htitle = histname + ";M02;#it{E}_{clus} (GeV); Particle type";
735  TH3* hM02VsParticleTypePeripheral = fHistManager.CreateTH3(histname.Data(), htitle.Data(), fNM02HistBins, fM02HistBins, fNPtHistBins, fPtHistBins, nParticleTypes, particleTypeBins);
736  SetParticleTypeLabels(hM02VsParticleTypePeripheral->GetZaxis());
737 
738  // Plot photon energy in photon-hadron overlap clusters (Centrality, Photon energy, M02)
739  histname = "ClusterHistogramsMC/hPhotonHadronPhotonEnergy";
740  htitle = histname + ";Centrality (%);M02;#it{E}_{photon} (GeV)";
742 
743  // Plot hadron energy in hadron-photon overlap clusters (Centrality, Photon energy, M02)
744  histname = "ClusterHistogramsMC/hHadronPhotonHadronEnergy";
745  htitle = histname + ";Centrality (%);M02;#it{E}_{hadron} (GeV)";
747 
748  if (fPlotJetHistograms) {
749 
750  // M02 vs. Energy vs. Particle type vs. Jet pT, for particles inside jets
751  Int_t dim = 0;
752  TString title[20];
753  Int_t nbins[20] = {0};
754  Double_t min[30] = {0.};
755  Double_t max[30] = {0.};
756  Double_t *binEdges[20] = {0};
757 
758  title[dim] = "M02";
759  nbins[dim] = fNM02HistBins;
760  binEdges[dim] = fM02HistBins;
761  min[dim] = fM02HistBins[0];
762  max[dim] = fM02HistBins[fNM02HistBins];
763  dim++;
764 
765  title[dim] = "#it{E}_{clus} (GeV)";
766  nbins[dim] = fNPtHistBins;
767  binEdges[dim] = fPtHistBins;
768  min[dim] = fPtHistBins[0];
769  max[dim] = fPtHistBins[fNPtHistBins];
770  dim++;
771 
772  title[dim] = "Contributor type";
773  nbins[dim] = nContributorTypes;
774  min[dim] = -0.5;
775  max[dim] = 7.5;
776  binEdges[dim] = contributorTypeBins;
777  dim++;
778 
779  title[dim] = "#it{p}_{T,jet}^{corr}";
780  nbins[dim] = nPtBins;
781  min[dim] = 0;
782  max[dim] = fMaxPt;
783  binEdges[dim] = GenerateFixedBinArray(nbins[dim], min[dim], max[dim]);
784  dim++;
785 
786  TString thnname = "JetPerformanceMC/hM02VsContributorTypeJets";
787  THnSparse* hn = fHistManager.CreateTHnSparse(thnname.Data(), thnname.Data(), dim, nbins, min, max);
788  for (Int_t i = 0; i < dim; i++) {
789  hn->GetAxis(i)->SetTitle(title[i]);
790  hn->SetBinEdges(i, binEdges[i]);
791  }
792 
793  // Particle composition inside each jet -- jet pT vs. particle type vs. particle number vs. particle pT sum
794  // (One entry per jet for each particle type)
795  dim = 0;
796 
797  title[dim] = "#it{p}_{T,jet}^{corr}";
798  nbins[dim] = nPtBins;
799  min[dim] = 0;
800  max[dim] = fMaxPt;
801  binEdges[dim] = GenerateFixedBinArray(nbins[dim], min[dim], max[dim]);
802  dim++;
803 
804  title[dim] = "Contributor type";
805  nbins[dim] = nContributorTypes;
806  min[dim] = -0.5;
807  max[dim] = 7.5;
808  binEdges[dim] = contributorTypeBins;
809  dim++;
810 
811  title[dim] = "N";
812  nbins[dim] = 30;
813  min[dim] = -0.5;
814  max[dim] = 29.5;
815  binEdges[dim] = GenerateFixedBinArray(nbins[dim], min[dim], max[dim]);
816  dim++;
817 
818  title[dim] = "#it{p}_{T,sum} (GeV)";
819  nbins[dim] = fNPtHistBins;
820  binEdges[dim] = fPtHistBins;
821  min[dim] = fPtHistBins[0];
822  max[dim] = fPtHistBins[fNPtHistBins];
823  dim++;
824 
825  thnname = "JetPerformanceMC/hJetComposition";
826  THnSparse* thn = fHistManager.CreateTHnSparse(thnname.Data(), thnname.Data(), dim, nbins, min, max);
827  for (Int_t i = 0; i < dim; i++) {
828  thn->GetAxis(i)->SetTitle(title[i]);
829  thn->SetBinEdges(i, binEdges[i]);
830  }
831 
832  // Hadronic calo energy in each jet
833 
834  // Jet pT vs. Summed energy of hadronic clusters without a matched track
835  histname = "JetPerformance/hHadCaloEnergyUnmatched";
836  htitle = histname + ";#it{p}_{T,jet} (GeV);#it{p}_{T,had} (GeV)";
837  fHistManager.CreateTH2(histname.Data(), htitle.Data(), fNPtHistBins, fPtHistBins, fNPtHistBins, fPtHistBins);
838 
839  // Jet pT vs. Summed energy of hadronic clusters with a matched track (before hadronic correction)
840  histname = "JetPerformance/hHadCaloEnergyMatchedNonlincorr";
841  htitle = histname + ";#it{p}_{T,jet} (GeV);#it{p}_{T,had} (GeV)";
842  fHistManager.CreateTH2(histname.Data(), htitle.Data(), fNPtHistBins, fPtHistBins, fNPtHistBins, fPtHistBins);
843 
844  // Jet pT vs. Summed energy of hadronic clusters with a matched track (after hadronic correction)
845  histname = "JetPerformance/hHadCaloEnergyMatchedHadCorr";
846  htitle = histname + ";#it{p}_{T,jet} (GeV);#it{p}_{T,had} (GeV)";
847  fHistManager.CreateTH2(histname.Data(), htitle.Data(), fNPtHistBins, fPtHistBins, fNPtHistBins, fPtHistBins);
848 
849  }
850 
851 }
852 
853 /*
854  * This function sets axis labels for particle type histograms.
855  */
857 {
858  axis->SetBinLabel(1, "SinglePhoton");
859  axis->SetBinLabel(2, "SingleElectron");
860  axis->SetBinLabel(3, "SingleChargedPion");
861  axis->SetBinLabel(4, "SingleProton");
862  axis->SetBinLabel(5, "SingleAntiProton");
863  axis->SetBinLabel(6, "SingleChargedKaon");
864  axis->SetBinLabel(7, "SingleK0L");
865  axis->SetBinLabel(8, "SingleNeutron");
866  axis->SetBinLabel(9, "SingleAntiNeutron");
867  axis->SetBinLabel(10, "SingleOther");
868  axis->SetBinLabel(11, "PhotonHadron");
869  axis->SetBinLabel(12, "HadronPhoton");
870  axis->SetBinLabel(13, "MergedPi0");
871  axis->SetBinLabel(14, "PhotonPhotonOther");
872  axis->SetBinLabel(15, "HadronHadron");
873  axis->SetBinLabel(16, "TwoContributorsOther");
874  axis->SetBinLabel(17, "MoreThanTwoContributors");
875 }
876 
877 /*
878  * This function allocates background subtraction histograms, if enabled.
879  * A set of histograms is allocated per each jet container.
880  */
882 {
883  TString histname;
884  TString title;
885 
886  AliJetContainer* jets = 0;
887  TIter nextJetColl(&fJetCollArray);
888  while ((jets = static_cast<AliJetContainer*>(nextJetColl()))) {
889 
890  // EMCal
891  histname = TString::Format("%s/BackgroundHistograms/hScaleFactorEMCal", jets->GetArrayName().Data());
892  title = histname + ";Centrality;Scale factor;counts";
893  fHistManager.CreateTH2(histname.Data(), title.Data(), 50, 0, 100, 100, 0, 5);
894 
895  histname = TString::Format("%s/BackgroundHistograms/hDeltaPtEMCal", jets->GetArrayName().Data());
896  title = histname + ";Centrality (%);#delta#it{p}_{T} (GeV/#it{c});counts";
897  fHistManager.CreateTH2(histname.Data(), title.Data(), 10, 0, 100, 400, -50, 150);
898 
899  histname = TString::Format("%s/BackgroundHistograms/hScaleFactorEMCalFid", jets->GetArrayName().Data());
900  title = histname + ";Centrality;Scale factor;counts";
901  fHistManager.CreateTH2(histname.Data(), title.Data(), 50, 0, 100, 100, 0, 5);
902 
903  // DCal
904  if (fPlotDCal) {
905  histname = TString::Format("%s/BackgroundHistograms/hScaleFactorDCal", jets->GetArrayName().Data());
906  title = histname + ";Centrality;Scale factor;counts";
907  fHistManager.CreateTH2(histname.Data(), title.Data(), 50, 0, 100, 100, 0, 5);
908 
909  histname = TString::Format("%s/BackgroundHistograms/hDeltaPtDCal", jets->GetArrayName().Data());
910  title = histname + ";Centrality (%);#delta#it{p}_{T} (GeV/#it{c});counts";
911  fHistManager.CreateTH2(histname.Data(), title.Data(), 10, 0, 100, 400, -50, 150);
912 
913  histname = TString::Format("%s/BackgroundHistograms/hScaleFactorDCalFid", jets->GetArrayName().Data());
914  title = histname + ";Centrality;Scale factor;counts";
915  fHistManager.CreateTH2(histname.Data(), title.Data(), 50, 0, 100, 100, 0, 5);
916  }
917 
918  }
919 }
920 
921 /*
922  * This function allocates the histograms for single jets, when the "simulated" trigger has been fired.
923  * A set of histograms is allocated per each jet container.
924  */
926 {
927  TString histname;
928  TString title;
929  Int_t nPtBins1 = TMath::CeilNint(fMaxPt-fMinPt);
930  Int_t nPtBins2 = TMath::CeilNint((fMaxPt-fMinPt)/2);
931  Int_t nPtBins5 = TMath::CeilNint((fMaxPt-fMinPt)/5);
932 
933  //----------------------------------------------
934  // Trigger patch histograms
935 
936  // patch eta vs. phi
937  histname = "TriggerSimHistograms/hEtaVsPhi";
938  title = histname + ";#eta_{patch} (rad);#phi_{patch} (rad)";
939  fHistManager.CreateTH2(histname.Data(), title.Data(), 140, -0.7, 0.7, 500, 1., 6.);
940 
941  // N patches
942  histname = "TriggerSimHistograms/hNPatches";
943  title = histname + ";#it{N}_{patches};type";
944  fHistManager.CreateTH2(histname.Data(), title.Data(), 200, 0, 200, 2, -0.5, 1.5);
945 
946  // patch E vs. centrality
947  histname = "TriggerSimHistograms/hPatchE";
948  title = histname + ";Centrality (%);#it{E}_{patch} (GeV)";
949  fHistManager.CreateTH2(histname.Data(), title.Data(), 50, 0, 100, nPtBins2, 0, fMaxPt);
950 
951  // patch median vs. Centrality
952  histname = "TriggerSimHistograms/hPatchMedianE";
953  title = histname + ";Centrality (%);#it{E}_{patch,med} (GeV);type";
954  fHistManager.CreateTH3(histname.Data(), title.Data(), 50, 0, 100, 100, 0, 50, 2, -0.5, 1.5);
955 
956  if (fDoTriggerResponse) {
957  histname = "TriggerSimHistograms/hMaxPatchResponseMatrix";
958  title = histname + ";#it{E}_{maxpatch,det} (GeV);#it{E}_{maxpatch,truth} (GeV)";
959  fHistManager.CreateTH2(histname.Data(), title.Data(), 400, 0, 200, 400, 0, 200);
960 
961  histname = "TriggerSimHistograms/hEMCalEnergyResponseMatrix";
962  title = histname + ";#Sigma#it{E}_{clus,EMCal} (GeV);#it{E}_{EMCal,truth} (GeV)";
963  fHistManager.CreateTH2(histname.Data(), title.Data(), 500, 0, 500, 500, 0, 500);
964  }
965 
966  //----------------------------------------------
967  // Jet histograms for "triggered" events
968  AliJetContainer* jets = 0;
969  TIter nextJetColl(&fJetCollArray);
970  while ((jets = static_cast<AliJetContainer*>(nextJetColl()))) {
971 
972  // Jet rejection reason
973  histname = TString::Format("%s/TriggerSimHistograms/hJetRejectionReason", jets->GetArrayName().Data());
974  title = histname + ";Rejection reason;#it{p}_{T,jet} (GeV/#it{c});counts";
975  TH2* hist = fHistManager.CreateTH2(histname.Data(), title.Data(), 32, 0, 32, nPtBins5, fMinPt, fMaxPt);
976  SetRejectionReasonLabels(hist->GetXaxis());
977 
978  // Rho vs. Centrality
979  if (!jets->GetRhoName().IsNull()) {
980  histname = TString::Format("%s/TriggerSimHistograms/hRhoVsCent", jets->GetArrayName().Data());
981  title = histname + ";Centrality (%);#rho (GeV/#it{c});counts";
982  fHistManager.CreateTH2(histname.Data(), title.Data(), 50, 0, 100, 100, 0, 500);
983  }
984 
985  // (Centrality, pT, NEF)
986  Int_t nbinsx = 20; Int_t minx = 0; Int_t maxx = 100;
987  Int_t nbinsy = nPtBins1; Int_t miny = fMinPt; Int_t maxy = fMaxPt;
988  Int_t nbinsz = 50; Int_t minz = 0; Int_t maxz = 1.;
989 
990  histname = TString::Format("%s/TriggerSimHistograms/hNEFVsPtEMCal", jets->GetArrayName().Data());
991  title = histname + ";Centrality (%);#it{p}_{T}^{corr} (GeV/#it{c});NEF";
992  fHistManager.CreateTH3(histname.Data(), title.Data(), nbinsx, minx, maxx, nbinsy, miny, maxy, nbinsz, minz, maxz);
993 
994  if (fPlotDCal) {
995  histname = TString::Format("%s/TriggerSimHistograms/hNEFVsPtDCal", jets->GetArrayName().Data());
996  title = histname + ";Centrality (%);#it{p}_{T}^{corr} (GeV/#it{c});NEF";
997  fHistManager.CreateTH3(histname.Data(), title.Data(), nbinsx, minx, maxx, nbinsy, miny, maxy, nbinsz, minz, maxz);
998  }
999 
1000  // pT-leading vs. pT
1001  histname = TString::Format("%s/TriggerSimHistograms/hPtLeadingVsPt", jets->GetArrayName().Data());
1002  title = histname + ";#it{p}_{T}^{corr} (GeV/#it{c});#it{p}_{T,particle}^{leading} (GeV/#it{c})";
1003  fHistManager.CreateTH2(histname.Data(), title.Data(), nPtBins2, fMinPt, fMaxPt, nPtBins2, fMinPt, fMaxPt);
1004 
1005  // A vs. pT
1006  histname = TString::Format("%s/TriggerSimHistograms/hAreaVsPt", jets->GetArrayName().Data());
1007  title = histname + ";#it{p}_{T}^{corr} (GeV/#it{c});#it{A}_{jet}";
1008  fHistManager.CreateTH2(histname.Data(), title.Data(), nPtBins2, fMinPt, fMaxPt, 50, 0, 0.5);
1009 
1010  // (Centrality, pT, z-leading (charged))
1011  nbinsx = 20; minx = 0; maxx = 100;
1012  nbinsy = nPtBins2; miny = fMinPt; maxy = fMaxPt;
1013  nbinsz = 50; minz = 0; maxz = 1.;
1014 
1015  histname = TString::Format("%s/TriggerSimHistograms/hZLeadingVsPtEMCal", jets->GetArrayName().Data());
1016  title = histname + ";Centrality (%);#it{p}_{T}^{corr} (GeV/#it{c});#it{z}_{leading}";
1017  fHistManager.CreateTH3(histname.Data(), title.Data(), nbinsx, minx, maxx, nbinsy, miny, maxy, nbinsz, minz, maxz);
1018 
1019  if (fPlotDCal) {
1020  histname = TString::Format("%s/TriggerSimHistograms/hZLeadingVsPtDCal", jets->GetArrayName().Data());
1021  title = histname + ";Centrality (%);#it{p}_{T}^{corr} (GeV/#it{c});#it{z}_{leading}";
1022  fHistManager.CreateTH3(histname.Data(), title.Data(), nbinsx, minx, maxx, nbinsy, miny, maxy, nbinsz, minz, maxz);
1023  }
1024 
1025  // z (charged) vs. pT
1026  nbinsx = 20; minx = 0; maxx = 100;
1027  nbinsy = nPtBins2; miny = fMinPt; maxy = fMaxPt;
1028  nbinsz = 50; minz = 0; maxz = 1.;
1029 
1030  histname = TString::Format("%s/TriggerSimHistograms/hZVsPtEMCal", jets->GetArrayName().Data());
1031  title = histname + ";Centrality (%);#it{p}_{T}^{corr} (GeV/#it{c});#it{z}";
1032  fHistManager.CreateTH3(histname.Data(), title.Data(), nbinsx, minx, maxx, nbinsy, miny, maxy, nbinsz, minz, maxz);
1033 
1034  if (fPlotDCal) {
1035  histname = TString::Format("%s/TriggerSimHistograms/hZVsPtDCal", jets->GetArrayName().Data());
1036  title = histname + ";Centrality (%);#it{p}_{T}^{corr} (GeV/#it{c});#it{z}";
1037  fHistManager.CreateTH3(histname.Data(), title.Data(), nbinsx, minx, maxx, nbinsy, miny, maxy, nbinsz, minz, maxz);
1038  }
1039 
1040  // (Centrality, pT, Nconst)
1041  nbinsx = 20; minx = 0; maxx = 100;
1042  nbinsy = nPtBins2; miny = fMinPt; maxy = fMaxPt;
1043  nbinsz = 50; minz = 0; maxz = fMaxPt;
1044 
1045  histname = TString::Format("%s/TriggerSimHistograms/hNConstVsPtEMCal", jets->GetArrayName().Data());
1046  title = histname + ";Centrality (%);#it{p}_{T}^{corr} (GeV/#it{c});No. of constituents";
1047  fHistManager.CreateTH3(histname.Data(), title.Data(), nbinsx, minx, maxx, nbinsy, miny, maxy, nbinsz, minz, maxz);
1048 
1049  if (fPlotDCal) {
1050  histname = TString::Format("%s/TriggerSimHistograms/hNConstVsPtDCal", jets->GetArrayName().Data());
1051  title = histname + ";Centrality (%);#it{p}_{T}^{corr} (GeV/#it{c});No. of constituents";
1052  fHistManager.CreateTH3(histname.Data(), title.Data(), nbinsx, minx, maxx, nbinsy, miny, maxy, nbinsz, minz, maxz);
1053  }
1054 
1055  }
1056 }
1057 
1058 /*
1059  * This function allocates histograms for matched truth-det jets in the case of embedding.
1060  */
1062 {
1063  TString histname;
1064  TString title;
1065  Int_t nPtBins1 = TMath::CeilNint(fMaxPt-fMinPt);
1066  Int_t nPtBinsTruth2 = TMath::CeilNint(fMaxPt/2);
1067 
1068  // Response matrix, (centrality, pT-truth, pT-det)
1069  Int_t nbinsx = 20; Int_t minx = 0; Int_t maxx = 100;
1070  Int_t nbinsy = fMaxPt; Int_t miny = 0; Int_t maxy = fMaxPt;
1071  Int_t nbinsz = nPtBins1; Int_t minz = fMinPt; Int_t maxz = fMaxPt;
1072 
1073  histname = "MatchedJetHistograms/hResponseMatrixEMCal";
1074  if (fForceBeamType != kpp) {
1075  title = histname + ";Centrality (%);#it{p}_{T}^{truth} (GeV/#it{c});#it{p}_{T,corr}^{det} (GeV/#it{c})";
1076  fHistManager.CreateTH3(histname.Data(), title.Data(), nbinsx, minx, maxx, nbinsy, miny, maxy, nbinsz, minz, maxz);
1077  }
1078  else {
1079  title = histname + ";#it{p}_{T,corr}^{det} (GeV/#it{c});#it{p}_{T}^{truth} (GeV/#it{c})";
1080  fHistManager.CreateTH2(histname.Data(), title.Data(), nbinsy, miny, maxy, nbinsy, miny, maxy);
1081  }
1082 
1083  if (fPlotDCal) {
1084  histname = "MatchedJetHistograms/hResponseMatrixDCal";
1085  if (fForceBeamType != kpp) {
1086  title = histname + ";Centrality (%);#it{p}_{T}^{truth} (GeV/#it{c});#it{p}_{T,corr}^{det} (GeV/#it{c})";
1087  fHistManager.CreateTH3(histname.Data(), title.Data(), nbinsx, minx, maxx, nbinsy, miny, maxy, nbinsz, minz, maxz);
1088  }
1089  else {
1090  title = histname + ";#it{p}_{T}^{truth} (GeV/#it{c});#it{p}_{T,corr}^{det} (GeV/#it{c})";
1091  fHistManager.CreateTH2(histname.Data(), title.Data(), nbinsy, miny, maxy, nbinsz, minz, maxz);
1092  }
1093  }
1094 
1095  // JES shift, (centrality, pT-truth, (pT-det - pT-truth) / pT-truth)
1096  nbinsx = 20; minx = 0; maxx = 100;
1097  nbinsy = nPtBinsTruth2; miny = 0; maxy = fMaxPt;
1098  nbinsz = 250; minz = -5.; maxz = 5.;
1099 
1100  histname = "MatchedJetHistograms/hJESshiftEMCal";
1101  if (fForceBeamType != kpp) {
1102  title = histname + ";Centrality (%);#it{p}_{T}^{truth} (GeV/#it{c});#frac{#it{p}_{T,corr}^{det} - #it{p}_{T}^{truth}}{#it{p}_{T}^{truth}}";
1103  fHistManager.CreateTH3(histname.Data(), title.Data(), nbinsx, minx, maxx, nbinsy, miny, maxy, nbinsz, minz, maxz);
1104  }
1105  else {
1106  title = histname + ";#it{p}_{T}^{truth} (GeV/#it{c});#frac{#it{p}_{T,corr}^{det} - #it{p}_{T}^{truth}}{#it{p}_{T}^{truth}}";
1107  fHistManager.CreateTH2(histname.Data(), title.Data(), nbinsy, miny, maxy, nbinsz, minz, maxz);
1108  }
1109 
1110  if (fPlotDCal) {
1111  histname = "MatchedJetHistograms/hJESshiftDCal";
1112  if (fForceBeamType != kpp) {
1113  title = histname + ";Centrality (%);#it{p}_{T}^{truth} (GeV/#it{c});#frac{#it{p}_{T,corr}^{det} - #it{p}_{T}^{truth}}{#it{p}_{T}^{truth}}";
1114  fHistManager.CreateTH3(histname.Data(), title.Data(), nbinsx, minx, maxx, nbinsy, miny, maxy, nbinsz, minz, maxz);
1115  }
1116  else {
1117  title = histname + ";#it{p}_{T}^{truth} (GeV/#it{c});#frac{#it{p}_{T,corr}^{det} - #it{p}_{T}^{truth}}{#it{p}_{T}^{truth}}";
1118  fHistManager.CreateTH2(histname.Data(), title.Data(), nbinsy, miny, maxy, nbinsz, minz, maxz);
1119  }
1120  }
1121 
1122  // NEF of det-level matched jets, (centrality, pT-truth, NEF)
1123  nbinsx = 20; minx = 0; maxx = 100;
1124  nbinsy = fMaxPt; miny = 0; maxy = fMaxPt;
1125  nbinsz = 50; minz = 0; maxz = 1.;
1126 
1127  histname = "MatchedJetHistograms/hNEFVsPt";
1128  if (fForceBeamType != kpp) {
1129  title = histname + ";Centrality (%);#it{p}_{T}^{truth} (GeV/#it{c});Calo energy fraction";
1130  fHistManager.CreateTH3(histname.Data(), title.Data(), nbinsx, minx, maxx, nbinsy, miny, maxy, nbinsz, minz, maxz);
1131  }
1132  else {
1133  title = histname + ";#it{p}_{T}^{truth} (GeV/#it{c});Calo energy fraction";
1134  fHistManager.CreateTH2(histname.Data(), title.Data(), nbinsy, miny, maxy, nbinsz, minz, maxz);
1135  }
1136 
1137  // z-leading (charged) of det-level matched jets, (centrality, pT-truth, z-leading)
1138  nbinsx = 20; minx = 0; maxx = 100;
1139  nbinsy = nPtBinsTruth2; miny = 0; maxy = fMaxPt;
1140  nbinsz = 50; minz = 0; maxz = 1.;
1141 
1142  histname = "MatchedJetHistograms/hZLeadingVsPt";
1143  if (fForceBeamType != kpp) {
1144  title = histname + ";Centrality (%);#it{p}_{T}^{truth} (GeV/#it{c});#it{z}_{leading}";
1145  fHistManager.CreateTH3(histname.Data(), title.Data(), nbinsx, minx, maxx, nbinsy, miny, maxy, nbinsz, minz, maxz);
1146  }
1147  else {
1148  title = histname + ";#it{p}_{T}^{truth} (GeV/#it{c});#it{z}_{leading}";
1149  fHistManager.CreateTH2(histname.Data(), title.Data(), nbinsy, miny, maxy, nbinsz, minz, maxz);
1150  }
1151 
1153 
1154  // Matching distance, (pT-det, pT-truth, deltaR)
1155  nbinsx = nPtBins1; minx = fMinPt; maxx = fMaxPt;
1156  nbinsy = nPtBinsTruth2; miny = 0; maxy = fMaxPt;
1157  nbinsz = 100; minz = 0; maxz = 4.;
1158  histname = "MatchedJetHistograms/hMatchingDistance";
1159  title = histname + ";#it{p}_{T}^{det} (GeV/#it{c});#it{p}_{T}^{truth} (GeV/#it{c});R";
1160  fHistManager.CreateTH3(histname.Data(), title.Data(), nbinsx, minx, maxx, nbinsy, miny, maxy, nbinsz, minz, maxz);
1161 
1162  // Jet matching QA (copied from AliAnalysisTaskEmcalJetHCorrelations.cxx)
1163  histname = "MatchedJetHistograms/fHistJetMatchingQA";
1164  title = histname;
1165  std::vector<std::string> binLabels = {"noMatch", "matchedJet", "uniqueMatch", "jetDistance", "passedAllCuts"};
1166  auto histMatchedJetCuts = fHistManager.CreateTH1(histname.Data(), title.Data(), binLabels.size(), 0, binLabels.size());
1167  // Set label names
1168  for (unsigned int i = 1; i <= binLabels.size(); i++) {
1169  histMatchedJetCuts->GetXaxis()->SetBinLabel(i, binLabels.at(i-1).c_str());
1170  }
1171  histMatchedJetCuts->GetYaxis()->SetTitle("Number of jets");
1172  }
1173 
1175 
1176  // (det pT, shared MC fraction, deltaR) of closest jets
1177  nbinsx = nPtBins1; minx = fMinPt; maxx = fMaxPt;
1178  nbinsy = 20; miny = 0; maxy = 1.;
1179  nbinsz = 100; minz = 0; maxz = 4.;
1180  histname = "MatchedJetHistograms/hMatchingDistanceVsMCFraction";
1181  title = histname + ";#it{p}_{T}^{det} (GeV/#it{c});MC fraction;#DeltaR";
1182  fHistManager.CreateTH3(histname.Data(), title.Data(), nbinsx, minx, maxx, nbinsy, miny, maxy, nbinsz, minz, maxz);
1183 
1184  // Jet matching QA (copied from AliAnalysisTaskEmcalJetHCorrelations.cxx)
1185  histname = "MatchedJetHistograms/fHistJetMatchingQA";
1186  title = histname;
1187  std::vector<std::string> binLabels = {"noMatch", "matchedJet", "sharedMomentumFraction", "partLevelMatchedJet", "jetDistance", "passedAllCuts"};
1188  auto histMatchedJetCuts = fHistManager.CreateTH1(histname.Data(), title.Data(), binLabels.size(), 0, binLabels.size());
1189  // Set label names
1190  for (unsigned int i = 1; i <= binLabels.size(); i++) {
1191  histMatchedJetCuts->GetXaxis()->SetBinLabel(i, binLabels.at(i-1).c_str());
1192  }
1193  histMatchedJetCuts->GetYaxis()->SetTitle("Number of jets");
1194  }
1195 
1196 }
1197 
1203 {
1204  // Configure base class to set fTriggerPatchInfo to array of trigger patches, each event
1205  // (Need to call this before base class ExecOnce)
1206  if (fDoTriggerSimulation) {
1207  this->SetCaloTriggerPatchInfoName("EmcalTriggers");
1208  }
1209 
1211 
1212  fNeedEmcalGeom = kTRUE;
1213 
1214  // Check if trigger patches are loaded
1215  if (fDoTriggerSimulation) {
1216  if (fTriggerPatchInfo) {
1217  TString objname(fTriggerPatchInfo->GetClass()->GetName());
1218  TClass cls(objname);
1219  if (!cls.InheritsFrom("AliEMCALTriggerPatchInfo")) {
1220  AliError(Form("%s: Objects of type %s in %s are not inherited from AliEMCALTriggerPatchInfo!",
1221  GetName(), cls.GetName(), "EmcalTriggers"));
1222  fTriggerPatchInfo = 0;
1223  }
1224  }
1225  if (!fTriggerPatchInfo) {
1226  AliError(Form("%s: Unable to get trigger patch container with name %s. Aborting", GetName(), "EmcalTriggers"));
1227  return;
1228  }
1229  }
1230 }
1231 
1236 
1237  // Get the downscaling factors for MB triggers (to be used to calculate trigger efficiency)
1238 
1240 
1241  // Get instance of the downscale factor helper class
1243  downscaleOCDB->SetRun(InputEvent()->GetRunNumber());
1244 
1245  // There are two possible min bias triggers for LHC15o
1246  TString triggerNameMB1 = "CINT7-B-NOPF-CENT";
1247  TString triggerNameMB2 = "CV0L7-B-NOPF-CENT";
1248  TString triggerNameJE = "CINT7EJ1-B-NOPF-CENTNOPMD";
1249 
1250  // Get the downscale factor for whichever MB trigger exists in the given run
1251  std::vector<TString> runtriggers = downscaleOCDB->GetTriggerClasses();
1252  Double_t downscalefactor;
1253  for (auto i : runtriggers) {
1254  if (i.EqualTo(triggerNameMB1) || i.EqualTo(triggerNameMB2)) {
1255  downscalefactor = downscaleOCDB->GetDownscaleFactorForTriggerClass(i.Data());
1256  break;
1257  }
1258  }
1259 
1260  // Store the inverse of the downscale factor, used later to weight the pT spectrum
1261  fMBUpscaleFactor = 1/downscalefactor;
1262 
1263  TString histname = "Trigger/hMBDownscaleFactor";
1264  fHistManager.FillTH1(histname.Data(), fMBUpscaleFactor);
1265 
1266  }
1267 
1268 }
1269 
1274 {
1275  if (fUseAliEventCuts) {
1276  if (!fEventCuts.AcceptEvent(InputEvent()))
1277  {
1278  PostData(1, fOutput);
1279  return kFALSE;
1280  }
1281  }
1282  else {
1284  }
1285  return kTRUE;
1286 }
1287 
1296 {
1297  TString histname;
1298  AliJetContainer* jetCont = 0;
1299  TIter next(&fJetCollArray);
1300  while ((jetCont = static_cast<AliJetContainer*>(next()))) {
1301  TString jetContName = jetCont->GetName();
1302 
1303  // Do a simple trigger simulation (if requested)
1304  if (fDoTriggerSimulation) {
1306  }
1307 
1308  }
1309 
1310  // Compute the full jet background scale factor and delta-pt
1311  if (fComputeBackground) {
1313  }
1314 
1315  // Only fill the embedding qa plots if:
1316  // - We are using the embedding helper
1317  // - The class has been initialized
1318  // - Both jet collections are available
1319  if (fEmbeddingQA.IsInitialized()) {
1321  }
1322 
1323  return kTRUE;
1324 }
1325 
1330 {
1331  TString histname;
1332 
1333  // Check if trigger patches are loaded
1334  if (fTriggerPatchInfo) {
1335  TString objname(fTriggerPatchInfo->GetClass()->GetName());
1336  TClass cls(objname);
1337  if (!cls.InheritsFrom("AliEMCALTriggerPatchInfo")) {
1338  AliError(Form("%s: Objects of type %s in %s are not inherited from AliEMCALTriggerPatchInfo!",
1339  GetName(), cls.GetName(), "EmcalTriggers"));
1340  fTriggerPatchInfo = 0;
1341  }
1342  }
1343  if (!fTriggerPatchInfo) {
1344  AliError(Form("%s: Unable to get trigger patch container with name %s. Aborting", GetName(), "EmcalTriggers"));
1345  return;
1346  }
1347 
1348  // Compute patches in EMCal, DCal (I want offline simple trigger patch, i.e. patch calculated using FEE energy)
1349  std::vector<Double_t> vecEMCal;
1350  std::vector<Double_t> vecDCal;
1351  Double_t maxPatchEnergy = 0.;
1352  Double_t patchE;
1353  AliEMCALTriggerPatchInfo *maxPatch = nullptr;
1354  for(auto p : *fTriggerPatchInfo){
1355  AliEMCALTriggerPatchInfo *recpatch = static_cast<AliEMCALTriggerPatchInfo *>(p);
1356  if (recpatch) {
1357 
1358  if(!recpatch->IsJetHighSimple()) continue;
1359 
1360  patchE = recpatch->GetPatchE();
1361 
1362  histname = "TriggerSimHistograms/hEtaVsPhi";
1363  fHistManager.FillTH2(histname.Data(), recpatch->GetEtaGeo(), recpatch->GetPhiGeo());
1364 
1365  histname = "TriggerSimHistograms/hPatchE";
1366  fHistManager.FillTH2(histname.Data(), fCent, patchE);
1367 
1368  if (recpatch->IsEMCal()) {
1369  vecEMCal.push_back(patchE);
1370  if (patchE > maxPatchEnergy) {
1371  maxPatchEnergy = patchE;
1372  maxPatch = recpatch;
1373  }
1374  } else {
1375  vecDCal.push_back(patchE);
1376  }
1377 
1378  }
1379  }
1380 
1381  // Compute the median in each calorimeter
1382  const Int_t nBkgPatchesEMCal = vecEMCal.size(); // 6*8;
1383  const Int_t nBkgPatchesDCal = vecDCal.size(); // 4*5;
1384  fMedianEMCal = TMath::Median(nBkgPatchesEMCal, &vecEMCal[0]); // point to array used internally by vector
1385  fMedianDCal = TMath::Median(nBkgPatchesDCal, &vecDCal[0]);
1386 
1387  histname = "TriggerSimHistograms/hPatchMedianE";
1388  fHistManager.FillTH3(histname.Data(), fCent, fMedianEMCal, kEMCal);
1389  fHistManager.FillTH3(histname.Data(), fCent, fMedianDCal, kDCal);
1390 
1391  histname = "TriggerSimHistograms/hNPatches";
1392  fHistManager.FillTH2(histname.Data(), nBkgPatchesEMCal, kEMCal);
1393  fHistManager.FillTH2(histname.Data(), nBkgPatchesDCal, kDCal);
1394 
1395  // Fill max patch response and total EMCal energy response, if requested
1397 
1398  Double_t phiMinEMCal = fGeom->GetArm1PhiMin() * TMath::DegToRad(); // 80
1399  Double_t phiMaxEMCal = fGeom->GetEMCALPhiMax() * TMath::DegToRad(); // ~188
1400  Double_t detEnergyEMCal = 0;
1401  Double_t truthEnergyEMCal = 0;
1402  Double_t truthEnergyPatch = 0;
1403 
1404  AliTLorentzVector part;
1405  Double_t partEta;
1406  Double_t partPhi;
1407  Double_t partE;
1408  for (auto mcpart:fGeneratorLevel->accepted_momentum()) {
1409  part.Clear();
1410  part = mcpart.first;
1411  partEta = part.Eta();
1412  partPhi = part.Phi_0_2pi();
1413  partE = part.E();
1414  if (maxPatch) {
1415  if (partPhi < maxPatch->GetPhiMax() && partPhi > maxPatch->GetPhiMin()) {
1416  if (partEta < maxPatch->GetEtaMax() && partEta > maxPatch->GetEtaMin()) {
1417  truthEnergyPatch += partE;
1418  }
1419  }
1420  }
1421  if (partPhi < phiMaxEMCal && partPhi > phiMinEMCal) {
1422  if (partEta < 0.7 && partEta > -0.7) {
1423  truthEnergyEMCal += partE;
1424  }
1425  }
1426  }
1427 
1428  AliTLorentzVector clusFourVec;
1429  const AliVCluster* clus;
1430  Double_t clusEta;
1431  Double_t clusPhi;
1433  for (auto clusIterator : clusters->accepted_momentum() ) {
1434  clusFourVec.Clear();
1435  clusFourVec = clusIterator.first;
1436  clusEta = clusFourVec.Eta();
1437  clusPhi = clusFourVec.Phi_0_2pi();
1438  clus = clusIterator.second;
1439  if (clus->IsEMCAL()) {
1440  if (clusPhi < phiMaxEMCal && clusPhi > phiMinEMCal) {
1441  if (clusEta < 0.7 && clusEta > -0.7) {
1442  detEnergyEMCal += clusFourVec.E();
1443  }
1444  }
1445  }
1446  }
1447 
1448  if (maxPatch) {
1449  histname = "TriggerSimHistograms/hMaxPatchResponseMatrix";
1450  fHistManager.FillTH2(histname.Data(), maxPatchEnergy, truthEnergyPatch);
1451  }
1452 
1453  histname = "TriggerSimHistograms/hEMCalEnergyResponseMatrix";
1454  fHistManager.FillTH2(histname.Data(), detEnergyEMCal, truthEnergyEMCal);
1455 
1456  }
1457 
1458  // Then compute background subtracted patches, by subtracting from each patch the median patch E from the opposite hemisphere
1459  // If a patch is above threshold, the event is "triggered"
1460  Bool_t fkEMCEJE = kFALSE;
1461  Double_t threshold = 20;
1462  for(auto p : *fTriggerPatchInfo){
1463  AliEMCALTriggerPatchInfo *recpatch = static_cast<AliEMCALTriggerPatchInfo *>(p);
1464  if (recpatch) {
1465 
1466  if(!recpatch->IsJetHighSimple()) continue;
1467 
1468  if (recpatch->IsEMCal()) {
1469  if ((recpatch->GetPatchE() - fMedianDCal) > threshold) {
1470  fkEMCEJE = kTRUE;
1471  break;
1472  }
1473  } else {
1474  if ((recpatch->GetPatchE() - fMedianEMCal) > threshold) {
1475  fkEMCEJE = kTRUE;
1476  break;
1477  }
1478  }
1479  }
1480  }
1481 
1482  if (fkEMCEJE) {
1484  }
1485 
1486 }
1487 
1495 {
1496 
1497  if (fPlotJetHistograms) {
1499  }
1500  if (fPlotClusterHistograms) {
1502  }
1505  }
1508  }
1509 
1510  return kTRUE;
1511 }
1512 
1518 {
1519  TString histname;
1520  AliJetContainer* jets = 0;
1521  TIter nextJetColl(&fJetCollArray);
1522  while ((jets = static_cast<AliJetContainer*>(nextJetColl()))) {
1523  TString jetContName = jets->GetName();
1524 
1526  if (!fDoClosureTest) {
1527  continue; // don't plot det-level histograms if embedding
1528  }
1529  }
1530 
1531  Double_t rhoVal = 0;
1532  if (jets->GetRhoParameter()) {
1533  rhoVal = jets->GetRhoVal();
1534  histname = TString::Format("%s/JetHistograms/hRhoVsCent", jets->GetArrayName().Data());
1535  fHistManager.FillTH2(histname.Data(), fCent, rhoVal);
1536  }
1537 
1538  for (auto jet : jets->all()) {
1539 
1540  Float_t ptLeading = jets->GetLeadingHadronPt(jet);
1541  Float_t corrPt = GetJetPt(jet, rhoVal);
1542 
1543  // compute jet acceptance type
1544  Double_t type = GetJetType(jet);
1545  if ( type != kEMCal ) {
1546  if ( type != kDCal || !fPlotDCal ) {
1547  continue;
1548  }
1549  }
1550 
1551  // (Centrality, Area, pT) (fill before area cut)
1552  histname = TString::Format("%s/JetHistograms/hAreaVsPt", jets->GetArrayName().Data());
1553  if (fForceBeamType != kpp) {
1554  fHistManager.FillTH3(histname.Data(), fCent, corrPt, jet->Area());
1555  }
1556  else {
1557  fHistManager.FillTH2(histname.Data(), corrPt, jet->Area());
1558  }
1559 
1560  // (Centrality, pT-leading, pT) (before leading hadron cuts)
1561  histname = TString::Format("%s/JetHistograms/hPtLeadingVsPt", jets->GetArrayName().Data());
1562  if (fForceBeamType != kpp) {
1563  fHistManager.FillTH3(histname.Data(), fCent, corrPt, ptLeading);
1564  }
1565  else {
1566  fHistManager.FillTH2(histname.Data(), corrPt, ptLeading);
1567  }
1568 
1569  // (Centrality, pT, z-leading (charged)) (before leading hadron cuts)
1570  if (type == kEMCal) {
1571  histname = TString::Format("%s/JetHistograms/hZLeadingVsPtEMCal", jets->GetArrayName().Data());
1572  }
1573  else if (type == kDCal) {
1574  histname = TString::Format("%s/JetHistograms/hZLeadingVsPtDCal", jets->GetArrayName().Data());
1575  }
1576  TLorentzVector leadPart;
1577  jets->GetLeadingHadronMomentum(leadPart, jet);
1578  Double_t z = GetParallelFraction(leadPart.Vect(), jet);
1579  if (z == 1 || (z > 1 && z - 1 < 1e-3)) z = 0.999; // so that it will contribute to the bin <1
1580  if (fForceBeamType != kpp) {
1581  fHistManager.FillTH3(histname, fCent, corrPt, z);
1582  }
1583  else {
1584  fHistManager.FillTH2(histname, corrPt, z);
1585  }
1586 
1587  // Rejection reason
1588  UInt_t rejectionReason = 0;
1589  if (!jets->AcceptJet(jet, rejectionReason)) {
1590  histname = TString::Format("%s/JetHistograms/hJetRejectionReason", jets->GetArrayName().Data());
1591  if (fForceBeamType != kpp) {
1592  fHistManager.FillTH3(histname.Data(), fCent, jets->GetRejectionReasonBitPosition(rejectionReason), corrPt);
1593  }
1594  else {
1595  fHistManager.FillTH2(histname.Data(), jets->GetRejectionReasonBitPosition(rejectionReason), corrPt);
1596  }
1597  continue;
1598  }
1599 
1600  // (Centrality, pT, NEF)
1601  if (type == kEMCal) {
1602  histname = TString::Format("%s/JetHistograms/hNEFVsPtEMCal", jets->GetArrayName().Data());
1603  }
1604  else if (type == kDCal) {
1605  histname = TString::Format("%s/JetHistograms/hNEFVsPtDCal", jets->GetArrayName().Data());
1606  }
1607  if (fForceBeamType != kpp) {
1608  fHistManager.FillTH3(histname, fCent, corrPt, jet->NEF());
1609  }
1610  else {
1611  fHistManager.FillTH2(histname, corrPt, jet->NEF());
1612  }
1613 
1614  // (Centrality, pT upscaled, calo type)
1615  if (fComputeMBDownscaling) {
1616  histname = TString::Format("%s/JetHistograms/hPtUpscaledMB", jets->GetArrayName().Data());
1617  fHistManager.FillTH3(histname.Data(), fCent, corrPt, type, fMBUpscaleFactor);
1618  }
1619 
1620  // (Centrality, pT, z (charged))
1621  if (type == kEMCal) {
1622  histname = TString::Format("%s/JetHistograms/hZVsPtEMCal", jets->GetArrayName().Data());
1623  }
1624  else if (type == kDCal) {
1625  histname = TString::Format("%s/JetHistograms/hZVsPtDCal", jets->GetArrayName().Data());
1626  }
1627  const AliVTrack* track;
1628  for (Int_t i=0; i<jet->GetNumberOfTracks(); i++) {
1629  track = static_cast<AliVTrack*>(jet->Track(i));
1630  z = track->Pt() / TMath::Abs(corrPt);
1631  if (fForceBeamType != kpp) {
1632  fHistManager.FillTH3(histname, fCent, corrPt, z);
1633  }
1634  else {
1635  fHistManager.FillTH2(histname, corrPt, z);
1636  }
1637  }
1638 
1639  // (Centrality, pT, Nconst)
1640  if (type == kEMCal) {
1641  histname = TString::Format("%s/JetHistograms/hNConstVsPtEMCal", jets->GetArrayName().Data());
1642  }
1643  else if (type == kDCal) {
1644  histname = TString::Format("%s/JetHistograms/hNConstVsPtDCal", jets->GetArrayName().Data());
1645  }
1646  if (fForceBeamType != kpp) {
1647  fHistManager.FillTH3(histname, fCent, corrPt, 1.*jet->GetNumberOfConstituents());
1648  }
1649  else {
1650  fHistManager.FillTH2(histname, corrPt, 1.*jet->GetNumberOfConstituents());
1651  }
1652 
1653  // (Centrality, pT) for eta<0 and eta>0
1654  if (fForceBeamType != kpp) {
1655  if (type == kEMCal) {
1656  if (jet->Eta() > 0) {
1657  histname = TString::Format("%s/JetHistograms/hEtaPosVsPtEMCal", jets->GetArrayName().Data());
1658  fHistManager.FillTH2(histname, fCent, corrPt);
1659  }
1660  else if (jet->Eta() < 0) {
1661  histname = TString::Format("%s/JetHistograms/hEtaNegVsPtEMCal", jets->GetArrayName().Data());
1662  fHistManager.FillTH2(histname, fCent, corrPt);
1663  }
1664  }
1665  }
1666 
1667  // (Centrality, jet pT, Enonlincorr - Ehadcorr)
1668  if (fForceBeamType != kpp) {
1669  Double_t deltaEhadcorr = 0;
1670  const AliVCluster* clus = nullptr;
1671  Int_t nClusters = jet->GetNumberOfClusters();
1672  for (Int_t iClus = 0; iClus < nClusters; iClus++) {
1673  clus = jet->Cluster(iClus);
1674  deltaEhadcorr += (clus->GetNonLinCorrEnergy() - clus->GetHadCorrEnergy());
1675  }
1676 
1677  histname = TString::Format("%s/JetHistograms/hDeltaEHadCorr", jets->GetArrayName().Data());
1678  fHistManager.FillTH3(histname, fCent, corrPt, deltaEhadcorr);
1679  }
1680 
1681 
1682  // (Median patch energy, calo type, jet pT, centrality)
1683  if (fDoTriggerSimulation) {
1684  histname = TString::Format("%s/JetHistograms/hMedPatchJet", jets->GetArrayName().Data());
1685  Double_t x[4] = {fMedianEMCal, kEMCal, corrPt, fCent};
1686  fHistManager.FillTHnSparse(histname, x);
1687  Double_t y[4] = {fMedianDCal, kDCal, corrPt, fCent};
1688  fHistManager.FillTHnSparse(histname, y);
1689  }
1690 
1691  } //jet loop
1692 
1693  }
1694 }
1695 
1696 /*
1697  * This function fills the histograms for the calorimeter performance study.
1698  */
1700 {
1701  TString histname;
1702 
1703  // Loop through clusters
1705  const AliVCluster* clus;
1706  for (auto it : clusters->accepted_momentum()) {
1707 
1708  clus = it.second;
1709  Double_t clusPhi = it.first.Phi_0_2pi();
1710  Double_t clusEta = it.first.Eta();
1711 
1712  // Include only EMCal/DCal clusters
1713  if (!clus->IsEMCAL()) {
1714  continue;
1715  }
1716 
1717  // Compute the sum of matched track momentum, and various track matching / hadronic corretion quantities
1718  Double_t trackPSum = 0;
1719  Int_t nTracksMatched = 0;
1720  const AliVTrack* track = nullptr;
1721  for (Int_t itrack=0; itrack < clus->GetNTracksMatched(); itrack++) {
1722  track = dynamic_cast<AliVTrack*>(clus->GetTrackMatched(itrack));
1723  if (!track) {
1724  continue;
1725  }
1726 
1727  Double_t trackPhi = TVector2::Phi_0_2pi(track->GetTrackPhiOnEMCal());
1728  Double_t trackEta = track->GetTrackEtaOnEMCal();
1729  Double_t deta = TMath::Abs(clusEta - trackEta);
1730  Double_t dphi = TMath::Abs(clusPhi - trackPhi);
1731 
1733  trackPSum += track->P();
1734  nTracksMatched++;
1735  }
1736  }
1737 
1738  Double_t EoverP = 0;
1739  if (trackPSum > 1e-3) {
1740  EoverP = clus->GetNonLinCorrEnergy() / trackPSum;
1741  }
1742 
1743  Double_t deltaE = clus->GetNonLinCorrEnergy() - clus->GetHadCorrEnergy();
1744  Double_t Rcorr = 0;
1745  if (trackPSum > 1e-3) {
1746  Rcorr = deltaE / trackPSum;
1747  }
1748  Double_t RcorrClus = deltaE / clus->GetNonLinCorrEnergy();
1749 
1752 
1753  // Fill M02 distribution (centrality, Eclus nonlincorr, M02)
1754  histname = "ClusterHistograms/hM02";
1755  fHistManager.FillTH3(histname, fCent, clus->GetNonLinCorrEnergy(), clus->GetM02());
1756 
1757  // Plot Ncell distribution for M02 > 0.4 or M02 < 0.4 (centrality, Eclus nonlincorr, Ncells)
1758  if (clus->GetM02() > 0.4) {
1759  histname = "ClusterHistograms/hNcellsM02G04";
1760  fHistManager.FillTH3(histname, fCent, clus->GetNonLinCorrEnergy(), clus->GetNCells());
1761  }
1762  if (clus->GetM02() > 0.1 && clus->GetM02() < 0.4) {
1763  histname = "ClusterHistograms/hNcellsM02L04";
1764  fHistManager.FillTH3(histname, fCent, clus->GetNonLinCorrEnergy(), clus->GetNCells());
1765  }
1766 
1769 
1770  // Plot matched track pT (centrality, Eclus nonlincorr, trackPsum)
1771  histname = "ClusterHistograms/hMatchedTrackPt";
1772  fHistManager.FillTH3(histname, fCent, clus->GetNonLinCorrEnergy(), trackPSum);
1773 
1774  if (clus->GetM02() > 0.4) {
1775  histname = "ClusterHistograms/hMatchedTrackPtM02G04";
1776  fHistManager.FillTH3(histname, fCent, clus->GetNonLinCorrEnergy(), trackPSum);
1777  }
1778  if (clus->GetM02() > 0.1 && clus->GetM02() < 0.4) {
1779  histname = "ClusterHistograms/hMatchedTrackPtM02L04";
1780  fHistManager.FillTH3(histname, fCent, clus->GetNonLinCorrEnergy(), trackPSum);
1781  }
1782 
1783  // Plot number of matched tracks (centrality, Eclus nonlincorr, N matches)
1784  histname = "ClusterHistograms/hMatchedTrackN";
1785  fHistManager.FillTH3(histname, fCent, clus->GetNonLinCorrEnergy(), nTracksMatched);
1786 
1787  if (clus->GetM02() > 0.4) {
1788  histname = "ClusterHistograms/hMatchedTrackNM02G04";
1789  fHistManager.FillTH3(histname, fCent, clus->GetNonLinCorrEnergy(), nTracksMatched);
1790  }
1791  if (clus->GetM02() > 0.1 && clus->GetM02() < 0.4) {
1792  histname = "ClusterHistograms/hMatchedTrackNM02L04";
1793  fHistManager.FillTH3(histname, fCent, clus->GetNonLinCorrEnergy(), nTracksMatched);
1794  }
1795 
1796  // Plot M02 distribution for clusters with matched tracks (centrality, Eclus nonlincorr, M02)
1797  histname = "ClusterHistograms/hM02Matched";
1798  if (nTracksMatched > 0) {
1799  fHistManager.FillTH3(histname, fCent, clus->GetNonLinCorrEnergy(), clus->GetM02());
1800  }
1801 
1802  // Plot M02 distribution for clusters without matched tracks (centrality, Eclus nonlincorr, M02)
1803  histname = "ClusterHistograms/hM02Unmatched";
1804  if (nTracksMatched == 0) {
1805  fHistManager.FillTH3(histname, fCent, clus->GetNonLinCorrEnergy(), clus->GetM02());
1806  }
1807 
1808  // Plot clus-track deltaEta if there is one matched track (deltaEta, Eclus, M02)
1809  if (nTracksMatched == 1) {
1810 
1811  const AliVTrack* track = dynamic_cast<AliVTrack*>(clus->GetTrackMatched(0));
1812  if (track) {
1813  Double_t trackEta = track->GetTrackEtaOnEMCal();
1814  Double_t deta = trackEta - clusEta;
1815 
1816  if (fCent > 0 && fCent < 10) {
1817  histname = "ClusterHistograms/hDeltaEtaCentral";
1818  fHistManager.FillTH3(histname.Data(), deta, clus->GetNonLinCorrEnergy(), clus->GetM02());
1819  }
1820 
1821  if (fCent > 50 && fCent < 90) {
1822  histname = "ClusterHistograms/hDeltaEtaPeripheral";
1823  fHistManager.FillTH3(histname.Data(), deta, clus->GetNonLinCorrEnergy(), clus->GetM02());
1824  }
1825  }
1826  }
1827 
1830 
1831  // Plot E/p vs. M02 for 0-10% and 50-90% (Eclus nonlincorr, Eclus nonlincorr / trackPsum, M02)
1832  if (fCent > 0 && fCent < 10) {
1833  histname = "ClusterHistograms/hEoverPM02Central";
1834  fHistManager.FillTH3(histname, clus->GetNonLinCorrEnergy(), EoverP, clus->GetM02());
1835  }
1836  if (fCent > 50 && fCent < 90) {
1837  histname = "ClusterHistograms/hEoverPM02Peripheral";
1838  fHistManager.FillTH3(histname, clus->GetNonLinCorrEnergy(), EoverP, clus->GetM02());
1839  }
1840 
1843 
1844  // Fill Rcorr distribution (centrality, trackPSum, Rcorr = (Enonlincorr - Ehadcorr) / trackPSum)
1845  histname = "ClusterHistograms/hRcorrVsCent";
1846  fHistManager.FillTH3(histname, fCent, trackPSum, Rcorr);
1847 
1848  // Fill Rcorr distribution for 0-10% centrality (Eclus nonlincorr, trackPSum, Rcorr)
1849  if (fCent > 0 && fCent < 10) {
1850  histname = "ClusterHistograms/hRcorr0-10";
1851  fHistManager.FillTH3(histname, clus->GetNonLinCorrEnergy(), trackPSum, Rcorr);
1852  }
1853 
1854  // Fill Rcorr distribution for 50-90% centrality (Eclus nonlincorr, trackPSum, Rcorr)
1855  if (fCent > 50 && fCent < 90) {
1856  histname = "ClusterHistograms/hRcorr50-90";
1857  fHistManager.FillTH3(histname, clus->GetNonLinCorrEnergy(), trackPSum, Rcorr);
1858  }
1859 
1860  // Fill also Rcorr-clus (centrality, trackPSum, Rcorr-clus = (Enonlincorr - Ehadcorr) / Enonlincorr )
1861  histname = "ClusterHistograms/hRcorrClusVsCent";
1862  fHistManager.FillTH3(histname, fCent, trackPSum, RcorrClus);
1863 
1864  // Fill Rcorr-clus for 0-10% centrality (Eclus nonlincorr, trackPSum, Rcorr-clus = (Enonlincorr - Ehadcorr) / Enonlincorr )
1865  if (fCent > 0 && fCent < 10) {
1866  histname = "ClusterHistograms/hRcorrClus0-10";
1867  fHistManager.FillTH3(histname, clus->GetNonLinCorrEnergy(), trackPSum, RcorrClus);
1868  }
1869 
1870  // Fill Rcorr-clus for 50-90% centrality (Eclus nonlincorr, trackPSum, Rcorr-clus = (Enonlincorr - Ehadcorr) / Enonlincorr )
1871  if (fCent > 50 && fCent < 90) {
1872  histname = "ClusterHistograms/hRcorrClus50-90";
1873  fHistManager.FillTH3(histname, clus->GetNonLinCorrEnergy(), trackPSum, RcorrClus);
1874  }
1875 
1876  }
1877 
1878  // Fill total track multiplicity
1879  histname = "ClusterHistograms/hTrackMultiplicity";
1880  AliTrackContainer* trackCont = dynamic_cast<AliTrackContainer*>(GetParticleContainer("tracks"));
1881  Int_t nTracks = trackCont->GetNAcceptedTracks();
1882  fHistManager.FillTH2(histname.Data(), nTracks, fCent);
1883 
1884 }
1885 
1886 /*
1887  * This function fills particle composition histograms for the calorimeter performance study in MC.
1888  */
1890 {
1891  // If MC, get the MC event
1892  const AliMCEvent* mcevent = nullptr;
1893  if (fGeneratorLevel) {
1894  mcevent = MCEvent();
1895  }
1896  else {
1897  return;
1898  }
1899 
1900  // Loop through clusters, and plot M02 for each particle type
1902 
1903  // Loop through jets, to fill various histograms
1904  if (fPlotJetHistograms) {
1906  }
1907 
1908 }
1909 
1910 /*
1911  * Loop through clusters, and plot M02 for each particle type
1912  */
1914 {
1915  TString histname;
1917  const AliVCluster* clus;
1918  std::vector<ContributorType> vecContributorTypes;
1919  std::vector<Int_t> vecContributorLabels;
1920  for (auto it : clusters->accepted_momentum()) {
1921 
1922  clus = it.second;
1923 
1924  // Include only EMCal/DCal clusters (reject PHOS clusters)
1925  if (!clus->IsEMCAL()) {
1926  continue;
1927  }
1928 
1929  // Loop through the cluster's contributors in order to classify its type
1930  ParticleType particleType = kNotDefined;
1931  ContributorType contributorType = kUndefined;
1932  const Int_t nLabels = clus->GetNLabels();
1933 
1934  // Create a vector to store the contributor types for PhysicalPrimary particles
1935  vecContributorTypes.clear();
1936  vecContributorLabels.clear();
1937  for (Int_t iLabel=0; iLabel<nLabels; iLabel++) {
1938 
1939  Int_t label = clus->GetLabels()[iLabel];
1940  if (TMath::Abs(label) > 0) { // if the particle has a truth-level match, the label is nonzero
1941  contributorType = GetContributorType(clus, mcevent, label);
1942  if (contributorType != kUndefined) {
1943  vecContributorTypes.push_back(contributorType);
1944  vecContributorLabels.push_back(label);
1945  }
1946  }
1947  }
1948 
1949  Int_t nLabelsPhysPrim = vecContributorTypes.size();
1950 
1951  if (nLabelsPhysPrim == 1) {
1952 
1953  contributorType = vecContributorTypes[0];
1954 
1955  if (contributorType == kPhoton) {
1956  particleType = kSinglePhoton;
1957  }
1958  else if (contributorType == kElectron) {
1959  particleType = kSingleElectron;
1960  }
1961  else if (contributorType == kChargedPion) {
1962  particleType = kSingleChargedPion;
1963  }
1964  else if (contributorType == kProton) {
1965  particleType = kSingleProton;
1966  }
1967  else if (contributorType == kAntiProton) {
1968  particleType = kSingleAntiProton;
1969  }
1970  else if (contributorType == kChargedKaon) {
1971  particleType = kSingleChargedKaon;
1972  }
1973  else if (contributorType == kK0L) {
1974  particleType = kSingleK0L;
1975  }
1976  else if (contributorType == kNeutron) {
1977  particleType = kSingleNeutron;
1978  }
1979  else if (contributorType == kAntiNeutron) {
1980  particleType = kSingleAntiNeutron;
1981  }
1982  else {
1983  particleType = kSingleOther;
1984  }
1985 
1986  }
1987  else if (nLabelsPhysPrim == 2) {
1988 
1989  // Get the contributor particle types
1990  ContributorType contributorType1 = vecContributorTypes[0];
1991  ContributorType contributorType2 = vecContributorTypes[1];
1992 
1993  // Get the fraction of cluster energy from each contributor
1994  //Double_t frac0 = clus->GetClusterMCEdepFraction(0);
1995  Double_t frac1 = clus->GetClusterMCEdepFraction(1);
1996 
1997  // Check whether the leading/subleading contributors are photons/hadrons
1998  Bool_t isHadron1 = IsHadron(contributorType1);
1999  Bool_t isHadron2 = IsHadron(contributorType2);
2000  Bool_t isPhoton1 = contributorType1 == kPhoton;
2001  Bool_t isPhoton2 = contributorType2 == kPhoton;
2002 
2003  if (isHadron1 && isHadron2) {
2004  particleType = kHadronHadron;
2005  }
2006  else if (isPhoton1 && isHadron2) {
2007  particleType = kPhotonHadron;
2008 
2009  // Plot cluster energy when subleading hadron is subtracted
2010  Double_t photonEnergy = clus->GetNonLinCorrEnergy() * (1 - frac1);
2011  histname = "ClusterHistogramsMC/hPhotonHadronPhotonEnergy";
2012  fHistManager.FillTH3(histname.Data(), fCent, clus->GetM02(), photonEnergy);
2013  }
2014  else if (isHadron1 && isPhoton2) {
2015  particleType = kHadronPhoton;
2016 
2017  // Plot cluster energy when subleading hadron is subtracted
2018  Double_t hadronEnergy = clus->GetNonLinCorrEnergy() * (1 - frac1);
2019  histname = "ClusterHistogramsMC/hHadronPhotonHadronEnergy";
2020  fHistManager.FillTH3(histname.Data(), fCent, clus->GetM02(), hadronEnergy);
2021  }
2022  else if (isPhoton1 && isPhoton2) {
2023 
2024  // By default, assume the two photons are not a merged pi0
2025  particleType = kPhotonPhotonOther;
2026 
2027  // Using the vector of accepted contributor labels, check whether the two photons are from the same pi0
2028  AliAODMCParticle *part1 = fGeneratorLevel->GetMCParticleWithLabel(vecContributorLabels[0]);
2029  AliAODMCParticle *part2 = fGeneratorLevel->GetMCParticleWithLabel(vecContributorLabels[1]);
2030  if (part1 && part2) {
2031  Int_t iMother1 = part1->GetMother();
2032  Int_t iMother2 = part2->GetMother();
2033  AliVParticle *mother1 = mcevent->GetTrack(iMother1);
2034  AliVParticle *mother2 = mcevent->GetTrack(iMother2);
2035 
2036  if (mother1 && mother2) {
2037  if ( (mother1->PdgCode() == 111) && (mother2->PdgCode() == 111) ) {
2038  if (iMother1 == iMother2) {
2039  particleType = kMergedPi0;
2040  }
2041  }
2042  }
2043  }
2044  }
2045  else {
2046  particleType = kTwoContributorsOther; // this includes partially contained conversion overlaps
2047  }
2048 
2049  }
2050  else if (nLabelsPhysPrim > 2) {
2051  particleType = kMoreThanTwoContributors;
2052  }
2053 
2054  // (M02, Eclus, part type)
2055  if (fCent > 0 && fCent < 10) {
2056  histname = "ClusterHistogramsMC/hM02VsParticleTypeCentral";
2057  fHistManager.FillTH3(histname, clus->GetM02(), clus->GetNonLinCorrEnergy(), particleType);
2058  }
2059  if (fCent > 50 && fCent < 90) {
2060  histname = "ClusterHistogramsMC/hM02VsParticleTypePeripheral";
2061  fHistManager.FillTH3(histname, clus->GetM02(), clus->GetNonLinCorrEnergy(), particleType);
2062  }
2063 
2064  }
2065 }
2066 
2067 /*
2068  * Loop through jets, to fill particle composition histograms.
2069  */
2071 {
2072  TString histname;
2073 
2074  const AliVCluster* clus;
2075  AliJetContainer* jets = GetJetContainer(0); // there is only a single, det-level jet finder here
2076  for (const auto jet : jets->accepted()) {
2077 
2078  Double_t jetPt = GetJetPt(jet, 0);
2079 
2080  // Keep track of hadronic calo energy in each jet
2081  Double_t hadCaloEnergyUnmatched = 0;
2082  Double_t hadCaloEnergyMatchedNonlincorr = 0;
2083  Double_t hadCaloEnergyMatchedHadCorr = 0;
2084 
2085  // Loop through clusters in each jet, to plot several histograms
2086  Int_t nClusters = jet->GetNumberOfClusters();
2087  for (Int_t iClus = 0; iClus < nClusters; iClus++) {
2088 
2089  clus = jet->Cluster(iClus);
2090 
2091  // Get the particle type of the cluster
2092  ContributorType contributorType = kUndefined;
2093  Int_t label = TMath::Abs(clus->GetLabel());
2094  if (label > 0) {
2095  contributorType = GetContributorType(clus, mcevent, label);
2096  }
2097 
2098  // Plot M02 for each particle type
2099  histname = "JetPerformanceMC/hM02VsContributorTypeJets";
2100  Double_t x[4] = {clus->GetM02(), clus->GetNonLinCorrEnergy(), contributorType, jetPt};
2101  fHistManager.FillTHnSparse(histname, x);
2102 
2103  // If the cluster is a hadron, sum its energy to compute the jet's hadronic calo energy
2104  Bool_t isHadron = IsHadron(contributorType);
2105  if (isHadron) {
2106  Bool_t hasMatchedTrack = (clus->GetNTracksMatched() > 0);
2107  //Bool_t hasMatchedTrack = ((clus->GetNonLinCorrEnergy() - clus->GetHadCorrEnergy()) > 1e-3);
2108  if (hasMatchedTrack) {
2109  hadCaloEnergyMatchedNonlincorr += clus->GetNonLinCorrEnergy();
2110  hadCaloEnergyMatchedHadCorr += clus->GetHadCorrEnergy();
2111  }
2112  else {
2113  hadCaloEnergyUnmatched += clus->GetNonLinCorrEnergy();
2114  }
2115  }
2116 
2117  }
2118 
2119  // Fill hadronic calo energy in each jet
2120 
2121  // (Jet pT, Summed energy of hadronic clusters without a matched track)
2122  histname = "JetPerformance/hHadCaloEnergyUnmatched";
2123  fHistManager.FillTH2(histname, jetPt, hadCaloEnergyUnmatched);
2124 
2125  // (Jet pT vs. Summed energy of hadronic clusters with a matched track (before hadronic correction))
2126  histname = "JetPerformance/hHadCaloEnergyMatchedNonlincorr";
2127  fHistManager.FillTH2(histname, jetPt, hadCaloEnergyMatchedNonlincorr);
2128 
2129  // (Jet pT vs. Summed energy of hadronic clusters with a matched track (after hadronic correction))
2130  histname = "JetPerformance/hHadCaloEnergyMatchedHadCorr";
2131  fHistManager.FillTH2(histname, jetPt, hadCaloEnergyMatchedHadCorr);
2132 
2133  // Loop through particle types, and plot jet composition for each particle type
2134  histname = "JetPerformanceMC/hJetComposition";
2135  for (Int_t type = 0; type < 8; type++) {
2136 
2137  ContributorType contributorType = kUndefined;
2138  Double_t nSum = 0;
2139  Double_t pTsum = 0;
2140 
2141  // Loop through clusters in jet, and add to sum if cluster matches particle type
2142  for (Int_t iClus = 0; iClus < nClusters; iClus++) {
2143 
2144  clus = jet->Cluster(iClus);
2145 
2146  Int_t label = TMath::Abs(clus->GetLabel());
2147  if (label > 0) {
2148  contributorType = GetContributorType(clus, mcevent, label);
2149  }
2150 
2151  if (type == contributorType) {
2152  nSum++;
2153  pTsum += clus->GetNonLinCorrEnergy();
2154  }
2155  }
2156 
2157  // Fill jet composition histogram
2158  Double_t x[4] = {jetPt, 1.*type, nSum, pTsum};
2159  fHistManager.FillTHnSparse(histname, x);
2160 
2161  }
2162  }
2163 }
2164 
2169 {
2170  // Loop over tracks and clusters in order to:
2171  // (1) Compute scale factor for full jets
2172  // (2) Compute delta-pT for full jets, with the random cone method
2173 
2174  // Define the acceptance boundaries for the TPC and EMCal/DCal/PHOS
2175  Double_t etaTPC = 0.9;
2176  Double_t etaEMCal = 0.7;
2177  Double_t etaMinDCal = 0.22;
2178  Double_t phiMinEMCal = fGeom->GetArm1PhiMin() * TMath::DegToRad(); // 80
2179  Double_t phiMaxEMCal = fGeom->GetEMCALPhiMax() * TMath::DegToRad(); // ~188
2180  Double_t phiMinDCal = fGeom->GetDCALPhiMin() * TMath::DegToRad(); // 260
2181  Double_t phiMaxDCal = fGeom->GetDCALPhiMax() * TMath::DegToRad(); // ~327 (1/3 SMs start at 320)
2182 
2183  Double_t accTPC = 2 * etaTPC * 2 * TMath::Pi();
2184  Double_t accEMCal = 2 * etaEMCal * (phiMaxEMCal - phiMinEMCal);
2185  Double_t accDCal = 2 * (etaEMCal - etaMinDCal) * (phiMaxDCal - phiMinDCal);
2186 
2187  // Loop over jet containers
2188  AliJetContainer* jetCont = 0;
2189  TIter nextJetColl(&fJetCollArray);
2190  while ((jetCont = static_cast<AliJetContainer*>(nextJetColl()))) {
2191 
2192  // Define fiducial acceptances, to be used to generate random cones, and for scale factor studies
2193  TRandom3* r = new TRandom3(0);
2194  Double_t jetR = jetCont->GetJetRadius();
2195  Double_t etaEMCalfid = etaEMCal - jetR;
2196  Double_t etaMinDCalfid = etaMinDCal + jetR;
2197  Double_t phiMinEMCalfid = phiMinEMCal + jetR;
2198  Double_t phiMaxEMCalfid = phiMaxEMCal - jetR;
2199  Double_t phiMinDCalfid = phiMinDCal + jetR;
2200  Double_t phiMaxDCalfid = phiMaxDCal - jetR;
2201  Double_t accEMCalfid = 2 * etaEMCalfid * (phiMaxEMCalfid - phiMinEMCalfid);
2202  Double_t accDCalfid = 2 * (etaEMCalfid - etaMinDCalfid) * (phiMaxDCalfid - phiMinDCalfid);
2203  if ( (etaEMCalfid - etaMinDCalfid) < 0) {
2204  accDCalfid = 0;
2205  }
2206 
2207  // Generate EMCal random cone eta-phi
2208  Double_t etaEMCalRC = r->Uniform(-etaEMCalfid, etaEMCalfid);
2209  Double_t phiEMCalRC = r->Uniform(phiMinEMCalfid, phiMaxEMCalfid);
2210 
2211  // Generate DCal random cone eta-phi
2212  Double_t etaDCalRC = r->Uniform(etaMinDCalfid, etaEMCalfid);
2213  Double_t sign = r->Uniform(-1., 1.);
2214  if (sign < 0) {
2215  etaDCalRC = -1*etaDCalRC;
2216  }
2217  Double_t phiDCalRC = r->Uniform(phiMinDCalfid, phiMaxDCalfid);
2218 
2219  // Initialize the various sums to 0
2220  Double_t trackPtSumTPC = 0;
2221  Double_t trackPtSumEMCal = 0;
2222  Double_t trackPtSumEMCalfid = 0;
2223  Double_t trackPtSumDCal = 0;
2224  Double_t trackPtSumDCalfid = 0;
2225  Double_t trackPtSumEMCalRC = 0;
2226  Double_t trackPtSumDCalRC = 0;
2227  Double_t clusPtSumEMCal = 0;
2228  Double_t clusPtSumEMCalfid = 0;
2229  Double_t clusPtSumDCal = 0;
2230  Double_t clusPtSumDCalfid = 0;
2231  Double_t clusPtSumEMCalRC = 0;
2232  Double_t clusPtSumDCalRC = 0;
2233 
2234  // Loop over tracks. Sum the track pT:
2235  // (1) in the entire TPC, (2) in the EMCal, (3) in the EMCal fiducial volume, (4) in the DCal, (5) in the DCal fiducial volume, (6) in the EMCal random cone, (7) in the DCal random cone
2236  // Note: Loops over all det-level track containers. For data there should be only one. For embedding, there should be signal+background tracks.
2237  AliParticleContainer * partCont = 0;
2238  AliTLorentzVector track;
2239  Double_t trackEta;
2240  Double_t trackPhi;
2241  Double_t trackPt;
2242  Double_t deltaR;
2243  TIter nextPartCont(&fParticleCollArray);
2244  while ((partCont = static_cast<AliParticleContainer*>(nextPartCont()))) {
2245 
2246  TString partContName = partCont->GetName();
2247  if (!partContName.CompareTo("tracks") || !partContName.CompareTo("thermalparticles")) {
2248 
2249  for (auto trackIterator : partCont->accepted_momentum() ) {
2250 
2251  track.Clear();
2252  track = trackIterator.first;
2253  trackEta = track.Eta();
2254  trackPhi = track.Phi_0_2pi();
2255  trackPt = track.Pt();
2256 
2257  // (1)
2258  if (TMath::Abs(trackEta) < etaTPC) {
2259  trackPtSumTPC += trackPt;
2260  }
2261 
2262  // (2)
2263  if (TMath::Abs(trackEta) < etaEMCal && trackPhi > phiMinEMCal && trackPhi < phiMaxEMCal) {
2264  trackPtSumEMCal += trackPt;
2265  }
2266 
2267  // (3)
2268  if (TMath::Abs(trackEta) < etaEMCalfid && trackPhi > phiMinEMCalfid && trackPhi < phiMaxEMCalfid) {
2269  trackPtSumEMCalfid += trackPt;
2270  }
2271 
2272  // (4)
2273  if (fPlotDCal) {
2274  if (TMath::Abs(trackEta) > etaMinDCal && TMath::Abs(trackEta) < etaEMCal && trackPhi > phiMinDCal && trackPhi < phiMaxDCal) {
2275  trackPtSumDCal += trackPt;
2276  }
2277  }
2278 
2279  // (5)
2280  if (fPlotDCal) {
2281  if (TMath::Abs(trackEta) > etaMinDCalfid && TMath::Abs(trackEta) < etaEMCalfid && trackPhi > phiMinDCalfid && trackPhi < phiMaxDCalfid) {
2282  trackPtSumDCalfid += trackPt;
2283  }
2284  }
2285 
2286  // (6)
2287  deltaR = GetDeltaR(&track, etaEMCalRC, phiEMCalRC);
2288  if (deltaR < jetR) {
2289  trackPtSumEMCalRC += trackPt;
2290  }
2291 
2292  // (7)
2293  if (fPlotDCal) {
2294  deltaR = GetDeltaR(&track, etaDCalRC, phiDCalRC);
2295  if (deltaR < jetR) {
2296  trackPtSumDCalRC += trackPt;
2297  }
2298  }
2299  }
2300  }
2301  }
2302 
2303  // Loop over clusters, if the jet container is for full jets. Sum the cluster ET:
2304  // (1) in the EMCal, (2) in the EMCal fiducial volume, (3) in the DCal, (4), in the DCal fiducial volume, (5) in the EMCal random cone, (6) in the DCal random cone
2305  if (jetCont->GetClusterContainer()) {
2307 
2308  AliTLorentzVector clus;
2309  Double_t clusEta;
2310  Double_t clusPhi;
2311  Double_t clusPt;
2312  for (auto clusIterator : clusCont->accepted_momentum() ) {
2313 
2314  clus.Clear();
2315  clus = clusIterator.first;
2316  clusEta = clus.Eta();
2317  clusPhi = clus.Phi_0_2pi();
2318  clusPt = clus.Pt();
2319 
2320  // (1)
2321  if (TMath::Abs(clusEta) < etaEMCal && clusPhi > phiMinEMCal && clusPhi < phiMaxEMCal) {
2322  clusPtSumEMCal += clusPt;
2323  }
2324 
2325  // (2)
2326  if (TMath::Abs(clusEta) < etaEMCalfid && clusPhi > phiMinEMCalfid && clusPhi < phiMaxEMCalfid) {
2327  clusPtSumEMCalfid += clusPt;
2328  }
2329 
2330  // (3)
2331  if (fPlotDCal) {
2332  if (TMath::Abs(clusEta) > etaMinDCal && TMath::Abs(clusEta) < etaEMCal && clusPhi > phiMinDCal && clusPhi < phiMaxDCal) {
2333  clusPtSumDCal += clusPt;
2334  }
2335  }
2336 
2337  // (4)
2338  if (fPlotDCal) {
2339  if (TMath::Abs(clusEta) > etaMinDCalfid && TMath::Abs(clusEta) < etaEMCalfid && clusPhi > phiMinDCalfid && clusPhi < phiMaxDCalfid) {
2340  clusPtSumDCalfid += clusPt;
2341  }
2342  }
2343 
2344  // (5)
2345  deltaR = GetDeltaR(&clus, etaEMCalRC, phiEMCalRC);
2346  if (deltaR < jetR) {
2347  clusPtSumEMCalRC += clusPt;
2348  }
2349 
2350  // (6)
2351  if (fPlotDCal) {
2352  deltaR = GetDeltaR(&clus, etaDCalRC, phiDCalRC);
2353  if (deltaR < jetR) {
2354  clusPtSumDCalRC += clusPt;
2355  }
2356  }
2357 
2358  }
2359  }
2360 
2361  // Compute the scale factor, as a function of centrality, for (1) EMCal, (2) EMCalfid, (3) DCal, (4) DCalfid
2362  // (1)
2363  Double_t numerator = (trackPtSumEMCal + clusPtSumEMCal) / accEMCal;
2364  Double_t denominator = trackPtSumTPC / accTPC;
2365  Double_t scaleFactor = numerator / denominator;
2366  TString histname = TString::Format("%s/BackgroundHistograms/hScaleFactorEMCal", jetCont->GetArrayName().Data());
2367  fHistManager.FillTH2(histname, fCent, scaleFactor);
2368 
2369  // (2)
2370  if (accEMCalfid > 1e-3) {
2371  numerator = (trackPtSumEMCalfid + clusPtSumEMCalfid) / accEMCalfid;
2372  scaleFactor = numerator / denominator;
2373  histname = TString::Format("%s/BackgroundHistograms/hScaleFactorEMCalFid", jetCont->GetArrayName().Data());
2374  fHistManager.FillTH2(histname, fCent, scaleFactor);
2375  }
2376 
2377  // (3)
2378  if (fPlotDCal) {
2379  numerator = (trackPtSumDCal + clusPtSumDCal) / accDCal;
2380  scaleFactor = numerator / denominator;
2381  histname = TString::Format("%s/BackgroundHistograms/hScaleFactorDCal", jetCont->GetArrayName().Data());
2382  fHistManager.FillTH2(histname, fCent, scaleFactor);
2383  }
2384 
2385  // (4)
2386  if (fPlotDCal) {
2387  if (accDCalfid > 1e-3) {
2388  numerator = (trackPtSumDCalfid + clusPtSumDCalfid) / accDCalfid;
2389  scaleFactor = numerator / denominator;
2390  histname = TString::Format("%s/BackgroundHistograms/hScaleFactorDCalFid", jetCont->GetArrayName().Data());
2391  fHistManager.FillTH2(histname, fCent, scaleFactor);
2392  }
2393  }
2394 
2395  // Compute delta pT, as a function of centrality
2396 
2397  // EMCal
2398  Double_t rho = jetCont->GetRhoVal();
2399  Double_t deltaPt = trackPtSumEMCalRC + clusPtSumEMCalRC - rho * TMath::Pi() * jetR * jetR;
2400  histname = TString::Format("%s/BackgroundHistograms/hDeltaPtEMCal", jetCont->GetArrayName().Data());
2401  fHistManager.FillTH2(histname, fCent, deltaPt);
2402 
2403  // DCal
2404  if (fPlotDCal) {
2405  if (accDCalfid > 1e-3) {
2406  deltaPt = trackPtSumDCalRC + clusPtSumDCalRC - rho * TMath::Pi() * jetR * jetR;
2407  histname = TString::Format("%s/BackgroundHistograms/hDeltaPtDCal", jetCont->GetArrayName().Data());
2408  fHistManager.FillTH2(histname, fCent, deltaPt);
2409  }
2410  }
2411 
2412  delete r;
2413 
2414  }
2415 
2416 }
2417 
2423 {
2424  TString histname;
2425  AliJetContainer* jets = 0;
2426  TIter nextJetColl(&fJetCollArray);
2427  while ((jets = static_cast<AliJetContainer*>(nextJetColl()))) {
2428  TString jetContName = jets->GetName();
2429 
2430  Double_t rhoVal = 0;
2431  if (jets->GetRhoParameter()) {
2432  rhoVal = jets->GetRhoVal();
2433  histname = TString::Format("%s/TriggerSimHistograms/hRhoVsCent", jets->GetArrayName().Data());
2434  fHistManager.FillTH2(histname.Data(), fCent, rhoVal);
2435  }
2436 
2437  for (auto jet : jets->all()) {
2438 
2439  Float_t ptLeading = jets->GetLeadingHadronPt(jet);
2440  Float_t corrPt = GetJetPt(jet, rhoVal);
2441 
2442  // A vs. pT (fill before area cut)
2443  histname = TString::Format("%s/TriggerSimHistograms/hAreaVsPt", jets->GetArrayName().Data());
2444  fHistManager.FillTH2(histname.Data(), corrPt, jet->Area());
2445 
2446  // Rejection reason
2447  UInt_t rejectionReason = 0;
2448  if (!jets->AcceptJet(jet, rejectionReason)) {
2449  histname = TString::Format("%s/TriggerSimHistograms/hJetRejectionReason", jets->GetArrayName().Data());
2450  fHistManager.FillTH2(histname.Data(), jets->GetRejectionReasonBitPosition(rejectionReason), jet->Pt());
2451  continue;
2452  }
2453 
2454  // compute jet acceptance type
2455  Double_t type = GetJetType(jet);
2456  if ( type != kEMCal ) {
2457  if ( type != kDCal || !fPlotDCal ) {
2458  continue;
2459  }
2460  }
2461 
2462  // (Centrality, pT, NEF, calo type)
2463  if (type == kEMCal) {
2464  histname = TString::Format("%s/TriggerSimHistograms/hNEFVsPtEMCal", jets->GetArrayName().Data());
2465  }
2466  else if (type == kDCal) {
2467  histname = TString::Format("%s/TriggerSimHistograms/hNEFVsPtDCal", jets->GetArrayName().Data());
2468  }
2469  fHistManager.FillTH3(histname, fCent, corrPt, jet->NEF());
2470 
2471  // pT-leading vs. pT
2472  histname = TString::Format("%s/TriggerSimHistograms/hPtLeadingVsPt", jets->GetArrayName().Data());
2473  fHistManager.FillTH2(histname.Data(), corrPt, ptLeading);
2474 
2475  // (Centrality, pT, z-leading (charged), calo type)
2476  if (type == kEMCal) {
2477  histname = TString::Format("%s/TriggerSimHistograms/hZLeadingVsPtEMCal", jets->GetArrayName().Data());
2478  }
2479  else if (type == kDCal) {
2480  histname = TString::Format("%s/TriggerSimHistograms/hZLeadingVsPtDCal", jets->GetArrayName().Data());
2481  }
2482  TLorentzVector leadPart;
2483  jets->GetLeadingHadronMomentum(leadPart, jet);
2484  Double_t z = GetParallelFraction(leadPart.Vect(), jet);
2485  if (z == 1 || (z > 1 && z - 1 < 1e-3)) z = 0.999; // so that it will contribute to the bin <1
2486  fHistManager.FillTH3(histname, fCent, corrPt, z);
2487 
2488  // (Centrality, pT, z (charged), calo type)
2489  if (type == kEMCal) {
2490  histname = TString::Format("%s/TriggerSimHistograms/hZVsPtEMCal", jets->GetArrayName().Data());
2491  }
2492  else if (type == kDCal) {
2493  histname = TString::Format("%s/TriggerSimHistograms/hZVsPtDCal", jets->GetArrayName().Data());
2494  }
2495  const AliVTrack* track;
2496  for (Int_t i=0; i<jet->GetNumberOfTracks(); i++) {
2497  track = static_cast<AliVTrack*>(jet->Track(i));
2498  z = track->Pt() / TMath::Abs(corrPt);
2499  fHistManager.FillTH3(histname, fCent, corrPt, z);
2500  }
2501 
2502  // (Centrality, pT, Nconst)
2503  if (type == kEMCal) {
2504  histname = TString::Format("%s/TriggerSimHistograms/hNConstVsPtEMCal", jets->GetArrayName().Data());
2505  }
2506  else if (type == kDCal) {
2507  histname = TString::Format("%s/TriggerSimHistograms/hNConstVsPtDCal", jets->GetArrayName().Data());
2508  }
2509  fHistManager.FillTH3(histname, fCent, corrPt, 1.*jet->GetNumberOfConstituents());
2510 
2511  } //jet loop
2512 
2513  }
2514 }
2515 
2525 {
2526  // Loop over all jets and fill the ClosestJet(), i.e. the matching candidate.
2527  // Note: Allow truth jets to be outside of EMCALfid or fail 5 GeV requirement, since these can still contribute accepted det-jets
2528  // (but for the jet reconstruction efficiency denominator, the criteria should be enforced).
2531  }
2532  else if (fDoJetMatchingMCFraction) {
2533  ComputeJetMatches(fDetJetContainer, fDetJetContainerPPIntermediate, kTRUE); // First match PbPb-det to pp-det
2534  ComputeJetMatches(fDetJetContainerPPIntermediate, fMCJetContainer, kFALSE); // Then match pp-det to pp-truth
2535  }
2536 
2537  // Loop through accepted det-level jets, and retrieve matching candidate.
2538  // It match passes criteria (i.e. matching distance, uniqueness, MC fraction), fill matching histos.
2539  Double_t rhoVal = 0;
2541  rhoVal = fDetJetContainer->GetRhoVal();
2542  }
2543  for (auto jet : fDetJetContainer->accepted()) {
2544 
2545  Float_t detPt = GetJetPt(jet, rhoVal);
2546 
2547  // Get the matched part-level jet
2548  const AliEmcalJet* matchedPartLevelJet = GetMatchedPartLevelJet(jet, detPt);
2549  if (!matchedPartLevelJet) {
2550  continue;
2551  }
2552  Float_t truthPt = matchedPartLevelJet->Pt();
2553 
2554  // compute jet acceptance type
2555  Double_t type = GetJetType(jet);
2556  if ( type != kEMCal ) {
2557  if ( type != kDCal || !fPlotDCal ) {
2558  continue;
2559  }
2560  }
2561 
2562  // Fill response matrix (centrality, pT-truth, pT-det)
2563  TString histname;
2564  if (type == kEMCal) {
2565  histname = "MatchedJetHistograms/hResponseMatrixEMCal";
2566  }
2567  else if (type == kDCal) {
2568  histname = "MatchedJetHistograms/hResponseMatrixDCal";
2569  }
2570  if (fForceBeamType != kpp) {
2571  fHistManager.FillTH3(histname, fCent, truthPt, detPt);
2572  }
2573  else {
2574  fHistManager.FillTH2(histname, detPt, truthPt);
2575  }
2576 
2577  // Fill JES shift (centrality, pT-truth, (pT-det - pT-truth) / pT-truth)
2578  if (type == kEMCal) {
2579  histname = "MatchedJetHistograms/hJESshiftEMCal";
2580  }
2581  else if (type == kDCal) {
2582  histname = "MatchedJetHistograms/hJESshiftDCal";
2583  }
2584  if (fForceBeamType != kpp) {
2585  fHistManager.FillTH3(histname, fCent, truthPt, (detPt-truthPt)/truthPt );
2586  }
2587  else {
2588  fHistManager.FillTH2(histname, truthPt, (detPt-truthPt)/truthPt );
2589  }
2590 
2591  // Fill NEF of det-level matched jets (centrality, pT-truth, NEF)
2592  histname = "MatchedJetHistograms/hNEFVsPt";
2593  if (fForceBeamType != kpp) {
2594  fHistManager.FillTH3(histname, fCent, truthPt, jet->NEF());
2595  }
2596  else {
2597  fHistManager.FillTH2(histname, truthPt, jet->NEF());
2598  }
2599 
2600  // Fill z-leading (charged) of det-level matched jets (centrality, pT-truth, z-leading)
2601  histname = "MatchedJetHistograms/hZLeadingVsPt";
2602  TLorentzVector leadPart;
2604  Double_t z = GetParallelFraction(leadPart.Vect(), jet);
2605  if (z == 1 || (z > 1 && z - 1 < 1e-3)) z = 0.999; // so that it will contribute to the bin <1
2606  if (fForceBeamType != kpp) {
2607  fHistManager.FillTH3(histname, fCent, truthPt, z);
2608  }
2609  else {
2610  fHistManager.FillTH2(histname, truthPt, z);
2611  }
2612 
2613  } //jet loop
2614 }
2615 
2616 /*
2617  * Loop over jets of two specified jet collections, and fill the ClosestJet(), i.e. the matching candidate.
2618  * The first collection always uses the container acceptance criteria.
2619  * The second collection can be configured to use the container acceptance criteria or not.
2620  */
2622 
2623  for (auto jet1 : jetCont1->all()) {
2624  jet1->ResetMatching();
2625  }
2626  for (auto jet2 : jetCont2->all()) {
2627  jet2->ResetMatching();
2628  }
2629 
2630  for (auto jet1 : jetCont1->accepted()) {
2631 
2632  if (bUseJetCont2Acceptance) {
2633  for (auto jet2 : jetCont2->accepted()) {
2634  SetJetClosestCandidate(jet1, jet2);
2635  }
2636  }
2637  else {
2638  for (auto jet2 : jetCont2->all()) {
2639  if (jet2->Pt() < fMCJetMinMatchingPt) {
2640  continue;
2641  }
2642  SetJetClosestCandidate(jet1, jet2);
2643  }
2644  }
2645  }
2646 }
2647 
2648 /*
2649  * Given two jets, set them as closest if they are closer than the current closest jets.
2650  */
2652 
2653  Double_t deltaR = jet1->DeltaR(jet2);
2654  if (deltaR > 0.) {
2655 
2656  if (deltaR < jet1->ClosestJetDistance()) {
2657  jet1->SetClosestJet(jet2, deltaR);
2658  }
2659 
2660  if (deltaR < jet2->ClosestJetDistance()) {
2661  jet2->SetClosestJet(jet1, deltaR);
2662  }
2663  }
2664 }
2665 
2666 /*
2667  * Return a pointer to the matched truth-level jet, if it passes the matching criteria.
2668  * For fDoJetMatchingGeometrical, this means (1) within R = fJetMatchingR, (2) unique match.
2669  * For fDoJetMatchingMCFraction, this means also shared MC fraction requirement. That is, if the jet (combined jet)
2670  * is matched to a pp det-level jet, which is matched to a pp truth-level jet -- it must satisfy
2671  * (1) The shared momentum fraction being larger than some minimum value fMinSharedMomentumFraction, and
2672  * (2) Their matched distance being below the max matching distance fMaxMatchedJetDistance
2673  *
2674  * @param[in] detJet det-level jet to be checked for a successful truth-level match.
2675  * @return Pointer to truth-level matched jet, if it exists. False otherwise.
2676 */
2678 
2679  // Track in histogram how many matches pass each distinct matching criteria
2680  TString histNameQA = "MatchedJetHistograms/fHistJetMatchingQA";
2681 
2682  // Geometrical matching case (pp, p-Pb)
2684 
2685  bool returnValue = false;
2686  const AliEmcalJet* partLevelJet = detJet->ClosestJet();
2687  if (partLevelJet) {
2688  fHistManager.FillTH1(histNameQA, "matchedJet");
2689  returnValue = true;
2690 
2691  // Check if match is unique
2692  if (partLevelJet->ClosestJet() != detJet) {
2693  returnValue = false;
2694  }
2695  else {
2696  fHistManager.FillTH1(histNameQA, "uniqueMatch");
2697 
2698  // Fill matching distance between unique matches, without imposing deltaR cut (centrality, pT-truth, R)
2699  TString histname = "MatchedJetHistograms/hMatchingDistance";
2700  fHistManager.FillTH3(histname, detJetPt, partLevelJet->Pt(), detJet->ClosestJetDistance());
2701  }
2702 
2703  // Check if the matching distance cut is passed
2704  double matchedJetDistance = detJet->ClosestJetDistance();
2705  if (matchedJetDistance > fJetMatchingR) {
2706  returnValue = false;
2707  }
2708  else {
2709  fHistManager.FillTH1(histNameQA, "jetDistance");
2710  }
2711 
2712  // Record all cuts passed
2713  if (returnValue == true) {
2714  fHistManager.FillTH1(histNameQA, "passedAllCuts");
2715  }
2716  }
2717  else {
2718  fHistManager.FillTH1(histNameQA, "noMatch");
2719  returnValue = false;
2720  }
2721 
2722  if (returnValue) {
2723  return partLevelJet;
2724  }
2725  }
2726 
2727  // Shared MC fraction matching case (Pb-Pb)
2728  else if (fDoJetMatchingMCFraction) { // This function is essentially copied from AliAnalysisTaskEmcalJetHCorrelations::CheckForMatchedJet
2729 
2730  bool returnValue = false;
2731  const AliEmcalJet* partLevelJet = nullptr;
2732 
2733  // First, check if combined jet has a pp det-level match assigned
2734  if (detJet->ClosestJet()) {
2735  fHistManager.FillTH1(histNameQA, "matchedJet");
2736  returnValue = true;
2737 
2738  // Check shared momentum fraction.
2739  double sharedFraction = fDetJetContainer->GetFractionSharedPt(detJet, nullptr);
2740  if (sharedFraction < fMinSharedMomentumFraction) {
2741  returnValue = false;
2742  }
2743  else {
2744  fHistManager.FillTH1(histNameQA, "sharedMomentumFraction");
2745  }
2746 
2747  // Check that the combined jet has a particle-level match
2748  AliEmcalJet * detLevelJetPP = detJet->ClosestJet();
2749  partLevelJet = detLevelJetPP->ClosestJet();
2750  if (!partLevelJet) {
2751  returnValue = false;
2752  }
2753  else {
2754  fHistManager.FillTH1(histNameQA, "partLevelMatchedJet");
2755  }
2756 
2757  // Check the matching distance between the combined and pp det-level jets
2758  double matchedJetDistance = detJet->ClosestJetDistance();
2759  if (matchedJetDistance > fJetMatchingR) {
2760  returnValue = false;
2761  }
2762  else {
2763  fHistManager.FillTH1(histNameQA, "jetDistance");
2764  }
2765 
2766  // Record all cuts passed
2767  if (returnValue == true) {
2768  fHistManager.FillTH1(histNameQA, "passedAllCuts");
2769  }
2770 
2771  // Fill (det pT, shared MC fraction, deltaR) of closest jets
2772  TString histname = "MatchedJetHistograms/hMatchingDistanceVsMCFraction";
2773  fHistManager.FillTH3(histname, detJetPt, sharedFraction, matchedJetDistance);
2774 
2775  }
2776  else {
2777  fHistManager.FillTH1(histNameQA, "noMatch");
2778  returnValue = false;
2779  }
2780 
2781  if (returnValue) {
2782  return partLevelJet;
2783  }
2784  }
2785  return 0;
2786 }
2787 
2788 /*
2789  * Compute the MC particle type of a given cluster contributor, using the MC particle container
2790  */
2792 {
2793  ContributorType contributorType = kUndefined;
2794 
2795  AliAODMCParticle *part = fGeneratorLevel->GetMCParticleWithLabel(label);
2796  if (part) {
2797 
2798  TString histname = "ClusterHistogramsMC/hClusterRejectionReasonMC";
2799  UInt_t rejectionReason = 0;
2800  if (!fGeneratorLevel->AcceptMCParticle(part, rejectionReason)) {
2801  fHistManager.FillTH2(histname, fGeneratorLevel->GetRejectionReasonBitPosition(rejectionReason), clus->GetNonLinCorrEnergy());
2802  return contributorType;
2803  }
2804 
2805  if (part->GetGeneratorIndex() == 0) { // generator index in cocktail
2806 
2807  // select charged pions, protons, kaons, electrons, muons
2808  Int_t pdg = part->PdgCode();
2809 
2810  if (pdg == 22) { // gamma 22
2811  contributorType = kPhoton;
2812  }
2813  else if (TMath::Abs(pdg) == 211) { // pi+ 211 (abs value ensures both particles and antiparticles are included)
2814  contributorType = kChargedPion;
2815  }
2816  else if (pdg == 2212) { // proton 2212
2817  contributorType = kProton;
2818  }
2819  else if (pdg == -2212) {
2820  contributorType = kAntiProton;
2821  }
2822  else if (TMath::Abs(pdg) == 321) { // K+ 321
2823  contributorType = kChargedKaon;
2824  }
2825  else if (pdg == 130) { // K0L 130
2826  contributorType = kK0L;
2827  }
2828  else if (pdg == 2112) { // neutron 2112
2829  contributorType = kNeutron;
2830  }
2831  else if (pdg == -2112) {
2832  contributorType = kAntiNeutron;
2833  }
2834  else if (TMath::Abs(pdg) == 11) { // e- 11
2835  contributorType = kElectron;
2836  }
2837  else if (TMath::Abs(pdg) == 13) { // mu- 13
2838  contributorType = kMuon;
2839  }
2840  else {
2841  contributorType = kOther;
2842  }
2843  }
2844  }
2845  return contributorType;
2846 }
2847 
2852 {
2853  Double_t deltaPhi = TMath::Abs(part->Phi_0_2pi() - phiRef);
2854  Double_t deltaEta = TMath::Abs(part->Eta() - etaRef);
2855  Double_t deltaR = TMath::Sqrt( deltaPhi*deltaPhi + deltaEta*deltaEta );
2856  return deltaR;
2857 }
2858 
2863 {
2864  UInt_t jetType = jet->GetJetAcceptanceType();
2865  Double_t type = -1;
2866  if (jetType & AliEmcalJet::kEMCALfid) {
2867  type = kEMCal;
2868  }
2869  else if (jetType & AliEmcalJet::kDCALonlyfid) {
2870  type = kDCal;
2871  }
2872 
2873  return type;
2874 }
2875 
2880 {
2881  Double_t pT = jet->Pt() - rho * jet->Area();
2882  return pT;
2883 }
2884 
2889 {
2890  return (contributor == kChargedPion) || (contributor == kProton) || (contributor == kAntiProton) || (contributor == kChargedKaon) || (contributor == kK0L) || (contributor == kNeutron) || (contributor == kAntiNeutron);
2891 }
2892 
2897  const char *ntracks,
2898  const char *nclusters,
2899  const char *nGenLev,
2900  const Double_t minTrPt,
2901  const Double_t minClPt,
2902  const char *suffix)
2903 {
2904  // Get the pointer to the existing analysis manager via the static access method.
2905  //==============================================================================
2906  AliAnalysisManager *mgr = AliAnalysisManager::GetAnalysisManager();
2907  if (!mgr)
2908  {
2909  ::Error("AddTaskEmcalJetPerformance", "No analysis manager to connect to.");
2910  return 0;
2911  }
2912 
2913  // Check the analysis type using the event handlers connected to the analysis manager.
2914  //==============================================================================
2915  AliVEventHandler* handler = mgr->GetInputEventHandler();
2916  if (!handler)
2917  {
2918  ::Error("AddTaskEmcalJetPerformance", "This task requires an input event handler");
2919  return 0;
2920  }
2921 
2922  enum EDataType_t {
2923  kUnknown,
2924  kESD,
2925  kAOD
2926  };
2927 
2928  EDataType_t dataType = kUnknown;
2929 
2930  if (handler->InheritsFrom("AliESDInputHandler")) {
2931  dataType = kESD;
2932  }
2933  else if (handler->InheritsFrom("AliAODInputHandler")) {
2934  dataType = kAOD;
2935  }
2936 
2937  //-------------------------------------------------------
2938  // Init the task and do settings
2939  //-------------------------------------------------------
2940 
2941  TString trackName(ntracks);
2942  TString clusName(nclusters);
2943 
2944  if (trackName == "usedefault") {
2945  if (dataType == kESD) {
2946  trackName = "Tracks";
2947  }
2948  else if (dataType == kAOD) {
2949  trackName = "tracks";
2950  }
2951  else {
2952  trackName = "";
2953  }
2954  }
2955 
2956  if (clusName == "usedefault") {
2957  if (dataType == kESD) {
2958  clusName = "CaloClusters";
2959  }
2960  else if (dataType == kAOD) {
2961  clusName = "caloClusters";
2962  }
2963  else {
2964  clusName = "";
2965  }
2966  }
2967 
2968  TString name("AliAnalysisTaskEmcalJetPerformance");
2969  if (!trackName.IsNull()) {
2970  name += "_";
2971  name += trackName;
2972  }
2973  if (!clusName.IsNull()) {
2974  name += "_";
2975  name += clusName;
2976  }
2977  if (strcmp(suffix,"") != 0) {
2978  name += "_";
2979  name += suffix;
2980  }
2981 
2983  // Configure jet performance task
2985 
2987  // Create track and cluster containers with the standard cuts
2988 
2989  AliParticleContainer* partCont = 0;
2990  if (trackName == "mcparticles") {
2991  AliMCParticleContainer* mcpartCont = new AliMCParticleContainer(trackName);
2992  partCont = mcpartCont;
2993  }
2994  else if (trackName == "tracks" || trackName == "Tracks") {
2995  AliTrackContainer* trackCont = new AliTrackContainer(trackName);
2996  partCont = trackCont;
2997  }
2998  if (partCont) partCont->SetParticlePtCut(minTrPt);
2999  if (partCont) task->AdoptParticleContainer(partCont);
3000 
3001  // Add the generator-level container, if specified
3002  if (nGenLev && strcmp(nGenLev,"")!=0) {
3003  AliMCParticleContainer* mcpartCont = task->AddMCParticleContainer(nGenLev);
3004  mcpartCont->SelectPhysicalPrimaries(kTRUE);
3005  mcpartCont->SetParticlePtCut(0);
3006  }
3007 
3008  // Add the cluster container
3009  AliClusterContainer* clusCont = 0;
3010  if (!clusName.IsNull()) {
3011  clusCont = new AliClusterContainer(clusName);
3012  clusCont->SetClusECut(0.);
3013  clusCont->SetClusPtCut(0.);
3014  }
3015  if (clusCont) task->AdoptClusterContainer(clusCont);
3016 
3017  //-------------------------------------------------------
3018  // Final settings, pass to manager and set the containers
3019  //-------------------------------------------------------
3020 
3021  mgr->AddTask(task);
3022 
3023  // Create containers for input/output
3024  AliAnalysisDataContainer *cinput1 = mgr->GetCommonInputContainer() ;
3025  TString contname(name);
3026  contname += "_histos";
3027  AliAnalysisDataContainer *coutput1 = mgr->CreateContainer(contname.Data(),
3028  TList::Class(),AliAnalysisManager::kOutputContainer,
3029  Form("%s", AliAnalysisManager::GetCommonFileName()));
3030  mgr->ConnectInput (task, 0, cinput1 );
3031  mgr->ConnectOutput (task, 1, coutput1 );
3032 
3033  return task;
3034 }
3035 
Int_t pdg
static AliAnalysisTaskEmcalJetPerformance * AddTaskEmcalJetPerformance(const char *ntracks="usedefault", const char *nclusters="usedefault", const char *nGenLev="mcparticles", const Double_t minTrPt=0.15, const Double_t minClPt=0.30, const char *suffix="")
Double_t Area() const
Definition: AliEmcalJet.h:130
TObjArray fClusterCollArray
cluster collection array
void SetParticlePtCut(Double_t cut)
Double_t GetRhoVal() const
const TString & GetRhoName() const
double Double_t
Definition: External.C:58
TList * fEventCutList
! Output list for event cut histograms
Double_t GetJetPt(const AliEmcalJet *jet, Double_t rho)
const char * title
Definition: MakeQAPdf.C:27
AliEmcalJet * ClosestJet() const
Definition: AliEmcalJet.h:327
Bool_t fUseAliEventCuts
Flag to use AliEventCuts (otherwise AliAnalysisTaskEmcal will be used)
AliJetContainer * GetJetContainer(Int_t i=0) const
Definition: External.C:244
Double_t ClosestJetDistance() const
Definition: AliEmcalJet.h:328
UInt_t fOffTrigger
offline trigger for event selection
void AdoptParticleContainer(AliParticleContainer *cont)
const AliClusterIterableMomentumContainer accepted_momentum() const
const AliMCParticleIterableMomentumContainer accepted_momentum() const
Int_t nbinsy
void FillParticleCompositionJetHistograms(const AliMCEvent *mcevent)
Container with name, TClonesArray and cuts for particles.
void FillTH2(const char *hname, double x, double y, double weight=1., Option_t *opt="")
Fill a 2D histogram within the container.
Declaration of class AliTLorentzVector.
static AliEmcalDownscaleFactorsOCDB * Instance()
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.
Declaration of class AliAnalysisTaskEmcalEmbeddingHelper.
AliMCParticleContainer * fGeneratorLevel
! generator level container
virtual Bool_t AcceptJet(Int_t i, UInt_t &rejectionReason) const
Double_t GetDownscaleFactorForTriggerClass(const TString &trigger) const
Container for particles within the EMCAL framework.
void SetCaloTriggerPatchInfoName(const char *n)
Bool_t fPlotClusterHistograms
Set whether to plot cluster histograms.
Double_t GetDeltaR(const AliTLorentzVector *part, Double_t etaRef, Double_t phiRef)
void AdoptClusterContainer(AliClusterContainer *cont)
TObjArray fParticleCollArray
particle/track collection array
virtual Bool_t AcceptMCParticle(const AliAODMCParticle *vp, UInt_t &rejectionReason) const
bool AddQAPlotsToList(TList *list)
Bool_t fDoJetMatchingGeometrical
Do geometrical matching between det-level and truth-level jet container.
AliParticleContainer * GetParticleContainer(Int_t i=0) const
Get particle container attached to this task.
const Int_t nPtBins
Bool_t fComputeMBDownscaling
Set whether to compute and plot MB downscaling factors.
void GetLeadingHadronMomentum(TLorentzVector &mom, const AliEmcalJet *jet) const
AliEmcalEmbeddingQA fEmbeddingQA
! QA hists for embedding (will only be added if embedding)
void FillTHnSparse(const char *name, const double *x, double weight=1., Option_t *opt="")
bool IsInitialized() 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.
int Int_t
Definition: External.C:63
AliJetContainer * fMCJetContainer
!Pointer to jet container of truth-level jets
Double_t fTrackMatchingDeltaEtaMax
Maximum delta-eta to consider a track to be matched to a cluster.
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
UInt_t GetJetAcceptanceType() const
Definition: AliEmcalJet.h:367
Double_t Phi_0_2pi() const
Double_t fJetMatchingR
Jet matching R threshold.
void ComputeJetMatches(AliJetContainer *jetCont1, AliJetContainer *jetCont2, Bool_t bUseJetCont2Acceptance)
Implementation of task to embed external events.
TH1 * CreateTH1(const char *name, const char *title, int nbins, double xmin, double xmax, Option_t *opt="")
Create a new TH1 within the container.
Double_t fTrackMatchingDeltaPhiMax
Maximum delta-phi to consider a track to be matched to a cluster.
Photon+Photon (not from same pi0) are the only contributors.
BeamType fForceBeamType
forced beam type
Bool_t fUseManualEventCuts
Flag to use manual event cuts.
AliClusterContainer * GetClusterContainer(Int_t i=0) const
Get cluster container attached to this task.
Double_t fCent
!event centrality
Photon (direct or decay) is the only contributor.
Two particles from merged pi0 are the only contributors.
Bool_t fPlotJetHistograms
Set whether to enable inclusive jet histograms.
Int_t fNEoverPBins
! number of variable E/p bins
AliMCParticleContainer * AddMCParticleContainer(const char *n)
Create new container for MC particles and attach it to the task.
static Double_t GetParallelFraction(AliVParticle *part1, AliVParticle *part2)
Calculates the fraction of momentum z of part 1 w.r.t. part 2 in the direction of part 2...
void FillTH1(const char *hname, double x, double weight=1., Option_t *opt="")
Fill a 1D histogram within the container.
AliEventCuts fEventCuts
event selection utility
Bool_t fDoJetMatchingMCFraction
Do MC-fraction based matching using PbPb det-level, pp det-level, and pp truth-level jet containers...
Bool_t fPlotMatchedJetHistograms
Set whether to plot matched jet histograms.
TObjArray fJetCollArray
jet collection array
Double_t DeltaR(const AliVParticle *part) const
void SetJetClosestCandidate(AliEmcalJet *jet1, AliEmcalJet *jet2)
Bool_t fkEMCEJE
! flag telling whether the event is "triggered" or not in "simulation"
AliRhoParameter * GetRhoParameter()
Double_t Pt() const
Definition: AliEmcalJet.h:109
const AliEmcalJet * GetMatchedPartLevelJet(const AliEmcalJet *jet, Double_t detJetPt)
Photon+Hadron are the only contributors, with Photon leading.
Int_t nbinsx
Double_t fMCJetMinMatchingPt
Min jet pT for MC jets being matched, for when container criteria is not applied. ...
virtual Bool_t IsEventSelected()
Performing event selection.
Bool_t IsHadron(const ContributorType contributor)
static Double_t * GenerateFixedBinArray(Int_t n, Double_t min, Double_t max)
Bool_t fComputeBackground
Set whether to enable study of background.
Float_t GetJetRadius() const
Double_t fMinSharedMomentumFraction
Minimum shared momentum (pp det-level track pT in combined jet) / (pp det-level track pT) ...
AliEmcalList * fOutput
!output list
AliJetContainer * fDetJetContainerPPIntermediate
!Pointer to jet container of intermediate pp det-level jets, if MC-fraction matching ...
ContributorType GetContributorType(const AliVCluster *clus, const AliMCEvent *mcevent, Int_t label)
Handler for downscale factors for various triggers obtained from the OCDB.
void SetClosestJet(AliEmcalJet *j, Double_t d)
Definition: AliEmcalJet.h:323
AliMCParticleContainer * GetMCParticleContainer(Int_t i=0) const
Int_t fNPtHistBins
! number of variable pt bins
Definition: External.C:220
Bool_t fPlotParticleCompositionHistograms
Set whether to plot jet composition histograms.
const AliParticleIterableMomentumContainer accepted_momentum() const
void SelectPhysicalPrimaries(Bool_t s)
Bool_t fDoClosureTest
Set whether to do thermal model closure test.
TClonesArray * fTriggerPatchInfo
!trigger patch info array
Bool_t fDoTriggerResponse
flag whether to compute max patch response, in case of MC
Base task in the EMCAL jet framework.
Represent a jet reconstructed using the EMCal jet framework.
Definition: AliEmcalJet.h:51
Int_t fNM02HistBins
! number of variable M02 bins
Double_t fMBUpscaleFactor
! inverse of downscale factor, for MB trigger
AliJetContainer * fDetJetContainer
!Pointer to jet container of det-level jets
Double_t GetFractionSharedPt(const AliEmcalJet *jet, AliParticleContainer *cont2=0x0) const
Int_t GetRunNumber(TString)
Definition: PlotMuonQA.C:2235
void SetClusPtCut(Double_t cut)
void SetRejectionReasonLabels(TAxis *axis)
void UserCreateOutputObjects()
Main initialization function on the worker.
Bool_t fDoTriggerSimulation
Set whether to perform a simple trigger simulation.
const Int_t nbins
const AliJetIterableContainer accepted() const
bool Bool_t
Definition: External.C:53
EDataType_t
Switch for the data type.
void SetClusECut(Double_t cut)
Double_t fMedianDCal
! median patch energy in DCal, per event
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.
Double_t fMedianEMCal
! median patch energy in EMCal, per event
Container structure for EMCAL clusters.
DCal fiducial acceptance (each eta, phi edge narrowed by jet R)
Definition: AliEmcalJet.h:74
Container for MC-true particles within the EMCAL framework.
virtual AliAODMCParticle * GetMCParticleWithLabel(Int_t lab) const
Bool_t fNeedEmcalGeom
whether or not the task needs the emcal geometry
EMCal fiducial acceptance (each eta, phi edge narrowed by jet R)
Definition: AliEmcalJet.h:70
Container for jet within the EMCAL jet framework.
Hadron+Photon are the only contributors, with Hadron leading.
void FillParticleCompositionClusterHistograms(const AliMCEvent *mcevent)
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.
Bool_t fPlotDCal
Set whether to enable several DCal-specific histograms.
Double_t fPlotJetMatchCandThresh
Threshold for jet R to count candidates, affects plotting only.
const AliJetIterableContainer all() const
static const AliAnalysisTaskEmcalEmbeddingHelper * GetInstance()