]> git.uio.no Git - u/mrichter/AliRoot.git/blame - PWGJE/EMCALJetTasks/AliJetResponseMaker.cxx
git://git.uio.no / u / mrichter / AliRoot.git / blame
? search:
summary | shortlog | log | commit | commitdiff | tree
blob | blame (incremental) | history | HEAD
Change binning
[u/mrichter/AliRoot.git] / PWGJE / EMCALJetTasks / AliJetResponseMaker.cxx
CommitLineData
7549c451 1// $Id$
2949a2ec 2//
3// Emcal jet response matrix maker task.
4//
5// Author: S. Aiola
6
7#include "AliJetResponseMaker.h"
8
2949a2ec 9#include <TClonesArray.h>
2949a2ec 10#include <TH2F.h>
43032ce2 11#include <THnSparse.h>
2949a2ec 12#include <TLorentzVector.h>
13
ca5c29fa 14#include "AliAnalysisManager.h"
2949a2ec 15#include "AliVCluster.h"
16#include "AliVTrack.h"
17#include "AliEmcalJet.h"
787a3c4f 18#include "AliLog.h"
cd6431de 19#include "AliRhoParameter.h"
5be3857d 20#include "AliNamedArrayI.h"
7cd832c7 21#include "AliJetContainer.h"
22#include "AliParticleContainer.h"
23#include "AliClusterContainer.h"
2949a2ec 24
25ClassImp(AliJetResponseMaker)
26
27//________________________________________________________________________
28AliJetResponseMaker::AliJetResponseMaker() :
9239b066 29 AliAnalysisTaskEmcalJet("AliJetResponseMaker", kTRUE),
cd6431de 30 fMatching(kNoMatching),
7f76e479 31 fMatchingPar1(0),
32 fMatchingPar2(0),
a7477843 33 fUseCellsToMatch(kFALSE),
8159f4cd 34 fMinJetMCPt(1),
43032ce2 35 fHistoType(0),
36 fDeltaPtAxis(0),
37 fDeltaEtaDeltaPhiAxis(0),
faf21244 38 fNEFAxis(0),
39 fZAxis(0),
7cd832c7 40 fIsJet1Rho(kFALSE),
41 fIsJet2Rho(kFALSE),
43032ce2 42 fHistJets1(0),
43 fHistJets2(0),
43032ce2 44 fHistMatching(0),
cd6431de 45 fHistJets1PhiEta(0),
46 fHistJets1PtArea(0),
47 fHistJets1CorrPtArea(0),
63fac07f 48 fHistJets1NEFvsPt(0),
49 fHistJets1CEFvsCEFPt(0),
50 fHistJets1ZvsPt(0),
cd6431de 51 fHistJets2PhiEta(0),
52 fHistJets2PtArea(0),
53 fHistJets2CorrPtArea(0),
63fac07f 54 fHistJets2NEFvsPt(0),
55 fHistJets2CEFvsCEFPt(0),
56 fHistJets2ZvsPt(0),
ca5c29fa 57 fHistCommonEnergy1vsJet1Pt(0),
58 fHistCommonEnergy2vsJet2Pt(0),
59 fHistDistancevsJet1Pt(0),
60 fHistDistancevsJet2Pt(0),
6a20534a 61 fHistDistancevsCommonEnergy1(0),
62 fHistDistancevsCommonEnergy2(0),
05077f28 63 fHistCommonEnergy1vsCommonEnergy2(0),
43032ce2 64 fHistDeltaEtaDeltaPhi(0),
c560b734 65 fHistJet2PtOverJet1PtvsJet2Pt(0),
7030f36f 66 fHistJet1PtOverJet2PtvsJet1Pt(0),
7030f36f 67 fHistDeltaPtvsJet1Pt(0),
cfc2ac24 68 fHistDeltaPtvsJet2Pt(0),
4c00a7af 69 fHistDeltaPtOverJet1PtvsJet1Pt(0),
70 fHistDeltaPtOverJet2PtvsJet2Pt(0),
05077f28 71 fHistDeltaPtvsDistance(0),
72 fHistDeltaPtvsCommonEnergy1(0),
73 fHistDeltaPtvsCommonEnergy2(0),
74 fHistDeltaPtvsArea1(0),
75 fHistDeltaPtvsArea2(0),
76 fHistDeltaPtvsDeltaArea(0),
77 fHistJet1PtvsJet2Pt(0),
4c00a7af 78 fHistDeltaCorrPtOverJet1CorrPtvsJet1CorrPt(0),
79 fHistDeltaCorrPtOverJet2CorrPtvsJet2CorrPt(0),
43032ce2 80 fHistDeltaCorrPtvsJet1CorrPt(0),
81 fHistDeltaCorrPtvsJet2CorrPt(0),
05077f28 82 fHistDeltaCorrPtvsDistance(0),
83 fHistDeltaCorrPtvsCommonEnergy1(0),
84 fHistDeltaCorrPtvsCommonEnergy2(0),
85 fHistDeltaCorrPtvsArea1(0),
86 fHistDeltaCorrPtvsArea2(0),
87 fHistDeltaCorrPtvsDeltaArea(0),
cd6431de 88 fHistJet1CorrPtvsJet2CorrPt(0),
4c00a7af 89 fHistDeltaMCPtOverJet1MCPtvsJet1MCPt(0),
90 fHistDeltaMCPtOverJet2PtvsJet2Pt(0),
43032ce2 91 fHistDeltaMCPtvsJet1MCPt(0),
92 fHistDeltaMCPtvsJet2Pt(0),
05077f28 93 fHistDeltaMCPtvsDistance(0),
94 fHistDeltaMCPtvsCommonEnergy1(0),
95 fHistDeltaMCPtvsCommonEnergy2(0),
96 fHistDeltaMCPtvsArea1(0),
97 fHistDeltaMCPtvsArea2(0),
98 fHistDeltaMCPtvsDeltaArea(0),
99 fHistJet1MCPtvsJet2Pt(0)
2949a2ec 100{
101 // Default constructor.
4643d2e8 102
a487deae 103 SetMakeGeneralHistograms(kTRUE);
2949a2ec 104}
105
106//________________________________________________________________________
107AliJetResponseMaker::AliJetResponseMaker(const char *name) :
9239b066 108 AliAnalysisTaskEmcalJet(name, kTRUE),
cd6431de 109 fMatching(kNoMatching),
7f76e479 110 fMatchingPar1(0),
111 fMatchingPar2(0),
a7477843 112 fUseCellsToMatch(kFALSE),
8159f4cd 113 fMinJetMCPt(1),
43032ce2 114 fHistoType(0),
115 fDeltaPtAxis(0),
116 fDeltaEtaDeltaPhiAxis(0),
faf21244 117 fNEFAxis(0),
118 fZAxis(0),
7cd832c7 119 fIsJet1Rho(kFALSE),
120 fIsJet2Rho(kFALSE),
43032ce2 121 fHistJets1(0),
122 fHistJets2(0),
43032ce2 123 fHistMatching(0),
cd6431de 124 fHistJets1PhiEta(0),
125 fHistJets1PtArea(0),
126 fHistJets1CorrPtArea(0),
63fac07f 127 fHistJets1NEFvsPt(0),
128 fHistJets1CEFvsCEFPt(0),
129 fHistJets1ZvsPt(0),
cd6431de 130 fHistJets2PhiEta(0),
131 fHistJets2PtArea(0),
132 fHistJets2CorrPtArea(0),
63fac07f 133 fHistJets2NEFvsPt(0),
134 fHistJets2CEFvsCEFPt(0),
135 fHistJets2ZvsPt(0),
ca5c29fa 136 fHistCommonEnergy1vsJet1Pt(0),
137 fHistCommonEnergy2vsJet2Pt(0),
138 fHistDistancevsJet1Pt(0),
139 fHistDistancevsJet2Pt(0),
6a20534a 140 fHistDistancevsCommonEnergy1(0),
141 fHistDistancevsCommonEnergy2(0),
05077f28 142 fHistCommonEnergy1vsCommonEnergy2(0),
43032ce2 143 fHistDeltaEtaDeltaPhi(0),
c560b734 144 fHistJet2PtOverJet1PtvsJet2Pt(0),
7030f36f 145 fHistJet1PtOverJet2PtvsJet1Pt(0),
7030f36f 146 fHistDeltaPtvsJet1Pt(0),
cfc2ac24 147 fHistDeltaPtvsJet2Pt(0),
4c00a7af 148 fHistDeltaPtOverJet1PtvsJet1Pt(0),
149 fHistDeltaPtOverJet2PtvsJet2Pt(0),
05077f28 150 fHistDeltaPtvsDistance(0),
151 fHistDeltaPtvsCommonEnergy1(0),
152 fHistDeltaPtvsCommonEnergy2(0),
153 fHistDeltaPtvsArea1(0),
154 fHistDeltaPtvsArea2(0),
155 fHistDeltaPtvsDeltaArea(0),
156 fHistJet1PtvsJet2Pt(0),
4c00a7af 157 fHistDeltaCorrPtOverJet1CorrPtvsJet1CorrPt(0),
158 fHistDeltaCorrPtOverJet2CorrPtvsJet2CorrPt(0),
43032ce2 159 fHistDeltaCorrPtvsJet1CorrPt(0),
160 fHistDeltaCorrPtvsJet2CorrPt(0),
05077f28 161 fHistDeltaCorrPtvsDistance(0),
162 fHistDeltaCorrPtvsCommonEnergy1(0),
163 fHistDeltaCorrPtvsCommonEnergy2(0),
164 fHistDeltaCorrPtvsArea1(0),
165 fHistDeltaCorrPtvsArea2(0),
166 fHistDeltaCorrPtvsDeltaArea(0),
cd6431de 167 fHistJet1CorrPtvsJet2CorrPt(0),
4c00a7af 168 fHistDeltaMCPtOverJet1MCPtvsJet1MCPt(0),
169 fHistDeltaMCPtOverJet2PtvsJet2Pt(0),
43032ce2 170 fHistDeltaMCPtvsJet1MCPt(0),
171 fHistDeltaMCPtvsJet2Pt(0),
05077f28 172 fHistDeltaMCPtvsDistance(0),
173 fHistDeltaMCPtvsCommonEnergy1(0),
174 fHistDeltaMCPtvsCommonEnergy2(0),
175 fHistDeltaMCPtvsArea1(0),
176 fHistDeltaMCPtvsArea2(0),
177 fHistDeltaMCPtvsDeltaArea(0),
178 fHistJet1MCPtvsJet2Pt(0)
2949a2ec 179{
180 // Standard constructor.
4643d2e8 181
a487deae 182 SetMakeGeneralHistograms(kTRUE);
2949a2ec 183}
184
185//________________________________________________________________________
186AliJetResponseMaker::~AliJetResponseMaker()
187{
188 // Destructor
189}
190
43032ce2 191
2949a2ec 192//________________________________________________________________________
43032ce2 193void AliJetResponseMaker::AllocateTH2()
2949a2ec 194{
43032ce2 195 // Allocate TH2 histograms.
05077f28 196
cfc2ac24 197 fHistJets1PhiEta = new TH2F("fHistJets1PhiEta", "fHistJets1PhiEta", 40, -1, 1, 40, 0, TMath::Pi()*2);
cd6431de 198 fHistJets1PhiEta->GetXaxis()->SetTitle("#eta");
199 fHistJets1PhiEta->GetYaxis()->SetTitle("#phi");
200 fOutput->Add(fHistJets1PhiEta);
2949a2ec 201
faf21244 202 fHistJets1PtArea = new TH2F("fHistJets1PtArea", "fHistJets1PtArea", fNbins/2, 0, 1.5, fNbins, fMinBinPt, fMaxBinPt);
cd6431de 203 fHistJets1PtArea->GetXaxis()->SetTitle("area");
204 fHistJets1PtArea->GetYaxis()->SetTitle("p_{T,1} (GeV/c)");
205 fHistJets1PtArea->GetZaxis()->SetTitle("counts");
206 fOutput->Add(fHistJets1PtArea);
43032ce2 207
325e2b7b 208 if (fIsJet1Rho) {
9adcb46d 209 fHistJets1CorrPtArea = new TH2F("fHistJets1CorrPtArea", "fHistJets1CorrPtArea",
faf21244 210 fNbins/2, 0, 1.5, 2*fNbins, -fMaxBinPt, fMaxBinPt);
05077f28 211 fHistJets1CorrPtArea->GetXaxis()->SetTitle("area");
212 fHistJets1CorrPtArea->GetYaxis()->SetTitle("p_{T,1}^{corr} (GeV/c)");
213 fHistJets1CorrPtArea->GetZaxis()->SetTitle("counts");
214 fOutput->Add(fHistJets1CorrPtArea);
cd6431de 215 }
216
63fac07f 217 fHistJets1ZvsPt = new TH2F("fHistJets1ZvsPt", "fHistJets1ZvsPt", 120, 0, 1.2, fNbins, fMinBinPt, fMaxBinPt);
218 fHistJets1ZvsPt->GetXaxis()->SetTitle("Z");
219 fHistJets1ZvsPt->GetYaxis()->SetTitle("p_{T,1} (GeV/c)");
220 fHistJets1ZvsPt->GetZaxis()->SetTitle("counts");
221 fOutput->Add(fHistJets1ZvsPt);
222
223 fHistJets1NEFvsPt = new TH2F("fHistJets1NEFvsPt", "fHistJets1NEFvsPt", 120, 0, 1.2, fNbins, fMinBinPt, fMaxBinPt);
224 fHistJets1NEFvsPt->GetXaxis()->SetTitle("NEF");
225 fHistJets1NEFvsPt->GetYaxis()->SetTitle("p_{T,1} (GeV/c)");
226 fHistJets1NEFvsPt->GetZaxis()->SetTitle("counts");
227 fOutput->Add(fHistJets1NEFvsPt);
228
229 fHistJets1CEFvsCEFPt = new TH2F("fHistJets1CEFvsCEFPt", "fHistJets1CEFvsCEFPt", 120, 0, 1.2, fNbins, fMinBinPt, fMaxBinPt);
230 fHistJets1CEFvsCEFPt->GetXaxis()->SetTitle("1-NEF");
231 fHistJets1CEFvsCEFPt->GetYaxis()->SetTitle("(1-NEF)*p_{T,1} (GeV/c)");
232 fHistJets1CEFvsCEFPt->GetZaxis()->SetTitle("counts");
233 fOutput->Add(fHistJets1CEFvsCEFPt);
234
05077f28 235 // Jets 2 histograms
cfc2ac24 236
43032ce2 237 fHistJets2PhiEta = new TH2F("fHistJets2PhiEta", "fHistJets2PhiEta", 40, -1, 1, 40, 0, TMath::Pi()*2);
238 fHistJets2PhiEta->GetXaxis()->SetTitle("#eta");
239 fHistJets2PhiEta->GetYaxis()->SetTitle("#phi");
240 fOutput->Add(fHistJets2PhiEta);
cd6431de 241
faf21244 242 fHistJets2PtArea = new TH2F("fHistJets2PtArea", "fHistJets2PtArea", fNbins/2, 0, 1.5, fNbins, fMinBinPt, fMaxBinPt);
43032ce2 243 fHistJets2PtArea->GetXaxis()->SetTitle("area");
244 fHistJets2PtArea->GetYaxis()->SetTitle("p_{T,2} (GeV/c)");
245 fHistJets2PtArea->GetZaxis()->SetTitle("counts");
246 fOutput->Add(fHistJets2PtArea);
247
325e2b7b 248 if (fIsJet2Rho) {
43032ce2 249 fHistJets2CorrPtArea = new TH2F("fHistJets2CorrPtArea", "fHistJets2CorrPtArea",
faf21244 250 fNbins/2, 0, 1.5, 2*fNbins, -fMaxBinPt, fMaxBinPt);
43032ce2 251 fHistJets2CorrPtArea->GetXaxis()->SetTitle("area");
252 fHistJets2CorrPtArea->GetYaxis()->SetTitle("p_{T,2}^{corr} (GeV/c)");
253 fHistJets2CorrPtArea->GetZaxis()->SetTitle("counts");
254 fOutput->Add(fHistJets2CorrPtArea);
255 }
ca5c29fa 256
63fac07f 257 fHistJets2ZvsPt = new TH2F("fHistJets2ZvsPt", "fHistJets2ZvsPt", 120, 0, 1.2, fNbins, fMinBinPt, fMaxBinPt);
258 fHistJets2ZvsPt->GetXaxis()->SetTitle("Z");
259 fHistJets2ZvsPt->GetYaxis()->SetTitle("p_{T,2} (GeV/c)");
260 fHistJets2ZvsPt->GetZaxis()->SetTitle("counts");
261 fOutput->Add(fHistJets2ZvsPt);
262
263 fHistJets2NEFvsPt = new TH2F("fHistJets2NEFvsPt", "fHistJets2NEFvsPt", 120, 0, 1.2, fNbins, fMinBinPt, fMaxBinPt);
264 fHistJets2NEFvsPt->GetXaxis()->SetTitle("NEF");
265 fHistJets2NEFvsPt->GetYaxis()->SetTitle("p_{T,2} (GeV/c)");
266 fHistJets2NEFvsPt->GetZaxis()->SetTitle("counts");
267 fOutput->Add(fHistJets2NEFvsPt);
268
269 fHistJets2CEFvsCEFPt = new TH2F("fHistJets2CEFvsCEFPt", "fHistJets2CEFvsCEFPt", 120, 0, 1.2, fNbins, fMinBinPt, fMaxBinPt);
270 fHistJets2CEFvsCEFPt->GetXaxis()->SetTitle("1-NEF");
271 fHistJets2CEFvsCEFPt->GetYaxis()->SetTitle("(1-NEF)*p_{T,2} (GeV/c)");
272 fHistJets2CEFvsCEFPt->GetZaxis()->SetTitle("counts");
273 fOutput->Add(fHistJets2CEFvsCEFPt);
274
05077f28 275 // Matching histograms
276
ca5c29fa 277 fHistCommonEnergy1vsJet1Pt = new TH2F("fHistCommonEnergy1vsJet1Pt", "fHistCommonEnergy1vsJet1Pt", fNbins/2, 0, 1.2, fNbins/2, fMinBinPt, fMaxBinPt);
278 fHistCommonEnergy1vsJet1Pt->GetXaxis()->SetTitle("Common energy 1 (%)");
279 fHistCommonEnergy1vsJet1Pt->GetYaxis()->SetTitle("p_{T,1}");
280 fHistCommonEnergy1vsJet1Pt->GetZaxis()->SetTitle("counts");
281 fOutput->Add(fHistCommonEnergy1vsJet1Pt);
282
283 fHistCommonEnergy2vsJet2Pt = new TH2F("fHistCommonEnergy2vsJet2Pt", "fHistCommonEnergy2vsJet2Pt", fNbins/2, 0, 1.2, fNbins/2, fMinBinPt, fMaxBinPt);
284 fHistCommonEnergy2vsJet2Pt->GetXaxis()->SetTitle("Common energy 2 (%)");
285 fHistCommonEnergy2vsJet2Pt->GetYaxis()->SetTitle("p_{T,2}");
286 fHistCommonEnergy2vsJet2Pt->GetZaxis()->SetTitle("counts");
287 fOutput->Add(fHistCommonEnergy2vsJet2Pt);
288
289 fHistDistancevsJet1Pt = new TH2F("fHistDistancevsJet1Pt", "fHistDistancevsJet1Pt", fNbins/2, 0, 1.2, fNbins/2, fMinBinPt, fMaxBinPt);
290 fHistDistancevsJet1Pt->GetXaxis()->SetTitle("Distance");
291 fHistDistancevsJet1Pt->GetYaxis()->SetTitle("p_{T,1}");
292 fHistDistancevsJet1Pt->GetZaxis()->SetTitle("counts");
293 fOutput->Add(fHistDistancevsJet1Pt);
294
295 fHistDistancevsJet2Pt = new TH2F("fHistDistancevsJet2Pt", "fHistDistancevsJet2Pt", fNbins/2, 0, 1.2, fNbins/2, fMinBinPt, fMaxBinPt);
296 fHistDistancevsJet2Pt->GetXaxis()->SetTitle("Distance");
297 fHistDistancevsJet2Pt->GetYaxis()->SetTitle("p_{T,2}");
298 fHistDistancevsJet2Pt->GetZaxis()->SetTitle("counts");
299 fOutput->Add(fHistDistancevsJet2Pt);
cd6431de 300
6a20534a 301 fHistDistancevsCommonEnergy1 = new TH2F("fHistDistancevsCommonEnergy1", "fHistDistancevsCommonEnergy1", fNbins/2, 0, 1.2, fNbins/2, 0, 1.2);
302 fHistDistancevsCommonEnergy1->GetXaxis()->SetTitle("Distance");
303 fHistDistancevsCommonEnergy1->GetYaxis()->SetTitle("Common energy 1 (%)");
304 fHistDistancevsCommonEnergy1->GetZaxis()->SetTitle("counts");
305 fOutput->Add(fHistDistancevsCommonEnergy1);
306
307 fHistDistancevsCommonEnergy2 = new TH2F("fHistDistancevsCommonEnergy2", "fHistDistancevsCommonEnergy2", fNbins/2, 0, 1.2, fNbins/2, 0, 1.2);
308 fHistDistancevsCommonEnergy2->GetXaxis()->SetTitle("Distance");
309 fHistDistancevsCommonEnergy2->GetYaxis()->SetTitle("Common energy 2 (%)");
310 fHistDistancevsCommonEnergy2->GetZaxis()->SetTitle("counts");
311 fOutput->Add(fHistDistancevsCommonEnergy2);
cfc2ac24 312
05077f28 313 fHistCommonEnergy1vsCommonEnergy2 = new TH2F("fHistCommonEnergy1vsCommonEnergy2", "fHistCommonEnergy1vsCommonEnergy2", fNbins/2, 0, 1.2, fNbins/2, 0, 1.2);
314 fHistCommonEnergy1vsCommonEnergy2->GetXaxis()->SetTitle("Common energy 1 (%)");
315 fHistCommonEnergy1vsCommonEnergy2->GetYaxis()->SetTitle("Common energy 2 (%)");
316 fHistCommonEnergy1vsCommonEnergy2->GetZaxis()->SetTitle("counts");
317 fOutput->Add(fHistCommonEnergy1vsCommonEnergy2);
318
43032ce2 319 fHistDeltaEtaDeltaPhi = new TH2F("fHistDeltaEtaDeltaPhi", "fHistDeltaEtaDeltaPhi", fNbins/4, -1, 1, fNbins/4, -TMath::Pi()/2, TMath::Pi()*3/2);
320 fHistDeltaEtaDeltaPhi->GetXaxis()->SetTitle("Common energy 1 (%)");
321 fHistDeltaEtaDeltaPhi->GetYaxis()->SetTitle("Common energy 2 (%)");
322 fHistDeltaEtaDeltaPhi->GetZaxis()->SetTitle("counts");
323 fOutput->Add(fHistDeltaEtaDeltaPhi);
324
7030f36f 325 fHistJet2PtOverJet1PtvsJet2Pt = new TH2F("fHistJet2PtOverJet1PtvsJet2Pt", "fHistJet2PtOverJet1PtvsJet2Pt", fNbins, fMinBinPt, fMaxBinPt, 300, 0, 1.5);
c560b734 326 fHistJet2PtOverJet1PtvsJet2Pt->GetXaxis()->SetTitle("p_{T,2}");
327 fHistJet2PtOverJet1PtvsJet2Pt->GetYaxis()->SetTitle("p_{T,2} / p_{T,1}");
328 fHistJet2PtOverJet1PtvsJet2Pt->GetZaxis()->SetTitle("counts");
329 fOutput->Add(fHistJet2PtOverJet1PtvsJet2Pt);
330
7030f36f 331 fHistJet1PtOverJet2PtvsJet1Pt = new TH2F("fHistJet1PtOverJet2PtvsJet1Pt", "fHistJet1PtOverJet2PtvsJet1Pt", fNbins, fMinBinPt, fMaxBinPt, 300, 0, 1.5);
332 fHistJet1PtOverJet2PtvsJet1Pt->GetXaxis()->SetTitle("p_{T,1}");
333 fHistJet1PtOverJet2PtvsJet1Pt->GetYaxis()->SetTitle("p_{T,1} / p_{T,2}");
334 fHistJet1PtOverJet2PtvsJet1Pt->GetZaxis()->SetTitle("counts");
335 fOutput->Add(fHistJet1PtOverJet2PtvsJet1Pt);
336
9adcb46d 337 fHistDeltaPtvsJet1Pt = new TH2F("fHistDeltaPtvsJet1Pt", "fHistDeltaPtvsJet1Pt",
338 fNbins, fMinBinPt, fMaxBinPt, 2*fNbins, -fMaxBinPt, fMaxBinPt);
7030f36f 339 fHistDeltaPtvsJet1Pt->GetXaxis()->SetTitle("p_{T,1}");
340 fHistDeltaPtvsJet1Pt->GetYaxis()->SetTitle("#deltap_{T} (GeV/c)");
341 fHistDeltaPtvsJet1Pt->GetZaxis()->SetTitle("counts");
342 fOutput->Add(fHistDeltaPtvsJet1Pt);
343
9adcb46d 344 fHistDeltaPtvsJet2Pt = new TH2F("fHistDeltaPtvsJet2Pt", "fHistDeltaPtvsJet2Pt",
345 fNbins, fMinBinPt, fMaxBinPt, 2*fNbins, -fMaxBinPt, fMaxBinPt);
cfc2ac24 346 fHistDeltaPtvsJet2Pt->GetXaxis()->SetTitle("p_{T,2}");
347 fHistDeltaPtvsJet2Pt->GetYaxis()->SetTitle("#deltap_{T} (GeV/c)");
348 fHistDeltaPtvsJet2Pt->GetZaxis()->SetTitle("counts");
349 fOutput->Add(fHistDeltaPtvsJet2Pt);
350
4c00a7af 351 fHistDeltaPtOverJet1PtvsJet1Pt = new TH2F("fHistDeltaPtOverJet1PtvsJet1Pt", "fHistDeltaPtOverJet1PtvsJet1Pt",
352 fNbins, fMinBinPt, fMaxBinPt, fNbins, -5, 5);
353 fHistDeltaPtOverJet1PtvsJet1Pt->GetXaxis()->SetTitle("p_{T,1}");
354 fHistDeltaPtOverJet1PtvsJet1Pt->GetYaxis()->SetTitle("#deltap_{T} / p_{T,1}");
355 fHistDeltaPtOverJet1PtvsJet1Pt->GetZaxis()->SetTitle("counts");
356 fOutput->Add(fHistDeltaPtOverJet1PtvsJet1Pt);
357
358 fHistDeltaPtOverJet2PtvsJet2Pt = new TH2F("fHistDeltaPtOverJet2PtvsJet2Pt", "fHistDeltaPtOverJet2PtvsJet2Pt",
359 fNbins, fMinBinPt, fMaxBinPt, fNbins, -5, 5);
360 fHistDeltaPtOverJet2PtvsJet2Pt->GetXaxis()->SetTitle("p_{T,2}");
361 fHistDeltaPtOverJet2PtvsJet2Pt->GetYaxis()->SetTitle("#deltap_{T} / p_{T,2}");
362 fHistDeltaPtOverJet2PtvsJet2Pt->GetZaxis()->SetTitle("counts");
363 fOutput->Add(fHistDeltaPtOverJet2PtvsJet2Pt);
364
9adcb46d 365 fHistDeltaPtvsDistance = new TH2F("fHistDeltaPtvsDistance", "fHistDeltaPtvsDistance",
366 fNbins/2, 0, 1.2, 2*fNbins, -fMaxBinPt, fMaxBinPt);
05077f28 367 fHistDeltaPtvsDistance->GetXaxis()->SetTitle("Distance");
368 fHistDeltaPtvsDistance->GetYaxis()->SetTitle("#deltap_{T} (GeV/c)");
369 fHistDeltaPtvsDistance->GetZaxis()->SetTitle("counts");
370 fOutput->Add(fHistDeltaPtvsDistance);
371
9adcb46d 372 fHistDeltaPtvsCommonEnergy1 = new TH2F("fHistDeltaPtvsCommonEnergy1", "fHistDeltaPtvsCommonEnergy1",
373 fNbins/2, 0, 1.2, 2*fNbins, -fMaxBinPt, fMaxBinPt);
05077f28 374 fHistDeltaPtvsCommonEnergy1->GetXaxis()->SetTitle("Common energy 1 (%)");
375 fHistDeltaPtvsCommonEnergy1->GetYaxis()->SetTitle("#deltap_{T} (GeV/c)");
376 fHistDeltaPtvsCommonEnergy1->GetZaxis()->SetTitle("counts");
377 fOutput->Add(fHistDeltaPtvsCommonEnergy1);
378
9adcb46d 379 fHistDeltaPtvsCommonEnergy2 = new TH2F("fHistDeltaPtvsCommonEnergy2", "fHistDeltaPtvsCommonEnergy2",
380 fNbins/2, 0, 1.2, 2*fNbins, -fMaxBinPt, fMaxBinPt);
05077f28 381 fHistDeltaPtvsCommonEnergy2->GetXaxis()->SetTitle("Common energy 2 (%)");
382 fHistDeltaPtvsCommonEnergy2->GetYaxis()->SetTitle("#deltap_{T} (GeV/c)");
383 fHistDeltaPtvsCommonEnergy2->GetZaxis()->SetTitle("counts");
384 fOutput->Add(fHistDeltaPtvsCommonEnergy2);
385
9adcb46d 386 fHistDeltaPtvsArea1 = new TH2F("fHistDeltaPtvsArea1", "fHistDeltaPtvsArea1",
faf21244 387 fNbins/2, 0, 1.5, 2*fNbins, -fMaxBinPt, fMaxBinPt);
05077f28 388 fHistDeltaPtvsArea1->GetXaxis()->SetTitle("A_{jet,1}");
389 fHistDeltaPtvsArea1->GetYaxis()->SetTitle("#deltap_{T} (GeV/c)");
390 fHistDeltaPtvsArea1->GetZaxis()->SetTitle("counts");
391 fOutput->Add(fHistDeltaPtvsArea1);
392
9adcb46d 393 fHistDeltaPtvsArea2 = new TH2F("fHistDeltaPtvsArea2", "fHistDeltaPtvsArea2",
faf21244 394 fNbins/2, 0, 1.5, 2*fNbins, -fMaxBinPt, fMaxBinPt);
05077f28 395 fHistDeltaPtvsArea2->GetXaxis()->SetTitle("A_{jet,2}");
396 fHistDeltaPtvsArea2->GetYaxis()->SetTitle("#deltap_{T} (GeV/c)");
397 fHistDeltaPtvsArea2->GetZaxis()->SetTitle("counts");
398 fOutput->Add(fHistDeltaPtvsArea2);
399
400 fHistDeltaPtvsDeltaArea = new TH2F("fHistDeltaPtvsDeltaArea", "fHistDeltaPtvsDeltaArea",
43032ce2 401 fNbins, -1.+1./fNbins, 1.+1./fNbins, 2*fNbins, -fMaxBinPt, fMaxBinPt);
05077f28 402 fHistDeltaPtvsDeltaArea->GetXaxis()->SetTitle("#deltaA_{jet}");
403 fHistDeltaPtvsDeltaArea->GetYaxis()->SetTitle("#deltap_{T} (GeV/c)");
404 fHistDeltaPtvsDeltaArea->GetZaxis()->SetTitle("counts");
405 fOutput->Add(fHistDeltaPtvsDeltaArea);
406
407 fHistJet1PtvsJet2Pt = new TH2F("fHistJet1PtvsJet2Pt", "fHistJet1PtvsJet2Pt", fNbins, fMinBinPt, fMaxBinPt, fNbins, fMinBinPt, fMaxBinPt);
408 fHistJet1PtvsJet2Pt->GetXaxis()->SetTitle("p_{T,1}");
409 fHistJet1PtvsJet2Pt->GetYaxis()->SetTitle("p_{T,2}");
410 fHistJet1PtvsJet2Pt->GetZaxis()->SetTitle("counts");
411 fOutput->Add(fHistJet1PtvsJet2Pt);
cd6431de 412
7cd832c7 413 if (fIsJet1Rho || fIsJet2Rho) {
414 if (!fIsJet1Rho)
5c045cde 415 fHistDeltaCorrPtvsJet1CorrPt = new TH2F("fHistDeltaCorrPtvsJet1CorrPt", "fHistDeltaCorrPtvsJet1CorrPt",
416 fNbins, fMinBinPt, fMaxBinPt, 2*fNbins, -fMaxBinPt, fMaxBinPt);
417 else
418 fHistDeltaCorrPtvsJet1CorrPt = new TH2F("fHistDeltaCorrPtvsJet1CorrPt", "fHistDeltaCorrPtvsJet1CorrPt",
419 2*fNbins, -fMaxBinPt, fMaxBinPt, 2*fNbins, -fMaxBinPt, fMaxBinPt);
420
4c00a7af 421 fHistDeltaCorrPtvsJet1CorrPt->GetXaxis()->SetTitle("p_{T,1}^{corr}");
422 fHistDeltaCorrPtvsJet1CorrPt->GetYaxis()->SetTitle("#deltap_{T}^{corr} (GeV/c)");
423 fHistDeltaCorrPtvsJet1CorrPt->GetZaxis()->SetTitle("counts");
424 fOutput->Add(fHistDeltaCorrPtvsJet1CorrPt);
425
7cd832c7 426 if (!fIsJet2Rho)
5c045cde 427 fHistDeltaCorrPtvsJet2CorrPt = new TH2F("fHistDeltaCorrPtvsJet2CorrPt", "fHistDeltaCorrPtvsJet2CorrPt",
428 fNbins, fMinBinPt, fMaxBinPt, 2*fNbins, -fMaxBinPt, fMaxBinPt);
429 else
430 fHistDeltaCorrPtvsJet2CorrPt = new TH2F("fHistDeltaCorrPtvsJet2CorrPt", "fHistDeltaCorrPtvsJet2CorrPt",
431 2*fNbins, -fMaxBinPt, fMaxBinPt, 2*fNbins, -fMaxBinPt, fMaxBinPt);
432
4c00a7af 433 fHistDeltaCorrPtvsJet2CorrPt->GetXaxis()->SetTitle("p_{T,2}^{corr}");
434 fHistDeltaCorrPtvsJet2CorrPt->GetYaxis()->SetTitle("#deltap_{T}^{corr} (GeV/c)");
435 fHistDeltaCorrPtvsJet2CorrPt->GetZaxis()->SetTitle("counts");
436 fOutput->Add(fHistDeltaCorrPtvsJet2CorrPt);
437
438 fHistDeltaCorrPtOverJet1CorrPtvsJet1CorrPt = new TH2F("fHistDeltaCorrPtOverJet1CorrPtvsJet1CorrPt", "fHistDeltaCorrPtOverJet1CorrPtvsJet1CorrPt",
439 2*fNbins, -fMaxBinPt, fMaxBinPt, fNbins, -5, 5);
440 fHistDeltaCorrPtOverJet1CorrPtvsJet1CorrPt->GetXaxis()->SetTitle("p_{T,1}^{corr}");
441 fHistDeltaCorrPtOverJet1CorrPtvsJet1CorrPt->GetYaxis()->SetTitle("#deltap_{T}^{corr} / p_{T,1}^{corr}");
442 fHistDeltaCorrPtOverJet1CorrPtvsJet1CorrPt->GetZaxis()->SetTitle("counts");
443 fOutput->Add(fHistDeltaCorrPtOverJet1CorrPtvsJet1CorrPt);
444
445 fHistDeltaCorrPtOverJet2CorrPtvsJet2CorrPt = new TH2F("fHistDeltaCorrPtOverJet2CorrPtvsJet2CorrPt", "fHistDeltaCorrPtOverJet2CorrPtvsJet2CorrPt",
446 2*fNbins, -fMaxBinPt, fMaxBinPt, fNbins, -5, 5);
447 fHistDeltaCorrPtOverJet2CorrPtvsJet2CorrPt->GetXaxis()->SetTitle("p_{T,2}^{corr}");
448 fHistDeltaCorrPtOverJet2CorrPtvsJet2CorrPt->GetYaxis()->SetTitle("#deltap_{T}^{corr} / p_{T,2}^{corr}");
449 fHistDeltaCorrPtOverJet2CorrPtvsJet2CorrPt->GetZaxis()->SetTitle("counts");
450 fOutput->Add(fHistDeltaCorrPtOverJet2CorrPtvsJet2CorrPt);
cfc2ac24 451
9adcb46d 452 fHistDeltaCorrPtvsDistance = new TH2F("fHistDeltaCorrPtvsDistance", "fHistDeltaCorrPtvsDistance",
453 fNbins/2, 0, 1.2, 2*fNbins, -fMaxBinPt, fMaxBinPt);
05077f28 454 fHistDeltaCorrPtvsDistance->GetXaxis()->SetTitle("Distance");
455 fHistDeltaCorrPtvsDistance->GetYaxis()->SetTitle("#deltap_{T}^{corr} (GeV/c)");
456 fHistDeltaCorrPtvsDistance->GetZaxis()->SetTitle("counts");
457 fOutput->Add(fHistDeltaCorrPtvsDistance);
458
9adcb46d 459 fHistDeltaCorrPtvsCommonEnergy1 = new TH2F("fHistDeltaCorrPtvsCommonEnergy1", "fHistDeltaCorrPtvsCommonEnergy1",
460 fNbins/2, 0, 1.2, 2*fNbins, -fMaxBinPt, fMaxBinPt);
05077f28 461 fHistDeltaCorrPtvsCommonEnergy1->GetXaxis()->SetTitle("Common energy 1 (%)");
462 fHistDeltaCorrPtvsCommonEnergy1->GetYaxis()->SetTitle("#deltap_{T}^{corr} (GeV/c)");
463 fHistDeltaCorrPtvsCommonEnergy1->GetZaxis()->SetTitle("counts");
464 fOutput->Add(fHistDeltaCorrPtvsCommonEnergy1);
465
9adcb46d 466 fHistDeltaCorrPtvsCommonEnergy2 = new TH2F("fHistDeltaCorrPtvsCommonEnergy2", "fHistDeltaCorrPtvsCommonEnergy2",
467 fNbins/2, 0, 1.2, 2*fNbins, -fMaxBinPt, fMaxBinPt);
05077f28 468 fHistDeltaCorrPtvsCommonEnergy2->GetXaxis()->SetTitle("Common energy 2 (%)");
469 fHistDeltaCorrPtvsCommonEnergy2->GetYaxis()->SetTitle("#deltap_{T}^{corr} (GeV/c)");
470 fHistDeltaCorrPtvsCommonEnergy2->GetZaxis()->SetTitle("counts");
471 fOutput->Add(fHistDeltaCorrPtvsCommonEnergy2);
472
9adcb46d 473 fHistDeltaCorrPtvsArea1 = new TH2F("fHistDeltaCorrPtvsArea1", "fHistDeltaCorrPtvsArea1",
faf21244 474 fNbins/2, 0, 1.5, 2*fNbins, -fMaxBinPt, fMaxBinPt);
05077f28 475 fHistDeltaCorrPtvsArea1->GetXaxis()->SetTitle("A_{jet,1}");
476 fHistDeltaCorrPtvsArea1->GetYaxis()->SetTitle("#deltap_{T}^{corr} (GeV/c)");
477 fHistDeltaCorrPtvsArea1->GetZaxis()->SetTitle("counts");
478 fOutput->Add(fHistDeltaCorrPtvsArea1);
479
9adcb46d 480 fHistDeltaCorrPtvsArea2 = new TH2F("fHistDeltaCorrPtvsArea2", "fHistDeltaCorrPtvsArea2",
faf21244 481 fNbins/2, 0, 1.5, 2*fNbins, -fMaxBinPt, fMaxBinPt);
05077f28 482 fHistDeltaCorrPtvsArea2->GetXaxis()->SetTitle("A_{jet,2}");
483 fHistDeltaCorrPtvsArea2->GetYaxis()->SetTitle("#deltap_{T}^{corr} (GeV/c)");
484 fHistDeltaCorrPtvsArea2->GetZaxis()->SetTitle("counts");
485 fOutput->Add(fHistDeltaCorrPtvsArea2);
486
487 fHistDeltaCorrPtvsDeltaArea = new TH2F("fHistDeltaCorrPtvsDeltaArea", "fHistDeltaCorrPtvsDeltaArea",
43032ce2 488 fNbins, -1.+1./fNbins, 1.+1./fNbins, 2*fNbins, -fMaxBinPt, fMaxBinPt);
05077f28 489 fHistDeltaCorrPtvsDeltaArea->GetXaxis()->SetTitle("#deltaA_{jet}");
490 fHistDeltaCorrPtvsDeltaArea->GetYaxis()->SetTitle("#deltap_{T}^{corr} (GeV/c)");
491 fHistDeltaCorrPtvsDeltaArea->GetZaxis()->SetTitle("counts");
492 fOutput->Add(fHistDeltaCorrPtvsDeltaArea);
493
7cd832c7 494 if (!fIsJet1Rho)
9adcb46d 495 fHistJet1CorrPtvsJet2CorrPt = new TH2F("fHistJet1CorrPtvsJet2CorrPt", "fHistJet1CorrPtvsJet2CorrPt",
496 fNbins, fMinBinPt, fMaxBinPt, 2*fNbins, -fMaxBinPt, fMaxBinPt);
7cd832c7 497 else if (!fIsJet2Rho)
9adcb46d 498 fHistJet1CorrPtvsJet2CorrPt = new TH2F("fHistJet1CorrPtvsJet2CorrPt", "fHistJet1CorrPtvsJet2CorrPt",
499 2*fNbins, -fMaxBinPt, fMaxBinPt, fNbins, fMinBinPt, fMaxBinPt);
cd6431de 500 else
9adcb46d 501 fHistJet1CorrPtvsJet2CorrPt = new TH2F("fHistJet1CorrPtvsJet2CorrPt", "fHistJet1CorrPtvsJet2CorrPt",
502 2*fNbins, -fMaxBinPt, fMaxBinPt,
503 2*fNbins, -fMaxBinPt, fMaxBinPt);
cd6431de 504 fHistJet1CorrPtvsJet2CorrPt->GetXaxis()->SetTitle("p_{T,1}^{corr}");
505 fHistJet1CorrPtvsJet2CorrPt->GetYaxis()->SetTitle("p_{T,2}^{corr}");
506 fHistJet1CorrPtvsJet2CorrPt->GetZaxis()->SetTitle("counts");
507 fOutput->Add(fHistJet1CorrPtvsJet2CorrPt);
508 }
509
4358e58a 510 if (fIsEmbedded) {
4c00a7af 511 fHistDeltaMCPtvsJet1MCPt = new TH2F("fHistDeltaMCPtvsJet1MCPt", "fHistDeltaMCPtvsJet1MCPt",
43032ce2 512 fNbins, fMinBinPt, fMaxBinPt, 2*fNbins, -fMaxBinPt, fMaxBinPt);
4c00a7af 513 fHistDeltaMCPtvsJet1MCPt->GetXaxis()->SetTitle("p_{T,1}^{MC}");
514 fHistDeltaMCPtvsJet1MCPt->GetYaxis()->SetTitle("#deltap_{T} (GeV/c)");
515 fHistDeltaMCPtvsJet1MCPt->GetZaxis()->SetTitle("counts");
516 fOutput->Add(fHistDeltaMCPtvsJet1MCPt);
05077f28 517
9adcb46d 518 fHistDeltaMCPtvsJet2Pt = new TH2F("fHistDeltaMCPtvsJet2Pt", "fHistDeltaMCPtvsJet2Pt",
519 fNbins, fMinBinPt, fMaxBinPt, 2*fNbins, -fMaxBinPt, fMaxBinPt);
05077f28 520 fHistDeltaMCPtvsJet2Pt->GetXaxis()->SetTitle("p_{T,2}");
521 fHistDeltaMCPtvsJet2Pt->GetYaxis()->SetTitle("#deltap_{T} (GeV/c)");
522 fHistDeltaMCPtvsJet2Pt->GetZaxis()->SetTitle("counts");
523 fOutput->Add(fHistDeltaMCPtvsJet2Pt);
524
4c00a7af 525 fHistDeltaMCPtOverJet1MCPtvsJet1MCPt = new TH2F("fHistDeltaMCPtOverJet1MCPtvsJet1MCPt", "fHistDeltaMCPtOverJet1MCPtvsJet1MCPt",
526 fNbins, fMinBinPt, fMaxBinPt, fNbins, -5, 5);
527 fHistDeltaMCPtOverJet1MCPtvsJet1MCPt->GetXaxis()->SetTitle("p_{T,1}^{MC}");
528 fHistDeltaMCPtOverJet1MCPtvsJet1MCPt->GetYaxis()->SetTitle("#deltap_{T} / p_{T,1}^{MC}");
529 fHistDeltaMCPtOverJet1MCPtvsJet1MCPt->GetZaxis()->SetTitle("counts");
530 fOutput->Add(fHistDeltaMCPtOverJet1MCPtvsJet1MCPt);
531
532 fHistDeltaMCPtOverJet2PtvsJet2Pt = new TH2F("fHistDeltaMCPtOverJet2PtvsJet2Pt", "fHistDeltaMCPtOverJet2PtvsJet2Pt",
533 fNbins, fMinBinPt, fMaxBinPt, fNbins, -5, 5);
534 fHistDeltaMCPtOverJet2PtvsJet2Pt->GetXaxis()->SetTitle("p_{T,2}");
535 fHistDeltaMCPtOverJet2PtvsJet2Pt->GetYaxis()->SetTitle("#deltap_{T} / p_{T,2}");
536 fHistDeltaMCPtOverJet2PtvsJet2Pt->GetZaxis()->SetTitle("counts");
537 fOutput->Add(fHistDeltaMCPtOverJet2PtvsJet2Pt);
538
9adcb46d 539 fHistDeltaMCPtvsDistance = new TH2F("fHistDeltaMCPtvsDistance", "fHistDeltaMCPtvsDistance",
540 fNbins/2, 0, 1.2, 2*fNbins, -fMaxBinPt, fMaxBinPt);
05077f28 541 fHistDeltaMCPtvsDistance->GetXaxis()->SetTitle("Distance");
542 fHistDeltaMCPtvsDistance->GetYaxis()->SetTitle("#deltap_{T} (GeV/c)");
543 fHistDeltaMCPtvsDistance->GetZaxis()->SetTitle("counts");
544 fOutput->Add(fHistDeltaMCPtvsDistance);
545
9adcb46d 546 fHistDeltaMCPtvsCommonEnergy1 = new TH2F("fHistDeltaMCPtvsCommonEnergy1", "fHistDeltaMCPtvsCommonEnergy1",
547 fNbins/2, 0, 1.2, 2*fNbins, -fMaxBinPt, fMaxBinPt);
05077f28 548 fHistDeltaMCPtvsCommonEnergy1->GetXaxis()->SetTitle("Common energy 1 (%)");
549 fHistDeltaMCPtvsCommonEnergy1->GetYaxis()->SetTitle("#deltap_{T} (GeV/c)");
550 fHistDeltaMCPtvsCommonEnergy1->GetZaxis()->SetTitle("counts");
551 fOutput->Add(fHistDeltaMCPtvsCommonEnergy1);
552
9adcb46d 553 fHistDeltaMCPtvsCommonEnergy2 = new TH2F("fHistDeltaMCPtvsCommonEnergy2", "fHistDeltaMCPtvsCommonEnergy2",
554 fNbins/2, 0, 1.2, 2*fNbins, -fMaxBinPt, fMaxBinPt);
05077f28 555 fHistDeltaMCPtvsCommonEnergy2->GetXaxis()->SetTitle("Common energy 2 (%)");
556 fHistDeltaMCPtvsCommonEnergy2->GetYaxis()->SetTitle("#deltap_{T} (GeV/c)");
557 fHistDeltaMCPtvsCommonEnergy2->GetZaxis()->SetTitle("counts");
558 fOutput->Add(fHistDeltaMCPtvsCommonEnergy2);
559
9adcb46d 560 fHistDeltaMCPtvsArea1 = new TH2F("fHistDeltaMCPtvsArea1", "fHistDeltaMCPtvsArea1",
faf21244 561 fNbins/2, 0, 1.5, 2*fNbins, -fMaxBinPt, fMaxBinPt);
05077f28 562 fHistDeltaMCPtvsArea1->GetXaxis()->SetTitle("A_{jet,1}");
563 fHistDeltaMCPtvsArea1->GetYaxis()->SetTitle("#deltap_{T} (GeV/c)");
564 fHistDeltaMCPtvsArea1->GetZaxis()->SetTitle("counts");
565 fOutput->Add(fHistDeltaMCPtvsArea1);
566
9adcb46d 567 fHistDeltaMCPtvsArea2 = new TH2F("fHistDeltaMCPtvsArea2", "fHistDeltaMCPtvsArea2",
faf21244 568 fNbins/2, 0, 1.5, 2*fNbins, -fMaxBinPt, fMaxBinPt);
05077f28 569 fHistDeltaMCPtvsArea2->GetXaxis()->SetTitle("A_{jet,2}");
570 fHistDeltaMCPtvsArea2->GetYaxis()->SetTitle("#deltap_{T} (GeV/c)");
571 fHistDeltaMCPtvsArea2->GetZaxis()->SetTitle("counts");
572 fOutput->Add(fHistDeltaMCPtvsArea2);
573
574 fHistDeltaMCPtvsDeltaArea = new TH2F("fHistDeltaMCPtvsDeltaArea", "fHistDeltaMCPtvsDeltaArea",
43032ce2 575 fNbins, -1.+1./fNbins, 1.+1./fNbins, 2*fNbins, -fMaxBinPt, fMaxBinPt);
05077f28 576 fHistDeltaMCPtvsDeltaArea->GetXaxis()->SetTitle("#deltaA_{jet}");
577 fHistDeltaMCPtvsDeltaArea->GetYaxis()->SetTitle("#deltap_{T} (GeV/c)");
578 fHistDeltaMCPtvsDeltaArea->GetZaxis()->SetTitle("counts");
579 fOutput->Add(fHistDeltaMCPtvsDeltaArea);
580
9adcb46d 581 fHistJet1MCPtvsJet2Pt = new TH2F("fHistJet1MCPtvsJet2Pt", "fHistJet1MCPtvsJet2Pt",
582 fNbins, fMinBinPt, fMaxBinPt, fNbins, fMinBinPt, fMaxBinPt);
4358e58a 583 fHistJet1MCPtvsJet2Pt->GetXaxis()->SetTitle("p_{T,1}^{MC}");
584 fHistJet1MCPtvsJet2Pt->GetYaxis()->SetTitle("p_{T,2}");
585 fHistJet1MCPtvsJet2Pt->GetZaxis()->SetTitle("counts");
586 fOutput->Add(fHistJet1MCPtvsJet2Pt);
587 }
43032ce2 588}
589
590//________________________________________________________________________
591void AliJetResponseMaker::AllocateTHnSparse()
592{
593 // Allocate THnSparse histograms.
594
69a6bda0 595 TString title[25]= {""};
596 Int_t nbins[25] = {0};
597 Double_t min[25] = {0.};
598 Double_t max[25] = {0.};
43032ce2 599 Int_t dim = 0;
600
601 title[dim] = "#phi";
602 nbins[dim] = fNbins/4;
603 min[dim] = 0;
604 max[dim] = 2*TMath::Pi()*(1 + 1./(nbins[dim]-1));
605 dim++;
606
607 title[dim] = "#eta";
608 nbins[dim] = fNbins/4;
609 min[dim] = -1;
610 max[dim] = 1;
611 dim++;
612
613 title[dim] = "p_{T}";
614 nbins[dim] = fNbins;
615 min[dim] = 0;
616 max[dim] = 250;
617 dim++;
618
619 title[dim] = "A_{jet}";
620 nbins[dim] = fNbins/4;
621 min[dim] = 0;
faf21244 622 max[dim] = 1.5;
43032ce2 623 dim++;
624
d3c0e12a 625 if (fNEFAxis) {
626 title[dim] = "NEF";
627 nbins[dim] = fNbins/4;
628 min[dim] = 0;
629 max[dim] = 1.2;
630 dim++;
631 }
43032ce2 632
d3c0e12a 633 if (fZAxis) {
634 title[dim] = "Z";
635 nbins[dim] = fNbins/4;
636 min[dim] = 0;
637 max[dim] = 1.2;
638 dim++;
639 }
640
641 title[dim] = "p_{T,particle}^{leading} (GeV/c)";
642 nbins[dim] = 120;
43032ce2 643 min[dim] = 0;
d3c0e12a 644 max[dim] = 120;
43032ce2 645 dim++;
646
647 Int_t dim1 = dim, dim2 = dim;
648
7cd832c7 649 if (fIsJet1Rho) {
43032ce2 650 title[dim1] = "p_{T}^{corr}";
651 nbins[dim1] = fNbins*2;
652 min[dim1] = -250;
653 max[dim1] = 250;
654 dim1++;
655 }
656
657 if (fIsEmbedded) {
658 title[dim1] = "p_{T}^{MC}";
659 nbins[dim1] = fNbins;
660 min[dim1] = 0;
661 max[dim1] = 250;
662 dim1++;
663 }
664
665 fHistJets1 = new THnSparseD("fHistJets1","fHistJets1",dim1,nbins,min,max);
666 for (Int_t i = 0; i < dim1; i++)
667 fHistJets1->GetAxis(i)->SetTitle(title[i]);
668 fOutput->Add(fHistJets1);
669
7cd832c7 670 if (fIsJet2Rho) {
43032ce2 671 title[dim2] = "p_{T}^{corr}";
672 nbins[dim2] = fNbins*2;
673 min[dim2] = -250;
674 max[dim2] = 250;
675 dim2++;
676 }
677
d3c0e12a 678 fHistJets2 = new THnSparseD("fHistJets2","fHistJets2",dim2,nbins,min,max);
43032ce2 679 for (Int_t i = 0; i < dim2; i++)
d3c0e12a 680 fHistJets2->GetAxis(i)->SetTitle(title[i]);
681 fOutput->Add(fHistJets2);
43032ce2 682
683 // Matching
684
685 dim = 0;
686
687 title[dim] = "p_{T,1}";
688 nbins[dim] = fNbins;
689 min[dim] = 0;
690 max[dim] = 250;
691 dim++;
692
693 title[dim] = "p_{T,2}";
694 nbins[dim] = fNbins;
695 min[dim] = 0;
696 max[dim] = 250;
697 dim++;
698
699 title[dim] = "A_{jet,1}";
700 nbins[dim] = fNbins/4;
701 min[dim] = 0;
faf21244 702 max[dim] = 1.5;
43032ce2 703 dim++;
704
705 title[dim] = "A_{jet,2}";
706 nbins[dim] = fNbins/4;
707 min[dim] = 0;
faf21244 708 max[dim] = 1.5;
43032ce2 709 dim++;
710
711 title[dim] = "distance";
712 nbins[dim] = fNbins/2;
713 min[dim] = 0;
714 max[dim] = 1.2;
715 dim++;
716
717 title[dim] = "CE1";
718 nbins[dim] = fNbins/2;
719 min[dim] = 0;
720 max[dim] = 1.2;
721 dim++;
722
723 title[dim] = "CE2";
724 nbins[dim] = fNbins/2;
725 min[dim] = 0;
726 max[dim] = 1.2;
727 dim++;
728
d3c0e12a 729 title[dim] = "p_{T,particle,1}^{leading} (GeV/c)";
730 nbins[dim] = 120;
731 min[dim] = 0;
732 max[dim] = 120;
733 dim++;
734
735 title[dim] = "p_{T,particle,2}^{leading} (GeV/c)";
736 nbins[dim] = 120;
737 min[dim] = 0;
738 max[dim] = 120;
739 dim++;
740
43032ce2 741 if (fDeltaPtAxis) {
742 title[dim] = "#deltaA_{jet}";
743 nbins[dim] = fNbins/2;
744 min[dim] = -1;
745 max[dim] = 1;
746 dim++;
747
748 title[dim] = "#deltap_{T}";
749 nbins[dim] = fNbins*2;
750 min[dim] = -250;
751 max[dim] = 250;
752 dim++;
753 }
7cd832c7 754 if (fIsJet1Rho) {
43032ce2 755 title[dim] = "p_{T,1}^{corr}";
756 nbins[dim] = fNbins*2;
757 min[dim] = -250;
758 max[dim] = 250;
759 dim++;
760 }
7cd832c7 761 if (fIsJet2Rho) {
43032ce2 762 title[dim] = "p_{T,2}^{corr}";
763 nbins[dim] = fNbins*2;
764 min[dim] = -250;
765 max[dim] = 250;
766 dim++;
767 }
7cd832c7 768 if (fDeltaPtAxis && (fIsJet1Rho || fIsJet2Rho)) {
43032ce2 769 title[dim] = "#deltap_{T}^{corr}";
770 nbins[dim] = fNbins*2;
771 min[dim] = -250;
772 max[dim] = 250;
773 dim++;
774 }
775 if (fDeltaEtaDeltaPhiAxis) {
776 title[dim] = "#delta#eta";
777 nbins[dim] = fNbins/2;
778 min[dim] = -1;
779 max[dim] = 1;
780 dim++;
781
782 title[dim] = "#delta#phi";
783 nbins[dim] = fNbins/2;
784 min[dim] = -TMath::Pi()/2;
785 max[dim] = TMath::Pi()*3/2;
786 dim++;
787 }
788 if (fIsEmbedded) {
789 title[dim] = "p_{T,1}^{MC}";
790 nbins[dim] = fNbins;
791 min[dim] = 0;
792 max[dim] = 250;
793 dim++;
794
795 if (fDeltaPtAxis) {
796 title[dim] = "#deltap_{T}^{MC}";
797 nbins[dim] = fNbins*2;
798 min[dim] = -250;
799 max[dim] = 250;
800 dim++;
801 }
802 }
faf21244 803
804 if (fNEFAxis) {
805 title[dim] = "NEF_{1}";
806 nbins[dim] = fNbins/4;
807 min[dim] = 0;
808 max[dim] = 1.2;
809 dim++;
810
811 title[dim] = "NEF_{2}";
812 nbins[dim] = fNbins/4;
813 min[dim] = 0;
814 max[dim] = 1.2;
815 dim++;
816 }
817
818 if (fZAxis) {
819 title[dim] = "Z_{1}";
820 nbins[dim] = fNbins/4;
821 min[dim] = 0;
822 max[dim] = 1.2;
823 dim++;
824
825 title[dim] = "Z_{2}";
826 nbins[dim] = fNbins/4;
827 min[dim] = 0;
828 max[dim] = 1.2;
829 dim++;
830 }
43032ce2 831
832 fHistMatching = new THnSparseD("fHistMatching","fHistMatching",dim,nbins,min,max);
833
834 for (Int_t i = 0; i < dim; i++)
835 fHistMatching->GetAxis(i)->SetTitle(title[i]);
836
837 fOutput->Add(fHistMatching);
838}
839
840//________________________________________________________________________
841void AliJetResponseMaker::UserCreateOutputObjects()
842{
843 // Create user objects.
844
9239b066 845 AliAnalysisTaskEmcalJet::UserCreateOutputObjects();
43032ce2 846
7cd832c7 847 AliJetContainer *jets1 = static_cast<AliJetContainer*>(fJetCollArray.At(0));
848 AliJetContainer *jets2 = static_cast<AliJetContainer*>(fJetCollArray.At(1));
849
850 if (!jets1 || !jets2) return;
851
852 if (jets1->GetRhoName().IsNull()) fIsJet1Rho = kFALSE;
853 else fIsJet1Rho = kTRUE;
854 if (jets2->GetRhoName().IsNull()) fIsJet2Rho = kFALSE;
855 else fIsJet2Rho = kTRUE;
43032ce2 856
857 if (fHistoType==0)
858 AllocateTH2();
859 else
860 AllocateTHnSparse();
4358e58a 861
99c67c1b 862 PostData(1, fOutput); // Post data for ALL output slots > 0 here, to get at least an empty histogram
2949a2ec 863}
864
43032ce2 865//________________________________________________________________________
d3c0e12a 866void AliJetResponseMaker::FillJetHisto(Double_t Phi, Double_t Eta, Double_t Pt, Double_t A, Double_t NEF, Double_t LeadingPt, Double_t CorrPt, Double_t MCPt, Int_t Set)
43032ce2 867{
868 if (fHistoType==1) {
869 THnSparse *histo = 0;
870 if (Set==1)
871 histo = fHistJets1;
872 else if (Set==2)
873 histo = fHistJets2;
43032ce2 874
875 if (!histo)
876 return;
877
878 Double_t contents[20]={0};
879
880 for (Int_t i = 0; i < histo->GetNdimensions(); i++) {
881 TString title(histo->GetAxis(i)->GetTitle());
882 if (title=="#phi")
883 contents[i] = Phi;
884 else if (title=="#eta")
885 contents[i] = Eta;
886 else if (title=="p_{T}")
887 contents[i] = Pt;
888 else if (title=="A_{jet}")
889 contents[i] = A;
890 else if (title=="NEF")
891 contents[i] = NEF;
892 else if (title=="Z")
d3c0e12a 893 contents[i] = LeadingPt/Pt;
43032ce2 894 else if (title=="p_{T}^{corr}")
895 contents[i] = CorrPt;
896 else if (title=="p_{T}^{MC}")
897 contents[i] = MCPt;
d3c0e12a 898 else if (title=="p_{T,particle}^{leading} (GeV/c)")
899 contents[i] = LeadingPt;
43032ce2 900 else
901 AliWarning(Form("Unable to fill dimension %s!",title.Data()));
902 }
903
904 histo->Fill(contents);
905 }
906 else {
907 if (Set == 1) {
908 fHistJets1PtArea->Fill(A, Pt);
909 fHistJets1PhiEta->Fill(Eta, Phi);
d3c0e12a 910 fHistJets1ZvsPt->Fill(LeadingPt/Pt, Pt);
43032ce2 911 fHistJets1NEFvsPt->Fill(NEF, Pt);
912 fHistJets1CEFvsCEFPt->Fill(1-NEF, (1-NEF)*Pt);
7cd832c7 913 if (fIsJet1Rho)
43032ce2 914 fHistJets1CorrPtArea->Fill(A, CorrPt);
915 }
916 else if (Set == 2) {
917 fHistJets2PtArea->Fill(A, Pt);
918 fHistJets2PhiEta->Fill(Eta, Phi);
d3c0e12a 919 fHistJets2ZvsPt->Fill(LeadingPt/Pt, Pt);
43032ce2 920 fHistJets2NEFvsPt->Fill(NEF, Pt);
921 fHistJets2CEFvsCEFPt->Fill(1-NEF, (1-NEF)*Pt);
7cd832c7 922 if (fIsJet2Rho)
e8751ce9 923 fHistJets2CorrPtArea->Fill(A, CorrPt);
43032ce2 924 }
925 }
926}
927
928//________________________________________________________________________
faf21244 929void AliJetResponseMaker::FillMatchingHistos(Double_t Pt1, Double_t Pt2, Double_t Eta1, Double_t Eta2, Double_t Phi1, Double_t Phi2,
930 Double_t A1, Double_t A2, Double_t d, Double_t CE1, Double_t CE2, Double_t CorrPt1, Double_t CorrPt2,
d3c0e12a 931 Double_t MCPt1, Double_t NEF1, Double_t NEF2, Double_t LeadingPt1, Double_t LeadingPt2)
43032ce2 932{
933 if (fHistoType==1) {
934 Double_t contents[20]={0};
935
936 for (Int_t i = 0; i < fHistMatching->GetNdimensions(); i++) {
937 TString title(fHistMatching->GetAxis(i)->GetTitle());
938 if (title=="p_{T,1}")
939 contents[i] = Pt1;
940 else if (title=="p_{T,2}")
941 contents[i] = Pt2;
942 else if (title=="A_{jet,1}")
943 contents[i] = A1;
944 else if (title=="A_{jet,2}")
945 contents[i] = A2;
946 else if (title=="distance")
947 contents[i] = d;
948 else if (title=="CE1")
949 contents[i] = CE1;
950 else if (title=="CE2")
951 contents[i] = CE2;
952 else if (title=="#deltaA_{jet}")
953 contents[i] = A1-A2;
954 else if (title=="#deltap_{T}")
955 contents[i] = Pt1-Pt2;
956 else if (title=="#delta#eta")
957 contents[i] = Eta1-Eta2;
958 else if (title=="#delta#phi")
959 contents[i] = Phi1-Phi2;
960 else if (title=="p_{T,1}^{corr}")
961 contents[i] = CorrPt1;
962 else if (title=="p_{T,2}^{corr}")
963 contents[i] = CorrPt2;
964 else if (title=="#deltap_{T}^{corr}")
965 contents[i] = CorrPt1-CorrPt2;
966 else if (title=="p_{T,1}^{MC}")
967 contents[i] = MCPt1;
968 else if (title=="#deltap_{T}^{MC}")
969 contents[i] = MCPt1-Pt2;
faf21244 970 else if (title=="NEF_{1}")
971 contents[i] = NEF1;
972 else if (title=="NEF_{2}")
973 contents[i] = NEF2;
974 else if (title=="Z_{1}")
d3c0e12a 975 contents[i] = LeadingPt1/Pt1;
faf21244 976 else if (title=="Z_{2}")
d3c0e12a 977 contents[i] = LeadingPt2/Pt2;
978 else if (title=="p_{T,particle,1}^{leading} (GeV/c)")
979 contents[i] = LeadingPt1;
980 else if (title=="p_{T,particle,2}^{leading} (GeV/c)")
981 contents[i] = LeadingPt2;
43032ce2 982 else
983 AliWarning(Form("Unable to fill dimension %s!",title.Data()));
984 }
985
986 fHistMatching->Fill(contents);
987 }
988 else {
989 fHistCommonEnergy1vsJet1Pt->Fill(CE1, Pt1);
990 fHistCommonEnergy2vsJet2Pt->Fill(CE2, Pt2);
991
992 fHistDistancevsJet1Pt->Fill(d, Pt1);
993 fHistDistancevsJet2Pt->Fill(d, Pt2);
994
995 fHistDistancevsCommonEnergy1->Fill(d, CE1);
996 fHistDistancevsCommonEnergy2->Fill(d, CE2);
997 fHistCommonEnergy1vsCommonEnergy2->Fill(CE1, CE2);
998
999 fHistDeltaEtaDeltaPhi->Fill(Eta1-Eta2,Phi1-Phi2);
1000
1001 fHistJet2PtOverJet1PtvsJet2Pt->Fill(Pt2, Pt2 / Pt1);
1002 fHistJet1PtOverJet2PtvsJet1Pt->Fill(Pt1, Pt1 / Pt2);
1003
1004 Double_t dpt = Pt1 - Pt2;
1005 fHistDeltaPtvsJet1Pt->Fill(Pt1, dpt);
1006 fHistDeltaPtvsJet2Pt->Fill(Pt2, dpt);
1007 fHistDeltaPtOverJet1PtvsJet1Pt->Fill(Pt1, dpt/Pt1);
1008 fHistDeltaPtOverJet2PtvsJet2Pt->Fill(Pt2, dpt/Pt2);
1009
1010 fHistDeltaPtvsDistance->Fill(d, dpt);
1011 fHistDeltaPtvsCommonEnergy1->Fill(CE1, dpt);
1012 fHistDeltaPtvsCommonEnergy2->Fill(CE2, dpt);
1013
1014 fHistDeltaPtvsArea1->Fill(A1, dpt);
1015 fHistDeltaPtvsArea2->Fill(A2, dpt);
1016
1017 Double_t darea = A1 - A2;
1018 fHistDeltaPtvsDeltaArea->Fill(darea, dpt);
1019
1020 fHistJet1PtvsJet2Pt->Fill(Pt1, Pt2);
1021
7cd832c7 1022 if (fIsJet1Rho || fIsJet2Rho) {
43032ce2 1023 Double_t dcorrpt = CorrPt1 - CorrPt2;
1024 fHistDeltaCorrPtvsJet1CorrPt->Fill(CorrPt1, dcorrpt);
1025 fHistDeltaCorrPtvsJet2CorrPt->Fill(CorrPt2, dcorrpt);
1026 fHistDeltaCorrPtOverJet1CorrPtvsJet1CorrPt->Fill(CorrPt1, dcorrpt/CorrPt1);
1027 fHistDeltaCorrPtOverJet2CorrPtvsJet2CorrPt->Fill(CorrPt2, dcorrpt/CorrPt2);
1028 fHistDeltaCorrPtvsDistance->Fill(d, dcorrpt);
1029 fHistDeltaCorrPtvsCommonEnergy1->Fill(CE1, dcorrpt);
1030 fHistDeltaCorrPtvsCommonEnergy2->Fill(CE2, dcorrpt);
1031 fHistDeltaCorrPtvsArea1->Fill(A1, dcorrpt);
1032 fHistDeltaCorrPtvsArea2->Fill(A2, dcorrpt);
1033 fHistDeltaCorrPtvsDeltaArea->Fill(darea, dcorrpt);
1034 fHistJet1CorrPtvsJet2CorrPt->Fill(CorrPt1, CorrPt2);
1035 }
1036
1037 if (fIsEmbedded) {
1038 Double_t dmcpt = MCPt1 - Pt2;
1039 fHistDeltaMCPtvsJet1MCPt->Fill(MCPt1, dmcpt);
1040 fHistDeltaMCPtvsJet2Pt->Fill(Pt2, dmcpt);
1041 fHistDeltaMCPtOverJet1MCPtvsJet1MCPt->Fill(MCPt1, dmcpt/MCPt1);
1042 fHistDeltaMCPtOverJet2PtvsJet2Pt->Fill(Pt2, dmcpt/Pt2);
1043 fHistDeltaMCPtvsDistance->Fill(d, dmcpt);
1044 fHistDeltaMCPtvsCommonEnergy1->Fill(CE1, dmcpt);
1045 fHistDeltaMCPtvsCommonEnergy2->Fill(CE2, dmcpt);
1046 fHistDeltaMCPtvsArea1->Fill(A1, dmcpt);
1047 fHistDeltaMCPtvsArea2->Fill(A2, dmcpt);
1048 fHistDeltaMCPtvsDeltaArea->Fill(darea, dmcpt);
1049 fHistJet1MCPtvsJet2Pt->Fill(MCPt1, Pt2);
1050 }
1051 }
1052}
1053
b3ccdfe2 1054//________________________________________________________________________
1055void AliJetResponseMaker::ExecOnce()
1056{
1057 // Execute once.
1058
9239b066 1059 AliAnalysisTaskEmcalJet::ExecOnce();
cfc2ac24 1060
7cd832c7 1061 AliJetContainer *jets1 = static_cast<AliJetContainer*>(fJetCollArray.At(0));
1062 AliJetContainer *jets2 = static_cast<AliJetContainer*>(fJetCollArray.At(1));
cd6431de 1063
7cd832c7 1064 if (!jets1 || !jets1->GetArray() || !jets2 || !jets2->GetArray()) return;
cd6431de 1065
7cd832c7 1066 if (fMatching == kMCLabel && (!jets2->GetIsParticleLevel() || jets1->GetIsParticleLevel())) {
1067 if (jets1->GetIsParticleLevel() == jets2->GetIsParticleLevel()) {
43032ce2 1068 AliWarning("Changing matching type from MC label to same collection...");
1069 fMatching = kSameCollections;
1070 }
1071 else {
1072 AliWarning("Changing matching type from MC label to geometrical...");
1073 fMatching = kGeometrical;
1074 }
1075 }
7cd832c7 1076 else if (fMatching == kSameCollections && jets1->GetIsParticleLevel() != jets2->GetIsParticleLevel()) {
1077 if (jets2->GetIsParticleLevel() && !jets1->GetIsParticleLevel()) {
43032ce2 1078 AliWarning("Changing matching type from same collection to MC label...");
1079 fMatching = kMCLabel;
1080 }
1081 else {
1082 AliWarning("Changing matching type from same collection to geometrical...");
1083 fMatching = kGeometrical;
1084 }
1085 }
b3ccdfe2 1086}
1087
1088//________________________________________________________________________
1089Bool_t AliJetResponseMaker::Run()
1090{
1091 // Find the closest jets
cfc2ac24 1092
f660c2d6 1093 if (fMatching == kNoMatching)
aa4d701c 1094 return kTRUE;
f660c2d6 1095 else
1096 return DoJetMatching();
cfc2ac24 1097}
aa4d701c 1098
cfc2ac24 1099//________________________________________________________________________
f660c2d6 1100Bool_t AliJetResponseMaker::DoJetMatching()
cfc2ac24 1101{
7cd832c7 1102 AliJetContainer *jets1 = static_cast<AliJetContainer*>(fJetCollArray.At(0));
1103 AliJetContainer *jets2 = static_cast<AliJetContainer*>(fJetCollArray.At(1));
cfc2ac24 1104
7cd832c7 1105 if (!jets1 || !jets1->GetArray() || !jets2 || !jets2->GetArray()) return kFALSE;
cfc2ac24 1106
7cd832c7 1107 DoJetLoop();
cfc2ac24 1108
7cd832c7 1109 AliEmcalJet* jet1 = 0;
cfc2ac24 1110
7cd832c7 1111 jets1->ResetCurrentID();
1112 while ((jet1 = jets1->GetNextJet())) {
cfc2ac24 1113
7cd832c7 1114 AliEmcalJet *jet2 = jet1->ClosestJet();
cfc2ac24 1115
7cd832c7 1116 if (!jet2) continue;
1117 if (jet2->ClosestJet() != jet1) continue;
1118 if (jet1->ClosestJetDistance() > fMatchingPar1 || jet2->ClosestJetDistance() > fMatchingPar2) continue;
cfc2ac24 1119
7cd832c7 1120 // Matched jet found
1121 jet1->SetMatchedToClosest(fMatching);
1122 jet2->SetMatchedToClosest(fMatching);
1123 AliDebug(2,Form("Found matching: jet1 pt = %f, eta = %f, phi = %f, jet2 pt = %f, eta = %f, phi = %f",
1124 jet1->Pt(), jet1->Eta(), jet1->Phi(),
1125 jet2->Pt(), jet2->Eta(), jet2->Phi()));
aa4d701c 1126 }
1127
b3ccdfe2 1128 return kTRUE;
1129}
1130
2949a2ec 1131//________________________________________________________________________
7cd832c7 1132void AliJetResponseMaker::DoJetLoop()
2949a2ec 1133{
44476be7 1134 // Do the jet loop.
1ee1b5b8 1135
7cd832c7 1136 AliJetContainer *jets1 = static_cast<AliJetContainer*>(fJetCollArray.At(0));
1137 AliJetContainer *jets2 = static_cast<AliJetContainer*>(fJetCollArray.At(1));
fde82e42 1138
7cd832c7 1139 if (!jets1 || !jets1->GetArray() || !jets2 || !jets2->GetArray()) return;
fde82e42 1140
7cd832c7 1141 AliEmcalJet* jet1 = 0;
1142 AliEmcalJet* jet2 = 0;
fde82e42 1143
7cd832c7 1144 jets2->ResetCurrentID();
1145 while ((jet2 = jets2->GetNextJet())) jet2->ResetMatching();
05077f28 1146
7cd832c7 1147 jets1->ResetCurrentID();
1148 while ((jet1 = jets1->GetNextJet())) {
2103dc6a 1149 jet1->ResetMatching();
7cd832c7 1150
1151 if (jet1->MCPt() < fMinJetMCPt) continue;
44476be7 1152
7cd832c7 1153 jets2->ResetCurrentID();
1154 while ((jet2 = jets2->GetNextJet())) {
7f76e479 1155 SetMatchingLevel(jet1, jet2, fMatching);
2103dc6a 1156 } // jet2 loop
2103dc6a 1157 } // jet1 loop
1ee1b5b8 1158}
1159
cd6431de 1160//________________________________________________________________________
7f76e479 1161void AliJetResponseMaker::GetGeometricalMatchingLevel(AliEmcalJet *jet1, AliEmcalJet *jet2, Double_t &d) const
cd6431de 1162{
7f76e479 1163 Double_t deta = jet2->Eta() - jet1->Eta();
1164 Double_t dphi = jet2->Phi() - jet1->Phi();
1165 d = TMath::Sqrt(deta * deta + dphi * dphi);
1166}
cd6431de 1167
7f76e479 1168//________________________________________________________________________
1169void AliJetResponseMaker::GetMCLabelMatchingLevel(AliEmcalJet *jet1, AliEmcalJet *jet2, Double_t &d1, Double_t &d2) const
1170{
7cd832c7 1171 AliJetContainer *jets1 = static_cast<AliJetContainer*>(fJetCollArray.At(0));
1172 AliJetContainer *jets2 = static_cast<AliJetContainer*>(fJetCollArray.At(1));
1173
1174 if (!jets1 || !jets1->GetArray() || !jets2 || !jets2->GetArray()) return;
1175
1176 AliParticleContainer *tracks1 = jets1->GetParticleContainer();
1177 AliClusterContainer *clusters1 = jets1->GetClusterContainer();
1178 AliParticleContainer *tracks2 = jets2->GetParticleContainer();
1179
7f76e479 1180 // d1 and d2 represent the matching level: 0 = maximum level of matching, 1 = the two jets are completely unrelated
1181 d1 = jet1->Pt();
1182 d2 = jet2->Pt();
1183 Double_t totalPt1 = d1; // the total pt of the reconstructed jet will be cleaned from the background
1184
7cd832c7 1185 // remove completely tracks that are not MC particles (label == 0)
1186 if (tracks1 && tracks1->GetArray()) {
7f76e479 1187 for (Int_t iTrack = 0; iTrack < jet1->GetNumberOfTracks(); iTrack++) {
7cd832c7 1188 AliVParticle *track = jet1->TrackAt(iTrack,tracks1->GetArray());
7f76e479 1189 if (!track) {
1190 AliWarning(Form("Could not find track %d!", iTrack));
1191 continue;
1192 }
7cd832c7 1193
787a3c4f 1194 Int_t MClabel = TMath::Abs(track->GetLabel());
68faaa50 1195 MClabel -= fMCLabelShift;
7cd832c7 1196 if (MClabel != 0) continue;
1197
1198 // this is not a MC particle; remove it completely
1199 AliDebug(3,Form("Track %d (pT = %f) is not a MC particle (MClabel = %d)!",iTrack,track->Pt(),MClabel));
1200 totalPt1 -= track->Pt();
1201 d1 -= track->Pt();
1202 }
1203 }
1204
1205 // remove completely clusters that are not MC particles (label == 0)
1206 if (fUseCellsToMatch && fCaloCells) {
1207 for (Int_t iClus = 0; iClus < jet1->GetNumberOfClusters(); iClus++) {
1208 AliVCluster *clus = jet1->ClusterAt(iClus,clusters1->GetArray());
1209 if (!clus) {
1210 AliWarning(Form("Could not find cluster %d!", iClus));
7f76e479 1211 continue;
1212 }
7cd832c7 1213 TLorentzVector part;
1214 clus->GetMomentum(part, const_cast<Double_t*>(fVertex));
1215
1216 for (Int_t iCell = 0; iCell < clus->GetNCells(); iCell++) {
1217 Int_t cellId = clus->GetCellAbsId(iCell);
1218 Double_t cellFrac = clus->GetCellAmplitudeFraction(iCell);
1219
1220 Int_t MClabel = TMath::Abs(fCaloCells->GetCellMCLabel(cellId));
68faaa50 1221 MClabel -= fMCLabelShift;
7cd832c7 1222 if (MClabel != 0) continue;
1223
1224 // this is not a MC particle; remove it completely
1225 AliDebug(3,Form("Cell %d (frac = %f) is not a MC particle (MClabel = %d)!",iCell,cellFrac,MClabel));
1226 totalPt1 -= part.Pt() * cellFrac;
1227 d1 -= part.Pt() * cellFrac;
1228 }
1229 }
1230 }
1231 else {
1232 for (Int_t iClus = 0; iClus < jet1->GetNumberOfClusters(); iClus++) {
1233 AliVCluster *clus = jet1->ClusterAt(iClus,clusters1->GetArray());
1234 if (!clus) {
1235 AliWarning(Form("Could not find cluster %d!", iClus));
1236 continue;
7f76e479 1237 }
7cd832c7 1238 TLorentzVector part;
1239 clus->GetMomentum(part, const_cast<Double_t*>(fVertex));
7f76e479 1240
7cd832c7 1241 Int_t MClabel = TMath::Abs(clus->GetLabel());
68faaa50 1242 MClabel -= fMCLabelShift;
7cd832c7 1243 if (MClabel != 0) continue;
1244
1245 // this is not a MC particle; remove it completely
1246 AliDebug(3,Form("Cluster %d (pT = %f) is not a MC particle (MClabel = %d)!",iClus,part.Pt(),MClabel));
1247 totalPt1 -= part.Pt();
1248 d1 -= part.Pt();
1249 }
1250 }
1251
1252 for (Int_t iTrack2 = 0; iTrack2 < jet2->GetNumberOfTracks(); iTrack2++) {
1253 Bool_t track2Found = kFALSE;
1254 Int_t index2 = jet2->TrackAt(iTrack2);
1255
1256 // now look for common particles in the track array
1257 for (Int_t iTrack = 0; iTrack < jet1->GetNumberOfTracks(); iTrack++) {
1258 AliVParticle *track = jet1->TrackAt(iTrack,tracks1->GetArray());
1259 if (!track) {
1260 AliWarning(Form("Could not find track %d!", iTrack));
1261 continue;
1262 }
1263 Int_t MClabel = TMath::Abs(track->GetLabel());
68faaa50 1264 MClabel -= fMCLabelShift;
7cd832c7 1265 if (MClabel <= 0) continue;
1266
1267 Int_t index = -1;
1268 index = tracks2->GetIndexFromLabel(MClabel);
7f76e479 1269 if (index < 0) {
1270 AliDebug(2,Form("Track %d (pT = %f) does not have an associated MC particle (MClabel = %d)!",iTrack,track->Pt(),MClabel));
1271 continue;
1272 }
7cd832c7 1273
1274 if (index2 != index) continue;
1275
1276 // found common particle
1277 d1 -= track->Pt();
7f76e479 1278
7cd832c7 1279 if (!track2Found) {
1280 AliVParticle *MCpart = tracks2->GetParticle(index2);
7f76e479 1281 AliDebug(3,Form("Track %d (pT = %f, eta = %f, phi = %f) is associated with the MC particle %d (pT = %f, eta = %f, phi = %f)!",
1282 iTrack,track->Pt(),track->Eta(),track->Phi(),MClabel,MCpart->Pt(),MCpart->Eta(),MCpart->Phi()));
1283 d2 -= MCpart->Pt();
7f76e479 1284 }
7cd832c7 1285
1286 track2Found = kTRUE;
cd6431de 1287 }
7cd832c7 1288
1289 // now look for common particles in the cluster array
1290 if (fUseCellsToMatch && fCaloCells) { // if the cell colection is available, look for cells with a matched MC particle
1291 for (Int_t iClus = 0; iClus < jet1->GetNumberOfClusters(); iClus++) {
1292 AliVCluster *clus = jet1->ClusterAt(iClus,fCaloClusters);
1293 if (!clus) {
1294 AliWarning(Form("Could not find cluster %d!", iClus));
1295 continue;
1296 }
1297 TLorentzVector part;
1298 clus->GetMomentum(part, const_cast<Double_t*>(fVertex));
1299
7f76e479 1300 for (Int_t iCell = 0; iCell < clus->GetNCells(); iCell++) {
1301 Int_t cellId = clus->GetCellAbsId(iCell);
1302 Double_t cellFrac = clus->GetCellAmplitudeFraction(iCell);
7cd832c7 1303
787a3c4f 1304 Int_t MClabel = TMath::Abs(fCaloCells->GetCellMCLabel(cellId));
68faaa50 1305 MClabel -= fMCLabelShift;
7cd832c7 1306 if (MClabel <= 0) continue;
7f76e479 1307
7cd832c7 1308 Int_t index1 = -1;
1309 index1 = tracks2->GetIndexFromLabel(MClabel);
1310 if (index1 < 0) {
7f76e479 1311 AliDebug(3,Form("Cell %d (frac = %f) does not have an associated MC particle (MClabel = %d)!",iCell,cellFrac,MClabel));
1312 continue;
1313 }
7cd832c7 1314
1315 if (index2 != index1) continue;
1316
1317 // found common particle
1318 d1 -= part.Pt() * cellFrac;
1319
1320 if (!track2Found) { // only if it is not already found among charged tracks (charged particles are most likely already found)
1321 AliVParticle *MCpart = tracks2->GetParticle(index2);
1322 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)!",
1323 iCell,iClus,part.Pt(),part.Eta(),part.Phi(),MClabel,MCpart->Pt(),MCpart->Eta(),MCpart->Phi()));
1324 d2 -= MCpart->Pt() * cellFrac;
7f76e479 1325 }
7cd832c7 1326
1327 track2Found = kTRUE;
cfc2ac24 1328 }
1329 }
7cd832c7 1330 }
1331 else { //otherwise look for the first contributor to the cluster, and if matched to a MC label remove it
1332 for (Int_t iClus = 0; iClus < jet1->GetNumberOfClusters(); iClus++) {
1333 AliVCluster *clus = jet1->ClusterAt(iClus,fCaloClusters);
1334 if (!clus) {
1335 AliWarning(Form("Could not find cluster %d!", iClus));
cfc2ac24 1336 continue;
1337 }
7cd832c7 1338 TLorentzVector part;
1339 clus->GetMomentum(part, const_cast<Double_t*>(fVertex));
1340
1341 Int_t MClabel = TMath::Abs(clus->GetLabel());
68faaa50 1342 MClabel -= fMCLabelShift;
7cd832c7 1343 if (MClabel <= 0) continue;
1344
1345 Int_t index = -1;
1346 index = tracks2->GetIndexFromLabel(MClabel);
1347
cfc2ac24 1348 if (index < 0) {
7f76e479 1349 AliDebug(3,Form("Cluster %d (pT = %f) does not have an associated MC particle (MClabel = %d)!",iClus,part.Pt(),MClabel));
cfc2ac24 1350 continue;
1351 }
7cd832c7 1352
1353 if (index2 != index) continue;
1354
1355 // found common particle
1356 d1 -= part.Pt();
1357
1358 if (!track2Found) { // only if it is not already found among charged tracks (charged particles are most likely already found)
1359 AliVParticle *MCpart = tracks2->GetParticle(index2);
1360 AliDebug(3,Form("Cluster %d (pT = %f, eta = %f, phi = %f) is associated with the MC particle %d (pT = %f, eta = %f, phi = %f)!",
1361 iClus,part.Pt(),part.Eta(),part.Phi(),MClabel,MCpart->Pt(),MCpart->Eta(),MCpart->Phi()));
1362
1363 d2 -= MCpart->Pt();
cfc2ac24 1364 }
7cd832c7 1365
1366 track2Found = kTRUE;
cfc2ac24 1367 }
7f76e479 1368 }
1369 }
05077f28 1370
1371 if (d1 < 0)
7f76e479 1372 d1 = 0;
05077f28 1373
1374 if (d2 < 0)
1375 d2 = 0;
1376
1377 if (totalPt1 < 1)
1378 d1 = -1;
7f76e479 1379 else
1380 d1 /= totalPt1;
f660c2d6 1381
05077f28 1382 if (jet2->Pt() < 1)
1383 d2 = -1;
7f76e479 1384 else
05077f28 1385 d2 /= jet2->Pt();
7f76e479 1386}
1387
1388//________________________________________________________________________
1389void AliJetResponseMaker::GetSameCollectionsMatchingLevel(AliEmcalJet *jet1, AliEmcalJet *jet2, Double_t &d1, Double_t &d2) const
1390{
7cd832c7 1391 AliJetContainer *jets1 = static_cast<AliJetContainer*>(fJetCollArray.At(0));
1392 AliJetContainer *jets2 = static_cast<AliJetContainer*>(fJetCollArray.At(1));
1393
1394 if (!jets1 || !jets1->GetArray() || !jets2 || !jets2->GetArray()) return;
1395
1396 AliParticleContainer *tracks1 = jets1->GetParticleContainer();
1397 AliClusterContainer *clusters1 = jets1->GetClusterContainer();
1398 AliParticleContainer *tracks2 = jets2->GetParticleContainer();
1399 AliClusterContainer *clusters2 = jets2->GetClusterContainer();
1400
7f76e479 1401 // d1 and d2 represent the matching level: 0 = maximum level of matching, 1 = the two jets are completely unrelated
1402 d1 = jet1->Pt();
1403 d2 = jet2->Pt();
1404
7cd832c7 1405 if (tracks1 && tracks2) {
7f76e479 1406
1407 for (Int_t iTrack2 = 0; iTrack2 < jet2->GetNumberOfTracks(); iTrack2++) {
1408 Int_t index2 = jet2->TrackAt(iTrack2);
7cd832c7 1409 for (Int_t iTrack1 = 0; iTrack1 < jet1->GetNumberOfTracks(); iTrack1++) {
1410 Int_t index1 = jet1->TrackAt(iTrack1);
1411 if (index2 == index1) { // found common particle
1412 AliVParticle *part1 = tracks1->GetParticle(index1);
1413 if (!part1) {
1414 AliWarning(Form("Could not find track %d!", index1));
f660c2d6 1415 continue;
1416 }
7cd832c7 1417 AliVParticle *part2 = tracks2->GetParticle(index2);
7f76e479 1418 if (!part2) {
1419 AliWarning(Form("Could not find track %d!", index2));
cfc2ac24 1420 continue;
1421 }
7f76e479 1422
7cd832c7 1423 d1 -= part1->Pt();
7f76e479 1424 d2 -= part2->Pt();
1425 break;
cfc2ac24 1426 }
1427 }
cfc2ac24 1428 }
7f76e479 1429
1430 }
1431
7cd832c7 1432 if (clusters1 && clusters2) {
7f76e479 1433
a7477843 1434 if (fUseCellsToMatch) {
1435 const Int_t nClus1 = jet1->GetNumberOfClusters();
1436
1437 Int_t ncells1[nClus1];
1438 UShort_t *cellsId1[nClus1];
1439 Double_t *cellsFrac1[nClus1];
1440 Int_t *sortedIndexes1[nClus1];
1441 Double_t ptClus1[nClus1];
1442 for (Int_t iClus1 = 0; iClus1 < nClus1; iClus1++) {
1443 Int_t index1 = jet1->ClusterAt(iClus1);
7cd832c7 1444 AliVCluster *clus1 = clusters1->GetCluster(index1);
a7477843 1445 if (!clus1) {
1446 AliWarning(Form("Could not find cluster %d!", index1));
9adcb46d 1447 ncells1[iClus1] = 0;
1448 cellsId1[iClus1] = 0;
1449 cellsFrac1[iClus1] = 0;
1450 sortedIndexes1[iClus1] = 0;
1451 ptClus1[iClus1] = 0;
a7477843 1452 continue;
1453 }
9adcb46d 1454 TLorentzVector part1;
1455 clus1->GetMomentum(part1, const_cast<Double_t*>(fVertex));
1456
a7477843 1457 ncells1[iClus1] = clus1->GetNCells();
1458 cellsId1[iClus1] = clus1->GetCellsAbsId();
1459 cellsFrac1[iClus1] = clus1->GetCellsAmplitudeFraction();
1460 sortedIndexes1[iClus1] = new Int_t[ncells1[iClus1]];
a7477843 1461 ptClus1[iClus1] = part1.Pt();
1462
1463 TMath::Sort(ncells1[iClus1], cellsId1[iClus1], sortedIndexes1[iClus1], kFALSE);
1464 }
1465
1466 const Int_t nClus2 = jet2->GetNumberOfClusters();
1467
1468 const Int_t maxNcells2 = 11520;
1469 Int_t sortedIndexes2[maxNcells2];
1470 for (Int_t iClus2 = 0; iClus2 < nClus2; iClus2++) {
1471 Int_t index2 = jet2->ClusterAt(iClus2);
7cd832c7 1472 AliVCluster *clus2 = clusters2->GetCluster(index2);
a7477843 1473 if (!clus2) {
1474 AliWarning(Form("Could not find cluster %d!", index2));
1475 continue;
1476 }
1477 Int_t ncells2 = clus2->GetNCells();
9adcb46d 1478 if (ncells2 >= maxNcells2) {
1479 AliError(Form("Number of cells in the cluster %d >= %d",ncells2,maxNcells2));
a7477843 1480 continue;
1481 }
1482 UShort_t *cellsId2 = clus2->GetCellsAbsId();
1483 Double_t *cellsFrac2 = clus2->GetCellsAmplitudeFraction();
1484
1485 TLorentzVector part2;
1486 clus2->GetMomentum(part2, const_cast<Double_t*>(fVertex));
1487 Double_t ptClus2 = part2.Pt();
1488
1489 TMath::Sort(ncells2, cellsId2, sortedIndexes2, kFALSE);
1490
1491 for (Int_t iClus1 = 0; iClus1 < nClus1; iClus1++) {
9adcb46d 1492 if (sortedIndexes1[iClus1] == 0)
1493 continue;
a7477843 1494 Int_t iCell1 = 0, iCell2 = 0;
a7477843 1495 while (iCell1 < ncells1[iClus1] && iCell2 < ncells2) {
1496 if (cellsId1[iClus1][sortedIndexes1[iClus1][iCell1]] == cellsId2[sortedIndexes2[iCell2]]) { // found a common cell
1497 d1 -= cellsFrac1[iClus1][sortedIndexes1[iClus1][iCell1]] * ptClus1[iClus1];
1498 d2 -= cellsFrac2[sortedIndexes2[iCell2]] * ptClus2;
1499 iCell1++;
1500 iCell2++;
a7477843 1501 }
1502 else if (cellsId1[iClus1][sortedIndexes1[iClus1][iCell1]] > cellsId2[sortedIndexes2[iCell2]]) {
1503 iCell2++;
1504 }
1505 else {
1506 iCell1++;
1507 }
cac0e10b 1508 }
a7477843 1509 }
1510 }
1511 }
1512 else {
1513 for (Int_t iClus2 = 0; iClus2 < jet2->GetNumberOfClusters(); iClus2++) {
1514 Int_t index2 = jet2->ClusterAt(iClus2);
7cd832c7 1515 for (Int_t iClus1 = 0; iClus1 < jet1->GetNumberOfClusters(); iClus1++) {
1516 Int_t index1 = jet1->ClusterAt(iClus1);
1517 if (index2 == index1) { // found common particle
1518 AliVCluster *clus1 = clusters1->GetCluster(index1);
1519 if (!clus1) {
1520 AliWarning(Form("Could not find cluster %d!", index1));
a7477843 1521 continue;
1522 }
7cd832c7 1523 AliVCluster *clus2 = clusters2->GetCluster(index2);
a7477843 1524 if (!clus2) {
1525 AliWarning(Form("Could not find cluster %d!", index2));
1526 continue;
1527 }
7cd832c7 1528 TLorentzVector part1, part2;
1529 clus1->GetMomentum(part1, const_cast<Double_t*>(fVertex));
a7477843 1530 clus2->GetMomentum(part2, const_cast<Double_t*>(fVertex));
1531
7cd832c7 1532 d1 -= part1.Pt();
a7477843 1533 d2 -= part2.Pt();
1534 break;
cac0e10b 1535 }
7f76e479 1536 }
7f76e479 1537 }
1538 }
7f76e479 1539 }
1540
05077f28 1541 if (d1 < 0)
1542 d1 = 0;
1543
1544 if (d2 < 0)
1545 d2 = 0;
1546
1547 if (jet1->Pt() > 0)
7f76e479 1548 d1 /= jet1->Pt();
1549 else
05077f28 1550 d1 = -1;
7f76e479 1551
05077f28 1552 if (jet2->Pt() > 0)
7f76e479 1553 d2 /= jet2->Pt();
1554 else
05077f28 1555 d2 = -1;
7f76e479 1556}
1557
1558//________________________________________________________________________
1559void AliJetResponseMaker::SetMatchingLevel(AliEmcalJet *jet1, AliEmcalJet *jet2, MatchingType matching)
1560{
1561 Double_t d1 = -1;
1562 Double_t d2 = -1;
1563
1564 switch (matching) {
1565 case kGeometrical:
1566 GetGeometricalMatchingLevel(jet1,jet2,d1);
1567 d2 = d1;
1568 break;
1569 case kMCLabel: // jet1 = detector level and jet2 = particle level!
1570 GetMCLabelMatchingLevel(jet1,jet2,d1,d2);
1571 break;
1572 case kSameCollections:
1573 GetSameCollectionsMatchingLevel(jet1,jet2,d1,d2);
cd6431de 1574 break;
1575 default:
1576 ;
1577 }
f660c2d6 1578
05077f28 1579 if (d1 >= 0) {
f660c2d6 1580
1581 if (d1 < jet1->ClosestJetDistance()) {
1582 jet1->SetSecondClosestJet(jet1->ClosestJet(), jet1->ClosestJetDistance());
1583 jet1->SetClosestJet(jet2, d1);
1584 }
1585 else if (d1 < jet1->SecondClosestJetDistance()) {
1586 jet1->SetSecondClosestJet(jet2, d1);
1587 }
1588 }
1589
05077f28 1590 if (d2 >= 0) {
f660c2d6 1591
1592 if (d2 < jet2->ClosestJetDistance()) {
1593 jet2->SetSecondClosestJet(jet2->ClosestJet(), jet2->ClosestJetDistance());
1594 jet2->SetClosestJet(jet1, d2);
1595 }
1596 else if (d2 < jet2->SecondClosestJetDistance()) {
1597 jet2->SetSecondClosestJet(jet1, d2);
1598 }
1599 }
cd6431de 1600}
1601
1ee1b5b8 1602//________________________________________________________________________
1603Bool_t AliJetResponseMaker::FillHistograms()
1604{
1605 // Fill histograms.
1606
7cd832c7 1607 AliJetContainer *jets1 = static_cast<AliJetContainer*>(fJetCollArray.At(0));
1608 AliJetContainer *jets2 = static_cast<AliJetContainer*>(fJetCollArray.At(1));
2949a2ec 1609
7cd832c7 1610 if (!jets1 || !jets1->GetArray() || !jets2 || !jets2->GetArray()) return kFALSE;
ca5c29fa 1611
7cd832c7 1612 AliEmcalJet* jet1 = 0;
1613 AliEmcalJet* jet2 = 0;
1614
1615 jets2->ResetCurrentID();
1616 while ((jet2 = jets2->GetNextJet())) {
d3c0e12a 1617
7cd832c7 1618 AliDebug(2,Form("Processing jet (2) %d", jets2->GetCurrentID()));
ca5c29fa 1619
d3c0e12a 1620 Double_t ptLeading2 = jets2->GetLeadingHadronPt(jet2);
1621 Double_t corrpt2 = jet2->Pt() - jets2->GetRhoVal() * jet2->Area();
7cd832c7 1622
7cd832c7 1623 if (jets2->AcceptJet(jet2))
d3c0e12a 1624 FillJetHisto(jet2->Phi(), jet2->Eta(), jet2->Pt(), jet2->Area(), jet2->NEF(), ptLeading2,
1625 corrpt2, jet2->MCPt(), 2);
7cd832c7 1626
1627 jet1 = jet2->MatchedJet();
1628
1629 if (!jet1) continue;
1630 if (!jets1->AcceptJet(jet1)) continue;
1631 if (jet1->MCPt() < fMinJetMCPt) continue;
1632
1633 Double_t d=-1, ce1=-1, ce2=-1;
1634 if (jet2->GetMatchingType() == kGeometrical) {
1635 if (jets2->GetIsParticleLevel() && !jets1->GetIsParticleLevel()) // the other way around is not supported
1636 GetMCLabelMatchingLevel(jet1, jet2, ce1, ce2);
1637 else if (jets1->GetIsParticleLevel() == jets2->GetIsParticleLevel())
1638 GetSameCollectionsMatchingLevel(jet1, jet2, ce1, ce2);
43032ce2 1639
7cd832c7 1640 d = jet2->ClosestJetDistance();
43032ce2 1641 }
7cd832c7 1642 else if (jet2->GetMatchingType() == kMCLabel || jet2->GetMatchingType() == kSameCollections) {
1643 GetGeometricalMatchingLevel(jet1, jet2, d);
43032ce2 1644
7cd832c7 1645 ce1 = jet1->ClosestJetDistance();
1646 ce2 = jet2->ClosestJetDistance();
2949a2ec 1647 }
1b76c28f 1648
d3c0e12a 1649 Double_t ptLeading1 = jets1->GetLeadingHadronPt(jet1);
7cd832c7 1650 Double_t corrpt1 = jet1->Pt() - jets1->GetRhoVal() * jet1->Area();
1b76c28f 1651
7cd832c7 1652 FillMatchingHistos(jet1->Pt(), jet2->Pt(), jet1->Eta(), jet2->Eta(), jet1->Phi(), jet2->Phi(),
1653 jet1->Area(), jet2->Area(), d, ce1, ce2, corrpt1, corrpt2, jet1->MCPt(),
d3c0e12a 1654 jet1->NEF(), jet2->NEF(), ptLeading1, ptLeading2);
7cd832c7 1655 }
2949a2ec 1656
7cd832c7 1657 jets1->ResetCurrentID();
1658 while ((jet1 = jets1->GetNextAcceptJet())) {
1659 if (jet1->MCPt() < fMinJetMCPt) continue;
1660 AliDebug(2,Form("Processing jet (1) %d", jets1->GetCurrentID()));
ca5c29fa 1661
d3c0e12a 1662 Double_t ptLeading1 = jets1->GetLeadingHadronPt(jet1);
1663 Double_t corrpt1 = jet1->Pt() - jets1->GetRhoVal() * jet1->Area();
1664
1665 FillJetHisto(jet1->Phi(), jet1->Eta(), jet1->Pt(), jet1->Area(), jet1->NEF(), ptLeading1,
1666 corrpt1, jet1->MCPt(), 1);
1b76c28f 1667 }
6fd5039f 1668 return kTRUE;
1b76c28f 1669}