AliPhysics  86c65ee (86c65ee)
 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 
21 #include <AliVCluster.h>
22 #include <AliVParticle.h>
23 #include <AliLog.h>
24 
25 #include "AliTLorentzVector.h"
26 #include "AliEmcalJet.h"
27 #include "AliRhoParameter.h"
28 #include "AliJetContainer.h"
29 #include "AliParticleContainer.h"
30 #include "AliClusterContainer.h"
31 
33 
37 
43  fHistManager(),
44  fDeltaPhiMin(0),
45  fTrigJetMinPt(0),
46  fAssJetMinPt(0)
47 {
48 }
49 
56  AliAnalysisTaskEmcalJet(name, kTRUE),
57  fHistManager(name),
58  fDeltaPhiMin(2),
59  fTrigJetMinPt(0),
60  fAssJetMinPt(0)
61 {
63 }
64 
69 {
70 }
71 
77 {
79 
85 
86  TIter next(fHistManager.GetListOfHistograms());
87  TObject* obj = 0;
88  while ((obj = next())) {
89  fOutput->Add(obj);
90  }
91 
92  PostData(1, fOutput); // Post data for ALL output slots > 0 here.
93 }
94 
95 /*
96  * This function allocates the histograms for basic EMCal cluster QA.
97  * A set of histograms (energy, eta, phi, number of cluster) is allocated
98  * per each cluster container and per each centrality bin.
99  */
101 {
102  TString histname;
103  TString histtitle;
104  TString groupname;
105  AliClusterContainer* clusCont = 0;
106  TIter next(&fClusterCollArray);
107  while ((clusCont = static_cast<AliClusterContainer*>(next()))) {
108  groupname = clusCont->GetName();
109  fHistManager.CreateHistoGroup(groupname);
110  for (Int_t cent = 0; cent < fNcentBins; cent++) {
111 
112  // Cluster histograms (PHOS+EMCal)
113 
114  histname = TString::Format("%s/histClusterEnergy_%d", groupname.Data(), cent);
115  histtitle = TString::Format("%s;#it{E}_{cluster} (GeV);counts", histname.Data());
116  fHistManager.CreateTH1(histname, histtitle, fNbins / 2, fMinBinPt, fMaxBinPt / 2);
117 
118  histname = TString::Format("%s/histClusterEtaPhi_%d", groupname.Data(), cent);
119  histtitle = TString::Format("%s;#it{#eta}_{cluster};#it{#phi}_{cluster};counts", histname.Data());
120  fHistManager.CreateTH2(histname, histtitle, fNbins / 6, -1, 1, fNbins / 2, 0, TMath::TwoPi());
121 
122  histname = TString::Format("%s/histNClusters_%d", groupname.Data(), cent);
123  histtitle = TString::Format("%s;number of clusters;events", histname.Data());
124  if (fForceBeamType != kpp) {
125  fHistManager.CreateTH1(histname, histtitle, 500, 0, 3000);
126  }
127  else {
128  fHistManager.CreateTH1(histname, histtitle, 200, 0, 200);
129  }
130 
131  // EMCal cluster histograms
132 
133  histname = TString::Format("%s/histEMCalClusterEnergy_%d", groupname.Data(), cent);
134  histtitle = TString::Format("%s;#it{E}_{cluster} (GeV);counts", histname.Data());
135  fHistManager.CreateTH1(histname, histtitle, fNbins / 2, fMinBinPt, fMaxBinPt / 2);
136 
137  histname = TString::Format("%s/histEMCalClusterEnergyExotic_%d", groupname.Data(), cent);
138  histtitle = TString::Format("%s;#it{E}_{cluster}^{exotic} (GeV);counts", histname.Data());
139  fHistManager.CreateTH1(histname, histtitle, fNbins / 2, fMinBinPt, fMaxBinPt / 2);
140 
141  histname = TString::Format("%s/histEMCalClusterNonLinCorrEnergy_%d", groupname.Data(), cent);
142  histtitle = TString::Format("%s;#it{E}_{cluster}^{non-lin.corr.} (GeV);counts", histname.Data());
143  fHistManager.CreateTH1(histname, histtitle, fNbins / 2, fMinBinPt, fMaxBinPt / 2);
144 
145  histname = TString::Format("%s/histEMCalClusterHadCorrEnergy_%d", groupname.Data(), cent);
146  histtitle = TString::Format("%s;#it{E}_{cluster}^{had.corr.} (GeV);counts", histname.Data());
147  fHistManager.CreateTH1(histname, histtitle, fNbins / 2, fMinBinPt, fMaxBinPt / 2);
148 
149  histname = TString::Format("%s/histEMCalClusterPhi_%d", groupname.Data(), cent);
150  histtitle = TString::Format("%s;#it{#phi}_{cluster};counts", histname.Data());
151  fHistManager.CreateTH1(histname, histtitle, fNbins / 2, 0, TMath::TwoPi());
152 
153  histname = TString::Format("%s/histEMCalClusterEta_%d", groupname.Data(), cent);
154  histtitle = TString::Format("%s;#it{#eta}_{cluster};counts", histname.Data());
155  fHistManager.CreateTH1(histname, histtitle, fNbins / 6, -1, 1);
156 
157  histname = TString::Format("%s/histEMCalNClusters_%d", groupname.Data(), cent);
158  histtitle = TString::Format("%s;number of clusters;events", histname.Data());
159  if (fForceBeamType != kpp) {
160  fHistManager.CreateTH1(histname, histtitle, 500, 0, 3000);
161  }
162  else {
163  fHistManager.CreateTH1(histname, histtitle, 200, 0, 200);
164  }
165 
166  // PHOS cluster histograms
167 
168  histname = TString::Format("%s/histPHOSClusterEnergy_%d", groupname.Data(), cent);
169  histtitle = TString::Format("%s;#it{E}_{cluster} (GeV);counts", histname.Data());
170  fHistManager.CreateTH1(histname, histtitle, fNbins / 2, fMinBinPt, fMaxBinPt / 2);
171 
172  histname = TString::Format("%s/histPHOSClusterPhi_%d", groupname.Data(), cent);
173  histtitle = TString::Format("%s;#it{#phi}_{cluster};counts", histname.Data());
174  fHistManager.CreateTH1(histname, histtitle, fNbins / 2, 0, TMath::TwoPi());
175 
176  histname = TString::Format("%s/histPHOSClusterEta_%d", groupname.Data(), cent);
177  histtitle = TString::Format("%s;#it{#eta}_{cluster};counts", histname.Data());
178  fHistManager.CreateTH1(histname, histtitle, fNbins / 6, -1, 1);
179 
180  histname = TString::Format("%s/histPHOSNClusters_%d", groupname.Data(), cent);
181  histtitle = TString::Format("%s;number of clusters;events", histname.Data());
182  if (fForceBeamType != kpp) {
183  fHistManager.CreateTH1(histname, histtitle, 500, 0, 3000);
184  }
185  else {
186  fHistManager.CreateTH1(histname, histtitle, 200, 0, 200);
187  }
188 
189  }
190  }
191 }
192 
193 /*
194  * This function allocates the histograms for basic EMCal QA.
195  * One 2D histogram with the cell energy spectra and the number of cells
196  * per event is allocated per each centrality bin.
197  */
199 {
200  TString histname;
201  TString histtitle;
202  TString groupname(fCaloCellsName);
203 
204  fHistManager.CreateHistoGroup(groupname);
205  for (Int_t cent = 0; cent < fNcentBins; cent++) {
206  histname = TString::Format("%s/histCellEnergyvsAbsId_%d", groupname.Data(), cent);
207  histtitle = TString::Format("%s;cell abs. ID;#it{E}_{cell} (GeV);counts", histname.Data());
208  fHistManager.CreateTH2(histname, histtitle, 20000, 0, 20000, fNbins / 2, fMinBinPt, fMaxBinPt / 2);
209 
210  histname = TString::Format("%s/histNCells_%d", groupname.Data(), cent);
211  histtitle = TString::Format("%s;number of cells;events", histname.Data());
212  if (fForceBeamType != kpp) {
213  fHistManager.CreateTH1(histname, histtitle, 500, 0, 6000);
214  }
215  else {
216  fHistManager.CreateTH1(histname, histtitle, 200, 0, 200);
217  }
218  }
219 }
220 
221 /*
222  * This function allocates the histograms for basic tracking QA.
223  * A set of histograms (pT, eta, phi, difference between kinematic properties
224  * at the vertex and at the EMCal surface, number of tracks) is allocated
225  * per each particle container and per each centrality bin.
226  */
228 {
229  TString histname;
230  TString histtitle;
231  TString groupname;
232  AliParticleContainer* partCont = 0;
233  TIter next(&fParticleCollArray);
234  while ((partCont = static_cast<AliParticleContainer*>(next()))) {
235  groupname = partCont->GetName();
236  fHistManager.CreateHistoGroup(groupname);
237  for (Int_t cent = 0; cent < fNcentBins; cent++) {
238  histname = TString::Format("%s/histTrackPt_%d", groupname.Data(), cent);
239  histtitle = TString::Format("%s;#it{p}_{T,track} (GeV/#it{c});counts", histname.Data());
240  fHistManager.CreateTH1(histname, histtitle, fNbins / 2, fMinBinPt, fMaxBinPt / 2);
241 
242  histname = TString::Format("%s/histTrackPhi_%d", groupname.Data(), cent);
243  histtitle = TString::Format("%s;#it{#phi}_{track};counts", histname.Data());
244  fHistManager.CreateTH1(histname, histtitle, fNbins / 2, 0, TMath::TwoPi());
245 
246  histname = TString::Format("%s/histTrackEta_%d", groupname.Data(), cent);
247  histtitle = TString::Format("%s;#it{#eta}_{track};counts", histname.Data());
248  fHistManager.CreateTH1(histname, histtitle, fNbins / 6, -1, 1);
249 
250  if (TClass(partCont->GetClassName()).InheritsFrom("AliVTrack")) {
251  histname = TString::Format("%s/fHistDeltaEtaPt_%d", groupname.Data(), cent);
252  histtitle = TString::Format("%s;#it{p}_{T,track}^{vertex} (GeV/#it{c});#it{#eta}_{track}^{vertex} - #it{#eta}_{track}^{EMCal};counts", histname.Data());
253  fHistManager.CreateTH2(histname, histtitle, fNbins / 2, fMinBinPt, fMaxBinPt, 50, -0.5, 0.5);
254 
255  histname = TString::Format("%s/fHistDeltaPhiPt_%d", groupname.Data(), cent);
256  histtitle = TString::Format("%s;#it{p}_{T,track}^{vertex} (GeV/#it{c});#it{#phi}_{track}^{vertex} - #it{#phi}_{track}^{EMCal};counts", histname.Data());
257  fHistManager.CreateTH2(histname, histtitle, fNbins / 2, fMinBinPt, fMaxBinPt, 200, -2, 2);
258 
259  histname = TString::Format("%s/fHistDeltaPtvsPt_%d", groupname.Data(), cent);
260  histtitle = TString::Format("%s;#it{p}_{T,track}^{vertex} (GeV/#it{c});#it{p}_{T,track}^{vertex} - #it{p}_{T,track}^{EMCal} (GeV/#it{c});counts", histname.Data());
261  fHistManager.CreateTH2(histname, histtitle, fNbins / 2, fMinBinPt, fMaxBinPt, fNbins / 2, -fMaxBinPt/2, fMaxBinPt/2);
262 
263  histname = TString::Format("%s/fHistEoverPvsP_%d", groupname.Data(), cent);
264  histtitle = TString::Format("%s;#it{P}_{track} (GeV/#it{c});#it{E}_{cluster} / #it{P}_{track} #it{c};counts", histname.Data());
265  fHistManager.CreateTH2(histname, histtitle, fNbins / 2, fMinBinPt, fMaxBinPt, fNbins / 2, 0, 4);
266  }
267 
268  histname = TString::Format("%s/histNTracks_%d", groupname.Data(), cent);
269  histtitle = TString::Format("%s;number of tracks;events", histname.Data());
270  if (fForceBeamType != kpp) {
271  fHistManager.CreateTH1(histname, histtitle, 500, 0, 5000);
272  }
273  else {
274  fHistManager.CreateTH1(histname, histtitle, 200, 0, 200);
275  }
276  }
277  }
278 }
279 
280 /*
281  * This function allocates the histograms for basic jet QA.
282  * A set of histograms (pT, eta, phi, area, number of jets, corrected pT) is allocated
283  * per each jet container and per each centrality bin.
284  */
286 {
287  TString histname;
288  TString histtitle;
289  TString groupname;
290  AliJetContainer* jetCont = 0;
291  TIter next(&fJetCollArray);
292  while ((jetCont = static_cast<AliJetContainer*>(next()))) {
293  groupname = jetCont->GetName();
294  fHistManager.CreateHistoGroup(groupname);
295  for (Int_t cent = 0; cent < fNcentBins; cent++) {
296  histname = TString::Format("%s/histJetPt_%d", groupname.Data(), cent);
297  histtitle = TString::Format("%s;#it{p}_{T,jet} (GeV/#it{c});counts", histname.Data());
298  fHistManager.CreateTH1(histname, histtitle, fNbins, fMinBinPt, fMaxBinPt);
299 
300  histname = TString::Format("%s/histJetArea_%d", groupname.Data(), cent);
301  histtitle = TString::Format("%s;#it{A}_{jet};counts", histname.Data());
302  fHistManager.CreateTH1(histname, histtitle, fNbins / 2, 0, 1.5);
303 
304  histname = TString::Format("%s/histJetPhi_%d", groupname.Data(), cent);
305  histtitle = TString::Format("%s;#it{#phi}_{jet};counts", histname.Data());
306  fHistManager.CreateTH1(histname, histtitle, fNbins / 2, 0, TMath::TwoPi());
307 
308  histname = TString::Format("%s/histJetEta_%d", groupname.Data(), cent);
309  histtitle = TString::Format("%s;#it{#eta}_{jet};counts", histname.Data());
310  fHistManager.CreateTH1(histname, histtitle, fNbins / 6, -1, 1);
311 
312  histname = TString::Format("%s/histJetEtaPhi_%d", groupname.Data(), cent);
313  histtitle = TString::Format("%s;#it{#eta}_{jet};#it{#phi}_{jet};counts", histname.Data());
314  fHistManager.CreateTH2(histname, histtitle, fNbins / 6, -1, 1, fNbins / 2, 0, TMath::TwoPi());
315 
316  histname = TString::Format("%s/histNJets_%d", groupname.Data(), cent);
317  histtitle = TString::Format("%s;number of jets;events", histname.Data());
318  if (fForceBeamType != kpp) {
319  fHistManager.CreateTH1(histname, histtitle, 500, 0, 500);
320  }
321  else {
322  fHistManager.CreateTH1(histname, histtitle, 100, 0, 100);
323  }
324 
325  if (!jetCont->GetRhoName().IsNull()) {
326  histname = TString::Format("%s/histJetCorrPt_%d", groupname.Data(), cent);
327  histtitle = TString::Format("%s;#it{p}_{T,jet}^{corr} (GeV/#it{c});counts", histname.Data());
328  fHistManager.CreateTH1(histname, histtitle, fNbins, -fMaxBinPt / 2, fMaxBinPt / 2);
329 
330  histname = TString::Format("%s/histJetRho_%d", groupname.Data(), cent);
331  histtitle = TString::Format("%s;{#rho} (GeV);counts", histname.Data());
332  fHistManager.CreateTH1(histname, histtitle, fNbins, 0, 500);
333  }
334  }
335  }
336 }
337 
338 /*
339  * This function allocates the histograms for the dijet analysis.
340  * A set of histograms is allocated per each jet container and per each centrality bin.
341  */
343 {
344  TString histname;
345  TString histtitle;
346  TString groupname;
347  AliJetContainer* jetCont = 0;
348  TIter next(&fJetCollArray);
349  while ((jetCont = static_cast<AliJetContainer*>(next()))) {
350  groupname = jetCont->GetName();
351  //fHistManager.CreateHistoGroup(groupname);
352  for (Int_t cent = 0; cent < fNcentBins; cent++) {
353 
354  // Loop over leading hadron requirement
355  for (Int_t k=0; k<2; k++) {
356 
357  // Loop over trigJetMinPt thresholds
358  for (Int_t i=0; i<4; i++) {
359 
360  // Loop over assJetMinPt thresholds
361  for (Int_t j=0; j<4; j++) {
362 
363  TString label = TString::Format("_%d_had%d_trig%d_ass%d", cent, k, i, j);
364 
365  // Some basic plots for di-jet pairs
366 
367  histname = TString::Format("%s/histDijetLeadingJetPt", groupname.Data()) += label;
368  histtitle = TString::Format("%s;Leading Jet p_{T} (GeV);counts", histname.Data());
369  fHistManager.CreateTH1(histname, histtitle, fNbins, fMinBinPt, fMaxBinPt);
370 
371  histname = TString::Format("%s/histDijetLeadingJetPtuncorr", groupname.Data()) += label;
372  histtitle = TString::Format("%s;Uncorrected Leading Jet p_{T} (GeV);counts", histname.Data());
373  fHistManager.CreateTH1(histname, histtitle, fNbins, fMinBinPt, fMaxBinPt);
374 
375  histname = TString::Format("%s/histDijetSubleadingJetPt", groupname.Data()) += label;
376  histtitle = TString::Format("%s;Subleading Jet p_{T} (GeV);counts", histname.Data());
377  fHistManager.CreateTH1(histname, histtitle, fNbins, fMinBinPt, fMaxBinPt);
378 
379  histname = TString::Format("%s/histDijetLeadingJetPhi", groupname.Data()) += label;
380  histtitle = TString::Format("%s;Leading Jet #phi;counts", histname.Data());
381  fHistManager.CreateTH1(histname, histtitle, 100, 0, TMath::TwoPi());
382 
383  histname = TString::Format("%s/histDijetSubleadingJetPhi", groupname.Data()) += label;
384  histtitle = TString::Format("%s;Subleading Jet #phi;counts", histname.Data());
385  fHistManager.CreateTH1(histname, histtitle, 100, 0, TMath::TwoPi());
386 
387  //histname = TString::Format("%s/histDijetLeadingJetEta_%d", groupname.Data(), cent);
388  //histtitle = TString::Format("%s;Leading Jet #eta;counts", histname.Data());
389  //fHistManager.CreateTH1(histname, histtitle, 100, -1, 1);
390 
391  //histname = TString::Format("%s/histDijetSubleadingJetEta_%d", groupname.Data(), cent);
392  //histtitle = TString::Format("%s;Subleading Jet #eta;counts", histname.Data());
393  //fHistManager.CreateTH1(histname, histtitle, 100, -1, 1);
394 
395  //histname = TString::Format("%s/histDijetDeltaEta_%d", groupname.Data(), cent);
396  //histtitle = TString::Format("%s;#Delta#eta;counts", histname.Data());
397  //fHistManager.CreateTH1(histname, histtitle, 100, -2, 2);
398 
399  // Some di-jet observables
400 
401  histname = TString::Format("%s/histDijetAJ", groupname.Data()) += label;
402  histtitle = TString::Format("%s;A_{J};counts", histname.Data());
403  fHistManager.CreateTH1(histname, histtitle, 100, 0, 1);
404 
405  histname = TString::Format("%s/histDijetxJ", groupname.Data()) += label;
406  histtitle = TString::Format("%s;x_{J};counts", histname.Data());
407  fHistManager.CreateTH1(histname, histtitle, 100, 0, 1);
408 
409  histname = TString::Format("%s/histDijetDeltaPhi", groupname.Data()) += label;
410  histtitle = TString::Format("%s;#Delta#phi;counts", histname.Data());
411  fHistManager.CreateTH1(histname, histtitle, 100, 0, 4);
412 
413  //histname = TString::Format("%s/histDijetkT_%d", groupname.Data(), cent);
414  //histtitle = TString::Format("%s;k_{Ty} (GeV);counts", histname.Data());
415  //fHistManager.CreateTH1(histname, histtitle, 100, 0, 100);
416 
417  //histname = TString::Format("%s/histDijetLeadingJetNTracks_%d", groupname.Data(), cent);
418  //histtitle = TString::Format("%s;Leading Jet N tracks;counts", histname.Data());
419  //fHistManager.CreateTH1(histname, histtitle, 100, 0, 100);
420 
421  //histname = TString::Format("%s/histDijetSubleadingJetNTracks_%d", groupname.Data(), cent);
422  //histtitle = TString::Format("%s;Subleading Jet N tracks;counts", histname.Data());
423  //fHistManager.CreateTH1(histname, histtitle, 100, 0, 100);
424 
425  //histname = TString::Format("%s/histDijetLeadingJetArea_%d", groupname.Data(), cent);
426  //histtitle = TString::Format("%s;Leading Jet Area;counts", histname.Data());
427  //fHistManager.CreateTH1(histname, histtitle, 100, 0, 3);
428 
429  //histname = TString::Format("%s/histDijetSubleadingJetArea_%d", groupname.Data(), cent);
430  //histtitle = TString::Format("%s;Subleading Jet Area;counts", histname.Data());
431  //fHistManager.CreateTH1(histname, histtitle, 100, 0, 3);
432 
433  // Leading jets that don't have acceptable associated jet
434  histname = TString::Format("%s/histUnmatchedLeadingJetPt", groupname.Data()) += label;
435  histtitle = TString::Format("%s;Leading Jet p_{T} (GeV);counts", histname.Data());
436  fHistManager.CreateTH1(histname, histtitle, fNbins, fMinBinPt, fMaxBinPt);
437 
438  // When leading jet doesn't have an acceptable associated jet, plot the subleading jet in the event
439  histname = TString::Format("%s/histUnmatchedSubleadingJetPt", groupname.Data()) += label;
440  histtitle = TString::Format("%s;Subleading Jet p_{T} (GeV);counts", histname.Data());
441  fHistManager.CreateTH1(histname, histtitle, fNbins, fMinBinPt, fMaxBinPt);
442 
443  }
444  }
445  }
446  }
447  }
448 }
449 
457 {
458  DoJetLoop();
459  DoTrackLoop();
460  DoClusterLoop();
461  DoCellLoop();
462 
463  return kTRUE;
464 }
465 
471 {
472  TString histname;
473  TString groupname;
474  AliJetContainer* jetCont = 0;
475  TIter next(&fJetCollArray);
476  while ((jetCont = static_cast<AliJetContainer*>(next()))) {
477  groupname = jetCont->GetName();
478  UInt_t count = 0;
479  for(auto jet : jetCont->accepted()) {
480  if (!jet) continue;
481  count++;
482 
483  histname = TString::Format("%s/histJetPt_%d", groupname.Data(), fCentBin);
484  fHistManager.FillTH1(histname, jet->Pt());
485 
486  histname = TString::Format("%s/histJetArea_%d", groupname.Data(), fCentBin);
487  fHistManager.FillTH1(histname, jet->Area());
488 
489  histname = TString::Format("%s/histJetPhi_%d", groupname.Data(), fCentBin);
490  fHistManager.FillTH1(histname, jet->Phi());
491 
492  histname = TString::Format("%s/histJetEta_%d", groupname.Data(), fCentBin);
493  fHistManager.FillTH1(histname, jet->Eta());
494 
495  histname = TString::Format("%s/histJetEtaPhi_%d", groupname.Data(), fCentBin);
496  fHistManager.FillTH1(histname, jet->Eta(), jet->Phi());
497 
498  if (jetCont->GetRhoParameter()) {
499  histname = TString::Format("%s/histJetCorrPt_%d", groupname.Data(), fCentBin);
500  fHistManager.FillTH1(histname, jet->Pt() - jetCont->GetRhoVal() * jet->Area());
501  }
502  }
503  histname = TString::Format("%s/histNJets_%d", groupname.Data(), fCentBin);
504  fHistManager.FillTH1(histname, count);
505 
506  if (jetCont->GetRhoParameter()) {
507  histname = TString::Format("%s/histJetRho_%d", groupname.Data(), fCentBin);
508  fHistManager.FillTH1(histname, jetCont->GetRhoVal());
509  }
510  }
511 }
512 
518 {
520 
521  TString histname;
522  TString groupname;
523  AliParticleContainer* partCont = 0;
524  TIter next(&fParticleCollArray);
525  while ((partCont = static_cast<AliParticleContainer*>(next()))) {
526  groupname = partCont->GetName();
527  UInt_t count = 0;
528  for(auto part : partCont->accepted()) {
529  if (!part) continue;
530  count++;
531 
532  histname = TString::Format("%s/histTrackPt_%d", groupname.Data(), fCentBin);
533  fHistManager.FillTH1(histname, part->Pt());
534 
535  histname = TString::Format("%s/histTrackPhi_%d", groupname.Data(), fCentBin);
536  fHistManager.FillTH1(histname, part->Phi());
537 
538  histname = TString::Format("%s/histTrackEta_%d", groupname.Data(), fCentBin);
539  fHistManager.FillTH1(histname, part->Eta());
540 
541  if (partCont->GetLoadedClass()->InheritsFrom("AliVTrack")) {
542  const AliVTrack* track = static_cast<const AliVTrack*>(part);
543 
544  histname = TString::Format("%s/fHistDeltaEtaPt_%d", groupname.Data(), fCentBin);
545  fHistManager.FillTH1(histname, track->Pt(), track->Eta() - track->GetTrackEtaOnEMCal());
546 
547  histname = TString::Format("%s/fHistDeltaPhiPt_%d", groupname.Data(), fCentBin);
548  fHistManager.FillTH1(histname, track->Pt(), track->Phi() - track->GetTrackPhiOnEMCal());
549 
550  histname = TString::Format("%s/fHistDeltaPtvsPt_%d", groupname.Data(), fCentBin);
551  fHistManager.FillTH1(histname, track->Pt(), track->Pt() - track->GetTrackPtOnEMCal());
552 
553  if (clusCont) {
554  Int_t iCluster = track->GetEMCALcluster();
555  if (iCluster >= 0) {
556  AliVCluster* cluster = clusCont->GetAcceptCluster(iCluster);
557  if (cluster) {
558  histname = TString::Format("%s/fHistEoverPvsP_%d", groupname.Data(), fCentBin);
559  fHistManager.FillTH2(histname, track->P(), cluster->GetNonLinCorrEnergy() / track->P());
560  }
561  }
562  }
563  }
564  }
565 
566  histname = TString::Format("%s/histNTracks_%d", groupname.Data(), fCentBin);
567  fHistManager.FillTH1(histname, count);
568  }
569 }
570 
576 {
577  TString histname;
578  TString groupname;
579  AliClusterContainer* clusCont = 0;
580  TIter next(&fClusterCollArray);
581  while ((clusCont = static_cast<AliClusterContainer*>(next()))) {
582  groupname = clusCont->GetName();
583 
584  for(auto cluster : clusCont->all()) {
585  if (!cluster) continue;
586 
587  if (cluster->GetIsExotic()) {
588  histname = TString::Format("%s/histEMCalClusterEnergyExotic_%d", groupname.Data(), fCentBin);
589  fHistManager.FillTH1(histname, cluster->E());
590  }
591  }
592 
593  UInt_t count = 0;
594  UInt_t countEMCal = 0;
595  UInt_t countPHOS = 0;
596  for(auto cluster : clusCont->accepted()) {
597  if (!cluster) continue;
598  count++;
599 
600  AliTLorentzVector nPart;
601  cluster->GetMomentum(nPart, fVertex);
602 
603  histname = TString::Format("%s/histClusterEnergy_%d", groupname.Data(), fCentBin);
604  fHistManager.FillTH1(histname, cluster->E());
605 
606  histname = TString::Format("%s/histClusterEtaPhi_%d", groupname.Data(),fCentBin);
607  fHistManager.FillTH2(histname, nPart.Eta(), nPart.Phi_0_2pi());
608 
609  if (cluster->IsEMCAL()) {
610 
611  countEMCal++;
612 
613  histname = TString::Format("%s/histEMCalClusterEnergy_%d", groupname.Data(), fCentBin);
614  fHistManager.FillTH1(histname, cluster->E());
615 
616  histname = TString::Format("%s/histEMCalClusterNonLinCorrEnergy_%d", groupname.Data(), fCentBin);
617  fHistManager.FillTH1(histname, cluster->GetNonLinCorrEnergy());
618 
619  histname = TString::Format("%s/histEMCalClusterHadCorrEnergy_%d", groupname.Data(), fCentBin);
620  fHistManager.FillTH1(histname, cluster->GetHadCorrEnergy());
621 
622  histname = TString::Format("%s/histEMCalClusterPhi_%d", groupname.Data(), fCentBin);
623  fHistManager.FillTH1(histname, nPart.Phi_0_2pi());
624 
625  histname = TString::Format("%s/histEMCalClusterEta_%d", groupname.Data(), fCentBin);
626  fHistManager.FillTH1(histname, nPart.Eta());
627 
628  } else if (cluster->IsPHOS()){
629 
630  countPHOS++;
631 
632  histname = TString::Format("%s/histPHOSClusterEnergy_%d", groupname.Data(), fCentBin);
633  fHistManager.FillTH1(histname, cluster->E());
634 
635  histname = TString::Format("%s/histPHOSClusterPhi_%d", groupname.Data(), fCentBin);
636  fHistManager.FillTH1(histname, nPart.Phi_0_2pi());
637 
638  histname = TString::Format("%s/histPHOSClusterEta_%d", groupname.Data(), fCentBin);
639  fHistManager.FillTH1(histname, nPart.Eta());
640 
641  }
642  }
643 
644  histname = TString::Format("%s/histNClusters_%d", groupname.Data(), fCentBin);
645  fHistManager.FillTH1(histname, count);
646 
647  histname = TString::Format("%s/histEMCalNClusters_%d", groupname.Data(), fCentBin);
648  fHistManager.FillTH1(histname, countEMCal);
649 
650  histname = TString::Format("%s/histPHOSNClusters_%d", groupname.Data(), fCentBin);
651  fHistManager.FillTH1(histname, countPHOS);
652 
653  }
654 }
655 
661 {
662  if (!fCaloCells) return;
663 
664  TString histname;
665 
666  const Short_t ncells = fCaloCells->GetNumberOfCells();
667 
668  histname = TString::Format("%s/histNCells_%d", fCaloCellsName.Data(), fCentBin);
669  fHistManager.FillTH1(histname, ncells);
670 
671  histname = TString::Format("%s/histCellEnergyvsAbsId_%d", fCaloCellsName.Data(), fCentBin);
672  for (Short_t pos = 0; pos < ncells; pos++) {
673  Double_t amp = fCaloCells->GetAmplitude(pos);
674  Short_t absId = fCaloCells->GetCellNumber(pos);
675 
676  fHistManager.FillTH2(histname, absId, amp);
677  }
678 }
679 
685 {
687 }
688 
697 {
698  TString histname;
699  TString groupname;
700 
701  AliJetContainer* jetCont = 0;
702  TIter next(&fJetCollArray);
703  while ((jetCont = static_cast<AliJetContainer*>(next()))) {
704  groupname = jetCont->GetName();
705 
706  // Get trigger jet
707  AliEmcalJet* trigJet = 0;
708  if (jetCont->GetRhoParameter())
709  trigJet = jetCont->GetLeadingJet("rho");
710  else
711  trigJet = jetCont->GetLeadingJet();
712 
713  if(!trigJet) continue;
714  Double_t trigJetPt = trigJet->Pt() - jetCont->GetRhoVal() * trigJet->Area();
715  Double_t trigJetEta = trigJet->Eta();
716  Double_t trigJetPhi = trigJet->Phi();
717 
718  // Iterate separately for each leading hadron cut
719  Double_t leadingHadCut;
720  for (Int_t k=0; k<2; k++) {
721 
722  // Skip the event if the leading jet doesn't satisfy the leading hadron threshold
723  leadingHadCut = 5*k;
724  if (jetCont->GetLeadingHadronPt(trigJet) < leadingHadCut) continue;
725 
726  // Code for a single set of pT thresolds
727  /*
728  if( trigJetPt < fTrigJetMinPt) continue;
729 
730  // Look for associated jet
731  AliEmcalJet *assJet = 0;
732  for(auto assJetCand : jetCont->accepted()) {
733  if (!assJetCand) continue;
734  Double_t assJetCandPt = assJetCand->Pt() - jetCont->GetRhoVal() * assJetCand->Area();
735  if ( assJetCandPt < fAssJetMinPt ) continue;
736  if ( TMath::Abs(trigJet->Phi() - assJetCand->Phi()) < fDeltaPhiMin ) continue;
737  if (assJet) {
738  Double_t assJetPt = assJet->Pt() - jetCont->GetRhoVal() * assJet->Area();
739  if ( assJetCandPt < assJetPt ) continue;
740  }
741 
742  assJet = assJetCand;
743  }
744  */
745 
746  // Loop through a variety of leading jet pT thresholds
747  Double_t trigJetMinPt;
748  for (Int_t i = 0; i < 4; i++) {
749 
750  // Set the trig jet pT threshold
751  trigJetMinPt = 35 + 5*i;
752  if( trigJetPt < trigJetMinPt) continue;
753 
754  // Loop through a variety of subleading jet pT thresholds
755  Double_t assJetMinPtFrac;
756  Double_t assJetMinPt;
757  for (Int_t j=0; j < 4; j++) {
758  if (j==0) assJetMinPtFrac = 0;
759  if (j==1) assJetMinPtFrac = 0.4;
760  if (j==2) assJetMinPtFrac = 0.5;
761  if (j==3) assJetMinPtFrac = 0.6;
762  Double_t assJetMinPt = trigJetMinPt * assJetMinPtFrac;
763 
764  // histogram label
765  TString label = TString::Format("_%d_had%d_trig%d_ass%d", fCentBin, k, i, j);
766 
767  // Find the subleading jet in the opposite hemisphere
768  AliEmcalJet *assJet = 0;
769  for(auto assJetCand : jetCont->accepted()) {
770  if (!assJetCand) continue;
771  Double_t assJetCandPt = assJetCand->Pt() - jetCont->GetRhoVal() * assJetCand->Area();
772  if ( assJetCandPt < assJetMinPt ) continue;
773  if ( TMath::Abs(trigJet->Phi() - assJetCand->Phi()) < fDeltaPhiMin ) continue;
774  if (assJet) {
775  Double_t assJetPt = assJet->Pt() - jetCont->GetRhoVal() * assJet->Area();
776  if ( assJetCandPt < assJetPt ) continue;
777  }
778 
779  assJet = assJetCand;
780  }
781 
782 
783  // If we find an acceptable associated jet, fill the di-jet histograms
784  if(assJet){
785 
786  // Define kinematic variables for the di-jet pair
787  Double_t assJetPt = assJet->Pt() - jetCont->GetRhoVal() * assJet->Area();
788  Double_t assJetPhi = assJet->Phi();
789  Double_t assJetEta = assJet->Eta();
790 
791  // Some basic plots for the found di-jet pairs
792 
793  histname = TString::Format("%s/histDijetLeadingJetPt", groupname.Data()) += label;
794  fHistManager.FillTH1(histname, trigJetPt);
795 
796  histname = TString::Format("%s/histDijetLeadingJetPtuncorr", groupname.Data()) += label;
797  fHistManager.FillTH1(histname, trigJet->Pt());
798 
799  histname = TString::Format("%s/histDijetLeadingJetPhi", groupname.Data()) += label;
800  fHistManager.FillTH1(histname, trigJetPhi);
801 
802  //histname = TString::Format("%s/histDijetLeadingJetEta_%d", groupname.Data(), fCentBin);
803  //fHistManager.FillTH1(histname, trigJetEta);
804 
805  histname = TString::Format("%s/histDijetSubleadingJetPt", groupname.Data()) += label;
806  fHistManager.FillTH1(histname, assJetPt);
807 
808  histname = TString::Format("%s/histDijetSubleadingJetPhi", groupname.Data()) += label;
809  fHistManager.FillTH1(histname, assJetPhi);
810 
811  //histname = TString::Format("%s/histDijetSubleadingJetEta_%d", groupname.Data(), fCentBin);
812  //fHistManager.FillTH1(histname, assJetEta);
813 
814  //Double_t deltaEta = trigJetEta - assJetEta;
815  //histname = TString::Format("%s/histDijetDeltaEta_%d", groupname.Data(), fCentBin);
816  //fHistManager.FillTH1(histname, deltaEta);
817 
818  // Some di-jet observables
819 
820  Double_t AJ = (trigJetPt - assJetPt)/(trigJetPt + assJetPt);
821  histname = TString::Format("%s/histDijetAJ", groupname.Data()) += label;
822  fHistManager.FillTH1(histname, AJ);
823 
824  Double_t xJ = assJetPt / trigJetPt;
825  histname = TString::Format("%s/histDijetxJ", groupname.Data()) += label;
826  fHistManager.FillTH1(histname, xJ);
827 
828  Double_t deltaPhi = TMath::Abs(trigJetPhi - assJetPhi);
829  histname = TString::Format("%s/histDijetDeltaPhi", groupname.Data()) += label;
830  fHistManager.FillTH1(histname, deltaPhi);
831 
832  //Double_t kT = TMath::Abs( trigJetPt * TMath::Sin(deltaPhi) );
833  //histname = TString::Format("%s/histDijetkT_%d", groupname.Data(), fCentBin);
834  //fHistManager.FillTH1(histname, kT);
835 
836  // Study internal structure of the leading, subleading jets
837 
838  //Double_t nTracksLeadingJet = trigJet->GetNumberOfTracks();
839  //histname = TString::Format("%s/histDijetLeadingJetNTracks_%d", groupname.Data(), fCentBin);
840  //fHistManager.FillTH1(histname, nTracksLeadingJet);
841 
842  //Double_t nTracksSubleadingJet = assJet->GetNumberOfTracks();
843  //histname = TString::Format("%s/histDijetSubleadingJetNTracks_%d", groupname.Data(), fCentBin);
844  //fHistManager.FillTH1(histname, nTracksSubleadingJet);
845 
846  //histname = TString::Format("%s/histDijetLeadingJetArea_%d", groupname.Data(), fCentBin);
847  //fHistManager.FillTH1(histname, trigJet->Area());
848 
849  //histname = TString::Format("%s/histDijetSubleadingJetArea_%d", groupname.Data(), fCentBin);
850  //fHistManager.FillTH1(histname, assJet->Area());
851 
852  // Study momentum balance
853  // TO DO: loop through tracks/clusters in event (as in above loops) using ClusterAt(i)/TrackAt(i)
854  // and group according to proximity to leading, subleading jet axes
855 
856  } else { // If we don't find an associated jet
857 
858  // Plot leading jets that did not find any associated jet
859  histname = TString::Format("%s/histUnmatchedLeadingJetPt", groupname.Data()) += label;
860  fHistManager.FillTH1(histname, trigJetPt);
861 
862  // Find subleading jet in event
863  AliEmcalJet *subleadingJet = 0;
864  for(auto subleadingJetCand : jetCont->accepted()) {
865  if (!subleadingJetCand) continue;
866  Double_t subleadingJetCandPt = subleadingJetCand->Pt() - jetCont->GetRhoVal() * subleadingJetCand->Area();
867  if (subleadingJetCandPt < trigJetPt-0.01) {
868  if (subleadingJet) {
869  Double_t subleadingJetPt = subleadingJet->Pt() - jetCont->GetRhoVal() * subleadingJet->Area();
870  if ( subleadingJetCandPt < subleadingJetPt ) continue;
871  }
872 
873  subleadingJet = subleadingJetCand;
874  }
875  }
876  if(subleadingJet) {
877  histname = TString::Format("%s/histUnmatchedSubleadingJetPt", groupname.Data()) += label;
878  fHistManager.FillTH1(histname, subleadingJet->Pt() - jetCont->GetRhoVal() * subleadingJet->Area());
879  }
880  }
881  }
882  }
883  }
884 
885  }
886 
887  return kTRUE;
888 }
889 
894 {
895 
896 }
Double_t Area() const
Definition: AliEmcalJet.h:117
THashList * CreateHistoGroup(const char *groupname)
TObjArray fClusterCollArray
cluster collection array
Double_t GetRhoVal() const
const TString & GetRhoName() const
double Double_t
Definition: External.C:58
Double_t Eta() const
Definition: AliEmcalJet.h:108
Double_t Phi() const
Definition: AliEmcalJet.h:104
void FillTH2(const char *hname, double x, double y, double weight=1., Option_t *opt="")
Declaration of class AliTLorentzVector.
Double_t fMinBinPt
min pt in histograms
Int_t fCentBin
!event centrality bin
Container for particles within the EMCAL framework.
TObjArray fParticleCollArray
particle/track collection array
const AliClusterIterableContainer all() const
AliEmcalJet * GetLeadingJet(const char *opt="")
TH2 * CreateTH2(const char *name, const char *title, int nbinsx, double xmin, double xmax, int nbinsy, double ymin, double ymax, Option_t *opt="")
int Int_t
Definition: External.C:63
unsigned int UInt_t
Definition: External.C:33
THashList * GetListOfHistograms() const
Definition: THistManager.h:504
Double_t GetLeadingHadronPt(const AliEmcalJet *jet) const
Declaration of class AliAnalysisTaskEmcalDijetImbalance.
Double_t Phi_0_2pi() const
TH1 * CreateTH1(const char *name, const char *title, int nbins, double xmin, double xmax, Option_t *opt="")
BeamType fForceBeamType
forced beam type
Int_t fNcentBins
how many centrality bins
AliClusterContainer * GetClusterContainer(Int_t i=0) const
Double_t fDeltaPhiMin
minimum delta phi between di-jets
AliVCluster * GetAcceptCluster(Int_t i) const
const AliClusterIterableContainer accepted() const
TString fCaloCellsName
name of calo cell collection
void FillTH1(const char *hname, double x, double weight=1., Option_t *opt="")
TObjArray fJetCollArray
jet collection array
short Short_t
Definition: External.C:23
AliVCaloCells * fCaloCells
!cells
AliRhoParameter * GetRhoParameter()
Double_t Pt() const
Definition: AliEmcalJet.h:96
AliEmcalList * fOutput
!output list
Double_t fMaxBinPt
max pt in histograms
ClassImp(AliAnalysisTaskCRC) AliAnalysisTaskCRC
Double_t fVertex[3]
!event vertex
void SetMakeGeneralHistograms(Bool_t g)
Base task in the EMCAL jet framework.
Implementation of a sample jet analysis task.
Represent a jet reconstructed using the EMCal jet framework.
Definition: AliEmcalJet.h:44
const AliParticleIterableContainer accepted() const
const char Option_t
Definition: External.C:48
const AliJetIterableContainer accepted() const
bool Bool_t
Definition: External.C:53
Container structure for EMCAL clusters.
Container for jet within the EMCAL jet framework.
Int_t fNbins
no. of pt bins