AliPhysics  781d0c7 (781d0c7)
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
AliJetResponseMaker.cxx
Go to the documentation of this file.
1 /*************************************************************************
2  * Copyright(c) 1998-2009, 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 // Emcal jet response matrix maker task.
17 //
18 // Author: Salvatore Aiola, Yale University, salvatore.aiola@cern.ch
19 
20 #include "AliJetResponseMaker.h"
21 
22 #include <TClonesArray.h>
23 #include <TH2F.h>
24 #include <THnSparse.h>
25 
26 #include "AliTLorentzVector.h"
27 #include "AliAnalysisManager.h"
28 #include "AliVCluster.h"
29 #include "AliVTrack.h"
30 #include "AliEmcalJet.h"
31 #include "AliLog.h"
32 #include "AliRhoParameter.h"
33 #include "AliNamedArrayI.h"
34 #include "AliJetContainer.h"
35 #include "AliParticleContainer.h"
36 #include "AliClusterContainer.h"
37 
39 
40 //________________________________________________________________________
43  fMatching(kNoMatching),
44  fMatchingPar1(0),
45  fMatchingPar2(0),
46  fUseCellsToMatch(kFALSE),
47  fMinJetMCPt(1),
48  fHistoType(0),
49  fDeltaPtAxis(0),
50  fDeltaEtaDeltaPhiAxis(0),
51  fNEFAxis(0),
52  fZAxis(0),
53  fFlavourZAxis(0),
54  fFlavourPtAxis(0),
55  fIsJet1Rho(kFALSE),
56  fIsJet2Rho(kFALSE),
57  fHistRejectionReason1(0),
58  fHistRejectionReason2(0),
59  fHistJets1(0),
60  fHistJets2(0),
61  fHistMatching(0),
62  fHistJets1PhiEta(0),
63  fHistJets1PtArea(0),
64  fHistJets1CorrPtArea(0),
65  fHistJets1NEFvsPt(0),
66  fHistJets1CEFvsCEFPt(0),
67  fHistJets1ZvsPt(0),
68  fHistJets2PhiEta(0),
69  fHistJets2PtArea(0),
70  fHistJets2CorrPtArea(0),
71  fHistJets2NEFvsPt(0),
72  fHistJets2CEFvsCEFPt(0),
73  fHistJets2ZvsPt(0),
74  fHistCommonEnergy1vsJet1Pt(0),
75  fHistCommonEnergy2vsJet2Pt(0),
76  fHistDistancevsJet1Pt(0),
77  fHistDistancevsJet2Pt(0),
78  fHistDistancevsCommonEnergy1(0),
79  fHistDistancevsCommonEnergy2(0),
80  fHistCommonEnergy1vsCommonEnergy2(0),
81  fHistDeltaEtaDeltaPhi(0),
82  fHistJet2PtOverJet1PtvsJet2Pt(0),
83  fHistJet1PtOverJet2PtvsJet1Pt(0),
84  fHistDeltaPtvsJet1Pt(0),
85  fHistDeltaPtvsJet2Pt(0),
86  fHistDeltaPtOverJet1PtvsJet1Pt(0),
87  fHistDeltaPtOverJet2PtvsJet2Pt(0),
88  fHistDeltaPtvsDistance(0),
89  fHistDeltaPtvsCommonEnergy1(0),
90  fHistDeltaPtvsCommonEnergy2(0),
91  fHistDeltaPtvsArea1(0),
92  fHistDeltaPtvsArea2(0),
93  fHistDeltaPtvsDeltaArea(0),
94  fHistJet1PtvsJet2Pt(0),
95  fHistDeltaCorrPtOverJet1CorrPtvsJet1CorrPt(0),
96  fHistDeltaCorrPtOverJet2CorrPtvsJet2CorrPt(0),
97  fHistDeltaCorrPtvsJet1CorrPt(0),
98  fHistDeltaCorrPtvsJet2CorrPt(0),
99  fHistDeltaCorrPtvsDistance(0),
100  fHistDeltaCorrPtvsCommonEnergy1(0),
101  fHistDeltaCorrPtvsCommonEnergy2(0),
102  fHistDeltaCorrPtvsArea1(0),
103  fHistDeltaCorrPtvsArea2(0),
104  fHistDeltaCorrPtvsDeltaArea(0),
105  fHistJet1CorrPtvsJet2CorrPt(0),
106  fHistDeltaMCPtOverJet1MCPtvsJet1MCPt(0),
107  fHistDeltaMCPtOverJet2PtvsJet2Pt(0),
108  fHistDeltaMCPtvsJet1MCPt(0),
109  fHistDeltaMCPtvsJet2Pt(0),
110  fHistDeltaMCPtvsDistance(0),
111  fHistDeltaMCPtvsCommonEnergy1(0),
112  fHistDeltaMCPtvsCommonEnergy2(0),
113  fHistDeltaMCPtvsArea1(0),
114  fHistDeltaMCPtvsArea2(0),
115  fHistDeltaMCPtvsDeltaArea(0),
116  fHistJet1MCPtvsJet2Pt(0)
117 {
118  // Default constructor.
119 
120  SetMakeGeneralHistograms(kTRUE);
121 }
122 
123 //________________________________________________________________________
125  AliAnalysisTaskEmcalJet(name, kTRUE),
126  fMatching(kNoMatching),
127  fMatchingPar1(0),
128  fMatchingPar2(0),
129  fUseCellsToMatch(kFALSE),
130  fMinJetMCPt(1),
131  fHistoType(0),
132  fDeltaPtAxis(0),
133  fDeltaEtaDeltaPhiAxis(0),
134  fNEFAxis(0),
135  fZAxis(0),
136  fFlavourZAxis(0),
137  fFlavourPtAxis(0),
138  fIsJet1Rho(kFALSE),
139  fIsJet2Rho(kFALSE),
140  fHistRejectionReason1(0),
141  fHistRejectionReason2(0),
142  fHistJets1(0),
143  fHistJets2(0),
144  fHistMatching(0),
145  fHistJets1PhiEta(0),
146  fHistJets1PtArea(0),
147  fHistJets1CorrPtArea(0),
148  fHistJets1NEFvsPt(0),
149  fHistJets1CEFvsCEFPt(0),
150  fHistJets1ZvsPt(0),
151  fHistJets2PhiEta(0),
152  fHistJets2PtArea(0),
153  fHistJets2CorrPtArea(0),
154  fHistJets2NEFvsPt(0),
155  fHistJets2CEFvsCEFPt(0),
156  fHistJets2ZvsPt(0),
157  fHistCommonEnergy1vsJet1Pt(0),
158  fHistCommonEnergy2vsJet2Pt(0),
159  fHistDistancevsJet1Pt(0),
160  fHistDistancevsJet2Pt(0),
161  fHistDistancevsCommonEnergy1(0),
162  fHistDistancevsCommonEnergy2(0),
163  fHistCommonEnergy1vsCommonEnergy2(0),
164  fHistDeltaEtaDeltaPhi(0),
165  fHistJet2PtOverJet1PtvsJet2Pt(0),
166  fHistJet1PtOverJet2PtvsJet1Pt(0),
167  fHistDeltaPtvsJet1Pt(0),
168  fHistDeltaPtvsJet2Pt(0),
169  fHistDeltaPtOverJet1PtvsJet1Pt(0),
170  fHistDeltaPtOverJet2PtvsJet2Pt(0),
171  fHistDeltaPtvsDistance(0),
172  fHistDeltaPtvsCommonEnergy1(0),
173  fHistDeltaPtvsCommonEnergy2(0),
174  fHistDeltaPtvsArea1(0),
175  fHistDeltaPtvsArea2(0),
176  fHistDeltaPtvsDeltaArea(0),
177  fHistJet1PtvsJet2Pt(0),
178  fHistDeltaCorrPtOverJet1CorrPtvsJet1CorrPt(0),
179  fHistDeltaCorrPtOverJet2CorrPtvsJet2CorrPt(0),
180  fHistDeltaCorrPtvsJet1CorrPt(0),
181  fHistDeltaCorrPtvsJet2CorrPt(0),
182  fHistDeltaCorrPtvsDistance(0),
183  fHistDeltaCorrPtvsCommonEnergy1(0),
184  fHistDeltaCorrPtvsCommonEnergy2(0),
185  fHistDeltaCorrPtvsArea1(0),
186  fHistDeltaCorrPtvsArea2(0),
187  fHistDeltaCorrPtvsDeltaArea(0),
188  fHistJet1CorrPtvsJet2CorrPt(0),
189  fHistDeltaMCPtOverJet1MCPtvsJet1MCPt(0),
190  fHistDeltaMCPtOverJet2PtvsJet2Pt(0),
191  fHistDeltaMCPtvsJet1MCPt(0),
192  fHistDeltaMCPtvsJet2Pt(0),
193  fHistDeltaMCPtvsDistance(0),
194  fHistDeltaMCPtvsCommonEnergy1(0),
195  fHistDeltaMCPtvsCommonEnergy2(0),
196  fHistDeltaMCPtvsArea1(0),
197  fHistDeltaMCPtvsArea2(0),
198  fHistDeltaMCPtvsDeltaArea(0),
199  fHistJet1MCPtvsJet2Pt(0)
200 {
201  // Standard constructor.
202 
204 }
205 
206 //________________________________________________________________________
208 {
209  // Destructor
210 }
211 
212 
213 //________________________________________________________________________
215 {
216  // Allocate TH2 histograms.
217 
218  fHistJets1PhiEta = new TH2F("fHistJets1PhiEta", "fHistJets1PhiEta", 40, -1, 1, 40, 0, TMath::Pi()*2);
219  fHistJets1PhiEta->GetXaxis()->SetTitle("#eta");
220  fHistJets1PhiEta->GetYaxis()->SetTitle("#phi");
222 
223  fHistJets1PtArea = new TH2F("fHistJets1PtArea", "fHistJets1PtArea", fNbins/2, 0, 1.5, fNbins, fMinBinPt, fMaxBinPt);
224  fHistJets1PtArea->GetXaxis()->SetTitle("area");
225  fHistJets1PtArea->GetYaxis()->SetTitle("p_{T,1} (GeV/c)");
226  fHistJets1PtArea->GetZaxis()->SetTitle("counts");
228 
229  if (fIsJet1Rho) {
230  fHistJets1CorrPtArea = new TH2F("fHistJets1CorrPtArea", "fHistJets1CorrPtArea",
231  fNbins/2, 0, 1.5, 2*fNbins, -fMaxBinPt, fMaxBinPt);
232  fHistJets1CorrPtArea->GetXaxis()->SetTitle("area");
233  fHistJets1CorrPtArea->GetYaxis()->SetTitle("p_{T,1}^{corr} (GeV/c)");
234  fHistJets1CorrPtArea->GetZaxis()->SetTitle("counts");
236  }
237 
238  fHistJets1ZvsPt = new TH2F("fHistJets1ZvsPt", "fHistJets1ZvsPt", 120, 0, 1.2, fNbins, fMinBinPt, fMaxBinPt);
239  fHistJets1ZvsPt->GetXaxis()->SetTitle("Z");
240  fHistJets1ZvsPt->GetYaxis()->SetTitle("p_{T,1} (GeV/c)");
241  fHistJets1ZvsPt->GetZaxis()->SetTitle("counts");
242  fOutput->Add(fHistJets1ZvsPt);
243 
244  fHistJets1NEFvsPt = new TH2F("fHistJets1NEFvsPt", "fHistJets1NEFvsPt", 120, 0, 1.2, fNbins, fMinBinPt, fMaxBinPt);
245  fHistJets1NEFvsPt->GetXaxis()->SetTitle("NEF");
246  fHistJets1NEFvsPt->GetYaxis()->SetTitle("p_{T,1} (GeV/c)");
247  fHistJets1NEFvsPt->GetZaxis()->SetTitle("counts");
249 
250  fHistJets1CEFvsCEFPt = new TH2F("fHistJets1CEFvsCEFPt", "fHistJets1CEFvsCEFPt", 120, 0, 1.2, fNbins, fMinBinPt, fMaxBinPt);
251  fHistJets1CEFvsCEFPt->GetXaxis()->SetTitle("1-NEF");
252  fHistJets1CEFvsCEFPt->GetYaxis()->SetTitle("(1-NEF)*p_{T,1} (GeV/c)");
253  fHistJets1CEFvsCEFPt->GetZaxis()->SetTitle("counts");
255 
256  // Jets 2 histograms
257 
258  fHistJets2PhiEta = new TH2F("fHistJets2PhiEta", "fHistJets2PhiEta", 40, -1, 1, 40, 0, TMath::Pi()*2);
259  fHistJets2PhiEta->GetXaxis()->SetTitle("#eta");
260  fHistJets2PhiEta->GetYaxis()->SetTitle("#phi");
262 
263  fHistJets2PtArea = new TH2F("fHistJets2PtArea", "fHistJets2PtArea", fNbins/2, 0, 1.5, fNbins, fMinBinPt, fMaxBinPt);
264  fHistJets2PtArea->GetXaxis()->SetTitle("area");
265  fHistJets2PtArea->GetYaxis()->SetTitle("p_{T,2} (GeV/c)");
266  fHistJets2PtArea->GetZaxis()->SetTitle("counts");
268 
269  if (fIsJet2Rho) {
270  fHistJets2CorrPtArea = new TH2F("fHistJets2CorrPtArea", "fHistJets2CorrPtArea",
271  fNbins/2, 0, 1.5, 2*fNbins, -fMaxBinPt, fMaxBinPt);
272  fHistJets2CorrPtArea->GetXaxis()->SetTitle("area");
273  fHistJets2CorrPtArea->GetYaxis()->SetTitle("p_{T,2}^{corr} (GeV/c)");
274  fHistJets2CorrPtArea->GetZaxis()->SetTitle("counts");
276  }
277 
278  fHistJets2ZvsPt = new TH2F("fHistJets2ZvsPt", "fHistJets2ZvsPt", 120, 0, 1.2, fNbins, fMinBinPt, fMaxBinPt);
279  fHistJets2ZvsPt->GetXaxis()->SetTitle("Z");
280  fHistJets2ZvsPt->GetYaxis()->SetTitle("p_{T,2} (GeV/c)");
281  fHistJets2ZvsPt->GetZaxis()->SetTitle("counts");
282  fOutput->Add(fHistJets2ZvsPt);
283 
284  fHistJets2NEFvsPt = new TH2F("fHistJets2NEFvsPt", "fHistJets2NEFvsPt", 120, 0, 1.2, fNbins, fMinBinPt, fMaxBinPt);
285  fHistJets2NEFvsPt->GetXaxis()->SetTitle("NEF");
286  fHistJets2NEFvsPt->GetYaxis()->SetTitle("p_{T,2} (GeV/c)");
287  fHistJets2NEFvsPt->GetZaxis()->SetTitle("counts");
289 
290  fHistJets2CEFvsCEFPt = new TH2F("fHistJets2CEFvsCEFPt", "fHistJets2CEFvsCEFPt", 120, 0, 1.2, fNbins, fMinBinPt, fMaxBinPt);
291  fHistJets2CEFvsCEFPt->GetXaxis()->SetTitle("1-NEF");
292  fHistJets2CEFvsCEFPt->GetYaxis()->SetTitle("(1-NEF)*p_{T,2} (GeV/c)");
293  fHistJets2CEFvsCEFPt->GetZaxis()->SetTitle("counts");
295 
296  // Matching histograms
297 
298  fHistCommonEnergy1vsJet1Pt = new TH2F("fHistCommonEnergy1vsJet1Pt", "fHistCommonEnergy1vsJet1Pt", fNbins/2, 0, 1.2, fNbins/2, fMinBinPt, fMaxBinPt);
299  fHistCommonEnergy1vsJet1Pt->GetXaxis()->SetTitle("Common energy 1 (%)");
300  fHistCommonEnergy1vsJet1Pt->GetYaxis()->SetTitle("p_{T,1}");
301  fHistCommonEnergy1vsJet1Pt->GetZaxis()->SetTitle("counts");
303 
304  fHistCommonEnergy2vsJet2Pt = new TH2F("fHistCommonEnergy2vsJet2Pt", "fHistCommonEnergy2vsJet2Pt", fNbins/2, 0, 1.2, fNbins/2, fMinBinPt, fMaxBinPt);
305  fHistCommonEnergy2vsJet2Pt->GetXaxis()->SetTitle("Common energy 2 (%)");
306  fHistCommonEnergy2vsJet2Pt->GetYaxis()->SetTitle("p_{T,2}");
307  fHistCommonEnergy2vsJet2Pt->GetZaxis()->SetTitle("counts");
309 
310  fHistDistancevsJet1Pt = new TH2F("fHistDistancevsJet1Pt", "fHistDistancevsJet1Pt", fNbins/2, 0, 1.2, fNbins/2, fMinBinPt, fMaxBinPt);
311  fHistDistancevsJet1Pt->GetXaxis()->SetTitle("Distance");
312  fHistDistancevsJet1Pt->GetYaxis()->SetTitle("p_{T,1}");
313  fHistDistancevsJet1Pt->GetZaxis()->SetTitle("counts");
315 
316  fHistDistancevsJet2Pt = new TH2F("fHistDistancevsJet2Pt", "fHistDistancevsJet2Pt", fNbins/2, 0, 1.2, fNbins/2, fMinBinPt, fMaxBinPt);
317  fHistDistancevsJet2Pt->GetXaxis()->SetTitle("Distance");
318  fHistDistancevsJet2Pt->GetYaxis()->SetTitle("p_{T,2}");
319  fHistDistancevsJet2Pt->GetZaxis()->SetTitle("counts");
321 
322  fHistDistancevsCommonEnergy1 = new TH2F("fHistDistancevsCommonEnergy1", "fHistDistancevsCommonEnergy1", fNbins/2, 0, 1.2, fNbins/2, 0, 1.2);
323  fHistDistancevsCommonEnergy1->GetXaxis()->SetTitle("Distance");
324  fHistDistancevsCommonEnergy1->GetYaxis()->SetTitle("Common energy 1 (%)");
325  fHistDistancevsCommonEnergy1->GetZaxis()->SetTitle("counts");
327 
328  fHistDistancevsCommonEnergy2 = new TH2F("fHistDistancevsCommonEnergy2", "fHistDistancevsCommonEnergy2", fNbins/2, 0, 1.2, fNbins/2, 0, 1.2);
329  fHistDistancevsCommonEnergy2->GetXaxis()->SetTitle("Distance");
330  fHistDistancevsCommonEnergy2->GetYaxis()->SetTitle("Common energy 2 (%)");
331  fHistDistancevsCommonEnergy2->GetZaxis()->SetTitle("counts");
333 
334  fHistCommonEnergy1vsCommonEnergy2 = new TH2F("fHistCommonEnergy1vsCommonEnergy2", "fHistCommonEnergy1vsCommonEnergy2", fNbins/2, 0, 1.2, fNbins/2, 0, 1.2);
335  fHistCommonEnergy1vsCommonEnergy2->GetXaxis()->SetTitle("Common energy 1 (%)");
336  fHistCommonEnergy1vsCommonEnergy2->GetYaxis()->SetTitle("Common energy 2 (%)");
337  fHistCommonEnergy1vsCommonEnergy2->GetZaxis()->SetTitle("counts");
339 
340  fHistDeltaEtaDeltaPhi = new TH2F("fHistDeltaEtaDeltaPhi", "fHistDeltaEtaDeltaPhi", fNbins/4, -1, 1, fNbins/4, -TMath::Pi()/2, TMath::Pi()*3/2);
341  fHistDeltaEtaDeltaPhi->GetXaxis()->SetTitle("Common energy 1 (%)");
342  fHistDeltaEtaDeltaPhi->GetYaxis()->SetTitle("Common energy 2 (%)");
343  fHistDeltaEtaDeltaPhi->GetZaxis()->SetTitle("counts");
345 
346  fHistJet2PtOverJet1PtvsJet2Pt = new TH2F("fHistJet2PtOverJet1PtvsJet2Pt", "fHistJet2PtOverJet1PtvsJet2Pt", fNbins, fMinBinPt, fMaxBinPt, 300, 0, 1.5);
347  fHistJet2PtOverJet1PtvsJet2Pt->GetXaxis()->SetTitle("p_{T,2}");
348  fHistJet2PtOverJet1PtvsJet2Pt->GetYaxis()->SetTitle("p_{T,2} / p_{T,1}");
349  fHistJet2PtOverJet1PtvsJet2Pt->GetZaxis()->SetTitle("counts");
351 
352  fHistJet1PtOverJet2PtvsJet1Pt = new TH2F("fHistJet1PtOverJet2PtvsJet1Pt", "fHistJet1PtOverJet2PtvsJet1Pt", fNbins, fMinBinPt, fMaxBinPt, 300, 0, 1.5);
353  fHistJet1PtOverJet2PtvsJet1Pt->GetXaxis()->SetTitle("p_{T,1}");
354  fHistJet1PtOverJet2PtvsJet1Pt->GetYaxis()->SetTitle("p_{T,1} / p_{T,2}");
355  fHistJet1PtOverJet2PtvsJet1Pt->GetZaxis()->SetTitle("counts");
357 
358  fHistDeltaPtvsJet1Pt = new TH2F("fHistDeltaPtvsJet1Pt", "fHistDeltaPtvsJet1Pt",
360  fHistDeltaPtvsJet1Pt->GetXaxis()->SetTitle("p_{T,1}");
361  fHistDeltaPtvsJet1Pt->GetYaxis()->SetTitle("#deltap_{T} (GeV/c)");
362  fHistDeltaPtvsJet1Pt->GetZaxis()->SetTitle("counts");
364 
365  fHistDeltaPtvsJet2Pt = new TH2F("fHistDeltaPtvsJet2Pt", "fHistDeltaPtvsJet2Pt",
367  fHistDeltaPtvsJet2Pt->GetXaxis()->SetTitle("p_{T,2}");
368  fHistDeltaPtvsJet2Pt->GetYaxis()->SetTitle("#deltap_{T} (GeV/c)");
369  fHistDeltaPtvsJet2Pt->GetZaxis()->SetTitle("counts");
371 
372  fHistDeltaPtOverJet1PtvsJet1Pt = new TH2F("fHistDeltaPtOverJet1PtvsJet1Pt", "fHistDeltaPtOverJet1PtvsJet1Pt",
373  fNbins, fMinBinPt, fMaxBinPt, fNbins, -5, 5);
374  fHistDeltaPtOverJet1PtvsJet1Pt->GetXaxis()->SetTitle("p_{T,1}");
375  fHistDeltaPtOverJet1PtvsJet1Pt->GetYaxis()->SetTitle("#deltap_{T} / p_{T,1}");
376  fHistDeltaPtOverJet1PtvsJet1Pt->GetZaxis()->SetTitle("counts");
378 
379  fHistDeltaPtOverJet2PtvsJet2Pt = new TH2F("fHistDeltaPtOverJet2PtvsJet2Pt", "fHistDeltaPtOverJet2PtvsJet2Pt",
380  fNbins, fMinBinPt, fMaxBinPt, fNbins, -5, 5);
381  fHistDeltaPtOverJet2PtvsJet2Pt->GetXaxis()->SetTitle("p_{T,2}");
382  fHistDeltaPtOverJet2PtvsJet2Pt->GetYaxis()->SetTitle("#deltap_{T} / p_{T,2}");
383  fHistDeltaPtOverJet2PtvsJet2Pt->GetZaxis()->SetTitle("counts");
385 
386  fHistDeltaPtvsDistance = new TH2F("fHistDeltaPtvsDistance", "fHistDeltaPtvsDistance",
387  fNbins/2, 0, 1.2, 2*fNbins, -fMaxBinPt, fMaxBinPt);
388  fHistDeltaPtvsDistance->GetXaxis()->SetTitle("Distance");
389  fHistDeltaPtvsDistance->GetYaxis()->SetTitle("#deltap_{T} (GeV/c)");
390  fHistDeltaPtvsDistance->GetZaxis()->SetTitle("counts");
392 
393  fHistDeltaPtvsCommonEnergy1 = new TH2F("fHistDeltaPtvsCommonEnergy1", "fHistDeltaPtvsCommonEnergy1",
394  fNbins/2, 0, 1.2, 2*fNbins, -fMaxBinPt, fMaxBinPt);
395  fHistDeltaPtvsCommonEnergy1->GetXaxis()->SetTitle("Common energy 1 (%)");
396  fHistDeltaPtvsCommonEnergy1->GetYaxis()->SetTitle("#deltap_{T} (GeV/c)");
397  fHistDeltaPtvsCommonEnergy1->GetZaxis()->SetTitle("counts");
399 
400  fHistDeltaPtvsCommonEnergy2 = new TH2F("fHistDeltaPtvsCommonEnergy2", "fHistDeltaPtvsCommonEnergy2",
401  fNbins/2, 0, 1.2, 2*fNbins, -fMaxBinPt, fMaxBinPt);
402  fHistDeltaPtvsCommonEnergy2->GetXaxis()->SetTitle("Common energy 2 (%)");
403  fHistDeltaPtvsCommonEnergy2->GetYaxis()->SetTitle("#deltap_{T} (GeV/c)");
404  fHistDeltaPtvsCommonEnergy2->GetZaxis()->SetTitle("counts");
406 
407  fHistDeltaPtvsArea1 = new TH2F("fHistDeltaPtvsArea1", "fHistDeltaPtvsArea1",
408  fNbins/2, 0, 1.5, 2*fNbins, -fMaxBinPt, fMaxBinPt);
409  fHistDeltaPtvsArea1->GetXaxis()->SetTitle("A_{jet,1}");
410  fHistDeltaPtvsArea1->GetYaxis()->SetTitle("#deltap_{T} (GeV/c)");
411  fHistDeltaPtvsArea1->GetZaxis()->SetTitle("counts");
413 
414  fHistDeltaPtvsArea2 = new TH2F("fHistDeltaPtvsArea2", "fHistDeltaPtvsArea2",
415  fNbins/2, 0, 1.5, 2*fNbins, -fMaxBinPt, fMaxBinPt);
416  fHistDeltaPtvsArea2->GetXaxis()->SetTitle("A_{jet,2}");
417  fHistDeltaPtvsArea2->GetYaxis()->SetTitle("#deltap_{T} (GeV/c)");
418  fHistDeltaPtvsArea2->GetZaxis()->SetTitle("counts");
420 
421  fHistDeltaPtvsDeltaArea = new TH2F("fHistDeltaPtvsDeltaArea", "fHistDeltaPtvsDeltaArea",
422  fNbins, -1.+1./fNbins, 1.+1./fNbins, 2*fNbins, -fMaxBinPt, fMaxBinPt);
423  fHistDeltaPtvsDeltaArea->GetXaxis()->SetTitle("#deltaA_{jet}");
424  fHistDeltaPtvsDeltaArea->GetYaxis()->SetTitle("#deltap_{T} (GeV/c)");
425  fHistDeltaPtvsDeltaArea->GetZaxis()->SetTitle("counts");
427 
428  fHistJet1PtvsJet2Pt = new TH2F("fHistJet1PtvsJet2Pt", "fHistJet1PtvsJet2Pt", fNbins, fMinBinPt, fMaxBinPt, fNbins, fMinBinPt, fMaxBinPt);
429  fHistJet1PtvsJet2Pt->GetXaxis()->SetTitle("p_{T,1}");
430  fHistJet1PtvsJet2Pt->GetYaxis()->SetTitle("p_{T,2}");
431  fHistJet1PtvsJet2Pt->GetZaxis()->SetTitle("counts");
433 
434  if (fIsJet1Rho || fIsJet2Rho) {
435  if (!fIsJet1Rho)
436  fHistDeltaCorrPtvsJet1CorrPt = new TH2F("fHistDeltaCorrPtvsJet1CorrPt", "fHistDeltaCorrPtvsJet1CorrPt",
438  else
439  fHistDeltaCorrPtvsJet1CorrPt = new TH2F("fHistDeltaCorrPtvsJet1CorrPt", "fHistDeltaCorrPtvsJet1CorrPt",
441 
442  fHistDeltaCorrPtvsJet1CorrPt->GetXaxis()->SetTitle("p_{T,1}^{corr}");
443  fHistDeltaCorrPtvsJet1CorrPt->GetYaxis()->SetTitle("#deltap_{T}^{corr} (GeV/c)");
444  fHistDeltaCorrPtvsJet1CorrPt->GetZaxis()->SetTitle("counts");
446 
447  if (!fIsJet2Rho)
448  fHistDeltaCorrPtvsJet2CorrPt = new TH2F("fHistDeltaCorrPtvsJet2CorrPt", "fHistDeltaCorrPtvsJet2CorrPt",
450  else
451  fHistDeltaCorrPtvsJet2CorrPt = new TH2F("fHistDeltaCorrPtvsJet2CorrPt", "fHistDeltaCorrPtvsJet2CorrPt",
453 
454  fHistDeltaCorrPtvsJet2CorrPt->GetXaxis()->SetTitle("p_{T,2}^{corr}");
455  fHistDeltaCorrPtvsJet2CorrPt->GetYaxis()->SetTitle("#deltap_{T}^{corr} (GeV/c)");
456  fHistDeltaCorrPtvsJet2CorrPt->GetZaxis()->SetTitle("counts");
458 
459  fHistDeltaCorrPtOverJet1CorrPtvsJet1CorrPt = new TH2F("fHistDeltaCorrPtOverJet1CorrPtvsJet1CorrPt", "fHistDeltaCorrPtOverJet1CorrPtvsJet1CorrPt",
460  2*fNbins, -fMaxBinPt, fMaxBinPt, fNbins, -5, 5);
461  fHistDeltaCorrPtOverJet1CorrPtvsJet1CorrPt->GetXaxis()->SetTitle("p_{T,1}^{corr}");
462  fHistDeltaCorrPtOverJet1CorrPtvsJet1CorrPt->GetYaxis()->SetTitle("#deltap_{T}^{corr} / p_{T,1}^{corr}");
463  fHistDeltaCorrPtOverJet1CorrPtvsJet1CorrPt->GetZaxis()->SetTitle("counts");
465 
466  fHistDeltaCorrPtOverJet2CorrPtvsJet2CorrPt = new TH2F("fHistDeltaCorrPtOverJet2CorrPtvsJet2CorrPt", "fHistDeltaCorrPtOverJet2CorrPtvsJet2CorrPt",
467  2*fNbins, -fMaxBinPt, fMaxBinPt, fNbins, -5, 5);
468  fHistDeltaCorrPtOverJet2CorrPtvsJet2CorrPt->GetXaxis()->SetTitle("p_{T,2}^{corr}");
469  fHistDeltaCorrPtOverJet2CorrPtvsJet2CorrPt->GetYaxis()->SetTitle("#deltap_{T}^{corr} / p_{T,2}^{corr}");
470  fHistDeltaCorrPtOverJet2CorrPtvsJet2CorrPt->GetZaxis()->SetTitle("counts");
472 
473  fHistDeltaCorrPtvsDistance = new TH2F("fHistDeltaCorrPtvsDistance", "fHistDeltaCorrPtvsDistance",
474  fNbins/2, 0, 1.2, 2*fNbins, -fMaxBinPt, fMaxBinPt);
475  fHistDeltaCorrPtvsDistance->GetXaxis()->SetTitle("Distance");
476  fHistDeltaCorrPtvsDistance->GetYaxis()->SetTitle("#deltap_{T}^{corr} (GeV/c)");
477  fHistDeltaCorrPtvsDistance->GetZaxis()->SetTitle("counts");
479 
480  fHistDeltaCorrPtvsCommonEnergy1 = new TH2F("fHistDeltaCorrPtvsCommonEnergy1", "fHistDeltaCorrPtvsCommonEnergy1",
481  fNbins/2, 0, 1.2, 2*fNbins, -fMaxBinPt, fMaxBinPt);
482  fHistDeltaCorrPtvsCommonEnergy1->GetXaxis()->SetTitle("Common energy 1 (%)");
483  fHistDeltaCorrPtvsCommonEnergy1->GetYaxis()->SetTitle("#deltap_{T}^{corr} (GeV/c)");
484  fHistDeltaCorrPtvsCommonEnergy1->GetZaxis()->SetTitle("counts");
486 
487  fHistDeltaCorrPtvsCommonEnergy2 = new TH2F("fHistDeltaCorrPtvsCommonEnergy2", "fHistDeltaCorrPtvsCommonEnergy2",
488  fNbins/2, 0, 1.2, 2*fNbins, -fMaxBinPt, fMaxBinPt);
489  fHistDeltaCorrPtvsCommonEnergy2->GetXaxis()->SetTitle("Common energy 2 (%)");
490  fHistDeltaCorrPtvsCommonEnergy2->GetYaxis()->SetTitle("#deltap_{T}^{corr} (GeV/c)");
491  fHistDeltaCorrPtvsCommonEnergy2->GetZaxis()->SetTitle("counts");
493 
494  fHistDeltaCorrPtvsArea1 = new TH2F("fHistDeltaCorrPtvsArea1", "fHistDeltaCorrPtvsArea1",
495  fNbins/2, 0, 1.5, 2*fNbins, -fMaxBinPt, fMaxBinPt);
496  fHistDeltaCorrPtvsArea1->GetXaxis()->SetTitle("A_{jet,1}");
497  fHistDeltaCorrPtvsArea1->GetYaxis()->SetTitle("#deltap_{T}^{corr} (GeV/c)");
498  fHistDeltaCorrPtvsArea1->GetZaxis()->SetTitle("counts");
500 
501  fHistDeltaCorrPtvsArea2 = new TH2F("fHistDeltaCorrPtvsArea2", "fHistDeltaCorrPtvsArea2",
502  fNbins/2, 0, 1.5, 2*fNbins, -fMaxBinPt, fMaxBinPt);
503  fHistDeltaCorrPtvsArea2->GetXaxis()->SetTitle("A_{jet,2}");
504  fHistDeltaCorrPtvsArea2->GetYaxis()->SetTitle("#deltap_{T}^{corr} (GeV/c)");
505  fHistDeltaCorrPtvsArea2->GetZaxis()->SetTitle("counts");
507 
508  fHistDeltaCorrPtvsDeltaArea = new TH2F("fHistDeltaCorrPtvsDeltaArea", "fHistDeltaCorrPtvsDeltaArea",
509  fNbins, -1.+1./fNbins, 1.+1./fNbins, 2*fNbins, -fMaxBinPt, fMaxBinPt);
510  fHistDeltaCorrPtvsDeltaArea->GetXaxis()->SetTitle("#deltaA_{jet}");
511  fHistDeltaCorrPtvsDeltaArea->GetYaxis()->SetTitle("#deltap_{T}^{corr} (GeV/c)");
512  fHistDeltaCorrPtvsDeltaArea->GetZaxis()->SetTitle("counts");
514 
515  if (!fIsJet1Rho)
516  fHistJet1CorrPtvsJet2CorrPt = new TH2F("fHistJet1CorrPtvsJet2CorrPt", "fHistJet1CorrPtvsJet2CorrPt",
518  else if (!fIsJet2Rho)
519  fHistJet1CorrPtvsJet2CorrPt = new TH2F("fHistJet1CorrPtvsJet2CorrPt", "fHistJet1CorrPtvsJet2CorrPt",
521  else
522  fHistJet1CorrPtvsJet2CorrPt = new TH2F("fHistJet1CorrPtvsJet2CorrPt", "fHistJet1CorrPtvsJet2CorrPt",
524  2*fNbins, -fMaxBinPt, fMaxBinPt);
525  fHistJet1CorrPtvsJet2CorrPt->GetXaxis()->SetTitle("p_{T,1}^{corr}");
526  fHistJet1CorrPtvsJet2CorrPt->GetYaxis()->SetTitle("p_{T,2}^{corr}");
527  fHistJet1CorrPtvsJet2CorrPt->GetZaxis()->SetTitle("counts");
529  }
530 
531  if (fIsEmbedded) {
532  fHistDeltaMCPtvsJet1MCPt = new TH2F("fHistDeltaMCPtvsJet1MCPt", "fHistDeltaMCPtvsJet1MCPt",
534  fHistDeltaMCPtvsJet1MCPt->GetXaxis()->SetTitle("p_{T,1}^{MC}");
535  fHistDeltaMCPtvsJet1MCPt->GetYaxis()->SetTitle("#deltap_{T} (GeV/c)");
536  fHistDeltaMCPtvsJet1MCPt->GetZaxis()->SetTitle("counts");
538 
539  fHistDeltaMCPtvsJet2Pt = new TH2F("fHistDeltaMCPtvsJet2Pt", "fHistDeltaMCPtvsJet2Pt",
541  fHistDeltaMCPtvsJet2Pt->GetXaxis()->SetTitle("p_{T,2}");
542  fHistDeltaMCPtvsJet2Pt->GetYaxis()->SetTitle("#deltap_{T} (GeV/c)");
543  fHistDeltaMCPtvsJet2Pt->GetZaxis()->SetTitle("counts");
545 
546  fHistDeltaMCPtOverJet1MCPtvsJet1MCPt = new TH2F("fHistDeltaMCPtOverJet1MCPtvsJet1MCPt", "fHistDeltaMCPtOverJet1MCPtvsJet1MCPt",
547  fNbins, fMinBinPt, fMaxBinPt, fNbins, -5, 5);
548  fHistDeltaMCPtOverJet1MCPtvsJet1MCPt->GetXaxis()->SetTitle("p_{T,1}^{MC}");
549  fHistDeltaMCPtOverJet1MCPtvsJet1MCPt->GetYaxis()->SetTitle("#deltap_{T} / p_{T,1}^{MC}");
550  fHistDeltaMCPtOverJet1MCPtvsJet1MCPt->GetZaxis()->SetTitle("counts");
552 
553  fHistDeltaMCPtOverJet2PtvsJet2Pt = new TH2F("fHistDeltaMCPtOverJet2PtvsJet2Pt", "fHistDeltaMCPtOverJet2PtvsJet2Pt",
554  fNbins, fMinBinPt, fMaxBinPt, fNbins, -5, 5);
555  fHistDeltaMCPtOverJet2PtvsJet2Pt->GetXaxis()->SetTitle("p_{T,2}");
556  fHistDeltaMCPtOverJet2PtvsJet2Pt->GetYaxis()->SetTitle("#deltap_{T} / p_{T,2}");
557  fHistDeltaMCPtOverJet2PtvsJet2Pt->GetZaxis()->SetTitle("counts");
559 
560  fHistDeltaMCPtvsDistance = new TH2F("fHistDeltaMCPtvsDistance", "fHistDeltaMCPtvsDistance",
561  fNbins/2, 0, 1.2, 2*fNbins, -fMaxBinPt, fMaxBinPt);
562  fHistDeltaMCPtvsDistance->GetXaxis()->SetTitle("Distance");
563  fHistDeltaMCPtvsDistance->GetYaxis()->SetTitle("#deltap_{T} (GeV/c)");
564  fHistDeltaMCPtvsDistance->GetZaxis()->SetTitle("counts");
566 
567  fHistDeltaMCPtvsCommonEnergy1 = new TH2F("fHistDeltaMCPtvsCommonEnergy1", "fHistDeltaMCPtvsCommonEnergy1",
568  fNbins/2, 0, 1.2, 2*fNbins, -fMaxBinPt, fMaxBinPt);
569  fHistDeltaMCPtvsCommonEnergy1->GetXaxis()->SetTitle("Common energy 1 (%)");
570  fHistDeltaMCPtvsCommonEnergy1->GetYaxis()->SetTitle("#deltap_{T} (GeV/c)");
571  fHistDeltaMCPtvsCommonEnergy1->GetZaxis()->SetTitle("counts");
573 
574  fHistDeltaMCPtvsCommonEnergy2 = new TH2F("fHistDeltaMCPtvsCommonEnergy2", "fHistDeltaMCPtvsCommonEnergy2",
575  fNbins/2, 0, 1.2, 2*fNbins, -fMaxBinPt, fMaxBinPt);
576  fHistDeltaMCPtvsCommonEnergy2->GetXaxis()->SetTitle("Common energy 2 (%)");
577  fHistDeltaMCPtvsCommonEnergy2->GetYaxis()->SetTitle("#deltap_{T} (GeV/c)");
578  fHistDeltaMCPtvsCommonEnergy2->GetZaxis()->SetTitle("counts");
580 
581  fHistDeltaMCPtvsArea1 = new TH2F("fHistDeltaMCPtvsArea1", "fHistDeltaMCPtvsArea1",
582  fNbins/2, 0, 1.5, 2*fNbins, -fMaxBinPt, fMaxBinPt);
583  fHistDeltaMCPtvsArea1->GetXaxis()->SetTitle("A_{jet,1}");
584  fHistDeltaMCPtvsArea1->GetYaxis()->SetTitle("#deltap_{T} (GeV/c)");
585  fHistDeltaMCPtvsArea1->GetZaxis()->SetTitle("counts");
587 
588  fHistDeltaMCPtvsArea2 = new TH2F("fHistDeltaMCPtvsArea2", "fHistDeltaMCPtvsArea2",
589  fNbins/2, 0, 1.5, 2*fNbins, -fMaxBinPt, fMaxBinPt);
590  fHistDeltaMCPtvsArea2->GetXaxis()->SetTitle("A_{jet,2}");
591  fHistDeltaMCPtvsArea2->GetYaxis()->SetTitle("#deltap_{T} (GeV/c)");
592  fHistDeltaMCPtvsArea2->GetZaxis()->SetTitle("counts");
594 
595  fHistDeltaMCPtvsDeltaArea = new TH2F("fHistDeltaMCPtvsDeltaArea", "fHistDeltaMCPtvsDeltaArea",
596  fNbins, -1.+1./fNbins, 1.+1./fNbins, 2*fNbins, -fMaxBinPt, fMaxBinPt);
597  fHistDeltaMCPtvsDeltaArea->GetXaxis()->SetTitle("#deltaA_{jet}");
598  fHistDeltaMCPtvsDeltaArea->GetYaxis()->SetTitle("#deltap_{T} (GeV/c)");
599  fHistDeltaMCPtvsDeltaArea->GetZaxis()->SetTitle("counts");
601 
602  fHistJet1MCPtvsJet2Pt = new TH2F("fHistJet1MCPtvsJet2Pt", "fHistJet1MCPtvsJet2Pt",
604  fHistJet1MCPtvsJet2Pt->GetXaxis()->SetTitle("p_{T,1}^{MC}");
605  fHistJet1MCPtvsJet2Pt->GetYaxis()->SetTitle("p_{T,2}");
606  fHistJet1MCPtvsJet2Pt->GetZaxis()->SetTitle("counts");
608  }
609 }
610 
611 //________________________________________________________________________
613 {
614  // Allocate THnSparse histograms.
615 
616  TString title[25]= {""};
617  Int_t nbins[25] = {0};
618  Double_t min[25] = {0.};
619  Double_t max[25] = {0.};
620  Int_t dim = 0;
621 
622  title[dim] = "#phi";
623  nbins[dim] = fNbins/4;
624  min[dim] = 0;
625  max[dim] = 2*TMath::Pi()*(1 + 1./(nbins[dim]-1));
626  dim++;
627 
628  title[dim] = "#eta";
629  nbins[dim] = fNbins/4;
630  min[dim] = -1;
631  max[dim] = 1;
632  dim++;
633 
634  title[dim] = "p_{T}";
635  nbins[dim] = fNbins;
636  min[dim] = 0;
637  max[dim] = 250;
638  dim++;
639 
640  title[dim] = "A_{jet}";
641  nbins[dim] = fNbins/4;
642  min[dim] = 0;
643  max[dim] = 1.5;
644  dim++;
645 
646  if (fNEFAxis) {
647  title[dim] = "NEF";
648  nbins[dim] = 120;
649  min[dim] = 0.0;
650  max[dim] = 1.2;
651  dim++;
652  }
653 
654  if (fZAxis) {
655  title[dim] = "Z";
656  nbins[dim] = 120;
657  min[dim] = 0.0;
658  max[dim] = 1.2;
659  dim++;
660  }
661 
662  if (fFlavourZAxis) {
663  title[dim] = "z_{flavour}";
664  nbins[dim] = 100;
665  min[dim] = 0.0;
666  max[dim] = 2.0;
667  dim++;
668  }
669 
670  if (fFlavourPtAxis) {
671  title[dim] = "p_{T}^{D}";
672  nbins[dim] = fNbins/2;
673  min[dim] = 0;
674  max[dim] = 125;
675  dim++;
676  }
677 
678  title[dim] = "p_{T,particle}^{leading} (GeV/c)";
679  nbins[dim] = 120;
680  min[dim] = 0;
681  max[dim] = 120;
682  dim++;
683 
684  Int_t dim1 = dim, dim2 = dim;
685 
686  if (fIsJet1Rho) {
687  title[dim1] = "p_{T}^{corr}";
688  nbins[dim1] = fNbins*2;
689  min[dim1] = -250;
690  max[dim1] = 250;
691  dim1++;
692  }
693 
694  if (fIsEmbedded) {
695  title[dim1] = "p_{T}^{MC}";
696  nbins[dim1] = fNbins;
697  min[dim1] = 0;
698  max[dim1] = 250;
699  dim1++;
700  }
701 
702  fHistJets1 = new THnSparseD("fHistJets1","fHistJets1",dim1,nbins,min,max);
703  for (Int_t i = 0; i < dim1; i++)
704  fHistJets1->GetAxis(i)->SetTitle(title[i]);
705  fOutput->Add(fHistJets1);
706 
707  if (fIsJet2Rho) {
708  title[dim2] = "p_{T}^{corr}";
709  nbins[dim2] = fNbins*2;
710  min[dim2] = -250;
711  max[dim2] = 250;
712  dim2++;
713  }
714 
715  fHistJets2 = new THnSparseD("fHistJets2","fHistJets2",dim2,nbins,min,max);
716  for (Int_t i = 0; i < dim2; i++)
717  fHistJets2->GetAxis(i)->SetTitle(title[i]);
718  fOutput->Add(fHistJets2);
719 
720  // Matching
721 
722  dim = 0;
723 
724  title[dim] = "p_{T,1}";
725  nbins[dim] = fNbins;
726  min[dim] = 0;
727  max[dim] = 250;
728  dim++;
729 
730  title[dim] = "p_{T,2}";
731  nbins[dim] = fNbins;
732  min[dim] = 0;
733  max[dim] = 250;
734  dim++;
735 
736  title[dim] = "A_{jet,1}";
737  nbins[dim] = fNbins/4;
738  min[dim] = 0;
739  max[dim] = 1.5;
740  dim++;
741 
742  title[dim] = "A_{jet,2}";
743  nbins[dim] = fNbins/4;
744  min[dim] = 0;
745  max[dim] = 1.5;
746  dim++;
747 
748  title[dim] = "distance";
749  nbins[dim] = fNbins/2;
750  min[dim] = 0;
751  max[dim] = 1.2;
752  dim++;
753 
754  title[dim] = "CE1";
755  nbins[dim] = fNbins/2;
756  min[dim] = 0;
757  max[dim] = 1.2;
758  dim++;
759 
760  title[dim] = "CE2";
761  nbins[dim] = fNbins/2;
762  min[dim] = 0;
763  max[dim] = 1.2;
764  dim++;
765 
766  title[dim] = "p_{T,particle,1}^{leading} (GeV/c)";
767  nbins[dim] = 120;
768  min[dim] = 0;
769  max[dim] = 120;
770  dim++;
771 
772  title[dim] = "p_{T,particle,2}^{leading} (GeV/c)";
773  nbins[dim] = 120;
774  min[dim] = 0;
775  max[dim] = 120;
776  dim++;
777 
778  if (fDeltaPtAxis) {
779  title[dim] = "#deltaA_{jet}";
780  nbins[dim] = fNbins/2;
781  min[dim] = -1;
782  max[dim] = 1;
783  dim++;
784 
785  title[dim] = "#deltap_{T}";
786  nbins[dim] = fNbins*2;
787  min[dim] = -250;
788  max[dim] = 250;
789  dim++;
790  }
791  if (fIsJet1Rho) {
792  title[dim] = "p_{T,1}^{corr}";
793  nbins[dim] = fNbins*2;
794  min[dim] = -250;
795  max[dim] = 250;
796  dim++;
797  }
798  if (fIsJet2Rho) {
799  title[dim] = "p_{T,2}^{corr}";
800  nbins[dim] = fNbins*2;
801  min[dim] = -250;
802  max[dim] = 250;
803  dim++;
804  }
805  if (fDeltaPtAxis && (fIsJet1Rho || fIsJet2Rho)) {
806  title[dim] = "#deltap_{T}^{corr}";
807  nbins[dim] = fNbins*2;
808  min[dim] = -250;
809  max[dim] = 250;
810  dim++;
811  }
812  if (fDeltaEtaDeltaPhiAxis) {
813  title[dim] = "#delta#eta";
814  nbins[dim] = fNbins/2;
815  min[dim] = -1;
816  max[dim] = 1;
817  dim++;
818 
819  title[dim] = "#delta#phi";
820  nbins[dim] = fNbins/2;
821  min[dim] = -TMath::Pi()/2;
822  max[dim] = TMath::Pi()*3/2;
823  dim++;
824  }
825  if (fIsEmbedded) {
826  title[dim] = "p_{T,1}^{MC}";
827  nbins[dim] = fNbins;
828  min[dim] = 0;
829  max[dim] = 250;
830  dim++;
831 
832  if (fDeltaPtAxis) {
833  title[dim] = "#deltap_{T}^{MC}";
834  nbins[dim] = fNbins*2;
835  min[dim] = -250;
836  max[dim] = 250;
837  dim++;
838  }
839  }
840 
841  if (fNEFAxis) {
842  title[dim] = "NEF_{1}";
843  nbins[dim] = 120;
844  min[dim] = 0.0;
845  max[dim] = 1.2;
846  dim++;
847 
848  title[dim] = "NEF_{2}";
849  nbins[dim] = 120;
850  min[dim] = 0.0;
851  max[dim] = 1.2;
852  dim++;
853  }
854 
855  if (fZAxis) {
856  title[dim] = "Z_{1}";
857  nbins[dim] = 120;
858  min[dim] = 0.0;
859  max[dim] = 1.2;
860  dim++;
861 
862  title[dim] = "Z_{2}";
863  nbins[dim] = 120;
864  min[dim] = 0.0;
865  max[dim] = 1.2;
866  dim++;
867  }
868 
869  if (fFlavourZAxis) {
870  title[dim] = "z_{flavour,1}";
871  nbins[dim] = 100;
872  min[dim] = 0.0;
873  max[dim] = 2.0;
874  dim++;
875 
876  title[dim] = "z_{flavour,2}";
877  nbins[dim] = 100;
878  min[dim] = 0.0;
879  max[dim] = 2.0;
880  dim++;
881  }
882 
883  if (fFlavourPtAxis) {
884  title[dim] = "p_{T,1}^{D}";
885  nbins[dim] = fNbins/2;
886  min[dim] = 0;
887  max[dim] = 125;
888  dim++;
889 
890  title[dim] = "p_{T,2}^{D}";
891  nbins[dim] = fNbins/2;
892  min[dim] = 0;
893  max[dim] = 125;
894  dim++;
895  }
896 
897  fHistMatching = new THnSparseD("fHistMatching","fHistMatching",dim,nbins,min,max);
898 
899  for (Int_t i = 0; i < dim; i++)
900  fHistMatching->GetAxis(i)->SetTitle(title[i]);
901 
902  fOutput->Add(fHistMatching);
903 }
904 
905 //________________________________________________________________________
907 {
908  // Create user objects.
909 
911 
912  AliJetContainer *jets1 = static_cast<AliJetContainer*>(fJetCollArray.At(0));
913  AliJetContainer *jets2 = static_cast<AliJetContainer*>(fJetCollArray.At(1));
914 
915  if (!jets1 || !jets2) return;
916 
917  if (jets1->GetRhoName().IsNull()) fIsJet1Rho = kFALSE;
918  else fIsJet1Rho = kTRUE;
919  if (jets2->GetRhoName().IsNull()) fIsJet2Rho = kFALSE;
920  else fIsJet2Rho = kTRUE;
921 
922  fHistRejectionReason1 = new TH2F("fHistRejectionReason1", "fHistRejectionReason1", 32, 0, 32, 100, 0, 250);
923  fHistRejectionReason1->GetXaxis()->SetTitle("Rejection reason");
924  fHistRejectionReason1->GetYaxis()->SetTitle("p_{T,jet} (GeV/c)");
925  fHistRejectionReason1->GetZaxis()->SetTitle("counts");
928 
929  fHistRejectionReason2 = new TH2F("fHistRejectionReason2", "fHistRejectionReason2", 32, 0, 32, 100, 0, 250);
930  fHistRejectionReason2->GetXaxis()->SetTitle("Rejection reason");
931  fHistRejectionReason2->GetYaxis()->SetTitle("p_{T,jet} (GeV/c)");
932  fHistRejectionReason2->GetZaxis()->SetTitle("counts");
935 
936  if (fHistoType==0)
937  AllocateTH2();
938  else
940 
941  PostData(1, fOutput); // Post data for ALL output slots > 0 here, to get at least an empty histogram
942 }
943 
944 //________________________________________________________________________
946 {
947  AliJetContainer* jets = GetJetContainer(Set-1);
948 
949  AliTLorentzVector leadPart;
950  jets->GetLeadingHadronMomentum(leadPart, jet);
951  Double_t zleading = GetParallelFraction(leadPart.Vect(), jet);
952  if (zleading == 1 || (zleading > 1 && zleading - 1 < 1e-3)) zleading = 0.999; // so that it will contribute to the bin 0.9-1 rather than 1-1.1
953 
954  Double_t corrpt = jet->Pt() - jets->GetRhoVal() * jet->Area();
955  Double_t zflavour = 0;
956  Double_t ptflavour = 0;
957  AliVParticle* hftrack = jet->GetFlavourTrack();
958  if (hftrack) {
959  zflavour = GetParallelFraction(hftrack, jet);
960  ptflavour = hftrack->Pt();
961 
962  if (zflavour == 1 || (zflavour > 1 && zflavour - 1 < 1e-3)) zflavour = 0.999; // so that it will contribute to the bin 0.9-1 rather than 1-1.1
963  }
964 
965  if (fHistoType==1) {
966  THnSparse *histo = 0;
967  if (Set==1) {
968  histo = fHistJets1;
969  }
970  else if (Set==2) {
971  histo = fHistJets2;
972  }
973 
974  if (!histo) return;
975 
976  Double_t contents[20]={0};
977 
978  for (Int_t i = 0; i < histo->GetNdimensions(); i++) {
979  TString title(histo->GetAxis(i)->GetTitle());
980  if (title=="#phi")
981  contents[i] = jet->Phi();
982  else if (title=="#eta")
983  contents[i] = jet->Eta();
984  else if (title=="p_{T}")
985  contents[i] = jet->Pt();
986  else if (title=="A_{jet}")
987  contents[i] = jet->Area();
988  else if (title=="NEF")
989  contents[i] = jet->NEF();
990  else if (title=="Z")
991  contents[i] = zleading;
992  else if (title=="p_{T}^{corr}")
993  contents[i] = corrpt;
994  else if (title=="p_{T}^{MC}")
995  contents[i] = jet->MCPt();
996  else if (title=="p_{T,particle}^{leading} (GeV/c)")
997  contents[i] = leadPart.Pt();
998  else if (title=="z_{flavour}")
999  contents[i] = zflavour;
1000  else if (title=="p_{T}^{D}")
1001  contents[i] = ptflavour;
1002  else
1003  AliWarning(Form("Unable to fill dimension %s!",title.Data()));
1004  }
1005 
1006  histo->Fill(contents);
1007  }
1008  else {
1009  if (Set == 1) {
1010  fHistJets1PtArea->Fill(jet->Area(), jet->Pt());
1011  fHistJets1PhiEta->Fill(jet->Eta(), jet->Phi());
1012  fHistJets1ZvsPt->Fill(zleading, jet->Pt());
1013  fHistJets1NEFvsPt->Fill(jet->NEF(), jet->Pt());
1014  fHistJets1CEFvsCEFPt->Fill(1-jet->NEF(), (1-jet->NEF())*jet->Pt());
1015  if (fIsJet1Rho)
1016  fHistJets1CorrPtArea->Fill(jet->Area(), corrpt);
1017  }
1018  else if (Set == 2) {
1019  fHistJets2PtArea->Fill(jet->Area(), jet->Pt());
1020  fHistJets2PhiEta->Fill(jet->Eta(), jet->Phi());
1021  fHistJets2ZvsPt->Fill(zleading, jet->Pt());
1022  fHistJets2NEFvsPt->Fill(jet->NEF(), jet->Pt());
1023  fHistJets2CEFvsCEFPt->Fill(1-jet->NEF(), (1-jet->NEF())*jet->Pt());
1024  if (fIsJet2Rho)
1025  fHistJets2CorrPtArea->Fill(jet->Area(), corrpt);
1026  }
1027  }
1028 }
1029 
1030 //________________________________________________________________________
1032 {
1033  AliJetContainer* jets1 = GetJetContainer(0);
1034  AliJetContainer* jets2 = GetJetContainer(1);
1035 
1036  /*
1037  Printf("Detector level jet");
1038  jet1->Print();
1039 
1040  jet1->PrintConstituents(jets1->GetParticleContainer() != 0 ? jets1->GetParticleContainer()->GetArray() : 0,
1041  jets1->GetClusterContainer() != 0 ? jets1->GetClusterContainer()->GetArray() : 0);
1042 
1043  Printf("Matched with particle level jet");
1044  jet2->Print();
1045  jet2->PrintConstituents(jets2->GetParticleContainer() != 0 ? jets2->GetParticleContainer()->GetArray() : 0,
1046  jets2->GetClusterContainer() != 0 ? jets2->GetClusterContainer()->GetArray() : 0);
1047  */
1048 
1049  AliTLorentzVector leadPart1;
1050  jets1->GetLeadingHadronMomentum(leadPart1, jet1);
1051  Double_t zleading1 = GetParallelFraction(leadPart1.Vect(), jet1);
1052  if (zleading1 == 1 || (zleading1 > 1 && zleading1 - 1 < 1e-3)) zleading1 = 0.999; // so that it will contribute to the bin 0.9-1 rather than 1-1.1
1053 
1054  Double_t corrpt1 = jet1->Pt() - jets1->GetRhoVal() * jet1->Area();
1055  Double_t zflavour1 = 0;
1056  Double_t ptflavour1 = 0;
1057  AliVParticle* hftrack1 = jet1->GetFlavourTrack();
1058  if (hftrack1) {
1059  zflavour1 = GetParallelFraction(hftrack1, jet1);
1060  ptflavour1 = hftrack1->Pt();
1061 
1062  if (zflavour1 == 1 || (zflavour1 > 1 && zflavour1 - 1 < 1e-3)) zflavour1 = 0.999; // so that it will contribute to the bin 0.9-1 rather than 1-1.1
1063  }
1064 
1065 
1066  AliTLorentzVector leadPart2;
1067  jets2->GetLeadingHadronMomentum(leadPart2, jet2);
1068  Double_t zleading2 = GetParallelFraction(leadPart2.Vect(), jet2);
1069  if (zleading2 == 1 || (zleading2 > 1 && zleading2 - 1 < 1e-3)) zleading2 = 0.999; // so that it will contribute to the bin 0.9-1 rather than 1-1.1
1070 
1071  Double_t corrpt2 = jet2->Pt() - jets2->GetRhoVal() * jet2->Area();
1072  Double_t zflavour2 = 0;
1073  Double_t ptflavour2 = 0;
1074  AliVParticle* hftrack2 = jet2->GetFlavourTrack();
1075  if (hftrack2) {
1076  zflavour2 = GetParallelFraction(hftrack2, jet2);
1077  ptflavour2 = hftrack2->Pt();
1078 
1079  if (zflavour2 == 1 || (zflavour2 > 1 && zflavour2 - 1 < 1e-3)) zflavour2 = 0.999; // so that it will contribute to the bin 0.9-1 rather than 1-1.1
1080  }
1081 
1082  if (fHistoType==1) {
1083  Double_t contents[20]={0};
1084 
1085  for (Int_t i = 0; i < fHistMatching->GetNdimensions(); i++) {
1086  TString title(fHistMatching->GetAxis(i)->GetTitle());
1087  if (title=="p_{T,1}")
1088  contents[i] = jet1->Pt();
1089  else if (title=="p_{T,2}")
1090  contents[i] = jet2->Pt();
1091  else if (title=="A_{jet,1}")
1092  contents[i] = jet1->Area();
1093  else if (title=="A_{jet,2}")
1094  contents[i] = jet2->Area();
1095  else if (title=="distance")
1096  contents[i] = d;
1097  else if (title=="CE1")
1098  contents[i] = CE1;
1099  else if (title=="CE2")
1100  contents[i] = CE2;
1101  else if (title=="#deltaA_{jet}")
1102  contents[i] = jet1->Area()-jet2->Area();
1103  else if (title=="#deltap_{T}")
1104  contents[i] = jet1->Pt()-jet2->Pt();
1105  else if (title=="#delta#eta")
1106  contents[i] = jet1->Eta()-jet2->Eta();
1107  else if (title=="#delta#phi")
1108  contents[i] = jet1->Phi()-jet2->Phi();
1109  else if (title=="p_{T,1}^{corr}")
1110  contents[i] = corrpt1;
1111  else if (title=="p_{T,2}^{corr}")
1112  contents[i] = corrpt2;
1113  else if (title=="#deltap_{T}^{corr}")
1114  contents[i] = corrpt1-corrpt2;
1115  else if (title=="p_{T,1}^{MC}")
1116  contents[i] = jet1->MCPt();
1117  else if (title=="#deltap_{T}^{MC}")
1118  contents[i] = jet1->MCPt()-jet2->Pt();
1119  else if (title=="NEF_{1}")
1120  contents[i] = jet1->NEF();
1121  else if (title=="NEF_{2}")
1122  contents[i] = jet2->NEF();
1123  else if (title=="Z_{1}")
1124  contents[i] = zleading1;
1125  else if (title=="Z_{2}")
1126  contents[i] = zleading2;
1127  else if (title=="p_{T,particle,1}^{leading} (GeV/c)")
1128  contents[i] = leadPart1.Pt();
1129  else if (title=="p_{T,particle,2}^{leading} (GeV/c)")
1130  contents[i] = leadPart2.Pt();
1131  else if (title=="z_{flavour,1}")
1132  contents[i] = zflavour1;
1133  else if (title=="z_{flavour,2}")
1134  contents[i] = zflavour2;
1135  else if (title=="p_{T,1}^{D}")
1136  contents[i] = ptflavour1;
1137  else if (title=="p_{T,2}^{D}")
1138  contents[i] = ptflavour2;
1139  else
1140  AliWarning(Form("Unable to fill dimension %s!",title.Data()));
1141  }
1142 
1143  fHistMatching->Fill(contents);
1144  }
1145  else {
1146  fHistCommonEnergy1vsJet1Pt->Fill(CE1, jet1->Pt());
1147  fHistCommonEnergy2vsJet2Pt->Fill(CE2, jet2->Pt());
1148 
1149  fHistDistancevsJet1Pt->Fill(d, jet1->Pt());
1150  fHistDistancevsJet2Pt->Fill(d, jet2->Pt());
1151 
1152  fHistDistancevsCommonEnergy1->Fill(d, CE1);
1153  fHistDistancevsCommonEnergy2->Fill(d, CE2);
1154  fHistCommonEnergy1vsCommonEnergy2->Fill(CE1, CE2);
1155 
1156  fHistDeltaEtaDeltaPhi->Fill(jet1->Eta()-jet2->Eta(),jet1->Phi()-jet2->Phi());
1157 
1158  fHistJet2PtOverJet1PtvsJet2Pt->Fill(jet2->Pt(), jet2->Pt() / jet1->Pt());
1159  fHistJet1PtOverJet2PtvsJet1Pt->Fill(jet1->Pt(), jet1->Pt() / jet2->Pt());
1160 
1161  Double_t dpt = jet1->Pt() - jet2->Pt();
1162  fHistDeltaPtvsJet1Pt->Fill(jet1->Pt(), dpt);
1163  fHistDeltaPtvsJet2Pt->Fill(jet2->Pt(), dpt);
1164  fHistDeltaPtOverJet1PtvsJet1Pt->Fill(jet1->Pt(), dpt/jet1->Pt());
1165  fHistDeltaPtOverJet2PtvsJet2Pt->Fill(jet2->Pt(), dpt/jet2->Pt());
1166 
1167  fHistDeltaPtvsDistance->Fill(d, dpt);
1168  fHistDeltaPtvsCommonEnergy1->Fill(CE1, dpt);
1169  fHistDeltaPtvsCommonEnergy2->Fill(CE2, dpt);
1170 
1171  fHistDeltaPtvsArea1->Fill(jet1->Area(), dpt);
1172  fHistDeltaPtvsArea2->Fill(jet2->Area(), dpt);
1173 
1174  Double_t darea = jet1->Area() - jet2->Area();
1175  fHistDeltaPtvsDeltaArea->Fill(darea, dpt);
1176 
1177  fHistJet1PtvsJet2Pt->Fill(jet1->Pt(), jet2->Pt());
1178 
1179  if (fIsJet1Rho || fIsJet2Rho) {
1180  Double_t dcorrpt = corrpt1 - corrpt2;
1181  fHistDeltaCorrPtvsJet1CorrPt->Fill(corrpt1, dcorrpt);
1182  fHistDeltaCorrPtvsJet2CorrPt->Fill(corrpt2, dcorrpt);
1183  fHistDeltaCorrPtOverJet1CorrPtvsJet1CorrPt->Fill(corrpt1, dcorrpt/corrpt1);
1184  fHistDeltaCorrPtOverJet2CorrPtvsJet2CorrPt->Fill(corrpt2, dcorrpt/corrpt2);
1185  fHistDeltaCorrPtvsDistance->Fill(d, dcorrpt);
1186  fHistDeltaCorrPtvsCommonEnergy1->Fill(CE1, dcorrpt);
1187  fHistDeltaCorrPtvsCommonEnergy2->Fill(CE2, dcorrpt);
1188  fHistDeltaCorrPtvsArea1->Fill(jet1->Area(), dcorrpt);
1189  fHistDeltaCorrPtvsArea2->Fill(jet2->Area(), dcorrpt);
1190  fHistDeltaCorrPtvsDeltaArea->Fill(darea, dcorrpt);
1191  fHistJet1CorrPtvsJet2CorrPt->Fill(corrpt1, corrpt2);
1192  }
1193 
1194  if (fIsEmbedded) {
1195  Double_t dmcpt = jet1->MCPt() - jet2->Pt();
1196  fHistDeltaMCPtvsJet1MCPt->Fill(jet1->MCPt(), dmcpt);
1197  fHistDeltaMCPtvsJet2Pt->Fill(jet2->MCPt(), dmcpt);
1198  fHistDeltaMCPtOverJet1MCPtvsJet1MCPt->Fill(jet1->MCPt(), dmcpt/jet1->MCPt());
1199  fHistDeltaMCPtOverJet2PtvsJet2Pt->Fill(jet2->Pt(), dmcpt/jet2->Pt());
1200  fHistDeltaMCPtvsDistance->Fill(d, dmcpt);
1201  fHistDeltaMCPtvsCommonEnergy1->Fill(CE1, dmcpt);
1202  fHistDeltaMCPtvsCommonEnergy2->Fill(CE2, dmcpt);
1203  fHistDeltaMCPtvsArea1->Fill(jet1->Area(), dmcpt);
1204  fHistDeltaMCPtvsArea2->Fill(jet2->Area(), dmcpt);
1205  fHistDeltaMCPtvsDeltaArea->Fill(darea, dmcpt);
1206  fHistJet1MCPtvsJet2Pt->Fill(jet1->MCPt(), jet2->Pt());
1207  }
1208  }
1209 }
1210 
1211 //________________________________________________________________________
1213 {
1214  // Execute once.
1215 
1217 
1218  AliJetContainer *jets1 = static_cast<AliJetContainer*>(fJetCollArray.At(0));
1219  AliJetContainer *jets2 = static_cast<AliJetContainer*>(fJetCollArray.At(1));
1220 
1221  if (!jets1 || !jets1->GetArray() || !jets2 || !jets2->GetArray()) return;
1222 
1223  if (fMatching == kMCLabel && (!jets2->GetIsParticleLevel() || jets1->GetIsParticleLevel())) {
1224  if (jets1->GetIsParticleLevel() == jets2->GetIsParticleLevel()) {
1225  AliWarning("Changing matching type from MC label to same collection...");
1227  }
1228  else {
1229  AliWarning("Changing matching type from MC label to geometrical...");
1231  }
1232  }
1233  else if (fMatching == kSameCollections && jets1->GetIsParticleLevel() != jets2->GetIsParticleLevel()) {
1234  if (jets2->GetIsParticleLevel() && !jets1->GetIsParticleLevel()) {
1235  AliWarning("Changing matching type from same collection to MC label...");
1236  fMatching = kMCLabel;
1237  }
1238  else {
1239  AliWarning("Changing matching type from same collection to geometrical...");
1241  }
1242  }
1243 }
1244 
1245 //________________________________________________________________________
1247 {
1248  // Find the closest jets
1249 
1250  if (fMatching == kNoMatching)
1251  return kTRUE;
1252  else
1253  return DoJetMatching();
1254 }
1255 
1256 //________________________________________________________________________
1258 {
1259  AliJetContainer *jets1 = static_cast<AliJetContainer*>(fJetCollArray.At(0));
1260  AliJetContainer *jets2 = static_cast<AliJetContainer*>(fJetCollArray.At(1));
1261 
1262  if (!jets1 || !jets1->GetArray() || !jets2 || !jets2->GetArray()) return kFALSE;
1263 
1264  DoJetLoop();
1265 
1266  AliEmcalJet* jet1 = 0;
1267 
1268  jets1->ResetCurrentID();
1269  while ((jet1 = jets1->GetNextJet())) {
1270 
1271  AliEmcalJet *jet2 = jet1->ClosestJet();
1272 
1273  if (!jet2) continue;
1274  if (jet2->ClosestJet() != jet1) continue;
1275  if (jet1->ClosestJetDistance() > fMatchingPar1 || jet2->ClosestJetDistance() > fMatchingPar2) continue;
1276 
1277  // Matched jet found
1280  AliDebug(2,Form("Found matching: jet1 pt = %f, eta = %f, phi = %f, jet2 pt = %f, eta = %f, phi = %f",
1281  jet1->Pt(), jet1->Eta(), jet1->Phi(),
1282  jet2->Pt(), jet2->Eta(), jet2->Phi()));
1283  }
1284 
1285  return kTRUE;
1286 }
1287 
1288 //________________________________________________________________________
1290 {
1291  // Do the jet loop.
1292 
1293  AliJetContainer *jets1 = static_cast<AliJetContainer*>(fJetCollArray.At(0));
1294  AliJetContainer *jets2 = static_cast<AliJetContainer*>(fJetCollArray.At(1));
1295 
1296  if (!jets1 || !jets1->GetArray() || !jets2 || !jets2->GetArray()) return;
1297 
1298  AliEmcalJet* jet1 = 0;
1299  AliEmcalJet* jet2 = 0;
1300 
1301  jets2->ResetCurrentID();
1302  while ((jet2 = jets2->GetNextJet())) jet2->ResetMatching();
1303 
1304  jets1->ResetCurrentID();
1305  while ((jet1 = jets1->GetNextJet())) {
1306  jet1->ResetMatching();
1307 
1308  if (jet1->MCPt() < fMinJetMCPt) continue;
1309 
1310  jets2->ResetCurrentID();
1311  while ((jet2 = jets2->GetNextJet())) {
1312  SetMatchingLevel(jet1, jet2, fMatching);
1313  } // jet2 loop
1314  } // jet1 loop
1315 }
1316 
1317 //________________________________________________________________________
1319 {
1320  d = jet1->DeltaR(jet2);
1321 }
1322 
1323 //________________________________________________________________________
1325 {
1326  AliJetContainer *jets1 = static_cast<AliJetContainer*>(fJetCollArray.At(0));
1327  AliJetContainer *jets2 = static_cast<AliJetContainer*>(fJetCollArray.At(1));
1328 
1329  if (!jets1 || !jets1->GetArray() || !jets2 || !jets2->GetArray()) return;
1330 
1331  AliParticleContainer *tracks1 = jets1->GetParticleContainer();
1332  AliClusterContainer *clusters1 = jets1->GetClusterContainer();
1333  AliParticleContainer *tracks2 = jets2->GetParticleContainer();
1334 
1335  // d1 and d2 represent the matching level: 0 = maximum level of matching, 1 = the two jets are completely unrelated
1336  d1 = jet1->Pt();
1337  d2 = jet2->Pt();
1338  Double_t totalPt1 = d1; // the total pt of the reconstructed jet will be cleaned from the background
1339 
1340  // remove completely tracks that are not MC particles (label == 0)
1341  if (tracks1 && tracks1->GetArray()) {
1342  for (Int_t iTrack = 0; iTrack < jet1->GetNumberOfTracks(); iTrack++) {
1343  AliVParticle *track = jet1->TrackAt(iTrack,tracks1->GetArray());
1344  if (!track) {
1345  AliWarning(Form("Could not find track %d!", iTrack));
1346  continue;
1347  }
1348 
1349  Int_t MClabel = TMath::Abs(track->GetLabel());
1350  MClabel -= fMCLabelShift;
1351  if (MClabel != 0) continue;
1352 
1353  // this is not a MC particle; remove it completely
1354  AliDebug(3,Form("Track %d (pT = %f) is not a MC particle (MClabel = %d)!",iTrack,track->Pt(),MClabel));
1355  totalPt1 -= track->Pt();
1356  d1 -= track->Pt();
1357  }
1358  }
1359 
1360  // remove completely clusters that are not MC particles (label == 0)
1361  if (fUseCellsToMatch && fCaloCells) {
1362  for (Int_t iClus = 0; iClus < jet1->GetNumberOfClusters(); iClus++) {
1363  AliVCluster *clus = jet1->ClusterAt(iClus,clusters1->GetArray());
1364  if (!clus) {
1365  AliWarning(Form("Could not find cluster %d!", iClus));
1366  continue;
1367  }
1368  AliTLorentzVector part;
1369  clus->GetMomentum(part, fVertex);
1370 
1371  for (Int_t iCell = 0; iCell < clus->GetNCells(); iCell++) {
1372  Int_t cellId = clus->GetCellAbsId(iCell);
1373  Double_t cellFrac = clus->GetCellAmplitudeFraction(iCell);
1374 
1375  Int_t MClabel = TMath::Abs(fCaloCells->GetCellMCLabel(cellId));
1376  MClabel -= fMCLabelShift;
1377  if (MClabel != 0) continue;
1378 
1379  // this is not a MC particle; remove it completely
1380  AliDebug(3,Form("Cell %d (frac = %f) is not a MC particle (MClabel = %d)!",iCell,cellFrac,MClabel));
1381  totalPt1 -= part.Pt() * cellFrac;
1382  d1 -= part.Pt() * cellFrac;
1383  }
1384  }
1385  }
1386  else {
1387  for (Int_t iClus = 0; iClus < jet1->GetNumberOfClusters(); iClus++) {
1388  AliVCluster *clus = jet1->ClusterAt(iClus,clusters1->GetArray());
1389  if (!clus) {
1390  AliWarning(Form("Could not find cluster %d!", iClus));
1391  continue;
1392  }
1393  TLorentzVector part;
1394  clus->GetMomentum(part, fVertex);
1395 
1396  Int_t MClabel = TMath::Abs(clus->GetLabel());
1397  MClabel -= fMCLabelShift;
1398  if (MClabel != 0) continue;
1399 
1400  // this is not a MC particle; remove it completely
1401  AliDebug(3,Form("Cluster %d (pT = %f) is not a MC particle (MClabel = %d)!",iClus,part.Pt(),MClabel));
1402  totalPt1 -= part.Pt();
1403  d1 -= part.Pt();
1404  }
1405  }
1406 
1407  for (Int_t iTrack2 = 0; iTrack2 < jet2->GetNumberOfTracks(); iTrack2++) {
1408  Bool_t track2Found = kFALSE;
1409  Int_t index2 = jet2->TrackAt(iTrack2);
1410 
1411  // now look for common particles in the track array
1412  for (Int_t iTrack = 0; iTrack < jet1->GetNumberOfTracks(); iTrack++) {
1413  AliVParticle *track = jet1->TrackAt(iTrack,tracks1->GetArray());
1414  if (!track) {
1415  AliWarning(Form("Could not find track %d!", iTrack));
1416  continue;
1417  }
1418  Int_t MClabel = TMath::Abs(track->GetLabel());
1419  MClabel -= fMCLabelShift;
1420  if (MClabel <= 0) continue;
1421 
1422  Int_t index = -1;
1423  index = tracks2->GetIndexFromLabel(MClabel);
1424  if (index < 0) {
1425  AliDebug(2,Form("Track %d (pT = %f) does not have an associated MC particle (MClabel = %d)!",iTrack,track->Pt(),MClabel));
1426  continue;
1427  }
1428 
1429  if (index2 != index) continue;
1430 
1431  // found common particle
1432  d1 -= track->Pt();
1433 
1434  if (!track2Found) {
1435  AliVParticle *MCpart = tracks2->GetParticle(index2);
1436  AliDebug(3,Form("Track %d (pT = %f, eta = %f, phi = %f) is associated with the MC particle %d (pT = %f, eta = %f, phi = %f)!",
1437  iTrack,track->Pt(),track->Eta(),track->Phi(),MClabel,MCpart->Pt(),MCpart->Eta(),MCpart->Phi()));
1438  d2 -= MCpart->Pt();
1439  }
1440 
1441  track2Found = kTRUE;
1442  }
1443 
1444  // now look for common particles in the cluster array
1445  if (fUseCellsToMatch && fCaloCells) { // if the cell colection is available, look for cells with a matched MC particle
1446  for (Int_t iClus = 0; iClus < jet1->GetNumberOfClusters(); iClus++) {
1447  AliVCluster *clus = jet1->ClusterAt(iClus,fCaloClusters);
1448  if (!clus) {
1449  AliWarning(Form("Could not find cluster %d!", iClus));
1450  continue;
1451  }
1452  AliTLorentzVector part;
1453  clus->GetMomentum(part, fVertex);
1454 
1455  for (Int_t iCell = 0; iCell < clus->GetNCells(); iCell++) {
1456  Int_t cellId = clus->GetCellAbsId(iCell);
1457  Double_t cellFrac = clus->GetCellAmplitudeFraction(iCell);
1458 
1459  Int_t MClabel = TMath::Abs(fCaloCells->GetCellMCLabel(cellId));
1460  MClabel -= fMCLabelShift;
1461  if (MClabel <= 0) continue;
1462 
1463  Int_t index1 = -1;
1464  index1 = tracks2->GetIndexFromLabel(MClabel);
1465  if (index1 < 0) {
1466  AliDebug(3,Form("Cell %d (frac = %f) does not have an associated MC particle (MClabel = %d)!",iCell,cellFrac,MClabel));
1467  continue;
1468  }
1469 
1470  if (index2 != index1) continue;
1471 
1472  // found common particle
1473  d1 -= part.Pt() * cellFrac;
1474 
1475  if (!track2Found) { // only if it is not already found among charged tracks (charged particles are most likely already found)
1476  AliVParticle *MCpart = tracks2->GetParticle(index2);
1477  AliDebug(3,Form("Cell %d belonging to cluster %d (pT = %f, eta = %f, phi = %f) is associated with the MC particle %d (pT = %f, eta = %f, phi = %f)!",
1478  iCell,iClus,part.Pt(),part.Eta(),part.Phi_0_2pi(),MClabel,MCpart->Pt(),MCpart->Eta(),MCpart->Phi()));
1479  d2 -= MCpart->Pt() * cellFrac;
1480  }
1481 
1482  track2Found = kTRUE;
1483  }
1484  }
1485  }
1486  else { //otherwise look for the first contributor to the cluster, and if matched to a MC label remove it
1487  for (Int_t iClus = 0; iClus < jet1->GetNumberOfClusters(); iClus++) {
1488  AliVCluster *clus = jet1->ClusterAt(iClus,fCaloClusters);
1489  if (!clus) {
1490  AliWarning(Form("Could not find cluster %d!", iClus));
1491  continue;
1492  }
1493  AliTLorentzVector part;
1494  clus->GetMomentum(part, fVertex);
1495 
1496  Int_t MClabel = TMath::Abs(clus->GetLabel());
1497  MClabel -= fMCLabelShift;
1498  if (MClabel <= 0) continue;
1499 
1500  Int_t index = -1;
1501  index = tracks2->GetIndexFromLabel(MClabel);
1502 
1503  if (index < 0) {
1504  AliDebug(3,Form("Cluster %d (pT = %f) does not have an associated MC particle (MClabel = %d)!",iClus,part.Pt(),MClabel));
1505  continue;
1506  }
1507 
1508  if (index2 != index) continue;
1509 
1510  // found common particle
1511  d1 -= part.Pt();
1512 
1513  if (!track2Found) { // only if it is not already found among charged tracks (charged particles are most likely already found)
1514  AliVParticle *MCpart = tracks2->GetParticle(index2);
1515  AliDebug(3,Form("Cluster %d (pT = %f, eta = %f, phi = %f) is associated with the MC particle %d (pT = %f, eta = %f, phi = %f)!",
1516  iClus,part.Pt(),part.Eta(),part.Phi_0_2pi(),MClabel,MCpart->Pt(),MCpart->Eta(),MCpart->Phi()));
1517 
1518  d2 -= MCpart->Pt();
1519  }
1520 
1521  track2Found = kTRUE;
1522  }
1523  }
1524  }
1525 
1526  if (d1 < 0)
1527  d1 = 0;
1528 
1529  if (d2 < 0)
1530  d2 = 0;
1531 
1532  if (totalPt1 < 1)
1533  d1 = -1;
1534  else
1535  d1 /= totalPt1;
1536 
1537  if (jet2->Pt() < 1)
1538  d2 = -1;
1539  else
1540  d2 /= jet2->Pt();
1541 }
1542 
1543 //________________________________________________________________________
1545 {
1546  AliJetContainer *jets1 = static_cast<AliJetContainer*>(fJetCollArray.At(0));
1547  AliJetContainer *jets2 = static_cast<AliJetContainer*>(fJetCollArray.At(1));
1548 
1549  if (!jets1 || !jets1->GetArray() || !jets2 || !jets2->GetArray()) return;
1550 
1551  AliParticleContainer *tracks1 = jets1->GetParticleContainer();
1552  AliClusterContainer *clusters1 = jets1->GetClusterContainer();
1553  AliParticleContainer *tracks2 = jets2->GetParticleContainer();
1554  AliClusterContainer *clusters2 = jets2->GetClusterContainer();
1555 
1556  // d1 and d2 represent the matching level: 0 = maximum level of matching, 1 = the two jets are completely unrelated
1557  d1 = jet1->Pt();
1558  d2 = jet2->Pt();
1559 
1560  if (tracks1 && tracks2) {
1561 
1562  for (Int_t iTrack2 = 0; iTrack2 < jet2->GetNumberOfTracks(); iTrack2++) {
1563  Int_t index2 = jet2->TrackAt(iTrack2);
1564  for (Int_t iTrack1 = 0; iTrack1 < jet1->GetNumberOfTracks(); iTrack1++) {
1565  Int_t index1 = jet1->TrackAt(iTrack1);
1566  if (index2 == index1) { // found common particle
1567  AliVParticle *part1 = tracks1->GetParticle(index1);
1568  if (!part1) {
1569  AliWarning(Form("Could not find track %d!", index1));
1570  continue;
1571  }
1572  AliVParticle *part2 = tracks2->GetParticle(index2);
1573  if (!part2) {
1574  AliWarning(Form("Could not find track %d!", index2));
1575  continue;
1576  }
1577 
1578  d1 -= part1->Pt();
1579  d2 -= part2->Pt();
1580  break;
1581  }
1582  }
1583  }
1584 
1585  }
1586 
1587  if (clusters1 && clusters2) {
1588 
1589  if (fUseCellsToMatch && fCaloCells) {
1590  // Note: this section of the code needs to be revised and tested heavily
1591  // While fixing some inconsistencies in AliAnalysisTaskEMCALClusterizeFast
1592  // some issues came up, e.g. memory leaks (fixed) and inconsistent use of
1593  // fCaloCells. In principle the two cluster collections may use cells
1594  // from different sources (embedding / non-embedding). This is not handled
1595  // correctly in the current version of this code.
1596  AliWarning("ATTENTION ATTENTION ATTENTION: this section of the AliJetResponseMaker code needs to be revised and tested before using it for physics!!!");
1597  const Int_t nClus1 = jet1->GetNumberOfClusters();
1598 
1599  Int_t ncells1[nClus1];
1600  UShort_t *cellsId1[nClus1];
1601  Double_t *cellsFrac1[nClus1];
1602  Double_t *cellsClusFrac1[nClus1];
1603  Int_t *sortedIndexes1[nClus1];
1604  Double_t ptClus1[nClus1];
1605  for (Int_t iClus1 = 0; iClus1 < nClus1; iClus1++) {
1606  Int_t index1 = jet1->ClusterAt(iClus1);
1607  AliVCluster *clus1 = clusters1->GetCluster(index1);
1608  if (!clus1) {
1609  AliWarning(Form("Could not find cluster %d!", index1));
1610  ncells1[iClus1] = 0;
1611  cellsId1[iClus1] = 0;
1612  cellsFrac1[iClus1] = 0;
1613  cellsClusFrac1[iClus1] = 0;
1614  sortedIndexes1[iClus1] = 0;
1615  ptClus1[iClus1] = 0;
1616  continue;
1617  }
1618  TLorentzVector part1;
1619  clus1->GetMomentum(part1, fVertex);
1620 
1621  ncells1[iClus1] = clus1->GetNCells();
1622  cellsId1[iClus1] = clus1->GetCellsAbsId();
1623  cellsFrac1[iClus1] = clus1->GetCellsAmplitudeFraction();
1624  cellsClusFrac1[iClus1] = new Double_t[ncells1[iClus1]];
1625  sortedIndexes1[iClus1] = new Int_t[ncells1[iClus1]];
1626  ptClus1[iClus1] = part1.Pt();
1627 
1628  for (Int_t iCell = 0; iCell < ncells1[iClus1]; iCell++) {
1629  cellsClusFrac1[iClus1][iCell] = fCaloCells->GetCellAmplitude(cellsId1[iClus1][iCell]) / clus1->E();
1630  }
1631 
1632  TMath::Sort(ncells1[iClus1], cellsId1[iClus1], sortedIndexes1[iClus1], kFALSE);
1633  }
1634 
1635  const Int_t nClus2 = jet2->GetNumberOfClusters();
1636 
1637  const Int_t maxNcells2 = 11520;
1638  Int_t sortedIndexes2[maxNcells2];
1639  for (Int_t iClus2 = 0; iClus2 < nClus2; iClus2++) {
1640  Int_t index2 = jet2->ClusterAt(iClus2);
1641  AliVCluster *clus2 = clusters2->GetCluster(index2);
1642  if (!clus2) {
1643  AliWarning(Form("Could not find cluster %d!", index2));
1644  continue;
1645  }
1646  Int_t ncells2 = clus2->GetNCells();
1647  if (ncells2 >= maxNcells2) {
1648  AliError(Form("Number of cells in the cluster %d >= %d",ncells2,maxNcells2));
1649  continue;
1650  }
1651  UShort_t *cellsId2 = clus2->GetCellsAbsId();
1652  Double_t *cellsFrac2 = clus2->GetCellsAmplitudeFraction();
1653  Double_t *cellsClusFrac2 = new Double_t[ncells2];
1654 
1655  for (Int_t iCell = 0; iCell < ncells2; iCell++) {
1656  cellsClusFrac2[iCell] = fCaloCells->GetCellAmplitude(cellsId2[iCell]) / clus2->E();
1657  }
1658 
1659  TLorentzVector part2;
1660  clus2->GetMomentum(part2, fVertex);
1661  Double_t ptClus2 = part2.Pt();
1662 
1663  TMath::Sort(ncells2, cellsId2, sortedIndexes2, kFALSE);
1664 
1665  for (Int_t iClus1 = 0; iClus1 < nClus1; iClus1++) {
1666  if (sortedIndexes1[iClus1] == 0)
1667  continue;
1668  Int_t iCell1 = 0, iCell2 = 0;
1669  while (iCell1 < ncells1[iClus1] && iCell2 < ncells2) {
1670  if (cellsId1[iClus1][sortedIndexes1[iClus1][iCell1]] == cellsId2[sortedIndexes2[iCell2]]) { // found a common cell
1671  d1 -= cellsFrac1[iClus1][sortedIndexes1[iClus1][iCell1]] * cellsClusFrac1[iClus1][sortedIndexes1[iClus1][iCell1]] * ptClus1[iClus1];
1672  d2 -= cellsFrac2[sortedIndexes2[iCell2]] * cellsClusFrac2[sortedIndexes2[iCell2]] * ptClus2;
1673  iCell1++;
1674  iCell2++;
1675  }
1676  else if (cellsId1[iClus1][sortedIndexes1[iClus1][iCell1]] > cellsId2[sortedIndexes2[iCell2]]) {
1677  iCell2++;
1678  }
1679  else {
1680  iCell1++;
1681  }
1682  }
1683  }
1684  delete[] cellsClusFrac2;
1685  }
1686  for (Int_t iClus1 = 0; iClus1 < nClus1; iClus1++) {
1687  delete[] cellsClusFrac1[iClus1];
1688  delete[] sortedIndexes1[iClus1];
1689  }
1690  }
1691  else {
1692  for (Int_t iClus2 = 0; iClus2 < jet2->GetNumberOfClusters(); iClus2++) {
1693  Int_t index2 = jet2->ClusterAt(iClus2);
1694  for (Int_t iClus1 = 0; iClus1 < jet1->GetNumberOfClusters(); iClus1++) {
1695  Int_t index1 = jet1->ClusterAt(iClus1);
1696  if (index2 == index1) { // found common particle
1697  AliVCluster *clus1 = clusters1->GetCluster(index1);
1698  if (!clus1) {
1699  AliWarning(Form("Could not find cluster %d!", index1));
1700  continue;
1701  }
1702  AliVCluster *clus2 = clusters2->GetCluster(index2);
1703  if (!clus2) {
1704  AliWarning(Form("Could not find cluster %d!", index2));
1705  continue;
1706  }
1707  TLorentzVector part1, part2;
1708  clus1->GetMomentum(part1, fVertex);
1709  clus2->GetMomentum(part2, fVertex);
1710 
1711  d1 -= part1.Pt();
1712  d2 -= part2.Pt();
1713  break;
1714  }
1715  }
1716  }
1717  }
1718  }
1719 
1720  if (d1 < 0)
1721  d1 = 0;
1722 
1723  if (d2 < 0)
1724  d2 = 0;
1725 
1726  if (jet1->Pt() > 0)
1727  d1 /= jet1->Pt();
1728  else
1729  d1 = -1;
1730 
1731  if (jet2->Pt() > 0)
1732  d2 /= jet2->Pt();
1733  else
1734  d2 = -1;
1735 }
1736 
1737 //________________________________________________________________________
1739 {
1740  Double_t d1 = -1;
1741  Double_t d2 = -1;
1742 
1743  switch (matching) {
1744  case kGeometrical:
1745  GetGeometricalMatchingLevel(jet1,jet2,d1);
1746  d2 = d1;
1747  break;
1748  case kMCLabel: // jet1 = detector level and jet2 = particle level!
1749  GetMCLabelMatchingLevel(jet1,jet2,d1,d2);
1750  break;
1751  case kSameCollections:
1752  GetSameCollectionsMatchingLevel(jet1,jet2,d1,d2);
1753  break;
1754  default:
1755  ;
1756  }
1757 
1758  if (d1 >= 0) {
1759 
1760  if (d1 < jet1->ClosestJetDistance()) {
1761  jet1->SetSecondClosestJet(jet1->ClosestJet(), jet1->ClosestJetDistance());
1762  jet1->SetClosestJet(jet2, d1);
1763  }
1764  else if (d1 < jet1->SecondClosestJetDistance()) {
1765  jet1->SetSecondClosestJet(jet2, d1);
1766  }
1767  }
1768 
1769  if (d2 >= 0) {
1770 
1771  if (d2 < jet2->ClosestJetDistance()) {
1772  jet2->SetSecondClosestJet(jet2->ClosestJet(), jet2->ClosestJetDistance());
1773  jet2->SetClosestJet(jet1, d2);
1774  }
1775  else if (d2 < jet2->SecondClosestJetDistance()) {
1776  jet2->SetSecondClosestJet(jet1, d2);
1777  }
1778  }
1779 }
1780 
1781 //________________________________________________________________________
1783 {
1784  // Fill histograms.
1785 
1786  AliJetContainer *jets1 = static_cast<AliJetContainer*>(fJetCollArray.At(0));
1787  AliJetContainer *jets2 = static_cast<AliJetContainer*>(fJetCollArray.At(1));
1788 
1789  if (!jets1 || !jets1->GetArray() || !jets2 || !jets2->GetArray()) return kFALSE;
1790 
1791  AliEmcalJet* jet1 = 0;
1792  AliEmcalJet* jet2 = 0;
1793 
1794  jets2->ResetCurrentID();
1795  while ((jet2 = jets2->GetNextJet())) {
1796 
1797  AliDebug(2,Form("Processing jet (2) %d", jets2->GetCurrentID()));
1798 
1799  if (jet2->Pt() < jets2->GetJetPtCut()) continue;
1800 
1801  UInt_t rejectionReason = 0;
1802  if (jets2->AcceptJet(jet2, rejectionReason))
1803  FillJetHisto(jet2, 2);
1804  else
1805  fHistRejectionReason2->Fill(jets2->GetRejectionReasonBitPosition(rejectionReason), jet2->Pt());
1806 
1807  jet1 = jet2->MatchedJet();
1808 
1809  if (!jet1) continue;
1810  rejectionReason = 0;
1811  if (!jets1->AcceptJet(jet1, rejectionReason)) continue;
1812  if (jet1->MCPt() < fMinJetMCPt) continue;
1813 
1814  Double_t d=-1, ce1=-1, ce2=-1;
1815  if (jet2->GetMatchingType() == kGeometrical) {
1816  if (jets2->GetIsParticleLevel() && !jets1->GetIsParticleLevel()) // the other way around is not supported
1817  GetMCLabelMatchingLevel(jet1, jet2, ce1, ce2);
1818  else if (jets1->GetIsParticleLevel() == jets2->GetIsParticleLevel())
1819  GetSameCollectionsMatchingLevel(jet1, jet2, ce1, ce2);
1820 
1821  d = jet2->ClosestJetDistance();
1822  }
1823  else if (jet2->GetMatchingType() == kMCLabel || jet2->GetMatchingType() == kSameCollections) {
1824  GetGeometricalMatchingLevel(jet1, jet2, d);
1825 
1826  ce1 = jet1->ClosestJetDistance();
1827  ce2 = jet2->ClosestJetDistance();
1828  }
1829 
1830  FillMatchingHistos(jet1, jet2, d, ce1, ce2);
1831  }
1832 
1833  jets1->ResetCurrentID();
1834  while ((jet1 = jets1->GetNextJet())) {
1835  UInt_t rejectionReason = 0;
1836  if (!jets1->AcceptJet(jet1, rejectionReason)) {
1837  fHistRejectionReason1->Fill(jets1->GetRejectionReasonBitPosition(rejectionReason), jet1->Pt());
1838  continue;
1839  }
1840  if (jet1->MCPt() < fMinJetMCPt) continue;
1841  AliDebug(2,Form("Processing jet (1) %d", jets1->GetCurrentID()));
1842 
1843  FillJetHisto(jet1, 1);
1844  }
1845  return kTRUE;
1846 }
void SetSecondClosestJet(AliEmcalJet *j, Double_t d)
Definition: AliEmcalJet.h:210
TH2 * fHistJets2PtArea
phi-eta distribution of jets 2
Double_t Area() const
Definition: AliEmcalJet.h:117
TH2 * fHistRejectionReason1
whether the jet2 collection has to be average subtracted
TH2 * fHistJets2CorrPtArea
inclusive jet pt vs. area histogram 2
Double_t GetRhoVal() const
const TString & GetRhoName() const
double Double_t
Definition: External.C:58
TH2 * fHistDeltaPtvsCommonEnergy1
delta pt between matched jets vs distance
TH2 * fHistJet2PtOverJet1PtvsJet2Pt
delta eta vs delta phi of matched jets
Definition: External.C:236
const char * title
Definition: MakeQAPdf.C:26
Double_t MCPt() const
Definition: AliEmcalJet.h:140
TH2 * fHistDeltaPtvsDeltaArea
delta pt between matched jets vs jet 2 area
AliEmcalJet * ClosestJet() const
Definition: AliEmcalJet.h:213
AliJetContainer * GetJetContainer(Int_t i=0) const
Double_t ClosestJetDistance() const
Definition: AliEmcalJet.h:214
Double_t Eta() const
Definition: AliEmcalJet.h:108
TH2 * fHistJets2PhiEta
Constituent Pt over Jet Pt ratio vs. jet pt 1.
TH2 * fHistDeltaCorrPtvsArea2
delta pt corr between matched jets vs jet 1 area
void FillJetHisto(AliEmcalJet *jet, Int_t Set)
Double_t Phi() const
Definition: AliEmcalJet.h:104
AliEmcalJet * MatchedJet() const
Definition: AliEmcalJet.h:217
Declaration of class AliTLorentzVector.
void GetSameCollectionsMatchingLevel(AliEmcalJet *jet1, AliEmcalJet *jet2, Double_t &d1, Double_t &d2) const
Double_t fMinBinPt
min pt in histograms
TH2 * fHistDeltaMCPtvsCommonEnergy2
jet 1 MC pt - jet2 pt vs common energy 1 (%)
AliClusterContainer * GetClusterContainer() const
TH2 * fHistJets1ZvsPt
Jet charged energy fraction vs. charged jet pt 1.
virtual Bool_t AcceptJet(Int_t i, UInt_t &rejectionReason) const
THnSparse * fHistMatching
jet2 THnSparse
TH2 * fHistDeltaEtaDeltaPhi
common energy 1 (%) vs common energy 2 (%)
ClassImp(AliJetResponseMaker) AliJetResponseMaker
TH2 * fHistDeltaMCPtvsCommonEnergy1
jet 1 MC pt - jet2 pt vs distance
TH2 * fHistDeltaCorrPtOverJet2CorrPtvsJet2CorrPt
delta pt corr / jet 1 corr pt between matched jets vs jet 1 corr pt
TH2 * fHistCommonEnergy1vsCommonEnergy2
distance vs common energy 2 (%)
TH2 * fHistDeltaPtvsCommonEnergy2
delta pt between matched jets vs common energy 1 (%)
TH2 * fHistDeltaMCPtvsDistance
jet 1 MC pt - jet2 pt vs jet 2 pt
AliVParticle * GetFlavourTrack(Int_t i=0) const
Container for particles within the EMCAL framework.
TH2 * fHistJets1NEFvsPt
inclusive jet pt vs. area histogram 1
TH2 * fHistJets2ZvsPt
Jet charged energy fraction vs. charged jet pt 2.
Bool_t fIsEmbedded
trigger, embedded signal
UShort_t GetNumberOfTracks() const
Definition: AliEmcalJet.h:126
Bool_t fIsJet2Rho
whether the jet1 collection has to be average subtracted
void ResetMatching()
TH2 * fHistDeltaPtvsJet1Pt
jet 1 pt over jet 2 pt vs jet 1 pt
TH2 * fHistDeltaCorrPtvsCommonEnergy2
delta pt corr between matched jets vs common energy 1 (%)
AliParticleContainer * GetParticleContainer() const
void GetLeadingHadronMomentum(TLorentzVector &mom, const AliEmcalJet *jet) const
TClonesArray * fCaloClusters
!clusters
TH2 * fHistDeltaPtvsJet2Pt
delta pt between matched jets vs jet 1 pt
Short_t ClusterAt(Int_t idx) const
Definition: AliEmcalJet.h:124
TH2 * fHistJets1PhiEta
matching THnSparse
int Int_t
Definition: External.C:63
UShort_t GetNumberOfClusters() const
Definition: AliEmcalJet.h:125
unsigned int UInt_t
Definition: External.C:33
THnSparse * fHistJets1
Rejection reason vs. jet pt.
UShort_t GetMatchingType() const
Definition: AliEmcalJet.h:218
TH2 * fHistJet1PtOverJet2PtvsJet1Pt
jet 2 pt over jet 1 pt vs jet 2 pt
Double_t Phi_0_2pi() const
TH2 * fHistDeltaCorrPtvsDistance
delta pt corr between matched jets vs jet 2 corr pt
virtual AliVParticle * GetParticle(Int_t i=-1) const
void FillMatchingHistos(AliEmcalJet *jet1, AliEmcalJet *jet2, Double_t d, Double_t CE1, Double_t CE2)
TH2 * fHistJet1PtvsJet2Pt
delta pt between matched jets vs delta area
TH2 * fHistDeltaCorrPtvsJet1CorrPt
delta pt corr / jet 2 corr pt between matched jets vs jet 2 corr pt
TH2 * fHistDeltaMCPtOverJet2PtvsJet2Pt
jet 1 MC pt - jet2 pt / jet 1 MC pt vs jet 1 pt
TH2 * fHistDistancevsCommonEnergy2
distance vs common energy 1 (%)
Double_t GetJetPtCut() const
TH2 * fHistRejectionReason2
Rejection reason vs. jet pt.
Int_t fMCLabelShift
if MC label > fMCLabelShift, MC label -= fMCLabelShift
TH2 * fHistDistancevsJet1Pt
common energy 2 (%) vs jet 2 pt
TH2 * fHistDeltaMCPtOverJet1MCPtvsJet1MCPt
correlation jet 1 corr pt vs jet 2 corr pt
static Double_t GetParallelFraction(AliVParticle *part1, AliVParticle *part2)
AliVCluster * GetCluster(Int_t i) const
TH2 * fHistDeltaPtOverJet1PtvsJet1Pt
delta pt between matched jets vs jet 2 pt
TObjArray fJetCollArray
jet collection array
Double_t DeltaR(const AliVParticle *part) const
TH2 * fHistDeltaCorrPtvsArea1
delta pt corr between matched jets vs common energy 2 (%)
AliVCaloCells * fCaloCells
!cells
TH2 * fHistCommonEnergy1vsJet1Pt
Constituent Pt over Jet Pt ratio vs. jet pt 2.
void SetMatchingLevel(AliEmcalJet *jet1, AliEmcalJet *jet2, MatchingType matching)
Double_t Pt() const
Definition: AliEmcalJet.h:96
TH2 * fHistJets1CorrPtArea
inclusive jet pt vs. area histogram 1
AliEmcalJet * GetNextJet()
AliEmcalList * fOutput
!output list
Double_t fMaxBinPt
max pt in histograms
TH2 * fHistDeltaMCPtvsJet1MCPt
jet 1 MC pt - jet2 pt / jet 2 pt vs jet 2 pt
void SetClosestJet(AliEmcalJet *j, Double_t d)
Definition: AliEmcalJet.h:209
TH2 * fHistDeltaMCPtvsArea2
jet 1 MC pt - jet2 pt vs jet 1 area
Short_t TrackAt(Int_t idx) const
Definition: AliEmcalJet.h:147
TH2 * fHistDeltaMCPtvsJet2Pt
jet 1 MC pt - jet2 pt vs jet 1 MC pt
TH2 * fHistDistancevsJet2Pt
distance vs jet 1 pt
TH2 * fHistDeltaPtvsArea1
delta pt between matched jets vs common energy 2 (%)
Double_t fVertex[3]
!event vertex
TH2 * fHistDeltaCorrPtvsJet2CorrPt
delta pt corr between matched jets vs jet 1 corr pt
void GetGeometricalMatchingLevel(AliEmcalJet *jet1, AliEmcalJet *jet2, Double_t &d) const
THnSparse * fHistJets2
jet1 THnSparse
TH2 * fHistDeltaPtOverJet2PtvsJet2Pt
delta pt / jet 1 pt between matched jets vs jet 1 pt
TH2 * fHistDeltaPtvsArea2
delta pt between matched jets vs jet 1 area
void SetMakeGeneralHistograms(Bool_t g)
TH2 * fHistJets2NEFvsPt
inclusive jet pt vs. area histogram 2
Base task in the EMCAL jet framework.
Represent a jet reconstructed using the EMCal jet framework.
Definition: AliEmcalJet.h:44
TH2 * fHistCommonEnergy2vsJet2Pt
common energy 1 (%) vs jet 1 pt
unsigned short UShort_t
Definition: External.C:28
TH2 * fHistDeltaPtvsDistance
delta pt / jet 2 pt between matched jets vs jet 2 pt
TH2 * fHistJet1MCPtvsJet2Pt
jet 1 MC pt - jet2 pt vs delta area
TH2 * fHistJets2CEFvsCEFPt
Jet neutral energy fraction vs. jet pt 2.
void SetRejectionReasonLabels(TAxis *axis)
TH2 * fHistDistancevsCommonEnergy1
distance vs jet 2 pt
TH2 * fHistJets1CEFvsCEFPt
Jet neutral energy fraction vs. jet pt 1.
const Int_t nbins
bool Bool_t
Definition: External.C:53
TH2 * fHistDeltaMCPtvsArea1
jet 1 MC pt - jet2 pt vs common energy 2 (%)
TH2 * fHistDeltaCorrPtOverJet1CorrPtvsJet1CorrPt
correlation jet 1 pt vs jet 2 pt
TH2 * fHistJets1PtArea
phi-eta distribution of jets 1
Double_t NEF() const
Definition: AliEmcalJet.h:135
TH2 * fHistDeltaMCPtvsDeltaArea
jet 1 MC pt - jet2 pt vs jet 2 area
Container structure for EMCAL clusters.
TH2 * fHistDeltaCorrPtvsCommonEnergy1
delta pt corr between matched jets vs distance
void GetMCLabelMatchingLevel(AliEmcalJet *jet1, AliEmcalJet *jet2, Double_t &d1, Double_t &d2) const
Container for jet within the EMCAL jet framework.
void SetMatchedToClosest(UShort_t m)
Definition: AliEmcalJet.h:211
Int_t fNbins
no. of pt bins
TH2 * fHistJet1CorrPtvsJet2CorrPt
delta pt corr between matched jets vs delta area
TH2 * fHistDeltaCorrPtvsDeltaArea
delta pt corr between matched jets vs jet 2 area