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