updates 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
7#include <TChain.h>
8#include <TClonesArray.h>
9#include <TH1F.h>
10#include <TH2F.h>
3000c095 11#include <TH3F.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"
26
27#include "AliAnalysisTaskSAQA.h"
28
29ClassImp(AliAnalysisTaskSAQA)
30
31//________________________________________________________________________
32AliAnalysisTaskSAQA::AliAnalysisTaskSAQA() :
99cfd012 33 AliAnalysisTaskEmcalJet("AliAnalysisTaskSAQA", kTRUE),
c3ba2d3d 34 fCellEnergyCut(0.1),
1f6fff78 35 fDoTrigger(kFALSE),
11d18b51 36 fRepropagateTracks(kFALSE),
c3ba2d3d 37 fTrgClusName("ClustersL1GAMMAFEE"),
38 fTrgClusters(0),
7c1d624c 39 fNclusters(0),
40 fNtracks(0),
41 fNjets(0),
c3ba2d3d 42 fHistCentrality(0),
c92284d5 43 fHistZVertex(0),
c3ba2d3d 44 fHistTracksCent(0),
45 fHistClusCent(0),
7c1d624c 46 fHistJetsCent(0),
090a0c3e 47 fHistClusTracks(0),
7c1d624c 48 fHistJetsParts(0),
be7b6e63 49 fHistCellsCent(0),
50 fHistCellsTracks(0),
c3ba2d3d 51 fHistMaxL1FastORCent(0),
52 fHistMaxL1ClusCent(0),
53 fHistMaxL1ThrCent(0),
54 fHistTracksPt(0),
55 fHistTrPhiEta(0),
6fd5039f 56 fHistTrEmcPhiEta(0),
079b4732 57 fHistTrPhiEtaNonProp(0),
3e951e32 58 fHistDeltaEtaPt(0),
59 fHistDeltaPhiPt(0),
079b4732 60 fHistDeltaEtaNewProp(0),
61 fHistDeltaPhiNewProp(0),
7cf4626b 62 fHistClusPhiEtaEnergy(0),
63 fHistNCellsEnergy(0),
090a0c3e 64 fHistClusTimeEnergy(0),
c3ba2d3d 65 fHistCellsEnergy(0),
66 fHistChVSneCells(0),
67 fHistChVSneClus(0),
68 fHistChVSneCorrCells(0)
69{
70 // Default constructor.
71
72 for (Int_t i = 0; i < 5; i++) {
73 fHistTrackPhi[i] = 0;
74 fHistTrackEta[i] = 0;
75 }
6fd5039f 76
aa4d701c 77 for (Int_t i = 0; i < 6; i++) {
78 fHistTrackPhiPt[i] = 0;
79 fHistTrackEtaPt[i] = 0;
80 }
81
6fd5039f 82 for (Int_t i = 0; i < 4; i++) {
58285fc6 83 fHistJetsPhiEta[i] = 0;
99cfd012 84 fHistJetsPtNonBias[i] = 0;
6fd5039f 85 fHistJetsPtTrack[i] = 0;
86 fHistJetsPtClus[i] = 0;
87 fHistJetsPt[i] = 0;
98750b70 88 fHistJetsPtArea[i] = 0;
99cfd012 89 fHistJetsPtAreaNonBias[i] = 0;
6fd5039f 90 }
c3ba2d3d 91}
92
93//________________________________________________________________________
94AliAnalysisTaskSAQA::AliAnalysisTaskSAQA(const char *name) :
99cfd012 95 AliAnalysisTaskEmcalJet(name, kTRUE),
c3ba2d3d 96 fCellEnergyCut(0.1),
1f6fff78 97 fDoTrigger(kFALSE),
11d18b51 98 fRepropagateTracks(kFALSE),
c3ba2d3d 99 fTrgClusName("ClustersL1GAMMAFEE"),
100 fTrgClusters(0),
7c1d624c 101 fNclusters(0),
102 fNtracks(0),
103 fNjets(0),
c3ba2d3d 104 fHistCentrality(0),
c92284d5 105 fHistZVertex(0),
c3ba2d3d 106 fHistTracksCent(0),
107 fHistClusCent(0),
7c1d624c 108 fHistJetsCent(0),
090a0c3e 109 fHistClusTracks(0),
7c1d624c 110 fHistJetsParts(0),
be7b6e63 111 fHistCellsCent(0),
112 fHistCellsTracks(0),
c3ba2d3d 113 fHistMaxL1FastORCent(0),
114 fHistMaxL1ClusCent(0),
115 fHistMaxL1ThrCent(0),
116 fHistTracksPt(0),
117 fHistTrPhiEta(0),
6fd5039f 118 fHistTrEmcPhiEta(0),
079b4732 119 fHistTrPhiEtaNonProp(0),
3e951e32 120 fHistDeltaEtaPt(0),
121 fHistDeltaPhiPt(0),
079b4732 122 fHistDeltaEtaNewProp(0),
123 fHistDeltaPhiNewProp(0),
7cf4626b 124 fHistClusPhiEtaEnergy(0),
125 fHistNCellsEnergy(0),
090a0c3e 126 fHistClusTimeEnergy(0),
c3ba2d3d 127 fHistCellsEnergy(0),
128 fHistChVSneCells(0),
129 fHistChVSneClus(0),
130 fHistChVSneCorrCells(0)
131{
132 // Standard constructor.
133
134 for (Int_t i = 0; i < 5; i++) {
135 fHistTrackPhi[i] = 0;
136 fHistTrackEta[i] = 0;
137 }
6fd5039f 138
aa4d701c 139 for (Int_t i = 0; i < 6; i++) {
140 fHistTrackPhiPt[i] = 0;
141 fHistTrackEtaPt[i] = 0;
142 }
143
6fd5039f 144 for (Int_t i = 0; i < 4; i++) {
58285fc6 145 fHistJetsPhiEta[i] = 0;
99cfd012 146 fHistJetsPtNonBias[i] = 0;
6fd5039f 147 fHistJetsPtTrack[i] = 0;
148 fHistJetsPtClus[i] = 0;
149 fHistJetsPt[i] = 0;
98750b70 150 fHistJetsPtArea[i] = 0;
99cfd012 151 fHistJetsPtAreaNonBias[i] = 0;
6fd5039f 152 }
c3ba2d3d 153}
154
155//________________________________________________________________________
156AliAnalysisTaskSAQA::~AliAnalysisTaskSAQA()
157{
158 // Destructor
159}
160
161//________________________________________________________________________
162void AliAnalysisTaskSAQA::UserCreateOutputObjects()
163{
164 // Create histograms
165
166 OpenFile(1);
167 fOutput = new TList();
6fd5039f 168 fOutput->SetOwner();
c3ba2d3d 169
a5621834 170 fHistCentrality = new TH1F("fHistCentrality","Event centrality distribution", 100, 0, 100);
c3ba2d3d 171 fHistCentrality->GetXaxis()->SetTitle("Centrality (%)");
172 fHistCentrality->GetYaxis()->SetTitle("counts");
173 fOutput->Add(fHistCentrality);
174
c92284d5 175 fHistZVertex = new TH1F("fHistZVertex","Z vertex position", 60, -30, 30);
176 fHistZVertex->GetXaxis()->SetTitle("z");
177 fHistZVertex->GetYaxis()->SetTitle("counts");
178 fOutput->Add(fHistZVertex);
be7b6e63 179
a5621834 180 fHistTracksCent = new TH2F("fHistTracksCent","Tracks vs. centrality", 100, 0, 100, fNbins, 0, 4000);
c3ba2d3d 181 fHistTracksCent->GetXaxis()->SetTitle("Centrality (%)");
182 fHistTracksCent->GetYaxis()->SetTitle("No. of tracks");
183 fOutput->Add(fHistTracksCent);
184
7c1d624c 185 fHistJetsCent = new TH2F("fHistJetsCent","Jets vs. centrality", 100, 0, 100, 60, 0, 60);
186 fHistJetsCent->GetXaxis()->SetTitle("Centrality (%)");
187 fHistJetsCent->GetYaxis()->SetTitle("No. of jets");
188 fOutput->Add(fHistJetsCent);
189
190 fHistJetsParts = new TH2F("fHistJetsParts","Jets vs. centrality", fNbins, 0, 6000, 60, 0, 60);
191 fHistJetsParts->GetXaxis()->SetTitle("No. of particles");
192 fHistJetsParts->GetYaxis()->SetTitle("No. of jets");
193 fOutput->Add(fHistJetsParts);
194
090a0c3e 195 if (fAnaType == kEMCAL || fAnaType == kEMCALOnly) {
a5621834 196 fHistClusCent = new TH2F("fHistClusCent","Clusters vs. centrality", 100, 0, 100, fNbins, 0, 2000);
6fd5039f 197 fHistClusCent->GetXaxis()->SetTitle("Centrality (%)");
198 fHistClusCent->GetYaxis()->SetTitle("No. of clusters");
199 fOutput->Add(fHistClusCent);
090a0c3e 200
201 fHistClusTracks = new TH2F("fHistClusTracks","Clusters vs. tracks", fNbins, 0, 4000, fNbins, 0, 2000);
202 fHistClusTracks->GetXaxis()->SetTitle("No. of tracks");
203 fHistClusTracks->GetYaxis()->SetTitle("No. of clusters");
204 fOutput->Add(fHistClusTracks);
205
be7b6e63 206 fHistCellsCent = new TH2F("fHistCellsCent","Cells vs. centrality", 100, 0, 100, fNbins, 0, 6000);
207 fHistCellsCent->GetXaxis()->SetTitle("Centrality (%)");
208 fHistCellsCent->GetYaxis()->SetTitle("No. of EMCal cells");
209 fOutput->Add(fHistCellsCent);
210
211 fHistCellsTracks = new TH2F("fHistCellsTracks","Cells vs. tracks", fNbins, 0, 4000, fNbins, 0, 6000);
212 fHistCellsTracks->GetXaxis()->SetTitle("No. of tracks");
213 fHistCellsTracks->GetYaxis()->SetTitle("No. of EMCal cells");
214 fOutput->Add(fHistCellsTracks);
215
781da0a3 216 fHistClusTimeEnergy = new TH2F("fHistClusTimeEnergy","Time vs. energy of clusters", fNbins, fMinBinPt, fMaxBinPt, fNbins, -1e-6, 1e-6);
090a0c3e 217 fHistClusTimeEnergy->GetXaxis()->SetTitle("Energy (GeV)");
218 fHistClusTimeEnergy->GetYaxis()->SetTitle("Time");
219 fOutput->Add(fHistClusTimeEnergy);
220
1f6fff78 221 if (fDoTrigger) {
222 fHistMaxL1FastORCent = new TH2F("fHistMaxL1FastORCent","fHistMaxL1ClusCent", 100, 0, 100, 250, 0, 250);
223 fHistMaxL1FastORCent->GetXaxis()->SetTitle("Centrality [%]");
224 fHistMaxL1FastORCent->GetYaxis()->SetTitle("Maximum L1 FastOR");
225 fOutput->Add(fHistMaxL1FastORCent);
226
227 fHistMaxL1ClusCent = new TH2F("fHistMaxL1ClusCent","fHistMaxL1ClusCent", 100, 0, 100, 250, 0, 250);
228 fHistMaxL1ClusCent->GetXaxis()->SetTitle("Centrality [%]");
229 fHistMaxL1ClusCent->GetYaxis()->SetTitle("Maximum L1 trigger cluster");
230 fOutput->Add(fHistMaxL1ClusCent);
231
232 fHistMaxL1ThrCent = new TH2F("fHistMaxL1ThrCent","fHistMaxL1ThrCent", 100, 0, 100, 250, 0, 250);
233 fHistMaxL1ThrCent->GetXaxis()->SetTitle("Centrality [%]");
234 fHistMaxL1ThrCent->GetYaxis()->SetTitle("Maximum L1 threshold");
235 fOutput->Add(fHistMaxL1ThrCent);
236 }
6fd5039f 237 }
238
239 fHistTracksPt = new TH1F("fHistTracksPt","p_{T} spectrum of reconstructed tracks", fNbins, fMinBinPt, fMaxBinPt);
240 fHistTracksPt->GetXaxis()->SetTitle("p_{T} [GeV/c]");
c3ba2d3d 241 fHistTracksPt->GetYaxis()->SetTitle("counts");
242 fOutput->Add(fHistTracksPt);
243
6fd5039f 244 fHistTrPhiEta = new TH2F("fHistTrPhiEta","Phi-Eta distribution of tracks", 80, -2, 2, 128, 0, 6.4);
c3ba2d3d 245 fHistTrPhiEta->GetXaxis()->SetTitle("#eta");
246 fHistTrPhiEta->GetYaxis()->SetTitle("#phi");
247 fOutput->Add(fHistTrPhiEta);
248
6fd5039f 249 fHistTrEmcPhiEta = new TH2F("fHistTrEmcPhiEta","Phi-Eta emcal distribution of tracks", 80, -2, 2, 128, 0, 6.4);
250 fHistTrEmcPhiEta->GetXaxis()->SetTitle("#eta");
251 fHistTrEmcPhiEta->GetYaxis()->SetTitle("#phi");
252 fOutput->Add(fHistTrEmcPhiEta);
253
079b4732 254 fHistTrPhiEtaNonProp = new TH2F("fHistTrPhiEtaNonProp","fHistTrPhiEtaNonProp", 80, -2, 2, 128, 0, 6.4);
255 fHistTrPhiEtaNonProp->GetXaxis()->SetTitle("#eta");
256 fHistTrPhiEtaNonProp->GetYaxis()->SetTitle("#phi");
257 fOutput->Add(fHistTrPhiEtaNonProp);
258
3e951e32 259 fHistDeltaEtaPt = new TH2F("fHistDeltaEtaPt","fHistDeltaEtaPt", fNbins, fMinBinPt, fMaxBinPt, 80, -0.5, 0.5);
260 fHistDeltaEtaPt->GetXaxis()->SetTitle("p_{T} [GeV/c]");
079b4732 261 fHistDeltaEtaPt->GetYaxis()->SetTitle("#delta#eta");
3e951e32 262 fOutput->Add(fHistDeltaEtaPt);
263
264 fHistDeltaPhiPt = new TH2F("fHistDeltaPhiPt","fHistDeltaPhiPt", fNbins, fMinBinPt, fMaxBinPt, 256, -1.6, 4.8);
265 fHistDeltaPhiPt->GetXaxis()->SetTitle("p_{T} [GeV/c]");
079b4732 266 fHistDeltaPhiPt->GetYaxis()->SetTitle("#delta#phi");
3e951e32 267 fOutput->Add(fHistDeltaPhiPt);
268
11d18b51 269 if (fRepropagateTracks) {
270 fHistDeltaEtaNewProp = new TH1F("fHistDeltaEtaNewProp","fHistDeltaEtaNewProp", 800, -0.4, 0.4);
271 fHistDeltaEtaNewProp->GetXaxis()->SetTitle("#delta#eta");
272 fHistDeltaEtaNewProp->GetYaxis()->SetTitle("counts");
273 fOutput->Add(fHistDeltaEtaNewProp);
274
275 fHistDeltaPhiNewProp = new TH1F("fHistDeltaPhiNewProp","fHistDeltaPhiNewProp", 2560, -1.6, 3.2);
276 fHistDeltaPhiNewProp->GetXaxis()->SetTitle("#delta#phi");
277 fHistDeltaPhiNewProp->GetYaxis()->SetTitle("counts");
278 fOutput->Add(fHistDeltaPhiNewProp);
279 }
079b4732 280
090a0c3e 281 if (fAnaType == kEMCAL || fAnaType == kEMCALOnly) {
7cf4626b 282 fHistClusPhiEtaEnergy = new TH3F("fHistClusPhiEtaEnergy","Phi-Eta-Energy distribution of clusters", fNbins, fMinBinPt, fMaxBinPt, 80, -2, 2, 128, 0, 6.4);
283 fHistClusPhiEtaEnergy->GetXaxis()->SetTitle("E [GeV]");
284 fHistClusPhiEtaEnergy->GetYaxis()->SetTitle("#eta");
285 fHistClusPhiEtaEnergy->GetZaxis()->SetTitle("#phi");
286 fOutput->Add(fHistClusPhiEtaEnergy);
287
288 fHistNCellsEnergy = new TH2F("fHistNCellsEnergy","Number of cells vs. energy of clusters", fNbins, fMinBinPt, fMaxBinPt, 30, 0, 30);
289 fHistNCellsEnergy->GetXaxis()->SetTitle("E [GeV]");
290 fHistNCellsEnergy->GetYaxis()->SetTitle("N_{cells}");
291 fOutput->Add(fHistNCellsEnergy);
6fd5039f 292 }
c3ba2d3d 293
090a0c3e 294 if (fAnaType == kEMCAL || fAnaType == kEMCALOnly) {
6fd5039f 295
296 fHistCellsEnergy = new TH1F("fHistCellsEnergy","Energy spectrum of cells", fNbins, fMinBinPt, fMaxBinPt);
297 fHistCellsEnergy->GetXaxis()->SetTitle("E [GeV]");
298 fHistCellsEnergy->GetYaxis()->SetTitle("counts");
299 fOutput->Add(fHistCellsEnergy);
300
98750b70 301 fHistChVSneCells = new TH2F("fHistChVSneCells","Charged energy vs. neutral (cells) energy",
be7b6e63 302 (Int_t)(fNbins * 2.5), fMinBinPt, fMaxBinPt * 2.5, (Int_t)(fNbins * 2.5), fMinBinPt, fMaxBinPt * 2.5);
6fd5039f 303 fHistChVSneCells->GetXaxis()->SetTitle("E [GeV]");
304 fHistChVSneCells->GetYaxis()->SetTitle("P [GeV/c]");
305 fOutput->Add(fHistChVSneCells);
306
98750b70 307 fHistChVSneClus = new TH2F("fHistChVSneClus","Charged energy vs. neutral (clusters) energy",
be7b6e63 308 (Int_t)(fNbins * 2.5), fMinBinPt, fMaxBinPt * 2.5, (Int_t)(fNbins * 2.5), fMinBinPt, fMaxBinPt * 2.5);
6fd5039f 309 fHistChVSneClus->GetXaxis()->SetTitle("E [GeV]");
310 fHistChVSneClus->GetYaxis()->SetTitle("P [GeV/c]");
311 fOutput->Add(fHistChVSneClus);
312
98750b70 313 fHistChVSneCorrCells = new TH2F("fHistChVSneCorrCells","Charged energy vs. neutral (corrected cells) energy",
be7b6e63 314 (Int_t)(fNbins * 2.5), fMinBinPt, fMaxBinPt * 2.5, (Int_t)(fNbins * 2.5), fMinBinPt , fMaxBinPt * 2.5);
6fd5039f 315 fHistChVSneCorrCells->GetXaxis()->SetTitle("E [GeV]");
316 fHistChVSneCorrCells->GetYaxis()->SetTitle("P [GeV/c]");
317 fOutput->Add(fHistChVSneCorrCells);
c3ba2d3d 318 }
c3ba2d3d 319
320 for (Int_t i = 0; i < 5; i++) {
321 TString histnamephi("fHistTrackPhi_");
322 histnamephi += i;
323 fHistTrackPhi[i] = new TH1F(histnamephi.Data(),histnamephi.Data(), 128, 0, 6.4);
324 fHistTrackPhi[i]->GetXaxis()->SetTitle("Phi");
325 fOutput->Add(fHistTrackPhi[i]);
326
327 TString histnameeta("fHistTrackEta_");
328 histnameeta += i;
329 fHistTrackEta[i] = new TH1F(histnameeta.Data(),histnameeta.Data(), 100, -2, 2);
330 fHistTrackEta[i]->GetXaxis()->SetTitle("Eta");
331 fOutput->Add(fHistTrackEta[i]);
332 }
aa4d701c 333
334 for (Int_t i = 0; i < 6; i++) {
335 TString histnamephipt("fHistTrackPhiPt_");
336 histnamephipt += i;
337 fHistTrackPhiPt[i] = new TH1F(histnamephipt.Data(),histnamephipt.Data(), 128, 0, 6.4);
338 fHistTrackPhiPt[i]->GetXaxis()->SetTitle("Phi");
339 fOutput->Add(fHistTrackPhiPt[i]);
340
341 TString histnameetapt("fHistTrackEtaPt_");
342 histnameetapt += i;
343 fHistTrackEtaPt[i] = new TH1F(histnameetapt.Data(),histnameetapt.Data(), 100, -2, 2);
344 fHistTrackEtaPt[i]->GetXaxis()->SetTitle("Eta");
345 fOutput->Add(fHistTrackEtaPt[i]);
346 }
c3ba2d3d 347
348 fHistTrackPhi[0]->SetLineColor(kRed);
349 fHistTrackEta[0]->SetLineColor(kRed);
350 fHistTrackPhi[1]->SetLineColor(kBlue);
351 fHistTrackEta[1]->SetLineColor(kBlue);
352 fHistTrackPhi[2]->SetLineColor(kGreen);
353 fHistTrackEta[2]->SetLineColor(kGreen);
354 fHistTrackPhi[3]->SetLineColor(kOrange);
355 fHistTrackEta[3]->SetLineColor(kOrange);
356 fHistTrackPhi[4]->SetLineColor(kBlack);
357 fHistTrackEta[4]->SetLineColor(kBlack);
358
e44e8726 359 if (!fJetsName.IsNull()) {
6fd5039f 360
e44e8726 361 TString histname;
362
363 for (Int_t i = 0; i < 4; i++) {
58285fc6 364 histname = "fHistJetsPhiEta_";
3000c095 365 histname += i;
58285fc6 366 fHistJetsPhiEta[i] = new TH2F(histname.Data(), histname.Data(), 80, -2, 2, 128, 0, 6.4);
367 fHistJetsPhiEta[i]->GetXaxis()->SetTitle("#eta");
368 fHistJetsPhiEta[i]->GetYaxis()->SetTitle("#phi");
369 fHistJetsPhiEta[i]->GetZaxis()->SetTitle("p_{T} [GeV/c]");
370 fOutput->Add(fHistJetsPhiEta[i]);
3000c095 371
e44e8726 372 histname = "fHistJetsPtNonBias_";
98750b70 373 histname += i;
be7b6e63 374 fHistJetsPtNonBias[i] = new TH1F(histname.Data(), histname.Data(), (Int_t)(fNbins * 2.5), fMinBinPt, fMaxBinPt * 2.5);
e44e8726 375 fHistJetsPtNonBias[i]->GetXaxis()->SetTitle("p_{T} [GeV/c]");
376 fHistJetsPtNonBias[i]->GetYaxis()->SetTitle("counts");
377 fOutput->Add(fHistJetsPtNonBias[i]);
378
379 histname = "fHistJetsPtTrack_";
380 histname += i;
be7b6e63 381 fHistJetsPtTrack[i] = new TH1F(histname.Data(), histname.Data(), (Int_t)(fNbins * 2.5), fMinBinPt, fMaxBinPt * 2.5);
e44e8726 382 fHistJetsPtTrack[i]->GetXaxis()->SetTitle("p_{T} [GeV/c]");
383 fHistJetsPtTrack[i]->GetYaxis()->SetTitle("counts");
384 fOutput->Add(fHistJetsPtTrack[i]);
385
090a0c3e 386 if (fAnaType == kEMCAL || fAnaType == kEMCALOnly) {
e44e8726 387 histname = "fHistJetsPtClus_";
388 histname += i;
be7b6e63 389 fHistJetsPtClus[i] = new TH1F(histname.Data(), histname.Data(), (Int_t)(fNbins * 2.5), fMinBinPt, fMaxBinPt * 2.5);
e44e8726 390 fHistJetsPtClus[i]->GetXaxis()->SetTitle("p_{T} [GeV/c]");
391 fHistJetsPtClus[i]->GetYaxis()->SetTitle("counts");
392 fOutput->Add(fHistJetsPtClus[i]);
393 }
394
395 histname = "fHistJetsPt_";
396 histname += i;
be7b6e63 397 fHistJetsPt[i] = new TH1F(histname.Data(), histname.Data(), (Int_t)(fNbins * 2.5), fMinBinPt, fMaxBinPt * 2.5);
e44e8726 398 fHistJetsPt[i]->GetXaxis()->SetTitle("p_{T} [GeV/c]");
399 fHistJetsPt[i]->GetYaxis()->SetTitle("counts");
400 fOutput->Add(fHistJetsPt[i]);
99cfd012 401
e44e8726 402 histname = "fHistJetsPtAreaNonBias_";
403 histname += i;
be7b6e63 404 fHistJetsPtAreaNonBias[i] = new TH2F(histname.Data(), histname.Data(), (Int_t)(fNbins * 2.5), fMinBinPt, fMaxBinPt * 2.5, 20, 0, fJetRadius * fJetRadius * TMath::Pi() * 1.5);
e44e8726 405 fHistJetsPtAreaNonBias[i]->GetXaxis()->SetTitle("p_{T} [GeV/c]");
406 fHistJetsPtAreaNonBias[i]->GetYaxis()->SetTitle("area");
407 fOutput->Add(fHistJetsPtAreaNonBias[i]);
408
409 histname = "fHistJetsPtArea_";
410 histname += i;
be7b6e63 411 fHistJetsPtArea[i] = new TH2F(histname.Data(), histname.Data(), (Int_t)(fNbins * 2.5), fMinBinPt, fMaxBinPt * 2.5, 20, 0, fJetRadius * fJetRadius * TMath::Pi() * 1.5);
e44e8726 412 fHistJetsPtArea[i]->GetXaxis()->SetTitle("p_{T} [GeV/c]");
413 fHistJetsPtArea[i]->GetYaxis()->SetTitle("area");
414 fOutput->Add(fHistJetsPtArea[i]);
415 }
6fd5039f 416 }
417
c3ba2d3d 418 PostData(1, fOutput); // Post data for ALL output slots >0 here, to get at least an empty histogram
419}
420
421//________________________________________________________________________
6fd5039f 422Bool_t AliAnalysisTaskSAQA::RetrieveEventObjects()
c3ba2d3d 423{
c3604199 424 // Retrieve event objects.
425
426 if (!AliAnalysisTaskEmcalJet::RetrieveEventObjects())
6fd5039f 427 return kFALSE;
c3ba2d3d 428
7c1d624c 429 fNclusters = 0;
430 fNtracks = 0;
431 fNjets = 0;
432
e44e8726 433 if (!fTrgClusName.IsNull() && fDoTrigger && !fTrgClusters) {
c3ba2d3d 434 fTrgClusters = dynamic_cast<TClonesArray*>(InputEvent()->FindListObject(fTrgClusName));
435 if (!fTrgClusters) {
e44e8726 436 AliError(Form("%s: Could not retrieve trigger clusters %s!", GetName(), fTrgClusName.Data()));
437 return kFALSE;
438 }
439 else {
440 TClass *cl = fTrgClusters->GetClass();
441 if (!cl->GetBaseClass("AliVCluster") && !cl->GetBaseClass("AliEmcalParticle")) {
442 AliError(Form("%s: Collection %s does not contain AliVCluster nor AliEmcalParticle objects!", GetName(), fTrgClusName.Data()));
443 fTrgClusters = 0;
444 return kFALSE;
445 }
c3ba2d3d 446 }
447 }
6fd5039f 448
449 return kTRUE;
c3ba2d3d 450}
451
7c1d624c 452
c3ba2d3d 453//________________________________________________________________________
6fd5039f 454Bool_t AliAnalysisTaskSAQA::FillHistograms()
c3ba2d3d 455{
c3604199 456 // Fill histograms.
457
c3ba2d3d 458 fHistCentrality->Fill(fCent);
c92284d5 459 fHistZVertex->Fill(fVertex[2]);
090a0c3e 460
c3ba2d3d 461 Float_t trackSum = DoTrackLoop();
462
c3ba2d3d 463 DoJetLoop();
464
7c1d624c 465 fHistTracksCent->Fill(fCent, fNtracks);
be7b6e63 466
7c1d624c 467 if (fAnaType == kEMCAL || fAnaType == kEMCALOnly) {
be7b6e63 468
6fd5039f 469 Float_t clusSum = DoClusterLoop();
c3ba2d3d 470
7c1d624c 471 fHistClusCent->Fill(fCent, fNclusters);
472 fHistClusTracks->Fill(fNtracks, fNclusters);
473
6fd5039f 474 Float_t cellSum = 0, cellCutSum = 0;
475
be7b6e63 476 Int_t ncells = DoCellLoop(cellSum, cellCutSum);
477
478 if (fTracks)
479 fHistCellsTracks->Fill(fTracks->GetEntriesFast(), ncells);
480
481 fHistCellsCent->Fill(fCent, ncells);
6fd5039f 482
483 fHistChVSneCells->Fill(cellSum, trackSum);
484 fHistChVSneClus->Fill(clusSum, trackSum);
485 fHistChVSneCorrCells->Fill(cellCutSum, trackSum);
c3ba2d3d 486
1f6fff78 487 if (fDoTrigger) {
488 Float_t maxTrgClus = DoTriggerClusLoop();
489 fHistMaxL1ClusCent->Fill(fCent, maxTrgClus);
6fd5039f 490
1f6fff78 491 Int_t maxL1amp = -1;
492 Int_t maxL1thr = -1;
6fd5039f 493
1f6fff78 494 DoTriggerPrimitives(maxL1amp, maxL1thr);
495
496 if (maxL1amp > -1)
497 fHistMaxL1FastORCent->Fill(fCent, maxL1amp);
498
499 if (maxL1thr > -1)
500 fHistMaxL1ThrCent->Fill(fCent, maxL1thr);
501 }
6fd5039f 502 }
c3ba2d3d 503
7c1d624c 504 fHistJetsCent->Fill(fCent, fNjets);
505 fHistJetsParts->Fill(fNtracks + fNclusters, fNjets);
506
6fd5039f 507 return kTRUE;
c3ba2d3d 508}
509
510//________________________________________________________________________
be7b6e63 511Int_t AliAnalysisTaskSAQA::DoCellLoop(Float_t &sum, Float_t &sum_cut)
c3ba2d3d 512{
c3604199 513 // Do cell loop.
514
c3ba2d3d 515 AliVCaloCells *cells = InputEvent()->GetEMCALCells();
516
517 if (!cells)
be7b6e63 518 return 0;
c3ba2d3d 519
c3604199 520 const Int_t ncells = cells->GetNumberOfCells();
c3ba2d3d 521
522 for (Int_t pos = 0; pos < ncells; pos++) {
c3ba2d3d 523 Float_t amp = cells->GetAmplitude(pos);
c3ba2d3d 524 fHistCellsEnergy->Fill(amp);
c3ba2d3d 525 sum += amp;
c3ba2d3d 526 if (amp < fCellEnergyCut)
527 continue;
c3ba2d3d 528 sum_cut += amp;
c3ba2d3d 529 }
be7b6e63 530
531 return ncells;
c3ba2d3d 532}
533
534//________________________________________________________________________
535Float_t AliAnalysisTaskSAQA::DoClusterLoop()
536{
c3604199 537 // Do cluster loop.
538
c3ba2d3d 539 if (!fCaloClusters)
540 return 0;
541
542 Float_t sum = 0;
543
544 // Cluster loop
545 Int_t nclusters = fCaloClusters->GetEntriesFast();
546
547 for (Int_t iClusters = 0; iClusters < nclusters; iClusters++) {
e44e8726 548 AliVCluster* cluster = static_cast<AliVCluster*>(fCaloClusters->At(iClusters));
c3ba2d3d 549 if (!cluster) {
550 AliError(Form("Could not receive cluster %d", iClusters));
551 continue;
552 }
c3ba2d3d 553
e44e8726 554 if (!AcceptCluster(cluster, kTRUE))
555 continue;
c3ba2d3d 556
c3ba2d3d 557 sum += cluster->E();
558
7cf4626b 559 TLorentzVector nPart;
560 cluster->GetMomentum(nPart, fVertex);
c3ba2d3d 561
7cf4626b 562 fHistClusPhiEtaEnergy->Fill(cluster->E(), nPart.Eta(), nPart.Phi());
563 fHistNCellsEnergy->Fill(cluster->E(), cluster->GetNCells());
090a0c3e 564
565 fHistClusTimeEnergy->Fill(cluster->E(), cluster->GetTOF());
7c1d624c 566
567 fNclusters++;
c3604199 568 }
c3ba2d3d 569
570 return sum;
571}
572
573//________________________________________________________________________
574Float_t AliAnalysisTaskSAQA::DoTrackLoop()
575{
c3604199 576 // Do track loop.
577
c3ba2d3d 578 if (!fTracks)
579 return 0;
580
581 Float_t sum = 0;
582
c3ba2d3d 583 Int_t ntracks = fTracks->GetEntriesFast();
584 Int_t nclusters = 0;
585 if (fCaloClusters)
586 nclusters = fCaloClusters->GetEntriesFast();
587
c3604199 588 for (Int_t i = 0; i < ntracks; i++) {
c3ba2d3d 589
e44e8726 590 AliVParticle* track = static_cast<AliVParticle*>(fTracks->At(i)); // pointer to reconstructed to track
1f6fff78 591
c3604199 592 if (!track) {
c3ba2d3d 593 AliError(Form("Could not retrieve track %d",i));
594 continue;
595 }
c3ba2d3d 596
1f6fff78 597 AliVTrack* vtrack = dynamic_cast<AliVTrack*>(track);
598
599 if (vtrack && !AcceptTrack(vtrack, kTRUE))
600 continue;
7c1d624c 601
602 fNtracks++;
1f6fff78 603
c3ba2d3d 604 fHistTracksPt->Fill(track->Pt());
605
606 sum += track->P();
607
608 Int_t label = track->GetLabel();
609
610 fHistTrPhiEta->Fill(track->Eta(), track->Phi());
611
612 fHistTrackEta[4]->Fill(track->Eta());
613 fHistTrackPhi[4]->Fill(track->Phi());
614
615 if (label >= 0 && label < 4) {
616 fHistTrackEta[label]->Fill(track->Eta());
617 fHistTrackPhi[label]->Fill(track->Phi());
618 }
619
aa4d701c 620 if (track->Pt() < 0.5) {
621 fHistTrackPhiPt[0]->Fill(track->Phi());
622 fHistTrackEtaPt[0]->Fill(track->Eta());
623 }
624 else if (track->Pt() < 1) {
625 fHistTrackPhiPt[1]->Fill(track->Phi());
626 fHistTrackEtaPt[1]->Fill(track->Eta());
627 }
628 else if (track->Pt() < 2) {
629 fHistTrackPhiPt[2]->Fill(track->Phi());
630 fHistTrackEtaPt[2]->Fill(track->Eta());
631 }
632 else if (track->Pt() < 3) {
633 fHistTrackPhiPt[3]->Fill(track->Phi());
634 fHistTrackEtaPt[3]->Fill(track->Eta());
635 }
636 else if (track->Pt() < 5) {
637 fHistTrackPhiPt[4]->Fill(track->Phi());
638 fHistTrackEtaPt[4]->Fill(track->Eta());
639 }
640 else {
641 fHistTrackPhiPt[5]->Fill(track->Phi());
642 fHistTrackEtaPt[5]->Fill(track->Eta());
643 }
644
1f6fff78 645 if (!vtrack)
c3ba2d3d 646 continue;
647
079b4732 648 if (vtrack->GetTrackEtaOnEMCal() == -999 || vtrack->GetTrackPhiOnEMCal() == -999)
649 fHistTrPhiEtaNonProp->Fill(vtrack->Eta(), vtrack->Phi());
650
1f6fff78 651 fHistTrEmcPhiEta->Fill(vtrack->GetTrackEtaOnEMCal(), vtrack->GetTrackPhiOnEMCal());
3e951e32 652 fHistDeltaEtaPt->Fill(vtrack->Pt(), vtrack->Eta() - vtrack->GetTrackEtaOnEMCal());
653 fHistDeltaPhiPt->Fill(vtrack->Pt(), vtrack->Phi() - vtrack->GetTrackPhiOnEMCal());
c92284d5 654
11d18b51 655 if (fRepropagateTracks && vtrack->GetTrackEtaOnEMCal() > -2) {
c92284d5 656 Float_t propeta = -999, propphi = -999;
657 PropagateTrack(vtrack, propeta, propphi);
658 fHistDeltaEtaNewProp->Fill(propeta - vtrack->GetTrackEtaOnEMCal());
659 fHistDeltaPhiNewProp->Fill(propphi - vtrack->GetTrackPhiOnEMCal());
660 }
c3ba2d3d 661 }
662
663 return sum;
664}
665
079b4732 666//____________________________________________________________________________
667void AliAnalysisTaskSAQA::PropagateTrack(AliVTrack *track, Float_t &eta, Float_t &phi)
668{
669 eta = -999;
670 phi = -999;
671
672 if (!track)
673 return;
674
c92284d5 675 if (track->Pt() == 0)
676 return;
677
079b4732 678 // init the magnetic field if not already on
679 if(!TGeoGlobalMagField::Instance()->GetField()) {
680 AliInfo("Init the magnetic field\n");
681 AliAODEvent* aodevent = dynamic_cast<AliAODEvent*>(InputEvent());
682 if (aodevent) {
683 Double_t curSol = 30000*aodevent->GetMagneticField()/5.00668;
684 Double_t curDip = 6000 *aodevent->GetMuonMagFieldScale();
685 AliMagF *field = AliMagF::CreateFieldMap(curSol,curDip);
686 TGeoGlobalMagField::Instance()->SetField(field);
687 }
688 }
689
690 Double_t cv[21];
691 for (Int_t i = 0; i < 21; i++) cv[i] = 0;
692
693 Double_t pos[3], mom[3];
694 track->GetXYZ(pos);
695 track->GetPxPyPz(mom);
696 AliExternalTrackParam *trackParam = new AliExternalTrackParam(pos, mom, cv, track->Charge());
697
c92284d5 698 if(!AliTrackerBase::PropagateTrackToBxByBz(trackParam, 430., 0.139, 20, kTRUE, 0.8, -1)) return;
079b4732 699 Double_t trkPos[3] = {0., 0., 0.};
700 if(!trackParam->GetXYZ(trkPos)) return;
701 TVector3 trkPosVec(trkPos[0], trkPos[1], trkPos[2]);
702 eta = trkPosVec.Eta();
703 phi = trkPosVec.Phi();
704 if(phi < 0)
705 phi += 2 * TMath::Pi();
706}
707
c3ba2d3d 708//________________________________________________________________________
709void AliAnalysisTaskSAQA::DoJetLoop()
710{
c3604199 711 // Do jet loop.
712
c3ba2d3d 713 if (!fJets)
714 return;
715
716 Int_t njets = fJets->GetEntriesFast();
717
718 for (Int_t ij = 0; ij < njets; ij++) {
719
e44e8726 720 AliEmcalJet* jet = static_cast<AliEmcalJet*>(fJets->At(ij));
c3ba2d3d 721
722 if (!jet) {
723 AliError(Form("Could not receive jet %d", ij));
724 continue;
725 }
726
6fd5039f 727 if (!AcceptJet(jet, kFALSE))
c3ba2d3d 728 continue;
6fd5039f 729
99cfd012 730 fHistJetsPtNonBias[fCentBin]->Fill(jet->Pt());
731 fHistJetsPtAreaNonBias[fCentBin]->Fill(jet->Pt(), jet->Area());
6fd5039f 732
7c1d624c 733 fNjets++;
734
6fd5039f 735 if (jet->MaxTrackPt() > fPtBiasJetTrack)
736 fHistJetsPtTrack[fCentBin]->Fill(jet->Pt());
737
738 if (fAnaType == kEMCAL && jet->MaxClusterPt() > fPtBiasJetClus)
739 fHistJetsPtClus[fCentBin]->Fill(jet->Pt());
c3ba2d3d 740
99cfd012 741 if (!AcceptBiasJet(jet))
742 continue;
743
744 fHistJetsPt[fCentBin]->Fill(jet->Pt());
745 fHistJetsPtArea[fCentBin]->Fill(jet->Pt(), jet->Area());
746
58285fc6 747 fHistJetsPhiEta[fCentBin]->Fill(jet->Eta(), jet->Phi());
c3ba2d3d 748 }
749}
750
751//________________________________________________________________________
752Float_t AliAnalysisTaskSAQA::DoTriggerClusLoop()
753{
c3604199 754 // Do trigger cluster loop.
755
c3ba2d3d 756 if (!fTrgClusters)
757 return 0;
758
759 Int_t ntrgclusters = fTrgClusters->GetEntriesFast();
760 Float_t maxTrgClus = 0;
761
762 for (Int_t iClusters = 0; iClusters < ntrgclusters; iClusters++) {
e44e8726 763 AliVCluster* cluster = static_cast<AliVCluster*>(fTrgClusters->At(iClusters));
c3ba2d3d 764 if (!cluster) {
765 AliError(Form("Could not receive cluster %d", iClusters));
766 continue;
767 }
768
769 if (!(cluster->IsEMCAL())) continue;
770
771 if (cluster->E() > maxTrgClus)
772 maxTrgClus = cluster->E();
773
774 }
775 return maxTrgClus;
776}
777
778//________________________________________________________________________
779void AliAnalysisTaskSAQA::DoTriggerPrimitives(Int_t &maxL1amp, Int_t &maxL1thr)
780{
c3604199 781 // Do trigger primitives loop.
782
c3ba2d3d 783 AliVCaloTrigger *triggers = InputEvent()->GetCaloTrigger("EMCAL");
784
785 if (!triggers || triggers->GetEntries() == 0)
786 return;
787
788 triggers->Reset();
789 Int_t L1amp = 0;
790 Int_t L1thr = 0;
791 maxL1amp = -1;
792 maxL1thr = -1;
793
794 while (triggers->Next()) {
795 triggers->GetL1TimeSum(L1amp);
796 if (maxL1amp < L1amp)
797 maxL1amp = L1amp;
798
799 triggers->GetL1Threshold(L1thr);
800 if (maxL1thr < L1thr)
801 maxL1thr = L1thr;
802 }
803}
804
805//________________________________________________________________________
806void AliAnalysisTaskSAQA::Terminate(Option_t *)
807{
808 // Called once at the end of the analysis.
809}