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