Adding missing include file
[u/mrichter/AliRoot.git] / PWG1 / VZERO / AliAnaVZEROPbPb.cxx
CommitLineData
a148e064 1#include "TChain.h"
2#include "TTree.h"
3#include "TH1F.h"
4#include "TH2F.h"
5
6#include "AliAnalysisTask.h"
7#include "AliAnalysisManager.h"
8
9#include "AliESDEvent.h"
10#include "AliESDInputHandler.h"
11#include "AliAnaVZEROPbPb.h"
12#include "AliMultiplicity.h"
13#include "AliESDtrackCuts.h"
14#include "AliESDUtils.h"
15#include "AliCentrality.h"
16
17// VZERO includes
18#include "AliESDVZERO.h"
19#include "AliESDVZEROfriend.h"
20
21ClassImp(AliAnaVZEROPbPb)
22
23
24AliAnaVZEROPbPb::AliAnaVZEROPbPb()
25 : AliAnalysisTaskSE("AliAnaVZEROPbPb"), fESD(0), fEsdV0(0), fOutputList(0), fClassesNames(0),
26 fNFlags(0),
27 fhAdcPMTNoTime(0),
28 fhAdcPMTWithTime(0),
29 fhTimePMT(0),
30 fhWidthPMT(0),
31 fhTimeCorr(0),
32 fhPmtMult(0),
33 fhV0ampl(0),
34 fhEvents(0),
35 fhVtxXYBB(0),
36 fhVtxZBB(0),
37 fhVtxXYBGA(0),
38 fhVtxZBGA(0),
39 fhVtxXYBGC(0),
40 fhVtxZBGC(0),
41 fhL2Triggers(0),
42fhOnlineCharge(0),
43fhRecoMult(0),
44fhV0vsSPDCentrality(0),
45fhTriggerBits(0),
46fhTotRecoMult(0),
47fhCentrality(0)
48{
49 // Constructor
50 for(Int_t i = 0; i < 2; ++i) {
51 fhAdcNoTime[i] = fhAdcWithTime[i] = fhWidth[i] = fhTime[i] = NULL;
52 fhAdcTime[i] = fhAdcWidth[i] = NULL;
53 }
54
55 // Define input and output slots here
56 // Input slot #0 works with a TChain
57 DefineInput(0, TChain::Class());
58 // Output slot #0 id reserved by the base class for AOD
59 // Output slot #1 writes into a TH1 container
60 DefineOutput(1, TList::Class());
61}
62
63//________________________________________________________________________
64AliAnaVZEROPbPb::AliAnaVZEROPbPb(const char *name)
65 : AliAnalysisTaskSE(name), fESD(0), fEsdV0(0), fOutputList(0),fNClasses(0),fClassesNames(0),
66 fNFlags(0),
67 fhAdcPMTNoTime(0),
68 fhAdcPMTWithTime(0),
69 fhTimePMT(0),
70 fhWidthPMT(0),
71 fhTimeCorr(0),
72 fhPmtMult(0),
73 fhV0ampl(0),
74 fhEvents(0),
75 fhVtxXYBB(0),
76 fhVtxZBB(0),
77 fhVtxXYBGA(0),
78 fhVtxZBGA(0),
79 fhVtxXYBGC(0),
80 fhVtxZBGC(0),
81 fhL2Triggers(0),
82fhOnlineCharge(0),
83fhRecoMult(0),
84fhV0vsSPDCentrality(0),
85fhTriggerBits(0),
86fhTotRecoMult(0),
87fhCentrality(0)
88{
89 // Constructor
90 for(Int_t i = 0; i < 2; ++i) {
91 fhAdcNoTime[i] = fhAdcWithTime[i] = fhWidth[i] = fhTime[i] = NULL;
92 fhAdcTime[i] = fhAdcWidth[i] = NULL;
93 }
94 // Define input and output slots here
95 // Input slot #0 works with a TChain
96 DefineInput(0, TChain::Class());
97 // Output slot #0 id reserved by the base class for AOD
98 // Output slot #1 writes into a TH1 container
99 DefineOutput(1, TList::Class());
100}
101//________________________________________________________________________
102void AliAnaVZEROPbPb::SetClassesNames(const Char_t * nameList){
103 TString names("AllClasses,");
104 names += nameList;
105 fClassesNames = names.Tokenize(",");
106 fNClasses = fClassesNames->GetEntriesFast();
107}
108//________________________________________________________________________
109TH1F * AliAnaVZEROPbPb::CreateHisto1D(const char* name, const char* title,Int_t nBins,
110 Double_t xMin, Double_t xMax,
111 const char* xLabel, const char* yLabel)
112{
113 // create a histogram
114 TH1F* result = new TH1F(name, title, nBins, xMin, xMax);
115 result->SetOption("hist");
116 if (xLabel) result->GetXaxis()->SetTitle(xLabel);
117 if (yLabel) result->GetYaxis()->SetTitle(yLabel);
118 result->SetMarkerStyle(kFullCircle);
119 return result;
120}
121
122//________________________________________________________________________
123TH2F * AliAnaVZEROPbPb::CreateHisto2D(const char* name, const char* title,Int_t nBinsX,
124 Double_t xMin, Double_t xMax,
125 Int_t nBinsY,
126 Double_t yMin, Double_t yMax,
127 const char* xLabel, const char* yLabel)
128{
129 // create a histogram
130 TH2F* result = new TH2F(name, title, nBinsX, xMin, xMax, nBinsY, yMin, yMax);
131 if (xLabel) result->GetXaxis()->SetTitle(xLabel);
132 if (yLabel) result->GetYaxis()->SetTitle(yLabel);
133 return result;
134}
135
136
137//________________________________________________________________________
138void AliAnaVZEROPbPb::UserCreateOutputObjects()
139{
140 // Create histograms
141 // Called once
142 if(fNClasses==0) {
143 AliFatal("No Classes Defined");
144 return;
145 }
146
147 fOutputList = new TList();
148 fOutputList->SetOwner(kTRUE);
149
150
151 CreateQAHistos();
152 CreateHistosPerL2Trigger();
153
154 PostData(1, fOutputList);
155 }
156//________________________________________________________________________
157void AliAnaVZEROPbPb::Init()
158{
159
160}
161
162//________________________________________________________________________
163void AliAnaVZEROPbPb::UserExec(Option_t *)
164{
165 // Main loop
166 // Called for each event
167
168 fESD = dynamic_cast<AliESDEvent*>(InputEvent());
169 if (!fESD) {
170 printf("ERROR: fESD not available\n");
171 return;
172 }
173
174 fEsdV0 = fESD->GetVZEROData();
175 if (!fEsdV0) {
176 Printf("ERROR: esd V0 not available");
177 return;
178 }
179
180 FillQAHistos();
181
182 FillPerL2TriggerHistos();
183
184
185
186 PostData(1, fOutputList);
187}
188//________________________________________________________________________
189void AliAnaVZEROPbPb::CreateHistosPerL2Trigger(){
190
191 fhOnlineCharge= new TH2F*[fNClasses];
192 fhCentrality= new TH1F*[fNClasses];
193 fhV0vsSPDCentrality= new TH2F*[fNClasses];
194 fhRecoMult= new TH2F*[fNClasses];
195 fhTotRecoMult= new TH1F*[fNClasses];
196 fhTriggerBits= new TH1F*[fNClasses];
197
198 fhL2Triggers = CreateHisto1D("hL2Triggers","L2 Triggers",fNClasses,0,fNClasses);
199 fOutputList->Add(fhL2Triggers);
200
201
202 TIter iter(fClassesNames);
203 TObjString* name;
204 Int_t iClass =0;
205 while((name = (TObjString*) iter.Next())){
206 fhL2Triggers->GetXaxis()->SetBinLabel(iClass+1,name->String().Data());
207
208 fhOnlineCharge[iClass] = CreateHisto2D(Form("hOnlineCharge_%s",name->String().Data()),Form("Online Charge for %s",name->String().Data()),100,0.,32000.,100,0.,32000,"V0A","V0C");
209 fOutputList->Add(fhOnlineCharge[iClass]);
210
211 fhCentrality[iClass] = CreateHisto1D(Form("hV0Centrality_%s",name->String().Data()),Form("V0 centrality for %s",name->String().Data()),100,0.,100.);
212 fOutputList->Add(fhCentrality[iClass]);
213
214 fhV0vsSPDCentrality[iClass] = CreateHisto2D(Form("hV0vsSPDCentrality_%s",name->String().Data()),Form("Centrality for %s",name->String().Data()),100,0.,100.,100,0.,100,"SPD Centrality (CL1)","V0 Centrality (V0M)");
215 fOutputList->Add(fhV0vsSPDCentrality[iClass]);
216
217 fhRecoMult[iClass] = CreateHisto2D(Form("hRecoMult_%s",name->String().Data()),Form("Reco Multiplicity for %s",name->String().Data()),100,0.,10000.,100,0.,20000,"V0A Mult","V0C Mult");
218 fOutputList->Add(fhRecoMult[iClass]);
219
220 fhTotRecoMult[iClass] = CreateHisto1D(Form("hTotRecoMult_%s",name->String().Data()),Form("Total Reco Multiplicity for %s",name->String().Data()),200,0.,20000.,"V0A + V0C Mult");
221 fOutputList->Add(fhTotRecoMult[iClass]);
222
223 fhTriggerBits[iClass] = CreateHisto1D(Form("hTriggerBits_%s",name->String().Data()),Form("Trigger Bits for %s",name->String().Data()),16,-0.5,15.5);
224 fhTriggerBits[iClass]->GetXaxis()->SetBinLabel(1,"BBA_AND_BBC");
225 fhTriggerBits[iClass]->GetXaxis()->SetBinLabel(2,"BBA_OR_BBC");
226 fhTriggerBits[iClass]->GetXaxis()->SetBinLabel(3,"BGA_AND_BBC");
227 fhTriggerBits[iClass]->GetXaxis()->SetBinLabel(4,"BGA");
228 fhTriggerBits[iClass]->GetXaxis()->SetBinLabel(5,"BBA_AND_BGC");
229 fhTriggerBits[iClass]->GetXaxis()->SetBinLabel(6,"BGC");
230 fhTriggerBits[iClass]->GetXaxis()->SetBinLabel(7,"CTA1_AND_CTC1");
231 fhTriggerBits[iClass]->GetXaxis()->SetBinLabel(8,"CTA1_OR_CTC1");
232 fhTriggerBits[iClass]->GetXaxis()->SetBinLabel(9,"CTA2_AND_CTC2");
233 fhTriggerBits[iClass]->GetXaxis()->SetBinLabel(10,"CTA1_OR_CTC1");
234 fhTriggerBits[iClass]->GetXaxis()->SetBinLabel(11,"MTA_AND_MTC");
235 fhTriggerBits[iClass]->GetXaxis()->SetBinLabel(12,"MTA_OR_MTC");
236 fhTriggerBits[iClass]->GetXaxis()->SetBinLabel(13,"BBA");
237 fhTriggerBits[iClass]->GetXaxis()->SetBinLabel(14,"BBC");
238 fhTriggerBits[iClass]->GetXaxis()->SetBinLabel(15,"BGA_OR_BGC");
239 fhTriggerBits[iClass]->GetXaxis()->SetBinLabel(16,"All True BG");
240 fOutputList->Add(fhTriggerBits[iClass]);
241
242 iClass++;
243 }
244
245}
246//________________________________________________________________________
247void AliAnaVZEROPbPb::CreateQAHistos(){
248
249 fNFlags = CreateHisto2D("hNFlags","BB Flags",33,-0.5,32.5,33,-0.5,32.5,"V0A","V0C");
250 fOutputList->Add(fNFlags);
251
252 for(int iSide = 0; iSide < 2; ++iSide){
253 TString side;
254 if(iSide) side = "V0A";
255 else side = "V0C";
256
257 fhAdcNoTime[iSide] = CreateHisto1D(Form("hAdcNoTime%s",side.Data()),Form("ADC (no Leading Time) %s",side.Data()),200,0,200,"ADC charge","Entries");
258 fOutputList->Add(fhAdcNoTime[iSide]);
259 fhAdcWithTime[iSide] = CreateHisto1D(Form("hAdcWithTime%s",side.Data()),Form("ADC (with Leading Time) %s",side.Data()),200,0,200,"ADC charge","Entries");
260 fOutputList->Add(fhAdcWithTime[iSide]);
261 fhTime[iSide] = CreateHisto1D(Form("htimepmt%s",side.Data()),Form("Time measured by TDC %s",side.Data()),400,-100,100,"Leading time (ns)","Entries");
262 fOutputList->Add(fhTime[iSide]);
263 fhWidth[iSide] = CreateHisto1D(Form("hwidth%s",side.Data()),Form("Signal width measured by TDC %s",side.Data()),128,0,200,"Signal width (ns)","Entries");
264 fOutputList->Add(fhWidth[iSide]);
265 fhAdcWidth[iSide] = CreateHisto2D(Form("hadcwidth%s",side.Data()),Form("Time width vs ADC %s",side.Data()),200,0,200,128,0,200,"ADC charge","Width (ns)");
266 fOutputList->Add(fhAdcWidth[iSide]);
267 fhAdcTime[iSide] = CreateHisto2D(Form("hAdcTime%s",side.Data()),Form("ADC vs Time %s",side.Data()),1000,-100,100,200,0,200,"Time (ns)","ADC charge");
268 fOutputList->Add(fhAdcTime[iSide]);
269 }
270
271 fhAdcPMTNoTime = CreateHisto2D("hadcpmtnotime","ADC vs PMT index (no leading time)",64,-0.5,63.5,200,0,200,"PMT index","ADC charge");
272 fhAdcPMTWithTime = CreateHisto2D("hadcpmtwithtime","ADC vs PMT index (with leading time)",64,-0.5,63.5,200,0,200,"PMT index","ADC charge");
273
274 fhTimePMT = CreateHisto2D("htimepmt","Time measured by TDC vs PMT index",64,-0.5,63.5,200,0,100,"PMT Index","Leading time (ns)");
275 fhWidthPMT = CreateHisto2D("hwidthpmt","Time width vs PMT index",64,-0.5,63.5,128,0,200,"PMT Index","Signal width (ns)");
276
277 fhTimeCorr = CreateHisto2D("htimecorr","Average time C side vs. A side",200,0,100,200,0,100,"Time V0A (ns)","Time V0C (ns");
278
279 fhV0ampl = CreateHisto1D("hV0ampl","V0 multiplicity in single channel (all V0 channels)",500,-0.5,499.5);
280
281 fhEvents = CreateHisto2D("hTriggerDecision","Trigger Decision",3,-0.5,2.5,3,-0.5,2.5,"V0A Decision","V0C Decision");
282 fhEvents->GetXaxis()->SetBinLabel(1,"Fake");
283 fhEvents->GetXaxis()->SetBinLabel(2,"BB");
284 fhEvents->GetXaxis()->SetBinLabel(3,"BG");
285 fhEvents->GetYaxis()->SetBinLabel(1,"Fake");
286 fhEvents->GetYaxis()->SetBinLabel(2,"BB");
287 fhEvents->GetYaxis()->SetBinLabel(3,"BG");
288
289 fhVtxXYBB = CreateHisto2D("fhVtxXYBB","XY SPD vertex (bb)",200,-2,2,200,-2,2);
290 fhVtxZBB = CreateHisto1D("fhVtxZBB","Z SPD vertex (bb)",400,-50,50);
291 fhVtxXYBGA = CreateHisto2D("fhVtxXYBGA","XY SPD vertex (bga)",200,-2,2,200,-2,2);
292 fhVtxZBGA = CreateHisto1D("fhVtxZBGA","Z SPD vertex (bga)",400,-50,50);
293 fhVtxXYBGC = CreateHisto2D("fhVtxXYBGC","XY SPD vertex (bgc)",200,-2,2,200,-2,2);
294 fhVtxZBGC = CreateHisto1D("fhVtxZBGC","Z SPD vertex (bgc)",400,-50,50);
295
296 fhPmtMult = CreateHisto2D("hV0CellMult","Number of fired PMTs (V0C vs V0A)",33,-0.5,32.5,33,-0.5,32.5,"# Cell (V0A)","# Cell (V0C)");
297 fOutputList->Add(fhPmtMult);
298
299
300 fOutputList->Add(fhAdcPMTNoTime);
301 fOutputList->Add(fhAdcPMTWithTime);
302
303 fOutputList->Add(fhTimePMT);
304 fOutputList->Add(fhWidthPMT);
305
306 fOutputList->Add(fhTimeCorr);
307 fOutputList->Add(fhV0ampl);
308
309 fOutputList->Add(fhEvents);
310
311 fOutputList->Add(fhVtxXYBB);
312 fOutputList->Add(fhVtxZBB);
313 fOutputList->Add(fhVtxXYBGA);
314 fOutputList->Add(fhVtxZBGA);
315 fOutputList->Add(fhVtxXYBGC);
316 fOutputList->Add(fhVtxZBGC);
317
318}
319//________________________________________________________________________
320Bool_t AliAnaVZEROPbPb::IsGoodEvent(){
321 Bool_t isSelected = (((AliInputEventHandler*)(AliAnalysisManager::GetAnalysisManager()->GetInputEventHandler()))->IsEventSelected() & AliVEvent::kMB);
322 if (!isSelected) return kFALSE;
323 const AliESDVertex *primaryVtx = fESD->GetPrimaryVertex();
324 if (!primaryVtx) return kFALSE;
325 if (!primaryVtx->GetStatus()) return kFALSE;
326 Double_t tPrimaryVtxPosition[3];
327 primaryVtx->GetXYZ(tPrimaryVtxPosition);
328 if (TMath::Abs(tPrimaryVtxPosition[2]) > 10.0) return kFALSE;
329
330 AliCentrality *centrality = fESD->GetCentrality();
331 if (centrality->GetQuality()) return kFALSE;
332
333 return kTRUE;
334}
335//________________________________________________________________________
336void AliAnaVZEROPbPb::FillPerL2TriggerHistos(){
337 TString trigStr(fESD->GetFiredTriggerClasses());
338
339// Bool_t isGoodEvt = IsGoodEvent();
340 TIter iter(fClassesNames);
341 TObjString* name;
342 Int_t iClass =0;
343 while((name = (TObjString*) iter.Next())){
344
345 if (!trigStr.Contains(name->String().Data()) && iClass>0) continue;
346
347 fhOnlineCharge[iClass]->Fill(fEsdV0->GetTriggerChargeA(),fEsdV0->GetTriggerChargeC());
348
349 fhL2Triggers->Fill(iClass);
350
351 fhRecoMult[iClass]->Fill(fEsdV0->GetMTotV0A(),fEsdV0->GetMTotV0C());
352 fhTotRecoMult[iClass]->Fill(fEsdV0->GetMTotV0A()+fEsdV0->GetMTotV0C());
353
354 for(int iTrig = 0; iTrig < 16; ++iTrig){
355 if(fEsdV0->GetTriggerBits() & (1<<iTrig)) fhTriggerBits[iClass]->Fill(iTrig);
356 }
357
358 AliCentrality *centrality = fESD->GetCentrality();
359 Float_t percentile = centrality->GetCentralityPercentile("V0M");
360 Float_t spdPercentile = centrality->GetCentralityPercentile("CL1");
361 if (spdPercentile < 0) spdPercentile = 0;
362 if (percentile < 0) percentile = 0;
363 fhCentrality[iClass]->Fill(percentile);
364 fhV0vsSPDCentrality[iClass]->Fill(spdPercentile,percentile);
365 iClass++;
366 }
367}
368//________________________________________________________________________
369void AliAnaVZEROPbPb::FillQAHistos(){
370 // Count V0 flags
371 Int_t nV0A = 0, nV0C = 0;
372 for(Int_t i = 0; i < 32; ++i) {
373 if (fEsdV0->GetBBFlag(i)) nV0C++;
374 if (fEsdV0->GetBBFlag(i+32)) nV0A++;
375 }
376 fNFlags->Fill((Float_t)nV0A,(Float_t)nV0C);
377
378 for (Int_t i=0; i<64; ++i) {
379 Int_t side = i/32;
380 if (fEsdV0->GetTime(i) < 1e-6) {
381 fhAdcNoTime[side]->Fill(fEsdV0->GetAdc(i));
382 fhAdcPMTNoTime->Fill(i,fEsdV0->GetAdc(i));
383 } else {
384 fhAdcWithTime[side]->Fill(fEsdV0->GetAdc(i));
385 fhAdcPMTWithTime->Fill(i,fEsdV0->GetAdc(i));
386 }
387
388 fhTime[side]->Fill(fEsdV0->GetTime(i));
389 fhWidth[side]->Fill(fEsdV0->GetWidth(i));
390 fhAdcWidth[side]->Fill(fEsdV0->GetAdc(i),fEsdV0->GetWidth(i));
391 fhAdcTime[side]->Fill(fEsdV0->GetTime(i),fEsdV0->GetAdc(i));
392
393 fhTimePMT->Fill(i,fEsdV0->GetTime(i));
394 fhWidthPMT->Fill(i,fEsdV0->GetWidth(i));
395
396 }
397
398 fhTimeCorr->Fill(fEsdV0->GetV0ATime(),fEsdV0->GetV0CTime());
399
400 AliESDVZERO::Decision flaga = fEsdV0->GetV0ADecision();
401 AliESDVZERO::Decision flagc = fEsdV0->GetV0CDecision();
402
403 fhEvents->Fill(flaga,flagc);
404
405 const AliESDVertex *vtx = fESD->GetPrimaryVertexSPD();
406
407 if (flaga <= 1 && flagc <=1) {
408 fhVtxXYBB->Fill(vtx->GetXv(),vtx->GetYv());
409 fhVtxZBB->Fill(vtx->GetZv());
410 }
411 else {
412 if (flaga == 2) {
413 fhVtxXYBGA->Fill(vtx->GetXv(),vtx->GetYv());
414 fhVtxZBGA->Fill(vtx->GetZv());
415 }
416 if (flagc == 2) {
417 fhVtxXYBGC->Fill(vtx->GetXv(),vtx->GetYv());
418 fhVtxZBGC->Fill(vtx->GetZv());
419 }
420 }
421
422 fhPmtMult->Fill(fEsdV0->GetNbPMV0A(),fEsdV0->GetNbPMV0C());
423 for(Int_t i = 0; i < 64; i++) {
424 fhV0ampl->Fill(fEsdV0->GetMultiplicity(i));
425 }
426
427}
428//________________________________________________________________________
429void AliAnaVZEROPbPb::Terminate(Option_t *)
430{
431 // Draw result to the screen
432 // Called once at the end of the query
433
434 fOutputList = dynamic_cast<TList*> (GetOutputData(1));
435 if (!fOutputList) {
436 printf("ERROR: Output list not available\n");
437 return;
438 }
439}