AliPhysics  58eb684 (58eb684)
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
AliAnalysisTaskEmcalJetPerformance.cxx
Go to the documentation of this file.
1 /**********************************************************************************
2 * Copyright (C) 2016, Copyright Holders of the ALICE Collaboration *
3 * All rights reserved. *
4 * *
5 * Redistribution and use in source and binary forms, with or without *
6 * modification, are permitted provided that the following conditions are met: *
7 * * Redistributions of source code must retain the above copyright *
8 * notice, this list of conditions and the following disclaimer. *
9 * * Redistributions in binary form must reproduce the above copyright *
10 * notice, this list of conditions and the following disclaimer in the *
11 * documentation and/or other materials provided with the distribution. *
12 * * Neither the name of the <organization> nor the *
13 * names of its contributors may be used to endorse or promote products *
14 * derived from this software without specific prior written permission. *
15 * *
16 * THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS "AS IS" AND *
17 * ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE IMPLIED *
18 * WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE ARE *
19 * DISCLAIMED. IN NO EVENT SHALL ALICE COLLABORATION BE LIABLE FOR ANY *
20 * DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES *
21 * (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; *
22 * LOSS OF USE, DATA, OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND *
23 * ON ANY THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT *
24 * (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS *
25 * SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE. *
26 * *********************************************************************************/
27 
28 #include <vector>
29 
30 #include <TClonesArray.h>
31 #include <TH1F.h>
32 #include <TH2F.h>
33 #include <TH3.h>
34 #include <TList.h>
35 #include <THnSparse.h>
36 #include <TRandom3.h>
37 #include <TGrid.h>
38 #include <TFile.h>
39 
40 #include <AliVCluster.h>
41 #include <AliVParticle.h>
42 #include <AliLog.h>
43 
44 #include "AliAnalysisManager.h"
45 #include <AliVEventHandler.h>
46 #include "AliTLorentzVector.h"
47 #include "AliEmcalJet.h"
48 #include "AliRhoParameter.h"
49 #include "AliJetContainer.h"
50 #include "AliParticleContainer.h"
51 #include "AliClusterContainer.h"
52 #include "AliEMCALGeometry.h"
53 #include "AliOADBContainer.h"
55 #include "AliEMCALTriggerPatchInfo.h"
57 #include "AliMCEvent.h"
58 
60 
64 
70  fPlotJetHistograms(kFALSE),
71  fPlotClusterHistograms(kFALSE),
72  fPlotParticleCompositionHistograms(kFALSE),
73  fComputeBackground(kFALSE),
74  fDoTriggerSimulation(kFALSE),
75  fPlotMatchedJetHistograms(kFALSE),
76  fComputeMBDownscaling(kFALSE),
77  fMaxPt(200),
78  fNEtaBins(40),
79  fNPhiBins(200),
80  fNCentHistBins(0),
81  fCentHistBins(0),
82  fNPtHistBins(0),
83  fPtHistBins(0),
84  fNM02HistBins(0),
85  fM02HistBins(0),
86  fNEoverPBins(0),
87  fEoverPBins(0),
88  fTrackMatchingDeltaEtaMax(0.015),
89  fTrackMatchingDeltaPhiMax(0.030),
90  fMBUpscaleFactor(1.),
91  fMedianEMCal(0.),
92  fMedianDCal(0.),
93  fkEMCEJE(kFALSE),
94  fEmbeddingQA(),
95  fUseAliEventCuts(kTRUE),
96  fEventCuts(0),
97  fEventCutList(0),
98  fUseManualEventCuts(kFALSE),
99  fGeneratorLevel(0),
100  fHistManager()
101 {
103 }
104 
111  AliAnalysisTaskEmcalJet(name, kTRUE),
112  fPlotJetHistograms(kFALSE),
113  fPlotClusterHistograms(kFALSE),
114  fPlotParticleCompositionHistograms(kFALSE),
115  fComputeBackground(kFALSE),
116  fDoTriggerSimulation(kFALSE),
117  fPlotMatchedJetHistograms(kFALSE),
118  fComputeMBDownscaling(kFALSE),
119  fMaxPt(200),
120  fNEtaBins(40),
121  fNPhiBins(200),
122  fNCentHistBins(0),
123  fCentHistBins(0),
124  fNPtHistBins(0),
125  fPtHistBins(0),
126  fNM02HistBins(0),
127  fM02HistBins(0),
128  fNEoverPBins(0),
129  fEoverPBins(0),
130  fTrackMatchingDeltaEtaMax(0.015),
131  fTrackMatchingDeltaPhiMax(0.030),
132  fMBUpscaleFactor(1.),
133  fMedianEMCal(0.),
134  fMedianDCal(0.),
135  fkEMCEJE(kFALSE),
136  fEmbeddingQA(),
137  fUseAliEventCuts(kTRUE),
138  fEventCuts(0),
139  fEventCutList(0),
140  fUseManualEventCuts(kFALSE),
141  fGeneratorLevel(0),
142  fHistManager(name)
143 {
145 }
146 
151 {
152 }
153 
158 {
159  fNCentHistBins = 4;
161  fCentHistBins[0] = 0;
162  fCentHistBins[1] = 10;
163  fCentHistBins[2] = 30;
164  fCentHistBins[3] = 50;
165  fCentHistBins[4] = 90;
166 
167  fNPtHistBins = 82;
170  GenerateFixedBinArray(7, 0.3, 1, fPtHistBins+6);
171  GenerateFixedBinArray(10, 1, 3, fPtHistBins+13);
172  GenerateFixedBinArray(14, 3, 10, fPtHistBins+23);
173  GenerateFixedBinArray(10, 10, 20, fPtHistBins+37);
174  GenerateFixedBinArray(15, 20, 50, fPtHistBins+47);
175  GenerateFixedBinArray(20, 50, 150, fPtHistBins+62);
176 
177  fNM02HistBins = 81;
179  GenerateFixedBinArray(35, 0, 0.7, fM02HistBins);
180  GenerateFixedBinArray(6, 0.7, 1., fM02HistBins+35);
181  GenerateFixedBinArray(20, 1., 3., fM02HistBins+41);
182  GenerateFixedBinArray(10, 3., 5., fM02HistBins+61);
183  GenerateFixedBinArray(10, 5., 10., fM02HistBins+71);
184 
185  fNEoverPBins = 47;
187  GenerateFixedBinArray(30, 0, 1.5, fEoverPBins);
188  GenerateFixedBinArray(10, 1.5, 3.5, fEoverPBins+30);
189  GenerateFixedBinArray(7, 3.5, 10.5, fEoverPBins+40);
190 }
191 
197 {
199 
200  // Intialize AliEventCuts
201  if (fUseAliEventCuts) {
202  fEventCutList = new TList();
203  fEventCutList ->SetOwner();
204  fEventCutList ->SetName("EventCutOutput");
205 
206  fEventCuts.OverrideAutomaticTriggerSelection(fOffTrigger);
207  if(fUseManualEventCuts==1)
208  {
209  fEventCuts.SetManualMode();
210  // Configure manual settings here
211  // ...
212  }
213  fEventCuts.AddQAplotsToList(fEventCutList);
214  fOutput->Add(fEventCutList);
215  }
216 
217  // Get the MC particle branch, in case it exists
218  fGeneratorLevel = GetMCParticleContainer("mcparticles");
219 
220  // Allocate histograms
221  if (fPlotJetHistograms) {
223  }
226  }
229  }
230  if (fComputeBackground) {
232  }
233  if (fDoTriggerSimulation) {
235  }
238  }
239 
240  // Initialize embedding QA
242  if (embeddingHelper) {
243  bool res = fEmbeddingQA.Initialize();
244  if (res) {
246  }
247  }
248 
249  TIter next(fHistManager.GetListOfHistograms());
250  TObject* obj = 0;
251  while ((obj = next())) {
252  fOutput->Add(obj);
253  }
254 
255  PostData(1, fOutput); // Post data for ALL output slots > 0 here.
256 }
257 
258 /*
259  * This function allocates the histograms for single jets.
260  * A set of histograms is allocated per each jet container.
261  */
263 {
264  TString histname;
265  TString title;
266 
267  Int_t nPtBins = TMath::CeilNint(fMaxPt/2);
268 
269  AliJetContainer* jets = 0;
270  TIter nextJetColl(&fJetCollArray);
271  while ((jets = static_cast<AliJetContainer*>(nextJetColl()))) {
272 
273  // Jet rejection reason
274  histname = TString::Format("%s/JetHistograms/hJetRejectionReason", jets->GetArrayName().Data());
275  title = histname + ";Rejection reason;#it{p}_{T,jet} (GeV/#it{c});counts";
276  TH2* hist = fHistManager.CreateTH2(histname.Data(), title.Data(), 32, 0, 32, 50, 0, fMaxPt);
277  SetRejectionReasonLabels(hist->GetXaxis());
278 
279  // Rho vs. Centrality
280  if (!jets->GetRhoName().IsNull()) {
281  histname = TString::Format("%s/JetHistograms/hRhoVsCent", jets->GetArrayName().Data());
282  title = histname + ";Centrality (%);#rho (GeV/#it{c});counts";
283  fHistManager.CreateTH2(histname.Data(), title.Data(), 50, 0, 100, 100, 0, 500);
284  }
285 
286  // (Centrality, pT, NEF)
287  Int_t nbinsx = 20; Int_t minx = 0; Int_t maxx = 100;
288  Int_t nbinsy = nPtBins; Int_t miny = 0; Int_t maxy = fMaxPt;
289  Int_t nbinsz = 50; Int_t minz = 0; Int_t maxz = 1.;
290 
291  histname = TString::Format("%s/JetHistograms/hNEFVsPtEMCal", jets->GetArrayName().Data());
292  title = histname + ";Centrality (%);#it{p}_{T}^{corr} (GeV/#it{c});NEF";
293  fHistManager.CreateTH3(histname.Data(), title.Data(), nbinsx, minx, maxx, nbinsy, miny, maxy, nbinsz, minz, maxz);
294 
295  histname = TString::Format("%s/JetHistograms/hNEFVsPtDCal", jets->GetArrayName().Data());
296  title = histname + ";Centrality (%);#it{p}_{T}^{corr} (GeV/#it{c});NEF";
297  fHistManager.CreateTH3(histname.Data(), title.Data(), nbinsx, minx, maxx, nbinsy, miny, maxy, nbinsz, minz, maxz);
298 
299  // (Centrality, pT upscaled, calo type)
300  if (fComputeMBDownscaling) {
301  histname = TString::Format("%s/JetHistograms/hPtUpscaledMB", jets->GetArrayName().Data());
302  title = histname + ";Centrality (%);#it{p}_{T}^{corr} (GeV/#it{c});type";
303  fHistManager.CreateTH3(histname.Data(), title.Data(), 20, 0, 100, nPtBins, 0, fMaxPt, 2, -0.5, 1.5, "s");
304  }
305 
306  // pT-leading vs. pT
307  histname = TString::Format("%s/JetHistograms/hPtLeadingVsPt", jets->GetArrayName().Data());
308  title = histname + ";#it{p}_{T}^{corr} (GeV/#it{c});#it{p}_{T,particle}^{leading} (GeV/#it{c})";
309  fHistManager.CreateTH2(histname.Data(), title.Data(), nPtBins, 0, fMaxPt, nPtBins, 0, fMaxPt);
310 
311  // A vs. pT
312  histname = TString::Format("%s/JetHistograms/hAreaVsPt", jets->GetArrayName().Data());
313  title = histname + ";#it{p}_{T}^{corr} (GeV/#it{c});#it{A}_{jet}";
314  fHistManager.CreateTH2(histname.Data(), title.Data(), nPtBins, 0, fMaxPt, fMaxPt/3, 0, 0.5);
315 
316  // (Centrality, pT, z-leading (charged))
317  nbinsx = 20; minx = 0; maxx = 100;
318  nbinsy = nPtBins; miny = 0; maxy = fMaxPt;
319  nbinsz = 50; minz = 0; maxz = 1.;
320 
321  histname = TString::Format("%s/JetHistograms/hZLeadingVsPtEMCal", jets->GetArrayName().Data());
322  title = histname + ";Centrality (%);#it{p}_{T}^{corr} (GeV/#it{c});#it{z}_{leading}";
323  fHistManager.CreateTH3(histname.Data(), title.Data(), nbinsx, minx, maxx, nbinsy, miny, maxy, nbinsz, minz, maxz);
324 
325  histname = TString::Format("%s/JetHistograms/hZLeadingVsPtDCal", jets->GetArrayName().Data());
326  title = histname + ";Centrality (%);#it{p}_{T}^{corr} (GeV/#it{c});#it{z}_{leading}";
327  fHistManager.CreateTH3(histname.Data(), title.Data(), nbinsx, minx, maxx, nbinsy, miny, maxy, nbinsz, minz, maxz);
328 
329  // (Centrality, pT, z (charged))
330  nbinsx = 20; minx = 0; maxx = 100;
331  nbinsy = nPtBins; miny = 0; maxy = fMaxPt;
332  nbinsz = 50; minz = 0; maxz = 1.;
333 
334  histname = TString::Format("%s/JetHistograms/hZVsPtEMCal", jets->GetArrayName().Data());
335  title = histname + ";Centrality (%);#it{p}_{T}^{corr} (GeV/#it{c});#it{z}";
336  fHistManager.CreateTH3(histname.Data(), title.Data(), nbinsx, minx, maxx, nbinsy, miny, maxy, nbinsz, minz, maxz);
337 
338  histname = TString::Format("%s/JetHistograms/hZVsPtDCal", jets->GetArrayName().Data());
339  title = histname + ";Centrality (%);#it{p}_{T}^{corr} (GeV/#it{c});#it{z}";
340  fHistManager.CreateTH3(histname.Data(), title.Data(), nbinsx, minx, maxx, nbinsy, miny, maxy, nbinsz, minz, maxz);
341 
342  // (Centrality, pT, Nconst, calo type)
343  nbinsx = 20; minx = 0; maxx = 100;
344  nbinsy = nPtBins; miny = 0; maxy = fMaxPt;
345  nbinsz = 50; minz = 0; maxz = fMaxPt;
346 
347  histname = TString::Format("%s/JetHistograms/hNConstVsPtEMCal", jets->GetArrayName().Data());
348  title = histname + ";Centrality (%);#it{p}_{T}^{corr} (GeV/#it{c});No. of constituents";
349  fHistManager.CreateTH3(histname.Data(), title.Data(), nbinsx, minx, maxx, nbinsy, miny, maxy, nbinsz, minz, maxz);
350 
351  histname = TString::Format("%s/JetHistograms/hNConstVsPtDCal", jets->GetArrayName().Data());
352  title = histname + ";Centrality (%);#it{p}_{T}^{corr} (GeV/#it{c});No. of constituents";
353  fHistManager.CreateTH3(histname.Data(), title.Data(), nbinsx, minx, maxx, nbinsy, miny, maxy, nbinsz, minz, maxz);
354 
355  // (Centrality, jet pT, Enonlincorr - Ehadcorr)
356  nbinsx = 20; minx = 0; maxx = 100;
357  nbinsy = nPtBins; miny = 0; maxy = fMaxPt;
358  nbinsz = nPtBins; minz = 0; maxz = fMaxPt;
359 
360  histname = TString::Format("%s/JetHistograms/hDeltaEHadCorr", jets->GetArrayName().Data());
361  title = histname + ";Centrality (%);#it{p}_{T}^{corr} (GeV/#it{c});#sum#it{E}_{nonlincorr} - #it{E}_{hadcorr}";
362  fHistManager.CreateTH3(histname.Data(), title.Data(), nbinsx, minx, maxx, nbinsy, miny, maxy, nbinsz, minz, maxz);
363 
364  // (Median patch energy, calo type, jet pT, centrality)
365  if (fDoTriggerSimulation) {
366  histname = TString::Format("%s/JetHistograms/hMedPatchJet", jets->GetArrayName().Data());
367  title = histname + ";#it{E}_{patch,med};type;#it{p}_{T}^{corr} (GeV/#it{c});Centrality (%)";
368  Int_t nbins5[4] = {100, 2, nPtBins, 50};
369  Double_t min5[4] = {0,-0.5, 0, 0};
370  Double_t max5[4] = {50,1.5, fMaxPt, 100};
371  fHistManager.CreateTHnSparse(histname.Data(), title.Data(), 4, nbins5, min5, max5);
372  }
373 
374  }
375 
376  // MB downscale factor histogram
377  if (fComputeMBDownscaling) {
378  histname = "Trigger/hMBDownscaleFactor";
379  title = histname + ";Downscale factor;counts";
380  fHistManager.CreateTH1(histname.Data(), title.Data(), 200, 0, 200);
381  }
382 
383 }
384 
385 /*
386  * This function allocates the histograms for the calorimeter performance study.
387  */
389 {
390  TString histname;
391  TString htitle;
392 
393  const Int_t nRcorrBins = 50;
394  Double_t *RcorrBins = GenerateFixedBinArray(nRcorrBins, 0., 1.);
395  const Int_t nCellBins = 30;
396  Double_t *cellBins = GenerateFixedBinArray(nCellBins, -0.5, 29.5);
397  const Int_t nMatchedTrackBins = 5;
398  Double_t *matchedTrackBins = GenerateFixedBinArray(nMatchedTrackBins, -0.5, 4.5);
399 
402 
403  // Plot M02 distribution (centrality, Eclus nonlincorr, M02)
404  histname = "ClusterHistograms/hM02";
405  htitle = histname + ";Centrality (%);#it{E}_{clus} (GeV); M02";
407 
408  // Plot Ncell distribution for M02 > 0.4 and 0.1 < M02 < 0.4 (centrality, Eclus nonlincorr, Ncells)
409  histname = "ClusterHistograms/hNcellsM02G04";
410  htitle = histname + ";Centrality (%);#it{E}_{clus} (GeV); Ncells";
411  fHistManager.CreateTH3(histname.Data(), htitle.Data(), fNCentHistBins, fCentHistBins, fNPtHistBins, fPtHistBins, nCellBins, cellBins);
412 
413  histname = "ClusterHistograms/hNcellsM02L04";
414  htitle = histname + ";Centrality (%);#it{E}_{clus} (GeV); Ncells";
415  fHistManager.CreateTH3(histname.Data(), htitle.Data(), fNCentHistBins, fCentHistBins, fNPtHistBins, fPtHistBins, nCellBins, cellBins);
416 
419 
420  // Plot matched track pT for all clusters, M02 > 0.4 clusters, and 0.1 < M02 < 0.4 clusters (centrality, Eclus nonlincorr, trackPsum)
421  histname = "ClusterHistograms/hMatchedTrackPt";
422  htitle = histname + ";Centrality (%);#it{E}_{clus} (GeV); #Sigma#it{p}_{track} (GeV/c)";
424 
425  histname = "ClusterHistograms/hMatchedTrackPtM02G04";
426  htitle = histname + ";Centrality (%);#it{E}_{clus} (GeV); #Sigma#it{p}_{track} (GeV/c)";
428 
429  histname = "ClusterHistograms/hMatchedTrackPtM02L04";
430  htitle = histname + ";Centrality (%);#it{E}_{clus} (GeV); #Sigma#it{p}_{track} (GeV/c)";
432 
433  // Plot number of matched tracks for all clusters, M02 > 0.4 clusters, and 0.1 < M02 < 0.4 clusters (centrality, Eclus nonlincorr, N matches)
434  histname = "ClusterHistograms/hMatchedTrackN";
435  htitle = histname + ";Centrality (%);#it{E}_{clus} (GeV); N_{tracks}";
436  fHistManager.CreateTH3(histname.Data(), htitle.Data(), fNCentHistBins, fCentHistBins, fNPtHistBins, fPtHistBins, nMatchedTrackBins, matchedTrackBins);
437 
438  histname = "ClusterHistograms/hMatchedTrackNM02G04";
439  htitle = histname + ";Centrality (%);#it{E}_{clus} (GeV); N_{tracks}";
440  fHistManager.CreateTH3(histname.Data(), htitle.Data(), fNCentHistBins, fCentHistBins, fNPtHistBins, fPtHistBins, nMatchedTrackBins, matchedTrackBins);
441 
442  histname = "ClusterHistograms/hMatchedTrackNM02L04";
443  htitle = histname + ";Centrality (%);#it{E}_{clus} (GeV); N_{tracks}";
444  fHistManager.CreateTH3(histname.Data(), htitle.Data(), fNCentHistBins, fCentHistBins, fNPtHistBins, fPtHistBins, nMatchedTrackBins, matchedTrackBins);
445 
446  // Plot M02 distribution for clusters with matched tracks (centrality, Eclus nonlincorr, M02)
447  histname = "ClusterHistograms/hM02Matched";
448  htitle = histname + ";Centrality (%);#it{E}_{clus} (GeV); M02";
450 
451  // Plot M02 distribution for clusters without matched tracks (centrality, Eclus nonlincorr, M02)
452  histname = "ClusterHistograms/hM02Unmatched";
453  htitle = histname + ";Centrality (%);#it{E}_{clus} (GeV); M02";
455 
458 
459  // Plot E/p vs. M02 for 0-10% and 50-90% (Eclus nonlincorr, Eclus nonlincorr / trackPsum, M02)
460  histname = "ClusterHistograms/hEoverPM02Central";
461  htitle = histname + ";#it{E}_{clus} (GeV); #it{E}_{clus} / #Sigma#it{p}_{track} (GeV); M02";
463 
464  histname = "ClusterHistograms/hEoverPM02Peripheral";
465  htitle = histname + ";#it{E}_{clus} (GeV); #it{E}_{clus} / #Sigma#it{p}_{track} (GeV); M02";
467 
470 
471  // Plot Rcorr distribution (centrality, trackPSum, Rcorr = (Enonlincorr - Ehadcorr) / trackPSum)
472  histname = "ClusterHistograms/hRcorrVsCent";
473  htitle = histname + ";Centrality (%);#Sigma#it{p}_{track} (GeV); R_{corr} = #frac{#DeltaE_{clus}}{#Sigmap_{track}}";
474  fHistManager.CreateTH3(histname.Data(), htitle.Data(), fNCentHistBins, fCentHistBins, fNPtHistBins, fPtHistBins, nRcorrBins, RcorrBins);
475 
476  // Plot Rcorr distribution for 0-10% centrality (Eclus nonlincorr, trackPSum, Rcorr)
477  histname = "ClusterHistograms/hRcorr0-10";
478  htitle = histname + ";#it{E}_{clus} (GeV);#Sigma#it{p}_{track} (GeV); R_{corr} = #frac{#DeltaE_{clus}}{#Sigmap_{track}}";
479  fHistManager.CreateTH3(histname.Data(), htitle.Data(), fNPtHistBins, fPtHistBins, fNPtHistBins, fPtHistBins, nRcorrBins, RcorrBins);
480 
481  // Plot Rcorr distribution for 50-90% centrality (Eclus nonlincorr, trackPSum, Rcorr)
482  histname = "ClusterHistograms/hRcorr50-90";
483  htitle = histname + ";#it{E}_{clus} (GeV);#Sigma#it{p}_{track} (GeV); R_{corr} = #frac{#DeltaE_{clus}}{#Sigmap_{track}}";
484  fHistManager.CreateTH3(histname.Data(), htitle.Data(), fNPtHistBins, fPtHistBins, fNPtHistBins, fPtHistBins, nRcorrBins, RcorrBins);
485 
486  // Plot also Rcorr-clus (centrality, trackPSum, Rcorr-clus = (Enonlincorr - Ehadcorr) / Enonlincorr )
487  histname = "ClusterHistograms/hRcorrClusVsCent";
488  htitle = histname + ";Centrality (%);#Sigma#it{p}_{track} (GeV); #frac{#DeltaE_{clus}}{E_{clus}}";
489  fHistManager.CreateTH3(histname.Data(), htitle.Data(), fNCentHistBins, fCentHistBins, fNPtHistBins, fPtHistBins, nRcorrBins, RcorrBins);
490 
491  // Rcorr-clus for 0-10% centrality (Eclus nonlincorr, trackPSum, Rcorr-clus = (Enonlincorr - Ehadcorr) / Enonlincorr )
492  histname = "ClusterHistograms/hRcorrClus0-10";
493  htitle = histname + ";#it{E}_{clus} (GeV);#Sigma#it{p}_{track} (GeV); #frac{#DeltaE_{clus}}{E_{clus}}";
494  fHistManager.CreateTH3(histname.Data(), htitle.Data(), fNPtHistBins, fPtHistBins, fNPtHistBins, fPtHistBins, nRcorrBins, RcorrBins);
495 
496  // Rcorr-clus for 50-90% centrality (Eclus nonlincorr, trackPSum, Rcorr-clus = (Enonlincorr - Ehadcorr) / Enonlincorr )
497  histname = "ClusterHistograms/hRcorrClus50-90";
498  htitle = histname + ";#it{E}_{clus} (GeV);#Sigma#it{p}_{track} (GeV); #frac{#DeltaE_{clus}}{E_{clus}}";
499  fHistManager.CreateTH3(histname.Data(), htitle.Data(), fNPtHistBins, fPtHistBins, fNPtHistBins, fPtHistBins, nRcorrBins, RcorrBins);
500 
501 }
502 
503 /*
504  * This function allocates the histograms for the jet composition study.
505  */
507 {
508  TString histname;
509  TString htitle;
510  Int_t nPtBins = TMath::CeilNint(fMaxPt/2);
511 
512  const Int_t nRejBins = 32;
513  Double_t* rejReasonBins = new Double_t[nRejBins+1];
514  GenerateFixedBinArray(nRejBins, 0, nRejBins, rejReasonBins);
515  const Int_t nContributorTypes = 11;
516  Double_t *contributorTypeBins = GenerateFixedBinArray(nContributorTypes, -0.5, 10.5);
517  const Int_t nParticleTypes = 16;
518  Double_t *particleTypeBins = GenerateFixedBinArray(nParticleTypes, -0.5, 15.5);
519 
520  AliEmcalContainer* cont = 0;
521  TIter nextClusColl(&fClusterCollArray);
522  while ((cont = static_cast<AliEmcalContainer*>(nextClusColl()))) {
523 
524  histname = "ClusterHistogramsMC/hClusterRejectionReasonMC";
525  htitle = histname + ";Rejection reason;#it{E}_{clus} (GeV/)";
526  TH2* histMC2 = fHistManager.CreateTH2(histname.Data(), htitle.Data(), nRejBins, rejReasonBins, fNPtHistBins, fPtHistBins);
527  SetRejectionReasonLabels(histMC2->GetXaxis());
528  }
529 
530  // M02 vs. Energy vs. Particle type
531  histname = "ClusterHistogramsMC/hM02VsParticleTypeCentral";
532  htitle = histname + ";M02;#it{E}_{clus} (GeV); Particle type";
533  TH3* hM02VsParticleTypeCentral = fHistManager.CreateTH3(histname.Data(), htitle.Data(), fNM02HistBins, fM02HistBins, fNPtHistBins, fPtHistBins, nParticleTypes, particleTypeBins);
534  SetParticleTypeLabels(hM02VsParticleTypeCentral->GetZaxis());
535 
536  histname = "ClusterHistogramsMC/hM02VsParticleTypePeripheral";
537  htitle = histname + ";M02;#it{E}_{clus} (GeV); Particle type";
538  TH3* hM02VsParticleTypePeripheral = fHistManager.CreateTH3(histname.Data(), htitle.Data(), fNM02HistBins, fM02HistBins, fNPtHistBins, fPtHistBins, nParticleTypes, particleTypeBins);
539  SetParticleTypeLabels(hM02VsParticleTypePeripheral->GetZaxis());
540 
541  if (fPlotJetHistograms) {
542 
543  // M02 vs. Energy vs. Particle type vs. Jet pT, for particles inside jets
544  Int_t dim = 0;
545  TString title[20];
546  Int_t nbins[20] = {0};
547  Double_t min[30] = {0.};
548  Double_t max[30] = {0.};
549  Double_t *binEdges[20] = {0};
550 
551  title[dim] = "M02";
552  nbins[dim] = fNM02HistBins;
553  binEdges[dim] = fM02HistBins;
554  min[dim] = fM02HistBins[0];
555  max[dim] = fM02HistBins[fNM02HistBins];
556  dim++;
557 
558  title[dim] = "#it{E}_{clus} (GeV)";
559  nbins[dim] = fNPtHistBins;
560  binEdges[dim] = fPtHistBins;
561  min[dim] = fPtHistBins[0];
562  max[dim] = fPtHistBins[fNPtHistBins];
563  dim++;
564 
565  title[dim] = "Contributor type";
566  nbins[dim] = nContributorTypes;
567  min[dim] = -0.5;
568  max[dim] = 7.5;
569  binEdges[dim] = contributorTypeBins;
570  dim++;
571 
572  title[dim] = "#it{p}_{T,jet}^{corr}";
573  nbins[dim] = nPtBins;
574  min[dim] = 0;
575  max[dim] = fMaxPt;
576  binEdges[dim] = GenerateFixedBinArray(nbins[dim], min[dim], max[dim]);
577  dim++;
578 
579  TString thnname = "JetPerformanceMC/hM02VsContributorTypeJets";
580  THnSparse* hn = fHistManager.CreateTHnSparse(thnname.Data(), thnname.Data(), dim, nbins, min, max);
581  for (Int_t i = 0; i < dim; i++) {
582  hn->GetAxis(i)->SetTitle(title[i]);
583  hn->SetBinEdges(i, binEdges[i]);
584  }
585 
586  // Particle composition inside each jet -- jet pT vs. particle type vs. particle number vs. particle pT sum
587  // (One entry per jet for each particle type)
588  dim = 0;
589 
590  title[dim] = "#it{p}_{T,jet}^{corr}";
591  nbins[dim] = nPtBins;
592  min[dim] = 0;
593  max[dim] = fMaxPt;
594  binEdges[dim] = GenerateFixedBinArray(nbins[dim], min[dim], max[dim]);
595  dim++;
596 
597  title[dim] = "Contributor type";
598  nbins[dim] = nContributorTypes;
599  min[dim] = -0.5;
600  max[dim] = 7.5;
601  binEdges[dim] = contributorTypeBins;
602  dim++;
603 
604  title[dim] = "N";
605  nbins[dim] = 30;
606  min[dim] = -0.5;
607  max[dim] = 29.5;
608  binEdges[dim] = GenerateFixedBinArray(nbins[dim], min[dim], max[dim]);
609  dim++;
610 
611  title[dim] = "#it{p}_{T,sum} (GeV)";
612  nbins[dim] = fNPtHistBins;
613  binEdges[dim] = fPtHistBins;
614  min[dim] = fPtHistBins[0];
615  max[dim] = fPtHistBins[fNPtHistBins];
616  dim++;
617 
618  thnname = "JetPerformanceMC/hJetComposition";
619  THnSparse* thn = fHistManager.CreateTHnSparse(thnname.Data(), thnname.Data(), dim, nbins, min, max);
620  for (Int_t i = 0; i < dim; i++) {
621  thn->GetAxis(i)->SetTitle(title[i]);
622  thn->SetBinEdges(i, binEdges[i]);
623  }
624 
625  // Hadronic calo energy in each jet
626 
627  // Jet pT vs. Summed energy of hadronic clusters without a matched track
628  histname = "JetPerformance/hHadCaloEnergyUnmatched";
629  htitle = histname + ";#it{p}_{T,jet} (GeV);#it{p}_{T,had} (GeV)";
630  fHistManager.CreateTH2(histname.Data(), htitle.Data(), fNPtHistBins, fPtHistBins, fNPtHistBins, fPtHistBins);
631 
632  // Jet pT vs. Summed energy of hadronic clusters with a matched track (before hadronic correction)
633  histname = "JetPerformance/hHadCaloEnergyMatchedNonlincorr";
634  htitle = histname + ";#it{p}_{T,jet} (GeV);#it{p}_{T,had} (GeV)";
635  fHistManager.CreateTH2(histname.Data(), htitle.Data(), fNPtHistBins, fPtHistBins, fNPtHistBins, fPtHistBins);
636 
637  // Jet pT vs. Summed energy of hadronic clusters with a matched track (after hadronic correction)
638  histname = "JetPerformance/hHadCaloEnergyMatchedHadCorr";
639  htitle = histname + ";#it{p}_{T,jet} (GeV);#it{p}_{T,had} (GeV)";
640  fHistManager.CreateTH2(histname.Data(), htitle.Data(), fNPtHistBins, fPtHistBins, fNPtHistBins, fPtHistBins);
641 
642  }
643 
644 }
645 
646 /*
647  * This function sets axis labels for particle type histograms.
648  */
650 {
651  axis->SetBinLabel(1, "SinglePhoton");
652  axis->SetBinLabel(2, "SingleElectron");
653  axis->SetBinLabel(3, "SingleChargedPion");
654  axis->SetBinLabel(4, "SingleProton");
655  axis->SetBinLabel(5, "SingleAntiProton");
656  axis->SetBinLabel(6, "SingleChargedKaon");
657  axis->SetBinLabel(7, "SingleK0L");
658  axis->SetBinLabel(8, "SingleNeutron");
659  axis->SetBinLabel(9, "SingleAntiNeutron");
660  axis->SetBinLabel(10, "SingleOther");
661  axis->SetBinLabel(11, "PhotonHadron");
662  axis->SetBinLabel(12, "MergedPi0");
663  axis->SetBinLabel(13, "PhotonPhotonOther");
664  axis->SetBinLabel(14, "HadronHadron");
665  axis->SetBinLabel(15, "TwoContributorsOther");
666  axis->SetBinLabel(16, "MoreThanTwoContributors");
667 }
668 
669 /*
670  * This function allocates background subtraction histograms, if enabled.
671  * A set of histograms is allocated per each jet container.
672  */
674 {
675  TString histname;
676  TString title;
677 
678  AliJetContainer* jets = 0;
679  TIter nextJetColl(&fJetCollArray);
680  while ((jets = static_cast<AliJetContainer*>(nextJetColl()))) {
681 
682  histname = TString::Format("%s/BackgroundHistograms/hScaleFactorEMCal", jets->GetArrayName().Data());
683  title = histname + ";Centrality;Scale factor;counts";
684  fHistManager.CreateTH2(histname.Data(), title.Data(), 50, 0, 100, 100, 0, 5);
685 
686  histname = TString::Format("%s/BackgroundHistograms/hDeltaPtEMCal", jets->GetArrayName().Data());
687  title = histname + ";Centrality (%);#delta#it{p}_{T} (GeV/#it{c});counts";
688  fHistManager.CreateTH2(histname.Data(), title.Data(), 10, 0, 100, 400, -50, 150);
689 
690  }
691 }
692 
693 /*
694  * This function allocates the histograms for single jets, when the "simulated" trigger has been fired.
695  * A set of histograms is allocated per each jet container.
696  */
698 {
699  TString histname;
700  TString title;
701  Int_t nPtBins = TMath::CeilNint(fMaxPt/2);
702 
703  //----------------------------------------------
704  // Trigger patch histograms
705 
706  // patch eta vs. phi
707  histname = "TriggerSimHistograms/hEtaVsPhi";
708  title = histname + ";#eta_{patch} (rad);#phi_{patch} (rad)";
709  fHistManager.CreateTH2(histname.Data(), title.Data(), 140, -0.7, 0.7, 500, 1., 6.);
710 
711  // N patches
712  histname = "TriggerSimHistograms/hNPatches";
713  title = histname + ";#it{N}_{patches};type";
714  fHistManager.CreateTH2(histname.Data(), title.Data(), 200, 0, 200, 2, -0.5, 1.5);
715 
716  // patch E vs. centrality
717  histname = "TriggerSimHistograms/hPatchE";
718  title = histname + ";Centrality (%);#it{E}_{patch} (GeV)";
719  fHistManager.CreateTH2(histname.Data(), title.Data(), 50, 0, 100, nPtBins, 0, fMaxPt);
720 
721  // patch median vs. Centrality
722  histname = "TriggerSimHistograms/hPatchMedianE";
723  title = histname + ";Centrality (%);#it{E}_{patch,med} (GeV);type";
724  fHistManager.CreateTH3(histname.Data(), title.Data(), 50, 0, 100, 100, 0, 50, 2, -0.5, 1.5);
725 
726  //----------------------------------------------
727  // Jet histograms for "triggered" events
728  AliJetContainer* jets = 0;
729  TIter nextJetColl(&fJetCollArray);
730  while ((jets = static_cast<AliJetContainer*>(nextJetColl()))) {
731 
732  // Jet rejection reason
733  histname = TString::Format("%s/TriggerSimHistograms/hJetRejectionReason", jets->GetArrayName().Data());
734  title = histname + ";Rejection reason;#it{p}_{T,jet} (GeV/#it{c});counts";
735  TH2* hist = fHistManager.CreateTH2(histname.Data(), title.Data(), 32, 0, 32, 50, 0, fMaxPt);
736  SetRejectionReasonLabels(hist->GetXaxis());
737 
738  // Rho vs. Centrality
739  if (!jets->GetRhoName().IsNull()) {
740  histname = TString::Format("%s/TriggerSimHistograms/hRhoVsCent", jets->GetArrayName().Data());
741  title = histname + ";Centrality (%);#rho (GeV/#it{c});counts";
742  fHistManager.CreateTH2(histname.Data(), title.Data(), 50, 0, 100, 100, 0, 500);
743  }
744 
745  // (Centrality, pT, NEF)
746  Int_t nbinsx = 20; Int_t minx = 0; Int_t maxx = 100;
747  Int_t nbinsy = nPtBins; Int_t miny = 0; Int_t maxy = fMaxPt;
748  Int_t nbinsz = 50; Int_t minz = 0; Int_t maxz = 1.;
749 
750  histname = TString::Format("%s/TriggerSimHistograms/hNEFVsPtEMCal", jets->GetArrayName().Data());
751  title = histname + ";Centrality (%);#it{p}_{T}^{corr} (GeV/#it{c});NEF";
752  fHistManager.CreateTH3(histname.Data(), title.Data(), nbinsx, minx, maxx, nbinsy, miny, maxy, nbinsz, minz, maxz);
753 
754  histname = TString::Format("%s/TriggerSimHistograms/hNEFVsPtDCal", jets->GetArrayName().Data());
755  title = histname + ";Centrality (%);#it{p}_{T}^{corr} (GeV/#it{c});NEF";
756  fHistManager.CreateTH3(histname.Data(), title.Data(), nbinsx, minx, maxx, nbinsy, miny, maxy, nbinsz, minz, maxz);
757 
758  // pT-leading vs. pT
759  histname = TString::Format("%s/TriggerSimHistograms/hPtLeadingVsPt", jets->GetArrayName().Data());
760  title = histname + ";#it{p}_{T}^{corr} (GeV/#it{c});#it{p}_{T,particle}^{leading} (GeV/#it{c})";
761  fHistManager.CreateTH2(histname.Data(), title.Data(), nPtBins, 0, fMaxPt, nPtBins, 0, fMaxPt);
762 
763  // A vs. pT
764  histname = TString::Format("%s/TriggerSimHistograms/hAreaVsPt", jets->GetArrayName().Data());
765  title = histname + ";#it{p}_{T}^{corr} (GeV/#it{c});#it{A}_{jet}";
766  fHistManager.CreateTH2(histname.Data(), title.Data(), nPtBins, 0, fMaxPt, fMaxPt/3, 0, 0.5);
767 
768  // (Centrality, pT, z-leading (charged))
769  nbinsx = 20; minx = 0; maxx = 100;
770  nbinsy = nPtBins; miny = 0; maxy = fMaxPt;
771  nbinsz = 50; minz = 0; maxz = 1.;
772 
773  histname = TString::Format("%s/TriggerSimHistograms/hZLeadingVsPtEMCal", jets->GetArrayName().Data());
774  title = histname + ";Centrality (%);#it{p}_{T}^{corr} (GeV/#it{c});#it{z}_{leading}";
775  fHistManager.CreateTH3(histname.Data(), title.Data(), nbinsx, minx, maxx, nbinsy, miny, maxy, nbinsz, minz, maxz);
776 
777  histname = TString::Format("%s/TriggerSimHistograms/hZLeadingVsPtDCal", jets->GetArrayName().Data());
778  title = histname + ";Centrality (%);#it{p}_{T}^{corr} (GeV/#it{c});#it{z}_{leading}";
779  fHistManager.CreateTH3(histname.Data(), title.Data(), nbinsx, minx, maxx, nbinsy, miny, maxy, nbinsz, minz, maxz);
780 
781  // z (charged) vs. pT
782  nbinsx = 20; minx = 0; maxx = 100;
783  nbinsy = nPtBins; miny = 0; maxy = fMaxPt;
784  nbinsz = 50; minz = 0; maxz = 1.;
785 
786  histname = TString::Format("%s/TriggerSimHistograms/hZVsPtEMCal", jets->GetArrayName().Data());
787  title = histname + ";Centrality (%);#it{p}_{T}^{corr} (GeV/#it{c});#it{z}";
788  fHistManager.CreateTH3(histname.Data(), title.Data(), nbinsx, minx, maxx, nbinsy, miny, maxy, nbinsz, minz, maxz);
789 
790  histname = TString::Format("%s/TriggerSimHistograms/hZVsPtDCal", jets->GetArrayName().Data());
791  title = histname + ";Centrality (%);#it{p}_{T}^{corr} (GeV/#it{c});#it{z}";
792  fHistManager.CreateTH3(histname.Data(), title.Data(), nbinsx, minx, maxx, nbinsy, miny, maxy, nbinsz, minz, maxz);
793 
794  // (Centrality, pT, Nconst)
795  nbinsx = 20; minx = 0; maxx = 100;
796  nbinsy = nPtBins; miny = 0; maxy = fMaxPt;
797  nbinsz = 50; minz = 0; maxz = fMaxPt;
798 
799  histname = TString::Format("%s/TriggerSimHistograms/hNConstVsPtEMCal", jets->GetArrayName().Data());
800  title = histname + ";Centrality (%);#it{p}_{T}^{corr} (GeV/#it{c});No. of constituents";
801  fHistManager.CreateTH3(histname.Data(), title.Data(), nbinsx, minx, maxx, nbinsy, miny, maxy, nbinsz, minz, maxz);
802 
803  histname = TString::Format("%s/TriggerSimHistograms/hNConstVsPtDCal", jets->GetArrayName().Data());
804  title = histname + ";Centrality (%);#it{p}_{T}^{corr} (GeV/#it{c});No. of constituents";
805  fHistManager.CreateTH3(histname.Data(), title.Data(), nbinsx, minx, maxx, nbinsy, miny, maxy, nbinsz, minz, maxz);
806 
807  }
808 }
809 
810 /*
811  * This function allocates histograms for matched truth-det jets in the case of embedding.
812  * The jet matching information must be previously filled by another task, such as AliJetResponseMaker.
813  */
815 {
816  TString histname;
817  TString title;
818  Int_t nPtBins = TMath::CeilNint(fMaxPt/2);
819 
820  // Response matrix, (centrality, pT-truth, pT-det)
821  Int_t nbinsx = 20; Int_t minx = 0; Int_t maxx = 100;
822  Int_t nbinsy = nPtBins; Int_t miny = 0; Int_t maxy = fMaxPt;
823  Int_t nbinsz = nPtBins; Int_t minz = 0; Int_t maxz = fMaxPt;
824 
825  histname = "MatchedJetHistograms/hResponseMatrixEMCal";
826  title = histname + ";Centrality (%);#it{p}_{T}^{truth} (GeV/#it{c});#it{p}_{T,corr}^{det} (GeV/#it{c})";
827  fHistManager.CreateTH3(histname.Data(), title.Data(), nbinsx, minx, maxx, nbinsy, miny, maxy, nbinsz, minz, maxz);
828 
829  histname = "MatchedJetHistograms/hResponseMatrixDCal";
830  title = histname + ";Centrality (%);#it{p}_{T}^{truth} (GeV/#it{c});#it{p}_{T,corr}^{det} (GeV/#it{c})";
831  fHistManager.CreateTH3(histname.Data(), title.Data(), nbinsx, minx, maxx, nbinsy, miny, maxy, nbinsz, minz, maxz);
832 
833  // JES shift, (centrality, pT-truth, (pT-det - pT-truth) / pT-truth)
834  nbinsx = 20; minx = 0; maxx = 100;
835  nbinsy = nPtBins; miny = 0; maxy = fMaxPt;
836  nbinsz = 250; minz = -5.; maxz = 5.;
837 
838  histname = "MatchedJetHistograms/hJESshiftEMCal";
839  title = histname + ";Centrality (%);#it{p}_{T}^{truth} (GeV/#it{c});#frac{#it{p}_{T,corr}^{det} - #it{p}_{T}^{truth}}{#it{p}_{T}^{truth}}";
840  fHistManager.CreateTH3(histname.Data(), title.Data(), nbinsx, minx, maxx, nbinsy, miny, maxy, nbinsz, minz, maxz);
841 
842  histname = "MatchedJetHistograms/hJESshiftDCal";
843  title = histname + ";Centrality (%);#it{p}_{T}^{truth} (GeV/#it{c});#frac{#it{p}_{T,corr}^{det} - #it{p}_{T}^{truth}}{#it{p}_{T}^{truth}}";
844  fHistManager.CreateTH3(histname.Data(), title.Data(), nbinsx, minx, maxx, nbinsy, miny, maxy, nbinsz, minz, maxz);
845 
846  // NEF of det-level matched jets, (centrality, pT-truth, NEF)
847  nbinsx = 20; minx = 0; maxx = 100;
848  nbinsy = nPtBins; miny = 0; maxy = fMaxPt;
849  nbinsz = 50; minz = 0; maxz = 1.;
850 
851  histname = "MatchedJetHistograms/hNEFVsPt";
852  title = histname + ";Centrality (%);#it{p}_{T,corr}^{det} (GeV/#it{c});Calo energy fraction";
853  fHistManager.CreateTH3(histname.Data(), title.Data(), nbinsx, minx, maxx, nbinsy, miny, maxy, nbinsz, minz, maxz);
854 
855  // z-leading (charged) of det-level matched jets, (centrality, pT-truth, z-leading)
856  nbinsx = 20; minx = 0; maxx = 100;
857  nbinsy = nPtBins; miny = 0; maxy = fMaxPt;
858  nbinsz = 50; minz = 0; maxz = 1.;
859 
860  histname = "MatchedJetHistograms/hZLeadingVsPt";
861  title = histname + ";Centrality (%);#it{p}_{T,corr}^{det} (GeV/#it{c});#it{z}_{leading}";
862  fHistManager.CreateTH3(histname.Data(), title.Data(), nbinsx, minx, maxx, nbinsy, miny, maxy, nbinsz, minz, maxz);
863 
864  // Matching distance, (centrality, pT-truth, R)
865  nbinsx = 20; minx = 0; maxx = 100;
866  nbinsy = nPtBins; miny = 0; maxy = fMaxPt;
867  nbinsz = 50; minz = 0; maxz = 1.;
868 
869  histname = "MatchedJetHistograms/hMatchingDistance";
870  title = histname + ";Centrality (%);#it{p}_{T}^{truth} (GeV/#it{c});R";
871  fHistManager.CreateTH3(histname.Data(), title.Data(), nbinsx, minx, maxx, nbinsy, miny, maxy, nbinsz, minz, maxz);
872 
873 }
874 
880 {
881  // Configure base class to set fTriggerPatchInfo to array of trigger patches, each event
882  // (Need to call this before base class ExecOnce)
883  if (fDoTriggerSimulation) {
884  this->SetCaloTriggerPatchInfoName("EmcalTriggers");
885  }
886 
888 
889  fNeedEmcalGeom = kTRUE;
890 
891  // Check if trigger patches are loaded
892  if (fDoTriggerSimulation) {
893  if (fTriggerPatchInfo) {
894  TString objname(fTriggerPatchInfo->GetClass()->GetName());
895  TClass cls(objname);
896  if (!cls.InheritsFrom("AliEMCALTriggerPatchInfo")) {
897  AliError(Form("%s: Objects of type %s in %s are not inherited from AliEMCALTriggerPatchInfo!",
898  GetName(), cls.GetName(), "EmcalTriggers"));
899  fTriggerPatchInfo = 0;
900  }
901  }
902  if (!fTriggerPatchInfo) {
903  AliError(Form("%s: Unable to get trigger patch container with name %s. Aborting", GetName(), "EmcalTriggers"));
904  return;
905  }
906  }
907 }
908 
913 
914  // Get the downscaling factors for MB triggers (to be used to calculate trigger efficiency)
915 
917 
918  // Get instance of the downscale factor helper class
920  downscaleOCDB->SetRun(InputEvent()->GetRunNumber());
921 
922  // There are two possible min bias triggers for LHC15o
923  TString triggerNameMB1 = "CINT7-B-NOPF-CENT";
924  TString triggerNameMB2 = "CV0L7-B-NOPF-CENT";
925  TString triggerNameJE = "CINT7EJ1-B-NOPF-CENTNOPMD";
926 
927  // Get the downscale factor for whichever MB trigger exists in the given run
928  std::vector<TString> runtriggers = downscaleOCDB->GetTriggerClasses();
929  Double_t downscalefactor;
930  for (auto i : runtriggers) {
931  if (i.EqualTo(triggerNameMB1) || i.EqualTo(triggerNameMB2)) {
932  downscalefactor = downscaleOCDB->GetDownscaleFactorForTriggerClass(i.Data());
933  break;
934  }
935  }
936 
937  // Store the inverse of the downscale factor, used later to weight the pT spectrum
938  fMBUpscaleFactor = 1/downscalefactor;
939 
940  TString histname = "Trigger/hMBDownscaleFactor";
941  fHistManager.FillTH1(histname.Data(), fMBUpscaleFactor);
942 
943  }
944 
945 }
946 
951 {
952  if (fUseAliEventCuts) {
953  if (!fEventCuts.AcceptEvent(InputEvent()))
954  {
955  PostData(1, fOutput);
956  return kFALSE;
957  }
958  }
959  else {
961  }
962  return kTRUE;
963 }
964 
973 {
974  TString histname;
975  AliJetContainer* jetCont = 0;
976  TIter next(&fJetCollArray);
977  while ((jetCont = static_cast<AliJetContainer*>(next()))) {
978  TString jetContName = jetCont->GetName();
979 
980  // Compute the full jet background scale factor and delta-pt
981  if (fComputeBackground) {
983  }
984 
985  // Do a simple trigger simulation (if requested)
986  if (fDoTriggerSimulation) {
988  }
989 
990  }
991 
992  // Only fill the embedding qa plots if:
993  // - We are using the embedding helper
994  // - The class has been initialized
995  // - Both jet collections are available
996  if (fEmbeddingQA.IsInitialized()) {
998  }
999 
1000  return kTRUE;
1001 }
1002 
1007 {
1008  TString histname;
1009 
1010  // Check if trigger patches are loaded
1011  if (fTriggerPatchInfo) {
1012  TString objname(fTriggerPatchInfo->GetClass()->GetName());
1013  TClass cls(objname);
1014  if (!cls.InheritsFrom("AliEMCALTriggerPatchInfo")) {
1015  AliError(Form("%s: Objects of type %s in %s are not inherited from AliEMCALTriggerPatchInfo!",
1016  GetName(), cls.GetName(), "EmcalTriggers"));
1017  fTriggerPatchInfo = 0;
1018  }
1019  }
1020  if (!fTriggerPatchInfo) {
1021  AliError(Form("%s: Unable to get trigger patch container with name %s. Aborting", GetName(), "EmcalTriggers"));
1022  return;
1023  }
1024 
1025  // Compute patches in EMCal, DCal (I want offline simple trigger patch, i.e. patch calculated using FEE energy)
1026  std::vector<Double_t> vecEMCal;
1027  std::vector<Double_t> vecDCal;
1028  for(auto p : *fTriggerPatchInfo){
1029  AliEMCALTriggerPatchInfo *recpatch = static_cast<AliEMCALTriggerPatchInfo *>(p);
1030  if (recpatch) {
1031 
1032  if(!recpatch->IsJetHighSimple()) continue;
1033 
1034  histname = "TriggerSimHistograms/hEtaVsPhi";
1035  fHistManager.FillTH2(histname.Data(), recpatch->GetEtaGeo(), recpatch->GetPhiGeo());
1036 
1037  histname = "TriggerSimHistograms/hPatchE";
1038  fHistManager.FillTH2(histname.Data(), fCent, recpatch->GetPatchE());
1039 
1040  if (recpatch->IsEMCal()) {
1041  vecEMCal.push_back(recpatch->GetPatchE());
1042  } else {
1043  vecDCal.push_back(recpatch->GetPatchE());
1044  }
1045 
1046  }
1047  }
1048 
1049  // Compute the median in each calorimeter
1050  const Int_t nBkgPatchesEMCal = vecEMCal.size(); // 6*8;
1051  const Int_t nBkgPatchesDCal = vecDCal.size(); // 4*5;
1052  fMedianEMCal = TMath::Median(nBkgPatchesEMCal, &vecEMCal[0]); // point to array used internally by vector
1053  fMedianDCal = TMath::Median(nBkgPatchesDCal, &vecDCal[0]);
1054 
1055  histname = "TriggerSimHistograms/hPatchMedianE";
1056  fHistManager.FillTH3(histname.Data(), fCent, fMedianEMCal, kEMCal);
1057  fHistManager.FillTH3(histname.Data(), fCent, fMedianDCal, kDCal);
1058 
1059  histname = "TriggerSimHistograms/hNPatches";
1060  fHistManager.FillTH2(histname.Data(), nBkgPatchesEMCal, kEMCal);
1061  fHistManager.FillTH2(histname.Data(), nBkgPatchesDCal, kDCal);
1062 
1063  // Then compute background subtracted patches, by subtracting from each patch the median patch E from the opposite hemisphere
1064  // If a patch is above threshold, the event is "triggered"
1065  Bool_t fkEMCEJE = kFALSE;
1066  Double_t threshold = 20;
1067  for(auto p : *fTriggerPatchInfo){
1068  AliEMCALTriggerPatchInfo *recpatch = static_cast<AliEMCALTriggerPatchInfo *>(p);
1069  if (recpatch) {
1070 
1071  if(!recpatch->IsJetHighSimple()) continue;
1072 
1073  if (recpatch->IsEMCal()) {
1074  if ((recpatch->GetPatchE() - fMedianDCal) > threshold) {
1075  fkEMCEJE = kTRUE;
1076  break;
1077  }
1078  } else {
1079  if ((recpatch->GetPatchE() - fMedianEMCal) > threshold) {
1080  fkEMCEJE = kTRUE;
1081  break;
1082  }
1083  }
1084  }
1085  }
1086 
1087  if (fkEMCEJE) {
1089  }
1090 
1091 }
1092 
1100 {
1101 
1102  if (fPlotJetHistograms) {
1104  }
1105  if (fPlotClusterHistograms) {
1107  }
1110  }
1113  }
1114 
1115  return kTRUE;
1116 }
1117 
1123 {
1124  TString histname;
1125  AliJetContainer* jets = 0;
1126  TIter nextJetColl(&fJetCollArray);
1127  while ((jets = static_cast<AliJetContainer*>(nextJetColl()))) {
1128  TString jetContName = jets->GetName();
1129 
1130  Double_t rhoVal = 0;
1131  if (jets->GetRhoParameter()) {
1132  rhoVal = jets->GetRhoVal();
1133  histname = TString::Format("%s/JetHistograms/hRhoVsCent", jets->GetArrayName().Data());
1134  fHistManager.FillTH2(histname.Data(), fCent, rhoVal);
1135  }
1136 
1137  for (auto jet : jets->all()) {
1138 
1139  Float_t ptLeading = jets->GetLeadingHadronPt(jet);
1140  Float_t corrPt = GetJetPt(jet, rhoVal);
1141 
1142  // A vs. pT (fill before area cut)
1143  histname = TString::Format("%s/JetHistograms/hAreaVsPt", jets->GetArrayName().Data());
1144  fHistManager.FillTH2(histname.Data(), corrPt, jet->Area());
1145 
1146 
1147  // Rejection reason
1148  UInt_t rejectionReason = 0;
1149  if (!jets->AcceptJet(jet, rejectionReason)) {
1150  histname = TString::Format("%s/JetHistograms/hJetRejectionReason", jets->GetArrayName().Data());
1151  fHistManager.FillTH2(histname.Data(), jets->GetRejectionReasonBitPosition(rejectionReason), jet->Pt());
1152  continue;
1153  }
1154 
1155  // compute jet acceptance type
1156  Double_t type = GetJetType(jet);
1157  if ( (type != kEMCal) && (type != kDCal) ) {
1158  continue;
1159  }
1160 
1161  // (Centrality, pT, NEF)
1162  if (type == kEMCal) {
1163  histname = TString::Format("%s/JetHistograms/hNEFVsPtEMCal", jets->GetArrayName().Data());
1164  }
1165  else if (type == kDCal) {
1166  histname = TString::Format("%s/JetHistograms/hNEFVsPtDCal", jets->GetArrayName().Data());
1167  }
1168  fHistManager.FillTH3(histname, fCent, corrPt, jet->NEF());
1169 
1170  // (Centrality, pT upscaled, calo type)
1171  if (fComputeMBDownscaling) {
1172  histname = TString::Format("%s/JetHistograms/hPtUpscaledMB", jets->GetArrayName().Data());
1173  fHistManager.FillTH3(histname.Data(), fCent, corrPt, type, fMBUpscaleFactor);
1174  }
1175 
1176  // pT-leading vs. pT
1177  histname = TString::Format("%s/JetHistograms/hPtLeadingVsPt", jets->GetArrayName().Data());
1178  fHistManager.FillTH2(histname.Data(), corrPt, ptLeading);
1179 
1180  // (Centrality, pT, z-leading (charged))
1181  if (type == kEMCal) {
1182  histname = TString::Format("%s/JetHistograms/hZLeadingVsPtEMCal", jets->GetArrayName().Data());
1183  }
1184  else if (type == kDCal) {
1185  histname = TString::Format("%s/JetHistograms/hZLeadingVsPtDCal", jets->GetArrayName().Data());
1186  }
1187  TLorentzVector leadPart;
1188  jets->GetLeadingHadronMomentum(leadPart, jet);
1189  Double_t z = GetParallelFraction(leadPart.Vect(), jet);
1190  if (z == 1 || (z > 1 && z - 1 < 1e-3)) z = 0.999; // so that it will contribute to the bin <1
1191  fHistManager.FillTH3(histname, fCent, corrPt, z);
1192 
1193  // (Centrality, pT, z (charged))
1194  if (type == kEMCal) {
1195  histname = TString::Format("%s/JetHistograms/hZVsPtEMCal", jets->GetArrayName().Data());
1196  }
1197  else if (type == kDCal) {
1198  histname = TString::Format("%s/JetHistograms/hZVsPtDCal", jets->GetArrayName().Data());
1199  }
1200  const AliVTrack* track;
1201  for (Int_t i=0; i<jet->GetNumberOfTracks(); i++) {
1202  track = static_cast<AliVTrack*>(jet->Track(i));
1203  z = track->Pt() / TMath::Abs(corrPt);
1204  fHistManager.FillTH3(histname, fCent, corrPt, z);
1205  }
1206 
1207  // (Centrality, pT, Nconst)
1208  if (type == kEMCal) {
1209  histname = TString::Format("%s/JetHistograms/hNConstVsPtEMCal", jets->GetArrayName().Data());
1210  }
1211  else if (type == kDCal) {
1212  histname = TString::Format("%s/JetHistograms/hNConstVsPtDCal", jets->GetArrayName().Data());
1213  }
1214  fHistManager.FillTH3(histname, fCent, corrPt, 1.*jet->GetNumberOfConstituents());
1215 
1216  // (Centrality, jet pT, Enonlincorr - Ehadcorr)
1217  Double_t deltaEhadcorr = 0;
1218  const AliVCluster* clus = nullptr;
1219  Int_t nClusters = jet->GetNumberOfClusters();
1220  for (Int_t iClus = 0; iClus < nClusters; iClus++) {
1221  clus = jet->Cluster(iClus);
1222  deltaEhadcorr += (clus->GetNonLinCorrEnergy() - clus->GetHadCorrEnergy());
1223  }
1224 
1225  histname = TString::Format("%s/JetHistograms/hDeltaEHadCorr", jets->GetArrayName().Data());
1226  fHistManager.FillTH3(histname, fCent, corrPt, deltaEhadcorr);
1227 
1228 
1229  // (Median patch energy, calo type, jet pT, centrality)
1230  if (fDoTriggerSimulation) {
1231  histname = TString::Format("%s/JetHistograms/hMedPatchJet", jets->GetArrayName().Data());
1232  Double_t x[4] = {fMedianEMCal, kEMCal, corrPt, fCent};
1233  fHistManager.FillTHnSparse(histname, x);
1234  Double_t y[4] = {fMedianDCal, kDCal, corrPt, fCent};
1235  fHistManager.FillTHnSparse(histname, y);
1236  }
1237 
1238  } //jet loop
1239 
1240  }
1241 }
1242 
1243 /*
1244  * This function fills the histograms for the calorimeter performance study.
1245  */
1247 {
1248  TString histname;
1249 
1250  // Loop through clusters
1252  const AliVCluster* clus;
1253  for (auto it : clusters->accepted_momentum()) {
1254 
1255  clus = it.second;
1256  Double_t clusPhi = it.first.Phi_0_2pi();
1257  Double_t clusEta = it.first.Eta();
1258 
1259  // Include only EMCal/DCal clusters
1260  if (!clus->IsEMCAL()) {
1261  continue;
1262  }
1263 
1264  // Compute the sum of matched track momentum, and various track matching / hadronic corretion quantities
1265  Double_t trackPSum = 0;
1266  Int_t nTracksMatched = 0;
1267  const AliVTrack* track = nullptr;
1268  for (Int_t itrack=0; itrack < clus->GetNTracksMatched(); itrack++) {
1269  track = dynamic_cast<AliVTrack*>(clus->GetTrackMatched(itrack));
1270  if (!track) {
1271  continue;
1272  }
1273 
1274  Double_t trackPhi = TVector2::Phi_0_2pi(track->GetTrackPhiOnEMCal());
1275  Double_t trackEta = track->GetTrackEtaOnEMCal();
1276  Double_t deta = TMath::Abs(clusPhi - trackPhi);
1277  Double_t dphi = TMath::Abs(clusEta - trackEta);
1278 
1280  trackPSum += track->P();
1281  nTracksMatched++;
1282  }
1283  }
1284 
1285  Double_t EoverP = 0;
1286  if (trackPSum > 1e-3) {
1287  EoverP = clus->GetNonLinCorrEnergy() / trackPSum;
1288  }
1289 
1290  Double_t deltaE = clus->GetNonLinCorrEnergy() - clus->GetHadCorrEnergy();
1291  Double_t Rcorr = 0;
1292  if (trackPSum > 1e-3) {
1293  Rcorr = deltaE / trackPSum;
1294  }
1295  Double_t RcorrClus = deltaE / clus->GetNonLinCorrEnergy();
1296 
1299 
1300  // Fill M02 distribution (centrality, Eclus nonlincorr, M02)
1301  histname = "ClusterHistograms/hM02";
1302  fHistManager.FillTH3(histname, fCent, clus->GetNonLinCorrEnergy(), clus->GetM02());
1303 
1304  // Plot Ncell distribution for M02 > 0.4 or M02 < 0.4 (centrality, Eclus nonlincorr, Ncells)
1305  if (clus->GetM02() > 0.4) {
1306  histname = "ClusterHistograms/hNcellsM02G04";
1307  fHistManager.FillTH3(histname, fCent, clus->GetNonLinCorrEnergy(), clus->GetNCells());
1308  }
1309  if (clus->GetM02() > 0.1 && clus->GetM02() < 0.4) {
1310  histname = "ClusterHistograms/hNcellsM02L04";
1311  fHistManager.FillTH3(histname, fCent, clus->GetNonLinCorrEnergy(), clus->GetNCells());
1312  }
1313 
1316 
1317  // Plot matched track pT (centrality, Eclus nonlincorr, trackPsum)
1318  histname = "ClusterHistograms/hMatchedTrackPt";
1319  fHistManager.FillTH3(histname, fCent, clus->GetNonLinCorrEnergy(), trackPSum);
1320 
1321  if (clus->GetM02() > 0.4) {
1322  histname = "ClusterHistograms/hMatchedTrackPtM02G04";
1323  fHistManager.FillTH3(histname, fCent, clus->GetNonLinCorrEnergy(), trackPSum);
1324  }
1325  if (clus->GetM02() > 0.1 && clus->GetM02() < 0.4) {
1326  histname = "ClusterHistograms/hMatchedTrackPtM02L04";
1327  fHistManager.FillTH3(histname, fCent, clus->GetNonLinCorrEnergy(), trackPSum);
1328  }
1329 
1330  // Plot number of matched tracks (centrality, Eclus nonlincorr, N matches)
1331  histname = "ClusterHistograms/hMatchedTrackN";
1332  fHistManager.FillTH3(histname, fCent, clus->GetNonLinCorrEnergy(), nTracksMatched);
1333 
1334  if (clus->GetM02() > 0.4) {
1335  histname = "ClusterHistograms/hMatchedTrackNM02G04";
1336  fHistManager.FillTH3(histname, fCent, clus->GetNonLinCorrEnergy(), nTracksMatched);
1337  }
1338  if (clus->GetM02() > 0.1 && clus->GetM02() < 0.4) {
1339  histname = "ClusterHistograms/hMatchedTrackNM02L04";
1340  fHistManager.FillTH3(histname, fCent, clus->GetNonLinCorrEnergy(), nTracksMatched);
1341  }
1342 
1343  // Plot M02 distribution for clusters with matched tracks (centrality, Eclus nonlincorr, M02)
1344  histname = "ClusterHistograms/hM02Matched";
1345  if (nTracksMatched > 0) {
1346  fHistManager.FillTH3(histname, fCent, clus->GetNonLinCorrEnergy(), clus->GetM02());
1347  }
1348 
1349  // Plot M02 distribution for clusters without matched tracks (centrality, Eclus nonlincorr, M02)
1350  histname = "ClusterHistograms/hM02Unmatched";
1351  if (nTracksMatched == 0) {
1352  fHistManager.FillTH3(histname, fCent, clus->GetNonLinCorrEnergy(), clus->GetM02());
1353  }
1354 
1357 
1358  // Plot E/p vs. M02 for 0-10% and 50-90% (Eclus nonlincorr, Eclus nonlincorr / trackPsum, M02)
1359  if (fCent < 10) {
1360  histname = "ClusterHistograms/hEoverPM02Central";
1361  fHistManager.FillTH3(histname, clus->GetNonLinCorrEnergy(), EoverP, clus->GetM02());
1362  }
1363  if (fCent > 50 && fCent < 90) {
1364  histname = "ClusterHistograms/hEoverPM02Peripheral";
1365  fHistManager.FillTH3(histname, clus->GetNonLinCorrEnergy(), EoverP, clus->GetM02());
1366  }
1367 
1370 
1371  // Fill Rcorr distribution (centrality, trackPSum, Rcorr = (Enonlincorr - Ehadcorr) / trackPSum)
1372  histname = "ClusterHistograms/hRcorrVsCent";
1373  fHistManager.FillTH3(histname, fCent, trackPSum, Rcorr);
1374 
1375  // Fill Rcorr distribution for 0-10% centrality (Eclus nonlincorr, trackPSum, Rcorr)
1376  if (fCent > 0 && fCent < 10) {
1377  histname = "ClusterHistograms/hRcorr0-10";
1378  fHistManager.FillTH3(histname, clus->GetNonLinCorrEnergy(), trackPSum, Rcorr);
1379  }
1380 
1381  // Fill Rcorr distribution for 50-90% centrality (Eclus nonlincorr, trackPSum, Rcorr)
1382  if (fCent > 50 && fCent < 90) {
1383  histname = "ClusterHistograms/hRcorr50-90";
1384  fHistManager.FillTH3(histname, clus->GetNonLinCorrEnergy(), trackPSum, Rcorr);
1385  }
1386 
1387  // Fill also Rcorr-clus (centrality, trackPSum, Rcorr-clus = (Enonlincorr - Ehadcorr) / Enonlincorr )
1388  histname = "ClusterHistograms/hRcorrClusVsCent";
1389  fHistManager.FillTH3(histname, fCent, trackPSum, RcorrClus);
1390 
1391  // Fill Rcorr-clus for 0-10% centrality (Eclus nonlincorr, trackPSum, Rcorr-clus = (Enonlincorr - Ehadcorr) / Enonlincorr )
1392  if (fCent > 0 && fCent < 10) {
1393  histname = "ClusterHistograms/hRcorrClus0-10";
1394  fHistManager.FillTH3(histname, clus->GetNonLinCorrEnergy(), trackPSum, RcorrClus);
1395  }
1396 
1397  // Fill Rcorr-clus for 50-90% centrality (Eclus nonlincorr, trackPSum, Rcorr-clus = (Enonlincorr - Ehadcorr) / Enonlincorr )
1398  if (fCent > 50 && fCent < 90) {
1399  histname = "ClusterHistograms/hRcorrClus50-90";
1400  fHistManager.FillTH3(histname, clus->GetNonLinCorrEnergy(), trackPSum, RcorrClus);
1401  }
1402 
1403  }
1404 
1405 }
1406 
1407 /*
1408  * This function fills particle composition histograms for the calorimeter performance study in MC.
1409  */
1411 {
1412  // If MC, get the MC event
1413  const AliMCEvent* mcevent = nullptr;
1414  if (fGeneratorLevel) {
1415  mcevent = MCEvent();
1416  }
1417  else {
1418  return;
1419  }
1420 
1421  // Loop through clusters, and plot M02 for each particle type
1423 
1424  // Loop through jets, to fill various histograms
1425  if (fPlotJetHistograms) {
1427  }
1428 
1429 }
1430 
1431 /*
1432  * Loop through clusters, and plot M02 for each particle type
1433  */
1435 {
1436  TString histname;
1438  const AliVCluster* clus;
1439  std::vector<ContributorType> vecContributorTypes;
1440  std::vector<Int_t> vecContributorLabels;
1441  for (auto it : clusters->accepted_momentum()) {
1442 
1443  clus = it.second;
1444 
1445  // Include only EMCal clusters (reject DCal and PHOS clusters)
1446  if (!clus->IsEMCAL()) {
1447  continue;
1448  }
1449  if (it.first.Phi_0_2pi() > fgkEMCalDCalPhiDivide) {
1450  continue;
1451  }
1452 
1453  // Loop through the cluster's contributors in order to classify its type
1454  ParticleType particleType = kNotDefined;
1455  ContributorType contributorType = kUndefined;
1456  const Int_t nLabels = clus->GetNLabels();
1457 
1458  // Create a vector to store the contributor types for PhysicalPrimary particles
1459  vecContributorTypes.clear();
1460  vecContributorLabels.clear();
1461  for (Int_t iLabel=0; iLabel<nLabels; iLabel++) {
1462 
1463  Int_t label = clus->GetLabels()[iLabel];
1464  if (TMath::Abs(label) > 0) { // if the particle has a truth-level match, the label is nonzero
1465  contributorType = GetContributorType(clus, mcevent, label);
1466  if (contributorType != kUndefined) {
1467  vecContributorTypes.push_back(contributorType);
1468  vecContributorLabels.push_back(label);
1469  }
1470  }
1471  }
1472 
1473  Int_t nLabelsPhysPrim = vecContributorTypes.size();
1474 
1475  if (nLabelsPhysPrim == 1) {
1476 
1477  contributorType = vecContributorTypes[0];
1478 
1479  if (contributorType == kPhoton) {
1480  particleType = kSinglePhoton;
1481  }
1482  else if (contributorType == kElectron) {
1483  particleType = kSingleElectron;
1484  }
1485  else if (contributorType == kChargedPion) {
1486  particleType = kSingleChargedPion;
1487  }
1488  else if (contributorType == kProton) {
1489  particleType = kSingleProton;
1490  }
1491  else if (contributorType == kAntiProton) {
1492  particleType = kSingleAntiProton;
1493  }
1494  else if (contributorType == kChargedKaon) {
1495  particleType = kSingleChargedKaon;
1496  }
1497  else if (contributorType == kK0L) {
1498  particleType = kSingleK0L;
1499  }
1500  else if (contributorType == kNeutron) {
1501  particleType = kSingleNeutron;
1502  }
1503  else if (contributorType == kAntiNeutron) {
1504  particleType = kSingleAntiNeutron;
1505  }
1506  else {
1507  particleType = kSingleOther;
1508  }
1509 
1510  }
1511  else if (nLabelsPhysPrim == 2) {
1512 
1513  ContributorType contributorType1 = vecContributorTypes[0];
1514  ContributorType contributorType2 = vecContributorTypes[1];
1515 
1516  Bool_t isHadron1 = IsHadron(contributorType1);
1517  Bool_t isHadron2 = IsHadron(contributorType2);
1518  Bool_t isPhoton1 = contributorType1 == kPhoton;
1519  Bool_t isPhoton2 = contributorType2 == kPhoton;
1520 
1521  if (isHadron1 && isHadron2) {
1522  particleType = kHadronHadron;
1523  }
1524  else if ((isHadron1 && isPhoton2) || (isHadron2 && isPhoton1)) {
1525  particleType = kPhotonHadron;
1526  }
1527  else if (isPhoton1 && isPhoton2) {
1528 
1529  // By default, assume the two photons are not a merged pi0
1530  particleType = kPhotonPhotonOther;
1531 
1532  // Using the vector of accepted contributor labels, check whether the two photons are from the same pi0
1533  AliAODMCParticle *part1 = fGeneratorLevel->GetMCParticleWithLabel(vecContributorLabels[0]);
1534  AliAODMCParticle *part2 = fGeneratorLevel->GetMCParticleWithLabel(vecContributorLabels[1]);
1535  if (part1 && part2) {
1536  Int_t iMother1 = part1->GetMother();
1537  Int_t iMother2 = part2->GetMother();
1538  AliVParticle *mother1 = mcevent->GetTrack(iMother1);
1539  AliVParticle *mother2 = mcevent->GetTrack(iMother2);
1540 
1541  if (mother1 && mother2) {
1542  if ( (mother1->PdgCode() == 111) && (mother2->PdgCode() == 111) ) {
1543  if (iMother1 == iMother2) {
1544  particleType = kMergedPi0;
1545  }
1546  }
1547  }
1548  }
1549  }
1550  else {
1551  particleType = kTwoContributorsOther; // this includes partially contained conversion overlaps
1552  }
1553 
1554  }
1555  else if (nLabelsPhysPrim > 2) {
1556  particleType = kMoreThanTwoContributors;
1557  }
1558 
1559  // (M02, Eclus, part type)
1560  if (fCent > 0 && fCent < 10) {
1561  histname = "ClusterHistogramsMC/hM02VsParticleTypeCentral";
1562  fHistManager.FillTH3(histname, clus->GetM02(), clus->GetNonLinCorrEnergy(), particleType);
1563  }
1564  if (fCent > 50 && fCent < 90) {
1565  histname = "ClusterHistogramsMC/hM02VsParticleTypePeripheral";
1566  fHistManager.FillTH3(histname, clus->GetM02(), clus->GetNonLinCorrEnergy(), particleType);
1567  }
1568 
1569  }
1570 }
1571 
1572 /*
1573  * Loop through jets, to fill particle composition histograms.
1574  */
1576 {
1577  TString histname;
1578 
1579  const AliVCluster* clus;
1580  AliJetContainer* jets = GetJetContainer(0); // there is only a single, det-level jet finder here
1581  for (const auto jet : jets->accepted()) {
1582 
1583  Double_t jetPt = GetJetPt(jet, 0);
1584 
1585  // Keep track of hadronic calo energy in each jet
1586  Double_t hadCaloEnergyUnmatched = 0;
1587  Double_t hadCaloEnergyMatchedNonlincorr = 0;
1588  Double_t hadCaloEnergyMatchedHadCorr = 0;
1589 
1590  // Loop through clusters in each jet, to plot several histograms
1591  Int_t nClusters = jet->GetNumberOfClusters();
1592  for (Int_t iClus = 0; iClus < nClusters; iClus++) {
1593 
1594  clus = jet->Cluster(iClus);
1595 
1596  // Get the particle type of the cluster
1597  ContributorType contributorType = kUndefined;
1598  Int_t label = TMath::Abs(clus->GetLabel());
1599  if (label > 0) {
1600  contributorType = GetContributorType(clus, mcevent, label);
1601  }
1602 
1603  // Plot M02 for each particle type
1604  histname = "JetPerformanceMC/hM02VsContributorTypeJets";
1605  Double_t x[4] = {clus->GetM02(), clus->GetNonLinCorrEnergy(), contributorType, jetPt};
1606  fHistManager.FillTHnSparse(histname, x);
1607 
1608  // If the cluster is a hadron, sum its energy to compute the jet's hadronic calo energy
1609  Bool_t isHadron = IsHadron(contributorType);
1610  if (isHadron) {
1611  Bool_t hasMatchedTrack = (clus->GetNTracksMatched() > 0);
1612  //Bool_t hasMatchedTrack = ((clus->GetNonLinCorrEnergy() - clus->GetHadCorrEnergy()) > 1e-3);
1613  if (hasMatchedTrack) {
1614  hadCaloEnergyMatchedNonlincorr += clus->GetNonLinCorrEnergy();
1615  hadCaloEnergyMatchedHadCorr += clus->GetHadCorrEnergy();
1616  }
1617  else {
1618  hadCaloEnergyUnmatched += clus->GetNonLinCorrEnergy();
1619  }
1620  }
1621 
1622  }
1623 
1624  // Fill hadronic calo energy in each jet
1625 
1626  // (Jet pT, Summed energy of hadronic clusters without a matched track)
1627  histname = "JetPerformance/hHadCaloEnergyUnmatched";
1628  fHistManager.FillTH2(histname, jetPt, hadCaloEnergyUnmatched);
1629 
1630  // (Jet pT vs. Summed energy of hadronic clusters with a matched track (before hadronic correction))
1631  histname = "JetPerformance/hHadCaloEnergyMatchedNonlincorr";
1632  fHistManager.FillTH2(histname, jetPt, hadCaloEnergyMatchedNonlincorr);
1633 
1634  // (Jet pT vs. Summed energy of hadronic clusters with a matched track (after hadronic correction))
1635  histname = "JetPerformance/hHadCaloEnergyMatchedHadCorr";
1636  fHistManager.FillTH2(histname, jetPt, hadCaloEnergyMatchedHadCorr);
1637 
1638  // Loop through particle types, and plot jet composition for each particle type
1639  histname = "JetPerformanceMC/hJetComposition";
1640  for (Int_t type = 0; type < 8; type++) {
1641 
1642  ContributorType contributorType = kUndefined;
1643  Double_t nSum = 0;
1644  Double_t pTsum = 0;
1645 
1646  // Loop through clusters in jet, and add to sum if cluster matches particle type
1647  for (Int_t iClus = 0; iClus < nClusters; iClus++) {
1648 
1649  clus = jet->Cluster(iClus);
1650 
1651  Int_t label = TMath::Abs(clus->GetLabel());
1652  if (label > 0) {
1653  contributorType = GetContributorType(clus, mcevent, label);
1654  }
1655 
1656  if (type == contributorType) {
1657  nSum++;
1658  pTsum += clus->GetNonLinCorrEnergy();
1659  }
1660  }
1661 
1662  // Fill jet composition histogram
1663  Double_t x[4] = {jetPt, 1.*type, nSum, pTsum};
1664  fHistManager.FillTHnSparse(histname, x);
1665 
1666  }
1667  }
1668 }
1669 
1674 {
1675  // Loop over tracks and clusters in order to:
1676  // (1) Compute scale factor for full jets
1677  // (2) Compute delta-pT for full jets, with the random cone method
1678 
1679  // Define the acceptance boundaries for the TPC and EMCal/DCal/PHOS
1680  Double_t etaTPC = 0.9;
1681  Double_t etaEMCal = 0.7;
1682  //Double_t etaMinDCal = 0.22;
1683  Double_t phiMinEMCal = fGeom->GetArm1PhiMin() * TMath::DegToRad(); // 80
1684  Double_t phiMaxEMCal = fGeom->GetEMCALPhiMax() * TMath::DegToRad(); // ~188
1685  //Double_t phiMinDCal = fGeom->GetDCALPhiMin() * TMath::DegToRad(); // 260
1686  //Double_t phiMaxDCal = fGeom->GetDCALPhiMax() * TMath::DegToRad(); // ~327 (1/3 SMs start at 320)
1687 
1688  Double_t accTPC = 2 * etaTPC * 2 * TMath::Pi();
1689  Double_t accEMCal = 2 * etaEMCal * (phiMaxEMCal - phiMinEMCal);
1690  //Double_t accDCalRegion = 2 * etaEMCal * (phiMaxDCal - phiMinDCal);
1691 
1692  // Loop over jet containers
1693  AliJetContainer* jetCont = 0;
1694  TIter nextJetColl(&fJetCollArray);
1695  while ((jetCont = static_cast<AliJetContainer*>(nextJetColl()))) {
1696 
1697  // Define fiducial acceptances, to be used to generate random cones
1698  TRandom3* r = new TRandom3(0);
1699  Double_t jetR = jetCont->GetJetRadius();
1700  Double_t etaEMCalfid = etaEMCal - jetR;
1701  Double_t phiMinEMCalfid = phiMinEMCal + jetR;
1702  Double_t phiMaxEMCalfid = phiMaxEMCal - jetR;
1703 
1704  // Generate EMCal random cone eta-phi
1705  Double_t etaEMCalRC = r->Uniform(-etaEMCalfid, etaEMCalfid);
1706  Double_t phiEMCalRC = r->Uniform(phiMinEMCalfid, phiMaxEMCalfid);
1707 
1708  // Initialize the various sums to 0
1709  Double_t trackPtSumTPC = 0;
1710  Double_t trackPtSumEMCal = 0;
1711  Double_t trackPtSumEMCalRC = 0;
1712  Double_t clusESumEMCal = 0;
1713  Double_t clusESumEMCalRC = 0;
1714 
1715  // Define a 2D vector (initialized to 0) to store the sum of track pT, and another for cluster ET
1716  std::vector<std::vector<Double_t>> trackPtSumDCalRC(fNEtaBins, std::vector<Double_t>(fNPhiBins));
1717  std::vector<std::vector<Double_t>> clusESumDCalRC(fNEtaBins, std::vector<Double_t>(fNPhiBins));
1718 
1719  // Loop over tracks. Sum the track pT:
1720  // (1) in the entire TPC, (2) in the EMCal, (3) in the EMCal random cone,
1721  AliTrackContainer* trackCont = dynamic_cast<AliTrackContainer*>(GetParticleContainer("tracks"));
1722  AliTLorentzVector track;
1723  Double_t trackEta;
1724  Double_t trackPhi;
1725  Double_t trackPt;
1726  Double_t deltaR;
1727  for (auto trackIterator : trackCont->accepted_momentum() ) {
1728 
1729  track.Clear();
1730  track = trackIterator.first;
1731  trackEta = track.Eta();
1732  trackPhi = track.Phi_0_2pi();
1733  trackPt = track.Pt();
1734 
1735  // (1)
1736  if (TMath::Abs(trackEta) < etaTPC) {
1737  trackPtSumTPC += trackPt;
1738  }
1739 
1740  // (2)
1741  if (TMath::Abs(trackEta) < etaEMCal && trackPhi > phiMinEMCal && trackPhi < phiMaxEMCal) {
1742  trackPtSumEMCal += trackPt;
1743  }
1744 
1745  // (3)
1746  deltaR = GetDeltaR(&track, etaEMCalRC, phiEMCalRC);
1747  if (deltaR < jetR) {
1748  trackPtSumEMCalRC += trackPt;
1749  }
1750 
1751  }
1752 
1753  // Loop over clusters. Sum the cluster ET:
1754  // (1) in the EMCal, (2) in the EMCal random cone
1756  AliTLorentzVector clus;
1757  Double_t clusEta;
1758  Double_t clusPhi;
1759  Double_t clusE;
1760  for (auto clusIterator : clusCont->accepted_momentum() ) {
1761 
1762  clus.Clear();
1763  clus = clusIterator.first;
1764  clusEta = clus.Eta();
1765  clusPhi = clus.Phi_0_2pi();
1766  clusE = clus.E();
1767 
1768  // (1)
1769  if (TMath::Abs(clusEta) < etaEMCal && clusPhi > phiMinEMCal && clusPhi < phiMaxEMCal) {
1770  clusESumEMCal += clusE;
1771  }
1772 
1773  // (2)
1774  deltaR = GetDeltaR(&clus, etaEMCalRC, phiEMCalRC);
1775  if (deltaR < jetR) {
1776  clusESumEMCalRC += clusE;
1777  }
1778 
1779  }
1780 
1781  // Compute the scale factor for EMCal, as a function of centrality
1782  Double_t numerator = (trackPtSumEMCal + clusESumEMCal) / accEMCal;
1783  Double_t denominator = trackPtSumTPC / accTPC;
1784  Double_t scaleFactor = numerator / denominator;
1785  TString histname = TString::Format("%s/BackgroundHistograms/hScaleFactorEMCal", jetCont->GetArrayName().Data());
1786  fHistManager.FillTH2(histname, fCent, scaleFactor);
1787 
1788  // Compute delta pT for EMCal, as a function of centrality
1789  Double_t rho = jetCont->GetRhoVal();
1790  Double_t deltaPt = trackPtSumEMCalRC + clusESumEMCalRC - rho * TMath::Pi() * jetR * jetR;
1791  histname = TString::Format("%s/BackgroundHistograms/hDeltaPtEMCal", jetCont->GetArrayName().Data());
1792  fHistManager.FillTH2(histname, fCent, deltaPt);
1793 
1794  delete r;
1795 
1796  }
1797 
1798 }
1799 
1805 {
1806  TString histname;
1807  AliJetContainer* jets = 0;
1808  TIter nextJetColl(&fJetCollArray);
1809  while ((jets = static_cast<AliJetContainer*>(nextJetColl()))) {
1810  TString jetContName = jets->GetName();
1811 
1812  Double_t rhoVal = 0;
1813  if (jets->GetRhoParameter()) {
1814  rhoVal = jets->GetRhoVal();
1815  histname = TString::Format("%s/TriggerSimHistograms/hRhoVsCent", jets->GetArrayName().Data());
1816  fHistManager.FillTH2(histname.Data(), fCent, rhoVal);
1817  }
1818 
1819  for (auto jet : jets->all()) {
1820 
1821  Float_t ptLeading = jets->GetLeadingHadronPt(jet);
1822  Float_t corrPt = GetJetPt(jet, rhoVal);
1823 
1824  // A vs. pT (fill before area cut)
1825  histname = TString::Format("%s/TriggerSimHistograms/hAreaVsPt", jets->GetArrayName().Data());
1826  fHistManager.FillTH2(histname.Data(), corrPt, jet->Area());
1827 
1828  // Rejection reason
1829  UInt_t rejectionReason = 0;
1830  if (!jets->AcceptJet(jet, rejectionReason)) {
1831  histname = TString::Format("%s/TriggerSimHistograms/hJetRejectionReason", jets->GetArrayName().Data());
1832  fHistManager.FillTH2(histname.Data(), jets->GetRejectionReasonBitPosition(rejectionReason), jet->Pt());
1833  continue;
1834  }
1835 
1836  // compute jet acceptance type
1837  Double_t type = GetJetType(jet);
1838  if ( (type != kEMCal) && (type != kDCal) ) {
1839  continue;
1840  }
1841 
1842  // (Centrality, pT, NEF, calo type)
1843  if (type == kEMCal) {
1844  histname = TString::Format("%s/TriggerSimHistograms/hNEFVsPtEMCal", jets->GetArrayName().Data());
1845  }
1846  else if (type == kDCal) {
1847  histname = TString::Format("%s/TriggerSimHistograms/hNEFVsPtDCal", jets->GetArrayName().Data());
1848  }
1849  fHistManager.FillTH3(histname, fCent, corrPt, jet->NEF());
1850 
1851  // pT-leading vs. pT
1852  histname = TString::Format("%s/TriggerSimHistograms/hPtLeadingVsPt", jets->GetArrayName().Data());
1853  fHistManager.FillTH2(histname.Data(), corrPt, ptLeading);
1854 
1855  // (Centrality, pT, z-leading (charged), calo type)
1856  if (type == kEMCal) {
1857  histname = TString::Format("%s/TriggerSimHistograms/hZLeadingVsPtEMCal", jets->GetArrayName().Data());
1858  }
1859  else if (type == kDCal) {
1860  histname = TString::Format("%s/TriggerSimHistograms/hZLeadingVsPtDCal", jets->GetArrayName().Data());
1861  }
1862  TLorentzVector leadPart;
1863  jets->GetLeadingHadronMomentum(leadPart, jet);
1864  Double_t z = GetParallelFraction(leadPart.Vect(), jet);
1865  if (z == 1 || (z > 1 && z - 1 < 1e-3)) z = 0.999; // so that it will contribute to the bin <1
1866  fHistManager.FillTH3(histname, fCent, corrPt, z);
1867 
1868  // (Centrality, pT, z (charged), calo type)
1869  if (type == kEMCal) {
1870  histname = TString::Format("%s/TriggerSimHistograms/hZVsPtEMCal", jets->GetArrayName().Data());
1871  }
1872  else if (type == kDCal) {
1873  histname = TString::Format("%s/TriggerSimHistograms/hZVsPtDCal", jets->GetArrayName().Data());
1874  }
1875  const AliVTrack* track;
1876  for (Int_t i=0; i<jet->GetNumberOfTracks(); i++) {
1877  track = static_cast<AliVTrack*>(jet->Track(i));
1878  z = track->Pt() / TMath::Abs(corrPt);
1879  fHistManager.FillTH3(histname, fCent, corrPt, z);
1880  }
1881 
1882  // (Centrality, pT, Nconst)
1883  if (type == kEMCal) {
1884  histname = TString::Format("%s/TriggerSimHistograms/hNConstVsPtEMCal", jets->GetArrayName().Data());
1885  }
1886  else if (type == kDCal) {
1887  histname = TString::Format("%s/TriggerSimHistograms/hNConstVsPtDCal", jets->GetArrayName().Data());
1888  }
1889  fHistManager.FillTH3(histname, fCent, corrPt, 1.*jet->GetNumberOfConstituents());
1890 
1891  } //jet loop
1892 
1893  }
1894 }
1895 
1901 {
1902  TString histname;
1903  AliJetContainer* jets = 0;
1904  const AliEmcalJet* matchedJet = nullptr;
1905  TIter nextJetColl(&fJetCollArray);
1906  while ((jets = static_cast<AliJetContainer*>(nextJetColl()))) {
1907  TString jetContName = jets->GetName();
1908 
1909  // Only loop over jets in the detector-level jet container
1910  if (jetContName.Contains("mcparticles")) {
1911  continue;
1912  }
1913 
1914  Double_t rhoVal = 0;
1915  if (jets->GetRhoParameter()) {
1916  rhoVal = jets->GetRhoVal();
1917  }
1918 
1919  for (auto jet : jets->accepted()) {
1920 
1921  // Get the matched jet, if it exists
1922  matchedJet = jet->MatchedJet();
1923  if (!matchedJet) {
1924  continue;
1925  }
1926 
1927  // compute jet acceptance type
1928  Double_t type = GetJetType(jet);
1929  if ( (type != kEMCal) && (type != kDCal) ) {
1930  continue;
1931  }
1932 
1933  Float_t detPt = GetJetPt(jet, rhoVal);
1934  Float_t truthPt = matchedJet->Pt();
1935 
1936  // Fill response matrix (centrality, pT-truth, pT-det)
1937  if (type == kEMCal) {
1938  histname = "MatchedJetHistograms/hResponseMatrixEMCal";
1939  }
1940  else if (type == kDCal) {
1941  histname = "MatchedJetHistograms/hResponseMatrixDCal";
1942  }
1943  fHistManager.FillTH3(histname, fCent, truthPt, detPt);
1944 
1945  // Fill JES shift (centrality, pT-truth, (pT-det - pT-truth) / pT-truth)
1946  if (type == kEMCal) {
1947  histname = "MatchedJetHistograms/hJESshiftEMCal";
1948  }
1949  else if (type == kDCal) {
1950  histname = "MatchedJetHistograms/hJESshiftDCal";
1951  }
1952  fHistManager.FillTH3(histname, fCent, truthPt, (detPt-truthPt)/truthPt );
1953 
1954  // Fill NEF of det-level matched jets (centrality, pT-truth, NEF)
1955  histname = "MatchedJetHistograms/hNEFVsPt";
1956  fHistManager.FillTH3(histname, fCent, truthPt, jet->NEF());
1957 
1958  // Fill z-leading (charged) of det-level matched jets (centrality, pT-truth, z-leading)
1959  histname = "MatchedJetHistograms/hZLeadingVsPt";
1960  TLorentzVector leadPart;
1961  jets->GetLeadingHadronMomentum(leadPart, jet);
1962  Double_t z = GetParallelFraction(leadPart.Vect(), jet);
1963  if (z == 1 || (z > 1 && z - 1 < 1e-3)) z = 0.999; // so that it will contribute to the bin <1
1964  fHistManager.FillTH3(histname, fCent, truthPt, z);
1965 
1966  // Fill matching distance (centrality, pT-truth, R)
1967  histname = "MatchedJetHistograms/hMatchingDistance";
1968  fHistManager.FillTH3(histname, fCent, truthPt, jet->ClosestJetDistance());
1969 
1970  } //jet loop
1971  }
1972 }
1973 
1974 /*
1975  * Compute the MC particle type of a given cluster contributor, using the MC particle container
1976  */
1978 {
1979  ContributorType contributorType = kUndefined;
1980 
1981  AliAODMCParticle *part = fGeneratorLevel->GetMCParticleWithLabel(label);
1982  if (part) {
1983 
1984  TString histname = "ClusterHistogramsMC/hClusterRejectionReasonMC";
1985  UInt_t rejectionReason = 0;
1986  if (!fGeneratorLevel->AcceptMCParticle(part, rejectionReason)) {
1987  fHistManager.FillTH2(histname, fGeneratorLevel->GetRejectionReasonBitPosition(rejectionReason), clus->GetNonLinCorrEnergy());
1988  return contributorType;
1989  }
1990 
1991  if (part->GetGeneratorIndex() == 0) { // generator index in cocktail
1992 
1993  // select charged pions, protons, kaons, electrons, muons
1994  Int_t pdg = part->PdgCode();
1995 
1996  if (pdg == 22) { // gamma 22
1997  contributorType = kPhoton;
1998  }
1999  else if (TMath::Abs(pdg) == 211) { // pi+ 211 (abs value ensures both particles and antiparticles are included)
2000  contributorType = kChargedPion;
2001  }
2002  else if (pdg == 2212) { // proton 2212
2003  contributorType = kProton;
2004  }
2005  else if (pdg == -2212) {
2006  contributorType = kAntiProton;
2007  }
2008  else if (TMath::Abs(pdg) == 321) { // K+ 321
2009  contributorType = kChargedKaon;
2010  }
2011  else if (pdg == 130) { // K0L 130
2012  contributorType = kK0L;
2013  }
2014  else if (pdg == 2112) { // neutron 2112
2015  contributorType = kNeutron;
2016  }
2017  else if (pdg == -2112) {
2018  contributorType = kAntiNeutron;
2019  }
2020  else if (TMath::Abs(pdg) == 11) { // e- 11
2021  contributorType = kElectron;
2022  }
2023  else if (TMath::Abs(pdg) == 13) { // mu- 13
2024  contributorType = kMuon;
2025  }
2026  else {
2027  contributorType = kOther;
2028  }
2029  }
2030  }
2031  return contributorType;
2032 }
2033 
2038 {
2039  Double_t deltaPhi = TMath::Abs(part->Phi_0_2pi() - phiRef);
2040  Double_t deltaEta = TMath::Abs(part->Eta() - etaRef);
2041  Double_t deltaR = TMath::Sqrt( deltaPhi*deltaPhi + deltaEta*deltaEta );
2042  return deltaR;
2043 }
2044 
2049 {
2050  UInt_t jetType = jet->GetJetAcceptanceType();
2051  Double_t type = -1;
2052  if (jetType & AliEmcalJet::kEMCAL) {
2053  type = kEMCal;
2054  }
2055  else if (jetType & AliEmcalJet::kDCALonly) {
2056  type = kDCal;
2057  }
2058 
2059  return type;
2060 }
2061 
2066 {
2067  Double_t pT = jet->Pt() - rho * jet->Area();
2068  return pT;
2069 }
2070 
2075 {
2076  return (contributor == kChargedPion) || (contributor == kProton) || (contributor == kAntiProton) || (contributor == kChargedKaon) || (contributor == kK0L) || (contributor == kNeutron) || (contributor == kAntiNeutron);
2077 }
2078 
2083  const char *ntracks,
2084  const char *nclusters,
2085  const char *nGenLev,
2086  const Double_t minTrPt,
2087  const Double_t minClPt,
2088  const char *suffix)
2089 {
2090  // Get the pointer to the existing analysis manager via the static access method.
2091  //==============================================================================
2092  AliAnalysisManager *mgr = AliAnalysisManager::GetAnalysisManager();
2093  if (!mgr)
2094  {
2095  ::Error("AddTaskEmcalJetPerformance", "No analysis manager to connect to.");
2096  return 0;
2097  }
2098 
2099  // Check the analysis type using the event handlers connected to the analysis manager.
2100  //==============================================================================
2101  AliVEventHandler* handler = mgr->GetInputEventHandler();
2102  if (!handler)
2103  {
2104  ::Error("AddTaskEmcalJetPerformance", "This task requires an input event handler");
2105  return 0;
2106  }
2107 
2108  enum EDataType_t {
2109  kUnknown,
2110  kESD,
2111  kAOD
2112  };
2113 
2114  EDataType_t dataType = kUnknown;
2115 
2116  if (handler->InheritsFrom("AliESDInputHandler")) {
2117  dataType = kESD;
2118  }
2119  else if (handler->InheritsFrom("AliAODInputHandler")) {
2120  dataType = kAOD;
2121  }
2122 
2123  //-------------------------------------------------------
2124  // Init the task and do settings
2125  //-------------------------------------------------------
2126 
2127  TString trackName(ntracks);
2128  TString clusName(nclusters);
2129 
2130  if (trackName == "usedefault") {
2131  if (dataType == kESD) {
2132  trackName = "Tracks";
2133  }
2134  else if (dataType == kAOD) {
2135  trackName = "tracks";
2136  }
2137  else {
2138  trackName = "";
2139  }
2140  }
2141 
2142  if (clusName == "usedefault") {
2143  if (dataType == kESD) {
2144  clusName = "CaloClusters";
2145  }
2146  else if (dataType == kAOD) {
2147  clusName = "caloClusters";
2148  }
2149  else {
2150  clusName = "";
2151  }
2152  }
2153 
2154  TString name("AliAnalysisTaskEmcalJetPerformance");
2155  if (!trackName.IsNull()) {
2156  name += "_";
2157  name += trackName;
2158  }
2159  if (!clusName.IsNull()) {
2160  name += "_";
2161  name += clusName;
2162  }
2163  if (strcmp(suffix,"") != 0) {
2164  name += "_";
2165  name += suffix;
2166  }
2167 
2169  // Configure di-jet task
2171 
2173  // Create track and cluster containers with the standard cuts
2174 
2175  AliParticleContainer* partCont = 0;
2176  if (trackName == "mcparticles") {
2177  AliMCParticleContainer* mcpartCont = new AliMCParticleContainer(trackName);
2178  partCont = mcpartCont;
2179  }
2180  else if (trackName == "tracks" || trackName == "Tracks") {
2181  AliTrackContainer* trackCont = new AliTrackContainer(trackName);
2182  partCont = trackCont;
2183  }
2184  if (partCont) partCont->SetParticlePtCut(minTrPt);
2185  if (partCont) task->AdoptParticleContainer(partCont);
2186 
2187  // Add the generator-level container, if specified
2188  if (nGenLev && strcmp(nGenLev,"")!=0) {
2189  AliMCParticleContainer* mcpartCont = task->AddMCParticleContainer(nGenLev);
2190  mcpartCont->SelectPhysicalPrimaries(kTRUE);
2191  mcpartCont->SetParticlePtCut(0);
2192  }
2193 
2194  // Add the cluster container
2195  AliClusterContainer* clusCont = 0;
2196  if (!clusName.IsNull()) {
2197  clusCont = new AliClusterContainer(clusName);
2198  clusCont->SetClusECut(0.);
2199  clusCont->SetClusPtCut(0.);
2200  }
2201  if (clusCont) task->AdoptClusterContainer(clusCont);
2202 
2203  //-------------------------------------------------------
2204  // Final settings, pass to manager and set the containers
2205  //-------------------------------------------------------
2206 
2207  mgr->AddTask(task);
2208 
2209  // Create containers for input/output
2210  AliAnalysisDataContainer *cinput1 = mgr->GetCommonInputContainer() ;
2211  TString contname(name);
2212  contname += "_histos";
2213  AliAnalysisDataContainer *coutput1 = mgr->CreateContainer(contname.Data(),
2214  TList::Class(),AliAnalysisManager::kOutputContainer,
2215  Form("%s", AliAnalysisManager::GetCommonFileName()));
2216  mgr->ConnectInput (task, 0, cinput1 );
2217  mgr->ConnectOutput (task, 1, coutput1 );
2218 
2219  return task;
2220 }
2221 
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
Bool_t fUseAliEventCuts
Flag to use AliEventCuts (otherwise AliAnalysisTaskEmcal will be used)
AliJetContainer * GetJetContainer(Int_t i=0) const
Definition: External.C:244
UInt_t fOffTrigger
offline trigger for event selection
void AdoptParticleContainer(AliParticleContainer *cont)
static AliEmcalDownscaleFactorsOCDB * Instance()
const AliClusterIterableMomentumContainer 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.
AliEmcalJet * MatchedJet() const
Definition: AliEmcalJet.h:331
Declaration of class AliTLorentzVector.
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)
virtual Bool_t AcceptMCParticle(const AliAODMCParticle *vp, UInt_t &rejectionReason) const
bool AddQAPlotsToList(TList *list)
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
Double_t fTrackMatchingDeltaEtaMax
Maximum delta-eta to consider a track to be matched to a cluster.
unsigned int UInt_t
Definition: External.C:33
THashList * GetListOfHistograms() const
Get the list of histograms.
Definition: THistManager.h:671
Double_t GetLeadingHadronPt(const AliEmcalJet *jet) const
float Float_t
Definition: External.C:68
AliEMCALGeometry * fGeom
!emcal geometry
UInt_t GetJetAcceptanceType() const
Definition: AliEmcalJet.h:367
Double_t Phi_0_2pi() const
EMCal acceptance.
Definition: AliEmcalJet.h:69
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.
Bool_t fUseManualEventCuts
Flag to use manual event cuts.
AliClusterContainer * GetClusterContainer(Int_t i=0) const
Get cluster container attached to this task.
std::vector< TString > GetTriggerClasses() const
Double_t fCent
!event centrality
Photon (direct or decay) is the only contributor.
Two particles from merged pi0 are the only contributors.
Bool_t fPlotJetHistograms
Set whether to enable inclusive jet histograms.
Int_t fNEoverPBins
! number of variable E/p bins
AliMCParticleContainer * AddMCParticleContainer(const char *n)
Create new container for MC particles and attach it to the task.
static Double_t GetParallelFraction(AliVParticle *part1, AliVParticle *part2)
Calculates the fraction of momentum z of part 1 w.r.t. part 2 in the direction of part 2...
void FillTH1(const char *hname, double x, double weight=1., Option_t *opt="")
Fill a 1D histogram within the container.
AliEventCuts fEventCuts
event selection utility
Bool_t fPlotMatchedJetHistograms
Set whether to plot matched jet histograms (must run ResponseMaker first)
TObjArray fJetCollArray
jet collection array
Bool_t fkEMCEJE
! flag telling whether the event is "triggered" or not in "simulation"
AliRhoParameter * GetRhoParameter()
Double_t Pt() const
Definition: AliEmcalJet.h:109
Int_t nbinsx
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.
AliEmcalList * fOutput
!output list
ContributorType GetContributorType(const AliVCluster *clus, const AliMCEvent *mcevent, Int_t label)
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.
Handler for downscale factors for various triggers obtained from the OCDB.
void SelectPhysicalPrimaries(Bool_t s)
TClonesArray * fTriggerPatchInfo
!trigger patch info array
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
const AliTrackIterableMomentumContainer accepted_momentum() const
Int_t GetRunNumber(TString)
Definition: PlotMuonQA.C:2235
void SetClusPtCut(Double_t cut)
void SetRejectionReasonLabels(TAxis *axis)
void UserCreateOutputObjects()
Main initialization function on the worker.
Bool_t fDoTriggerSimulation
Set whether to perform a simple trigger simulation.
const Int_t nbins
const AliJetIterableContainer accepted() const
bool Bool_t
Definition: External.C:53
EDataType_t
Switch for the data type.
DCal acceptance – spans ONLY DCal (no PHOS or gap)
Definition: AliEmcalJet.h:73
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.
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
Container for jet within the EMCAL jet framework.
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.
const AliJetIterableContainer all() const
static const AliAnalysisTaskEmcalEmbeddingHelper * GetInstance()
static Double_t fgkEMCalDCalPhiDivide
phi value used to distinguish between DCal and EMCal