]> git.uio.no Git - u/mrichter/AliRoot.git/blame - PWGJE/EMCALJetTasks/UserTasks/AliAnalysisTaskSAQA.cxx
Merge branch 'master' into TPCdev
[u/mrichter/AliRoot.git] / PWGJE / EMCALJetTasks / UserTasks / AliAnalysisTaskSAQA.cxx
CommitLineData
6e8d91c9 1// $Id$
c3ba2d3d 2//
297edd60 3// General QA task.
c3ba2d3d 4//
297edd60 5// Author: S.Aiola
c3ba2d3d 6
c3ba2d3d 7#include <TClonesArray.h>
8#include <TH1F.h>
9#include <TH2F.h>
3000c095 10#include <TH3F.h>
43032ce2 11#include <THnSparse.h>
c3ba2d3d 12#include <TList.h>
13#include <TLorentzVector.h>
14
15#include "AliAnalysisManager.h"
16#include "AliCentrality.h"
17#include "AliVCluster.h"
18#include "AliVParticle.h"
19#include "AliVTrack.h"
20#include "AliEmcalJet.h"
21#include "AliVEventHandler.h"
079b4732 22#include "AliAODEvent.h"
23#include "AliExternalTrackParam.h"
24#include "AliTrackerBase.h"
c3ba2d3d 25#include "AliLog.h"
59f16b27 26#include "AliEMCALGeometry.h"
27#include "AliEMCALGeoParams.h"
b16bb001 28#include "AliPicoTrack.h"
43032ce2 29#include "AliVVZERO.h"
30#include "AliESDUtils.h"
c3ba2d3d 31
32#include "AliAnalysisTaskSAQA.h"
33
34ClassImp(AliAnalysisTaskSAQA)
35
36//________________________________________________________________________
37AliAnalysisTaskSAQA::AliAnalysisTaskSAQA() :
9239b066 38 AliAnalysisTaskEmcalJet("AliAnalysisTaskSAQA", kTRUE),
5a025875 39 fCellEnergyCut(0.05),
1f9c287f 40 fParticleLevel(kFALSE),
5be3857d 41 fIsMC(kFALSE),
43032ce2 42 fCentMethod2(""),
43 fCentMethod3(""),
44 fDoV0QA(0),
3fe08cdb 45 fDoEPQA(0),
7cd832c7 46 fDoLeadingObjectPosition(0),
26516bb5 47 fMaxCellsInCluster(30),
43032ce2 48 fCent2(0),
49 fCent3(0),
50 fVZERO(0),
51 fV0ATotMult(0),
52 fV0CTotMult(0),
f30ab63b 53 fHistEventQA(0)
c3ba2d3d 54{
55 // Default constructor.
56
6fd5039f 57 for (Int_t i = 0; i < 4; i++) {
6f18d73a 58 for (Int_t j = 0; j < 4; j++) fHistTrPhiEtaPt[i][j] = 0;
f30ab63b 59
60 fHistTrNegativeLabels[i] = 0;
61 fHistTrZeroLabels[i] = 0;
43032ce2 62 fHistTrPhiEtaZeroLab[i] = 0;
63 fHistTrPtZeroLab[i] = 0;
56bd3193 64 fHistTrEmcPhiEta[i] = 0;
65 fHistTrEmcPt[i] = 0;
66 fHistTrPhiEtaNonProp[i] = 0;
43032ce2 67 fHistTrPtNonProp[i] = 0;
56bd3193 68 fHistDeltaEtaPt[i] = 0;
69 fHistDeltaPhiPt[i] = 0;
e2b76953 70 fHistDeltaPtvsPt[i] = 0;
f30ab63b 71
a487deae 72 fHistClusPhiEtaEnergy[i] = 0;
e2b76953 73 fHistClusDeltaPhiEPEnergy[i] = 0;
f30ab63b 74 fHistNCellsEnergy[i] = 0;
75 fHistFcrossEnergy[i] = 0;
76 fHistClusTimeEnergy[i] = 0;
8d3d1996 77 fHistClusMCEnergyFraction[i] = 0;
f30ab63b 78
79 fHistCellsAbsIdEnergy[i] = 0;
80
43032ce2 81 fHistJetsPhiEta[i] = 0;
98750b70 82 fHistJetsPtArea[i] = 0;
6fd5039f 83 }
a487deae 84
85 SetMakeGeneralHistograms(kTRUE);
c3ba2d3d 86}
87
88//________________________________________________________________________
89AliAnalysisTaskSAQA::AliAnalysisTaskSAQA(const char *name) :
9239b066 90 AliAnalysisTaskEmcalJet(name, kTRUE),
5a025875 91 fCellEnergyCut(0.05),
1f9c287f 92 fParticleLevel(kFALSE),
5be3857d 93 fIsMC(kFALSE),
43032ce2 94 fCentMethod2(""),
95 fCentMethod3(""),
96 fDoV0QA(0),
3fe08cdb 97 fDoEPQA(0),
7cd832c7 98 fDoLeadingObjectPosition(0),
26516bb5 99 fMaxCellsInCluster(30),
43032ce2 100 fCent2(0),
101 fCent3(0),
102 fVZERO(0),
103 fV0ATotMult(0),
104 fV0CTotMult(0),
f30ab63b 105 fHistEventQA(0)
c3ba2d3d 106{
107 // Standard constructor.
108
6fd5039f 109 for (Int_t i = 0; i < 4; i++) {
6f18d73a 110 for (Int_t j = 0; j < 4; j++) fHistTrPhiEtaPt[i][j] = 0;
f30ab63b 111
112 fHistTrNegativeLabels[i] = 0;
113 fHistTrZeroLabels[i] = 0;
43032ce2 114 fHistTrPhiEtaZeroLab[i] = 0;
115 fHistTrPtZeroLab[i] = 0;
56bd3193 116 fHistTrEmcPhiEta[i] = 0;
117 fHistTrEmcPt[i] = 0;
118 fHistTrPhiEtaNonProp[i] = 0;
43032ce2 119 fHistTrPtNonProp[i] = 0;
56bd3193 120 fHistDeltaEtaPt[i] = 0;
121 fHistDeltaPhiPt[i] = 0;
e2b76953 122 fHistDeltaPtvsPt[i] = 0;
f30ab63b 123
a487deae 124 fHistClusPhiEtaEnergy[i] = 0;
e2b76953 125 fHistClusDeltaPhiEPEnergy[i] = 0;
f30ab63b 126 fHistNCellsEnergy[i] = 0;
127 fHistFcrossEnergy[i] = 0;
128 fHistClusTimeEnergy[i] = 0;
8d3d1996 129 fHistClusMCEnergyFraction[i] = 0;
f30ab63b 130
131 fHistCellsAbsIdEnergy[i] = 0;
132
43032ce2 133 fHistJetsPhiEta[i] = 0;
98750b70 134 fHistJetsPtArea[i] = 0;
6fd5039f 135 }
a487deae 136
137 SetMakeGeneralHistograms(kTRUE);
c3ba2d3d 138}
139
140//________________________________________________________________________
141AliAnalysisTaskSAQA::~AliAnalysisTaskSAQA()
142{
143 // Destructor
144}
145
146//________________________________________________________________________
147void AliAnalysisTaskSAQA::UserCreateOutputObjects()
148{
149 // Create histograms
150
9239b066 151 AliAnalysisTaskEmcalJet::UserCreateOutputObjects();
c3ba2d3d 152
f30ab63b 153 TString histname;
154
6421eeb0 155 if (fParticleCollArray.GetEntriesFast()>0) {
4358e58a 156 if (!fParticleLevel && fIsMC) {
f30ab63b 157 for (Int_t i = 0; i < fNcentBins; i++) {
158 histname = Form("fHistTrNegativeLabels_%d",i);
159 fHistTrNegativeLabels[i] = new TH1F(histname,histname, 500, 0, 1);
160 fHistTrNegativeLabels[i]->GetXaxis()->SetTitle("% of negative labels");
161 fHistTrNegativeLabels[i]->GetYaxis()->SetTitle("counts");
162 fOutput->Add(fHistTrNegativeLabels[i]);
163
164 histname = Form("fHistTrZeroLabels_%d",i);
165 fHistTrZeroLabels[i] = new TH1F(histname,histname, 500, 0, 1);
166 fHistTrZeroLabels[i]->GetXaxis()->SetTitle("% of negative labels");
167 fHistTrZeroLabels[i]->GetYaxis()->SetTitle("counts");
168 fOutput->Add(fHistTrZeroLabels[i]);
169 }
4358e58a 170 }
171
393ec2ad 172 Int_t nlabels = 4;
173 if (fParticleLevel)
174 nlabels = 1;
175
176 for (Int_t i = 0; i < fNcentBins; i++) {
177 for (Int_t j = 0; j < nlabels; j++) {
6f18d73a 178 histname = Form("fHistTrPhiEtaPt_%d_%d",i,j);
7cd832c7 179 fHistTrPhiEtaPt[i][j] = new TH3F(histname,histname, 100, -1, 1, 101, 0, TMath::Pi() * 2.02, fNbins, fMinBinPt, fMaxBinPt);
6f18d73a 180 fHistTrPhiEtaPt[i][j]->GetXaxis()->SetTitle("#eta");
181 fHistTrPhiEtaPt[i][j]->GetYaxis()->SetTitle("#phi");
182 fHistTrPhiEtaPt[i][j]->GetZaxis()->SetTitle("p_{T} (GeV/c)");
183 fOutput->Add(fHistTrPhiEtaPt[i][j]);
184 }
c3ba2d3d 185
56bd3193 186 if (!fParticleLevel) {
187 if (fIsMC) {
43032ce2 188 histname = Form("fHistTrPhiEtaZeroLab_%d",i);
7cd832c7 189 fHistTrPhiEtaZeroLab[i] = new TH2F(histname,histname, 100, -1, 1, 101, 0, TMath::Pi() * 2.02);
43032ce2 190 fHistTrPhiEtaZeroLab[i]->GetXaxis()->SetTitle("#eta");
191 fHistTrPhiEtaZeroLab[i]->GetYaxis()->SetTitle("#phi");
192 fHistTrPhiEtaZeroLab[i]->GetZaxis()->SetTitle("counts");
193 fOutput->Add(fHistTrPhiEtaZeroLab[i]);
194
195 histname = Form("fHistTrPtZeroLab_%d",i);
196 fHistTrPtZeroLab[i] = new TH1F(histname,histname, fNbins, fMinBinPt, fMaxBinPt);
197 fHistTrPtZeroLab[i]->GetZaxis()->SetTitle("p_{T} (GeV/c)");
198 fHistTrPtZeroLab[i]->GetYaxis()->SetTitle("counts");
199 fOutput->Add(fHistTrPtZeroLab[i]);
56bd3193 200 }
201
202 histname = Form("fHistTrEmcPhiEta_%d",i);
7cd832c7 203 fHistTrEmcPhiEta[i] = new TH2F(histname,histname, 100, -1, 1, 101, 0, TMath::Pi() * 2.02);
56bd3193 204 fHistTrEmcPhiEta[i]->GetXaxis()->SetTitle("#eta");
205 fHistTrEmcPhiEta[i]->GetYaxis()->SetTitle("#phi");
206 fOutput->Add(fHistTrEmcPhiEta[i]);
207
208 histname = Form("fHistTrEmcPt_%d",i);
209 fHistTrEmcPt[i] = new TH1F(histname,histname, fNbins, fMinBinPt, fMaxBinPt);
210 fHistTrEmcPt[i]->GetXaxis()->SetTitle("p_{T} (GeV/c)");
211 fHistTrEmcPt[i]->GetYaxis()->SetTitle("counts");
212 fOutput->Add(fHistTrEmcPt[i]);
213
214 histname = Form("fHistTrPhiEtaNonProp_%d",i);
7cd832c7 215 fHistTrPhiEtaNonProp[i] = new TH2F(histname,histname, 100, -1, 1, 101, 0, TMath::Pi() * 2.02);
56bd3193 216 fHistTrPhiEtaNonProp[i]->GetXaxis()->SetTitle("#eta");
217 fHistTrPhiEtaNonProp[i]->GetYaxis()->SetTitle("#phi");
218 fOutput->Add(fHistTrPhiEtaNonProp[i]);
43032ce2 219
220 histname = Form("fHistTrPtNonProp_%d",i);
221 fHistTrPtNonProp[i] = new TH1F(histname,histname, fNbins, fMinBinPt, fMaxBinPt);
222 fHistTrPtNonProp[i]->GetXaxis()->SetTitle("p_{T} (GeV/c)");
223 fHistTrPtNonProp[i]->GetYaxis()->SetTitle("counts");
224 fOutput->Add(fHistTrPtNonProp[i]);
56bd3193 225
226 histname = Form("fHistDeltaEtaPt_%d",i);
227 fHistDeltaEtaPt[i] = new TH2F(histname,histname, fNbins, fMinBinPt, fMaxBinPt, 50, -0.5, 0.5);
228 fHistDeltaEtaPt[i]->GetXaxis()->SetTitle("p_{T} (GeV/c)");
229 fHistDeltaEtaPt[i]->GetYaxis()->SetTitle("#delta#eta");
230 fOutput->Add(fHistDeltaEtaPt[i]);
231
232 histname = Form("fHistDeltaPhiPt_%d",i);
233 fHistDeltaPhiPt[i] = new TH2F(histname,histname, fNbins, fMinBinPt, fMaxBinPt, 200, -2, 2);
234 fHistDeltaPhiPt[i]->GetXaxis()->SetTitle("p_{T} (GeV/c)");
235 fHistDeltaPhiPt[i]->GetYaxis()->SetTitle("#delta#phi");
236 fOutput->Add(fHistDeltaPhiPt[i]);
237
e2b76953 238 histname = Form("fHistDeltaPtvsPt_%d",i);
239 fHistDeltaPtvsPt[i] = new TH2F(histname,histname, fNbins, fMinBinPt, fMaxBinPt, fNbins, -fMaxBinPt/2, fMaxBinPt/2);
240 fHistDeltaPtvsPt[i]->GetXaxis()->SetTitle("p_{T} (GeV/c)");
241 fHistDeltaPtvsPt[i]->GetYaxis()->SetTitle("#deltap_{T} (GeV/c)");
242 fHistDeltaPtvsPt[i]->GetZaxis()->SetTitle("counts");
243 fOutput->Add(fHistDeltaPtvsPt[i]);
56bd3193 244 }
2103dc6a 245 }
a487deae 246 }
7c1d624c 247
6421eeb0 248 if (fClusterCollArray.GetEntriesFast()>0) {
393ec2ad 249 for (Int_t i = 0; i < fNcentBins; i++) {
a487deae 250 histname = "fHistClusPhiEtaEnergy_";
251 histname += i;
7cd832c7 252 fHistClusPhiEtaEnergy[i] = new TH3F(histname, histname, 100, -1, 1, 101, 0, TMath::Pi() * 2.02, fNbins, fMinBinPt, fMaxBinPt);
a487deae 253 fHistClusPhiEtaEnergy[i]->GetXaxis()->SetTitle("#eta");
254 fHistClusPhiEtaEnergy[i]->GetYaxis()->SetTitle("#phi");
59f16b27 255 fHistClusPhiEtaEnergy[i]->GetZaxis()->SetTitle("E_{cluster} (GeV)");
a487deae 256 fOutput->Add(fHistClusPhiEtaEnergy[i]);
8d3d1996 257
e2b76953 258 histname = "fHistClusDeltaPhiEPEnergy_";
259 histname += i;
260 fHistClusDeltaPhiEPEnergy[i] = new TH2F(histname, histname, fNbins, fMinBinPt, fMaxBinPt, 100, 0, TMath::Pi());
261 fHistClusDeltaPhiEPEnergy[i]->GetXaxis()->SetTitle("E_{cluster} (GeV)");
262 fHistClusDeltaPhiEPEnergy[i]->GetYaxis()->SetTitle("#phi_{cluster} - #psi_{RP}");
263 fOutput->Add(fHistClusDeltaPhiEPEnergy[i]);
264
8d3d1996 265 if (fIsEmbedded) {
266 histname = "fHistClusMCEnergyFraction_";
267 histname += i;
268 fHistClusMCEnergyFraction[i] = new TH1F(histname, histname, fNbins, 0, 1.2);
269 fHistClusMCEnergyFraction[i]->GetXaxis()->SetTitle("MC fraction");
270 fHistClusMCEnergyFraction[i]->GetYaxis()->SetTitle("counts");
271 fOutput->Add(fHistClusMCEnergyFraction[i]);
272 }
a487deae 273
f30ab63b 274 histname = "fHistClusTimeEnergy_";
275 histname += i;
276 fHistClusTimeEnergy[i] = new TH2F(histname,histname, fNbins, fMinBinPt, fMaxBinPt, fNbins, -1e-6, 1e-6);
277 fHistClusTimeEnergy[i]->GetXaxis()->SetTitle("E_{cluster} (GeV)");
278 fHistClusTimeEnergy[i]->GetYaxis()->SetTitle("Time");
279 fOutput->Add(fHistClusTimeEnergy[i]);
280
281 Int_t nbins = fMaxCellsInCluster;
282 while (nbins > fNbins) nbins /= 2;
283 histname = "fHistNCellsEnergy_";
284 histname += i;
285 fHistNCellsEnergy[i] = new TH2F(histname,histname, fNbins, fMinBinPt, fMaxBinPt, nbins, -0.5, fMaxCellsInCluster-0.5);
286 fHistNCellsEnergy[i]->GetXaxis()->SetTitle("E_{cluster} (GeV)");
287 fHistNCellsEnergy[i]->GetYaxis()->SetTitle("N_{cells}");
288 fOutput->Add(fHistNCellsEnergy[i]);
289
290 histname = "fHistFcrossEnergy_";
291 histname += i;
292 fHistFcrossEnergy[i] = new TH2F(histname,histname, fNbins, fMinBinPt, fMaxBinPt, 200, -3.5, 1.5);
293 fHistFcrossEnergy[i]->GetXaxis()->SetTitle("E_{cluster} (GeV)");
294 fHistFcrossEnergy[i]->GetYaxis()->SetTitle("F_{cross}");
295 fOutput->Add(fHistFcrossEnergy[i]);
296
297 histname = "fHistCellsAbsIdEnergy_";
298 histname += i;
299 fHistCellsAbsIdEnergy[i] = new TH2F(histname,histname, 11600,0,11599,(Int_t)(fNbins / 2), fMinBinPt, fMaxBinPt / 2);
300 fHistCellsAbsIdEnergy[i]->GetXaxis()->SetTitle("cell abs. Id");
301 fHistCellsAbsIdEnergy[i]->GetYaxis()->SetTitle("E_{cluster} (GeV)");
302 fHistCellsAbsIdEnergy[i]->GetZaxis()->SetTitle("counts");
303 fOutput->Add(fHistCellsAbsIdEnergy[i]);
304 }
c3ba2d3d 305 }
a487deae 306
6421eeb0 307 if (fJetCollArray.GetEntriesFast()>0) {
393ec2ad 308 for (Int_t i = 0; i < fNcentBins; i++) {
43032ce2 309 histname = "fHistJetsPhiEta_";
e44e8726 310 histname += i;
7cd832c7 311 fHistJetsPhiEta[i] = new TH2F(histname.Data(), histname.Data(), 100, -1, 1, 101, 0, TMath::Pi() * 2.02);
43032ce2 312 fHistJetsPhiEta[i]->GetXaxis()->SetTitle("#eta");
313 fHistJetsPhiEta[i]->GetYaxis()->SetTitle("#phi");
314 fHistJetsPhiEta[i]->GetZaxis()->SetTitle("counts");
315 fOutput->Add(fHistJetsPhiEta[i]);
e44e8726 316
317 histname = "fHistJetsPtArea_";
318 histname += i;
6421eeb0 319 fHistJetsPtArea[i] = new TH2F(histname.Data(), histname.Data(), (Int_t)(fNbins * 2.5), fMinBinPt, fMaxBinPt * 2.5, 50, 0, 1.5);
a487deae 320 fHistJetsPtArea[i]->GetXaxis()->SetTitle("p_{T} (GeV/c)");
e44e8726 321 fHistJetsPtArea[i]->GetYaxis()->SetTitle("area");
322 fOutput->Add(fHistJetsPtArea[i]);
323 }
6fd5039f 324 }
43032ce2 325
326 Int_t dim = 0;
26516bb5 327 TString title[20];
328 Int_t nbins[20] = {0};
329 Double_t min[20] = {0};
330 Double_t max[20] = {0};
43032ce2 331
9239b066 332 if (fForceBeamType != AliAnalysisTaskEmcal::kpp) {
43032ce2 333 title[dim] = "Centrality %";
334 nbins[dim] = 101;
335 min[dim] = 0;
336 max[dim] = 101;
337 dim++;
338
339 if (!fCentMethod2.IsNull()) {
340 title[dim] = Form("Centrality %s %%", fCentMethod2.Data());
341 nbins[dim] = 101;
342 min[dim] = 0;
343 max[dim] = 101;
344 dim++;
345 }
346
347 if (!fCentMethod3.IsNull()) {
348 title[dim] = Form("Centrality %s %%", fCentMethod3.Data());
349 nbins[dim] = 101;
350 min[dim] = 0;
351 max[dim] = 101;
352 dim++;
353 }
354
26516bb5 355 if (fDoV0QA==1) {
43032ce2 356 title[dim] = "V0A total multiplicity";
357 nbins[dim] = 200;
358 min[dim] = 0;
359 max[dim] = 20000;
360 dim++;
361
362 title[dim] = "V0C total multiplicity";
363 nbins[dim] = 200;
364 min[dim] = 0;
365 max[dim] = 20000;
366 dim++;
367 }
26516bb5 368 else if (fDoV0QA==2) {
369 title[dim] = "V0A+V0C total multiplicity";
370 nbins[dim] = 300;
371 min[dim] = 0;
372 max[dim] = 30000;
373 dim++;
374 }
43032ce2 375
376 if (!fRhoName.IsNull()) {
377 title[dim] = "#rho (GeV/c)";
378 nbins[dim] = fNbins*4;
379 min[dim] = 0;
380 max[dim] = 400;
381 dim++;
382 }
3fe08cdb 383
384 if (fDoEPQA) {
385 title[dim] = "#psi_{RP}";
386 nbins[dim] = 200;
387 min[dim] = -TMath::Pi();
388 max[dim] = TMath::Pi();
389 dim++;
390 }
43032ce2 391 }
392
6421eeb0 393 if (fParticleCollArray.GetEntriesFast()>0) {
43032ce2 394 title[dim] = "No. of tracks";
26516bb5 395 nbins[dim] = 3000;
43032ce2 396 min[dim] = -0.5;
397 max[dim] = 6000-0.5;
398 dim++;
26516bb5 399
400 title[dim] = "p_{T,track}^{leading} (GeV/c)";
401 nbins[dim] = fNbins;
402 min[dim] = fMinBinPt;
403 max[dim] = fMaxBinPt;
404 dim++;
7cd832c7 405
406 if (fDoLeadingObjectPosition) {
407 title[dim] = "#eta_{track}^{leading}";
408 nbins[dim] = 100;
409 min[dim] = -1;
410 max[dim] = 1;
411 dim++;
412
413 title[dim] = "#phi_{track}^{leading}";
414 nbins[dim] = 101;
415 min[dim] = 0;
416 max[dim] = TMath::Pi() * 2.02;
417 dim++;
418 }
43032ce2 419 }
420
6421eeb0 421 if (fClusterCollArray.GetEntriesFast()>0) {
43032ce2 422 title[dim] = "No. of clusters";
26516bb5 423 nbins[dim] = 2000;
43032ce2 424 min[dim] = 0;
425 max[dim] = 4000-0.5;
426 dim++;
26516bb5 427
428 title[dim] = "E_{cluster}^{leading} (GeV)";
429 nbins[dim] = fNbins;
430 min[dim] = fMinBinPt;
431 max[dim] = fMaxBinPt;
432 dim++;
7cd832c7 433
434 if (fDoLeadingObjectPosition) {
435 title[dim] = "#eta_{cluster}^{leading}";
436 nbins[dim] = 100;
437 min[dim] = -1;
438 max[dim] = 1;
439 dim++;
440
441 title[dim] = "#phi_{cluster}^{leading}";
442 nbins[dim] = 101;
443 min[dim] = 0;
444 max[dim] = TMath::Pi() * 2.02;
445 dim++;
446 }
43032ce2 447 }
448
449 if (!fCaloCellsName.IsNull()) {
450 title[dim] = "No. of cells";
26516bb5 451 nbins[dim] = 3000;
43032ce2 452 min[dim] = 0;
453 max[dim] = 6000-0.5;
454 dim++;
455 }
456
6421eeb0 457 if (fJetCollArray.GetEntriesFast()>0) {
43032ce2 458 title[dim] = "No. of jets";
459 nbins[dim] = 200;
460 min[dim] = 0;
461 max[dim] = 200-0.5;
462 dim++;
26516bb5 463
464 title[dim] = "p_{T,jet}^{leading} (GeV/c)";
465 nbins[dim] = fNbins;
466 min[dim] = fMinBinPt;
467 max[dim] = fMaxBinPt;
468 dim++;
7cd832c7 469
470 if (fDoLeadingObjectPosition) {
471 title[dim] = "#eta_{jet}^{leading}";
472 nbins[dim] = 100;
473 min[dim] = -1;
474 max[dim] = 1;
475 dim++;
476
477 title[dim] = "#phi_{jet}^{leading}";
478 nbins[dim] = 101;
479 min[dim] = 0;
480 max[dim] = TMath::Pi() * 2.02;
481 dim++;
482 }
43032ce2 483 }
484
485 fHistEventQA = new THnSparseF("fHistEventQA","fHistEventQA",dim,nbins,min,max);
486 for (Int_t i = 0; i < dim; i++)
487 fHistEventQA->GetAxis(i)->SetTitle(title[i]);
488 fOutput->Add(fHistEventQA);
6fd5039f 489
c3ba2d3d 490 PostData(1, fOutput); // Post data for ALL output slots >0 here, to get at least an empty histogram
491}
492
43032ce2 493//________________________________________________________________________
494void AliAnalysisTaskSAQA::ExecOnce()
495{
9239b066 496 AliAnalysisTaskEmcalJet::ExecOnce();
43032ce2 497
498 if (fDoV0QA) {
499 fVZERO = InputEvent()->GetVZEROData();
500 if (!fVZERO) {
501 AliError("AliVVZERO not available");
502 }
503 }
504}
505
c3ba2d3d 506//________________________________________________________________________
6fd5039f 507Bool_t AliAnalysisTaskSAQA::RetrieveEventObjects()
c3ba2d3d 508{
c3604199 509 // Retrieve event objects.
510
9239b066 511 if (!AliAnalysisTaskEmcalJet::RetrieveEventObjects())
6fd5039f 512 return kFALSE;
c3ba2d3d 513
43032ce2 514 if (!fCentMethod2.IsNull() || !fCentMethod3.IsNull()) {
515 if (fBeamType == kAA || fBeamType == kpA ) {
516 AliCentrality *aliCent = InputEvent()->GetCentrality();
517 if (aliCent) {
518 if (!fCentMethod2.IsNull())
519 fCent2 = aliCent->GetCentralityPercentile(fCentMethod2);
520 if (!fCentMethod3.IsNull())
521 fCent3 = aliCent->GetCentralityPercentile(fCentMethod3);
522 }
523 }
524 }
525
526 if (fVZERO) {
527 fV0ATotMult = AliESDUtils::GetCorrV0A(fVZERO->GetMTotV0A(),fVertex[2]);
528 fV0CTotMult = AliESDUtils::GetCorrV0C(fVZERO->GetMTotV0C(),fVertex[2]);
529 }
7c1d624c 530
6fd5039f 531 return kTRUE;
c3ba2d3d 532}
533
7c1d624c 534
c3ba2d3d 535//________________________________________________________________________
6fd5039f 536Bool_t AliAnalysisTaskSAQA::FillHistograms()
c3ba2d3d 537{
c3604199 538 // Fill histograms.
539
a487deae 540 Float_t trackSum = 0;
43032ce2 541 Float_t clusSum = 0;
542 Float_t cellSum = 0;
090a0c3e 543
43032ce2 544 Int_t ntracks = 0;
545 Int_t nclusters = 0;
546 Int_t ncells = 0;
547 Int_t njets = 0;
26516bb5 548
7cd832c7 549 Float_t leadingClusE = 0;
550 Float_t leadingClusEta = 0;
551 Float_t leadingClusPhi = 0;
552
553 Float_t leadingTrackPt = 0;
554 Float_t leadingTrackEta = 0;
555 Float_t leadingTrackPhi = 0;
556
557 Float_t leadingJetPt = 0;
558 Float_t leadingJetEta = 0;
559 Float_t leadingJetPhi = 0;
43032ce2 560
a487deae 561 if (fTracks) {
7cd832c7 562 AliVParticle *leadingTrack = 0;
563
26516bb5 564 ntracks = DoTrackLoop(trackSum, leadingTrack);
43032ce2 565 AliDebug(2,Form("%d tracks found in the event", ntracks));
7cd832c7 566
567 if (leadingTrack) {
568 leadingTrackPt = leadingTrack->Pt();
569 leadingTrackEta = leadingTrack->Eta();
570 leadingTrackPhi = leadingTrack->Phi();
571 }
a487deae 572 }
be7b6e63 573
a487deae 574 if (fCaloClusters) {
7cd832c7 575 AliVCluster *leadingClus = 0;
576
26516bb5 577 nclusters = DoClusterLoop(clusSum, leadingClus);
43032ce2 578 AliDebug(2,Form("%d clusters found in the event", nclusters));
c3ba2d3d 579
7cd832c7 580 if (leadingClus) {
581 TLorentzVector leadingClusVect;
582 leadingClus->GetMomentum(leadingClusVect, fVertex);
583 leadingClusE = leadingClus->E();
584 leadingClusEta = leadingClusVect.Eta();
585 leadingClusPhi = leadingClusVect.Phi();
586 }
43032ce2 587 }
588
589 if (fCaloCells) {
5a025875 590 ncells = DoCellLoop(cellSum);
43032ce2 591 AliDebug(2,Form("%d cells found in the event", ncells));
a487deae 592 }
c3ba2d3d 593
a487deae 594 if (fJets) {
7cd832c7 595 AliEmcalJet *leadingJet = 0;
596
26516bb5 597 njets = DoJetLoop(leadingJet);
43032ce2 598 AliDebug(2,Form("%d jets found in the event", njets));
7cd832c7 599
600 if (leadingJet) {
601 leadingJetPt = leadingJet->Pt();
602 leadingJetEta = leadingJet->Eta();
603 leadingJetPhi = leadingJet->Phi();
604 }
6fd5039f 605 }
c3ba2d3d 606
7cd832c7 607 FillEventQAHisto(fCent, fCent2, fCent3, fV0ATotMult, fV0CTotMult, fEPV0, fRhoVal,
608 ntracks, nclusters, ncells, njets,
609 leadingTrackPt, leadingTrackEta, leadingTrackPhi,
610 leadingClusE, leadingClusEta, leadingClusPhi,
611 leadingJetPt, leadingJetEta, leadingJetPhi);
43032ce2 612
6fd5039f 613 return kTRUE;
c3ba2d3d 614}
615
43032ce2 616//________________________________________________________________________
26516bb5 617void AliAnalysisTaskSAQA::FillEventQAHisto(Float_t cent, Float_t cent2, Float_t cent3, Float_t v0a, Float_t v0c,
618 Float_t ep, Float_t rho, Int_t ntracks, Int_t nclusters, Int_t ncells, Int_t njets,
7cd832c7 619 Float_t maxTrackPt, Float_t maxTrackEta, Float_t maxTrackPhi,
620 Float_t maxClusterE, Float_t maxClusterEta, Float_t maxClusterPhi,
621 Float_t maxJetPt, Float_t maxJetEta, Float_t maxJetPhi)
43032ce2 622{
26516bb5 623 Double_t contents[20]={0};
43032ce2 624
625 for (Int_t i = 0; i < fHistEventQA->GetNdimensions(); i++) {
626 TString title(fHistEventQA->GetAxis(i)->GetTitle());
627 if (title=="Centrality %")
628 contents[i] = cent;
629 else if (title==Form("Centrality %s %%", fCentMethod2.Data()))
630 contents[i] = cent2;
631 else if (title==Form("Centrality %s %%", fCentMethod3.Data()))
632 contents[i] = cent3;
633 else if (title=="V0A total multiplicity")
634 contents[i] = v0a;
635 else if (title=="V0C total multiplicity")
636 contents[i] = v0c;
26516bb5 637 else if (title=="V0A+V0C total multiplicity")
638 contents[i] = v0a+v0c;
3fe08cdb 639 else if (title=="#psi_{RP}")
640 contents[i] = ep;
43032ce2 641 else if (title=="#rho (GeV/c)")
642 contents[i] = rho;
643 else if (title=="No. of tracks")
644 contents[i] = ntracks;
645 else if (title=="No. of clusters")
646 contents[i] = nclusters;
647 else if (title=="No. of cells")
648 contents[i] = ncells;
649 else if (title=="No. of jets")
650 contents[i] = njets;
26516bb5 651 else if (title=="p_{T,track}^{leading} (GeV/c)")
7cd832c7 652 contents[i] = maxTrackPt;
653 else if (title=="#eta_{track}^{leading}")
654 contents[i] = maxTrackEta;
655 else if (title=="#phi_{track}^{leading}")
656 contents[i] = maxTrackPhi;
26516bb5 657 else if (title=="E_{cluster}^{leading} (GeV)")
7cd832c7 658 contents[i] = maxClusterE;
659 else if (title=="#eta_{cluster}^{leading}")
660 contents[i] = maxClusterEta;
661 else if (title=="#phi_{cluster}^{leading}")
662 contents[i] = maxClusterPhi;
26516bb5 663 else if (title=="p_{T,jet}^{leading} (GeV/c)")
7cd832c7 664 contents[i] = maxJetPt;
665 else if (title=="#eta_{jet}^{leading}")
666 contents[i] = maxJetEta;
667 else if (title=="#phi_{jet}^{leading}")
668 contents[i] = maxJetPhi;
43032ce2 669 else
670 AliWarning(Form("Unable to fill dimension %s!",title.Data()));
671 }
672
673 fHistEventQA->Fill(contents);
674}
675
c3ba2d3d 676//________________________________________________________________________
5a025875 677Int_t AliAnalysisTaskSAQA::DoCellLoop(Float_t &sum)
c3ba2d3d 678{
c3604199 679 // Do cell loop.
680
c3ba2d3d 681 AliVCaloCells *cells = InputEvent()->GetEMCALCells();
682
683 if (!cells)
be7b6e63 684 return 0;
c3ba2d3d 685
c3604199 686 const Int_t ncells = cells->GetNumberOfCells();
5a025875 687 Int_t nAccCells = 0;
c3ba2d3d 688
689 for (Int_t pos = 0; pos < ncells; pos++) {
f483218e 690 Float_t amp = cells->GetAmplitude(pos);
691 Int_t absId = cells->GetCellNumber(pos);
5a025875 692
c3ba2d3d 693 if (amp < fCellEnergyCut)
694 continue;
be7b6e63 695
5a025875 696 fHistCellsAbsIdEnergy[fCentBin]->Fill(absId,amp);
697 nAccCells++;
698 sum += amp;
699 }
c3ba2d3d 700
5a025875 701 return nAccCells;
e2b76953 702}
703
59f16b27 704//________________________________________________________________________
705Double_t AliAnalysisTaskSAQA::GetFcross(AliVCluster *cluster, AliVCaloCells *cells)
706{
707 Int_t AbsIdseed = -1;
708 Double_t Eseed = 0;
709 for (Int_t i = 0; i < cluster->GetNCells(); i++) {
710 if (cells->GetCellAmplitude(cluster->GetCellAbsId(i)) > AbsIdseed) {
711 Eseed = cells->GetCellAmplitude(cluster->GetCellAbsId(i));
712 AbsIdseed = cluster->GetCellAbsId(i);
713 }
714 }
715
716 if (Eseed < 1e-9)
717 return 100;
718
719 Int_t imod = -1, iphi =-1, ieta=-1,iTower = -1, iIphi = -1, iIeta = -1;
720 fGeom->GetCellIndex(AbsIdseed,imod,iTower,iIphi,iIeta);
721 fGeom->GetCellPhiEtaIndexInSModule(imod,iTower,iIphi,iIeta,iphi,ieta);
722
723 //Get close cells index and energy, not in corners
724
725 Int_t absID1 = -1;
726 Int_t absID2 = -1;
727
728 if (iphi < AliEMCALGeoParams::fgkEMCALRows-1) absID1 = fGeom->GetAbsCellIdFromCellIndexes(imod, iphi+1, ieta);
729 if (iphi > 0) absID2 = fGeom->GetAbsCellIdFromCellIndexes(imod, iphi-1, ieta);
730
731 // In case of cell in eta = 0 border, depending on SM shift the cross cell index
732
733 Int_t absID3 = -1;
734 Int_t absID4 = -1;
735
736 if (ieta == AliEMCALGeoParams::fgkEMCALCols-1 && !(imod%2)) {
737 absID3 = fGeom->GetAbsCellIdFromCellIndexes(imod+1, iphi, 0);
738 absID4 = fGeom->GetAbsCellIdFromCellIndexes(imod, iphi, ieta-1);
739 }
740 else if (ieta == 0 && imod%2) {
741 absID3 = fGeom->GetAbsCellIdFromCellIndexes(imod, iphi, ieta+1);
742 absID4 = fGeom->GetAbsCellIdFromCellIndexes(imod-1, iphi, AliEMCALGeoParams::fgkEMCALCols-1);
743 }
744 else {
745 if (ieta < AliEMCALGeoParams::fgkEMCALCols-1)
746 absID3 = fGeom->GetAbsCellIdFromCellIndexes(imod, iphi, ieta+1);
747 if (ieta > 0)
748 absID4 = fGeom->GetAbsCellIdFromCellIndexes(imod, iphi, ieta-1);
749 }
750
751 Double_t ecell1 = cells->GetCellAmplitude(absID1);
752 Double_t ecell2 = cells->GetCellAmplitude(absID2);
753 Double_t ecell3 = cells->GetCellAmplitude(absID3);
754 Double_t ecell4 = cells->GetCellAmplitude(absID4);
755
756 Double_t Ecross = ecell1 + ecell2 + ecell3 + ecell4;
757
758 Double_t Fcross = 1 - Ecross/Eseed;
759
760 return Fcross;
761}
762
c3ba2d3d 763//________________________________________________________________________
7cd832c7 764Int_t AliAnalysisTaskSAQA::DoClusterLoop(Float_t &sum, AliVCluster* &leading)
c3ba2d3d 765{
c3604199 766 // Do cluster loop.
767
c3ba2d3d 768 if (!fCaloClusters)
769 return 0;
770
43032ce2 771 Int_t nAccClusters = 0;
772
59f16b27 773 AliVCaloCells *cells = InputEvent()->GetEMCALCells();
774
43032ce2 775 sum = 0;
26516bb5 776 leading = 0;
c3ba2d3d 777
778 // Cluster loop
43032ce2 779 const Int_t nclusters = fCaloClusters->GetEntriesFast();
c3ba2d3d 780
781 for (Int_t iClusters = 0; iClusters < nclusters; iClusters++) {
e44e8726 782 AliVCluster* cluster = static_cast<AliVCluster*>(fCaloClusters->At(iClusters));
c3ba2d3d 783 if (!cluster) {
784 AliError(Form("Could not receive cluster %d", iClusters));
785 continue;
786 }
c3ba2d3d 787
f660c2d6 788 if (!AcceptCluster(cluster))
e44e8726 789 continue;
c3ba2d3d 790
c3ba2d3d 791 sum += cluster->E();
792
7cd832c7 793 if (!leading || leading->E() < cluster->E()) leading = cluster;
26516bb5 794
7cf4626b 795 TLorentzVector nPart;
796 cluster->GetMomentum(nPart, fVertex);
c3ba2d3d 797
a487deae 798 fHistClusPhiEtaEnergy[fCentBin]->Fill(nPart.Eta(), nPart.Phi(), cluster->E());
e2b76953 799
800 Double_t ep = nPart.Phi() - fEPV0;
801 while (ep < 0) ep += TMath::Pi();
802 while (ep >= TMath::Pi()) ep -= TMath::Pi();
803 fHistClusDeltaPhiEPEnergy[fCentBin]->Fill(cluster->E(), ep);
804
f30ab63b 805 fHistNCellsEnergy[fCentBin]->Fill(cluster->E(), cluster->GetNCells());
090a0c3e 806
f30ab63b 807 fHistClusTimeEnergy[fCentBin]->Fill(cluster->E(), cluster->GetTOF());
7c1d624c 808
f30ab63b 809 if (cells) fHistFcrossEnergy[fCentBin]->Fill(cluster->E(), GetFcross(cluster, cells));
59f16b27 810
8d3d1996 811 if (fHistClusMCEnergyFraction[fCentBin])
812 fHistClusMCEnergyFraction[fCentBin]->Fill(cluster->GetMCEnergyFraction());
1b47ed10 813
43032ce2 814 nAccClusters++;
c3604199 815 }
c3ba2d3d 816
43032ce2 817 return nAccClusters;
c3ba2d3d 818}
819
820//________________________________________________________________________
7cd832c7 821Int_t AliAnalysisTaskSAQA::DoTrackLoop(Float_t &sum, AliVParticle* &leading)
c3ba2d3d 822{
c3604199 823 // Do track loop.
c3ba2d3d 824 if (!fTracks)
825 return 0;
826
43032ce2 827 Int_t nAccTracks = 0;
828
829 sum = 0;
26516bb5 830 leading = 0;
c3ba2d3d 831
43032ce2 832 const Int_t ntracks = fTracks->GetEntriesFast();
4358e58a 833 Int_t neg = 0;
834 Int_t zero = 0;
c3ba2d3d 835
c3604199 836 for (Int_t i = 0; i < ntracks; i++) {
e44e8726 837 AliVParticle* track = static_cast<AliVParticle*>(fTracks->At(i)); // pointer to reconstructed to track
1f6fff78 838
c3604199 839 if (!track) {
c3ba2d3d 840 AliError(Form("Could not retrieve track %d",i));
841 continue;
842 }
c3ba2d3d 843
5be3857d 844 if (!AcceptTrack(track))
1f6fff78 845 continue;
7c1d624c 846
43032ce2 847 nAccTracks++;
c3ba2d3d 848
849 sum += track->P();
850
7cd832c7 851 if (!leading || leading->Pt() < track->Pt()) leading = track;
26516bb5 852
393ec2ad 853 if (fParticleLevel) {
854 fHistTrPhiEtaPt[fCentBin][0]->Fill(track->Eta(), track->Phi(), track->Pt());
855 }
856 else {
857 fHistTrPhiEtaPt[fCentBin][3]->Fill(track->Eta(), track->Phi(), track->Pt());
4358e58a 858 if (track->GetLabel() == 0) {
859 zero++;
43032ce2 860 if (fHistTrPhiEtaZeroLab[fCentBin]) {
861 fHistTrPhiEtaZeroLab[fCentBin]->Fill(track->Eta(), track->Phi());
862 fHistTrPtZeroLab[fCentBin]->Fill(track->Pt());
863 }
4358e58a 864 }
865
866 if (track->GetLabel() < 0)
867 neg++;
787a3c4f 868
b16bb001 869 Int_t type = 0;
393ec2ad 870
b16bb001 871 AliPicoTrack* ptrack = dynamic_cast<AliPicoTrack*>(track);
872 if (ptrack)
873 type = ptrack->GetTrackType();
874
875 if (type >= 0 && type < 3)
876 fHistTrPhiEtaPt[fCentBin][type]->Fill(track->Eta(), track->Phi(), track->Pt());
393ec2ad 877 else
5be3857d 878 AliDebug(2,Form("%s: track type %d not recognized!", GetName(), type));
393ec2ad 879 }
aa4d701c 880
5be3857d 881 AliVTrack* vtrack = dynamic_cast<AliVTrack*>(track);
882
1f6fff78 883 if (!vtrack)
c3ba2d3d 884 continue;
885
43032ce2 886 if ((vtrack->GetTrackEtaOnEMCal() == -999 || vtrack->GetTrackPhiOnEMCal() == -999) && fHistTrPhiEtaNonProp[fCentBin]) {
56bd3193 887 fHistTrPhiEtaNonProp[fCentBin]->Fill(vtrack->Eta(), vtrack->Phi());
43032ce2 888 fHistTrPtNonProp[fCentBin]->Fill(vtrack->Pt());
889 }
56bd3193 890
891 if (fHistTrEmcPhiEta[fCentBin])
892 fHistTrEmcPhiEta[fCentBin]->Fill(vtrack->GetTrackEtaOnEMCal(), vtrack->GetTrackPhiOnEMCal());
893 if (fHistTrEmcPt[fCentBin])
894 fHistTrEmcPt[fCentBin]->Fill(vtrack->GetTrackPtOnEMCal());
895 if (fHistDeltaEtaPt[fCentBin])
896 fHistDeltaEtaPt[fCentBin]->Fill(vtrack->Pt(), vtrack->Eta() - vtrack->GetTrackEtaOnEMCal());
897 if (fHistDeltaPhiPt[fCentBin])
898 fHistDeltaPhiPt[fCentBin]->Fill(vtrack->Pt(), vtrack->Phi() - vtrack->GetTrackPhiOnEMCal());
e2b76953 899 if (fHistDeltaPtvsPt[fCentBin])
900 fHistDeltaPtvsPt[fCentBin]->Fill(vtrack->Pt(), vtrack->Pt() - vtrack->GetTrackPtOnEMCal());
c3ba2d3d 901 }
4358e58a 902
f30ab63b 903 if (fHistTrNegativeLabels[fCentBin])
904 fHistTrNegativeLabels[fCentBin]->Fill(1. * neg / ntracks);
4358e58a 905
f30ab63b 906 if (fHistTrZeroLabels[fCentBin])
907 fHistTrZeroLabels[fCentBin]->Fill(1. * zero / ntracks);
4358e58a 908
43032ce2 909 return nAccTracks;
c3ba2d3d 910}
911
912//________________________________________________________________________
7cd832c7 913Int_t AliAnalysisTaskSAQA::DoJetLoop(AliEmcalJet* &leading)
c3ba2d3d 914{
c3604199 915 // Do jet loop.
916
c3ba2d3d 917 if (!fJets)
43032ce2 918 return 0;
919
920 Int_t nAccJets = 0;
c3ba2d3d 921
26516bb5 922 leading = 0;
923
c3ba2d3d 924 Int_t njets = fJets->GetEntriesFast();
925
926 for (Int_t ij = 0; ij < njets; ij++) {
927
e44e8726 928 AliEmcalJet* jet = static_cast<AliEmcalJet*>(fJets->At(ij));
c3ba2d3d 929
930 if (!jet) {
931 AliError(Form("Could not receive jet %d", ij));
932 continue;
933 }
934
a487deae 935 if (!AcceptJet(jet))
c3ba2d3d 936 continue;
6fd5039f 937
7cd832c7 938 if (!leading || leading->Pt() < jet->Pt()) leading = jet;
26516bb5 939
43032ce2 940 nAccJets++;
7c1d624c 941
43032ce2 942 fHistJetsPhiEta[fCentBin]->Fill(jet->Eta(), jet->Phi());
99cfd012 943 fHistJetsPtArea[fCentBin]->Fill(jet->Pt(), jet->Area());
c3ba2d3d 944 }
c3ba2d3d 945
43032ce2 946 return nAccJets;
c3ba2d3d 947}