AliPhysics  4e47bdd (4e47bdd)
 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(250),
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(250),
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  Int_t minPtBin = -fMaxPt/2 + 25;
178  Int_t maxPtBin = fMaxPt/2 + 25;
179 
180  AliJetContainer* jets = 0;
181  TIter nextJetColl(&fJetCollArray);
182  while ((jets = static_cast<AliJetContainer*>(nextJetColl()))) {
183 
184  // Jet rejection reason
185  histname = TString::Format("%s/JetHistograms/hJetRejectionReason", jets->GetArrayName().Data());
186  title = histname + ";Rejection reason;#it{p}_{T,jet} (GeV/#it{c});counts";
187  TH2* hist = fHistManager.CreateTH2(histname.Data(), title.Data(), 32, 0, 32, 50, 0, 250);
188  SetRejectionReasonLabels(hist->GetXaxis());
189 
190  // Rho vs. Centrality
191  if (!jets->GetRhoName().IsNull()) {
192  histname = TString::Format("%s/JetHistograms/hRhoVsCent", jets->GetArrayName().Data());
193  title = histname + ";Centrality (%);#rho (GeV/#it{c});counts";
194  fHistManager.CreateTH2(histname.Data(), title.Data(), 101, 0, 101, 100, 0, 500);
195  }
196 
197  // Centrality vs. pT
198  histname = TString::Format("%s/JetHistograms/hCentVsPt", jets->GetArrayName().Data());
199  title = histname + ";#it{p}_{T}^{corr} (GeV/#it{c});Centrality (%);counts";
200  fHistManager.CreateTH2(histname.Data(), title.Data(), nPtBins, minPtBin, maxPtBin, 10, 0, 100);
201 
202  // pT vs. eta vs. phi
203  histname = TString::Format("%s/JetHistograms/hPtVsEtaVsPhi", jets->GetArrayName().Data());
204  title = histname + ";#eta_{jet} (rad);#phi_{jet} (rad);#it{p}_{T}^{corr} (GeV/#it{c})";
205  fHistManager.CreateTH3(histname.Data(), title.Data(), 50, -0.5, 0.5, 101, 0, TMath::Pi() * 2.02, 75, 0, maxPtBin);
206 
207  // pT upscaled
208  histname = TString::Format("%s/JetHistograms/hPtUpscaledMB", jets->GetArrayName().Data());
209  title = histname + ";#it{p}_{T}^{corr} (GeV/#it{c})";
210  fHistManager.CreateTH1(histname.Data(), title.Data(), 75, 0, maxPtBin, "s");
211 
212  // pT-leading vs. pT
213  histname = TString::Format("%s/JetHistograms/hPtLeadingVsPt", jets->GetArrayName().Data());
214  title = histname + ";#it{p}_{T}^{corr} (GeV/#it{c});#it{p}_{T,particle}^{leading} (GeV/#it{c})";
215  fHistManager.CreateTH2(histname.Data(), title.Data(), nPtBins, minPtBin, maxPtBin, 150, 0, 150);
216 
217  // A vs. pT
218  histname = TString::Format("%s/JetHistograms/hAreaVsPt", jets->GetArrayName().Data());
219  title = histname + ";#it{p}_{T}^{corr} (GeV/#it{c});#it{A}_{jet}";
220  fHistManager.CreateTH2(histname.Data(), title.Data(), nPtBins, minPtBin, maxPtBin, fMaxPt/3, 0, 1.5);
221 
222  // NEF vs. pT
223  histname = TString::Format("%s/JetHistograms/hNEFVsPt", jets->GetArrayName().Data());
224  title = histname + ";#it{p}_{T}^{corr} (GeV/#it{c});NEF";
225  fHistManager.CreateTH2(histname.Data(), title.Data(), nPtBins, minPtBin, maxPtBin, fMaxPt/5, 0, 1.0);
226 
227  // z-leading vs. pT
228  histname = TString::Format("%s/JetHistograms/hZLeadingVsPt", jets->GetArrayName().Data());
229  title = histname + ";#it{p}_{T}^{corr} (GeV/#it{c});#it{z}_{leading}";
230  fHistManager.CreateTH2(histname.Data(), title.Data(), nPtBins, minPtBin, maxPtBin, fMaxPt/5, 0, 1.0);
231 
232  // Nconst vs. pT
233  histname = TString::Format("%s/JetHistograms/hNConstVsPt", jets->GetArrayName().Data());
234  title = histname + ";#it{p}_{T}^{corr} (GeV/#it{c});No. of constituents";
235  fHistManager.CreateTH2(histname.Data(), title.Data(), nPtBins, minPtBin, maxPtBin, fMaxPt/5, 0, fMaxPt);
236 
237  // Allocate background subtraction histograms, if enabled
238  if (fComputeBackground) {
239 
240  histname = TString::Format("%s/BackgroundHistograms/hScaleFactorEMCal", jets->GetArrayName().Data());
241  title = histname + ";Centrality;Scale factor;counts";
242  fHistManager.CreateTH2(histname.Data(), title.Data(), 50, 0, 100, 100, 0, 5);
243 
244  histname = TString::Format("%s/BackgroundHistograms/hDeltaPtEMCal", jets->GetArrayName().Data());
245  title = histname + ";#delta#it{p}_{T} (GeV/#it{c});counts";
246  fHistManager.CreateTH1(histname.Data(), title.Data(), 400, -100, 100);
247 
248  histname = TString::Format("%s/BackgroundHistograms/hScaleFactorDCalRegion", jets->GetArrayName().Data());
249  title = histname + ";Centrality;Scale factor;Eta bin";
250  fHistManager.CreateTH3(histname.Data(), title.Data(), 50, 0, 100, 200, 0, 10, 28, -0.5, 27.5);
251 
252  histname = TString::Format("%s/BackgroundHistograms/hDeltaPtDCalRegion", jets->GetArrayName().Data());
253  title = histname + ";#delta#it{p}_{T} (GeV/#it{c});Eta bin;counts";
254  fHistManager.CreateTH2(histname.Data(), title.Data(), 400, -100, 100, 20, -0.5, 19.5);
255 
256  }
257 
258  }
259 
260  // MB downscale factor histogram
261  histname = "Trigger/hMBDownscaleFactor";
262  title = histname + ";Downscale factor;counts";
263  TH1* hist = fHistManager.CreateTH1(histname.Data(), title.Data(), 200, 0, 200);
264 
265 }
266 
267 /*
268  * This function allocates the histograms for dijet candidates, i.e. dijet pairs with the leading jet
269  * passing the trigger condition, but no condition on the subleading jet. In particular this histogram is
270  * designed to study the kinematic selection of dijets.
271  */
273 {
274  // Allocate dijet THnSparse
275  AliJetContainer* jets = 0;
276  TIter nextJetColl(&fJetCollArray);
277  while ((jets = static_cast<AliJetContainer*>(nextJetColl()))) {
278  TString axisTitle[30]= {""};
279  Int_t nbins[30] = {0};
280  Double_t min[30] = {0.};
281  Double_t max[30] = {0.};
282  Double_t *binEdges[20] = {0};
283  Int_t dim = 0;
284 
285  if (fForceBeamType != kpp) {
286  axisTitle[dim] = "Centrality (%)";
287  nbins[dim] = fNCentHistBins;
288  binEdges[dim] = fCentHistBins;
289  min[dim] = fCentHistBins[0];
290  max[dim] = fCentHistBins[fNCentHistBins];
291  dim++;
292  }
293 
294  axisTitle[dim] = "LeadingHadronRequired";
295  nbins[dim] = 2;
296  min[dim] = -0.5;
297  max[dim] = 1.5;
298  binEdges[dim] = GenerateFixedBinArray(nbins[dim], min[dim], max[dim]);
299  dim++;
300 
301  axisTitle[dim] = "#it{p}_{T,trig jet} (GeV/#it{c})";
302  nbins[dim] = fMaxPt/3;
303  min[dim] = -fMaxPt/2 + 25;
304  max[dim] = fMaxPt/2 + 25;
305  binEdges[dim] = GenerateFixedBinArray(nbins[dim], min[dim], max[dim]);
306  dim++;
307 
308  axisTitle[dim] = "#it{p}_{T,ass jet} (GeV/#it{c})";
309  nbins[dim] = fMaxPt/3;
310  min[dim] = -fMaxPt/2 + 25;
311  max[dim] = fMaxPt/2 + 25;
312  binEdges[dim] = GenerateFixedBinArray(nbins[dim], min[dim], max[dim]);
313  dim++;
314 
315  axisTitle[dim] = "#phi_{trig jet}";
316  nbins[dim] = fMaxPt/3;
317  min[dim] = 0;
318  max[dim] = TMath::TwoPi();
319  binEdges[dim] = GenerateFixedBinArray(nbins[dim], min[dim], max[dim]);
320  dim++;
321 
322  axisTitle[dim] = "#phi_{ass jet}";
323  nbins[dim] = fMaxPt/3;
324  min[dim] = 0;
325  max[dim] = TMath::TwoPi();
326  binEdges[dim] = GenerateFixedBinArray(nbins[dim], min[dim], max[dim]);
327  dim++;
328 
329  axisTitle[dim] = "#eta_{trig jet}";
330  nbins[dim] = fMaxPt/3;
331  min[dim] = -1;
332  max[dim] = 1;
333  binEdges[dim] = GenerateFixedBinArray(nbins[dim], min[dim], max[dim]);
334  dim++;
335 
336  axisTitle[dim] = "#eta_{ass jet}";
337  nbins[dim] = fMaxPt/3;
338  min[dim] = -1;
339  max[dim] = 1;
340  binEdges[dim] = GenerateFixedBinArray(nbins[dim], min[dim], max[dim]);
341  dim++;
342 
343  TString thnname = TString::Format("%s/DijetCandObservables", jets->GetArrayName().Data());
344  THnSparse* hn = fHistManager.CreateTHnSparse(thnname.Data(), thnname.Data(), dim, nbins, min, max);
345  for (Int_t i = 0; i < dim; i++) {
346  hn->GetAxis(i)->SetTitle(axisTitle[i]);
347  hn->SetBinEdges(i, binEdges[i]);
348  }
349  }
350 }
351 
352 /*
353  * This function allocates the histograms for accepted dijets.
354  * The purpose is to study in detail the imbalance properties of the accepted dijets.
355  */
357 {
358  // Allocate dijet imbalance THnSparse
359  AliJetContainer* jets = 0;
360  TIter nextJetColl(&fJetCollArray);
361  while ((jets = static_cast<AliJetContainer*>(nextJetColl()))) {
362  TString axisTitle[30]= {""};
363  Int_t nbins[30] = {0};
364  Double_t min[30] = {0.};
365  Double_t max[30] = {0.};
366  Double_t *binEdges[20] = {0};
367  Int_t dim = 0;
368 
369  if (fForceBeamType != kpp) {
370  axisTitle[dim] = "Centrality (%)";
371  nbins[dim] = fNCentHistBins;
372  binEdges[dim] = fCentHistBins;
373  min[dim] = fCentHistBins[0];
374  max[dim] = fCentHistBins[fNCentHistBins];
375  dim++;
376  }
377 
378  axisTitle[dim] = "#Delta#phi";
379  nbins[dim] = 100;
380  min[dim] = 0;
381  max[dim] = 4;
382  binEdges[dim] = GenerateFixedBinArray(nbins[dim], min[dim], max[dim]);
383  dim++;
384 
385  axisTitle[dim] = "#Delta#eta";
386  nbins[dim] = 100;
387  min[dim] = -2;
388  max[dim] = 2;
389  binEdges[dim] = GenerateFixedBinArray(nbins[dim], min[dim], max[dim]);
390  dim++;
391 
392  axisTitle[dim] = "A_{J}";
393  nbins[dim] = 100;
394  min[dim] = 0;
395  max[dim] = 1;
396  binEdges[dim] = GenerateFixedBinArray(nbins[dim], min[dim], max[dim]);
397  dim++;
398 
399  axisTitle[dim] = "x_{J}";
400  nbins[dim] = 100;
401  min[dim] = 0;
402  max[dim] = 1;
403  binEdges[dim] = GenerateFixedBinArray(nbins[dim], min[dim], max[dim]);
404  dim++;
405 
406  axisTitle[dim] = "k_{Ty} (GeV)";
407  nbins[dim] = 100;
408  min[dim] = 0;
409  max[dim] = 100;
410  binEdges[dim] = GenerateFixedBinArray(nbins[dim], min[dim], max[dim]);
411  dim++;
412 
413  axisTitle[dim] = "N_{tracks, trig jet}";
414  nbins[dim] = fMaxPt/5;
415  min[dim] = 0;
416  max[dim] = 100;
417  binEdges[dim] = GenerateFixedBinArray(nbins[dim], min[dim], max[dim]);
418  dim++;
419 
420  axisTitle[dim] = "N_{tracks, ass jet}";
421  nbins[dim] = fMaxPt/5;
422  min[dim] = 0;
423  max[dim] = 100;
424  binEdges[dim] = GenerateFixedBinArray(nbins[dim], min[dim], max[dim]);
425  dim++;
426 
427  TString thnname = TString::Format("%s/DijetImbalanceObservables", jets->GetArrayName().Data());
428  THnSparse* hn = fHistManager.CreateTHnSparse(thnname.Data(), thnname.Data(), dim, nbins, min, max);
429  for (Int_t i = 0; i < dim; i++) {
430  hn->GetAxis(i)->SetTitle(axisTitle[i]);
431  hn->SetBinEdges(i, binEdges[i]);
432  }
433  }
434 }
435 
436 /*
437  * This function allocates the histograms for the momentum balance study.
438  */
440 {
441  AliJetContainer* jets = 0;
442  TIter nextJetColl(&fJetCollArray);
443  while ((jets = static_cast<AliJetContainer*>(nextJetColl()))) {
444 
445  // Allocate THnSparse
446  TString axisTitle[30]= {""};
447  Int_t nbins[30] = {0};
448  Double_t min[30] = {0.};
449  Double_t max[30] = {0.};
450  Double_t *binEdges[20] = {0};
451  Int_t dim = 0;
452 
453  axisTitle[dim] = "A_{J}";
454  nbins[dim] = 100;
455  min[dim] = 0;
456  max[dim] = 1;
457  binEdges[dim] = GenerateFixedBinArray(nbins[dim], min[dim], max[dim]);
458  dim++;
459 
460  axisTitle[dim] = "#Delta#phi";
461  nbins[dim] = 100;
462  min[dim] = -4;
463  max[dim] = 4;
464  binEdges[dim] = GenerateFixedBinArray(nbins[dim], min[dim], max[dim]);
465  dim++;
466 
467  axisTitle[dim] = "#it{p}_{T,particle} (GeV/#it{c})";
468  nbins[dim] = 9;
469  Double_t* pTParticleBins = new Double_t[nbins[dim]+1];
470  GenerateFixedBinArray(1, 0.15, 0.3, pTParticleBins);
471  GenerateFixedBinArray(1, 0.3, 0.5, pTParticleBins+1);
472  GenerateFixedBinArray(1, 0.5, 1, pTParticleBins+2);
473  GenerateFixedBinArray(2, 1, 5, pTParticleBins+3);
474  GenerateFixedBinArray(3, 5, 20, pTParticleBins+5);
475  GenerateFixedBinArray(1, 20, 150, pTParticleBins+8);
476  min[dim] = 0;
477  max[dim] = pTParticleBins[nbins[dim]];
478  binEdges[dim] = pTParticleBins;
479  dim++;
480 
481  axisTitle[dim] = "#it{p}_{T}#parallel (GeV/#it{c})";
482  nbins[dim] = 80;
483  Double_t* pTParallelBins = new Double_t[nbins[dim]+1];
484  GenerateFixedBinArray(20, 0, 2, pTParallelBins);
485  GenerateFixedBinArray(16, 2, 10, pTParallelBins+20);
486  GenerateFixedBinArray(10, 10, 20, pTParallelBins+36);
487  GenerateFixedBinArray(10, 20, 40, pTParallelBins+46);
488  GenerateFixedBinArray(24, 40, 150, pTParallelBins+56);
489  min[dim] = 0;
490  max[dim] = pTParallelBins[nbins[dim]];
491  binEdges[dim] = pTParallelBins;
492  dim++;
493 
494  TString thnname = TString::Format("%s/MomentumBalance", jets->GetArrayName().Data());
495  THnSparse* hn = fHistManager.CreateTHnSparse(thnname.Data(), thnname.Data(), dim, nbins, min, max);
496  for (Int_t i = 0; i < dim; i++) {
497  hn->GetAxis(i)->SetTitle(axisTitle[i]);
498  hn->SetBinEdges(i, binEdges[i]);
499  }
500  }
501 }
502 
503 /*
504  * This function allocates the histograms for the constituent dijet study.
505  */
507 {
508  // Allocate geometrical matching THnSparse
509  TString axisTitle[30]= {""};
510  Int_t nbins[30] = {0};
511  Double_t min[30] = {0.};
512  Double_t max[30] = {0.};
513  Double_t *binEdges[20] = {0};
514  Int_t dim = 0;
515 
516  if (fForceBeamType != kpp) {
517  axisTitle[dim] = "Centrality (%)";
518  nbins[dim] = fNCentHistBins;
519  binEdges[dim] = fCentHistBins;
520  min[dim] = fCentHistBins[0];
521  max[dim] = fCentHistBins[fNCentHistBins];
522  dim++;
523  }
524 
525  axisTitle[dim] = "isSwitched";
526  nbins[dim] = 2;
527  min[dim] = -0.5;
528  max[dim] = 1.5;
529  binEdges[dim] = GenerateFixedBinArray(nbins[dim], min[dim], max[dim]);
530  dim++;
531 
532  axisTitle[dim] = "#DeltaR_{trig}";
533  nbins[dim] = 50;
534  min[dim] = 0;
535  max[dim] = 0.5;
536  binEdges[dim] = GenerateFixedBinArray(nbins[dim], min[dim], max[dim]);
537  dim++;
538 
539  axisTitle[dim] = "#DeltaR_{ass}";
540  nbins[dim] = 50;
541  min[dim] = 0;
542  max[dim] = 0.5;
543  binEdges[dim] = GenerateFixedBinArray(nbins[dim], min[dim], max[dim]);
544  dim++;
545 
546  axisTitle[dim] = "trig #it{p}_{T,low-thresh} - #it{p}_{T,hard-core}";
547  nbins[dim] = 100;
548  min[dim] = -50;
549  max[dim] = 50;
550  binEdges[dim] = GenerateFixedBinArray(nbins[dim], min[dim], max[dim]);
551  dim++;
552 
553  axisTitle[dim] = "ass #it{p}_{T,low-thresh} - #it{p}_{T,hard-core}";
554  nbins[dim] = 100;
555  min[dim] = -50;
556  max[dim] = 50;
557  binEdges[dim] = GenerateFixedBinArray(nbins[dim], min[dim], max[dim]);
558  dim++;
559 
560  axisTitle[dim] = "A_{J} low-threshold";
561  nbins[dim] = 100;
562  min[dim] = 0;
563  max[dim] = 1;
564  binEdges[dim] = GenerateFixedBinArray(nbins[dim], min[dim], max[dim]);
565  dim++;
566 
567  axisTitle[dim] = "A_{J} hard-core";
568  nbins[dim] = 100;
569  min[dim] = 0;
570  max[dim] = 1;
571  binEdges[dim] = GenerateFixedBinArray(nbins[dim], min[dim], max[dim]);
572  dim++;
573 
574  TString thnname = "GeometricalMatching";
575  THnSparse* hn = fHistManager.CreateTHnSparse(thnname.Data(), thnname.Data(), dim, nbins, min, max);
576  for (Int_t i = 0; i < dim; i++) {
577  hn->GetAxis(i)->SetTitle(axisTitle[i]);
578  hn->SetBinEdges(i, binEdges[i]);
579  }
580 
581  // Allocate other histograms
582  TString histname;
583  TString title;
584  histname = "GeometricalMatchingEfficiency";
585  title = histname + ";isMatched;counts";
586  TH1* hist = fHistManager.CreateTH1(histname.Data(), title.Data(), 2, -0.5, 1.5);
587 }
588 
594 {
596 
597  fNeedEmcalGeom = kTRUE;
598 
599  AliInfo(Form("Trigger jet threshold = %f, Associated jet threshold = %f", fMinTrigJetPt, fMinAssJetPt));
600  AliInfo(Form("Leading hadron threshold (for dijet leading jet): %f GeV", fDijetLeadingHadronPt));
601  AliInfo(Form("Momentum balance study: %d", fDoMomentumBalance));
602  AliInfo(Form("Geometrical matching study: %d", fDoGeometricalMatching));
603 
604 }
605 
610 
611  // Get the downscaling factors for MB triggers (to be used to calculate trigger efficiency)
612 
613  // Get instance of the downscale factor helper class
615  downscaleOCDB->SetRun(InputEvent()->GetRunNumber());
616 
617  // There are two possible min bias triggers for LHC15o
618  TString triggerNameMB1 = "CINT7-B-NOPF-CENT";
619  TString triggerNameMB2 = "CV0L7-B-NOPF-CENT";
620  TString triggerNameJE = "CINT7EJ1-B-NOPF-CENTNOPMD";
621 
622  // Get the downscale factor for whichever MB trigger exists in the given run
623  std::vector<TString> runtriggers = downscaleOCDB->GetTriggerClasses();
624  Double_t downscalefactor;
625  for (auto i : runtriggers) {
626  if (i.EqualTo(triggerNameMB1) || i.EqualTo(triggerNameMB2)) {
627  downscalefactor = downscaleOCDB->GetDownscaleFactorForTriggerClass(i.Data());
628  break;
629  }
630  }
631 
632  // Store the inverse of the downscale factor, used later to weight the pT spectrum
633  fMBUpscaleFactor = 1/downscalefactor;
634 
635  TString histname = "Trigger/hMBDownscaleFactor";
636  fHistManager.FillTH1(histname.Data(), fMBUpscaleFactor);
637 
638 }
639 
644 {
645  if (!fEventCuts.AcceptEvent(InputEvent()))
646  {
647  PostData(1, fOutput);
648  return kFALSE;
649  }
650  return kTRUE;
651 }
652 
661 {
662  TString histname;
663  AliJetContainer* jetCont = 0;
664  TIter next(&fJetCollArray);
665  while ((jetCont = static_cast<AliJetContainer*>(next()))) {
666  TString jetContName = jetCont->GetName();
667  if (jetContName.Contains("HardCore")) continue;
668 
669  //-----------------------------------------------------------------------------
670  // Find the leading di-jet candidate in each event, and if it satisfies the
671  // trig jet pT threshold, then fill di-jet candidate histogram (regardless of ass jet).
672  // The idea is to study the kinematic selections in post-processing.
673 
674  // Loop over leading hadron cut or not
675  for (Int_t leadingHadronCutType=0; leadingHadronCutType<2; leadingHadronCutType++) {
676 
677  // Find the dijet candidate of the event and store its info in struct fDijet
678  FindDijet(jetCont, leadingHadronCutType);
679 
680  // If we find a dijet candidate (i.e. acceptable trig jet; ass jet accepted or not), fill the di-jet candidate histogram
682  TString histname = TString::Format("%s/DijetCandObservables", jetCont->GetArrayName().Data());
683  FillDijetCandHistograms(histname);
684  }
685 
686  }
687 
688  //---------------------------------------------------------------------------------------------------
689  // Now, study the accepted dijet selection -- specified by the trig/ass jet pT conditions
690 
691  // Find the dijet candidate of the event and store its info in struct fDijet
692  FindDijet(jetCont, 0);
693 
694  // If we find an accepted dijet, fill the dijet imbalance histogram
696  TString histname = TString::Format("%s/DijetImbalanceObservables", jetCont->GetArrayName().Data());
698  }
699 
700  // If we find an accepted dijet, perform momentum-balance study (if requested)
702  histname = TString::Format("%s/MomentumBalance", jetCont->GetArrayName().Data());
703  DoMomentumBalance(histname);
704  }
705 
706  }
707 
708  //---------------------------------------------------------------------------
709  // Do the constituent threshold and geometrical matching study (if requested)
712 
713  return kTRUE;
714 }
715 
724 {
725  fDijet.clear();
726  fDijet.leadingHadronCutType = leadingHadronCutBin;
727 
728  // Get trigger jet
729  AliEmcalJet* trigJet = 0;
730  if (jetCont->GetRhoParameter())
731  trigJet = jetCont->GetLeadingJet("rho");
732  else
733  trigJet = jetCont->GetLeadingJet();
734  if(!trigJet) return;
735 
736  // Skip the event if the leading jet doesn't satisfy the pT threshold
737  Double_t trigJetPt = GetJetPt(jetCont, trigJet);
738  if ( trigJetPt < fMinTrigJetPt ) return;
739 
740  // Skip the event if the leading jet doesn't satisfy the leading hadron threshold
741  if (jetCont->GetLeadingHadronPt(trigJet) < fDijetLeadingHadronPt*leadingHadronCutBin) return;
742 
743  // Fill the dijet struct
744  fDijet.trigJet = trigJet;
745  fDijet.trigJetPt = trigJetPt;
746  fDijet.trigJetEta = trigJet->Eta();
747  fDijet.trigJetPhi = trigJet->Phi();
748 
749  // Find the subleading jet in the opposite hemisphere
750  AliEmcalJet *assJet = 0;
751  for(auto assJetCand : jetCont->accepted()) {
752  if (!assJetCand) continue;
753  Double_t assJetCandPt = GetJetPt(jetCont, assJetCand);
754  if ( TMath::Abs(trigJet->Phi() - assJetCand->Phi()) < fDeltaPhiMin ) continue;
755  if (assJet) {
756  Double_t assJetPt = GetJetPt(jetCont, assJet);
757  if ( assJetCandPt < assJetPt ) continue;
758  }
759  assJet = assJetCand;
760  }
761  if (!assJet) return;
762 
763  // Fill the dijet struct
764  fDijet.assJet = assJet;
765  fDijet.assJetPt = GetJetPt(jetCont, assJet);
766  fDijet.assJetPhi = assJet->Phi();
767  fDijet.assJetEta = assJet->Eta();
769 
770  fDijet.deltaPhi = TMath::Abs(trigJet->Phi() - assJet->Phi());
771  fDijet.deltaEta = trigJet->Eta() - assJet->Eta();
774  fDijet.kTy = TMath::Abs( fDijet.trigJetPt * TMath::Sin(fDijet.deltaPhi) );
775 }
776 
782 {
783 
784  AliTrackContainer* trackCont = dynamic_cast<AliTrackContainer*>(GetParticleContainer("tracks"));
785 
786  AliVTrack* track;
787  for (auto trackIterator : trackCont->accepted_momentum() ) {
788 
789  track = trackIterator.second;
790 
791  // Compute the delta phi between the track and its nearest jet (of the two jets in the dijet),
792  // as well as its pT-parallel projection onto the nearest jet's axis.
793 
794  Double_t trackPt = track->Pt();
795  Double_t trackPhi = track->Phi();
796  Double_t trigJetPhi = fDijet.trigJet->Phi();
797  Double_t assJetPhi = fDijet.assJet->Phi();
798 
799  Double_t deltaPhiTrigJet = TMath::Abs(trackPhi - trigJetPhi);
800  Double_t deltaPhiAssJet = TMath::Abs(trackPhi - assJetPhi);
801  Bool_t isNearside = deltaPhiTrigJet < deltaPhiAssJet;
802 
803  Double_t deltaPhi;
804  Double_t balancePt;
805  if (isNearside) {
806  deltaPhi = trackPhi - trigJetPhi;
807  balancePt = trackPt * TMath::Cos(deltaPhi);
808  }
809  else {
810  deltaPhi = trackPhi - assJetPhi;
811  balancePt = -trackPt * TMath::Cos(deltaPhi);
812  }
813 
814  FillMomentumBalanceHistograms(histname, deltaPhi, trackPt, balancePt);
815 
816  }
817 }
818 
823 {
824  // Get jet container with minimum constituent pT,E thresholds
825  TString jetContAllName = Form("Jet_AKTFullR0%d0_tracks_pT0150_caloClusters_E0300_pt_scheme", (int) (fMatchingJetR*10) );
826  AliJetContainer* jetContAll = GetJetContainer(jetContAllName.Data());
827 
828  // Get jet container with X GeV constituent pT,E thresholds
829  Int_t trackThreshold = (int) (fTrackConstituentThreshold*1000); // in MeV
830  Int_t clusThreshold = (int) (fClusterConstituentThreshold*1000); // in MeV
831  TString jetContHardCoreName = Form("JetHardCore_AKTFullR0%d0_tracks_pT%d_caloClusters_E%d_pt_scheme", (int) (fMatchingJetR*10), trackThreshold, clusThreshold);
832  AliJetContainer* jetContHardCore = GetJetContainer(jetContHardCoreName.Data());
833 
834  // Find the di-jet in the hard-core jet sample, then find the matching di-jet and fill histograms
835  FindDijet(jetContHardCore, 0);
836  if (fDijet.isAccepted) {
837  FindMatchingDijet(jetContAll);
839  }
840 }
841 
847 {
849 
850  // Loop over jets and find leading jet within R of fDijet.trigJet
851  AliEmcalJet *matchingTrigJet = 0;
852  for(auto matchingTrigJetCand : jetCont->accepted()) {
853  if (!matchingTrigJetCand) continue;
854  if ( GetDeltaR(matchingTrigJetCand, fDijet.trigJet) > fMatchingJetR ) continue;
855  if (matchingTrigJet) {
856  if ( GetJetPt(jetCont, matchingTrigJetCand) < GetJetPt(jetCont, matchingTrigJet) ) continue;
857  }
858  matchingTrigJet = matchingTrigJetCand;
859  }
860  if (!matchingTrigJet) return;
861 
862  // Loop over jets and find leading jet within R of fDijet.assJet
863  AliEmcalJet *matchingAssJet = 0;
864  for(auto matchingAssJetCand : jetCont->accepted()) {
865  if (!matchingAssJetCand) continue;
866  if ( GetDeltaR(matchingAssJetCand, fDijet.assJet) > fMatchingJetR ) continue;
867  if (matchingAssJet) {
868  if ( GetJetPt(jetCont, matchingAssJetCand) < GetJetPt(jetCont, matchingAssJet) ) continue;
869  }
870  matchingAssJet = matchingAssJetCand;
871  }
872 
873  // Determine which matching jet is the leading jet (i.e. allow them to flip)
874  if (matchingAssJet) {
875  AliEmcalJet* trigJet = matchingTrigJet;
876  AliEmcalJet* assJet = matchingAssJet;
877  if ( GetJetPt(jetCont, matchingTrigJet) < GetJetPt(jetCont, matchingAssJet) ) {
878  AliEmcalJet* trigJet = matchingAssJet;
879  AliEmcalJet* assJet = matchingTrigJet;
880  }
881 
882  // Fill the dijet struct
883  fMatchingDijet.trigJet = trigJet;
884  fMatchingDijet.trigJetPt = GetJetPt(jetCont, trigJet);
885  fMatchingDijet.trigJetEta = trigJet->Eta();
886  fMatchingDijet.trigJetPhi = trigJet->Phi();
887 
888  fMatchingDijet.assJet = assJet;
889  fMatchingDijet.assJetPt = GetJetPt(jetCont, assJet);
890  fMatchingDijet.assJetPhi = assJet->Phi();
891  fMatchingDijet.assJetEta = assJet->Eta();
893 
894  fMatchingDijet.deltaPhi = TMath::Abs(trigJet->Phi() - assJet->Phi());
895  fMatchingDijet.deltaEta = trigJet->Eta() - assJet->Eta();
898  fMatchingDijet.kTy = TMath::Abs( fMatchingDijet.trigJetPt * TMath::Sin(fMatchingDijet.deltaPhi) );
899  }
900 }
901 
906 {
907  // Loop over tracks and clusters in order to:
908  // (1) Compute scale factor for full jets
909  // (2) Compute delta-pT for full jets, with the random cone method
910  // For both the scale factor and delta-pT, we compute only one histogram each for EMCal.
911  // But for DCal, we bin in eta, in order to study DCal vs. PHOS vs. gap
912 
913  // Define the acceptance boundaries for the TPC and EMCal/DCal/PHOS
914  Double_t etaTPC = 0.9;
915  Double_t etaEMCal = 0.7;
916  Double_t etaMinDCal = 0.22;
917  Double_t etaMaxPHOS = 0.13;
918  Double_t phiMinEMCal = fGeom->GetArm1PhiMin() * TMath::DegToRad(); // 80
919  Double_t phiMaxEMCal = fGeom->GetEMCALPhiMax() * TMath::DegToRad(); // ~188
920  Double_t phiMinDCal = fGeom->GetDCALPhiMin() * TMath::DegToRad(); // 260
921  Double_t phiMaxDCal = fGeom->GetDCALPhiMax() * TMath::DegToRad(); // ~327 (1/3 SMs start at 320)
922  Double_t phiMinPHOS = 250 * TMath::DegToRad();
923  Double_t phiMaxPHOS = 320 * TMath::DegToRad();
924 
925  Double_t accTPC = 2 * etaTPC * 2 * TMath::Pi();
926  Double_t accEMCal = 2 * etaEMCal * (phiMaxEMCal - phiMinEMCal);
927  Double_t accDCalRegion = 2 * etaEMCal * (phiMaxDCal - phiMinDCal);
928 
929  // Define fiducial acceptances, to be used to generate random cones
930  TRandom3* r = new TRandom3(0);
931  Double_t jetR = jetCont->GetJetRadius();
932  Double_t etaEMCalfid = etaEMCal - jetR;
933  Double_t phiMinEMCalfid = phiMinEMCal + jetR;
934  Double_t phiMaxEMCalfid = phiMaxEMCal - jetR;
935  Double_t phiMinDCalRegionfid = phiMinDCal + jetR;
936  Double_t phiMaxDCalRegionfid = phiMaxDCal - jetR;
937 
938  // Generate EMCal random cone eta-phi
939  Double_t etaEMCalRC = r->Uniform(-etaEMCalfid, etaEMCalfid);
940  Double_t phiEMCalRC = r->Uniform(phiMinEMCalfid, phiMaxEMCalfid);
941 
942  // Generate DCalRegion random cone eta-phi inside each eta slice (same phi used for all)
943  Double_t etaStep = 0.05;
944  const Int_t nEtaBinsSF = 28; // 2 * 0.7 / 0.05
945  const Int_t nEtaBinsRC = 20; // 2 * 0.5 / 0.05
946 
947  Double_t phiDCalRC = r->Uniform(phiMinDCalRegionfid, phiMaxDCalRegionfid);
948  Double_t etaDCalRC[nEtaBinsRC];
949  Double_t etaMin;
950  Double_t etaMax;
951  for (Int_t bin=0; bin < nEtaBinsRC; bin++) {
952  etaMin = -etaEMCalfid + bin*etaStep;
953  etaMax = etaMin + etaStep;
954  etaDCalRC[bin] = r->Uniform(etaMin, etaMax);
955  }
956 
957  // Initialize the various sums to 0
958  Double_t trackPtSumTPC = 0;
959  Double_t trackPtSumEMCal = 0;
960  Double_t trackPtSumEMCalRC = 0;
961  Double_t clusESumEMCal = 0;
962  Double_t clusESumEMCalRC = 0;
963  Double_t trackPtSumDCal[nEtaBinsSF] = {0.};
964  Double_t trackPtSumDCalRC[nEtaBinsRC] = {0.};
965  Double_t clusESumDCal[nEtaBinsSF] = {0.};
966  Double_t clusESumDCalRC[nEtaBinsRC] = {0.};
967 
968  // Loop over tracks. Sum the track pT:
969  // (1) in the entire TPC, (2) in the EMCal, (3) in the EMCal random cone,
970  // (4) in the DCalRegion at each eta, (5) in the DCalRegion random cone at each eta
971  AliTrackContainer* trackCont = dynamic_cast<AliTrackContainer*>(GetParticleContainer("tracks"));
972  AliTLorentzVector track;
973  Double_t trackEta;
974  Double_t trackPhi;
975  Double_t trackPt;
976  Double_t deltaR;
977  for (auto trackIterator : trackCont->accepted_momentum() ) {
978 
979  track.Clear();
980  track = trackIterator.first;
981  trackEta = track.Eta();
982  trackPhi = track.Phi_0_2pi();
983  trackPt = track.Pt();
984 
985  // (1)
986  if (TMath::Abs(trackEta) < etaTPC) {
987  trackPtSumTPC += trackPt;
988  }
989 
990  // (2)
991  if (TMath::Abs(trackEta) < etaEMCal && trackPhi > phiMinEMCal && trackPhi < phiMaxEMCal) {
992  trackPtSumEMCal += trackPt;
993  }
994 
995  // (3)
996  deltaR = GetDeltaR(&track, etaEMCalRC, phiEMCalRC);
997  if (deltaR < jetR) {
998  trackPtSumEMCalRC += trackPt;
999  }
1000 
1001  // (4)
1002  for (Int_t bin=0; bin < nEtaBinsSF; bin++) {
1003  if (trackPhi > phiMinDCal && trackPhi < phiMaxDCal) {
1004  etaMin = -etaEMCal+ bin*etaStep;
1005  etaMax = etaMin + etaStep;
1006  if (trackEta > etaMin && trackEta < etaMax) {
1007  trackPtSumDCal[bin] += trackPt;
1008  }
1009  }
1010  }
1011 
1012  // (5)
1013  for (Int_t bin=0; bin < nEtaBinsRC; bin++) {
1014  deltaR = GetDeltaR(&track, etaDCalRC[bin], phiDCalRC);
1015  if (deltaR < jetR) {
1016  trackPtSumDCalRC[bin] += trackPt;
1017  }
1018  }
1019 
1020  }
1021 
1022  // Loop over clusters. Sum the cluster ET:
1023  // (1) in the EMCal, (2) in the EMCal random cone, (3) in the DCalRegion at each eta,
1024  // (4) in the DCalRegion random cone at each eta
1026  AliTLorentzVector clus;
1027  Double_t clusEta;
1028  Double_t clusPhi;
1029  Double_t clusE;
1030  for (auto clusIterator : clusCont->accepted_momentum() ) {
1031 
1032  clus.Clear();
1033  clus = clusIterator.first;
1034  clusEta = clus.Eta();
1035  clusPhi = clus.Phi_0_2pi();
1036  clusE = clus.E();
1037 
1038  // (1)
1039  if (TMath::Abs(clusEta) < etaEMCal && clusPhi > phiMinEMCal && clusPhi < phiMaxEMCal) {
1040  clusESumEMCal += clusE;
1041  }
1042 
1043  // (2)
1044  deltaR = GetDeltaR(&clus, etaEMCalRC, phiEMCalRC);
1045  if (deltaR < jetR) {
1046  clusESumEMCalRC += clusE;
1047  }
1048 
1049  // (3)
1050  for (Int_t bin=0; bin < nEtaBinsSF; bin++) {
1051  if (clusPhi > phiMinDCal && clusPhi < phiMaxDCal) {
1052  etaMin = -etaEMCal+ bin*etaStep;
1053  etaMax = etaMin + etaStep;
1054  if (clusEta > etaMin && clusEta < etaMax) {
1055  clusESumDCal[bin] += clusE;
1056  }
1057  }
1058  }
1059 
1060  // (4)
1061  for (Int_t bin=0; bin < nEtaBinsRC; bin++) {
1062  deltaR = GetDeltaR(&clus, etaDCalRC[bin], phiDCalRC);
1063  if (deltaR < jetR) {
1064  clusESumDCalRC[bin] += clusE;
1065  }
1066  }
1067 
1068  }
1069 
1070  // Compute the scale factor for EMCal
1071  Double_t numerator = (trackPtSumEMCal + clusESumEMCal) / accEMCal;
1072  Double_t denominator = trackPtSumTPC / accTPC;
1073  Double_t scaleFactor = numerator / denominator;
1074  TString histname = TString::Format("%s/BackgroundHistograms/hScaleFactorEMCal", jetCont->GetArrayName().Data());
1075  fHistManager.FillTH2(histname, fCent, scaleFactor);
1076 
1077  // Compute the scale factor for DCalRegion
1078  Double_t accDCalRegionBin = accDCalRegion / nEtaBinsSF;
1079  for (Int_t bin=0; bin < nEtaBinsSF; bin++) {
1080  numerator = (trackPtSumDCal[bin] + clusESumDCal[bin]) / accDCalRegionBin;
1081  scaleFactor = numerator / denominator;
1082  histname = TString::Format("%s/BackgroundHistograms/hScaleFactorDCalRegion", jetCont->GetArrayName().Data());
1083  fHistManager.FillTH3(histname, fCent, scaleFactor, bin);
1084  }
1085 
1086  // Compute delta pT for EMCal
1087  Double_t rho = jetCont->GetRhoVal();
1088  Double_t deltaPt = trackPtSumEMCalRC + clusESumEMCalRC - rho * TMath::Pi() * jetR * jetR;
1089  histname = TString::Format("%s/BackgroundHistograms/hDeltaPtEMCal", jetCont->GetArrayName().Data());
1090  fHistManager.FillTH1(histname, deltaPt);
1091 
1092  // Compute delta pT for DCalRegion
1093  for (Int_t bin=0; bin < nEtaBinsRC; bin++) {
1094  deltaPt = trackPtSumDCalRC[bin] + clusESumDCalRC[bin] - rho * TMath::Pi() * jetR * jetR;
1095  histname = TString::Format("%s/BackgroundHistograms/hDeltaPtDCalRegion", jetCont->GetArrayName().Data());
1096  fHistManager.FillTH2(histname, deltaPt, bin);
1097  }
1098 
1099 }
1100 
1105 {
1106  Double_t pT = jet->Pt() - jetCont->GetRhoVal() * jet->Area();
1107  TString jetContName = jetCont->GetName();
1108  if (jetContName.Contains("HardCore")) pT = jet->Pt();
1109  return pT;
1110 }
1111 
1116 {
1117  Double_t deltaPhi = TMath::Abs(jet1->Phi() - jet2->Phi());
1118  Double_t deltaEta = TMath::Abs(jet1->Eta() - jet2->Eta());
1119  Double_t deltaR = TMath::Sqrt( deltaPhi*deltaPhi + deltaEta*deltaEta );
1120  return deltaR;
1121 }
1122 
1127 {
1128  Double_t deltaPhi = TMath::Abs(part->Phi_0_2pi() - phiRef);
1129  Double_t deltaEta = TMath::Abs(part->Eta() - etaRef);
1130  Double_t deltaR = TMath::Sqrt( deltaPhi*deltaPhi + deltaEta*deltaEta );
1131  return deltaR;
1132 }
1133 
1141 {
1143 
1144  return kTRUE;
1145 }
1146 
1152 {
1153  TString histname;
1154  AliJetContainer* jets = 0;
1155  TIter nextJetColl(&fJetCollArray);
1156  while ((jets = static_cast<AliJetContainer*>(nextJetColl()))) {
1157  TString jetContName = jets->GetName();
1158  if (jetContName.Contains("HardCore")) continue;
1159  Double_t rhoVal = 0;
1160  if (jets->GetRhoParameter()) {
1161  rhoVal = jets->GetRhoVal();
1162  histname = TString::Format("%s/JetHistograms/hRhoVsCent", jets->GetArrayName().Data());
1163  fHistManager.FillTH2(histname.Data(), fCent, rhoVal);
1164  }
1165 
1166  for (auto jet : jets->all()) {
1167 
1168  Float_t ptLeading = jets->GetLeadingHadronPt(jet);
1169  Float_t corrPt = jet->Pt() - rhoVal * jet->Area();
1170 
1171  // A vs. pT (fill before area cut)
1172  histname = TString::Format("%s/JetHistograms/hAreaVsPt", jets->GetArrayName().Data());
1173  fHistManager.FillTH2(histname.Data(), corrPt, jet->Area());
1174 
1175 
1176  // Rejection reason
1177  UInt_t rejectionReason = 0;
1178  if (!jets->AcceptJet(jet, rejectionReason)) {
1179  histname = TString::Format("%s/JetHistograms/hJetRejectionReason", jets->GetArrayName().Data());
1180  fHistManager.FillTH2(histname.Data(), jets->GetRejectionReasonBitPosition(rejectionReason), jet->Pt());
1181  continue;
1182  }
1183 
1184  // Centrality vs. pT
1185  histname = TString::Format("%s/JetHistograms/hCentVsPt", jets->GetArrayName().Data());
1186  fHistManager.FillTH2(histname.Data(), corrPt, fCent);
1187 
1188  // pT vs. eta vs. phi
1189  histname = TString::Format("%s/JetHistograms/hPtVsEtaVsPhi", jets->GetArrayName().Data());
1190  fHistManager.FillTH3(histname.Data(), jet->Eta(), jet->Phi_0_2pi(), corrPt);
1191 
1192  // pT un-downscaled
1193  histname = TString::Format("%s/JetHistograms/hPtUpscaledMB", jets->GetArrayName().Data());
1194  fHistManager.FillTH1(histname.Data(), corrPt, fMBUpscaleFactor);
1195 
1196  // pT-leading vs. pT
1197  histname = TString::Format("%s/JetHistograms/hPtLeadingVsPt", jets->GetArrayName().Data());
1198  fHistManager.FillTH2(histname.Data(), corrPt, ptLeading);
1199 
1200  // NEF vs. pT
1201  histname = TString::Format("%s/JetHistograms/hNEFVsPt", jets->GetArrayName().Data());
1202  fHistManager.FillTH2(histname.Data(), corrPt, jet->NEF());
1203 
1204  // z-leading vs. pT
1205  TLorentzVector leadPart;
1206  jets->GetLeadingHadronMomentum(leadPart, jet);
1207  Double_t z = GetParallelFraction(leadPart.Vect(), jet);
1208  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
1209  histname = TString::Format("%s/JetHistograms/hZLeadingVsPt", jets->GetArrayName().Data());
1210  fHistManager.FillTH2(histname.Data(), corrPt, z);
1211 
1212  // Nconst vs. pT
1213  histname = TString::Format("%s/JetHistograms/hNConstVsPt", jets->GetArrayName().Data());
1214  fHistManager.FillTH2(histname.Data(), corrPt, jet->GetNumberOfConstituents());
1215 
1216 
1217  } //jet loop
1218 
1219  //---------------------------------------------------------------------------
1220  // Do study of background (if requested)
1222  }
1223 }
1224 
1229 {
1230  Double_t contents[30]={0};
1231  THnSparse* histJetObservables = static_cast<THnSparse*>(fHistManager.FindObject(histname));
1232  if (!histJetObservables) return;
1233  for (Int_t n = 0; n < histJetObservables->GetNdimensions(); n++) {
1234  TString title(histJetObservables->GetAxis(n)->GetTitle());
1235  if (title=="Centrality (%)")
1236  contents[n] = fCent;
1237  else if (title=="LeadingHadronRequired")
1238  contents[n] = fDijet.leadingHadronCutType;
1239  else if (title=="#it{p}_{T,trig jet} (GeV/#it{c})")
1240  contents[n] = fDijet.trigJetPt;
1241  else if (title=="#it{p}_{T,ass jet} (GeV/#it{c})")
1242  contents[n] = fDijet.assJetPt;
1243  else if (title=="#phi_{trig jet}")
1244  contents[n] = fDijet.trigJetPhi;
1245  else if (title=="#phi_{ass jet}")
1246  contents[n] = fDijet.assJetPhi;
1247  else if (title=="#eta_{trig jet}")
1248  contents[n] = fDijet.trigJetEta;
1249  else if (title=="#eta_{ass jet}")
1250  contents[n] = fDijet.assJetEta;
1251  else
1252  AliWarning(Form("Unable to fill dimension %s!",title.Data()));
1253  }
1254  histJetObservables->Fill(contents);
1255 }
1256 
1261 {
1262  Double_t contents[30]={0};
1263  THnSparse* histJetObservables = static_cast<THnSparse*>(fHistManager.FindObject(histname));
1264  if (!histJetObservables) return;
1265  for (Int_t n = 0; n < histJetObservables->GetNdimensions(); n++) {
1266  TString title(histJetObservables->GetAxis(n)->GetTitle());
1267  if (title=="Centrality (%)")
1268  contents[n] = fCent;
1269  else if (title=="#Delta#phi")
1270  contents[n] = fDijet.deltaPhi;
1271  else if (title=="#Delta#eta")
1272  contents[n] = fDijet.deltaEta;
1273  else if (title=="A_{J}")
1274  contents[n] = fDijet.AJ;
1275  else if (title=="x_{J}")
1276  contents[n] = fDijet.xJ;
1277  else if (title=="k_{Ty} (GeV)")
1278  contents[n] = fDijet.kTy;
1279  else if (title=="N_{tracks, trig jet}")
1280  contents[n] = fDijet.trigJet->GetNumberOfTracks();
1281  else if (title=="N_{tracks, ass jet}")
1282  contents[n] = fDijet.assJet->GetNumberOfTracks();
1283  else
1284  AliWarning(Form("Unable to fill dimension %s!",title.Data()));
1285  }
1286  histJetObservables->Fill(contents);
1287 }
1288 
1293 {
1294  Double_t contents[30]={0};
1295  THnSparse* histJetObservables = static_cast<THnSparse*>(fHistManager.FindObject(histname));
1296  if (!histJetObservables) return;
1297  for (Int_t n = 0; n < histJetObservables->GetNdimensions(); n++) {
1298  TString title(histJetObservables->GetAxis(n)->GetTitle());
1299  if (title=="A_{J}")
1300  contents[n] = fDijet.AJ;
1301  else if (title=="#Delta#phi")
1302  contents[n] = deltaPhi;
1303  else if (title=="#it{p}_{T,particle} (GeV/#it{c})")
1304  contents[n] = trackPt;
1305  else if (title=="#it{p}_{T}#parallel (GeV/#it{c})")
1306  contents[n] = balancePt;
1307  else
1308  AliWarning(Form("Unable to fill dimension %s!",title.Data()));
1309  }
1310  histJetObservables->Fill(contents);
1311 }
1312 
1317 {
1318  // Matching efficiency histogram
1319  TString histname = "GeometricalMatchingEfficiency";
1320 
1321  // If we have a matching di-jet, fill the geometrical matching histogram
1322  if (fMatchingDijet.assJet) {
1323  fHistManager.FillTH1(histname.Data(), 1);
1324 
1327  Bool_t isSwitched = trigDeltaR > fMatchingJetR;
1328 
1329  TString thnname = "GeometricalMatching";
1330  Double_t contents[30]={0};
1331  THnSparse* histJetObservables = static_cast<THnSparse*>(fHistManager.FindObject(thnname.Data()));
1332  if (!histJetObservables) return;
1333  for (Int_t n = 0; n < histJetObservables->GetNdimensions(); n++) {
1334  TString title(histJetObservables->GetAxis(n)->GetTitle());
1335  if (title=="Centrality (%)")
1336  contents[n] = fCent;
1337  else if (title=="isSwitched")
1338  contents[n] = isSwitched;
1339  else if (title=="#DeltaR_{trig}")
1340  contents[n] = trigDeltaR;
1341  else if (title=="#DeltaR_{ass}")
1342  contents[n] = assDeltaR;
1343  else if (title=="trig #it{p}_{T,low-thresh} - #it{p}_{T,hard-core}")
1344  contents[n] = fMatchingDijet.trigJetPt - fDijet.trigJetPt;
1345  else if (title=="ass #it{p}_{T,low-thresh} - #it{p}_{T,hard-core}")
1346  contents[n] = fMatchingDijet.assJetPt - fDijet.assJetPt;
1347  else if (title=="A_{J} low-threshold")
1348  contents[n] = fMatchingDijet.AJ;
1349  else if (title=="A_{J} hard-core")
1350  contents[n] = fDijet.AJ;
1351  else
1352  AliWarning(Form("Unable to fill dimension %s!",title.Data()));
1353  }
1354  histJetObservables->Fill(contents);
1355  }
1356  else {
1357  fHistManager.FillTH1(histname.Data(), 0.);
1358  }
1359 }
1360 
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