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