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