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