]> git.uio.no Git - u/mrichter/AliRoot.git/blame - PWGJE/EMCALJetTasks/AliAnalysisTaskDeltaPt.cxx
git://git.uio.no / u / mrichter / AliRoot.git / blame
? search:
summary | shortlog | log | commit | commitdiff | tree
blob | blame (incremental) | history | HEAD
RelVal: small error message fix
[u/mrichter/AliRoot.git] / PWGJE / EMCALJetTasks / AliAnalysisTaskDeltaPt.cxx
CommitLineData
a487deae 1// $Id$
2//
3// Jet deltaPt task.
4//
5// Author: S.Aiola
6
a487deae 7#include <TClonesArray.h>
8#include <TH1F.h>
9#include <TH2F.h>
ca5c29fa 10#include <TH3F.h>
a487deae 11#include <TList.h>
12#include <TLorentzVector.h>
13#include <TRandom3.h>
a487deae 14
a487deae 15#include "AliVCluster.h"
16#include "AliVParticle.h"
17#include "AliVTrack.h"
18#include "AliEmcalJet.h"
a487deae 19#include "AliRhoParameter.h"
20#include "AliLog.h"
6421eeb0 21#include "AliJetContainer.h"
22#include "AliParticleContainer.h"
23#include "AliClusterContainer.h"
a487deae 24
25#include "AliAnalysisTaskDeltaPt.h"
26
27ClassImp(AliAnalysisTaskDeltaPt)
28
29//________________________________________________________________________
30AliAnalysisTaskDeltaPt::AliAnalysisTaskDeltaPt() :
9239b066 31 AliAnalysisTaskEmcalJet("AliAnalysisTaskDeltaPt", kTRUE),
f660c2d6 32 fMCJetPtThreshold(1),
a487deae 33 fMinRC2LJ(-1),
a487deae 34 fRCperEvent(-1),
6421eeb0 35 fConeRadius(0.2),
36 fConeMinEta(-0.9),
37 fConeMaxEta(0.9),
38 fConeMinPhi(0),
39 fConeMaxPhi(TMath::Pi()*2),
40 fJetsCont(0),
41 fTracksCont(0),
42 fCaloClustersCont(0),
43 fEmbJetsCont(0),
44 fEmbTracksCont(0),
45 fEmbCaloClustersCont(0),
46 fRandTracksCont(0),
47 fRandCaloClustersCont(0),
56f370b0 48 fHistRCPhiEta(0),
49 fHistRCPt(0),
50 fHistRCPtExLJ(0),
51 fHistRCPtExPartialLJ(0),
52 fHistRCPtRand(0),
53 fHistRhoVSRCPt(0),
faf21244 54 fHistDeltaPtRCvsEP(0),
56f370b0 55 fHistDeltaPtRCExLJ(0),
56 fHistDeltaPtRCExPartialLJ(0),
57 fHistDeltaPtRCRand(0),
56f370b0 58 fHistEmbJetsPtArea(0),
59 fHistEmbJetsCorrPtArea(0),
60 fHistEmbPartPtvsJetPt(0),
61 fHistEmbPartPtvsJetCorrPt(0),
62 fHistJetPtvsJetCorrPt(0),
63 fHistDistLeadPart2JetAxis(0),
64 fHistEmbBkgArea(0),
65 fHistRhoVSEmbBkg(0),
66 fHistDeltaPtEmbArea(0),
faf21244 67 fHistDeltaPtEmbvsEP(0),
a487deae 68 fHistRCPtExLJVSDPhiLJ(0),
56f370b0 69 fHistRCPtExPartialLJVSDPhiLJ(0),
6a20534a 70 fHistEmbJetsPhiEta(0),
71 fHistLeadPartPhiEta(0)
a487deae 72{
73 // Default constructor.
74
ed062c95 75 fHistRCPt = 0;
76 fHistRCPtExLJ = 0;
77 fHistRCPtExPartialLJ = 0;
78 fHistRCPtRand = 0;
79 fHistRhoVSRCPt = 0;
80 fHistDeltaPtRCvsEP = 0;
81 fHistDeltaPtRCExLJ = 0;
82 fHistDeltaPtRCExPartialLJ = 0;
83 fHistDeltaPtRCRand = 0;
84 fHistEmbJetsPtArea = 0;
85 fHistEmbJetsCorrPtArea = 0;
86 fHistEmbPartPtvsJetPt = 0;
87 fHistEmbPartPtvsJetCorrPt = 0;
88 fHistJetPtvsJetCorrPt = 0;
89 fHistDistLeadPart2JetAxis = 0;
90 fHistEmbBkgArea = 0;
91 fHistRhoVSEmbBkg = 0;
92 fHistDeltaPtEmbArea = 0;
93 fHistDeltaPtEmbvsEP = 0;
a487deae 94
95 SetMakeGeneralHistograms(kTRUE);
96}
97
98//________________________________________________________________________
99AliAnalysisTaskDeltaPt::AliAnalysisTaskDeltaPt(const char *name) :
9239b066 100 AliAnalysisTaskEmcalJet(name, kTRUE),
f660c2d6 101 fMCJetPtThreshold(1),
a487deae 102 fMinRC2LJ(-1),
a487deae 103 fRCperEvent(-1),
6421eeb0 104 fConeRadius(0.2),
105 fConeMinEta(-0.9),
106 fConeMaxEta(0.9),
107 fConeMinPhi(0),
108 fConeMaxPhi(TMath::Pi()*2),
109 fJetsCont(0),
110 fTracksCont(0),
111 fCaloClustersCont(0),
112 fEmbJetsCont(0),
113 fEmbTracksCont(0),
114 fEmbCaloClustersCont(0),
115 fRandTracksCont(0),
116 fRandCaloClustersCont(0),
56f370b0 117 fHistRCPhiEta(0),
118 fHistRCPt(0),
119 fHistRCPtExLJ(0),
120 fHistRCPtExPartialLJ(0),
121 fHistRCPtRand(0),
122 fHistRhoVSRCPt(0),
faf21244 123 fHistDeltaPtRCvsEP(0),
56f370b0 124 fHistDeltaPtRCExLJ(0),
125 fHistDeltaPtRCExPartialLJ(0),
126 fHistDeltaPtRCRand(0),
56f370b0 127 fHistEmbJetsPtArea(0),
128 fHistEmbJetsCorrPtArea(0),
129 fHistEmbPartPtvsJetPt(0),
130 fHistEmbPartPtvsJetCorrPt(0),
131 fHistJetPtvsJetCorrPt(0),
132 fHistDistLeadPart2JetAxis(0),
133 fHistEmbBkgArea(0),
134 fHistRhoVSEmbBkg(0),
135 fHistDeltaPtEmbArea(0),
faf21244 136 fHistDeltaPtEmbvsEP(0),
a487deae 137 fHistRCPtExLJVSDPhiLJ(0),
56f370b0 138 fHistRCPtExPartialLJVSDPhiLJ(0),
6a20534a 139 fHistEmbJetsPhiEta(0),
140 fHistLeadPartPhiEta(0)
a487deae 141{
142 // Standard constructor.
143
ed062c95 144 fHistRCPt = 0;
145 fHistRCPtExLJ = 0;
146 fHistRCPtExPartialLJ = 0;
147 fHistRCPtRand = 0;
148 fHistRhoVSRCPt = 0;
149 fHistDeltaPtRCvsEP = 0;
150 fHistDeltaPtRCExLJ = 0;
151 fHistDeltaPtRCExPartialLJ = 0;
152 fHistDeltaPtRCRand = 0;
153 fHistEmbJetsPtArea = 0;
154 fHistEmbJetsCorrPtArea = 0;
155 fHistEmbPartPtvsJetPt = 0;
156 fHistEmbPartPtvsJetCorrPt = 0;
157 fHistJetPtvsJetCorrPt = 0;
158 fHistDistLeadPart2JetAxis = 0;
159 fHistEmbBkgArea = 0;
160 fHistRhoVSEmbBkg = 0;
161 fHistDeltaPtEmbArea = 0;
162 fHistDeltaPtEmbvsEP = 0;
163
164 SetMakeGeneralHistograms(kTRUE);
165}
166
167//________________________________________________________________________
168void AliAnalysisTaskDeltaPt::AllocateHistogramArrays()
169{
56f370b0 170 fHistRCPt = new TH1*[fNcentBins];
171 fHistRCPtExLJ = new TH1*[fNcentBins];
172 fHistRCPtExPartialLJ = new TH1*[fNcentBins];
173 fHistRCPtRand = new TH1*[fNcentBins];
174 fHistRhoVSRCPt = new TH2*[fNcentBins];
faf21244 175 fHistDeltaPtRCvsEP = new TH2*[fNcentBins];
56f370b0 176 fHistDeltaPtRCExLJ = new TH1*[fNcentBins];
177 fHistDeltaPtRCExPartialLJ = new TH1*[fNcentBins];
178 fHistDeltaPtRCRand = new TH1*[fNcentBins];
56f370b0 179 fHistEmbJetsPtArea = new TH3*[fNcentBins];
180 fHistEmbJetsCorrPtArea = new TH3*[fNcentBins];
181 fHistEmbPartPtvsJetPt = new TH2*[fNcentBins];
182 fHistEmbPartPtvsJetCorrPt = new TH2*[fNcentBins];
183 fHistJetPtvsJetCorrPt = new TH2*[fNcentBins];
184 fHistDistLeadPart2JetAxis = new TH1*[fNcentBins];
185 fHistEmbBkgArea = new TH2*[fNcentBins];
186 fHistRhoVSEmbBkg = new TH2*[fNcentBins];
187 fHistDeltaPtEmbArea = new TH2*[fNcentBins];
faf21244 188 fHistDeltaPtEmbvsEP = new TH2*[fNcentBins];
56f370b0 189
190 for (Int_t i = 0; i < fNcentBins; i++) {
a487deae 191 fHistRCPt[i] = 0;
192 fHistRCPtExLJ[i] = 0;
56f370b0 193 fHistRCPtExPartialLJ[i] = 0;
194 fHistRCPtRand[i] = 0;
59f16b27 195 fHistRhoVSRCPt[i] = 0;
faf21244 196 fHistDeltaPtRCvsEP[i] = 0;
a487deae 197 fHistDeltaPtRCExLJ[i] = 0;
56f370b0 198 fHistDeltaPtRCExPartialLJ[i] = 0;
a487deae 199 fHistDeltaPtRCRand[i] = 0;
1f9c287f 200 fHistEmbJetsPtArea[i] = 0;
201 fHistEmbJetsCorrPtArea[i] = 0;
6a20534a 202 fHistEmbPartPtvsJetPt[i] = 0;
ca5c29fa 203 fHistEmbPartPtvsJetCorrPt[i] = 0;
2103dc6a 204 fHistJetPtvsJetCorrPt[i] = 0;
6a20534a 205 fHistDistLeadPart2JetAxis[i] = 0;
1f9c287f 206 fHistEmbBkgArea[i] = 0;
59f16b27 207 fHistRhoVSEmbBkg[i] = 0;
1f9c287f 208 fHistDeltaPtEmbArea[i] = 0;
faf21244 209 fHistDeltaPtEmbvsEP[i] = 0;
a487deae 210 }
a487deae 211}
212
a487deae 213//________________________________________________________________________
214void AliAnalysisTaskDeltaPt::UserCreateOutputObjects()
215{
216 // Create user output.
217
9239b066 218 AliAnalysisTaskEmcalJet::UserCreateOutputObjects();
a487deae 219
ed062c95 220 AllocateHistogramArrays();
221
6421eeb0 222 fJetsCont = GetJetContainer("Jets");
223 fTracksCont = GetParticleContainer("Tracks");
224 fCaloClustersCont = GetClusterContainer("CaloClusters");
225 fEmbJetsCont = GetJetContainer("EmbJets");
226 fEmbTracksCont = GetParticleContainer("EmbTracks");
227 fEmbCaloClustersCont = GetClusterContainer("EmbCaloClusters");
228 fRandTracksCont = GetParticleContainer("RandTracks");
229 fRandCaloClustersCont = GetClusterContainer("RandCaloClusters");
230
231 if (fTracksCont || fCaloClustersCont) {
5ce8ae64 232 fHistRCPhiEta = new TH2F("fHistRCPhiEta","fHistRCPhiEta", 100, -1, 1, 201, 0, TMath::Pi() * 2.01);
6a20534a 233 fHistRCPhiEta->GetXaxis()->SetTitle("#eta");
234 fHistRCPhiEta->GetYaxis()->SetTitle("#phi");
235 fOutput->Add(fHistRCPhiEta);
236
6421eeb0 237 if (fJetsCont) {
6a20534a 238 fHistRCPtExLJVSDPhiLJ = new TH2F("fHistRCPtExLJVSDPhiLJ","fHistRCPtExLJVSDPhiLJ", fNbins, fMinBinPt, fMaxBinPt, 128, -1.6, 4.8);
239 fHistRCPtExLJVSDPhiLJ->GetXaxis()->SetTitle("#it{p}_{T} (GeV/#it{c})");
240 fHistRCPtExLJVSDPhiLJ->GetYaxis()->SetTitle("#Delta#phi");
241 fOutput->Add(fHistRCPtExLJVSDPhiLJ);
56f370b0 242
243 fHistRCPtExPartialLJVSDPhiLJ = new TH2F("fHistRCPtExPartialLJVSDPhiLJ","fHistRCPtExPartialLJVSDPhiLJ", fNbins, fMinBinPt, fMaxBinPt, 128, -1.6, 4.8);
244 fHistRCPtExPartialLJVSDPhiLJ->GetXaxis()->SetTitle("#it{p}_{T} (GeV/#it{c})");
245 fHistRCPtExPartialLJVSDPhiLJ->GetYaxis()->SetTitle("#Delta#phi");
246 fOutput->Add(fHistRCPtExPartialLJVSDPhiLJ);
6a20534a 247 }
248 }
249
6421eeb0 250 if (fEmbJetsCont) {
5ce8ae64 251 fHistEmbJetsPhiEta = new TH2F("fHistEmbJetsPhiEta","fHistEmbJetsPhiEta", 100, -1, 1, 201, 0, TMath::Pi() * 2.01);
6a20534a 252 fHistEmbJetsPhiEta->GetXaxis()->SetTitle("#eta");
253 fHistEmbJetsPhiEta->GetYaxis()->SetTitle("#phi");
254 fOutput->Add(fHistEmbJetsPhiEta);
a487deae 255
5ce8ae64 256 fHistLeadPartPhiEta = new TH2F("fHistLeadPartPhiEta","fHistLeadPartPhiEta", 100, -1, 1, 201, 0, TMath::Pi() * 2.01);
6a20534a 257 fHistLeadPartPhiEta->GetXaxis()->SetTitle("#eta");
258 fHistLeadPartPhiEta->GetYaxis()->SetTitle("#phi");
259 fOutput->Add(fHistLeadPartPhiEta);
a487deae 260 }
261
262 TString histname;
263
ca5c29fa 264 const Int_t nbinsZ = 12;
265 Float_t binsZ[nbinsZ+1] = {0,1,2,3,4,5,6,7,8,9,10,20,1000};
266
267 Float_t *binsPt = GenerateFixedBinArray(fNbins, fMinBinPt, fMaxBinPt);
268 Float_t *binsCorrPt = GenerateFixedBinArray(fNbins*2, -fMaxBinPt, fMaxBinPt);
9adcb46d 269 Float_t *binsArea = GenerateFixedBinArray(50, 0, 2);
ca5c29fa 270
6a20534a 271 for (Int_t i = 0; i < fNcentBins; i++) {
6421eeb0 272 if (fTracksCont || fCaloClustersCont) {
6a20534a 273 histname = "fHistRCPt_";
274 histname += i;
275 fHistRCPt[i] = new TH1F(histname.Data(), histname.Data(), fNbins, fMinBinPt, fMaxBinPt * 2);
276 fHistRCPt[i]->GetXaxis()->SetTitle("#it{p}_{T} (GeV/#it{c})");
277 fHistRCPt[i]->GetYaxis()->SetTitle("counts");
278 fOutput->Add(fHistRCPt[i]);
279
280 histname = "fHistRhoVSRCPt_";
281 histname += i;
282 fHistRhoVSRCPt[i] = new TH2F(histname.Data(), histname.Data(), fNbins, fMinBinPt, fMaxBinPt, fNbins, fMinBinPt, fMaxBinPt);
283 fHistRhoVSRCPt[i]->GetXaxis()->SetTitle("A#rho (GeV/#it{c})");
284 fHistRhoVSRCPt[i]->GetYaxis()->SetTitle("#it{p}_{T} (GeV/#it{c})");
285 fOutput->Add(fHistRhoVSRCPt[i]);
286
faf21244 287 histname = "fHistDeltaPtRCvsEP_";
6a20534a 288 histname += i;
9f6dc7fd 289 fHistDeltaPtRCvsEP[i] = new TH2F(histname.Data(), histname.Data(), 101, 0, TMath::Pi()*1.01, fNbins * 2, -fMaxBinPt, fMaxBinPt);
290 fHistDeltaPtRCvsEP[i]->GetXaxis()->SetTitle("#phi_{RC} - #psi_{RP}");
faf21244 291 fHistDeltaPtRCvsEP[i]->GetYaxis()->SetTitle("#delta#it{p}_{T}^{RC} (GeV/#it{c})");
292 fHistDeltaPtRCvsEP[i]->GetZaxis()->SetTitle("counts");
293 fOutput->Add(fHistDeltaPtRCvsEP[i]);
6a20534a 294
6421eeb0 295 if (fJetsCont) {
6a20534a 296 histname = "fHistRCPtExLJ_";
297 histname += i;
298 fHistRCPtExLJ[i] = new TH1F(histname.Data(), histname.Data(), fNbins, fMinBinPt, fMaxBinPt * 2);
299 fHistRCPtExLJ[i]->GetXaxis()->SetTitle("#it{p}_{T}^{RC} (GeV/#it{c})");
300 fHistRCPtExLJ[i]->GetYaxis()->SetTitle("counts");
301 fOutput->Add(fHistRCPtExLJ[i]);
302
303 histname = "fHistDeltaPtRCExLJ_";
304 histname += i;
305 fHistDeltaPtRCExLJ[i] = new TH1F(histname.Data(), histname.Data(), fNbins * 2, -fMaxBinPt, fMaxBinPt);
306 fHistDeltaPtRCExLJ[i]->GetXaxis()->SetTitle("#delta#it{p}_{T}^{RC} (GeV/#it{c})");
307 fHistDeltaPtRCExLJ[i]->GetYaxis()->SetTitle("counts");
308 fOutput->Add(fHistDeltaPtRCExLJ[i]);
56f370b0 309
310 histname = "fHistRCPtExPartialLJ_";
311 histname += i;
312 fHistRCPtExPartialLJ[i] = new TH1F(histname.Data(), histname.Data(), fNbins, fMinBinPt, fMaxBinPt * 2);
313 fHistRCPtExPartialLJ[i]->GetXaxis()->SetTitle("#it{p}_{T}^{RC} (GeV/#it{c})");
314 fHistRCPtExPartialLJ[i]->GetYaxis()->SetTitle("counts");
315 fOutput->Add(fHistRCPtExPartialLJ[i]);
316
317 histname = "fHistDeltaPtRCExPartialLJ_";
318 histname += i;
319 fHistDeltaPtRCExPartialLJ[i] = new TH1F(histname.Data(), histname.Data(), fNbins * 2, -fMaxBinPt, fMaxBinPt);
320 fHistDeltaPtRCExPartialLJ[i]->GetXaxis()->SetTitle("#delta#it{p}_{T}^{RC} (GeV/#it{c})");
321 fHistDeltaPtRCExPartialLJ[i]->GetYaxis()->SetTitle("counts");
322 fOutput->Add(fHistDeltaPtRCExPartialLJ[i]);
6a20534a 323 }
324 }
325
6421eeb0 326 if (fRandTracksCont || fRandCaloClustersCont) {
6a20534a 327 histname = "fHistRCPtRand_";
328 histname += i;
329 fHistRCPtRand[i] = new TH1F(histname.Data(), histname.Data(), fNbins, fMinBinPt, fMaxBinPt * 2);
330 fHistRCPtRand[i]->GetXaxis()->SetTitle("#it{p}_{T}^{RC} (GeV/#it{c})");
331 fHistRCPtRand[i]->GetYaxis()->SetTitle("counts");
332 fOutput->Add(fHistRCPtRand[i]);
333
334 histname = "fHistDeltaPtRCRand_";
335 histname += i;
336 fHistDeltaPtRCRand[i] = new TH1F(histname.Data(), histname.Data(), fNbins * 2, -fMaxBinPt, fMaxBinPt);
337 fHistDeltaPtRCRand[i]->GetXaxis()->SetTitle("#delta#it{p}_{T}^{RC} (GeV/#it{c})");
338 fHistDeltaPtRCRand[i]->GetYaxis()->SetTitle("counts");
339 fOutput->Add(fHistDeltaPtRCRand[i]);
340 }
341
6421eeb0 342 if (fEmbJetsCont) {
1f9c287f 343 histname = "fHistEmbJetsPtArea_";
a487deae 344 histname += i;
9adcb46d 345 fHistEmbJetsPtArea[i] = new TH3F(histname.Data(), histname.Data(), 50, binsArea, fNbins, binsPt, nbinsZ, binsZ);
1f9c287f 346 fHistEmbJetsPtArea[i]->GetXaxis()->SetTitle("area");
6a20534a 347 fHistEmbJetsPtArea[i]->GetYaxis()->SetTitle("#it{p}_{T,jet}^{emb,raw} (GeV/#it{c})");
1f9c287f 348 fOutput->Add(fHistEmbJetsPtArea[i]);
a487deae 349
1f9c287f 350 histname = "fHistEmbJetsCorrPtArea_";
a487deae 351 histname += i;
9adcb46d 352 fHistEmbJetsCorrPtArea[i] = new TH3F(histname.Data(), histname.Data(), 50, binsArea, fNbins * 2, binsCorrPt, nbinsZ, binsZ);
1f9c287f 353 fHistEmbJetsCorrPtArea[i]->GetXaxis()->SetTitle("area");
6a20534a 354 fHistEmbJetsCorrPtArea[i]->GetYaxis()->SetTitle("#it{p}_{T,jet}^{emb,corr} (GeV/#it{c})");
1f9c287f 355 fOutput->Add(fHistEmbJetsCorrPtArea[i]);
a487deae 356
6a20534a 357 histname = "fHistEmbPartPtvsJetPt_";
358 histname += i;
359 fHistEmbPartPtvsJetPt[i] = new TH2F(histname.Data(), histname.Data(), fNbins, fMinBinPt, fMaxBinPt, fNbins, fMinBinPt, fMaxBinPt);
360 fHistEmbPartPtvsJetPt[i]->GetXaxis()->SetTitle("#sum#it{p}_{T,const}^{emb} (GeV/#it{c})");
361 fHistEmbPartPtvsJetPt[i]->GetYaxis()->SetTitle("#it{p}_{T,jet}^{emb} (GeV/#it{c})");
362 fHistEmbPartPtvsJetPt[i]->GetZaxis()->SetTitle("counts");
363 fOutput->Add(fHistEmbPartPtvsJetPt[i]);
364
ca5c29fa 365 histname = "fHistEmbPartPtvsJetCorrPt_";
a487deae 366 histname += i;
9adcb46d 367 fHistEmbPartPtvsJetCorrPt[i] = new TH2F(histname.Data(), histname.Data(),
368 fNbins, fMinBinPt, fMaxBinPt, fNbins*2, -fMaxBinPt, fMaxBinPt);
ca5c29fa 369 fHistEmbPartPtvsJetCorrPt[i]->GetXaxis()->SetTitle("#sum#it{p}_{T,const}^{emb} (GeV/#it{c})");
370 fHistEmbPartPtvsJetCorrPt[i]->GetYaxis()->SetTitle("#it{p}_{T,jet}^{emb} - A#rho (GeV/#it{c})");
371 fHistEmbPartPtvsJetCorrPt[i]->GetZaxis()->SetTitle("counts");
372 fOutput->Add(fHistEmbPartPtvsJetCorrPt[i]);
a487deae 373
2103dc6a 374 histname = "fHistJetPtvsJetCorrPt_";
375 histname += i;
9adcb46d 376 fHistJetPtvsJetCorrPt[i] = new TH2F(histname.Data(), histname.Data(),
377 fNbins, fMinBinPt, fMaxBinPt, fNbins*2, -fMaxBinPt, fMaxBinPt);
2103dc6a 378 fHistJetPtvsJetCorrPt[i]->GetXaxis()->SetTitle("#it{p}_{T,jet}^{emb} (GeV/#it{c})");
379 fHistJetPtvsJetCorrPt[i]->GetYaxis()->SetTitle("#it{p}_{T,jet}^{emb} - A#rho (GeV/#it{c})");
380 fHistJetPtvsJetCorrPt[i]->GetZaxis()->SetTitle("counts");
381 fOutput->Add(fHistJetPtvsJetCorrPt[i]);
382
6a20534a 383 histname = "fHistDistLeadPart2JetAxis_";
a487deae 384 histname += i;
6a20534a 385 fHistDistLeadPart2JetAxis[i] = new TH1F(histname.Data(), histname.Data(), 50, 0, 0.5);
386 fHistDistLeadPart2JetAxis[i]->GetXaxis()->SetTitle("distance");
387 fHistDistLeadPart2JetAxis[i]->GetYaxis()->SetTitle("counts");
388 fOutput->Add(fHistDistLeadPart2JetAxis[i]);
a487deae 389
1f9c287f 390 histname = "fHistEmbBkgArea_";
391 histname += i;
9adcb46d 392 fHistEmbBkgArea[i] = new TH2F(histname.Data(), histname.Data(), 50, 0, 2, fNbins, fMinBinPt, fMaxBinPt);
1f9c287f 393 fHistEmbBkgArea[i]->GetXaxis()->SetTitle("area");
6a20534a 394 fHistEmbBkgArea[i]->GetYaxis()->SetTitle("#it{p}_{T,jet}^{emb} - #sum#it{p}_{T,const}^{emb} (GeV/#it{c})");
1f9c287f 395 fOutput->Add(fHistEmbBkgArea[i]);
396
59f16b27 397 histname = "fHistRhoVSEmbBkg_";
398 histname += i;
1f9c287f 399 fHistRhoVSEmbBkg[i] = new TH2F(histname.Data(), histname.Data(), fNbins, fMinBinPt, fMaxBinPt, fNbins, fMinBinPt, fMaxBinPt);
400 fHistRhoVSEmbBkg[i]->GetXaxis()->SetTitle("A#rho (GeV/#it{c})");
6a20534a 401 fHistRhoVSEmbBkg[i]->GetYaxis()->SetTitle("#it{p}_{T,jet}^{emb} - #sum#it{p}_{T,const}^{emb} (GeV/#it{c})");
59f16b27 402 fOutput->Add(fHistRhoVSEmbBkg[i]);
a487deae 403
1f9c287f 404 histname = "fHistDeltaPtEmbArea_";
a487deae 405 histname += i;
9adcb46d 406 fHistDeltaPtEmbArea[i] = new TH2F(histname.Data(), histname.Data(),
faf21244 407 50, 0, 2, fNbins * 2, -fMaxBinPt, fMaxBinPt);
1f9c287f 408 fHistDeltaPtEmbArea[i]->GetXaxis()->SetTitle("area");
409 fHistDeltaPtEmbArea[i]->GetYaxis()->SetTitle("#delta#it{p}_{T}^{emb} (GeV/#it{c})");
faf21244 410 fHistDeltaPtEmbArea[i]->GetZaxis()->SetTitle("counts");
1f9c287f 411 fOutput->Add(fHistDeltaPtEmbArea[i]);
faf21244 412
413 histname = "fHistDeltaPtEmbvsEP_";
414 histname += i;
7090b9be 415 fHistDeltaPtEmbvsEP[i] = new TH2F(histname.Data(), histname.Data(), 101, 0, TMath::Pi()*1.01, fNbins * 2, -fMaxBinPt, fMaxBinPt);
faf21244 416 fHistDeltaPtEmbvsEP[i]->GetXaxis()->SetTitle("#phi_{jet} - #Psi_{EP}");
417 fHistDeltaPtEmbvsEP[i]->GetYaxis()->SetTitle("#delta#it{p}_{T}^{emb} (GeV/#it{c})");
418 fHistDeltaPtEmbvsEP[i]->GetZaxis()->SetTitle("counts");
419 fOutput->Add(fHistDeltaPtEmbvsEP[i]);
a487deae 420 }
421 }
422
cac0e10b 423 delete[] binsPt;
424 delete[] binsCorrPt;
425 delete[] binsArea;
426
a487deae 427 PostData(1, fOutput); // Post data for ALL output slots >0 here, to get at least an empty histogram
428}
429
430//________________________________________________________________________
431Bool_t AliAnalysisTaskDeltaPt::FillHistograms()
432{
433 // Fill histograms.
434
a487deae 435 // ************
436 // Random cones
437 // _________________________________
438
6421eeb0 439 const Float_t rcArea = fConeRadius * fConeRadius * TMath::Pi();
6a20534a 440 Float_t RCpt = 0;
441 Float_t RCeta = 0;
442 Float_t RCphi = 0;
a487deae 443
6421eeb0 444 if (fTracksCont || fCaloClustersCont) {
6a20534a 445
446 for (Int_t i = 0; i < fRCperEvent; i++) {
447 // Simple random cones
448 RCpt = 0;
449 RCeta = 0;
450 RCphi = 0;
6421eeb0 451 GetRandomCone(RCpt, RCeta, RCphi, fTracksCont, fCaloClustersCont, 0);
6a20534a 452 if (RCpt > 0) {
453 fHistRCPhiEta->Fill(RCeta, RCphi);
454 fHistRhoVSRCPt[fCentBin]->Fill(fRhoVal * rcArea, RCpt);
455
456 fHistRCPt[fCentBin]->Fill(RCpt);
9f6dc7fd 457
458 Double_t ep = RCphi - fEPV0;
459 while (ep < 0) ep += TMath::Pi();
460 while (ep >= TMath::Pi()) ep -= TMath::Pi();
461
462 fHistDeltaPtRCvsEP[fCentBin]->Fill(ep, RCpt - rcArea * fRhoVal);
6a20534a 463 }
6421eeb0 464
465 if (fJetsCont) {
6a20534a 466
467 // Random cones far from leading jet
6421eeb0 468 AliEmcalJet* jet = fJetsCont->GetLeadingJet("rho");
6a20534a 469
470 RCpt = 0;
471 RCeta = 0;
472 RCphi = 0;
6421eeb0 473 GetRandomCone(RCpt, RCeta, RCphi, fTracksCont, fCaloClustersCont, jet);
6a20534a 474 if (RCpt > 0) {
475 if (jet) {
476 Float_t dphi = RCphi - jet->Phi();
477 if (dphi > 4.8) dphi -= TMath::Pi() * 2;
478 if (dphi < -1.6) dphi += TMath::Pi() * 2;
479 fHistRCPtExLJVSDPhiLJ->Fill(RCpt, dphi);
480 }
481 fHistRCPtExLJ[fCentBin]->Fill(RCpt);
482 fHistDeltaPtRCExLJ[fCentBin]->Fill(RCpt - rcArea * fRhoVal);
483 }
56f370b0 484
485 //partial exclusion
486 if(fBeamType == kpA) {
487
488 RCpt = 0;
489 RCeta = 0;
490 RCphi = 0;
6421eeb0 491 GetRandomCone(RCpt, RCeta, RCphi, fTracksCont, fCaloClustersCont, jet, kTRUE);
56f370b0 492
493 if (RCpt > 0) {
494 if (jet) {
495 Float_t dphi = RCphi - jet->Phi();
496 if (dphi > 4.8) dphi -= TMath::Pi() * 2;
497 if (dphi < -1.6) dphi += TMath::Pi() * 2;
498 fHistRCPtExPartialLJVSDPhiLJ->Fill(RCpt, dphi);
499 }
500 fHistRCPtExPartialLJ[fCentBin]->Fill(RCpt);
501 fHistDeltaPtRCExPartialLJ[fCentBin]->Fill(RCpt - rcArea * fRhoVal);
502 }
503 }
a487deae 504 }
a487deae 505 }
6a20534a 506 }
507
508 // Random cones with randomized particles
6421eeb0 509 if (fRandTracksCont || fRandCaloClustersCont) {
a487deae 510 RCpt = 0;
511 RCeta = 0;
512 RCphi = 0;
6421eeb0 513 GetRandomCone(RCpt, RCeta, RCphi, fRandTracksCont, fRandCaloClustersCont, 0);
a487deae 514 if (RCpt > 0) {
1f9c287f 515 fHistRCPtRand[fCentBin]->Fill(RCpt);
a487deae 516 fHistDeltaPtRCRand[fCentBin]->Fill(RCpt - rcArea * fRhoVal);
517 }
518 }
6a20534a 519
a487deae 520 // ************
521 // Embedding
522 // _________________________________
523
6421eeb0 524 if (fEmbJetsCont) {
6a20534a 525
6421eeb0 526 AliEmcalJet *embJet = NextEmbeddedJet(kTRUE);
6a20534a 527
528 while (embJet != 0) {
6421eeb0 529 TLorentzVector mom;
530 fEmbJetsCont->GetLeadingHadronMomentum(mom,embJet);
6a20534a 531
6421eeb0 532 Double_t distLeading2Jet = TMath::Sqrt((embJet->Eta() - mom.Eta()) * (embJet->Eta() - mom.Eta()) + (embJet->Phi() - mom.Phi()) * (embJet->Phi() - mom.Phi()));
6a20534a 533
534 fHistEmbPartPtvsJetPt[fCentBin]->Fill(embJet->MCPt(), embJet->Pt());
ca5c29fa 535 fHistEmbPartPtvsJetCorrPt[fCentBin]->Fill(embJet->MCPt(), embJet->Pt() - embJet->Area() * fRhoVal);
6421eeb0 536 fHistLeadPartPhiEta->Fill(mom.Eta(), mom.Phi());
6a20534a 537 fHistDistLeadPart2JetAxis[fCentBin]->Fill(distLeading2Jet);
538
6421eeb0 539 fHistEmbJetsPtArea[fCentBin]->Fill(embJet->Area(), embJet->Pt(), mom.Pt());
540 fHistEmbJetsCorrPtArea[fCentBin]->Fill(embJet->Area(), embJet->Pt() - fRhoVal * embJet->Area(), mom.Pt());
6a20534a 541 fHistEmbJetsPhiEta->Fill(embJet->Eta(), embJet->Phi());
2103dc6a 542 fHistJetPtvsJetCorrPt[fCentBin]->Fill(embJet->Pt(), embJet->Pt() - fRhoVal * embJet->Area());
6a20534a 543
544 fHistEmbBkgArea[fCentBin]->Fill(embJet->Area(), embJet->Pt() - embJet->MCPt());
545 fHistRhoVSEmbBkg[fCentBin]->Fill(fRhoVal * embJet->Area(), embJet->Pt() - embJet->MCPt());
546 fHistDeltaPtEmbArea[fCentBin]->Fill(embJet->Area(), embJet->Pt() - embJet->Area() * fRhoVal - embJet->MCPt());
7090b9be 547
548 Double_t ep = embJet->Phi() - fEPV0;
549 while (ep < 0) ep += TMath::Pi();
550 while (ep >= TMath::Pi()) ep -= TMath::Pi();
551
552 fHistDeltaPtEmbvsEP[fCentBin]->Fill(ep, embJet->Pt() - embJet->Area() * fRhoVal - embJet->MCPt());
f660c2d6 553
6a20534a 554 embJet = NextEmbeddedJet();
a487deae 555 }
a487deae 556 }
557
558 return kTRUE;
559}
560
561//________________________________________________________________________
6421eeb0 562AliEmcalJet* AliAnalysisTaskDeltaPt::NextEmbeddedJet(Bool_t reset)
a487deae 563{
6421eeb0 564 // Get the next accepted embedded jet.
f660c2d6 565
6421eeb0 566 Int_t i = reset ? 0 : -1;
a487deae 567
6421eeb0 568 AliEmcalJet* jet = fEmbJetsCont->GetNextAcceptJet(i);
569 while (jet && jet->MCPt() < fMCJetPtThreshold) jet = fEmbJetsCont->GetNextAcceptJet();
f660c2d6 570
6421eeb0 571 return jet;
a487deae 572}
573
574//________________________________________________________________________
575void AliAnalysisTaskDeltaPt::GetRandomCone(Float_t &pt, Float_t &eta, Float_t &phi,
6421eeb0 576 AliParticleContainer* tracks, AliClusterContainer* clusters,
577 AliEmcalJet *jet, Bool_t bPartialExclusion) const
a487deae 578{
579 // Get rigid cone.
580
6421eeb0 581 eta = -999;
582 phi = -999;
a487deae 583 pt = 0;
584
585 if (!tracks && !clusters)
586 return;
587
588 Float_t LJeta = 999;
589 Float_t LJphi = 999;
590
591 if (jet) {
592 LJeta = jet->Eta();
593 LJphi = jet->Phi();
594 }
595
6421eeb0 596 Float_t maxEta = fConeMaxEta;
597 Float_t minEta = fConeMinEta;
598 Float_t maxPhi = fConeMaxPhi;
599 Float_t minPhi = fConeMinPhi;
a487deae 600
601 if (maxPhi > TMath::Pi() * 2) maxPhi = TMath::Pi() * 2;
602 if (minPhi < 0) minPhi = 0;
603
604 Float_t dLJ = 0;
605 Int_t repeats = 0;
56f370b0 606 Bool_t reject = kTRUE;
a487deae 607 do {
608 eta = gRandom->Rndm() * (maxEta - minEta) + minEta;
609 phi = gRandom->Rndm() * (maxPhi - minPhi) + minPhi;
610 dLJ = TMath::Sqrt((LJeta - eta) * (LJeta - eta) + (LJphi - phi) * (LJphi - phi));
56f370b0 611
612 if(bPartialExclusion) {
613 reject = kFALSE;
614
615 TRandom3 rnd;
616 rnd.SetSeed(0);
617
618 Double_t ncoll = GetNColl();
619
620 Double_t prob = 0.;
621 if(ncoll>0)
622 prob = 1./ncoll;
623
624 if(rnd.Rndm()<=prob) reject = kTRUE; //reject cone
625 }
626
a487deae 627 repeats++;
56f370b0 628 } while (dLJ < fMinRC2LJ && repeats < 999 && reject);
a487deae 629
630 if (repeats == 999) {
631 AliWarning(Form("%s: Could not get random cone!", GetName()));
632 return;
633 }
634
635 if (clusters) {
6421eeb0 636 AliVCluster* cluster = clusters->GetNextAcceptCluster(0);
637 while (cluster) {
a487deae 638 TLorentzVector nPart;
639 cluster->GetMomentum(nPart, const_cast<Double_t*>(fVertex));
a487deae 640
6421eeb0 641 Float_t cluseta = nPart.Eta();
642 Float_t clusphi = nPart.Phi();
643
644 if (TMath::Abs(clusphi - phi) > TMath::Abs(clusphi - phi + 2 * TMath::Pi()))
645 clusphi += 2 * TMath::Pi();
646 if (TMath::Abs(clusphi - phi) > TMath::Abs(clusphi - phi - 2 * TMath::Pi()))
647 clusphi -= 2 * TMath::Pi();
648
649 Float_t d = TMath::Sqrt((cluseta - eta) * (cluseta - eta) + (clusphi - phi) * (clusphi - phi));
650 if (d <= fConeRadius)
a487deae 651 pt += nPart.Pt();
6421eeb0 652
653 cluster = clusters->GetNextAcceptCluster();
a487deae 654 }
655 }
656
657 if (tracks) {
6421eeb0 658 AliVParticle* track = tracks->GetNextAcceptParticle(0);
659 while(track) {
a487deae 660 Float_t tracketa = track->Eta();
661 Float_t trackphi = track->Phi();
662
663 if (TMath::Abs(trackphi - phi) > TMath::Abs(trackphi - phi + 2 * TMath::Pi()))
664 trackphi += 2 * TMath::Pi();
665 if (TMath::Abs(trackphi - phi) > TMath::Abs(trackphi - phi - 2 * TMath::Pi()))
666 trackphi -= 2 * TMath::Pi();
667
668 Float_t d = TMath::Sqrt((tracketa - eta) * (tracketa - eta) + (trackphi - phi) * (trackphi - phi));
6421eeb0 669 if (d <= fConeRadius)
a487deae 670 pt += track->Pt();
6421eeb0 671
672 track = tracks->GetNextAcceptParticle();
a487deae 673 }
674 }
675}
676
6421eeb0 677//________________________________________________________________________
678void AliAnalysisTaskDeltaPt::SetConeEtaPhiEMCAL()
679{
680 // Set default cuts for full cones
681
682 SetConeEtaLimits(-0.7+fConeRadius,0.7-fConeRadius);
683 SetConePhiLimits(1.4+fConeRadius,TMath::Pi()-fConeRadius);
684}
685
686//________________________________________________________________________
687void AliAnalysisTaskDeltaPt::SetConeEtaPhiTPC()
688{
689 // Set default cuts for charged cones
690
691 SetConeEtaLimits(-0.9+fConeRadius, 0.9-fConeRadius);
692 SetConePhiLimits(-10, 10);
693}
694
a487deae 695//________________________________________________________________________
696void AliAnalysisTaskDeltaPt::ExecOnce()
697{
698 // Initialize the analysis.
699
9239b066 700 AliAnalysisTaskEmcalJet::ExecOnce();
a487deae 701
6421eeb0 702 if (fTracksCont && fTracksCont->GetArray() == 0) fTracksCont = 0;
703 if (fCaloClustersCont && fCaloClustersCont->GetArray() == 0) fCaloClustersCont = 0;
704 if (fEmbTracksCont && fEmbTracksCont->GetArray() == 0) fEmbTracksCont = 0;
705 if (fEmbCaloClustersCont && fEmbCaloClustersCont->GetArray() == 0) fEmbCaloClustersCont = 0;
706 if (fRandTracksCont && fRandTracksCont->GetArray() == 0) fRandTracksCont = 0;
707 if (fRandCaloClustersCont && fRandCaloClustersCont->GetArray() == 0) fRandCaloClustersCont = 0;
7cd832c7 708 if (fJetsCont && fJetsCont->GetArray() == 0) fJetsCont = 0;
709 if (fEmbJetsCont && fEmbJetsCont->GetArray() == 0) fEmbJetsCont = 0;
a487deae 710
a487deae 711 if (fRCperEvent < 0) {
6421eeb0 712 Double_t area = (fConeMaxEta - fConeMinEta) * (fConeMaxPhi - fConeMinPhi);
713 Double_t rcArea = TMath::Pi() * fConeRadius * fConeRadius;
714 fRCperEvent = TMath::FloorNint(area / rcArea - 0.5);
a487deae 715 if (fRCperEvent == 0)
716 fRCperEvent = 1;
717 }
718
719 if (fMinRC2LJ < 0)
6421eeb0 720 fMinRC2LJ = fConeRadius * 1.5;
a487deae 721
6421eeb0 722 const Float_t maxDist = TMath::Max(fConeMaxPhi - fConeMinPhi, fConeMaxEta - fConeMinEta) / 2;
a487deae 723 if (fMinRC2LJ > maxDist) {
724 AliWarning(Form("The parameter fMinRC2LJ = %f is too large for the considered acceptance. "
725 "Will use fMinRC2LJ = %f", fMinRC2LJ, maxDist));
726 fMinRC2LJ = maxDist;
727 }
728}
729
56f370b0 730//________________________________________________________________________
731Double_t AliAnalysisTaskDeltaPt::GetNColl() const {
732 // Get NColl - returns value of corresponding bin
733 // only works for pA
734 // values taken from V0A slicing https://twiki.cern.ch/twiki/bin/viewauth/ALICE/PACentStudies#Tables_with_centrality_bins_from
735
736 if(fBeamType == kpA) {
737
738 const Int_t nNCollBins = 7;
739
740 Double_t centMin[nNCollBins] = {0.,5.,10.,20.,40.,60.,80.};
741 Double_t centMax[nNCollBins] = {5.,10.,20.,40.,60.,80.,100.};
742
743 Double_t nColl[nNCollBins] = {14.7,13.,11.7,9.38,6.49,3.96,1.52};
744
745 for(Int_t i = 0; i<nNCollBins; i++) {
746 if(fCent>=centMin[i] && fCent<centMax[i])
747 return nColl[i];
748 }
749
750 return -1.;
751 }
752 else {
753 AliWarning(Form("%s: Only works for pA analysis. Returning -1",GetName()));
754 return -1.;
755 }
a487deae 756}