]> git.uio.no Git - u/mrichter/AliRoot.git/blame - PWGJE/EMCALJetTasks/UserTasks/AliAnalysisTaskSAQA.cxx
From Salvatore:
[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() :
99cfd012 38 AliAnalysisTaskEmcalJet("AliAnalysisTaskSAQA", kTRUE),
c3ba2d3d 39 fCellEnergyCut(0.1),
1f9c287f 40 fParticleLevel(kFALSE),
5be3857d 41 fIsMC(kFALSE),
43032ce2 42 fCentMethod2(""),
43 fCentMethod3(""),
44 fDoV0QA(0),
45 fCent2(0),
46 fCent3(0),
47 fVZERO(0),
48 fV0ATotMult(0),
49 fV0CTotMult(0),
50 fHistEventQA(0),
4358e58a 51 fHistTrNegativeLabels(0),
52 fHistTrZeroLabels(0),
7cf4626b 53 fHistNCellsEnergy(0),
59f16b27 54 fHistFcrossEnergy(0),
090a0c3e 55 fHistClusTimeEnergy(0),
f483218e 56 fHistCellsAbsIdEnergy(0),
c3ba2d3d 57 fHistChVSneCells(0),
58 fHistChVSneClus(0),
59 fHistChVSneCorrCells(0)
60{
61 // Default constructor.
62
6fd5039f 63 for (Int_t i = 0; i < 4; i++) {
6f18d73a 64 for (Int_t j = 0; j < 4; j++) fHistTrPhiEtaPt[i][j] = 0;
43032ce2 65 fHistTrPhiEtaZeroLab[i] = 0;
66 fHistTrPtZeroLab[i] = 0;
56bd3193 67 fHistTrEmcPhiEta[i] = 0;
68 fHistTrEmcPt[i] = 0;
69 fHistTrPhiEtaNonProp[i] = 0;
43032ce2 70 fHistTrPtNonProp[i] = 0;
56bd3193 71 fHistDeltaEtaPt[i] = 0;
72 fHistDeltaPhiPt[i] = 0;
73 fHistDeltaPtvsPt[i] = 0;
a487deae 74 fHistClusPhiEtaEnergy[i] = 0;
8d3d1996 75 fHistClusMCEnergyFraction[i] = 0;
43032ce2 76 fHistJetsPhiEta[i] = 0;
98750b70 77 fHistJetsPtArea[i] = 0;
6fd5039f 78 }
a487deae 79
80 SetMakeGeneralHistograms(kTRUE);
c3ba2d3d 81}
82
83//________________________________________________________________________
84AliAnalysisTaskSAQA::AliAnalysisTaskSAQA(const char *name) :
99cfd012 85 AliAnalysisTaskEmcalJet(name, kTRUE),
c3ba2d3d 86 fCellEnergyCut(0.1),
1f9c287f 87 fParticleLevel(kFALSE),
5be3857d 88 fIsMC(kFALSE),
43032ce2 89 fCentMethod2(""),
90 fCentMethod3(""),
91 fDoV0QA(0),
92 fCent2(0),
93 fCent3(0),
94 fVZERO(0),
95 fV0ATotMult(0),
96 fV0CTotMult(0),
97 fHistEventQA(0),
4358e58a 98 fHistTrNegativeLabels(0),
99 fHistTrZeroLabels(0),
7cf4626b 100 fHistNCellsEnergy(0),
59f16b27 101 fHistFcrossEnergy(0),
090a0c3e 102 fHistClusTimeEnergy(0),
f483218e 103 fHistCellsAbsIdEnergy(0),
c3ba2d3d 104 fHistChVSneCells(0),
105 fHistChVSneClus(0),
106 fHistChVSneCorrCells(0)
107{
108 // Standard constructor.
109
6fd5039f 110 for (Int_t i = 0; i < 4; i++) {
6f18d73a 111 for (Int_t j = 0; j < 4; j++) fHistTrPhiEtaPt[i][j] = 0;
43032ce2 112 fHistTrPhiEtaZeroLab[i] = 0;
113 fHistTrPtZeroLab[i] = 0;
56bd3193 114 fHistTrEmcPhiEta[i] = 0;
115 fHistTrEmcPt[i] = 0;
116 fHistTrPhiEtaNonProp[i] = 0;
43032ce2 117 fHistTrPtNonProp[i] = 0;
56bd3193 118 fHistDeltaEtaPt[i] = 0;
119 fHistDeltaPhiPt[i] = 0;
120 fHistDeltaPtvsPt[i] = 0;
a487deae 121 fHistClusPhiEtaEnergy[i] = 0;
8d3d1996 122 fHistClusMCEnergyFraction[i] = 0;
43032ce2 123 fHistJetsPhiEta[i] = 0;
98750b70 124 fHistJetsPtArea[i] = 0;
6fd5039f 125 }
a487deae 126
127 SetMakeGeneralHistograms(kTRUE);
c3ba2d3d 128}
129
130//________________________________________________________________________
131AliAnalysisTaskSAQA::~AliAnalysisTaskSAQA()
132{
133 // Destructor
134}
135
136//________________________________________________________________________
137void AliAnalysisTaskSAQA::UserCreateOutputObjects()
138{
139 // Create histograms
140
a487deae 141 AliAnalysisTaskEmcalJet::UserCreateOutputObjects();
c3ba2d3d 142
a487deae 143 if (!fTracksName.IsNull()) {
4358e58a 144 if (!fParticleLevel && fIsMC) {
145 fHistTrNegativeLabels = new TH1F("fHistTrNegativeLabels","fHistTrNegativeLabels", 500, 0, 1);
146 fHistTrNegativeLabels->GetXaxis()->SetTitle("% of negative labels");
147 fHistTrNegativeLabels->GetYaxis()->SetTitle("counts");
148 fOutput->Add(fHistTrNegativeLabels);
149
150 fHistTrZeroLabels = new TH1F("fHistTrZeroLabels","fHistTrZeroLabels", 500, 0, 1);
151 fHistTrZeroLabels->GetXaxis()->SetTitle("% of negative labels");
152 fHistTrZeroLabels->GetYaxis()->SetTitle("counts");
153 fOutput->Add(fHistTrZeroLabels);
154 }
155
a487deae 156 TString histname;
be7b6e63 157
393ec2ad 158 Int_t nlabels = 4;
159 if (fParticleLevel)
160 nlabels = 1;
161
162 for (Int_t i = 0; i < fNcentBins; i++) {
163 for (Int_t j = 0; j < nlabels; j++) {
6f18d73a 164 histname = Form("fHistTrPhiEtaPt_%d_%d",i,j);
165 fHistTrPhiEtaPt[i][j] = new TH3F(histname,histname, 100, -1, 1, 201, 0, TMath::Pi() * 2.01, fNbins, fMinBinPt, fMaxBinPt);
166 fHistTrPhiEtaPt[i][j]->GetXaxis()->SetTitle("#eta");
167 fHistTrPhiEtaPt[i][j]->GetYaxis()->SetTitle("#phi");
168 fHistTrPhiEtaPt[i][j]->GetZaxis()->SetTitle("p_{T} (GeV/c)");
169 fOutput->Add(fHistTrPhiEtaPt[i][j]);
170 }
c3ba2d3d 171
56bd3193 172 if (!fParticleLevel) {
173 if (fIsMC) {
43032ce2 174 histname = Form("fHistTrPhiEtaZeroLab_%d",i);
175 fHistTrPhiEtaZeroLab[i] = new TH2F(histname,histname, 100, -1, 1, 201, 0, TMath::Pi() * 2.01);
176 fHistTrPhiEtaZeroLab[i]->GetXaxis()->SetTitle("#eta");
177 fHistTrPhiEtaZeroLab[i]->GetYaxis()->SetTitle("#phi");
178 fHistTrPhiEtaZeroLab[i]->GetZaxis()->SetTitle("counts");
179 fOutput->Add(fHistTrPhiEtaZeroLab[i]);
180
181 histname = Form("fHistTrPtZeroLab_%d",i);
182 fHistTrPtZeroLab[i] = new TH1F(histname,histname, fNbins, fMinBinPt, fMaxBinPt);
183 fHistTrPtZeroLab[i]->GetZaxis()->SetTitle("p_{T} (GeV/c)");
184 fHistTrPtZeroLab[i]->GetYaxis()->SetTitle("counts");
185 fOutput->Add(fHistTrPtZeroLab[i]);
56bd3193 186 }
187
188 histname = Form("fHistTrEmcPhiEta_%d",i);
189 fHistTrEmcPhiEta[i] = new TH2F(histname,histname, 100, -1, 1, 201, 0, TMath::Pi() * 2.01);
190 fHistTrEmcPhiEta[i]->GetXaxis()->SetTitle("#eta");
191 fHistTrEmcPhiEta[i]->GetYaxis()->SetTitle("#phi");
192 fOutput->Add(fHistTrEmcPhiEta[i]);
193
194 histname = Form("fHistTrEmcPt_%d",i);
195 fHistTrEmcPt[i] = new TH1F(histname,histname, fNbins, fMinBinPt, fMaxBinPt);
196 fHistTrEmcPt[i]->GetXaxis()->SetTitle("p_{T} (GeV/c)");
197 fHistTrEmcPt[i]->GetYaxis()->SetTitle("counts");
198 fOutput->Add(fHistTrEmcPt[i]);
199
200 histname = Form("fHistTrPhiEtaNonProp_%d",i);
201 fHistTrPhiEtaNonProp[i] = new TH2F(histname,histname, 100, -1, 1, 201, 0, TMath::Pi() * 2.01);
202 fHistTrPhiEtaNonProp[i]->GetXaxis()->SetTitle("#eta");
203 fHistTrPhiEtaNonProp[i]->GetYaxis()->SetTitle("#phi");
204 fOutput->Add(fHistTrPhiEtaNonProp[i]);
43032ce2 205
206 histname = Form("fHistTrPtNonProp_%d",i);
207 fHistTrPtNonProp[i] = new TH1F(histname,histname, fNbins, fMinBinPt, fMaxBinPt);
208 fHistTrPtNonProp[i]->GetXaxis()->SetTitle("p_{T} (GeV/c)");
209 fHistTrPtNonProp[i]->GetYaxis()->SetTitle("counts");
210 fOutput->Add(fHistTrPtNonProp[i]);
56bd3193 211
212 histname = Form("fHistDeltaEtaPt_%d",i);
213 fHistDeltaEtaPt[i] = new TH2F(histname,histname, fNbins, fMinBinPt, fMaxBinPt, 50, -0.5, 0.5);
214 fHistDeltaEtaPt[i]->GetXaxis()->SetTitle("p_{T} (GeV/c)");
215 fHistDeltaEtaPt[i]->GetYaxis()->SetTitle("#delta#eta");
216 fOutput->Add(fHistDeltaEtaPt[i]);
217
218 histname = Form("fHistDeltaPhiPt_%d",i);
219 fHistDeltaPhiPt[i] = new TH2F(histname,histname, fNbins, fMinBinPt, fMaxBinPt, 200, -2, 2);
220 fHistDeltaPhiPt[i]->GetXaxis()->SetTitle("p_{T} (GeV/c)");
221 fHistDeltaPhiPt[i]->GetYaxis()->SetTitle("#delta#phi");
222 fOutput->Add(fHistDeltaPhiPt[i]);
223
224 histname = Form("fHistDeltaPtvsPt_%d",i);
225 fHistDeltaPtvsPt[i] = new TH2F(histname,histname, fNbins, fMinBinPt, fMaxBinPt, fNbins, -fMaxBinPt/2, fMaxBinPt/2);
226 fHistDeltaPtvsPt[i]->GetXaxis()->SetTitle("p_{T} (GeV/c)");
227 fHistDeltaPtvsPt[i]->GetYaxis()->SetTitle("#deltap_{T} (GeV/c)");
228 fOutput->Add(fHistDeltaPtvsPt[i]);
229 }
2103dc6a 230 }
a487deae 231 }
7c1d624c 232
a487deae 233 if (!fCaloName.IsNull()) {
a487deae 234 TString histname;
235
393ec2ad 236 for (Int_t i = 0; i < fNcentBins; i++) {
a487deae 237 histname = "fHistClusPhiEtaEnergy_";
238 histname += i;
239 fHistClusPhiEtaEnergy[i] = new TH3F(histname, histname, 100, -1.2, 1.2, 201, 0, TMath::Pi() * 2.01, fNbins, fMinBinPt, fMaxBinPt);
240 fHistClusPhiEtaEnergy[i]->GetXaxis()->SetTitle("#eta");
241 fHistClusPhiEtaEnergy[i]->GetYaxis()->SetTitle("#phi");
59f16b27 242 fHistClusPhiEtaEnergy[i]->GetZaxis()->SetTitle("E_{cluster} (GeV)");
a487deae 243 fOutput->Add(fHistClusPhiEtaEnergy[i]);
8d3d1996 244
245 if (fIsEmbedded) {
246 histname = "fHistClusMCEnergyFraction_";
247 histname += i;
248 fHistClusMCEnergyFraction[i] = new TH1F(histname, histname, fNbins, 0, 1.2);
249 fHistClusMCEnergyFraction[i]->GetXaxis()->SetTitle("MC fraction");
250 fHistClusMCEnergyFraction[i]->GetYaxis()->SetTitle("counts");
251 fOutput->Add(fHistClusMCEnergyFraction[i]);
252 }
a487deae 253 }
254
781da0a3 255 fHistClusTimeEnergy = new TH2F("fHistClusTimeEnergy","Time vs. energy of clusters", fNbins, fMinBinPt, fMaxBinPt, fNbins, -1e-6, 1e-6);
59f16b27 256 fHistClusTimeEnergy->GetXaxis()->SetTitle("E_{cluster} (GeV)");
090a0c3e 257 fHistClusTimeEnergy->GetYaxis()->SetTitle("Time");
258 fOutput->Add(fHistClusTimeEnergy);
259
7cf4626b 260 fHistNCellsEnergy = new TH2F("fHistNCellsEnergy","Number of cells vs. energy of clusters", fNbins, fMinBinPt, fMaxBinPt, 30, 0, 30);
59f16b27 261 fHistNCellsEnergy->GetXaxis()->SetTitle("E_{cluster} (GeV)");
7cf4626b 262 fHistNCellsEnergy->GetYaxis()->SetTitle("N_{cells}");
a487deae 263 fOutput->Add(fHistNCellsEnergy);
59f16b27 264
265 fHistFcrossEnergy = new TH2F("fHistFcrossEnergy","fHistFcrossEnergy", fNbins, fMinBinPt, fMaxBinPt, 200, -3.5, 1.5);
266 fHistFcrossEnergy->GetXaxis()->SetTitle("E_{cluster} (GeV)");
267 fHistFcrossEnergy->GetYaxis()->SetTitle("F_{cross}");
268 fOutput->Add(fHistFcrossEnergy);
a487deae 269
f483218e 270 fHistCellsAbsIdEnergy = new TH2F("fHistCellsAbsIdEnergy","fHistCellsAbsIdEnergy", 11600,0,11599,(Int_t)(fNbins / 2), fMinBinPt, fMaxBinPt / 2);
271 fHistCellsAbsIdEnergy->GetXaxis()->SetTitle("cell abs. Id");
59f16b27 272 fHistCellsAbsIdEnergy->GetYaxis()->SetTitle("E_{cluster} (GeV)");
f483218e 273 fHistCellsAbsIdEnergy->GetZaxis()->SetTitle("counts");
274 fOutput->Add(fHistCellsAbsIdEnergy);
6fd5039f 275
98750b70 276 fHistChVSneCells = new TH2F("fHistChVSneCells","Charged energy vs. neutral (cells) energy",
be7b6e63 277 (Int_t)(fNbins * 2.5), fMinBinPt, fMaxBinPt * 2.5, (Int_t)(fNbins * 2.5), fMinBinPt, fMaxBinPt * 2.5);
a487deae 278 fHistChVSneCells->GetXaxis()->SetTitle("Energy (GeV)");
279 fHistChVSneCells->GetYaxis()->SetTitle("Momentum (GeV/c)");
6fd5039f 280 fOutput->Add(fHistChVSneCells);
281
98750b70 282 fHistChVSneClus = new TH2F("fHistChVSneClus","Charged energy vs. neutral (clusters) energy",
be7b6e63 283 (Int_t)(fNbins * 2.5), fMinBinPt, fMaxBinPt * 2.5, (Int_t)(fNbins * 2.5), fMinBinPt, fMaxBinPt * 2.5);
a487deae 284 fHistChVSneClus->GetXaxis()->SetTitle("Energy (GeV)");
285 fHistChVSneClus->GetYaxis()->SetTitle("Momentum (GeV/c)");
6fd5039f 286 fOutput->Add(fHistChVSneClus);
287
98750b70 288 fHistChVSneCorrCells = new TH2F("fHistChVSneCorrCells","Charged energy vs. neutral (corrected cells) energy",
be7b6e63 289 (Int_t)(fNbins * 2.5), fMinBinPt, fMaxBinPt * 2.5, (Int_t)(fNbins * 2.5), fMinBinPt , fMaxBinPt * 2.5);
a487deae 290 fHistChVSneCorrCells->GetXaxis()->SetTitle("Energy (GeV)");
291 fHistChVSneCorrCells->GetYaxis()->SetTitle("Momentum (GeV/c)");
6fd5039f 292 fOutput->Add(fHistChVSneCorrCells);
c3ba2d3d 293 }
a487deae 294
e44e8726 295 if (!fJetsName.IsNull()) {
6fd5039f 296
e44e8726 297 TString histname;
298
393ec2ad 299 for (Int_t i = 0; i < fNcentBins; i++) {
43032ce2 300 histname = "fHistJetsPhiEta_";
e44e8726 301 histname += i;
43032ce2 302 fHistJetsPhiEta[i] = new TH2F(histname.Data(), histname.Data(), 100, -1.2, 1.2, 201, 0, TMath::Pi() * 2.01);
303 fHistJetsPhiEta[i]->GetXaxis()->SetTitle("#eta");
304 fHistJetsPhiEta[i]->GetYaxis()->SetTitle("#phi");
305 fHistJetsPhiEta[i]->GetZaxis()->SetTitle("counts");
306 fOutput->Add(fHistJetsPhiEta[i]);
e44e8726 307
308 histname = "fHistJetsPtArea_";
309 histname += i;
a487deae 310 fHistJetsPtArea[i] = new TH2F(histname.Data(), histname.Data(), (Int_t)(fNbins * 2.5), fMinBinPt, fMaxBinPt * 2.5, 30, 0, fJetRadius * fJetRadius * TMath::Pi() * 3);
311 fHistJetsPtArea[i]->GetXaxis()->SetTitle("p_{T} (GeV/c)");
e44e8726 312 fHistJetsPtArea[i]->GetYaxis()->SetTitle("area");
313 fOutput->Add(fHistJetsPtArea[i]);
314 }
6fd5039f 315 }
43032ce2 316
317 Int_t dim = 0;
318 TString title[10];
319 Int_t nbins[10] = {0};
320 Double_t min[10] = {0};
321 Double_t max[10] = {0};
322
323 if (fForceBeamType != AliAnalysisTaskEmcal::kpp) {
324 title[dim] = "Centrality %";
325 nbins[dim] = 101;
326 min[dim] = 0;
327 max[dim] = 101;
328 dim++;
329
330 if (!fCentMethod2.IsNull()) {
331 title[dim] = Form("Centrality %s %%", fCentMethod2.Data());
332 nbins[dim] = 101;
333 min[dim] = 0;
334 max[dim] = 101;
335 dim++;
336 }
337
338 if (!fCentMethod3.IsNull()) {
339 title[dim] = Form("Centrality %s %%", fCentMethod3.Data());
340 nbins[dim] = 101;
341 min[dim] = 0;
342 max[dim] = 101;
343 dim++;
344 }
345
346 if (fDoV0QA) {
347 title[dim] = "V0A total multiplicity";
348 nbins[dim] = 200;
349 min[dim] = 0;
350 max[dim] = 20000;
351 dim++;
352
353 title[dim] = "V0C total multiplicity";
354 nbins[dim] = 200;
355 min[dim] = 0;
356 max[dim] = 20000;
357 dim++;
358 }
359
360 if (!fRhoName.IsNull()) {
361 title[dim] = "#rho (GeV/c)";
362 nbins[dim] = fNbins*4;
363 min[dim] = 0;
364 max[dim] = 400;
365 dim++;
366 }
367 }
368
369 if (!fTracksName.IsNull()) {
370 title[dim] = "No. of tracks";
371 nbins[dim] = 6000;
372 min[dim] = -0.5;
373 max[dim] = 6000-0.5;
374 dim++;
375 }
376
377 if (!fCaloName.IsNull()) {
378 title[dim] = "No. of clusters";
379 nbins[dim] = 4000;
380 min[dim] = 0;
381 max[dim] = 4000-0.5;
382 dim++;
383 }
384
385 if (!fCaloCellsName.IsNull()) {
386 title[dim] = "No. of cells";
387 nbins[dim] = 6000;
388 min[dim] = 0;
389 max[dim] = 6000-0.5;
390 dim++;
391 }
392
393 if (!fJetsName.IsNull()) {
394 title[dim] = "No. of jets";
395 nbins[dim] = 200;
396 min[dim] = 0;
397 max[dim] = 200-0.5;
398 dim++;
399 }
400
401 fHistEventQA = new THnSparseF("fHistEventQA","fHistEventQA",dim,nbins,min,max);
402 for (Int_t i = 0; i < dim; i++)
403 fHistEventQA->GetAxis(i)->SetTitle(title[i]);
404 fOutput->Add(fHistEventQA);
6fd5039f 405
c3ba2d3d 406 PostData(1, fOutput); // Post data for ALL output slots >0 here, to get at least an empty histogram
407}
408
43032ce2 409//________________________________________________________________________
410void AliAnalysisTaskSAQA::ExecOnce()
411{
412 AliAnalysisTaskEmcalJet::ExecOnce();
413
414 if (fDoV0QA) {
415 fVZERO = InputEvent()->GetVZEROData();
416 if (!fVZERO) {
417 AliError("AliVVZERO not available");
418 }
419 }
420}
421
c3ba2d3d 422//________________________________________________________________________
6fd5039f 423Bool_t AliAnalysisTaskSAQA::RetrieveEventObjects()
c3ba2d3d 424{
c3604199 425 // Retrieve event objects.
426
427 if (!AliAnalysisTaskEmcalJet::RetrieveEventObjects())
6fd5039f 428 return kFALSE;
c3ba2d3d 429
43032ce2 430 if (!fCentMethod2.IsNull() || !fCentMethod3.IsNull()) {
431 if (fBeamType == kAA || fBeamType == kpA ) {
432 AliCentrality *aliCent = InputEvent()->GetCentrality();
433 if (aliCent) {
434 if (!fCentMethod2.IsNull())
435 fCent2 = aliCent->GetCentralityPercentile(fCentMethod2);
436 if (!fCentMethod3.IsNull())
437 fCent3 = aliCent->GetCentralityPercentile(fCentMethod3);
438 }
439 }
440 }
441
442 if (fVZERO) {
443 fV0ATotMult = AliESDUtils::GetCorrV0A(fVZERO->GetMTotV0A(),fVertex[2]);
444 fV0CTotMult = AliESDUtils::GetCorrV0C(fVZERO->GetMTotV0C(),fVertex[2]);
445 }
7c1d624c 446
6fd5039f 447 return kTRUE;
c3ba2d3d 448}
449
7c1d624c 450
c3ba2d3d 451//________________________________________________________________________
6fd5039f 452Bool_t AliAnalysisTaskSAQA::FillHistograms()
c3ba2d3d 453{
c3604199 454 // Fill histograms.
455
a487deae 456 Float_t trackSum = 0;
43032ce2 457 Float_t clusSum = 0;
458 Float_t cellSum = 0;
459 Float_t cellCutSum = 0;
090a0c3e 460
43032ce2 461 Int_t ntracks = 0;
462 Int_t nclusters = 0;
463 Int_t ncells = 0;
464 Int_t njets = 0;
465
a487deae 466 if (fTracks) {
43032ce2 467 ntracks = DoTrackLoop(trackSum);
468 AliDebug(2,Form("%d tracks found in the event", ntracks));
a487deae 469 }
be7b6e63 470
a487deae 471 if (fCaloClusters) {
43032ce2 472 nclusters = DoClusterLoop(clusSum);
473 AliDebug(2,Form("%d clusters found in the event", nclusters));
c3ba2d3d 474
43032ce2 475 fHistChVSneClus->Fill(clusSum, trackSum);
476 }
477
478 if (fCaloCells) {
479 ncells = DoCellLoop(cellSum, cellCutSum);
480 AliDebug(2,Form("%d cells found in the event", ncells));
6fd5039f 481
482 fHistChVSneCells->Fill(cellSum, trackSum);
6fd5039f 483 fHistChVSneCorrCells->Fill(cellCutSum, trackSum);
a487deae 484 }
c3ba2d3d 485
a487deae 486 if (fJets) {
43032ce2 487 njets = DoJetLoop();
488 AliDebug(2,Form("%d jets found in the event", njets));
6fd5039f 489 }
c3ba2d3d 490
43032ce2 491 FillEventQAHisto(fCent, fCent2, fCent3, fV0ATotMult, fV0CTotMult, fRhoVal, ntracks, nclusters, ncells, njets);
492
6fd5039f 493 return kTRUE;
c3ba2d3d 494}
495
43032ce2 496//________________________________________________________________________
497void AliAnalysisTaskSAQA::FillEventQAHisto(Float_t cent, Float_t cent2, Float_t cent3, Float_t v0a, Float_t v0c, Float_t rho, Int_t ntracks, Int_t nclusters, Int_t ncells, Int_t njets)
498{
499 Double_t contents[10]={0};
500
501 for (Int_t i = 0; i < fHistEventQA->GetNdimensions(); i++) {
502 TString title(fHistEventQA->GetAxis(i)->GetTitle());
503 if (title=="Centrality %")
504 contents[i] = cent;
505 else if (title==Form("Centrality %s %%", fCentMethod2.Data()))
506 contents[i] = cent2;
507 else if (title==Form("Centrality %s %%", fCentMethod3.Data()))
508 contents[i] = cent3;
509 else if (title=="V0A total multiplicity")
510 contents[i] = v0a;
511 else if (title=="V0C total multiplicity")
512 contents[i] = v0c;
513 else if (title=="#rho (GeV/c)")
514 contents[i] = rho;
515 else if (title=="No. of tracks")
516 contents[i] = ntracks;
517 else if (title=="No. of clusters")
518 contents[i] = nclusters;
519 else if (title=="No. of cells")
520 contents[i] = ncells;
521 else if (title=="No. of jets")
522 contents[i] = njets;
523 else
524 AliWarning(Form("Unable to fill dimension %s!",title.Data()));
525 }
526
527 fHistEventQA->Fill(contents);
528}
529
c3ba2d3d 530//________________________________________________________________________
be7b6e63 531Int_t AliAnalysisTaskSAQA::DoCellLoop(Float_t &sum, Float_t &sum_cut)
c3ba2d3d 532{
c3604199 533 // Do cell loop.
534
c3ba2d3d 535 AliVCaloCells *cells = InputEvent()->GetEMCALCells();
536
537 if (!cells)
be7b6e63 538 return 0;
c3ba2d3d 539
c3604199 540 const Int_t ncells = cells->GetNumberOfCells();
c3ba2d3d 541
542 for (Int_t pos = 0; pos < ncells; pos++) {
f483218e 543 Float_t amp = cells->GetAmplitude(pos);
544 Int_t absId = cells->GetCellNumber(pos);
545 fHistCellsAbsIdEnergy->Fill(absId,amp);
c3ba2d3d 546 sum += amp;
c3ba2d3d 547 if (amp < fCellEnergyCut)
548 continue;
c3ba2d3d 549 sum_cut += amp;
c3ba2d3d 550 }
be7b6e63 551
552 return ncells;
c3ba2d3d 553}
554
59f16b27 555//________________________________________________________________________
556Double_t AliAnalysisTaskSAQA::GetFcross(AliVCluster *cluster, AliVCaloCells *cells)
557{
558 Int_t AbsIdseed = -1;
559 Double_t Eseed = 0;
560 for (Int_t i = 0; i < cluster->GetNCells(); i++) {
561 if (cells->GetCellAmplitude(cluster->GetCellAbsId(i)) > AbsIdseed) {
562 Eseed = cells->GetCellAmplitude(cluster->GetCellAbsId(i));
563 AbsIdseed = cluster->GetCellAbsId(i);
564 }
565 }
566
567 if (Eseed < 1e-9)
568 return 100;
569
570 Int_t imod = -1, iphi =-1, ieta=-1,iTower = -1, iIphi = -1, iIeta = -1;
571 fGeom->GetCellIndex(AbsIdseed,imod,iTower,iIphi,iIeta);
572 fGeom->GetCellPhiEtaIndexInSModule(imod,iTower,iIphi,iIeta,iphi,ieta);
573
574 //Get close cells index and energy, not in corners
575
576 Int_t absID1 = -1;
577 Int_t absID2 = -1;
578
579 if (iphi < AliEMCALGeoParams::fgkEMCALRows-1) absID1 = fGeom->GetAbsCellIdFromCellIndexes(imod, iphi+1, ieta);
580 if (iphi > 0) absID2 = fGeom->GetAbsCellIdFromCellIndexes(imod, iphi-1, ieta);
581
582 // In case of cell in eta = 0 border, depending on SM shift the cross cell index
583
584 Int_t absID3 = -1;
585 Int_t absID4 = -1;
586
587 if (ieta == AliEMCALGeoParams::fgkEMCALCols-1 && !(imod%2)) {
588 absID3 = fGeom->GetAbsCellIdFromCellIndexes(imod+1, iphi, 0);
589 absID4 = fGeom->GetAbsCellIdFromCellIndexes(imod, iphi, ieta-1);
590 }
591 else if (ieta == 0 && imod%2) {
592 absID3 = fGeom->GetAbsCellIdFromCellIndexes(imod, iphi, ieta+1);
593 absID4 = fGeom->GetAbsCellIdFromCellIndexes(imod-1, iphi, AliEMCALGeoParams::fgkEMCALCols-1);
594 }
595 else {
596 if (ieta < AliEMCALGeoParams::fgkEMCALCols-1)
597 absID3 = fGeom->GetAbsCellIdFromCellIndexes(imod, iphi, ieta+1);
598 if (ieta > 0)
599 absID4 = fGeom->GetAbsCellIdFromCellIndexes(imod, iphi, ieta-1);
600 }
601
602 Double_t ecell1 = cells->GetCellAmplitude(absID1);
603 Double_t ecell2 = cells->GetCellAmplitude(absID2);
604 Double_t ecell3 = cells->GetCellAmplitude(absID3);
605 Double_t ecell4 = cells->GetCellAmplitude(absID4);
606
607 Double_t Ecross = ecell1 + ecell2 + ecell3 + ecell4;
608
609 Double_t Fcross = 1 - Ecross/Eseed;
610
611 return Fcross;
612}
613
c3ba2d3d 614//________________________________________________________________________
43032ce2 615Int_t AliAnalysisTaskSAQA::DoClusterLoop(Float_t &sum)
c3ba2d3d 616{
c3604199 617 // Do cluster loop.
618
c3ba2d3d 619 if (!fCaloClusters)
620 return 0;
621
43032ce2 622 Int_t nAccClusters = 0;
623
59f16b27 624 AliVCaloCells *cells = InputEvent()->GetEMCALCells();
625
43032ce2 626 sum = 0;
c3ba2d3d 627
628 // Cluster loop
43032ce2 629 const Int_t nclusters = fCaloClusters->GetEntriesFast();
c3ba2d3d 630
631 for (Int_t iClusters = 0; iClusters < nclusters; iClusters++) {
e44e8726 632 AliVCluster* cluster = static_cast<AliVCluster*>(fCaloClusters->At(iClusters));
c3ba2d3d 633 if (!cluster) {
634 AliError(Form("Could not receive cluster %d", iClusters));
635 continue;
636 }
c3ba2d3d 637
f660c2d6 638 if (!AcceptCluster(cluster))
e44e8726 639 continue;
c3ba2d3d 640
c3ba2d3d 641 sum += cluster->E();
642
7cf4626b 643 TLorentzVector nPart;
644 cluster->GetMomentum(nPart, fVertex);
c3ba2d3d 645
a487deae 646 fHistClusPhiEtaEnergy[fCentBin]->Fill(nPart.Eta(), nPart.Phi(), cluster->E());
7cf4626b 647 fHistNCellsEnergy->Fill(cluster->E(), cluster->GetNCells());
090a0c3e 648
649 fHistClusTimeEnergy->Fill(cluster->E(), cluster->GetTOF());
7c1d624c 650
59f16b27 651 if (cells)
652 fHistFcrossEnergy->Fill(cluster->E(), GetFcross(cluster, cells));
653
8d3d1996 654 if (fHistClusMCEnergyFraction[fCentBin])
655 fHistClusMCEnergyFraction[fCentBin]->Fill(cluster->GetMCEnergyFraction());
1b47ed10 656
43032ce2 657 nAccClusters++;
c3604199 658 }
c3ba2d3d 659
43032ce2 660 return nAccClusters;
c3ba2d3d 661}
662
663//________________________________________________________________________
43032ce2 664Int_t AliAnalysisTaskSAQA::DoTrackLoop(Float_t &sum)
c3ba2d3d 665{
c3604199 666 // Do track loop.
667
c3ba2d3d 668 if (!fTracks)
669 return 0;
670
43032ce2 671 Int_t nAccTracks = 0;
672
673 sum = 0;
c3ba2d3d 674
43032ce2 675 const Int_t ntracks = fTracks->GetEntriesFast();
4358e58a 676 Int_t neg = 0;
677 Int_t zero = 0;
c3ba2d3d 678
c3604199 679 for (Int_t i = 0; i < ntracks; i++) {
c3ba2d3d 680
e44e8726 681 AliVParticle* track = static_cast<AliVParticle*>(fTracks->At(i)); // pointer to reconstructed to track
1f6fff78 682
c3604199 683 if (!track) {
c3ba2d3d 684 AliError(Form("Could not retrieve track %d",i));
685 continue;
686 }
c3ba2d3d 687
5be3857d 688 if (!AcceptTrack(track))
1f6fff78 689 continue;
7c1d624c 690
43032ce2 691 nAccTracks++;
c3ba2d3d 692
693 sum += track->P();
694
393ec2ad 695 if (fParticleLevel) {
696 fHistTrPhiEtaPt[fCentBin][0]->Fill(track->Eta(), track->Phi(), track->Pt());
697 }
698 else {
699 fHistTrPhiEtaPt[fCentBin][3]->Fill(track->Eta(), track->Phi(), track->Pt());
4358e58a 700 if (track->GetLabel() == 0) {
701 zero++;
43032ce2 702 if (fHistTrPhiEtaZeroLab[fCentBin]) {
703 fHistTrPhiEtaZeroLab[fCentBin]->Fill(track->Eta(), track->Phi());
704 fHistTrPtZeroLab[fCentBin]->Fill(track->Pt());
705 }
4358e58a 706 }
707
708 if (track->GetLabel() < 0)
709 neg++;
787a3c4f 710
b16bb001 711 Int_t type = 0;
393ec2ad 712
b16bb001 713 AliPicoTrack* ptrack = dynamic_cast<AliPicoTrack*>(track);
714 if (ptrack)
715 type = ptrack->GetTrackType();
716
717 if (type >= 0 && type < 3)
718 fHistTrPhiEtaPt[fCentBin][type]->Fill(track->Eta(), track->Phi(), track->Pt());
393ec2ad 719 else
5be3857d 720 AliDebug(2,Form("%s: track type %d not recognized!", GetName(), type));
393ec2ad 721 }
aa4d701c 722
5be3857d 723 AliVTrack* vtrack = dynamic_cast<AliVTrack*>(track);
724
1f6fff78 725 if (!vtrack)
c3ba2d3d 726 continue;
727
43032ce2 728 if ((vtrack->GetTrackEtaOnEMCal() == -999 || vtrack->GetTrackPhiOnEMCal() == -999) && fHistTrPhiEtaNonProp[fCentBin]) {
56bd3193 729 fHistTrPhiEtaNonProp[fCentBin]->Fill(vtrack->Eta(), vtrack->Phi());
43032ce2 730 fHistTrPtNonProp[fCentBin]->Fill(vtrack->Pt());
731 }
56bd3193 732
733 if (fHistTrEmcPhiEta[fCentBin])
734 fHistTrEmcPhiEta[fCentBin]->Fill(vtrack->GetTrackEtaOnEMCal(), vtrack->GetTrackPhiOnEMCal());
735 if (fHistTrEmcPt[fCentBin])
736 fHistTrEmcPt[fCentBin]->Fill(vtrack->GetTrackPtOnEMCal());
737 if (fHistDeltaEtaPt[fCentBin])
738 fHistDeltaEtaPt[fCentBin]->Fill(vtrack->Pt(), vtrack->Eta() - vtrack->GetTrackEtaOnEMCal());
739 if (fHistDeltaPhiPt[fCentBin])
740 fHistDeltaPhiPt[fCentBin]->Fill(vtrack->Pt(), vtrack->Phi() - vtrack->GetTrackPhiOnEMCal());
741 if (fHistDeltaPtvsPt[fCentBin])
742 fHistDeltaPtvsPt[fCentBin]->Fill(vtrack->Pt(), vtrack->Pt() - vtrack->GetTrackPtOnEMCal());
c3ba2d3d 743 }
4358e58a 744
745 if (fHistTrNegativeLabels)
746 fHistTrNegativeLabels->Fill(1. * neg / ntracks);
747
748 if (fHistTrZeroLabels)
749 fHistTrZeroLabels->Fill(1. * zero / ntracks);
750
43032ce2 751 return nAccTracks;
c3ba2d3d 752}
753
754//________________________________________________________________________
43032ce2 755Int_t AliAnalysisTaskSAQA::DoJetLoop()
c3ba2d3d 756{
c3604199 757 // Do jet loop.
758
c3ba2d3d 759 if (!fJets)
43032ce2 760 return 0;
761
762 Int_t nAccJets = 0;
c3ba2d3d 763
764 Int_t njets = fJets->GetEntriesFast();
765
766 for (Int_t ij = 0; ij < njets; ij++) {
767
e44e8726 768 AliEmcalJet* jet = static_cast<AliEmcalJet*>(fJets->At(ij));
c3ba2d3d 769
770 if (!jet) {
771 AliError(Form("Could not receive jet %d", ij));
772 continue;
773 }
774
a487deae 775 if (!AcceptJet(jet))
c3ba2d3d 776 continue;
6fd5039f 777
43032ce2 778 nAccJets++;
7c1d624c 779
43032ce2 780 fHistJetsPhiEta[fCentBin]->Fill(jet->Eta(), jet->Phi());
99cfd012 781 fHistJetsPtArea[fCentBin]->Fill(jet->Pt(), jet->Area());
c3ba2d3d 782 }
c3ba2d3d 783
43032ce2 784 return nAccJets;
c3ba2d3d 785}