AliPhysics  fceccc5 (fceccc5)
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
AliAnalysisTaskEmcalDijetImbalance.cxx
Go to the documentation of this file.
1 /**************************************************************************
2  * Copyright(c) 1998-2016, ALICE Experiment at CERN, All rights reserved. *
3  * *
4  * Author: The ALICE Off-line Project. *
5  * Contributors are mentioned in the code where appropriate. *
6  * *
7  * Permission to use, copy, modify and distribute this software and its *
8  * documentation strictly for non-commercial purposes is hereby granted *
9  * without fee, provided that the above copyright notice appears in all *
10  * copies and that both the copyright notice and this permission notice *
11  * appear in the supporting documentation. The authors make no claims *
12  * about the suitability of this software for any purpose. It is *
13  * provided "as is" without express or implied warranty. *
14  **************************************************************************/
15 
16 #include <TClonesArray.h>
17 #include <TH1F.h>
18 #include <TH2F.h>
19 #include <TList.h>
20 #include <THnSparse.h>
21 #include <TRandom3.h>
22 
23 #include <AliVCluster.h>
24 #include <AliVParticle.h>
25 #include <AliLog.h>
26 
27 #include "AliTLorentzVector.h"
28 #include "AliEmcalJet.h"
29 #include "AliRhoParameter.h"
30 #include "AliJetContainer.h"
31 #include "AliParticleContainer.h"
32 #include "AliClusterContainer.h"
33 #include "AliEMCALGeometry.h"
35 
37 
41 
47  fHistManager(),
48  fEventCuts(0),
49  fEventCutList(0),
50  fUseManualEventCuts(kFALSE),
51  fDeltaPhiMin(0),
52  fMinTrigJetPt(0),
53  fMinAssJetPt(0),
54  fDijetLeadingHadronPt(0),
55  fMaxPt(200),
56  fNCentHistBins(0),
57  fCentHistBins(0),
58  fPlotJetHistograms(kFALSE),
59  fPlotDijetCandHistograms(kFALSE),
60  fPlotDijetImbalanceHistograms(kFALSE),
61  fComputeBackground(kFALSE),
62  fDoMomentumBalance(kFALSE),
63  fDoGeometricalMatching(kFALSE),
64  fMatchingJetR(0.2),
65  fTrackConstituentThreshold(0),
66  fClusterConstituentThreshold(0),
67  fMBUpscaleFactor(1.)
68 {
72 }
73 
80  AliAnalysisTaskEmcalJet(name, kTRUE),
81  fHistManager(name),
82  fEventCuts(0),
83  fEventCutList(0),
84  fUseManualEventCuts(kFALSE),
85  fDeltaPhiMin(0),
86  fMinTrigJetPt(0),
87  fMinAssJetPt(0),
88  fDijetLeadingHadronPt(0),
89  fMaxPt(200),
90  fNCentHistBins(0),
91  fCentHistBins(0),
92  fPlotJetHistograms(kFALSE),
93  fPlotDijetCandHistograms(kFALSE),
94  fPlotDijetImbalanceHistograms(kFALSE),
95  fComputeBackground(kFALSE),
96  fDoMomentumBalance(kFALSE),
97  fDoGeometricalMatching(kFALSE),
98  fMatchingJetR(0.2),
99  fTrackConstituentThreshold(0),
100  fClusterConstituentThreshold(0),
101  fMBUpscaleFactor(1.)
102 {
104  Dijet_t fDijet;
106 }
107 
112 {
113 }
114 
119 {
120  fNCentHistBins = 4;
122  fCentHistBins[0] = 0;
123  fCentHistBins[1] = 10;
124  fCentHistBins[2] = 30;
125  fCentHistBins[3] = 50;
126  fCentHistBins[4] = 90;
127 }
128 
134 {
136 
142 
143  TIter next(fHistManager.GetListOfHistograms());
144  TObject* obj = 0;
145  while ((obj = next())) {
146  fOutput->Add(obj);
147  }
148 
149  // Intialize AliEventCuts
150  fEventCutList = new TList();
151  fEventCutList ->SetOwner();
152  fEventCutList ->SetName("EventCutOutput");
153 
154  fEventCuts.OverrideAutomaticTriggerSelection(fOffTrigger);
155  if(fUseManualEventCuts==1)
156  {
157  fEventCuts.SetManualMode();
158  // Configure manual settings here
159  // ...
160  }
161  fEventCuts.AddQAplotsToList(fEventCutList);
162  fOutput->Add(fEventCutList);
163 
164  PostData(1, fOutput); // Post data for ALL output slots > 0 here.
165 }
166 
167 /*
168  * This function allocates the histograms for single jets.
169  * A set of histograms is allocated per each jet container.
170  */
172 {
173  TString histname;
174  TString title;
175 
176  Int_t nPtBins = TMath::CeilNint(fMaxPt/2);
177 
178  AliJetContainer* jets = 0;
179  TIter nextJetColl(&fJetCollArray);
180  while ((jets = static_cast<AliJetContainer*>(nextJetColl()))) {
181 
182  // Jet rejection reason
183  histname = TString::Format("%s/JetHistograms/hJetRejectionReason", jets->GetArrayName().Data());
184  title = histname + ";Rejection reason;#it{p}_{T,jet} (GeV/#it{c});counts";
185  TH2* hist = fHistManager.CreateTH2(histname.Data(), title.Data(), 32, 0, 32, 50, 0, 250);
186  SetRejectionReasonLabels(hist->GetXaxis());
187 
188  // Rho vs. Centrality
189  if (!jets->GetRhoName().IsNull()) {
190  histname = TString::Format("%s/JetHistograms/hRhoVsCent", jets->GetArrayName().Data());
191  title = histname + ";Centrality (%);#rho (GeV/#it{c});counts";
192  fHistManager.CreateTH2(histname.Data(), title.Data(), 101, 0, 101, 100, 0, 500);
193  }
194 
195  // Centrality vs. pT
196  histname = TString::Format("%s/JetHistograms/hCentVsPt", jets->GetArrayName().Data());
197  title = histname + ";#it{p}_{T}^{corr} (GeV/#it{c});Centrality (%);counts";
198  fHistManager.CreateTH2(histname.Data(), title.Data(), nPtBins, 0, fMaxPt, 10, 0, 100);
199 
200  // pT vs. eta vs. phi
201  histname = TString::Format("%s/JetHistograms/hPtVsEtaVsPhi", jets->GetArrayName().Data());
202  title = histname + ";#eta_{jet} (rad);#phi_{jet} (rad);#it{p}_{T}^{corr} (GeV/#it{c})";
203  fHistManager.CreateTH3(histname.Data(), title.Data(), 50, -0.5, 0.5, 101, 0, TMath::Pi() * 2.02, nPtBins, 0, fMaxPt);
204 
205  // pT upscaled
206  histname = TString::Format("%s/JetHistograms/hPtUpscaledMB", jets->GetArrayName().Data());
207  title = histname + ";#it{p}_{T}^{corr} (GeV/#it{c})";
208  fHistManager.CreateTH1(histname.Data(), title.Data(), nPtBins, 0, fMaxPt, "s");
209 
210  // pT-leading vs. pT
211  histname = TString::Format("%s/JetHistograms/hPtLeadingVsPt", jets->GetArrayName().Data());
212  title = histname + ";#it{p}_{T}^{corr} (GeV/#it{c});#it{p}_{T,particle}^{leading} (GeV/#it{c})";
213  fHistManager.CreateTH2(histname.Data(), title.Data(), nPtBins, 0, fMaxPt, nPtBins, 0, fMaxPt);
214 
215  // A vs. pT
216  histname = TString::Format("%s/JetHistograms/hAreaVsPt", jets->GetArrayName().Data());
217  title = histname + ";#it{p}_{T}^{corr} (GeV/#it{c});#it{A}_{jet}";
218  fHistManager.CreateTH2(histname.Data(), title.Data(), nPtBins, 0, fMaxPt, fMaxPt/3, 0, 1.5);
219 
220  // NEF vs. pT
221  histname = TString::Format("%s/JetHistograms/hNEFVsPt", jets->GetArrayName().Data());
222  title = histname + ";#it{p}_{T}^{corr} (GeV/#it{c});NEF;";
223  fHistManager.CreateTH2(histname.Data(), title.Data(), nPtBins, 0, fMaxPt, fMaxPt/5, 0, 1.0);
224 
225  // z-leading (charged) vs. pT
226  histname = TString::Format("%s/JetHistograms/hZLeadingVsPt", jets->GetArrayName().Data());
227  title = histname + ";#it{p}_{T}^{corr} (GeV/#it{c});#it{z}_{leading}";
228  fHistManager.CreateTH2(histname.Data(), title.Data(), nPtBins, 0, fMaxPt, fMaxPt/5, 0, 1.0);
229 
230  // z (charged) vs. pT
231  histname = TString::Format("%s/JetHistograms/hZVsPt", jets->GetArrayName().Data());
232  title = histname + ";#it{p}_{T}^{corr} (GeV/#it{c});#it{z}";
233  fHistManager.CreateTH2(histname.Data(), title.Data(), nPtBins, 0, fMaxPt, fMaxPt/5, 0, 1.0);
234 
235  // Nconst vs. pT
236  histname = TString::Format("%s/JetHistograms/hNConstVsPt", jets->GetArrayName().Data());
237  title = histname + ";#it{p}_{T}^{corr} (GeV/#it{c});No. of constituents";
238  fHistManager.CreateTH2(histname.Data(), title.Data(), nPtBins, 0, fMaxPt, fMaxPt/5, 0, fMaxPt);
239 
240  // Allocate background subtraction histograms, if enabled
241  if (fComputeBackground) {
242 
243  histname = TString::Format("%s/BackgroundHistograms/hScaleFactorEMCal", jets->GetArrayName().Data());
244  title = histname + ";Centrality;Scale factor;counts";
245  fHistManager.CreateTH2(histname.Data(), title.Data(), 50, 0, 100, 100, 0, 5);
246 
247  histname = TString::Format("%s/BackgroundHistograms/hDeltaPtEMCal", jets->GetArrayName().Data());
248  title = histname + ";#delta#it{p}_{T} (GeV/#it{c});counts";
249  fHistManager.CreateTH1(histname.Data(), title.Data(), 400, -100, 100);
250 
251  histname = TString::Format("%s/BackgroundHistograms/hScaleFactorDCalRegion", jets->GetArrayName().Data());
252  title = histname + ";Centrality;Scale factor;Eta bin";
253  fHistManager.CreateTH3(histname.Data(), title.Data(), 50, 0, 100, 200, 0, 10, 28, -0.5, 27.5);
254 
255  histname = TString::Format("%s/BackgroundHistograms/hDeltaPtDCalRegion", jets->GetArrayName().Data());
256  title = histname + ";#delta#it{p}_{T} (GeV/#it{c});Eta bin;counts";
257  fHistManager.CreateTH2(histname.Data(), title.Data(), 400, -100, 100, 20, -0.5, 19.5);
258 
259  }
260 
261  }
262 
263  // MB downscale factor histogram
264  histname = "Trigger/hMBDownscaleFactor";
265  title = histname + ";Downscale factor;counts";
266  TH1* hist = fHistManager.CreateTH1(histname.Data(), title.Data(), 200, 0, 200);
267 
268 }
269 
270 /*
271  * This function allocates the histograms for dijet candidates, i.e. dijet pairs with the leading jet
272  * passing the trigger condition, but no condition on the subleading jet. In particular this histogram is
273  * designed to study the kinematic selection of dijets.
274  */
276 {
277  // Allocate dijet THnSparse
278  AliJetContainer* jets = 0;
279  TIter nextJetColl(&fJetCollArray);
280  while ((jets = static_cast<AliJetContainer*>(nextJetColl()))) {
281  TString axisTitle[30]= {""};
282  Int_t nbins[30] = {0};
283  Double_t min[30] = {0.};
284  Double_t max[30] = {0.};
285  Double_t *binEdges[20] = {0};
286  Int_t dim = 0;
287 
288  if (fForceBeamType != kpp) {
289  axisTitle[dim] = "Centrality (%)";
290  nbins[dim] = fNCentHistBins;
291  binEdges[dim] = fCentHistBins;
292  min[dim] = fCentHistBins[0];
293  max[dim] = fCentHistBins[fNCentHistBins];
294  dim++;
295  }
296 
297  axisTitle[dim] = "LeadingHadronRequired";
298  nbins[dim] = 2;
299  min[dim] = -0.5;
300  max[dim] = 1.5;
301  binEdges[dim] = GenerateFixedBinArray(nbins[dim], min[dim], max[dim]);
302  dim++;
303 
304  axisTitle[dim] = "#it{p}_{T,trig jet} (GeV/#it{c})";
305  nbins[dim] = TMath::CeilNint(fMaxPt/2);
306  min[dim] = 0;
307  max[dim] = fMaxPt;
308  binEdges[dim] = GenerateFixedBinArray(nbins[dim], min[dim], max[dim]);
309  dim++;
310 
311  axisTitle[dim] = "#it{p}_{T,ass jet} (GeV/#it{c})";
312  nbins[dim] = TMath::CeilNint(fMaxPt/2);
313  min[dim] = 0;
314  max[dim] = fMaxPt;
315  binEdges[dim] = GenerateFixedBinArray(nbins[dim], min[dim], max[dim]);
316  dim++;
317 
318  axisTitle[dim] = "#phi_{trig jet}";
319  nbins[dim] = TMath::CeilNint(fMaxPt/2);
320  min[dim] = 0;
321  max[dim] = TMath::TwoPi();
322  binEdges[dim] = GenerateFixedBinArray(nbins[dim], min[dim], max[dim]);
323  dim++;
324 
325  axisTitle[dim] = "#phi_{ass jet}";
326  nbins[dim] = TMath::CeilNint(fMaxPt/2);
327  min[dim] = 0;
328  max[dim] = TMath::TwoPi();
329  binEdges[dim] = GenerateFixedBinArray(nbins[dim], min[dim], max[dim]);
330  dim++;
331 
332  axisTitle[dim] = "#eta_{trig jet}";
333  nbins[dim] = TMath::CeilNint(fMaxPt/2);
334  min[dim] = -1;
335  max[dim] = 1;
336  binEdges[dim] = GenerateFixedBinArray(nbins[dim], min[dim], max[dim]);
337  dim++;
338 
339  axisTitle[dim] = "#eta_{ass jet}";
340  nbins[dim] = TMath::CeilNint(fMaxPt/2);
341  min[dim] = -1;
342  max[dim] = 1;
343  binEdges[dim] = GenerateFixedBinArray(nbins[dim], min[dim], max[dim]);
344  dim++;
345 
346  TString thnname = TString::Format("%s/DijetCandObservables", jets->GetArrayName().Data());
347  THnSparse* hn = fHistManager.CreateTHnSparse(thnname.Data(), thnname.Data(), dim, nbins, min, max);
348  for (Int_t i = 0; i < dim; i++) {
349  hn->GetAxis(i)->SetTitle(axisTitle[i]);
350  hn->SetBinEdges(i, binEdges[i]);
351  }
352  }
353 }
354 
355 /*
356  * This function allocates the histograms for accepted dijets.
357  * The purpose is to study in detail the imbalance properties of the accepted dijets.
358  */
360 {
361  // Allocate dijet imbalance THnSparse
362  AliJetContainer* jets = 0;
363  TIter nextJetColl(&fJetCollArray);
364  while ((jets = static_cast<AliJetContainer*>(nextJetColl()))) {
365  TString axisTitle[30]= {""};
366  Int_t nbins[30] = {0};
367  Double_t min[30] = {0.};
368  Double_t max[30] = {0.};
369  Double_t *binEdges[20] = {0};
370  Int_t dim = 0;
371 
372  if (fForceBeamType != kpp) {
373  axisTitle[dim] = "Centrality (%)";
374  nbins[dim] = fNCentHistBins;
375  binEdges[dim] = fCentHistBins;
376  min[dim] = fCentHistBins[0];
377  max[dim] = fCentHistBins[fNCentHistBins];
378  dim++;
379  }
380 
381  axisTitle[dim] = "#Delta#phi";
382  nbins[dim] = 100;
383  min[dim] = 0;
384  max[dim] = 4;
385  binEdges[dim] = GenerateFixedBinArray(nbins[dim], min[dim], max[dim]);
386  dim++;
387 
388  axisTitle[dim] = "#Delta#eta";
389  nbins[dim] = 100;
390  min[dim] = -2;
391  max[dim] = 2;
392  binEdges[dim] = GenerateFixedBinArray(nbins[dim], min[dim], max[dim]);
393  dim++;
394 
395  axisTitle[dim] = "A_{J}";
396  nbins[dim] = 100;
397  min[dim] = 0;
398  max[dim] = 1;
399  binEdges[dim] = GenerateFixedBinArray(nbins[dim], min[dim], max[dim]);
400  dim++;
401 
402  axisTitle[dim] = "x_{J}";
403  nbins[dim] = 100;
404  min[dim] = 0;
405  max[dim] = 1;
406  binEdges[dim] = GenerateFixedBinArray(nbins[dim], min[dim], max[dim]);
407  dim++;
408 
409  axisTitle[dim] = "k_{Ty} (GeV)";
410  nbins[dim] = 100;
411  min[dim] = 0;
412  max[dim] = 100;
413  binEdges[dim] = GenerateFixedBinArray(nbins[dim], min[dim], max[dim]);
414  dim++;
415 
416  axisTitle[dim] = "N_{tracks, trig jet}";
417  nbins[dim] = fMaxPt/5;
418  min[dim] = 0;
419  max[dim] = 100;
420  binEdges[dim] = GenerateFixedBinArray(nbins[dim], min[dim], max[dim]);
421  dim++;
422 
423  axisTitle[dim] = "N_{tracks, ass jet}";
424  nbins[dim] = fMaxPt/5;
425  min[dim] = 0;
426  max[dim] = 100;
427  binEdges[dim] = GenerateFixedBinArray(nbins[dim], min[dim], max[dim]);
428  dim++;
429 
430  TString thnname = TString::Format("%s/DijetImbalanceObservables", jets->GetArrayName().Data());
431  THnSparse* hn = fHistManager.CreateTHnSparse(thnname.Data(), thnname.Data(), dim, nbins, min, max);
432  for (Int_t i = 0; i < dim; i++) {
433  hn->GetAxis(i)->SetTitle(axisTitle[i]);
434  hn->SetBinEdges(i, binEdges[i]);
435  }
436  }
437 }
438 
439 /*
440  * This function allocates the histograms for the momentum balance study.
441  */
443 {
444  AliJetContainer* jets = 0;
445  TIter nextJetColl(&fJetCollArray);
446  while ((jets = static_cast<AliJetContainer*>(nextJetColl()))) {
447 
448  // Allocate THnSparse
449  TString axisTitle[30]= {""};
450  Int_t nbins[30] = {0};
451  Double_t min[30] = {0.};
452  Double_t max[30] = {0.};
453  Double_t *binEdges[20] = {0};
454  Int_t dim = 0;
455 
456  axisTitle[dim] = "A_{J}";
457  nbins[dim] = 100;
458  min[dim] = 0;
459  max[dim] = 1;
460  binEdges[dim] = GenerateFixedBinArray(nbins[dim], min[dim], max[dim]);
461  dim++;
462 
463  axisTitle[dim] = "#Delta#phi";
464  nbins[dim] = 100;
465  min[dim] = -4;
466  max[dim] = 4;
467  binEdges[dim] = GenerateFixedBinArray(nbins[dim], min[dim], max[dim]);
468  dim++;
469 
470  axisTitle[dim] = "#it{p}_{T,particle} (GeV/#it{c})";
471  nbins[dim] = 9;
472  Double_t* pTParticleBins = new Double_t[nbins[dim]+1];
473  GenerateFixedBinArray(1, 0.15, 0.3, pTParticleBins);
474  GenerateFixedBinArray(1, 0.3, 0.5, pTParticleBins+1);
475  GenerateFixedBinArray(1, 0.5, 1, pTParticleBins+2);
476  GenerateFixedBinArray(2, 1, 5, pTParticleBins+3);
477  GenerateFixedBinArray(3, 5, 20, pTParticleBins+5);
478  GenerateFixedBinArray(1, 20, 150, pTParticleBins+8);
479  min[dim] = 0;
480  max[dim] = pTParticleBins[nbins[dim]];
481  binEdges[dim] = pTParticleBins;
482  dim++;
483 
484  axisTitle[dim] = "#it{p}_{T}#parallel (GeV/#it{c})";
485  nbins[dim] = 80;
486  Double_t* pTParallelBins = new Double_t[nbins[dim]+1];
487  GenerateFixedBinArray(20, 0, 2, pTParallelBins);
488  GenerateFixedBinArray(16, 2, 10, pTParallelBins+20);
489  GenerateFixedBinArray(10, 10, 20, pTParallelBins+36);
490  GenerateFixedBinArray(10, 20, 40, pTParallelBins+46);
491  GenerateFixedBinArray(24, 40, 150, pTParallelBins+56);
492  min[dim] = 0;
493  max[dim] = pTParallelBins[nbins[dim]];
494  binEdges[dim] = pTParallelBins;
495  dim++;
496 
497  TString thnname = TString::Format("%s/MomentumBalance", jets->GetArrayName().Data());
498  THnSparse* hn = fHistManager.CreateTHnSparse(thnname.Data(), thnname.Data(), dim, nbins, min, max);
499  for (Int_t i = 0; i < dim; i++) {
500  hn->GetAxis(i)->SetTitle(axisTitle[i]);
501  hn->SetBinEdges(i, binEdges[i]);
502  }
503  }
504 }
505 
506 /*
507  * This function allocates the histograms for the constituent dijet study.
508  */
510 {
511  // Allocate geometrical matching THnSparse
512  TString axisTitle[30]= {""};
513  Int_t nbins[30] = {0};
514  Double_t min[30] = {0.};
515  Double_t max[30] = {0.};
516  Double_t *binEdges[20] = {0};
517  Int_t dim = 0;
518 
519  if (fForceBeamType != kpp) {
520  axisTitle[dim] = "Centrality (%)";
521  nbins[dim] = fNCentHistBins;
522  binEdges[dim] = fCentHistBins;
523  min[dim] = fCentHistBins[0];
524  max[dim] = fCentHistBins[fNCentHistBins];
525  dim++;
526  }
527 
528  axisTitle[dim] = "isSwitched";
529  nbins[dim] = 2;
530  min[dim] = -0.5;
531  max[dim] = 1.5;
532  binEdges[dim] = GenerateFixedBinArray(nbins[dim], min[dim], max[dim]);
533  dim++;
534 
535  axisTitle[dim] = "#DeltaR_{trig}";
536  nbins[dim] = 50;
537  min[dim] = 0;
538  max[dim] = 0.5;
539  binEdges[dim] = GenerateFixedBinArray(nbins[dim], min[dim], max[dim]);
540  dim++;
541 
542  axisTitle[dim] = "#DeltaR_{ass}";
543  nbins[dim] = 50;
544  min[dim] = 0;
545  max[dim] = 0.5;
546  binEdges[dim] = GenerateFixedBinArray(nbins[dim], min[dim], max[dim]);
547  dim++;
548 
549  axisTitle[dim] = "trig #it{p}_{T,low-thresh} - #it{p}_{T,hard-core}";
550  nbins[dim] = 100;
551  min[dim] = -50;
552  max[dim] = 50;
553  binEdges[dim] = GenerateFixedBinArray(nbins[dim], min[dim], max[dim]);
554  dim++;
555 
556  axisTitle[dim] = "ass #it{p}_{T,low-thresh} - #it{p}_{T,hard-core}";
557  nbins[dim] = 100;
558  min[dim] = -50;
559  max[dim] = 50;
560  binEdges[dim] = GenerateFixedBinArray(nbins[dim], min[dim], max[dim]);
561  dim++;
562 
563  axisTitle[dim] = "A_{J} low-threshold";
564  nbins[dim] = 100;
565  min[dim] = 0;
566  max[dim] = 1;
567  binEdges[dim] = GenerateFixedBinArray(nbins[dim], min[dim], max[dim]);
568  dim++;
569 
570  axisTitle[dim] = "A_{J} hard-core";
571  nbins[dim] = 100;
572  min[dim] = 0;
573  max[dim] = 1;
574  binEdges[dim] = GenerateFixedBinArray(nbins[dim], min[dim], max[dim]);
575  dim++;
576 
577  TString thnname = "GeometricalMatching";
578  THnSparse* hn = fHistManager.CreateTHnSparse(thnname.Data(), thnname.Data(), dim, nbins, min, max);
579  for (Int_t i = 0; i < dim; i++) {
580  hn->GetAxis(i)->SetTitle(axisTitle[i]);
581  hn->SetBinEdges(i, binEdges[i]);
582  }
583 
584  // Allocate other histograms
585  TString histname;
586  TString title;
587  histname = "GeometricalMatchingEfficiency";
588  title = histname + ";isMatched;counts";
589  TH1* hist = fHistManager.CreateTH1(histname.Data(), title.Data(), 2, -0.5, 1.5);
590 }
591 
597 {
599 
600  fNeedEmcalGeom = kTRUE;
601 
602  AliInfo(Form("Trigger jet threshold = %f, Associated jet threshold = %f", fMinTrigJetPt, fMinAssJetPt));
603  AliInfo(Form("Leading hadron threshold (for dijet leading jet): %f GeV", fDijetLeadingHadronPt));
604  AliInfo(Form("Momentum balance study: %d", fDoMomentumBalance));
605  AliInfo(Form("Geometrical matching study: %d", fDoGeometricalMatching));
606 
607 }
608 
613 
614  // Get the downscaling factors for MB triggers (to be used to calculate trigger efficiency)
615 
616  if (fPlotJetHistograms) {
617 
618  // Get instance of the downscale factor helper class
620  downscaleOCDB->SetRun(InputEvent()->GetRunNumber());
621 
622  // There are two possible min bias triggers for LHC15o
623  TString triggerNameMB1 = "CINT7-B-NOPF-CENT";
624  TString triggerNameMB2 = "CV0L7-B-NOPF-CENT";
625  TString triggerNameJE = "CINT7EJ1-B-NOPF-CENTNOPMD";
626 
627  // Get the downscale factor for whichever MB trigger exists in the given run
628  std::vector<TString> runtriggers = downscaleOCDB->GetTriggerClasses();
629  Double_t downscalefactor;
630  for (auto i : runtriggers) {
631  if (i.EqualTo(triggerNameMB1) || i.EqualTo(triggerNameMB2)) {
632  downscalefactor = downscaleOCDB->GetDownscaleFactorForTriggerClass(i.Data());
633  break;
634  }
635  }
636 
637  // Store the inverse of the downscale factor, used later to weight the pT spectrum
638  fMBUpscaleFactor = 1/downscalefactor;
639 
640  TString histname = "Trigger/hMBDownscaleFactor";
641  fHistManager.FillTH1(histname.Data(), fMBUpscaleFactor);
642 
643  }
644 
645 }
646 
651 {
652  if (!fEventCuts.AcceptEvent(InputEvent()))
653  {
654  PostData(1, fOutput);
655  return kFALSE;
656  }
657  return kTRUE;
658 }
659 
668 {
669  TString histname;
670  AliJetContainer* jetCont = 0;
671  TIter next(&fJetCollArray);
672  while ((jetCont = static_cast<AliJetContainer*>(next()))) {
673  TString jetContName = jetCont->GetName();
674  if (jetContName.Contains("HardCore")) continue;
675 
676  //-----------------------------------------------------------------------------
677  // Find the leading di-jet candidate in each event, and if it satisfies the
678  // trig jet pT threshold, then fill di-jet candidate histogram (regardless of ass jet).
679  // The idea is to study the kinematic selections in post-processing.
680 
681  // Loop over leading hadron cut or not
682  for (Int_t leadingHadronCutType=0; leadingHadronCutType<2; leadingHadronCutType++) {
683 
684  // Find the dijet candidate of the event and store its info in struct fDijet
685  FindDijet(jetCont, leadingHadronCutType);
686 
687  // If we find a dijet candidate (i.e. acceptable trig jet; ass jet accepted or not), fill the di-jet candidate histogram
689  TString histname = TString::Format("%s/DijetCandObservables", jetCont->GetArrayName().Data());
690  FillDijetCandHistograms(histname);
691  }
692 
693  }
694 
695  //---------------------------------------------------------------------------------------------------
696  // Now, study the accepted dijet selection -- specified by the trig/ass jet pT conditions
697 
698  // Find the dijet candidate of the event and store its info in struct fDijet
699  FindDijet(jetCont, 0);
700 
701  // If we find an accepted dijet, fill the dijet imbalance histogram
703  TString histname = TString::Format("%s/DijetImbalanceObservables", jetCont->GetArrayName().Data());
705  }
706 
707  // If we find an accepted dijet, perform momentum-balance study (if requested)
709  histname = TString::Format("%s/MomentumBalance", jetCont->GetArrayName().Data());
710  DoMomentumBalance(histname);
711  }
712 
713  }
714 
715  //---------------------------------------------------------------------------
716  // Do the constituent threshold and geometrical matching study (if requested)
719 
720  return kTRUE;
721 }
722 
731 {
732  fDijet.clear();
733  fDijet.leadingHadronCutType = leadingHadronCutBin;
734 
735  // Get trigger jet
736  AliEmcalJet* trigJet = 0;
737  if (jetCont->GetRhoParameter())
738  trigJet = jetCont->GetLeadingJet("rho");
739  else
740  trigJet = jetCont->GetLeadingJet();
741  if(!trigJet) return;
742 
743  // Skip the event if the leading jet doesn't satisfy the pT threshold
744  Double_t trigJetPt = GetJetPt(jetCont, trigJet);
745  if ( trigJetPt < fMinTrigJetPt ) return;
746 
747  // Skip the event if the leading jet doesn't satisfy the leading hadron threshold
748  if (jetCont->GetLeadingHadronPt(trigJet) < fDijetLeadingHadronPt*leadingHadronCutBin) return;
749 
750  // Fill the dijet struct
751  fDijet.trigJet = trigJet;
752  fDijet.trigJetPt = trigJetPt;
753  fDijet.trigJetEta = trigJet->Eta();
754  fDijet.trigJetPhi = trigJet->Phi();
755 
756  // Find the subleading jet in the opposite hemisphere
757  AliEmcalJet *assJet = 0;
758  for(auto assJetCand : jetCont->accepted()) {
759  if (!assJetCand) continue;
760  Double_t assJetCandPt = GetJetPt(jetCont, assJetCand);
761  if ( TMath::Abs(trigJet->Phi() - assJetCand->Phi()) < fDeltaPhiMin ) continue;
762  if (assJet) {
763  Double_t assJetPt = GetJetPt(jetCont, assJet);
764  if ( assJetCandPt < assJetPt ) continue;
765  }
766  assJet = assJetCand;
767  }
768  if (!assJet) return;
769 
770  // Fill the dijet struct
771  fDijet.assJet = assJet;
772  fDijet.assJetPt = GetJetPt(jetCont, assJet);
773  fDijet.assJetPhi = assJet->Phi();
774  fDijet.assJetEta = assJet->Eta();
776 
777  fDijet.deltaPhi = TMath::Abs(trigJet->Phi() - assJet->Phi());
778  fDijet.deltaEta = trigJet->Eta() - assJet->Eta();
781  fDijet.kTy = TMath::Abs( fDijet.trigJetPt * TMath::Sin(fDijet.deltaPhi) );
782 }
783 
789 {
790 
791  AliTrackContainer* trackCont = dynamic_cast<AliTrackContainer*>(GetParticleContainer("tracks"));
792 
793  AliVTrack* track;
794  for (auto trackIterator : trackCont->accepted_momentum() ) {
795 
796  track = trackIterator.second;
797 
798  // Compute the delta phi between the track and its nearest jet (of the two jets in the dijet),
799  // as well as its pT-parallel projection onto the nearest jet's axis.
800 
801  Double_t trackPt = track->Pt();
802  Double_t trackPhi = track->Phi();
803  Double_t trigJetPhi = fDijet.trigJet->Phi();
804  Double_t assJetPhi = fDijet.assJet->Phi();
805 
806  Double_t deltaPhiTrigJet = TMath::Abs(trackPhi - trigJetPhi);
807  Double_t deltaPhiAssJet = TMath::Abs(trackPhi - assJetPhi);
808  Bool_t isNearside = deltaPhiTrigJet < deltaPhiAssJet;
809 
810  Double_t deltaPhi;
811  Double_t balancePt;
812  if (isNearside) {
813  deltaPhi = trackPhi - trigJetPhi;
814  balancePt = trackPt * TMath::Cos(deltaPhi);
815  }
816  else {
817  deltaPhi = trackPhi - assJetPhi;
818  balancePt = -trackPt * TMath::Cos(deltaPhi);
819  }
820 
821  FillMomentumBalanceHistograms(histname, deltaPhi, trackPt, balancePt);
822 
823  }
824 }
825 
830 {
831  // Get jet container with minimum constituent pT,E thresholds
832  TString jetContAllName = Form("Jet_AKTFullR0%d0_tracks_pT0150_caloClusters_E0300_pt_scheme", (int) (fMatchingJetR*10) );
833  AliJetContainer* jetContAll = GetJetContainer(jetContAllName.Data());
834 
835  // Get jet container with X GeV constituent pT,E thresholds
836  Int_t trackThreshold = (int) (fTrackConstituentThreshold*1000); // in MeV
837  Int_t clusThreshold = (int) (fClusterConstituentThreshold*1000); // in MeV
838  TString jetContHardCoreName = Form("JetHardCore_AKTFullR0%d0_tracks_pT%d_caloClusters_E%d_pt_scheme", (int) (fMatchingJetR*10), trackThreshold, clusThreshold);
839  AliJetContainer* jetContHardCore = GetJetContainer(jetContHardCoreName.Data());
840 
841  // Find the di-jet in the hard-core jet sample, then find the matching di-jet and fill histograms
842  FindDijet(jetContHardCore, 0);
843  if (fDijet.isAccepted) {
844  FindMatchingDijet(jetContAll);
846  }
847 }
848 
854 {
856 
857  // Loop over jets and find leading jet within R of fDijet.trigJet
858  AliEmcalJet *matchingTrigJet = 0;
859  for(auto matchingTrigJetCand : jetCont->accepted()) {
860  if (!matchingTrigJetCand) continue;
861  if ( GetDeltaR(matchingTrigJetCand, fDijet.trigJet) > fMatchingJetR ) continue;
862  if (matchingTrigJet) {
863  if ( GetJetPt(jetCont, matchingTrigJetCand) < GetJetPt(jetCont, matchingTrigJet) ) continue;
864  }
865  matchingTrigJet = matchingTrigJetCand;
866  }
867  if (!matchingTrigJet) return;
868 
869  // Loop over jets and find leading jet within R of fDijet.assJet
870  AliEmcalJet *matchingAssJet = 0;
871  for(auto matchingAssJetCand : jetCont->accepted()) {
872  if (!matchingAssJetCand) continue;
873  if ( GetDeltaR(matchingAssJetCand, fDijet.assJet) > fMatchingJetR ) continue;
874  if (matchingAssJet) {
875  if ( GetJetPt(jetCont, matchingAssJetCand) < GetJetPt(jetCont, matchingAssJet) ) continue;
876  }
877  matchingAssJet = matchingAssJetCand;
878  }
879 
880  // Determine which matching jet is the leading jet (i.e. allow them to flip)
881  if (matchingAssJet) {
882  AliEmcalJet* trigJet = matchingTrigJet;
883  AliEmcalJet* assJet = matchingAssJet;
884  if ( GetJetPt(jetCont, matchingTrigJet) < GetJetPt(jetCont, matchingAssJet) ) {
885  AliEmcalJet* trigJet = matchingAssJet;
886  AliEmcalJet* assJet = matchingTrigJet;
887  }
888 
889  // Fill the dijet struct
890  fMatchingDijet.trigJet = trigJet;
891  fMatchingDijet.trigJetPt = GetJetPt(jetCont, trigJet);
892  fMatchingDijet.trigJetEta = trigJet->Eta();
893  fMatchingDijet.trigJetPhi = trigJet->Phi();
894 
895  fMatchingDijet.assJet = assJet;
896  fMatchingDijet.assJetPt = GetJetPt(jetCont, assJet);
897  fMatchingDijet.assJetPhi = assJet->Phi();
898  fMatchingDijet.assJetEta = assJet->Eta();
900 
901  fMatchingDijet.deltaPhi = TMath::Abs(trigJet->Phi() - assJet->Phi());
902  fMatchingDijet.deltaEta = trigJet->Eta() - assJet->Eta();
905  fMatchingDijet.kTy = TMath::Abs( fMatchingDijet.trigJetPt * TMath::Sin(fMatchingDijet.deltaPhi) );
906  }
907 }
908 
913 {
914  // Loop over tracks and clusters in order to:
915  // (1) Compute scale factor for full jets
916  // (2) Compute delta-pT for full jets, with the random cone method
917  // For both the scale factor and delta-pT, we compute only one histogram each for EMCal.
918  // But for DCal, we bin in eta, in order to study DCal vs. PHOS vs. gap
919 
920  // Define the acceptance boundaries for the TPC and EMCal/DCal/PHOS
921  Double_t etaTPC = 0.9;
922  Double_t etaEMCal = 0.7;
923  Double_t etaMinDCal = 0.22;
924  Double_t etaMaxPHOS = 0.13;
925  Double_t phiMinEMCal = fGeom->GetArm1PhiMin() * TMath::DegToRad(); // 80
926  Double_t phiMaxEMCal = fGeom->GetEMCALPhiMax() * TMath::DegToRad(); // ~188
927  Double_t phiMinDCal = fGeom->GetDCALPhiMin() * TMath::DegToRad(); // 260
928  Double_t phiMaxDCal = fGeom->GetDCALPhiMax() * TMath::DegToRad(); // ~327 (1/3 SMs start at 320)
929  Double_t phiMinPHOS = 250 * TMath::DegToRad();
930  Double_t phiMaxPHOS = 320 * TMath::DegToRad();
931 
932  Double_t accTPC = 2 * etaTPC * 2 * TMath::Pi();
933  Double_t accEMCal = 2 * etaEMCal * (phiMaxEMCal - phiMinEMCal);
934  Double_t accDCalRegion = 2 * etaEMCal * (phiMaxDCal - phiMinDCal);
935 
936  // Define fiducial acceptances, to be used to generate random cones
937  TRandom3* r = new TRandom3(0);
938  Double_t jetR = jetCont->GetJetRadius();
939  Double_t etaEMCalfid = etaEMCal - jetR;
940  Double_t phiMinEMCalfid = phiMinEMCal + jetR;
941  Double_t phiMaxEMCalfid = phiMaxEMCal - jetR;
942  Double_t phiMinDCalRegionfid = phiMinDCal + jetR;
943  Double_t phiMaxDCalRegionfid = phiMaxDCal - jetR;
944 
945  // Generate EMCal random cone eta-phi
946  Double_t etaEMCalRC = r->Uniform(-etaEMCalfid, etaEMCalfid);
947  Double_t phiEMCalRC = r->Uniform(phiMinEMCalfid, phiMaxEMCalfid);
948 
949  // Generate DCalRegion random cone eta-phi inside each eta slice (same phi used for all)
950  Double_t etaStep = 0.05;
951  const Int_t nEtaBinsSF = 28; // 2 * 0.7 / 0.05
952  const Int_t nEtaBinsRC = 20; // 2 * 0.5 / 0.05
953 
954  Double_t phiDCalRC = r->Uniform(phiMinDCalRegionfid, phiMaxDCalRegionfid);
955  Double_t etaDCalRC[nEtaBinsRC];
956  Double_t etaMin;
957  Double_t etaMax;
958  for (Int_t bin=0; bin < nEtaBinsRC; bin++) {
959  etaMin = -etaEMCalfid + bin*etaStep;
960  etaMax = etaMin + etaStep;
961  etaDCalRC[bin] = r->Uniform(etaMin, etaMax);
962  }
963 
964  // Initialize the various sums to 0
965  Double_t trackPtSumTPC = 0;
966  Double_t trackPtSumEMCal = 0;
967  Double_t trackPtSumEMCalRC = 0;
968  Double_t clusESumEMCal = 0;
969  Double_t clusESumEMCalRC = 0;
970  Double_t trackPtSumDCal[nEtaBinsSF] = {0.};
971  Double_t trackPtSumDCalRC[nEtaBinsRC] = {0.};
972  Double_t clusESumDCal[nEtaBinsSF] = {0.};
973  Double_t clusESumDCalRC[nEtaBinsRC] = {0.};
974 
975  // Loop over tracks. Sum the track pT:
976  // (1) in the entire TPC, (2) in the EMCal, (3) in the EMCal random cone,
977  // (4) in the DCalRegion at each eta, (5) in the DCalRegion random cone at each eta
978  AliTrackContainer* trackCont = dynamic_cast<AliTrackContainer*>(GetParticleContainer("tracks"));
979  AliTLorentzVector track;
980  Double_t trackEta;
981  Double_t trackPhi;
982  Double_t trackPt;
983  Double_t deltaR;
984  for (auto trackIterator : trackCont->accepted_momentum() ) {
985 
986  track.Clear();
987  track = trackIterator.first;
988  trackEta = track.Eta();
989  trackPhi = track.Phi_0_2pi();
990  trackPt = track.Pt();
991 
992  // (1)
993  if (TMath::Abs(trackEta) < etaTPC) {
994  trackPtSumTPC += trackPt;
995  }
996 
997  // (2)
998  if (TMath::Abs(trackEta) < etaEMCal && trackPhi > phiMinEMCal && trackPhi < phiMaxEMCal) {
999  trackPtSumEMCal += trackPt;
1000  }
1001 
1002  // (3)
1003  deltaR = GetDeltaR(&track, etaEMCalRC, phiEMCalRC);
1004  if (deltaR < jetR) {
1005  trackPtSumEMCalRC += trackPt;
1006  }
1007 
1008  // (4)
1009  for (Int_t bin=0; bin < nEtaBinsSF; bin++) {
1010  if (trackPhi > phiMinDCal && trackPhi < phiMaxDCal) {
1011  etaMin = -etaEMCal+ bin*etaStep;
1012  etaMax = etaMin + etaStep;
1013  if (trackEta > etaMin && trackEta < etaMax) {
1014  trackPtSumDCal[bin] += trackPt;
1015  }
1016  }
1017  }
1018 
1019  // (5)
1020  for (Int_t bin=0; bin < nEtaBinsRC; bin++) {
1021  deltaR = GetDeltaR(&track, etaDCalRC[bin], phiDCalRC);
1022  if (deltaR < jetR) {
1023  trackPtSumDCalRC[bin] += trackPt;
1024  }
1025  }
1026 
1027  }
1028 
1029  // Loop over clusters. Sum the cluster ET:
1030  // (1) in the EMCal, (2) in the EMCal random cone, (3) in the DCalRegion at each eta,
1031  // (4) in the DCalRegion random cone at each eta
1033  AliTLorentzVector clus;
1034  Double_t clusEta;
1035  Double_t clusPhi;
1036  Double_t clusE;
1037  for (auto clusIterator : clusCont->accepted_momentum() ) {
1038 
1039  clus.Clear();
1040  clus = clusIterator.first;
1041  clusEta = clus.Eta();
1042  clusPhi = clus.Phi_0_2pi();
1043  clusE = clus.E();
1044 
1045  // (1)
1046  if (TMath::Abs(clusEta) < etaEMCal && clusPhi > phiMinEMCal && clusPhi < phiMaxEMCal) {
1047  clusESumEMCal += clusE;
1048  }
1049 
1050  // (2)
1051  deltaR = GetDeltaR(&clus, etaEMCalRC, phiEMCalRC);
1052  if (deltaR < jetR) {
1053  clusESumEMCalRC += clusE;
1054  }
1055 
1056  // (3)
1057  for (Int_t bin=0; bin < nEtaBinsSF; bin++) {
1058  if (clusPhi > phiMinDCal && clusPhi < phiMaxDCal) {
1059  etaMin = -etaEMCal+ bin*etaStep;
1060  etaMax = etaMin + etaStep;
1061  if (clusEta > etaMin && clusEta < etaMax) {
1062  clusESumDCal[bin] += clusE;
1063  }
1064  }
1065  }
1066 
1067  // (4)
1068  for (Int_t bin=0; bin < nEtaBinsRC; bin++) {
1069  deltaR = GetDeltaR(&clus, etaDCalRC[bin], phiDCalRC);
1070  if (deltaR < jetR) {
1071  clusESumDCalRC[bin] += clusE;
1072  }
1073  }
1074 
1075  }
1076 
1077  // Compute the scale factor for EMCal
1078  Double_t numerator = (trackPtSumEMCal + clusESumEMCal) / accEMCal;
1079  Double_t denominator = trackPtSumTPC / accTPC;
1080  Double_t scaleFactor = numerator / denominator;
1081  TString histname = TString::Format("%s/BackgroundHistograms/hScaleFactorEMCal", jetCont->GetArrayName().Data());
1082  fHistManager.FillTH2(histname, fCent, scaleFactor);
1083 
1084  // Compute the scale factor for DCalRegion
1085  Double_t accDCalRegionBin = accDCalRegion / nEtaBinsSF;
1086  for (Int_t bin=0; bin < nEtaBinsSF; bin++) {
1087  numerator = (trackPtSumDCal[bin] + clusESumDCal[bin]) / accDCalRegionBin;
1088  scaleFactor = numerator / denominator;
1089  histname = TString::Format("%s/BackgroundHistograms/hScaleFactorDCalRegion", jetCont->GetArrayName().Data());
1090  fHistManager.FillTH3(histname, fCent, scaleFactor, bin);
1091  }
1092 
1093  // Compute delta pT for EMCal
1094  Double_t rho = jetCont->GetRhoVal();
1095  Double_t deltaPt = trackPtSumEMCalRC + clusESumEMCalRC - rho * TMath::Pi() * jetR * jetR;
1096  histname = TString::Format("%s/BackgroundHistograms/hDeltaPtEMCal", jetCont->GetArrayName().Data());
1097  fHistManager.FillTH1(histname, deltaPt);
1098 
1099  // Compute delta pT for DCalRegion
1100  for (Int_t bin=0; bin < nEtaBinsRC; bin++) {
1101  deltaPt = trackPtSumDCalRC[bin] + clusESumDCalRC[bin] - rho * TMath::Pi() * jetR * jetR;
1102  histname = TString::Format("%s/BackgroundHistograms/hDeltaPtDCalRegion", jetCont->GetArrayName().Data());
1103  fHistManager.FillTH2(histname, deltaPt, bin);
1104  }
1105 
1106 }
1107 
1112 {
1113  Double_t pT = jet->Pt() - jetCont->GetRhoVal() * jet->Area();
1114  TString jetContName = jetCont->GetName();
1115  if (jetContName.Contains("HardCore")) pT = jet->Pt();
1116  return pT;
1117 }
1118 
1123 {
1124  Double_t deltaPhi = TMath::Abs(jet1->Phi() - jet2->Phi());
1125  Double_t deltaEta = TMath::Abs(jet1->Eta() - jet2->Eta());
1126  Double_t deltaR = TMath::Sqrt( deltaPhi*deltaPhi + deltaEta*deltaEta );
1127  return deltaR;
1128 }
1129 
1134 {
1135  Double_t deltaPhi = TMath::Abs(part->Phi_0_2pi() - phiRef);
1136  Double_t deltaEta = TMath::Abs(part->Eta() - etaRef);
1137  Double_t deltaR = TMath::Sqrt( deltaPhi*deltaPhi + deltaEta*deltaEta );
1138  return deltaR;
1139 }
1140 
1148 {
1150 
1151  return kTRUE;
1152 }
1153 
1159 {
1160  TString histname;
1161  AliJetContainer* jets = 0;
1162  TIter nextJetColl(&fJetCollArray);
1163  while ((jets = static_cast<AliJetContainer*>(nextJetColl()))) {
1164  TString jetContName = jets->GetName();
1165  if (jetContName.Contains("HardCore")) continue;
1166  Double_t rhoVal = 0;
1167  if (jets->GetRhoParameter()) {
1168  rhoVal = jets->GetRhoVal();
1169  histname = TString::Format("%s/JetHistograms/hRhoVsCent", jets->GetArrayName().Data());
1170  fHistManager.FillTH2(histname.Data(), fCent, rhoVal);
1171  }
1172 
1173  for (auto jet : jets->all()) {
1174 
1175  Float_t ptLeading = jets->GetLeadingHadronPt(jet);
1176  Float_t corrPt = jet->Pt() - rhoVal * jet->Area();
1177 
1178  // A vs. pT (fill before area cut)
1179  histname = TString::Format("%s/JetHistograms/hAreaVsPt", jets->GetArrayName().Data());
1180  fHistManager.FillTH2(histname.Data(), corrPt, jet->Area());
1181 
1182 
1183  // Rejection reason
1184  UInt_t rejectionReason = 0;
1185  if (!jets->AcceptJet(jet, rejectionReason)) {
1186  histname = TString::Format("%s/JetHistograms/hJetRejectionReason", jets->GetArrayName().Data());
1187  fHistManager.FillTH2(histname.Data(), jets->GetRejectionReasonBitPosition(rejectionReason), jet->Pt());
1188  continue;
1189  }
1190 
1191  // Centrality vs. pT
1192  histname = TString::Format("%s/JetHistograms/hCentVsPt", jets->GetArrayName().Data());
1193  fHistManager.FillTH2(histname.Data(), corrPt, fCent);
1194 
1195  // pT vs. eta vs. phi
1196  histname = TString::Format("%s/JetHistograms/hPtVsEtaVsPhi", jets->GetArrayName().Data());
1197  fHistManager.FillTH3(histname.Data(), jet->Eta(), jet->Phi_0_2pi(), corrPt);
1198 
1199  // pT un-downscaled
1200  histname = TString::Format("%s/JetHistograms/hPtUpscaledMB", jets->GetArrayName().Data());
1201  fHistManager.FillTH1(histname.Data(), corrPt, fMBUpscaleFactor);
1202 
1203  // pT-leading vs. pT
1204  histname = TString::Format("%s/JetHistograms/hPtLeadingVsPt", jets->GetArrayName().Data());
1205  fHistManager.FillTH2(histname.Data(), corrPt, ptLeading);
1206 
1207  // NEF vs. pT
1208  histname = TString::Format("%s/JetHistograms/hNEFVsPt", jets->GetArrayName().Data());
1209  fHistManager.FillTH2(histname.Data(), corrPt, jet->NEF());
1210 
1211  // z-leading (charged) vs. pT
1212  TLorentzVector leadPart;
1213  jets->GetLeadingHadronMomentum(leadPart, jet);
1214  Double_t z = GetParallelFraction(leadPart.Vect(), jet);
1215  if (z == 1 || (z > 1 && z - 1 < 1e-3)) z = 0.999; // so that it will contribute to the bin 0.9-1 rather than 1-1.1
1216  histname = TString::Format("%s/JetHistograms/hZLeadingVsPt", jets->GetArrayName().Data());
1217  fHistManager.FillTH2(histname.Data(), corrPt, z);
1218 
1219  // z (charged) vs. pT
1220  histname = TString::Format("%s/JetHistograms/hZVsPt", jets->GetArrayName().Data());
1221  AliVTrack* track;
1222  for (Int_t i=0; i<jet->GetNumberOfTracks(); i++) {
1223  track = static_cast<AliVTrack*>(jet->Track(i));
1224  z = track->Pt() / TMath::Abs(corrPt);
1225  fHistManager.FillTH2(histname.Data(), corrPt, z);
1226  }
1227 
1228  // Nconst vs. pT
1229  histname = TString::Format("%s/JetHistograms/hNConstVsPt", jets->GetArrayName().Data());
1230  fHistManager.FillTH2(histname.Data(), corrPt, jet->GetNumberOfConstituents());
1231 
1232  } //jet loop
1233 
1234  //---------------------------------------------------------------------------
1235  // Do study of background (if requested)
1237  }
1238 }
1239 
1244 {
1245  Double_t contents[30]={0};
1246  THnSparse* histJetObservables = static_cast<THnSparse*>(fHistManager.FindObject(histname));
1247  if (!histJetObservables) return;
1248  for (Int_t n = 0; n < histJetObservables->GetNdimensions(); n++) {
1249  TString title(histJetObservables->GetAxis(n)->GetTitle());
1250  if (title=="Centrality (%)")
1251  contents[n] = fCent;
1252  else if (title=="LeadingHadronRequired")
1253  contents[n] = fDijet.leadingHadronCutType;
1254  else if (title=="#it{p}_{T,trig jet} (GeV/#it{c})")
1255  contents[n] = fDijet.trigJetPt;
1256  else if (title=="#it{p}_{T,ass jet} (GeV/#it{c})")
1257  contents[n] = fDijet.assJetPt;
1258  else if (title=="#phi_{trig jet}")
1259  contents[n] = fDijet.trigJetPhi;
1260  else if (title=="#phi_{ass jet}")
1261  contents[n] = fDijet.assJetPhi;
1262  else if (title=="#eta_{trig jet}")
1263  contents[n] = fDijet.trigJetEta;
1264  else if (title=="#eta_{ass jet}")
1265  contents[n] = fDijet.assJetEta;
1266  else
1267  AliWarning(Form("Unable to fill dimension %s!",title.Data()));
1268  }
1269  histJetObservables->Fill(contents);
1270 }
1271 
1276 {
1277  Double_t contents[30]={0};
1278  THnSparse* histJetObservables = static_cast<THnSparse*>(fHistManager.FindObject(histname));
1279  if (!histJetObservables) return;
1280  for (Int_t n = 0; n < histJetObservables->GetNdimensions(); n++) {
1281  TString title(histJetObservables->GetAxis(n)->GetTitle());
1282  if (title=="Centrality (%)")
1283  contents[n] = fCent;
1284  else if (title=="#Delta#phi")
1285  contents[n] = fDijet.deltaPhi;
1286  else if (title=="#Delta#eta")
1287  contents[n] = fDijet.deltaEta;
1288  else if (title=="A_{J}")
1289  contents[n] = fDijet.AJ;
1290  else if (title=="x_{J}")
1291  contents[n] = fDijet.xJ;
1292  else if (title=="k_{Ty} (GeV)")
1293  contents[n] = fDijet.kTy;
1294  else if (title=="N_{tracks, trig jet}")
1295  contents[n] = fDijet.trigJet->GetNumberOfTracks();
1296  else if (title=="N_{tracks, ass jet}")
1297  contents[n] = fDijet.assJet->GetNumberOfTracks();
1298  else
1299  AliWarning(Form("Unable to fill dimension %s!",title.Data()));
1300  }
1301  histJetObservables->Fill(contents);
1302 }
1303 
1308 {
1309  Double_t contents[30]={0};
1310  THnSparse* histJetObservables = static_cast<THnSparse*>(fHistManager.FindObject(histname));
1311  if (!histJetObservables) return;
1312  for (Int_t n = 0; n < histJetObservables->GetNdimensions(); n++) {
1313  TString title(histJetObservables->GetAxis(n)->GetTitle());
1314  if (title=="A_{J}")
1315  contents[n] = fDijet.AJ;
1316  else if (title=="#Delta#phi")
1317  contents[n] = deltaPhi;
1318  else if (title=="#it{p}_{T,particle} (GeV/#it{c})")
1319  contents[n] = trackPt;
1320  else if (title=="#it{p}_{T}#parallel (GeV/#it{c})")
1321  contents[n] = balancePt;
1322  else
1323  AliWarning(Form("Unable to fill dimension %s!",title.Data()));
1324  }
1325  histJetObservables->Fill(contents);
1326 }
1327 
1332 {
1333  // Matching efficiency histogram
1334  TString histname = "GeometricalMatchingEfficiency";
1335 
1336  // If we have a matching di-jet, fill the geometrical matching histogram
1337  if (fMatchingDijet.assJet) {
1338  fHistManager.FillTH1(histname.Data(), 1);
1339 
1342  Bool_t isSwitched = trigDeltaR > fMatchingJetR;
1343 
1344  TString thnname = "GeometricalMatching";
1345  Double_t contents[30]={0};
1346  THnSparse* histJetObservables = static_cast<THnSparse*>(fHistManager.FindObject(thnname.Data()));
1347  if (!histJetObservables) return;
1348  for (Int_t n = 0; n < histJetObservables->GetNdimensions(); n++) {
1349  TString title(histJetObservables->GetAxis(n)->GetTitle());
1350  if (title=="Centrality (%)")
1351  contents[n] = fCent;
1352  else if (title=="isSwitched")
1353  contents[n] = isSwitched;
1354  else if (title=="#DeltaR_{trig}")
1355  contents[n] = trigDeltaR;
1356  else if (title=="#DeltaR_{ass}")
1357  contents[n] = assDeltaR;
1358  else if (title=="trig #it{p}_{T,low-thresh} - #it{p}_{T,hard-core}")
1359  contents[n] = fMatchingDijet.trigJetPt - fDijet.trigJetPt;
1360  else if (title=="ass #it{p}_{T,low-thresh} - #it{p}_{T,hard-core}")
1361  contents[n] = fMatchingDijet.assJetPt - fDijet.assJetPt;
1362  else if (title=="A_{J} low-threshold")
1363  contents[n] = fMatchingDijet.AJ;
1364  else if (title=="A_{J} hard-core")
1365  contents[n] = fDijet.AJ;
1366  else
1367  AliWarning(Form("Unable to fill dimension %s!",title.Data()));
1368  }
1369  histJetObservables->Fill(contents);
1370  }
1371  else {
1372  fHistManager.FillTH1(histname.Data(), 0.);
1373  }
1374 }
1375 
Double_t Area() const
Definition: AliEmcalJet.h:123
Double_t GetRhoVal() const
const TString & GetRhoName() const
double Double_t
Definition: External.C:58
const char * title
Definition: MakeQAPdf.C:26
AliJetContainer * GetJetContainer(Int_t i=0) const
Bool_t fUseManualEventCuts
Flag to use manual event cuts.
UInt_t fOffTrigger
offline trigger for event selection
Double_t Eta() const
Definition: AliEmcalJet.h:114
static AliEmcalDownscaleFactorsOCDB * Instance()
const AliClusterIterableMomentumContainer accepted_momentum() const
Double_t Phi() const
Definition: AliEmcalJet.h:110
Container with name, TClonesArray and cuts for particles.
Bool_t fPlotJetHistograms
Set whether to enable inclusive jet histograms.
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.
TList * fEventCutList
! Output list for event cut histograms
void FillTH3(const char *hname, double x, double y, double z, double weight=1., Option_t *opt="")
Fill a 3D histogram within the container.
Double_t fDijetLeadingHadronPt
leading hadron pT threshold for leading jet in dijet
virtual Bool_t AcceptJet(Int_t i, UInt_t &rejectionReason) const
Double_t GetDownscaleFactorForTriggerClass(const TString &trigger) const
Double_t fMinAssJetPt
subleading jet min pT in a dijet pair, for it to be accepted
UShort_t GetNumberOfTracks() const
Definition: AliEmcalJet.h:132
Bool_t fDoMomentumBalance
Set whether to enable momentum balance study.
void FindDijet(AliJetContainer *jetCont, Int_t leadingHadronCutBin)
Bool_t fPlotDijetImbalanceHistograms
Set whether to enable dijet imbalance histograms.
AliEventCuts fEventCuts
event selection utility
AliParticleContainer * GetParticleContainer(Int_t i=0) const
Get particle container attached to this task.
const Int_t nPtBins
AliEmcalJet * GetLeadingJet(const char *opt="")
void GetLeadingHadronMomentum(TLorentzVector &mom, const AliEmcalJet *jet) const
Dijet_t fDijet
! dijet candidate (per event)
TH2 * CreateTH2(const char *name, const char *title, int nbinsx, double xmin, double xmax, int nbinsy, double ymin, double ymax, Option_t *opt="")
Create a new TH2 within the container.
TObject * FindObject(const char *name) const
Find an object inside the container.
int Int_t
Definition: External.C:63
unsigned int UInt_t
Definition: External.C:33
THashList * GetListOfHistograms() const
Get the list of histograms.
Definition: THistManager.h:658
Double_t GetDeltaR(AliEmcalJet *jet1, AliEmcalJet *jet2)
Double_t GetLeadingHadronPt(const AliEmcalJet *jet) const
float Float_t
Definition: External.C:68
Double_t fClusterConstituentThreshold
constituent threshold for matching study
Di-jet imbalance analysis with full jets.
AliEMCALGeometry * fGeom
!emcal geometry
Double_t Phi_0_2pi() const
TH1 * CreateTH1(const char *name, const char *title, int nbins, double xmin, double xmax, Option_t *opt="")
Create a new TH1 within the container.
BeamType fForceBeamType
forced beam type
AliClusterContainer * GetClusterContainer(Int_t i=0) const
Get cluster container attached to this task.
Double_t fDeltaPhiMin
minimum delta phi between di-jets
Double_t fMBUpscaleFactor
! inverse of downscale factor, for MB trigger
Bool_t fComputeBackground
Set whether to enable study of background.
std::vector< TString > GetTriggerClasses() const
Double_t fCent
!event centrality
Double_t GetJetPt(AliJetContainer *jetCont, AliEmcalJet *jet)
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 FillMomentumBalanceHistograms(TString histname, Double_t deltaPhi, Double_t trackPt, Double_t balancePt)
void FillTH1(const char *hname, double x, double weight=1., Option_t *opt="")
Fill a 1D histogram within the container.
TObjArray fJetCollArray
jet collection array
AliRhoParameter * GetRhoParameter()
Double_t Pt() const
Definition: AliEmcalJet.h:102
static Double_t * GenerateFixedBinArray(Int_t n, Double_t min, Double_t max)
Float_t GetJetRadius() const
AliEmcalList * fOutput
!output list
Definition: External.C:220
ClassImp(AliAnalysisTaskCRC) AliAnalysisTaskCRC
Handler for downscale factors for various triggers obtained from the OCDB.
Base task in the EMCAL jet framework.
Bool_t fDoGeometricalMatching
Set whether to enable constituent study with geometrical matching.
Represent a jet reconstructed using the EMCal jet framework.
Definition: AliEmcalJet.h:44
const AliTrackIterableMomentumContainer accepted_momentum() const
Int_t GetRunNumber(TString)
Definition: PlotMuonQA.C:2235
Dijet_t fMatchingDijet
! low-threshold matching dijet, for matching study
void SetRejectionReasonLabels(TAxis *axis)
void UserCreateOutputObjects()
Main initialization function on the worker.
const Int_t nbins
const AliJetIterableContainer accepted() const
bool Bool_t
Definition: External.C:53
Double_t fTrackConstituentThreshold
constituent threshold for matching study
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 fMinTrigJetPt
leading jet min pT in a dijet pair
Container structure for EMCAL clusters.
Bool_t fNeedEmcalGeom
whether or not the task needs the emcal geometry
Container for jet within the EMCAL jet framework.
Definition: External.C:196
Bool_t fPlotDijetCandHistograms
Set whether to enable dijet pair histograms.
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