]> git.uio.no Git - u/mrichter/AliRoot.git/blame - PWGJE/EMCALJetTasks/UserTasks/AliAnalysisTaskEmcalJetv2QA.cxx
git://git.uio.no / u / mrichter / AliRoot.git / blame
? search:
summary | shortlog | log | commit | commitdiff | tree
blob | blame (incremental) | history | HEAD
Merge branch 'feature-movesplit'
[u/mrichter/AliRoot.git] / PWGJE / EMCALJetTasks / UserTasks / AliAnalysisTaskEmcalJetv2QA.cxx
CommitLineData
51a8c8fd 1// Jet v2 task using QA method, based on jet sample task (S.Aiola).
2//
3// Authors: Jason Mueller (CERN summer student 2014) & Alice Ohlson
4
5
6#include <TClonesArray.h>
7#include <TH1F.h>
8#include <TH2F.h>
9#include <TH3F.h>
10#include <TList.h>
11#include <TLorentzVector.h>
12#include <TF1.h>
13
14#include "AliVCluster.h"
15#include "AliAODCaloCluster.h"
16#include "AliESDCaloCluster.h"
17#include "AliVTrack.h"
18#include "AliEmcalJet.h"
19#include "AliRhoParameter.h"
20#include "AliLog.h"
21#include "AliJetContainer.h"
22#include "AliParticleContainer.h"
23#include "AliClusterContainer.h"
24#include "AliPicoTrack.h"
25
26#include "AliAnalysisTaskEmcalJetv2QA.h"
27
2b010c34 28using std::cout;
29using std::endl;
30
51a8c8fd 31ClassImp(AliAnalysisTaskEmcalJetv2QA)
32
33//________________________________________________________________________
34AliAnalysisTaskEmcalJetv2QA::AliAnalysisTaskEmcalJetv2QA() :
35AliAnalysisTaskEmcalJet("AliAnalysisTaskEmcalJetv2QA", kTRUE),
36 nCentBins(0),
de2a3fac 37 nCentBins1(1),
51a8c8fd 38 centBins(0x0),
de2a3fac 39 nJetPtBins(0),
40 nJetPtBins1(1),
51a8c8fd 41 jetPtBins(0x0),
42 fJetv2(0x0),
43 doPtWeight(kFALSE),
44 fHistTracksPt(0),
45 fHistClustersPt(0),
46 fHistLeadingJetPt(0),
47 fHistLeadingJetPtCorr(0),
48 fHistJetsPhiEta(0),
49 fHistJetsPtArea(0),
50 fHistJetsPtLeadHad(0),
51 fHistJetsCorrPtArea(0),
52 fHistPtDEtaDPhiTrackClus(0),
53 fHistPtDEtaDPhiClusTrack(0),
54 fDPhiJet(0),
55 fDPhiJetPythia(0),
56 fDPhiEP(0),
57 hGx(0),
58 hGy2(0),
59 hGxGy2(0),
60 hGy4(0),
61 hGx2(0),
62 hGx2Gy2(0),
63 hGxGy4(0),
64 hGy6(0),
65 hGx2Gy4(0),
66 hGxGy6(0),
67 hGy8(0),
68 hGy(0),
69 hN(0),
70 htv2std(0),
71 htjv2std(0),
72 htj2v2std(0),
73 hV0jv2std(0),
74 htdPsi(0),
75 htjdPsi(0),
76 htj2dPsi(0),
77 hV0jdPsi(0),
78 hAx(0),
79 hAxDijet(0),
80 hQx(0),
81 hQy(0),
82 hEventData(0),
83 hNTracks(0),
84 hNTracksCent(0),
85 hGxTracks(0),
86 hGyTracks(0),
87 hGy2Tracks(0),
88 hGxGy2Tracks(0),
89 hGy4Tracks(0),
90 htEPRes(0),
91 htjEPRes(0),
92 htj2EPRes(0),
93 fJetsCont(0),
94 fTracksCont(0),
95 fCaloClustersCont(0)
96{
97 // Default constructor.
98
51a8c8fd 99 SetMakeGeneralHistograms(kTRUE);
100}
101
102//________________________________________________________________________
103AliAnalysisTaskEmcalJetv2QA::AliAnalysisTaskEmcalJetv2QA(const char *name) :
104 AliAnalysisTaskEmcalJet(name, kTRUE),
105 nCentBins(0),
de2a3fac 106 nCentBins1(1),
51a8c8fd 107 centBins(0x0),
de2a3fac 108 nJetPtBins(0),
109 nJetPtBins1(1),
51a8c8fd 110 jetPtBins(0x0),
111 fJetv2(0x0),
112 doPtWeight(kFALSE),
113 fHistTracksPt(0),
114 fHistClustersPt(0),
115 fHistLeadingJetPt(0),
116 fHistLeadingJetPtCorr(0),
117 fHistJetsPhiEta(0),
118 fHistJetsPtArea(0),
119 fHistJetsPtLeadHad(0),
120 fHistJetsCorrPtArea(0),
121 fHistPtDEtaDPhiTrackClus(0),
122 fHistPtDEtaDPhiClusTrack(0),
123 fDPhiJet(0),
124 fDPhiJetPythia(0),
125 fDPhiEP(0),
126 hGx(0),
127 hGy2(0),
128 hGxGy2(0),
129 hGy4(0),
130 hGx2(0),
131 hGx2Gy2(0),
132 hGxGy4(0),
133 hGy6(0),
134 hGx2Gy4(0),
135 hGxGy6(0),
136 hGy8(0),
137 hGy(0),
138 hN(0),
139 htv2std(0),
140 htjv2std(0),
141 htj2v2std(0),
142 hV0jv2std(0),
143 htdPsi(0),
144 htjdPsi(0),
145 htj2dPsi(0),
146 hV0jdPsi(0),
147 hAx(0),
148 hAxDijet(0),
149 hQx(0),
150 hQy(0),
151 hEventData(0),
152 hNTracks(0),
153 hNTracksCent(0),
154 hGxTracks(0),
155 hGyTracks(0),
156 hGy2Tracks(0),
157 hGxGy2Tracks(0),
158 hGy4Tracks(0),
159 htEPRes(0),
160 htjEPRes(0),
161 htj2EPRes(0),
162 fJetsCont(0),
163 fTracksCont(0),
164 fCaloClustersCont(0)
165{
166 // Standard constructor.
167
de2a3fac 168 // default binning
169 Double_t centBinsTemp[7] = {0,5,10,20,30,40,50};
170 SetCentBins(6,centBinsTemp);
171 Double_t jetPtBinsTemp[6] = {40.,50.,70.,90.,120.,200.};
172 SetJetPtBins(5,jetPtBinsTemp);
51a8c8fd 173
174 SetMakeGeneralHistograms(kTRUE);
175}
176
177//________________________________________________________________________
178AliAnalysisTaskEmcalJetv2QA::~AliAnalysisTaskEmcalJetv2QA()
179{
180 // Destructor.
181}
182
183//________________________________________________________________________
184void AliAnalysisTaskEmcalJetv2QA::UserCreateOutputObjects()
185{
186 // Create user output.
187
188 AliAnalysisTaskEmcalJet::UserCreateOutputObjects();
189
190 fJetsCont = GetJetContainer(0);
191 if(fJetsCont) { //get particles and clusters connected to jets
192 fTracksCont = fJetsCont->GetParticleContainer();
193 fCaloClustersCont = fJetsCont->GetClusterContainer();
194 } else { //no jets, just analysis tracks and clusters
195 fTracksCont = GetParticleContainer(0);
196 fCaloClustersCont = GetClusterContainer(0);
197 }
198
51a8c8fd 199 fTracksCont->SetClassName("AliVTrack");
200 fCaloClustersCont->SetClassName("AliVCluster");
201
202 fHistTracksPt = new TH1F("fHistTracksPt","fHistTracksPt", fNbins / 2, fMinBinPt, fMaxBinPt / 2);
203 fHistTracksPt->GetXaxis()->SetTitle("p_{T,track} (GeV/c)");
204 fHistTracksPt->GetYaxis()->SetTitle("counts");
205 fOutput->Add(fHistTracksPt);
206
207 fHistClustersPt = new TH1F("fHistClustersPt", "fHistClustersPt", fNbins / 2, fMinBinPt, fMaxBinPt / 2);
208 fHistClustersPt->GetXaxis()->SetTitle("p_{T,clus} (GeV/c)");
209 fHistClustersPt->GetYaxis()->SetTitle("counts");
210 fOutput->Add(fHistClustersPt);
211
212 fHistLeadingJetPt = new TH1F("fHistLeadingJetPt", "fHistLeadingJetPt", fNbins, fMinBinPt, fMaxBinPt);
213 fHistLeadingJetPt->GetXaxis()->SetTitle("p_{T}^{raw} (GeV/c)");
214 fHistLeadingJetPt->GetYaxis()->SetTitle("counts");
215 fOutput->Add(fHistLeadingJetPt);
216
217 fHistLeadingJetPtCorr = new TH1F("fHistLeadingJetPtCorr", "fHistLeadingJetPtCorr", fNbins, fMinBinPt, fMaxBinPt);
218 fHistLeadingJetPtCorr->GetXaxis()->SetTitle("p_{T}^{corr} (GeV/c)");
219 fHistLeadingJetPtCorr->GetYaxis()->SetTitle("counts");
220 fOutput->Add(fHistLeadingJetPtCorr);
221
222 fHistJetsPhiEta = new TH2F("fHistJetsPhiEta", "fHistJetsPhiEta", 50, -1, 1, 101, 0, TMath::Pi()*2 + TMath::Pi()/200);
223 fHistJetsPhiEta->GetXaxis()->SetTitle("#eta");
224 fHistJetsPhiEta->GetYaxis()->SetTitle("#phi");
225 fOutput->Add(fHistJetsPhiEta);
226
227 fHistJetsPtArea = new TH2F("fHistJetsPtArea", "fHistJetsPtArea", fNbins, fMinBinPt, fMaxBinPt, 30, 0, 3);
228 fHistJetsPtArea->GetXaxis()->SetTitle("p_{T}^{raw} (GeV/c)");
229 fHistJetsPtArea->GetYaxis()->SetTitle("area");
230 fOutput->Add(fHistJetsPtArea);
231
232 fHistJetsPtLeadHad = new TH2F("fHistJetsPtLeadHad", "fHistJetsPtLeadHad", fNbins, fMinBinPt, fMaxBinPt, fNbins / 2, fMinBinPt, fMaxBinPt / 2);
233 fHistJetsPtLeadHad->GetXaxis()->SetTitle("p_{T}^{raw} (GeV/c)");
234 fHistJetsPtLeadHad->GetYaxis()->SetTitle("p_{T,lead} (GeV/c)");
235 fHistJetsPtLeadHad->GetZaxis()->SetTitle("counts");
236 fOutput->Add(fHistJetsPtLeadHad);
237
238 fHistJetsCorrPtArea = new TH2F("fHistJetsCorrPtArea", "fHistJetsCorrPtArea", fNbins*2, -fMaxBinPt, fMaxBinPt, 30, 0, 3);
239 fHistJetsCorrPtArea->GetXaxis()->SetTitle("p_{T}^{corr} [GeV/c]");
240 fHistJetsCorrPtArea->GetYaxis()->SetTitle("area");
241 fOutput->Add(fHistJetsCorrPtArea);
242
243 fHistPtDEtaDPhiTrackClus = new TH3F("fHistPtDEtaDPhiTrackClus","fHistPtDEtaDPhiTrackClus;#it{p}_{T}^{track};#Delta#eta;#Delta#varphi",100,0.,100.,100,-0.1,0.1,100,-0.1,0.1);
244 fOutput->Add(fHistPtDEtaDPhiTrackClus);
245
246 fHistPtDEtaDPhiClusTrack = new TH3F("fHistPtDEtaDPhiClusTrack","fHistPtDEtaDPhiClusTrack;#it{p}_{T}^{clus};#Delta#eta;#Delta#varphi",100,0.,100.,100,-0.1,0.1,100,-0.1,0.1);
247 fOutput->Add(fHistPtDEtaDPhiClusTrack);
248
249 fDPhiJet = new TH1F("fDPhiJet","fDPhiJet",90, -TMath::Pi()/3, 5*TMath::Pi()/3);
250 fDPhiJet->GetXaxis()->SetTitle("#Delta#varphi");
251 fDPhiJet->GetYaxis()->SetTitle("counts");
252 fOutput->Add(fDPhiJet);
253
254 fDPhiJetPythia = new TH1F("fDPhiJetPythia","fDPhiJetPythia",90, -TMath::Pi()/3, 5*TMath::Pi()/3);
255 fDPhiJetPythia->GetXaxis()->SetTitle("#Delta#varphi");
256 fDPhiJetPythia->GetYaxis()->SetTitle("counts");
257 fOutput->Add(fDPhiJetPythia);
258
259 fDPhiEP = new TH1F("fDPhiEP","fDPhiEP",90, 0, 2*TMath::Pi());
260 fDPhiEP->GetXaxis()->SetTitle("#Delta#varphi");
261 fDPhiEP->GetYaxis()->SetTitle("counts");
262 fOutput->Add(fDPhiEP);
263
264 hGx = new TH2D("hGx", "Gx v Centrality v JetPt", nCentBins, centBins, nJetPtBins, jetPtBins);
265 hGx->GetXaxis()->SetTitle("Centrality (%)");
266 hGx->GetYaxis()->SetTitle("Leading Jet pT (GeV/c)");
267 fOutput->Add(hGx);
268
269 hGy2 = new TH2D("hGy2", "Gy2 v Centrality v JetPt", nCentBins, centBins, nJetPtBins, jetPtBins);
270 hGy2->GetXaxis()->SetTitle("Centrality (%)");
271 hGy2->GetYaxis()->SetTitle("Leading Jet pT (GeV/c)");
272 fOutput->Add(hGy2);
273
274 hGxGy2 = new TH2D("hGxGy2", "GxGy2 v Centrality v JetPt", nCentBins, centBins, nJetPtBins, jetPtBins);
275 hGxGy2->GetXaxis()->SetTitle("Centrality (%)");
276 hGxGy2->GetYaxis()->SetTitle("Leading Jet pT (GeV/c)");
277 fOutput->Add(hGxGy2);
278
279 hGy4 = new TH2D("hGy4", "Gy4 v Centrality v JetPt", nCentBins, centBins, nJetPtBins, jetPtBins);
280 hGy4->GetXaxis()->SetTitle("Centrality (%)");
281 hGy4->GetYaxis()->SetTitle("Leading Jet pT (GeV/c)");
282 fOutput->Add(hGy4);
283
284 hGx2 = new TH2D("hGx2", "Gx2 v Centrality v JetPt", nCentBins, centBins, nJetPtBins, jetPtBins);
285 hGx2->GetXaxis()->SetTitle("Centrality (%)");
286 hGx2->GetYaxis()->SetTitle("Leading Jet pT (GeV/c)");
287 fOutput->Add(hGx2);
288
289 hGx2Gy2 = new TH2D("hGx2Gy2", "Gx2Gy2 v Centrality v JetPt", nCentBins, centBins, nJetPtBins, jetPtBins);
290 hGx2Gy2->GetXaxis()->SetTitle("Centrality (%)");
291 hGx2Gy2->GetYaxis()->SetTitle("Leading Jet pT (GeV/c)");
292 fOutput->Add(hGx2Gy2);
293
294 hGxGy4 = new TH2D("hGxGy4", "GxGy4 v Centrality v JetPt", nCentBins, centBins, nJetPtBins, jetPtBins);
295 hGxGy4->GetXaxis()->SetTitle("Centrality (%)");
296 hGxGy4->GetYaxis()->SetTitle("Leading Jet pT (GeV/c)");
297 fOutput->Add(hGxGy4);
298
299 hGy6 = new TH2D("hGy6", "Gy6 v Centrality v JetPt", nCentBins, centBins, nJetPtBins, jetPtBins);
300 hGy6->GetXaxis()->SetTitle("Centrality (%)");
301 hGy6->GetYaxis()->SetTitle("Leading Jet pT (GeV/c)");
302 fOutput->Add(hGy6);
303
304 hGx2Gy4 = new TH2D("hGx2Gy4", "Gx2Gy4 v Centrality v JetPt", nCentBins, centBins, nJetPtBins, jetPtBins);
305 hGx2Gy4->GetXaxis()->SetTitle("Centrality (%)");
306 hGx2Gy4->GetYaxis()->SetTitle("Leading Jet pT (GeV/c)");
307 fOutput->Add(hGx2Gy4);
308
309 hGxGy6 = new TH2D("hGxGy6", "GxGy6 v Centrality v JetPt", nCentBins, centBins, nJetPtBins, jetPtBins);
310 hGxGy6->GetXaxis()->SetTitle("Centrality (%)");
311 hGxGy6->GetYaxis()->SetTitle("Leading Jet pT (GeV/c)");
312 fOutput->Add(hGxGy6);
313
314 hGy8 = new TH2D("hGy8", "Gy8 v Centrality v JetPt", nCentBins, centBins, nJetPtBins, jetPtBins);
315 hGy8->GetXaxis()->SetTitle("Centrality (%)");
316 hGy8->GetYaxis()->SetTitle("Leading Jet pT (GeV/c)");
317 fOutput->Add(hGy8);
318
319 hGy = new TH2D("hGy", "Gy v Centrality v JetPt", nCentBins, centBins, nJetPtBins, jetPtBins);
320 hGy->GetXaxis()->SetTitle("Centrality (%)");
321 hGy->GetYaxis()->SetTitle("Leading Jet pT (GeV/c)");
322 fOutput->Add(hGy);
323
324 hN = new TH2D("hN", "N v Centrality v JetPt", nCentBins, centBins, nJetPtBins, jetPtBins);
325 hN->GetXaxis()->SetTitle("Centrality (%)");
326 hN->GetYaxis()->SetTitle("Leading Jet pT (GeV/c)");
327 fOutput->Add(hN);
328
329 htv2std = new TH2D("htv2std", "v2std v Centrality v JetPt w/o jets", nCentBins, centBins, nJetPtBins, jetPtBins);
330 htv2std->GetXaxis()->SetTitle("Centrality (%)");
331 htv2std->GetYaxis()->SetTitle("Leading Jet pT (GeV/c)");
332 fOutput->Add(htv2std);
333
334 htjv2std = new TH2D("htjv2std", "v2std v Centrality v JetPt w/ jets", nCentBins, centBins, nJetPtBins, jetPtBins);
335 htjv2std->GetXaxis()->SetTitle("Centrality (%)");
336 htjv2std->GetYaxis()->SetTitle("Leading Jet pT (GeV/c)");
337 fOutput->Add(htjv2std);
338
339 htj2v2std = new TH2D("htj2v2std", "v2std v Centrality v JetPt w/ trackPt < 2 GeV", nCentBins, centBins, nJetPtBins, jetPtBins);
340 htj2v2std->GetXaxis()->SetTitle("Centrality (%)");
341 htj2v2std->GetYaxis()->SetTitle("Leading Jet pT (GeV/c)");
342 fOutput->Add(htj2v2std);
343
344 hV0jv2std = new TH2D("hV0jv2std", "v2std v Centrality v JetPt", nCentBins, centBins, nJetPtBins, jetPtBins);
345 hV0jv2std->GetXaxis()->SetTitle("Centrality (%)");
346 hV0jv2std->GetYaxis()->SetTitle("Leading Jet pT (GeV/c)");
347 fOutput->Add(hV0jv2std);
348
349 Int_t ndpsibins = 100;
350 Double_t dpsibins[101];
351 for(Int_t t = 0; t < 101; t++) dpsibins[t] = TMath::Pi()*t/50.;
352
353 htdPsi = new TH3D("htdPsi", "JetAxis - EventPlane w/o jets", nCentBins, centBins, nJetPtBins, jetPtBins, ndpsibins, dpsibins);
354 htdPsi->GetZaxis()->SetTitle("#Psi_{jet} - #Psi_{EP}");
355 fOutput->Add(htdPsi);
356
357 htjdPsi = new TH3D("htjdPsi", "JetAxis - EventPlane w/ jets", nCentBins, centBins, nJetPtBins, jetPtBins, ndpsibins, dpsibins);
358 htjdPsi->GetZaxis()->SetTitle("#Psi_{jet} - #Psi_{EP}");
359 fOutput->Add(htjdPsi);
360
361 htj2dPsi = new TH3D("htj2dPsi", "JetAxis - EventPlane w/ trackPt < 2 GeV", nCentBins, centBins, nJetPtBins, jetPtBins, ndpsibins, dpsibins);
362 htj2dPsi->GetZaxis()->SetTitle("#Psi_{jet} - #Psi_{EP}");
363 fOutput->Add(htj2dPsi);
364
365 hV0jdPsi = new TH3D("hV0jdPsi", "JetAxis - EventPlane", nCentBins, centBins, nJetPtBins, jetPtBins, ndpsibins, dpsibins);
366 hV0jdPsi->GetZaxis()->SetTitle("#Psi_{jet} - #Psi_{EP}");
367 fOutput->Add(hV0jdPsi);
368
369 hQx = new TH2D("hQx", "Qx Distribution in EP frame", 100, -0.3, 0.3, nCentBins, centBins);
370 hQx->GetXaxis()->SetTitle("Qx");
371 hQx->GetYaxis()->SetTitle("Centrality (%)");
372 fOutput->Add(hQx);
373
374 hQy = new TH2D("hQy", "Qy Distribution in EP frame", 100, -0.3, 0.3, nCentBins, centBins);
375 hQy->GetXaxis()->SetTitle("Qy");
376 hQy->GetYaxis()->SetTitle("Centrality (%)");
377 fOutput->Add(hQy);
378
379 hAx = new TH2D("hAx", "Ax Distribution in EP frame w/o Dijets", 100, -35, 70, nJetPtBins, jetPtBins);
380 hAx->GetXaxis()->SetTitle("Ax");
381 hAx->GetYaxis()->SetTitle("Leading Jet pT (GeV/c)");
382 fOutput->Add(hAx);
383
384 hAxDijet = new TH2D("hAxDijet", "Ax Distribution in EP frame w/ Dijets", 100, -35, 70, nJetPtBins, jetPtBins);
385 hAxDijet->GetXaxis()->SetTitle("Ax");
386 hAxDijet->GetYaxis()->SetTitle("Leading Jet pT (GeV/c)");
387 fOutput->Add(hAxDijet);
388
389 hEventData = new TH1F("hEventData","Events Kept and Discarded", 9, 0, 9);
390 hEventData->GetYaxis()->SetTitle("counts");
391 fOutput->Add(hEventData);
392
393 hNTracks = new TH2F("hNTracks","Number of Tracks Per Event", 100, 0, 3000, 3, 0, 3);
394 hNTracks->GetXaxis()->SetTitle("# tracks");
395 fOutput->Add(hNTracks);
396
397 hNTracksCent = new TH2D("hNTracksCent", "NTracks by centrality", 100, 0, 3000, (Int_t)centBins[nCentBins], centBins[0], centBins[nCentBins]);
398 hNTracksCent->GetXaxis()->SetTitle("# tracks");
399 hNTracksCent->GetYaxis()->SetTitle("Centrality (%)");
400 fOutput->Add(hNTracksCent);
401
402 hGxTracks = new TH2D("hGxTracks", "Gx by NTracks", 200, -200, 200, 100, 0, 3000);
403 hGxTracks->GetXaxis()->SetTitle("Gx");
404 hGxTracks->GetYaxis()->SetTitle("# tracks");
405 fOutput->Add(hGxTracks);
406
407 hGyTracks = new TH2D("hGyTracks", "Gy by NTracks", 200, -200, 200, 100, 0, 3000);
408 hGyTracks->GetXaxis()->SetTitle("Gy");
409 hGyTracks->GetYaxis()->SetTitle("# tracks");
410 fOutput->Add(hGyTracks);
411
412 hGy2Tracks = new TH2D("hGy2Tracks", "Gy2 by NTracks", 100, 0, 20000, 100, 0, 3000);
413 hGy2Tracks->GetXaxis()->SetTitle("Gy2");
414 hGy2Tracks->GetYaxis()->SetTitle("# tracks");
415 fOutput->Add(hGy2Tracks);
416
417 hGxGy2Tracks = new TH2D("hGxGy2Tracks", "GxGy2 by NTracks", 100, -100000, 100000, 100, 0, 3000);
418 hGxGy2Tracks->GetXaxis()->SetTitle("GxGy2");
419 hGxGy2Tracks->GetYaxis()->SetTitle("# tracks");
420 fOutput->Add(hGxGy2Tracks);
421
422 hGy4Tracks = new TH2D("hGy4Tracks", "Gy4 by NTracks", 100, 0, 100000000, 100, 0, 3000);
423 hGy4Tracks->GetXaxis()->SetTitle("Gy4");
424 hGy4Tracks->GetYaxis()->SetTitle("# tracks");
425 fOutput->Add(hGy4Tracks);
426
427 htEPRes = new TH2D("htEPRes", "EP Resolution w/o Jets", nCentBins, centBins, nJetPtBins, jetPtBins);
428 htEPRes->GetXaxis()->SetTitle("Centrality (%)");
429 htEPRes->GetYaxis()->SetTitle("Leading Jet pT (GeV/c)");
430 fOutput->Add(htEPRes);
431
432 htjEPRes = new TH2D("htjEPRes", "EP Resolution w/ Jets", nCentBins, centBins, nJetPtBins, jetPtBins);
433 htjEPRes->GetXaxis()->SetTitle("Centrality (%)");
434 htjEPRes->GetYaxis()->SetTitle("Leading Jet pT (GeV/c)");
435 fOutput->Add(htjEPRes);
436
437 htj2EPRes = new TH2D("htj2EPRes", "EP Resolution w/ trackPT < 2 GeV", nCentBins, centBins, nJetPtBins, jetPtBins);
438 htj2EPRes->GetXaxis()->SetTitle("Centrality (%)");
439 htj2EPRes->GetYaxis()->SetTitle("Leading Jet pT (GeV/c)");
440 fOutput->Add(htj2EPRes);
441
442 PostData(1, fOutput); // Post data for ALL output slots > 0 here.
443}
444
445//________________________________________________________________________
446Bool_t AliAnalysisTaskEmcalJetv2QA::FillHistograms()
447{
448 // Fill histograms.
449
450 return kTRUE;
451}
452
453
454//________________________________________________________________________
455void AliAnalysisTaskEmcalJetv2QA::ExecOnce() {
456
457 AliAnalysisTaskEmcalJet::ExecOnce();
458
459 if (fJetsCont && fJetsCont->GetArray() == 0) fJetsCont = 0;
460 if (fTracksCont && fTracksCont->GetArray() == 0) fTracksCont = 0;
461 if (fCaloClustersCont && fCaloClustersCont->GetArray() == 0) fCaloClustersCont = 0;
462
463}
464
465//________________________________________________________________________
466Bool_t AliAnalysisTaskEmcalJetv2QA::Run() // this part loops over each event
467{
468 // Run analysis code here, if needed. It will be executed before FillHistograms().
469
470 Double_t jetPhi = -999; // angle of leading jet
471 Double_t jetPt = -999; // pt of leading jet
472 Double_t jetArea = -999;
473 Double_t trackPt = -999;
474 Double_t phi = -999; // track phi
475 Double_t dPhi = -999; // track phi - jet axis
476 Double_t dPhiQA = -999; // track phi - EP
477 Double_t tSin = 0; // used for std EP calc
478 Double_t tCos = 0; // used for std EP calc
479 Double_t jSin = 0; // used for std EP calc
480 Double_t jCos = 0; // used for std EP calc
481 Double_t tSin2 = 0; // used for std EP calc with trackPt < 2 GeV
482 Double_t tCos2 = 0; // used for std EP calc with trackPt < 2 GeV
483 Double_t qx = 0; // used for Qx distribution
484 Double_t qy = 0; // used for Qy distribution
485 Double_t ax = 0; // used for Ax distribution
486 Double_t tEP = -999; // EP w/o jets
487 Double_t tjEP = -999; // EP w/ jets
488 Double_t tjEP2 = -999; // EP w/ jets
489 Double_t dPsi = -999; // jet axis - EP
490 Double_t gx = 0; // used for G moment calc
491 Double_t gy = 0; // used for G moment calc
492 Int_t isDijet = 0; // if 0, no dijet. if 1, dijet.
493 Int_t nTracksBkgnd = 0; // used to keep track of number of tracks in background
494 Int_t nTracksJet =0; // used to keep track of number of tracks in jets
495
496 if(fJetsCont)
497 {
498
499 AliEmcalJet *jettest = fJetsCont->GetNextAcceptJet(0);
500 while(jettest) {
501
502 fHistJetsPtArea->Fill(jettest->Pt(), jettest->Area());
503 fHistJetsPhiEta->Fill(jettest->Eta(), jettest->Phi());
504
505 Float_t ptLeading = fJetsCont->GetLeadingHadronPt(jettest);
506 fHistJetsPtLeadHad->Fill(jettest->Pt(), ptLeading);
507
508 if (fHistJetsCorrPtArea) {
509 Float_t corrPt = jettest->Pt() - fJetsCont->GetRhoVal() * jettest->Area();
510 fHistJetsCorrPtArea->Fill(corrPt, jettest->Area());
511 }
512 jettest = fJetsCont->GetNextAcceptJet();
513 }
514
515
516 AliEmcalJet *jet = fJetsCont->GetLeadingJet();
517 if(jet)
518 {
519 jetPhi = jet->Phi(); // get leading jet phi (jet axis)
520 jetPt = jet->Pt(); // get leading jet pT to filter out low pT jet events
521 jetArea = jet->Area();
522 }
523 }
524
525 // event cuts
526 if(!fTracksCont)
527 {
528 hEventData->Fill("!fTracksCont",1);
529 return kFALSE;
530 }
531 if(!fJetsCont)
532 {
533 hEventData->Fill("!fJetsCont",1);
534 return kFALSE;
535 }
536 if(jetPt == -999)
537 {
538 hEventData->Fill("jetPt=-999",1);
539 return kFALSE;
540 }
541 if(jetPt < jetPtBins[0])
542 {
543 hEventData->Fill("leadingJetPt<jetPtMin",1);
544 return kFALSE;
545 }
546 if(jetPt > jetPtBins[nJetPtBins])
547 {
548 hEventData->Fill("leadingJetPt>jetPtMax",1);
549 return kFALSE;
550 }
551 if(fCent < centBins[0])
552 {
553 hEventData->Fill("cent<centMin",1);
554 return kFALSE;
555 }
556 if(fCent > centBins[nCentBins])
557 {
558 hEventData->Fill("cent>centMax",1);
559 return kFALSE;
560 }
561 hEventData->Fill("good event",1);
562
563 AliEmcalJet *dijet = fJetsCont->GetNextAcceptJet(0); // check for dijet events
564 while(dijet)
565 {
de2a3fac 566 if(dijet->Pt() > jetPt*2./3. && fabs(jetPhi-dijet->Phi()-TMath::Pi()) < 0.4) // loop over jets with pT>50 and exclude leading jet and check that angular separation is < 0.4
51a8c8fd 567 isDijet = 1;
568 dijet = fJetsCont->GetNextAcceptJet();
569 }
570
571 if (fCaloClustersCont)
572 {
573 AliVCluster *cluster = fCaloClustersCont->GetNextAcceptCluster(0);
574 while(cluster) {
575 TLorentzVector nPart;
576 cluster->GetMomentum(nPart, fVertex);
577 fHistClustersPt->Fill(nPart.Pt());
578 cluster = fCaloClustersCont->GetNextAcceptCluster();
579 }
580 }
581
582 fHistLeadingJetPt->Fill(jetPt);
583 fHistLeadingJetPtCorr->Fill(jetPt-fJetsCont->GetRhoVal()*jetArea);
584
585 AliVTrack *track = static_cast<AliVTrack*>(fTracksCont->GetNextAcceptParticle(0));
586 while(track)
587 { // loop over all particles (including jet tracks)
588 trackPt = track->Pt();
589 fHistTracksPt->Fill(trackPt);
590 phi = track->Phi(); // get track phi
591
592 dPhi = phi-jetPhi; // get track phi - jet axis
593 if(dPhi < 0) dPhi += TMath::TwoPi();
594 if(dPhi > TMath::TwoPi()) dPhi -= TMath::TwoPi();
595
596 dPhiQA = phi-fEPV0; // get track phi - EP
597 if(dPhiQA < 0) dPhiQA += TMath::TwoPi();
598 if(dPhiQA > TMath::TwoPi()) dPhiQA -= TMath::TwoPi();
599 fDPhiEP->Fill(dPhiQA);
600
601 // fill jet-hadron correlation just to check if track labels make sense...
602 if(track->Pt()>1.)
603 {
604 Double_t dphiJet = dPhi;
605 if(dphiJet > 5*TMath::Pi()/3) dphiJet -= 2*TMath::Pi();
606 if(dphiJet < -TMath::Pi()/3) dphiJet += 2*TMath::Pi();
607 if(track->GetLabel() == 0)
608 fDPhiJet->Fill(dphiJet);
609 else
610 fDPhiJetPythia->Fill(dphiJet);
611 }
612
613 Double_t weight = 1.;
614 if(doPtWeight) weight = trackPt;
615
616 gx += weight*cos(2*dPhi);
617 gy += weight*sin(2*dPhi);
618
619 if(track->GetLabel() == 0)
620 { // sum for std EP method
621 tSin += weight*sin(2*phi); // bkgnd has label = 0
622 tCos += weight*cos(2*phi);
623 qx += weight*cos(2*dPhiQA);
624 qy += weight*sin(2*dPhiQA);
625 nTracksBkgnd += 1;
626 }
627 else
628 {
629 jSin += weight*sin(2*phi); // jets have label =/= 0
630 jCos += weight*cos(2*phi);
631 ax += weight*cos(2*dPhi);
632 nTracksJet += 1;
633 }
634
635 if(track->Pt() < 2.)
636 { // sum for std EP method w/ trackPt < 2 GeV
637 tSin2 += weight*sin(2*phi);
638 tCos2 += weight*cos(2*phi);
639 }
640
641
642 track = static_cast<AliVTrack*>(fTracksCont->GetNextAcceptParticle()); // increment to next track
643 } // close loop over particles
644
645 hNTracks->Fill(nTracksBkgnd,"Bkgnd Tracks",1);
646 hNTracks->Fill(nTracksJet,"Jet Tracks",1);
647 hNTracks->Fill(nTracksBkgnd+nTracksJet,"Total Tracks",1);
648 hNTracksCent->Fill(nTracksBkgnd+nTracksJet,fCent,1);
649
650 if(nTracksBkgnd == 0)
651 {
652 hEventData->Fill("no tracks",1);
653 hEventData->Fill("good event",-1);
654 return kFALSE;
655 }
656
657 Double_t v2weight = 1+2*fJetv2*cos(2*(jetPhi-fEPV0)); // set v2 weight for event
658
659 hGx->Fill(fCent, jetPt, v2weight*gx); // fill histograms for QA method
660 hGy2->Fill(fCent, jetPt, v2weight*gy*gy);
661 hGxGy2->Fill(fCent, jetPt, v2weight*gx*gy*gy);
662 hGy4->Fill(fCent, jetPt, v2weight*gy*gy*gy*gy);
663 hGx2->Fill(fCent, jetPt, v2weight*gx*gx);
664 hGx2Gy2->Fill(fCent, jetPt, v2weight*gx*gx*gy*gy);
665 hGxGy4->Fill(fCent, jetPt, v2weight*gx*gy*gy*gy*gy);
666 hGy6->Fill(fCent, jetPt, v2weight*gy*gy*gy*gy*gy*gy);
667 hGx2Gy4->Fill(fCent, jetPt, v2weight*gx*gx*gy*gy*gy*gy);
668 hGxGy6->Fill(fCent, jetPt, v2weight*gx*gy*gy*gy*gy*gy*gy);
669 hGy8->Fill(fCent, jetPt, v2weight*gy*gy*gy*gy*gy*gy*gy*gy);
670 hGy->Fill(fCent, jetPt, v2weight*gy);
671 hN->Fill(fCent, jetPt, v2weight);
672
673 hGxTracks->Fill(gx,nTracksBkgnd+nTracksJet);
674 hGyTracks->Fill(gy,nTracksBkgnd+nTracksJet);
675 hGy2Tracks->Fill(gy*gy,nTracksBkgnd+nTracksJet);
676 hGxGy2Tracks->Fill(gx*gy*gy,nTracksBkgnd+nTracksJet);
677 hGy4Tracks->Fill(gy*gy*gy*gy,nTracksBkgnd+nTracksJet);
678 hQx->Fill(qx/nTracksBkgnd,fCent);
679 hQy->Fill(qy/nTracksBkgnd,fCent);
680
681 if(isDijet == 0)
682 hAx->Fill(ax,jetPt);
683 if(isDijet == 1)
684 hAxDijet->Fill(ax,jetPt);
685
686 tEP = 0.5*atan2(tSin,tCos); // calculate EP w/o jets
687 tjEP = 0.5*atan2((tSin+jSin),(tCos+jCos)); // calculate EP w/ jets
688 tjEP2 = 0.5*atan2(tSin2,tCos2); // calculate EP w/ trackPt < 2 GeV
689
690 htEPRes->Fill(fCent, jetPt, v2weight*cos(2*(tEP-fEPV0)));
691 htjEPRes->Fill(fCent, jetPt, v2weight*cos(2*(tjEP-fEPV0)));
692 htj2EPRes->Fill(fCent, jetPt, v2weight*cos(2*(tjEP2-fEPV0)));
693
694
695 dPsi = jetPhi-tEP;
696 if(dPsi < 0) dPsi += TMath::TwoPi();
697 if(dPsi > TMath::TwoPi()) dPsi -= TMath::TwoPi();
698
699 htv2std->Fill(fCent, jetPt, v2weight*cos(2*dPsi)); // fill histogram with v2 data w/o jets
700 htdPsi->Fill(fCent,jetPt,dPsi); // fill histogram with jet axis - EP w/o jets
701
702 dPsi = jetPhi-tjEP;
703 if(dPsi < 0) dPsi += TMath::TwoPi();
704 if(dPsi > TMath::TwoPi()) dPsi -= TMath::TwoPi();
705
706 htjv2std->Fill(fCent, jetPt, v2weight*cos(2*dPsi)); // fill histogram with v2 data w/ jets
707 htjdPsi->Fill(fCent,jetPt,dPsi); // fill histogram with jet axis - EP w/ jets
708
709 dPsi = jetPhi-tjEP2;
710 if(dPsi < 0) dPsi += TMath::TwoPi();
711 if(dPsi > TMath::TwoPi()) dPsi -= TMath::TwoPi();
712
713 htj2v2std->Fill(fCent, jetPt, v2weight*cos(2*dPsi)); // fill histogram with v2 data w/ trackPt < 2 GeV
714 htj2dPsi->Fill(fCent,jetPt,dPsi); // fill histogram with jet axis - EP w/ trackPt < 2 GeV
715
716 dPsi = jetPhi-fEPV0;
717 if(dPsi < 0) dPsi += TMath::TwoPi();
718 if(dPsi > TMath::TwoPi()) dPsi -= TMath::TwoPi();
719
720 hV0jv2std->Fill(fCent, jetPt, v2weight*cos(2*dPsi)); // fill histogram with v2 data
721 hV0jdPsi->Fill(fCent,jetPt,dPsi); // fill histogram with jet axis - EPV0
722
723
724 return kTRUE; // If return kFALSE FillHistogram() will NOT be executed.
725}
726
727//________________________________________________________________________
728void AliAnalysisTaskEmcalJetv2QA::Terminate(Option_t *)
729{
730 // Called once at the end of the analysis.
731 if(centBins) delete [] centBins;
732 if(jetPtBins) delete [] jetPtBins;
733}
de2a3fac 734
735void AliAnalysisTaskEmcalJetv2QA::SetCentBins(Int_t n, Double_t* bins)
736{
737 if(centBins) delete [] centBins;
738 nCentBins = n;
739 nCentBins1 = n+1;
740 centBins = new Double_t[nCentBins+1];
741 for(Int_t i = 0; i < nCentBins+1; i++)
742 centBins[i]=bins[i];
743 cout << endl << "Setting " << nCentBins << " centrality bins: " << endl;
744 for(Int_t i = 0; i < nCentBins+1; i++)
745 cout << centBins[i] << " ";
746 cout << endl << endl;
747}
748
749void AliAnalysisTaskEmcalJetv2QA::SetJetPtBins(Int_t n, Double_t* bins)
750{
751 if(jetPtBins) delete [] jetPtBins;
752 nJetPtBins = n;
753 nJetPtBins1 = n+1;
754 jetPtBins = new Double_t[nJetPtBins+1];
755 for(Int_t i = 0; i < nJetPtBins+1; i++)
756 jetPtBins[i]=bins[i];
757 cout << endl << "Setting " << nJetPtBins << " jet pt bins: " << endl;
758 for(Int_t i = 0; i < nJetPtBins+1; i++)
759 cout << jetPtBins[i] << " ";
760 cout << endl << endl;
761}