]> git.uio.no Git - u/mrichter/AliRoot.git/blob - ANALYSIS/AliCentralitySelectionTask.cxx
Analysis code updated
[u/mrichter/AliRoot.git] / ANALYSIS / AliCentralitySelectionTask.cxx
1 /**************************************************************************
2  * Copyright(c) 1998-2008, ALICE Experiment at CERN, All rights reserved. *
3  *                                                                        *
4  * Author: The ALICE Off-line Project.                                    *
5  * Contributors are mentioned in the code where appropriate.              *
6  *                                                                        *
7  * Permission to use, copy, modify and distribute this software and its   *
8  * documentation strictly for non-commercial purposes is hereby granted   *
9  * without fee, provided that the above copyright notice appears in all   *
10  * copies and that both the copyright notice and this permission notice   *
11  * appear in the supporting documentation. The authors make no claims     *
12  * about the suitability of this software for any purpose. It is          *
13  * provided "as is" without express or implied warranty.                  *
14  **************************************************************************/
15
16 //*****************************************************
17 //   Class AliCentralitySelectionTask
18 //   Class to analyze determine centrality            
19 //   author: Alberica Toia
20 //*****************************************************
21
22 #include "AliCentralitySelectionTask.h"
23
24 #include <TTree.h>
25 #include <TList.h>
26 #include <TH1F.h>
27 #include <TH2F.h>
28 #include <TF1.h>
29 #include <TProfile.h>
30 #include <TFile.h>
31 #include <TObjString.h>
32 #include <TString.h>
33 #include <TCanvas.h>
34 #include <TROOT.h>
35 #include <TDirectory.h>
36 #include <TSystem.h>
37 #include <iostream>
38
39 #include "AliAnalysisManager.h"
40 #include "AliHeader.h"
41 #include "AliVEvent.h"
42 #include "AliESD.h"
43 #include "AliESDEvent.h"
44 #include "AliESDHeader.h"
45 #include "AliESDInputHandler.h"
46 #include "AliESDZDC.h"
47 #include "AliESDFMD.h"
48 #include "AliESDTZERO.h"
49 #include "AliESDVZERO.h"
50 #include "AliESDtrackCuts.h"
51 #include "AliESDVertex.h"
52 #include "AliCentrality.h"
53 #include "AliOADBCentrality.h"
54 #include "AliOADBContainer.h"
55 #include "AliMultiplicity.h"
56 #include "AliAODHandler.h"
57 #include "AliAODHeader.h"
58 #include "AliAODEvent.h"
59 #include "AliAODVertex.h"
60 #include "AliAODVZERO.h"
61 #include "AliAODTracklets.h"
62 #include "AliAODMCHeader.h"
63 #include "AliMCEventHandler.h"
64 #include "AliMCEvent.h"
65 #include "AliAODMCParticle.h"
66 #include "AliMCParticle.h"
67 #include "AliStack.h"
68 #include "AliAnalysisTaskSE.h"
69 #include "AliGenEventHeader.h"
70 #include "AliGenHijingEventHeader.h"
71 #include "AliGenDPMjetEventHeader.h"
72 #include "AliGenCocktailEventHeader.h"
73 #include "AliPhysicsSelectionTask.h"
74 #include "AliPhysicsSelection.h"
75 #include "AliBackgroundSelection.h"
76 #include "AliESDUtils.h"
77
78 ClassImp(AliCentralitySelectionTask)
79
80
81 //________________________________________________________________________
82 AliCentralitySelectionTask::AliCentralitySelectionTask():
83 AliAnalysisTaskSE(),
84   fAnalysisInput("ESD"),
85   fIsMCInput(kFALSE),
86   fCurrentRun(-1),
87   fUseScaling(0),
88   fUseCleaning(0),
89   fFillHistos(0),
90   fV0MScaleFactor(0),
91   fSPDScaleFactor(0),
92   fTPCScaleFactor(0),
93   fV0MScaleFactorMC(0),
94   fV0MSPDOutlierPar0(0),  
95   fV0MSPDOutlierPar1(0),  
96   fV0MTPCOutlierPar0(0),  
97   fV0MTPCOutlierPar1(0),  
98   fV0MSPDSigmaOutlierPar0(0),  
99   fV0MSPDSigmaOutlierPar1(0),  
100   fV0MSPDSigmaOutlierPar2(0),  
101   fV0MTPCSigmaOutlierPar0(0),  
102   fV0MTPCSigmaOutlierPar1(0),  
103   fV0MTPCSigmaOutlierPar2(0),   
104   fV0MZDCOutlierPar0(0),            
105   fV0MZDCOutlierPar1(0),            
106   fV0MZDCEcalOutlierPar0(0),   
107   fV0MZDCEcalOutlierPar1(0),   
108   fTrackCuts(0),
109   fEsdTrackCuts(0),
110   fEsdTrackCutsExtra1(0),
111   fEsdTrackCutsExtra2(0),
112   fZVCut(10),
113   fOutliersCut(5),
114   fQuality(999),
115   fIsSelected(0),
116   fMSL(0),
117   fMSH(0),
118   fMUL(0),
119   fMLL(0),
120   fEJE(0),
121   fEGA(0),
122   fPHS(0),
123   fMB(0),
124   fCVHN(0),
125   fCVLN(0),
126   fCVHNbit(0),
127   fCVLNbit(0),
128   fCCENT(0),
129   fCSEMI(0),
130   fCCENTbit(0),
131   fCSEMIbit(0),
132   fCentV0M(0),
133   fCentV0A(0),
134   fCentV0C(0),
135   fCentFMD(0),
136   fCentTRK(0),
137   fCentTKL(0),
138   fCentCL0(0),
139   fCentCL1(0),
140   fCentCND(0),
141   fCentNPA(0),
142   fCentZNA(0),
143   fCentV0MvsFMD(0),
144   fCentTKLvsV0M(0),
145   fCentZEMvsZDC(0),
146   fHtempV0M(0),
147   fHtempV0A(0),
148   fHtempV0C(0),
149   fHtempFMD(0),
150   fHtempTRK(0),
151   fHtempTKL(0),
152   fHtempCL0(0),
153   fHtempCL1(0),
154   fHtempCND(0),
155   fHtempNPA(0),
156   fHtempZNA(0),
157   fHtempV0MvsFMD(0),
158   fHtempTKLvsV0M(0),
159   fHtempZEMvsZDC(0),
160   fOutputList(0),
161   fHOutCentV0M    (0),
162   fHOutCentV0A    (0),
163   fHOutCentV0C    (0),
164   fHOutCentV0MCVHN(0),
165   fHOutCentV0MCVLN(0),
166   fHOutCentV0MCVHNinMB(0),
167   fHOutCentV0MCVLNinMB(0),
168   fHOutCentV0MCCENT(0),
169   fHOutCentV0MCSEMI(0),
170   fHOutCentV0MCCENTinMB(0),
171   fHOutCentV0MCSEMIinMB(0),
172   fHOutCentV0MMSL(0),
173   fHOutCentV0MMSH(0),
174   fHOutCentV0MMUL(0),
175   fHOutCentV0MMLL(0),
176   fHOutCentV0MEJE(0),
177   fHOutCentV0MEGA(0),
178   fHOutCentV0MPHS(0),
179   fHOutCentV0MMSLinMB(0),
180   fHOutCentV0MMSHinMB(0),
181   fHOutCentV0MMULinMB(0),
182   fHOutCentV0MMLLinMB(0),
183   fHOutCentV0MEJEinMB(0),
184   fHOutCentV0MEGAinMB(0),
185   fHOutCentV0MPHSinMB(0),
186   fHOutCentFMD     (0),
187   fHOutCentTRK     (0),
188   fHOutCentTKL     (0),
189   fHOutCentCL0     (0),
190   fHOutCentCL1     (0),
191   fHOutCentCND     (0),
192   fHOutCentNPA     (0),
193   fHOutCentZNA     (0),
194   fHOutCentV0MvsFMD(0),
195   fHOutCentTKLvsV0M(0),
196   fHOutCentZEMvsZDC(0),
197   fHOutCentV0MvsCentCL1(0),
198   fHOutCentV0MvsCentTRK(0),
199   fHOutCentTRKvsCentCL1(0),
200   fHOutCentV0MvsCentZDC(0),
201   fHOutCentV0AvsCentV0C(0),
202   fHOutCentV0AvsCentTRK(0),
203   fHOutCentV0AvsCentCND(0),
204   fHOutCentV0AvsCentCL1(0),
205   fHOutCentV0CvsCentTRK(0),
206   fHOutCentV0CvsCentCND(0),
207   fHOutCentV0CvsCentCL1(0),
208   fHOutCentNPAvsCentV0A(0),
209   fHOutCentNPAvsCentV0C(0),
210   fHOutCentNPAvsCentTRK(0),
211   fHOutCentNPAvsCentCND(0),
212   fHOutCentNPAvsCentCL1(0),
213   fHOutCentZNAvsCentV0A(0),
214   fHOutCentZNAvsCentV0C(0),
215   fHOutCentZNAvsCentTRK(0),
216   fHOutCentZNAvsCentCND(0),
217   fHOutCentZNAvsCentCL1(0),
218   fHOutMultV0AC(0),
219   fHOutMultV0M(0),
220   fHOutMultV0A(0),
221   fHOutMultV0C(0),
222   fHOutMultV0Mnc(0),
223   fHOutMultV0Anc(0),
224   fHOutMultV0Cnc(0),
225   fHOutMultV0O(0),
226   fHOutMultV0Cells(0),
227   fHOutMultFMD(0),
228   fHOutMultTRK(0),
229   fHOutMultTKL(0),
230   fHOutMultCL0(0),
231   fHOutMultCL1(0),
232   fHOutMultCND(0),
233   fHOutMultNPA(0),
234   fHOutMultZNA(0),
235   fHOutMultV0MvsZDN(0),
236   fHOutMultZEMvsZDN(0),
237   fHOutMultV0MvsZDC(0),
238   fHOutMultZEMvsZDC(0),
239   fHOutMultZEMvsZDCw(0),
240   fHOutMultV0MvsCL1(0),
241   fHOutMultV0MvsTRK(0),
242   fHOutMultTRKvsCL1(0),
243   fHOutMultV0MvsV0O(0),
244   fHOutMultV0OvsCL1(0),
245   fHOutMultV0OvsTRK(0),
246   fHOutMultCL1vsTKL(0),
247   fHOutCentV0Mqual1(0),
248   fHOutCentTRKqual1(0),
249   fHOutCentCL1qual1(0),
250   fHOutMultV0MvsCL1qual1(0),
251   fHOutMultV0MvsTRKqual1(0),
252   fHOutMultTRKvsCL1qual1(0),
253   fHOutCentV0Mqual2(0),
254   fHOutCentTRKqual2(0),
255   fHOutCentCL1qual2(0),
256   fHOutMultV0MvsCL1qual2(0),
257   fHOutMultV0MvsTRKqual2(0),
258   fHOutMultTRKvsCL1qual2(0),
259   fHOutQuality(0),
260   fHOutVertex(0),
261   fHOutVertexT0(0)
262 {   
263   // Default constructor
264   AliInfo("Centrality Selection enabled.");
265
266   fUseScaling=kTRUE;
267   fUseCleaning=kTRUE;
268   fFillHistos=kFALSE;
269   fBranchNames="ESD:AliESDRun.,AliESDHeader.,AliESDZDC.,AliESDFMD.,AliESDVZERO.,AliESDTZERO."
270     ",SPDVertex.,TPCVertex.,PrimaryVertex.,AliMultiplicity.,Tracks ";
271 }   
272
273 //________________________________________________________________________
274 AliCentralitySelectionTask::AliCentralitySelectionTask(const char *name):
275   AliAnalysisTaskSE(name),
276   fAnalysisInput("ESD"),
277   fIsMCInput(kFALSE),
278   fCurrentRun(-1),
279   fUseScaling(0),
280   fUseCleaning(0),
281   fFillHistos(0),
282   fV0MScaleFactor(0),
283   fSPDScaleFactor(0),
284   fTPCScaleFactor(0),
285   fV0MScaleFactorMC(0),
286   fV0MSPDOutlierPar0(0),  
287   fV0MSPDOutlierPar1(0),  
288   fV0MTPCOutlierPar0(0),  
289   fV0MTPCOutlierPar1(0),  
290   fV0MSPDSigmaOutlierPar0(0),  
291   fV0MSPDSigmaOutlierPar1(0),  
292   fV0MSPDSigmaOutlierPar2(0),  
293   fV0MTPCSigmaOutlierPar0(0),  
294   fV0MTPCSigmaOutlierPar1(0),  
295   fV0MTPCSigmaOutlierPar2(0),   
296   fV0MZDCOutlierPar0(0),            
297   fV0MZDCOutlierPar1(0),            
298   fV0MZDCEcalOutlierPar0(0),   
299   fV0MZDCEcalOutlierPar1(0),   
300   fTrackCuts(0),
301   fEsdTrackCuts(0),
302   fEsdTrackCutsExtra1(0),
303   fEsdTrackCutsExtra2(0),
304   fZVCut(10),
305   fOutliersCut(5),
306   fQuality(999),
307   fIsSelected(0),
308   fMSL(0),
309   fMSH(0),
310   fMUL(0),
311   fMLL(0),
312   fEJE(0),
313   fEGA(0),
314   fPHS(0),
315   fMB(0),
316   fCVHN(0),
317   fCVLN(0),
318   fCVHNbit(0),
319   fCVLNbit(0),
320   fCCENT(0),
321   fCSEMI(0),
322   fCCENTbit(0),
323   fCSEMIbit(0),
324   fCentV0M(0),
325   fCentV0A(0),
326   fCentV0C(0),
327   fCentFMD(0),
328   fCentTRK(0),
329   fCentTKL(0),
330   fCentCL0(0),
331   fCentCL1(0),
332   fCentCND(0),
333   fCentNPA(0),
334   fCentZNA(0),
335   fCentV0MvsFMD(0),
336   fCentTKLvsV0M(0),
337   fCentZEMvsZDC(0),
338   fHtempV0M(0),
339   fHtempV0A(0),
340   fHtempV0C(0),
341   fHtempFMD(0),
342   fHtempTRK(0),
343   fHtempTKL(0),
344   fHtempCL0(0),
345   fHtempCL1(0),
346   fHtempCND(0),
347   fHtempNPA(0),
348   fHtempZNA(0),
349   fHtempV0MvsFMD(0),
350   fHtempTKLvsV0M(0),
351   fHtempZEMvsZDC(0),
352   fOutputList(0),
353   fHOutCentV0M    (0),
354   fHOutCentV0A    (0),
355   fHOutCentV0C    (0),
356   fHOutCentV0MCVHN(0),
357   fHOutCentV0MCVLN(0),
358   fHOutCentV0MCVHNinMB(0),
359   fHOutCentV0MCVLNinMB(0),
360   fHOutCentV0MCCENT(0),
361   fHOutCentV0MCSEMI(0),
362   fHOutCentV0MCCENTinMB(0),
363   fHOutCentV0MCSEMIinMB(0),
364   fHOutCentV0MMSL(0),
365   fHOutCentV0MMSH(0),
366   fHOutCentV0MMUL(0),
367   fHOutCentV0MMLL(0),
368   fHOutCentV0MEJE(0),
369   fHOutCentV0MEGA(0),
370   fHOutCentV0MPHS(0),
371   fHOutCentV0MMSLinMB(0),
372   fHOutCentV0MMSHinMB(0),
373   fHOutCentV0MMULinMB(0),
374   fHOutCentV0MMLLinMB(0),
375   fHOutCentV0MEJEinMB(0),
376   fHOutCentV0MEGAinMB(0),
377   fHOutCentV0MPHSinMB(0),
378   fHOutCentFMD     (0),
379   fHOutCentTRK     (0),
380   fHOutCentTKL     (0),
381   fHOutCentCL0     (0),
382   fHOutCentCL1     (0),
383   fHOutCentCND     (0),
384   fHOutCentNPA     (0),
385   fHOutCentZNA     (0),
386   fHOutCentV0MvsFMD(0),
387   fHOutCentTKLvsV0M(0),
388   fHOutCentZEMvsZDC(0),
389   fHOutCentV0MvsCentCL1(0),
390   fHOutCentV0MvsCentTRK(0),
391   fHOutCentTRKvsCentCL1(0),
392   fHOutCentV0MvsCentZDC(0),
393   fHOutCentV0AvsCentV0C(0),
394   fHOutCentV0AvsCentTRK(0),
395   fHOutCentV0AvsCentCND(0),
396   fHOutCentV0AvsCentCL1(0),
397   fHOutCentV0CvsCentTRK(0),
398   fHOutCentV0CvsCentCND(0),
399   fHOutCentV0CvsCentCL1(0),
400   fHOutCentNPAvsCentV0A(0),
401   fHOutCentNPAvsCentV0C(0),
402   fHOutCentNPAvsCentTRK(0),
403   fHOutCentNPAvsCentCND(0),
404   fHOutCentNPAvsCentCL1(0),
405   fHOutCentZNAvsCentV0A(0),
406   fHOutCentZNAvsCentV0C(0),
407   fHOutCentZNAvsCentTRK(0),
408   fHOutCentZNAvsCentCND(0),
409   fHOutCentZNAvsCentCL1(0),
410   fHOutMultV0AC(0),
411   fHOutMultV0M(0),
412   fHOutMultV0A(0),
413   fHOutMultV0C(0),
414   fHOutMultV0Mnc(0),
415   fHOutMultV0Anc(0),
416   fHOutMultV0Cnc(0),
417   fHOutMultV0O(0),
418   fHOutMultV0Cells(0),
419   fHOutMultFMD(0),
420   fHOutMultTRK(0),
421   fHOutMultTKL(0),
422   fHOutMultCL0(0),
423   fHOutMultCL1(0),
424   fHOutMultCND(0),
425   fHOutMultNPA(0),
426   fHOutMultZNA(0),
427   fHOutMultV0MvsZDN(0),
428   fHOutMultZEMvsZDN(0),
429   fHOutMultV0MvsZDC(0),
430   fHOutMultZEMvsZDC(0),
431   fHOutMultZEMvsZDCw(0),
432   fHOutMultV0MvsCL1(0),
433   fHOutMultV0MvsTRK(0),
434   fHOutMultTRKvsCL1(0),
435   fHOutMultV0MvsV0O(0),
436   fHOutMultV0OvsCL1(0),
437   fHOutMultV0OvsTRK(0),
438   fHOutMultCL1vsTKL(0),
439   fHOutCentV0Mqual1(0),
440   fHOutCentTRKqual1(0),
441   fHOutCentCL1qual1(0),
442   fHOutMultV0MvsCL1qual1(0),
443   fHOutMultV0MvsTRKqual1(0),
444   fHOutMultTRKvsCL1qual1(0),
445   fHOutCentV0Mqual2(0),
446   fHOutCentTRKqual2(0),
447   fHOutCentCL1qual2(0),
448   fHOutMultV0MvsCL1qual2(0),
449   fHOutMultV0MvsTRKqual2(0),
450   fHOutMultTRKvsCL1qual2(0),
451   fHOutQuality(0),
452   fHOutVertex(0),
453   fHOutVertexT0(0)
454 {
455   // Default constructor
456   AliInfo("Centrality Selection enabled.");
457   //DefineOutput(1, TList::Class());
458   fUseScaling=kTRUE;
459   fUseCleaning=kTRUE;
460   fFillHistos=kFALSE;
461   fBranchNames="ESD:AliESDRun.,AliESDHeader.,AliESDZDC.,AliESDFMD.,AliESDVZERO.,AliESDTZERO."
462     ",SPDVertex.,TPCVertex.,PrimaryVertex.,AliMultiplicity.,Tracks ";
463 }
464
465 //________________________________________________________________________
466 AliCentralitySelectionTask& AliCentralitySelectionTask::operator=(const AliCentralitySelectionTask& c)
467 {
468   // Assignment operator
469   if (this!=&c) {
470     AliAnalysisTaskSE::operator=(c);
471   }
472   return *this;
473 }
474
475 //________________________________________________________________________
476 AliCentralitySelectionTask::AliCentralitySelectionTask(const AliCentralitySelectionTask& ana):
477   AliAnalysisTaskSE(ana),
478   fAnalysisInput(ana.fAnalysisInput),
479   fIsMCInput(ana.fIsMCInput),
480   fCurrentRun(ana.fCurrentRun),
481   fUseScaling(ana.fUseScaling),
482   fUseCleaning(ana.fUseCleaning),
483   fFillHistos(ana.fFillHistos),
484   fV0MScaleFactor(ana.fV0MScaleFactor),
485   fSPDScaleFactor(ana.fSPDScaleFactor),
486   fTPCScaleFactor(ana.fTPCScaleFactor),
487   fV0MScaleFactorMC(ana.fV0MScaleFactorMC),
488   fV0MSPDOutlierPar0(ana.fV0MSPDOutlierPar0),  
489   fV0MSPDOutlierPar1(ana.fV0MSPDOutlierPar1),  
490   fV0MTPCOutlierPar0(ana.fV0MTPCOutlierPar0),  
491   fV0MTPCOutlierPar1(ana.fV0MTPCOutlierPar1),  
492   fV0MSPDSigmaOutlierPar0(ana.fV0MSPDSigmaOutlierPar0),  
493   fV0MSPDSigmaOutlierPar1(ana.fV0MSPDSigmaOutlierPar1),  
494   fV0MSPDSigmaOutlierPar2(ana.fV0MSPDSigmaOutlierPar2),  
495   fV0MTPCSigmaOutlierPar0(ana.fV0MTPCSigmaOutlierPar0),  
496   fV0MTPCSigmaOutlierPar1(ana.fV0MTPCSigmaOutlierPar1),  
497   fV0MTPCSigmaOutlierPar2(ana.fV0MTPCSigmaOutlierPar2), 
498   fV0MZDCOutlierPar0(ana.fV0MZDCOutlierPar0),       
499   fV0MZDCOutlierPar1(ana.fV0MZDCOutlierPar1),       
500   fV0MZDCEcalOutlierPar0(ana.fV0MZDCEcalOutlierPar0),   
501   fV0MZDCEcalOutlierPar1(ana.fV0MZDCEcalOutlierPar1),   
502   fTrackCuts(ana.fTrackCuts),
503   fEsdTrackCuts(ana.fEsdTrackCuts),
504   fEsdTrackCutsExtra1(ana.fEsdTrackCutsExtra1),
505   fEsdTrackCutsExtra2(ana.fEsdTrackCutsExtra2),
506   fZVCut(ana.fZVCut),
507   fOutliersCut(ana.fOutliersCut),
508   fQuality(ana.fQuality),
509   fIsSelected(ana.fIsSelected),
510   fMSL(ana.fMSL),
511   fMSH(ana.fMSH),
512   fMUL(ana.fMUL),
513   fMLL(ana.fMLL),
514   fEJE(ana.fEJE),
515   fEGA(ana.fEGA),
516   fPHS(ana.fPHS),
517   fMB(ana.fMB),
518   fCVHN(ana.fCVHN),
519   fCVLN(ana.fCVLN),
520   fCVHNbit(ana.fCVHNbit),
521   fCVLNbit(ana.fCVLNbit),
522   fCCENT(ana.fCCENT),
523   fCSEMI(ana.fCSEMI),
524   fCCENTbit(ana.fCCENTbit),
525   fCSEMIbit(ana.fCSEMIbit),
526   fCentV0M(ana.fCentV0M),
527   fCentV0A(ana.fCentV0A),
528   fCentV0C(ana.fCentV0C),
529   fCentFMD(ana.fCentFMD),
530   fCentTRK(ana.fCentTRK),
531   fCentTKL(ana.fCentTKL),
532   fCentCL0(ana.fCentCL0),
533   fCentCL1(ana.fCentCL1),
534   fCentCND(ana.fCentCND),
535   fCentNPA(ana.fCentNPA),
536   fCentZNA(ana.fCentZNA),
537   fCentV0MvsFMD(ana.fCentV0MvsFMD),
538   fCentTKLvsV0M(ana.fCentTKLvsV0M),
539   fCentZEMvsZDC(ana.fCentZEMvsZDC),
540   fHtempV0M(ana.fHtempV0M),
541   fHtempV0A(ana.fHtempV0A),
542   fHtempV0C(ana.fHtempV0C),
543   fHtempFMD(ana.fHtempFMD),
544   fHtempTRK(ana.fHtempTRK),
545   fHtempTKL(ana.fHtempTKL),
546   fHtempCL0(ana.fHtempCL0),
547   fHtempCL1(ana.fHtempCL1),
548   fHtempCND(ana.fHtempCND),
549   fHtempNPA(ana.fHtempNPA),
550   fHtempZNA(ana.fHtempZNA),
551   fHtempV0MvsFMD(ana.fHtempV0MvsFMD),
552   fHtempTKLvsV0M(ana.fHtempTKLvsV0M),
553   fHtempZEMvsZDC(ana.fHtempZEMvsZDC),
554   fOutputList(ana.fOutputList),
555   fHOutCentV0M    (ana.fHOutCentV0M    ),
556   fHOutCentV0A    (ana.fHOutCentV0A    ),
557   fHOutCentV0C    (ana.fHOutCentV0C    ),
558   fHOutCentV0MCVHN(ana.fHOutCentV0MCVHN),
559   fHOutCentV0MCVLN(ana.fHOutCentV0MCVLN),
560   fHOutCentV0MCVHNinMB(ana.fHOutCentV0MCVHNinMB),
561   fHOutCentV0MCVLNinMB(ana.fHOutCentV0MCVLNinMB),
562   fHOutCentV0MCCENT(ana.fHOutCentV0MCCENT),
563   fHOutCentV0MCSEMI(ana.fHOutCentV0MCSEMI),
564   fHOutCentV0MCCENTinMB(ana.fHOutCentV0MCCENTinMB),
565   fHOutCentV0MCSEMIinMB(ana.fHOutCentV0MCSEMIinMB),
566   fHOutCentV0MMSL(ana.fHOutCentV0MMSL),
567   fHOutCentV0MMSH(ana.fHOutCentV0MMSH),
568   fHOutCentV0MMUL(ana.fHOutCentV0MMUL),
569   fHOutCentV0MMLL(ana.fHOutCentV0MMLL),
570   fHOutCentV0MEJE(ana.fHOutCentV0MEJE),
571   fHOutCentV0MEGA(ana.fHOutCentV0MEGA),
572   fHOutCentV0MPHS(ana.fHOutCentV0MPHS),
573   fHOutCentV0MMSLinMB(ana.fHOutCentV0MMSLinMB),
574   fHOutCentV0MMSHinMB(ana.fHOutCentV0MMSHinMB),
575   fHOutCentV0MMULinMB(ana.fHOutCentV0MMULinMB),
576   fHOutCentV0MMLLinMB(ana.fHOutCentV0MMLLinMB),
577   fHOutCentV0MEJEinMB(ana.fHOutCentV0MEJEinMB),
578   fHOutCentV0MEGAinMB(ana.fHOutCentV0MEGAinMB),
579   fHOutCentV0MPHSinMB(ana.fHOutCentV0MPHSinMB),
580   fHOutCentFMD     (ana.fHOutCentFMD     ),
581   fHOutCentTRK     (ana.fHOutCentTRK     ),
582   fHOutCentTKL     (ana.fHOutCentTKL     ),
583   fHOutCentCL0     (ana.fHOutCentCL0     ),
584   fHOutCentCL1     (ana.fHOutCentCL1     ),
585   fHOutCentCND     (ana.fHOutCentCND     ),
586   fHOutCentNPA     (ana.fHOutCentNPA     ),
587   fHOutCentZNA     (ana.fHOutCentZNA     ),
588   fHOutCentV0MvsFMD(ana.fHOutCentV0MvsFMD),
589   fHOutCentTKLvsV0M(ana.fHOutCentTKLvsV0M),
590   fHOutCentZEMvsZDC(ana.fHOutCentZEMvsZDC),
591   fHOutCentV0MvsCentCL1(ana.fHOutCentV0MvsCentCL1),
592   fHOutCentV0MvsCentTRK(ana.fHOutCentV0MvsCentTRK),
593   fHOutCentTRKvsCentCL1(ana.fHOutCentTRKvsCentCL1),
594   fHOutCentV0MvsCentZDC(ana.fHOutCentV0MvsCentZDC),
595   fHOutCentV0AvsCentV0C(ana.fHOutCentV0AvsCentV0C),
596   fHOutCentV0AvsCentTRK(ana.fHOutCentV0AvsCentTRK),
597   fHOutCentV0AvsCentCND(ana.fHOutCentV0AvsCentCND),
598   fHOutCentV0AvsCentCL1(ana.fHOutCentV0AvsCentCL1),
599   fHOutCentV0CvsCentTRK(ana.fHOutCentV0CvsCentTRK),
600   fHOutCentV0CvsCentCND(ana.fHOutCentV0CvsCentCND),
601   fHOutCentV0CvsCentCL1(ana.fHOutCentV0CvsCentCL1),
602   fHOutCentNPAvsCentV0A(ana.fHOutCentNPAvsCentV0A),
603   fHOutCentNPAvsCentV0C(ana.fHOutCentNPAvsCentV0C),
604   fHOutCentNPAvsCentTRK(ana.fHOutCentNPAvsCentTRK),
605   fHOutCentNPAvsCentCND(ana.fHOutCentNPAvsCentCND),
606   fHOutCentNPAvsCentCL1(ana.fHOutCentNPAvsCentCL1),
607   fHOutCentZNAvsCentV0A(ana.fHOutCentZNAvsCentV0A),
608   fHOutCentZNAvsCentV0C(ana.fHOutCentZNAvsCentV0C),
609   fHOutCentZNAvsCentTRK(ana.fHOutCentZNAvsCentTRK),
610   fHOutCentZNAvsCentCND(ana.fHOutCentZNAvsCentCND),
611   fHOutCentZNAvsCentCL1(ana.fHOutCentZNAvsCentCL1),
612   fHOutMultV0AC(ana.fHOutMultV0AC),
613   fHOutMultV0M(ana.fHOutMultV0M),
614   fHOutMultV0A(ana.fHOutMultV0A),
615   fHOutMultV0C(ana.fHOutMultV0C),
616   fHOutMultV0Mnc(ana.fHOutMultV0Mnc),
617   fHOutMultV0Anc(ana.fHOutMultV0Anc),
618   fHOutMultV0Cnc(ana.fHOutMultV0Cnc),
619   fHOutMultV0O(ana.fHOutMultV0O),
620   fHOutMultV0Cells(ana.fHOutMultV0Cells),
621   fHOutMultFMD(ana.fHOutMultFMD),
622   fHOutMultTRK(ana.fHOutMultTRK),
623   fHOutMultTKL(ana.fHOutMultTKL),
624   fHOutMultCL0(ana.fHOutMultCL0),
625   fHOutMultCL1(ana.fHOutMultCL1),
626   fHOutMultCND(ana.fHOutMultCND),
627   fHOutMultNPA(ana.fHOutMultNPA),
628   fHOutMultZNA(ana.fHOutMultZNA),
629   fHOutMultV0MvsZDN(ana.fHOutMultV0MvsZDN),
630   fHOutMultZEMvsZDN(ana.fHOutMultZEMvsZDN),
631   fHOutMultV0MvsZDC(ana.fHOutMultV0MvsZDC),
632   fHOutMultZEMvsZDC(ana.fHOutMultZEMvsZDC),
633   fHOutMultZEMvsZDCw(ana.fHOutMultZEMvsZDCw),
634   fHOutMultV0MvsCL1(ana.fHOutMultV0MvsCL1),
635   fHOutMultV0MvsTRK(ana.fHOutMultV0MvsTRK),
636   fHOutMultTRKvsCL1(ana.fHOutMultTRKvsCL1),
637   fHOutMultV0MvsV0O(ana.fHOutMultV0MvsV0O),
638   fHOutMultV0OvsCL1(ana.fHOutMultV0OvsCL1),
639   fHOutMultV0OvsTRK(ana.fHOutMultV0OvsTRK),
640   fHOutMultCL1vsTKL(ana.fHOutMultCL1vsTKL),
641   fHOutCentV0Mqual1(ana.fHOutCentV0Mqual1),
642   fHOutCentTRKqual1(ana.fHOutCentTRKqual1),
643   fHOutCentCL1qual1(ana.fHOutCentCL1qual1),
644   fHOutMultV0MvsCL1qual1(ana.fHOutMultV0MvsCL1qual1),
645   fHOutMultV0MvsTRKqual1(ana.fHOutMultV0MvsTRKqual1),
646   fHOutMultTRKvsCL1qual1(ana.fHOutMultTRKvsCL1qual1),
647   fHOutCentV0Mqual2(ana.fHOutCentV0Mqual2),
648   fHOutCentTRKqual2(ana.fHOutCentTRKqual2),
649   fHOutCentCL1qual2(ana.fHOutCentCL1qual2),
650   fHOutMultV0MvsCL1qual2(ana.fHOutMultV0MvsCL1qual2),
651   fHOutMultV0MvsTRKqual2(ana.fHOutMultV0MvsTRKqual2),
652   fHOutMultTRKvsCL1qual2(ana.fHOutMultTRKvsCL1qual2),
653   fHOutQuality(ana.fHOutQuality),
654   fHOutVertex(ana.fHOutVertex),
655   fHOutVertexT0(ana.fHOutVertexT0)
656 {
657   // Copy Constructor   
658
659 }
660
661 //________________________________________________________________________
662 AliCentralitySelectionTask::~AliCentralitySelectionTask()
663 {
664   // Destructor  
665   if (fOutputList && !AliAnalysisManager::GetAnalysisManager()->IsProofMode()) delete fOutputList;
666   if (fTrackCuts) delete fTrackCuts;
667   if (fEsdTrackCuts) delete fEsdTrackCuts;
668   if (fEsdTrackCutsExtra1) delete fEsdTrackCutsExtra1;
669   if (fEsdTrackCutsExtra2) delete fEsdTrackCutsExtra2;
670 }  
671
672 //________________________________________________________________________
673 void AliCentralitySelectionTask::UserCreateOutputObjects()
674 {  
675   // Create the output containers
676   if(fDebug>1) printf("AnalysisCentralitySelectionTask::UserCreateOutputObjects() \n");
677   AliLog::SetClassDebugLevel("AliCentralitySelectionTask", AliLog::kInfo);
678
679   if (fFillHistos) {    
680     fOutputList = new TList();
681     fOutputList->SetOwner();
682     fHOutCentV0M     = new TH1F("fHOutCentV0M","fHOutCentV0M; Centrality V0",505,0,101);
683     fHOutCentV0A    = new TH1F("fHOutCentV0A","fHOutCentV0A; Centrality V0A",505,0,101);
684     fHOutCentV0C    = new TH1F("fHOutCentV0C","fHOutCentV0C; Centrality V0C",505,0,101);
685     fHOutCentV0MCVHN= new TH1F("fHOutCentV0M_CVHN","fHOutCentV0M_CVHN; Centrality V0",505,0,101);
686     fHOutCentV0MCVLN= new TH1F("fHOutCentV0M_CVLN","fHOutCentV0M_CVLN; Centrality V0",505,0,101);
687     fHOutCentV0MCVHNinMB= new TH1F("fHOutCentV0M_CVHNinMB","fHOutCentV0M_CVHN; Centrality V0",505,0,101);
688     fHOutCentV0MCVLNinMB= new TH1F("fHOutCentV0M_CVLNinMB","fHOutCentV0M_CVLN; Centrality V0",505,0,101);
689     fHOutCentV0MCCENT= new TH1F("fHOutCentV0M_CCENT","fHOutCentV0M_CCENT; Centrality V0",505,0,101);
690     fHOutCentV0MCSEMI= new TH1F("fHOutCentV0M_CSEMI","fHOutCentV0M_CSEMI; Centrality V0",505,0,101);
691     fHOutCentV0MCCENTinMB= new TH1F("fHOutCentV0M_CCENTinMB","fHOutCentV0M_CCENT; Centrality V0",505,0,101);
692     fHOutCentV0MCSEMIinMB= new TH1F("fHOutCentV0M_CSEMIinMB","fHOutCentV0M_CSEMI; Centrality V0",505,0,101);
693     fHOutCentV0MMSL= new TH1F("fHOutCentV0M_MSL","fHOutCentV0M_MSL; Centrality V0",505,0,101);
694     fHOutCentV0MMSH= new TH1F("fHOutCentV0M_MSH","fHOutCentV0M_MSH; Centrality V0",505,0,101);
695     fHOutCentV0MMUL= new TH1F("fHOutCentV0M_MUL","fHOutCentV0M_MUL; Centrality V0",505,0,101);
696     fHOutCentV0MMLL= new TH1F("fHOutCentV0M_MLL","fHOutCentV0M_MLL; Centrality V0",505,0,101);
697     fHOutCentV0MEJE= new TH1F("fHOutCentV0M_EJE","fHOutCentV0M_EJE; Centrality V0",505,0,101);
698     fHOutCentV0MEGA= new TH1F("fHOutCentV0M_EGA","fHOutCentV0M_EGA; Centrality V0",505,0,101);
699     fHOutCentV0MPHS= new TH1F("fHOutCentV0M_PHS","fHOutCentV0M_PHS; Centrality V0",505,0,101);
700     fHOutCentV0MMSLinMB= new TH1F("fHOutCentV0M_MSLinMB","fHOutCentV0M_MSLinMB; Centrality V0",505,0,101);
701     fHOutCentV0MMSHinMB= new TH1F("fHOutCentV0M_MSHinMB","fHOutCentV0M_MSHinMB; Centrality V0",505,0,101);
702     fHOutCentV0MMULinMB= new TH1F("fHOutCentV0M_MULinMB","fHOutCentV0M_MULinMB; Centrality V0",505,0,101);
703     fHOutCentV0MMLLinMB= new TH1F("fHOutCentV0M_MLLinMB","fHOutCentV0M_MLLinMB; Centrality V0",505,0,101);
704     fHOutCentV0MEJEinMB= new TH1F("fHOutCentV0M_EJEinMB","fHOutCentV0M_EJEinMB; Centrality V0",505,0,101);
705     fHOutCentV0MEGAinMB= new TH1F("fHOutCentV0M_EGAinMB","fHOutCentV0M_EGAinMB; Centrality V0",505,0,101);
706     fHOutCentV0MPHSinMB= new TH1F("fHOutCentV0M_PHSinMB","fHOutCentV0M_PHSinMB; Centrality V0",505,0,101);
707     fHOutCentFMD     = new TH1F("fHOutCentFMD","fHOutCentFMD; Centrality FMD",505,0,101);
708     fHOutCentTRK     = new TH1F("fHOutCentTRK","fHOutCentTRK; Centrality TPC",505,0,101);
709     fHOutCentTKL     = new TH1F("fHOutCentTKL","fHOutCentTKL; Centrality tracklets",505,0,101);
710     fHOutCentCL0     = new TH1F("fHOutCentCL0","fHOutCentCL0; Centrality SPD inner",505,0,101);
711     fHOutCentCL1     = new TH1F("fHOutCentCL1","fHOutCentCL1; Centrality SPD outer",505,0,101);
712     fHOutCentCND     = new TH1F("fHOutCentCND","fHOutCentCND; Centrality candle",505,0,101);
713     fHOutCentNPA     = new TH1F("fHOutCentNPA","fHOutCentNPA; Centrality Npart",505,0,101);
714     fHOutCentZNA     = new TH1F("fHOutCentZNA","fHOutCentZNA; Centrality ZNA",505,0,101);
715     fHOutCentV0MvsFMD= new TH1F("fHOutCentV0MvsFMD","fHOutCentV0MvsFMD; Centrality V0 vs FMD",505,0,101);
716     fHOutCentTKLvsV0M= new TH1F("fHOutCentTKLvsV0M","fHOutCentTKLvsV0M; Centrality tracklets vs V0",505,0,101);
717     fHOutCentZEMvsZDC= new TH1F("fHOutCentZEMvsZDC","fHOutCentZEMvsZDC; Centrality ZEM vs ZDC",505,0,101);
718     fHOutCentV0MvsCentCL1= new TH2F("fHOutCentV0MvsCentCL1","fHOutCentV0MvsCentCL1; Cent V0; Cent SPD",505,0,101,505,0,101);
719     fHOutCentV0MvsCentTRK= new TH2F("fHOutCentV0MvsCentTRK","fHOutCentV0MvsCentTRK; Cent V0; Cent TPC",505,0,101,505,0,101);
720     fHOutCentTRKvsCentCL1= new TH2F("fHOutCentTRKvsCentCL1","fHOutCentTRKvsCentCL1; Cent TPC; Cent SPD",505,0,101,505,0,101);
721     fHOutCentV0MvsCentZDC= new TH2F("fHOutCentV0MvsCentZDC","fHOutCentV0MvsCentZDC; Cent V0; Cent ZDC",505,0,101,505,0,101);
722     fHOutCentV0AvsCentV0C= new TH2F("fHOutCentV0AvsCentV0C","fHOutCentV0AvsCentV0C; Cent V0A; Cent V0C;", 505,0,101,505,0,101);
723     fHOutCentV0AvsCentTRK= new TH2F("fHOutCentV0AvsCentTRK","fHOutCentV0AvsCentTRK; Cent V0A; Cent TRK;", 505,0,101,505,0,101);
724     fHOutCentV0AvsCentCND= new TH2F("fHOutCentV0AvsCentCND","fHOutCentV0AvsCentCND; Cent V0A; Cent CND;", 505,0,101,505,0,101);
725     fHOutCentV0AvsCentCL1= new TH2F("fHOutCentV0AvsCentCL1","fHOutCentV0AvsCentCL1; Cent V0A; Cent CL1;", 505,0,101,505,0,101);
726     fHOutCentV0CvsCentTRK= new TH2F("fHOutCentV0CvsCentTRK","fHOutCentV0CvsCentTRK; Cent V0C; Cent TRK;", 505,0,101,505,0,101);
727     fHOutCentV0CvsCentCND= new TH2F("fHOutCentV0CvsCentCND","fHOutCentV0CvsCentCND; Cent V0C; Cent CND;", 505,0,101,505,0,101);
728     fHOutCentV0CvsCentCL1= new TH2F("fHOutCentV0CvsCentCL1","fHOutCentV0CvsCentCL1; Cent V0C; Cent CL1;", 505,0,101,505,0,101);
729     fHOutCentNPAvsCentV0A= new TH2F("fHOutCentNPAvsCentV0A","fHOutCentNPAvsCentV0A; Cent NPA; Cent V0A;", 505,0,101,505,0,101);
730     fHOutCentNPAvsCentV0C= new TH2F("fHOutCentNPAvsCentV0C","fHOutCentNPAvsCentV0C; Cent NPA; Cent V0C;", 505,0,101,505,0,101);
731     fHOutCentNPAvsCentTRK= new TH2F("fHOutCentNPAvsCentTRK","fHOutCentNPAvsCentTRK; Cent NPA; Cent TRK;", 505,0,101,505,0,101);
732     fHOutCentNPAvsCentCND= new TH2F("fHOutCentNPAvsCentCND","fHOutCentNPAvsCentCND; Cent NPA; Cent CND;", 505,0,101,505,0,101);
733     fHOutCentNPAvsCentCL1= new TH2F("fHOutCentNPAvsCentCL1","fHOutCentNPAvsCentCL1; Cent NPA; Cent CL1;", 505,0,101,505,0,101);
734     fHOutCentZNAvsCentV0A= new TH2F("fHOutCentZNAvsCentV0A","fHOutCentZNAvsCentV0A; Cent ZNA; Cent V0A;", 505,0,101,505,0,101);
735     fHOutCentZNAvsCentV0C= new TH2F("fHOutCentZNAvsCentV0C","fHOutCentZNAvsCentV0C; Cent ZNA; Cent V0C;", 505,0,101,505,0,101);
736     fHOutCentZNAvsCentTRK= new TH2F("fHOutCentZNAvsCentTRK","fHOutCentZNAvsCentTRK; Cent ZNA; Cent TRK;", 505,0,101,505,0,101);
737     fHOutCentZNAvsCentCND= new TH2F("fHOutCentZNAvsCentCND","fHOutCentZNAvsCentCND; Cent ZNA; Cent CND;", 505,0,101,505,0,101);
738     fHOutCentZNAvsCentCL1= new TH2F("fHOutCentZNAvsCentCL1","fHOutCentZNAvsCentCL1; Cent ZNA; Cent CL1;", 505,0,101,505,0,101);
739
740     fHOutMultV0AC = new TH2F("fHOutMultV0AC","fHOutMultV0AC; Multiplicity V0A; Multiplicity V0C",1000,0,1000,1000,0,1000);
741     fHOutMultV0M  = new TH1F("fHOutMultV0M","fHOutMultV0M; Multiplicity V0",25000,0,25000);
742     fHOutMultV0A  = new TH1F("fHOutMultV0A","fHOutMultV0A; Multiplicity V0",25000,0,25000);
743     fHOutMultV0C  = new TH1F("fHOutMultV0C","fHOutMultV0C; Multiplicity V0",25000,0,25000);
744     fHOutMultV0Mnc= new TH1F("fHOutMultV0Mnc","fHOutMultV0Mnc; Multiplicity V0",25000,0,25000);
745     fHOutMultV0Anc= new TH1F("fHOutMultV0Anc","fHOutMultV0Anc; Multiplicity V0",25000,0,25000);
746     fHOutMultV0Cnc= new TH1F("fHOutMultV0Cnc","fHOutMultV0Cnc; Multiplicity V0",25000,0,25000);
747     fHOutMultV0O  = new TH1F("fHOutMultV0O","fHOutMultV0O; Multiplicity V0",40000,0,40000);
748     fHOutMultV0Cells = new TH2F("fHOutMultV0Cells","fHOutMultV0Cells",33,-0.5,32.5,33,-0.5,32.5); 
749     fHOutMultFMD = new TH1F("fHOutMultFMD","fHOutMultFMD; Multiplicity FMD",24000,0,24000);
750     fHOutMultTRK = new TH1F("fHOutMultTRK","fHOutMultTRK; Multiplicity TPC",4000,0,4000);
751     fHOutMultTKL = new TH1F("fHOutMultTKL","fHOutMultTKL; Multiplicity tracklets",5000,0,5000);
752     fHOutMultCL0 = new TH1F("fHOutMultCL0","fHOutMultCL0; Multiplicity SPD inner",7000,0,7000);
753     fHOutMultCL1 = new TH1F("fHOutMultCL1","fHOutMultCL1; Multiplicity SPD outer",7000,0,7000);
754     fHOutMultCND = new TH1F("fHOutMultCND","fHOutMultCND; Multiplicity candle",4000,0,4000);
755     fHOutMultNPA = new TH1F("fHOutMultNPA","fHOutMultNPA; Nparticipants",450,0,450);
756     fHOutMultZNA = new TH1F("fHOutMultZNA","fHOutMultZNA; ZNA Energy",500,0,2000);
757
758     fHOutMultV0MvsZDN = new TH2F("fHOutMultV0MvsZDN","fHOutMultV0MvsZDN; Multiplicity V0; Energy ZDC-N",500,0,30000,500,0,180000);
759     fHOutMultZEMvsZDN = new TH2F("fHOutMultZEMvsZDN","fHOutMultZEMvsZDN; Energy ZEM; Energy ZDC-N",500,0,2500,500,0,180000);
760     fHOutMultV0MvsZDC = new TH2F("fHOutMultV0MvsZDC","fHOutMultV0MvsZDC; Multiplicity V0; Energy ZDC",500,0,30000,500,0,200000);
761     fHOutMultZEMvsZDC = new TH2F("fHOutMultZEMvsZDC","fHOutMultZEMvsZDC; Energy ZEM; Energy ZDC",500,0,2500,500,0,200000);
762     fHOutMultZEMvsZDCw = new TH2F("fHOutMultZEMvsZDCw","fHOutMultZEMvsZDCw; Energy ZEM; Energy ZDC (weigthed with V0 percentile)",500,0,2500,500,0,200000);
763     fHOutMultV0MvsCL1 = new TH2F("fHOutMultV0MvsCL1","fHOutMultV0MvsCL1; Multiplicity V0; Multiplicity SPD outer",2500,0,30000,700,0,7000);
764     fHOutMultV0MvsTRK = new TH2F("fHOutMultV0MvsTRK","fHOutMultV0MvsTRK; Multiplicity V0; Multiplicity TPC",2500,0,30000,400,0,4000);
765     fHOutMultTRKvsCL1 = new TH2F("fHOutMultTRKvsCL1","fHOutMultTRKvsCL1; Multiplicity TPC; Multiplicity SPD outer",400,0,4000,700,0,7000);
766     fHOutMultV0MvsV0O = new TH2F("fHOutMultV0MvsV0O","fHOutMultV0MvsV0O; Multiplicity V0; Multiplicity V0 Online",500,0,30000,500,0,40000);
767     fHOutMultV0OvsCL1 = new TH2F("fHOutMultV0OvsCL1","fHOutMultV0OvsCL1; Multiplicity V0; Multiplicity SPD outer",500,0,40000,700,0,7000);
768     fHOutMultV0OvsTRK = new TH2F("fHOutMultV0OvsTRK","fHOutMultV0OvsTRK; Multiplicity V0; Multiplicity TPC",500,0,40000,400,0,4000);
769     fHOutMultV0MvsV0O = new TH2F("fHOutMultV0MvsV0O","fHOutMultV0MvsV0O; Multiplicity V0; Multiplicity V0 Online",500,0,30000,500,0,30000);
770     fHOutMultV0OvsCL1 = new TH2F("fHOutMultV0OvsCL1","fHOutMultV0OvsCL1; Multiplicity V0; Multiplicity SPD outer",2500,0,30000,700,0,7000);
771     fHOutMultV0OvsTRK = new TH2F("fHOutMultV0OvsTRK","fHOutMultV0OvsTRK; Multiplicity V0; Multiplicity TPC",2500,0,30000,400,0,4000);
772     fHOutMultCL1vsTKL = new TH2F ("fHOutMultCL1vsTKL","fHOutMultCL1vsTKL; Multiplicity SPD outer; Multiplicity tracklets",700,0,7000,700,0,7000);
773     
774     fHOutCentV0Mqual1 = new TH1F("fHOutCentV0M_qual1","fHOutCentV0M_qual1; Centrality V0",505,0,101);
775     fHOutCentTRKqual1 = new TH1F("fHOutCentTRK_qual1","fHOutCentTRK_qual1; Centrality TPC",505,0,101);
776     fHOutCentCL1qual1 = new TH1F("fHOutCentCL1_qual1","fHOutCentCL1_qual1; Centrality SPD outer",505,0,101);
777     fHOutMultV0MvsCL1qual1 = new TH2F("fHOutMultV0MvsCL1_qual1","fHOutMultV0MvsCL1_qual1; Multiplicity V0; Multiplicity SPD outer",2500,0,25000,700,0,7000);
778     fHOutMultV0MvsTRKqual1 = new TH2F("fHOutMultV0MvsTRK_qual1","fHOutMultV0MvsTRK_qual1; Multiplicity V0; Multiplicity TPC",2500,0,25000,400,0,4000);
779     fHOutMultTRKvsCL1qual1 = new TH2F("fHOutMultTRKvsCL1_qual1","fHOutMultTRKvsCL1_qual1; Multiplicity TPC; Multiplicity SPD outer",400,0,4000,700,0,7000);
780     
781     fHOutCentV0Mqual2 = new TH1F("fHOutCentV0M_qual2","fHOutCentV0M_qual2; Centrality V0",505,0,101);
782     fHOutCentTRKqual2 = new TH1F("fHOutCentTRK_qual2","fHOutCentTRK_qual2; Centrality TPC",505,0,101);
783     fHOutCentCL1qual2 = new TH1F("fHOutCentCL1_qual2","fHOutCentCL1_qual2; Centrality SPD outer",505,0,101);
784     fHOutMultV0MvsCL1qual2 = new TH2F("fHOutMultV0MvsCL1_qual2","fHOutMultV0MvsCL1_qual2; Multiplicity V0; Multiplicity SPD outer",2500,0,25000,700,0,7000);
785     fHOutMultV0MvsTRKqual2 = new TH2F("fHOutMultV0MvsTRK_qual2","fHOutMultV0MvsTRK_qual2; Multiplicity V0; Multiplicity TPC",2500,0,25000,400,0,4000);
786     fHOutMultTRKvsCL1qual2 = new TH2F("fHOutMultTRKvsCL1_qual2","fHOutMultTRKvsCL1_qual2; Multiplicity TPC; Multiplicity SPD outer",400,0,4000,700,0,7000);
787     
788     fHOutQuality = new TH1F("fHOutQuality", "fHOutQuality", 100,-0.5,99.5);
789     fHOutVertex  = new TH1F("fHOutVertex", "fHOutVertex", 100,-20,20);
790     fHOutVertexT0  = new TH1F("fHOutVertexT0", "fHOutVertexT0", 100,-20,20);
791   
792     fOutputList->Add(  fHOutCentV0M    );
793     fOutputList->Add(  fHOutCentV0A    );
794     fOutputList->Add(  fHOutCentV0C    );
795     fOutputList->Add(  fHOutCentV0MCVHN);
796     fOutputList->Add(  fHOutCentV0MCVLN);
797     fOutputList->Add(  fHOutCentV0MCVHNinMB);
798     fOutputList->Add(  fHOutCentV0MCVLNinMB);
799     fOutputList->Add(  fHOutCentV0MCCENT);
800     fOutputList->Add(  fHOutCentV0MCSEMI);
801     fOutputList->Add(  fHOutCentV0MCCENTinMB);
802     fOutputList->Add(  fHOutCentV0MCSEMIinMB);
803     fOutputList->Add(  fHOutCentV0MMSL    );
804     fOutputList->Add(  fHOutCentV0MMSH    );
805     fOutputList->Add(  fHOutCentV0MMUL    );
806     fOutputList->Add(  fHOutCentV0MMLL    );
807     fOutputList->Add(  fHOutCentV0MEJE    );
808     fOutputList->Add(  fHOutCentV0MEGA    );
809     fOutputList->Add(  fHOutCentV0MPHS   );
810     fOutputList->Add(  fHOutCentV0MMSLinMB);
811     fOutputList->Add(  fHOutCentV0MMSHinMB);
812     fOutputList->Add(  fHOutCentV0MMULinMB);
813     fOutputList->Add(  fHOutCentV0MMLLinMB);
814     fOutputList->Add(  fHOutCentV0MEJEinMB);
815     fOutputList->Add(  fHOutCentV0MEGAinMB);
816     fOutputList->Add(  fHOutCentV0MPHSinMB);
817     fOutputList->Add(  fHOutCentFMD     );
818     fOutputList->Add(  fHOutCentTRK     );
819     fOutputList->Add(  fHOutCentTKL     );
820     fOutputList->Add(  fHOutCentCL0     );
821     fOutputList->Add(  fHOutCentCL1     );
822     fOutputList->Add(  fHOutCentCND     );
823     fOutputList->Add(  fHOutCentNPA     );
824     fOutputList->Add(  fHOutCentZNA     );
825     fOutputList->Add(  fHOutCentV0MvsFMD);
826     fOutputList->Add(  fHOutCentTKLvsV0M);
827     fOutputList->Add(  fHOutCentZEMvsZDC);
828     fOutputList->Add(  fHOutCentV0MvsCentCL1);
829     fOutputList->Add(  fHOutCentV0MvsCentTRK);
830     fOutputList->Add(  fHOutCentTRKvsCentCL1);
831     fOutputList->Add(  fHOutCentV0MvsCentZDC);
832     fOutputList->Add(  fHOutCentV0AvsCentV0C);
833     fOutputList->Add(  fHOutCentV0AvsCentTRK);
834     fOutputList->Add(  fHOutCentV0AvsCentCND);
835     fOutputList->Add(  fHOutCentV0AvsCentCL1);
836     fOutputList->Add(  fHOutCentV0CvsCentTRK);
837     fOutputList->Add(  fHOutCentV0CvsCentCND);
838     fOutputList->Add(  fHOutCentV0CvsCentCL1);
839     fOutputList->Add(  fHOutCentNPAvsCentV0A);
840     fOutputList->Add(  fHOutCentNPAvsCentV0C);
841     fOutputList->Add(  fHOutCentNPAvsCentTRK);
842     fOutputList->Add(  fHOutCentNPAvsCentCND);
843     fOutputList->Add(  fHOutCentNPAvsCentCL1);
844     fOutputList->Add(  fHOutCentZNAvsCentV0A);
845     fOutputList->Add(  fHOutCentZNAvsCentV0C);
846     fOutputList->Add(  fHOutCentZNAvsCentTRK);
847     fOutputList->Add(  fHOutCentZNAvsCentCND);
848     fOutputList->Add(  fHOutCentZNAvsCentCL1);
849
850     fOutputList->Add(  fHOutMultV0AC); 
851     fOutputList->Add(  fHOutMultV0M); 
852     fOutputList->Add(  fHOutMultV0A); 
853     fOutputList->Add(  fHOutMultV0C); 
854     fOutputList->Add(  fHOutMultV0Mnc); 
855     fOutputList->Add(  fHOutMultV0Anc); 
856     fOutputList->Add(  fHOutMultV0Cnc); 
857     fOutputList->Add(  fHOutMultV0O);
858     fOutputList->Add(  fHOutMultV0Cells) ;   
859     fOutputList->Add(  fHOutMultFMD); 
860     fOutputList->Add(  fHOutMultTRK); 
861     fOutputList->Add(  fHOutMultTKL); 
862     fOutputList->Add(  fHOutMultCL0); 
863     fOutputList->Add(  fHOutMultCL1); 
864     fOutputList->Add(  fHOutMultCND); 
865     fOutputList->Add(  fHOutMultNPA); 
866     fOutputList->Add(  fHOutMultZNA); 
867     fOutputList->Add(  fHOutMultV0MvsZDN);
868     fOutputList->Add(  fHOutMultZEMvsZDN);
869     fOutputList->Add(  fHOutMultV0MvsZDC);
870     fOutputList->Add(  fHOutMultZEMvsZDC);
871     fOutputList->Add(  fHOutMultZEMvsZDCw);
872     fOutputList->Add(  fHOutMultV0MvsCL1);
873     fOutputList->Add(  fHOutMultV0MvsTRK);
874     fOutputList->Add(  fHOutMultTRKvsCL1);
875     fOutputList->Add(  fHOutMultV0MvsV0O);
876     fOutputList->Add(  fHOutMultV0OvsCL1);
877     fOutputList->Add(  fHOutMultV0OvsTRK);
878     fOutputList->Add(  fHOutMultCL1vsTKL);
879     fOutputList->Add(  fHOutCentV0Mqual1 );
880     fOutputList->Add(  fHOutCentTRKqual1 );
881     fOutputList->Add(  fHOutCentCL1qual1 );                   
882     fOutputList->Add(  fHOutMultV0MvsCL1qual1);
883     fOutputList->Add(  fHOutMultV0MvsTRKqual1);
884     fOutputList->Add(  fHOutMultTRKvsCL1qual1);
885     fOutputList->Add(  fHOutCentV0Mqual2 );
886     fOutputList->Add(  fHOutCentTRKqual2 );
887     fOutputList->Add(  fHOutCentCL1qual2 );
888     fOutputList->Add(  fHOutMultV0MvsCL1qual2);
889     fOutputList->Add(  fHOutMultV0MvsTRKqual2);
890     fOutputList->Add(  fHOutMultTRKvsCL1qual2);
891     fOutputList->Add(  fHOutQuality );
892     fOutputList->Add(  fHOutVertex );
893     fOutputList->Add(  fHOutVertexT0 );
894   
895     PostData(1, fOutputList); 
896   }
897   
898   fTrackCuts = AliESDtrackCuts::GetStandardTPCOnlyTrackCuts();
899   fEsdTrackCuts = AliESDtrackCuts::GetStandardITSTPCTrackCuts2011();
900   fEsdTrackCuts->SetClusterRequirementITS(AliESDtrackCuts::kSPD,AliESDtrackCuts::kOff);
901   // Add SPD requirement
902   fEsdTrackCutsExtra1 = new AliESDtrackCuts("SPD", "Require 1 cluster in SPD");
903   fEsdTrackCutsExtra1->SetClusterRequirementITS(AliESDtrackCuts::kSPD,AliESDtrackCuts::kAny);
904   // Add SDD requirement
905   fEsdTrackCutsExtra2 = new AliESDtrackCuts("SDD", "Require 1 cluster in first layer SDD");
906   fEsdTrackCutsExtra2->SetClusterRequirementITS(AliESDtrackCuts::kSDD,AliESDtrackCuts::kFirst);  
907 }
908
909 //________________________________________________________________________
910 void AliCentralitySelectionTask::UserExec(Option_t */*option*/)
911
912   // Execute analysis for current event:
913   if(fDebug>1) printf(" **** AliCentralitySelectionTask::UserExec() \n");
914
915   Float_t  zncEnergy = 0.;          //  ZNC Energy
916   Float_t  zpcEnergy = 0.;          //  ZPC Energy
917   Float_t  znaEnergy = 0.;          //  ZNA Energy
918   Float_t  zpaEnergy = 0.;          //  ZPA Energy
919   Float_t  zem1Energy = 0.;         //  ZEM1 Energy
920   Float_t  zem2Energy = 0.;         //  ZEM2 Energy
921   Bool_t   zdcEnergyCal = kFALSE;   // if zdc is calibrated (in pass2)
922   Double_t znaTower = 0.;           // common PMT of ZNA 
923   Bool_t   znaFired = kFALSE;
924
925   Int_t    nTracks = 0;             //  no. tracks
926   Int_t    nTracklets = 0;          //  no. tracklets
927   Int_t    nClusters[6] = {0};      //  no. clusters on 6 ITS layers
928   Int_t    nChips[2];               //  no. chips on 2 SPD layers
929   Float_t  spdCorr =0;              //  corrected spd2 multiplicity
930   Int_t    multCND = 0;             //  no. tracks (candle condition)
931
932   Float_t  multV0A  = 0;            //  multiplicity from V0 reco side A
933   Float_t  multV0C  = 0;            //  multiplicity from V0 reco side C
934   Float_t  multV0ACorr  = 0;            //  multiplicity from V0 reco side A
935   Float_t  multV0CCorr  = 0;            //  multiplicity from V0 reco side C
936   Short_t  multV0AOnline  = 0;      //  multiplicity from V0 reco side A
937   Short_t  multV0COnline  = 0;      //  multiplicity from V0 reco side C
938   Float_t  v0Corr = 0;              // corrected V0 multiplicity (used for MC)
939   Int_t nV0A = 0;
940   Int_t nV0C = 0;
941
942   Float_t  multFMDA = 0;            //  multiplicity from FMD on detector A
943   Float_t  multFMDC = 0;            //  multiplicity from FMD on detector C
944
945   Float_t zvtx =0;                  // z-vertex SPD
946   Int_t zvtxNcont =0;               // contributors to z-vertex SPD
947
948   Float_t zvtxT0 =0;                // z-vertex T0
949
950   Int_t Npart =0;                   // N. of participants (true MC)
951
952   AliCentrality *esdCent = 0;
953
954   if(fAnalysisInput.CompareTo("ESD")==0){
955   
956     AliVEvent* event = InputEvent();
957     AliESDEvent* esd = dynamic_cast<AliESDEvent*>(event);
958     if (!esd) {
959       AliError("No ESD Event");
960       return;
961     }
962   
963     LoadBranches();
964     
965     if (SetupRun(esd)<0) {
966       AliError("Centrality File not available for this run");
967       return;
968     }
969     
970     esdCent = esd->GetCentrality();
971
972     // ***** Vertex Info
973     const AliESDVertex* vtxESD = esd->GetPrimaryVertexSPD();
974     zvtx        = vtxESD->GetZ(); 
975     zvtxNcont   = vtxESD->GetNContributors();
976
977     // ***** V0 info    
978     AliESDVZERO* esdV0 = esd->GetVZEROData();
979     if (!esdV0)
980       {
981         AliError("AliESDVZERO not available");
982         return;
983       }
984     multV0A=esdV0->GetMTotV0A();
985     multV0C=esdV0->GetMTotV0C();
986
987     multV0ACorr = AliESDUtils::GetCorrV0A(multV0A,zvtx);    
988     multV0CCorr = AliESDUtils::GetCorrV0C(multV0C,zvtx);    
989
990     v0Corr = multV0A+multV0C;
991
992     multV0AOnline=esdV0->GetTriggerChargeA(); 
993     multV0COnline=esdV0->GetTriggerChargeC(); 
994
995     // Count V0 flags
996     for(Int_t i = 0; i < 32; ++i) {
997       if (esdV0->GetBBFlag(i)) nV0C++;
998       if (esdV0->GetBBFlag(i+32)) nV0A++;
999     }
1000     
1001     // ***** T0 info    
1002     const AliESDTZERO* esdT0 = esd->GetESDTZERO();
1003     if (!esdT0)
1004       {
1005         AliError("AliESDTZERO not available");
1006         return;
1007       }
1008     Int_t trig=esdT0->GetT0Trig();
1009     Bool_t kT0BB = kFALSE;    
1010     if(trig&1) kT0BB=kTRUE;
1011     zvtxT0=esdT0->GetT0zVertex();
1012
1013
1014     // ***** Trigger info    
1015     fIsSelected = ((esdV0->GetV0ADecision()==1) && (esdV0->GetV0CDecision()==1));
1016     TString trigStr(esd->GetFiredTriggerClasses());
1017     
1018     fMB=kFALSE;
1019     fCVHN=kFALSE; fCVLN=kFALSE; fCCENT=kFALSE; fCSEMI=kFALSE; 
1020     fMSL=kFALSE;  fMSH=kFALSE;  fMUL=kFALSE;   fMLL=kFALSE;
1021     fEJE=kFALSE;  fEGA=kFALSE;  fPHS=kFALSE;
1022
1023     if ( (trigStr.Contains("-B-")) &&  (trigStr.Contains("CPBI")) && (fIsSelected)) 
1024       fMB=kTRUE;
1025     if ( (trigStr.Contains("-B-")) &&  (trigStr.Contains("CVHN")) && (fIsSelected)) 
1026       fCVHN=kTRUE;
1027     if ( (trigStr.Contains("-B-")) &&  (trigStr.Contains("CVLN")) && (fIsSelected))
1028       fCVLN=kTRUE;
1029     if ( (trigStr.Contains("-B-")) &&  (trigStr.Contains("CCENT")) && (fIsSelected)) 
1030       fCCENT=kTRUE;
1031     if ( (trigStr.Contains("-B-")) &&  (trigStr.Contains("CSEMI")) && (fIsSelected))
1032       fCSEMI=kTRUE;
1033     
1034     if ( (trigStr.Contains("-B-")) &&  (trigStr.Contains("CPBI1MSL")) && (fIsSelected))
1035       fMSL=kTRUE;
1036     if ( (trigStr.Contains("-B-")) &&  (trigStr.Contains("CPBI1MSH")) && (fIsSelected))
1037       fMSH=kTRUE;
1038     if ( (trigStr.Contains("-B-")) &&  (trigStr.Contains("CPBI1MUL")) && (fIsSelected))
1039       fMUL=kTRUE;
1040     if ( (trigStr.Contains("-B-")) &&  (trigStr.Contains("CPBI1MLL")) && (fIsSelected))
1041       fMLL=kTRUE;
1042     if ( (trigStr.Contains("-B-")) &&  (trigStr.Contains("CPBI1EJE")) && (fIsSelected))
1043       fEJE=kTRUE;
1044     if ( (trigStr.Contains("-B-")) &&  (trigStr.Contains("CPBI1EGA")) && (fIsSelected))
1045       fEGA=kTRUE;
1046     if ( (trigStr.Contains("-B-")) &&  (trigStr.Contains("CPBI1PHS")) && (fIsSelected))
1047       fPHS=kTRUE;
1048
1049     fCVHNbit=kFALSE;    fCVLNbit=kFALSE;       fCCENTbit=kFALSE;    fCSEMIbit=kFALSE; 
1050     if (esdV0->GetTriggerBits() & (1<<8)) 
1051       fCVHNbit=kTRUE;
1052     if (esdV0->GetTriggerBits() & (1<<6)) 
1053       fCVLNbit=kTRUE;
1054     
1055     if (kT0BB && fCVHNbit)
1056       fCCENTbit=kTRUE;
1057     if (kT0BB && fCVLNbit)
1058       fCSEMIbit=kTRUE;
1059
1060     
1061     // ***** CB info (tracklets, clusters, chips)
1062     //nTracks    = event->GetNumberOfTracks();     
1063     nTracks    = fTrackCuts ? (Short_t)fTrackCuts->GetReferenceMultiplicity(esd,kTRUE):-1;
1064
1065     Short_t nTrTPCcandle = 0;
1066     for (Int_t iTracks = 0; iTracks < esd->GetNumberOfTracks(); iTracks++) {
1067
1068       AliESDtrack* track = esd->GetTrack(iTracks);
1069       if (!track) continue;
1070       
1071       if (! fEsdTrackCuts->IsSelected(track) )continue;
1072       
1073       if (fEsdTrackCutsExtra1 && fEsdTrackCutsExtra2 &&
1074           !fEsdTrackCutsExtra1->IsSelected(track) &&
1075           !fEsdTrackCutsExtra2->IsSelected(track)) continue;
1076       
1077       if (track->Pt() > 0.4 && TMath::Abs(track->Eta()) < 0.9)  nTrTPCcandle++;
1078     } 
1079     multCND = nTrTPCcandle;
1080
1081     const AliMultiplicity *mult = esd->GetMultiplicity();
1082     nTracklets = mult->GetNumberOfTracklets();
1083
1084     for(Int_t ilay=0; ilay<6; ilay++){
1085       nClusters[ilay] = mult->GetNumberOfITSClusters(ilay);
1086     }
1087
1088     for(Int_t ilay=0; ilay<2; ilay++){
1089       nChips[ilay] = mult->GetNumberOfFiredChips(ilay);
1090     }
1091
1092     spdCorr = AliESDUtils::GetCorrSPD2(nClusters[1],zvtx);    
1093
1094     // ***** FMD info
1095     AliESDFMD *fmd = esd->GetFMDData();
1096     Float_t totalMultA = 0;
1097     Float_t totalMultC = 0;
1098     const Float_t fFMDLowCut = 0.4;
1099   
1100     for(UShort_t det=1;det<=3;det++) {
1101       Int_t nRings = (det==1 ? 1 : 2);
1102       for (UShort_t ir = 0; ir < nRings; ir++) {          
1103         Char_t   ring = (ir == 0 ? 'I' : 'O');
1104         UShort_t nsec = (ir == 0 ? 20  : 40);
1105         UShort_t nstr = (ir == 0 ? 512 : 256);
1106         for(UShort_t sec =0; sec < nsec;  sec++)  {
1107           for(UShort_t strip = 0; strip < nstr; strip++) {
1108             
1109             Float_t fmdMult = fmd->Multiplicity(det,ring,sec,strip);
1110             if(fmdMult == 0 || fmdMult == AliESDFMD::kInvalidMult) continue;
1111             
1112             Float_t nParticles=0;
1113             
1114             if(fmdMult > fFMDLowCut) {
1115               nParticles = 1.;
1116             }
1117             
1118             if (det<3) totalMultA = totalMultA + nParticles;
1119             else totalMultC = totalMultC + nParticles;
1120             
1121           }
1122         }
1123       }
1124     }
1125     multFMDA = totalMultA;
1126     multFMDC = totalMultC;
1127   
1128     // ***** ZDC info
1129     AliESDZDC *esdZDC = esd->GetESDZDC();
1130     zdcEnergyCal = esdZDC->AliESDZDC::TestBit(AliESDZDC::kEnergyCalibratedSignal);
1131     if (zdcEnergyCal) {      
1132       zncEnergy = (Float_t) (esdZDC->GetZDCN1Energy());
1133       zpcEnergy = (Float_t) (esdZDC->GetZDCP1Energy());
1134       znaEnergy = (Float_t) (esdZDC->GetZDCN2Energy());
1135       zpaEnergy = (Float_t) (esdZDC->GetZDCP2Energy());
1136     } else {
1137       zncEnergy = (Float_t) (esdZDC->GetZDCN1Energy())/8.;
1138       zpcEnergy = (Float_t) (esdZDC->GetZDCP1Energy())/8.;
1139       znaEnergy = (Float_t) (esdZDC->GetZDCN2Energy())/8.;
1140       zpaEnergy = (Float_t) (esdZDC->GetZDCP2Energy())/8.;
1141     }
1142     zem1Energy = (Float_t) (esdZDC->GetZDCEMEnergy(0))/8.;
1143     zem2Energy = (Float_t) (esdZDC->GetZDCEMEnergy(1))/8.;
1144
1145     const Double_t *ZNAtower = esdZDC->GetZN2TowerEnergy(); 
1146     znaTower = ZNAtower[0];
1147
1148     for (Int_t j = 0; j < 4; ++j) {
1149       if (esdZDC->GetZDCTDCData(12,j) != 0) {
1150         znaFired = kTRUE;
1151       }
1152     } 
1153   
1154     // ***** MC info
1155     AliAnalysisManager* anMan = AliAnalysisManager::GetAnalysisManager();
1156     AliMCEventHandler* eventHandler = (AliMCEventHandler*)anMan->GetMCtruthEventHandler();
1157     AliStack*    stack=0;
1158     AliMCEvent*  mcEvent=0;
1159     if (fIsMCInput && eventHandler && (mcEvent=eventHandler->MCEvent()) && (stack=mcEvent->Stack())) {
1160       
1161       AliGenHijingEventHeader* hHijing=0;
1162       AliGenDPMjetEventHeader* dpmHeader=0;
1163       
1164       AliGenEventHeader* mcGenH = mcEvent->GenEventHeader();
1165       if (mcGenH->InheritsFrom(AliGenHijingEventHeader::Class())) {
1166         hHijing = (AliGenHijingEventHeader*)mcGenH;
1167         if(hHijing) Npart = hHijing->ProjectileParticipants()+hHijing->TargetParticipants();
1168       }
1169       else if (mcGenH->InheritsFrom(AliGenCocktailEventHeader::Class())) {
1170         TList* headers = ((AliGenCocktailEventHeader*)mcGenH)->GetHeaders();
1171         hHijing = dynamic_cast<AliGenHijingEventHeader*>(headers->FindObject("Hijing"));
1172         if(hHijing) Npart = hHijing->ProjectileParticipants()+hHijing->TargetParticipants();
1173       }
1174       else if (mcGenH->InheritsFrom(AliGenDPMjetEventHeader::Class())) {
1175         dpmHeader = (AliGenDPMjetEventHeader*)mcGenH;
1176         if(dpmHeader) Npart = dpmHeader->ProjectileParticipants()+ dpmHeader->TargetParticipants();
1177       }
1178     }
1179   }   
1180
1181   else if(fAnalysisInput.CompareTo("AOD")==0){
1182     //AliAODEvent *aod =  dynamic_cast<AliAODEvent*> (InputEvent());
1183     // to be implemented
1184     printf("  AOD analysis not yet implemented!!!\n\n");
1185     return;
1186   } 
1187
1188   // ***** Scaling for MC
1189   if (fIsMCInput) {
1190     fUseScaling=kFALSE;
1191     v0Corr  = Short_t((multV0A+multV0C)  * fV0MScaleFactorMC);
1192     multV0A  = multV0A * fV0MScaleFactorMC;
1193     multV0C  = multV0C * fV0MScaleFactorMC;
1194   }
1195   // ***** Scaling for Data 
1196   if (fUseScaling) {
1197     v0Corr  = Short_t(v0Corr /   fV0MScaleFactor);
1198     spdCorr = spdCorr / fSPDScaleFactor;
1199     nTracks = Int_t(nTracks /   fTPCScaleFactor);
1200   }
1201
1202   // ***** Centrality Selection
1203   if(fHtempV0M) fCentV0M = fHtempV0M->GetBinContent(fHtempV0M->FindBin((v0Corr)));
1204   if(fHtempV0A) fCentV0A = fHtempV0A->GetBinContent(fHtempV0A->FindBin((multV0ACorr)));
1205   if(fHtempV0C) fCentV0C = fHtempV0C->GetBinContent(fHtempV0C->FindBin((multV0CCorr)));
1206   if(fHtempFMD) fCentFMD = fHtempFMD->GetBinContent(fHtempFMD->FindBin((multFMDA+multFMDC)));
1207   if(fHtempTRK) fCentTRK = fHtempTRK->GetBinContent(fHtempTRK->FindBin(nTracks));
1208   if(fHtempTKL) fCentTKL = fHtempTKL->GetBinContent(fHtempTKL->FindBin(nTracklets));
1209   if(fHtempCL0) fCentCL0 = fHtempCL0->GetBinContent(fHtempCL0->FindBin(nClusters[0]));
1210   if(fHtempCL1) fCentCL1 = fHtempCL1->GetBinContent(fHtempCL1->FindBin(spdCorr));
1211   if(fHtempCND) fCentCND = fHtempCND->GetBinContent(fHtempCND->FindBin(multCND));
1212   if(fHtempNPA) fCentNPA = fHtempNPA->GetBinContent(fHtempNPA->FindBin(Npart));
1213   if(fHtempZNA) {
1214     if(znaFired) fCentZNA = fHtempZNA->GetBinContent(fHtempZNA->FindBin(znaTower));
1215     else fCentZNA = 101;
1216   }
1217   if(fHtempV0MvsFMD) fCentV0MvsFMD = fHtempV0MvsFMD->GetBinContent(fHtempV0MvsFMD->FindBin((multV0A+multV0C)));
1218   if(fHtempTKLvsV0M) fCentTKLvsV0M = fHtempTKLvsV0M->GetBinContent(fHtempTKLvsV0M->FindBin(nTracklets));
1219   if(fHtempZEMvsZDC) fCentZEMvsZDC = fHtempZEMvsZDC->GetBinContent(fHtempZEMvsZDC->FindBin(zem1Energy+zem2Energy,zncEnergy+znaEnergy+zpcEnergy+zpaEnergy));
1220
1221
1222   // ***** Cleaning
1223   if (fUseCleaning) {
1224     fQuality=0;
1225     
1226     // ***** vertex
1227     if (TMath::Abs(zvtx)>fZVCut || zvtxNcont<1) fQuality += 1;   
1228
1229     // ***** outliers, skip in case of MC input
1230     if (!fIsMCInput) {
1231       // **** V0 vs SPD
1232       if (IsOutlierV0MSPD(spdCorr, v0Corr, int(fCentV0M))) fQuality  += 2;
1233       // ***** V0 vs TPC
1234       if (IsOutlierV0MTPC(nTracks, v0Corr, int(fCentV0M))) fQuality  += 4;
1235       // ***** V0 vs ZDC
1236       if (IsOutlierV0MZDC((zncEnergy+znaEnergy+zpcEnergy+zpaEnergy), v0Corr) &&
1237           (zdcEnergyCal==kFALSE) ) fQuality  += 8;
1238       if (IsOutlierV0MZDCECal((zncEnergy+znaEnergy+zpcEnergy+zpaEnergy), v0Corr) &&
1239           (zdcEnergyCal==kTRUE) ) fQuality  += 8;
1240     }
1241   } else {
1242     fQuality = 0;
1243   }
1244
1245       
1246   if (esdCent) {
1247     esdCent->SetQuality(fQuality);
1248     esdCent->SetCentralityV0M(fCentV0M);
1249     esdCent->SetCentralityV0A(fCentV0A);
1250     esdCent->SetCentralityV0C(fCentV0C);
1251     esdCent->SetCentralityFMD(fCentFMD);
1252     esdCent->SetCentralityTRK(fCentTRK);
1253     esdCent->SetCentralityTKL(fCentTKL);
1254     esdCent->SetCentralityCL0(fCentCL0);
1255     esdCent->SetCentralityCL1(fCentCL1);
1256     esdCent->SetCentralityCND(fCentCND);
1257     esdCent->SetCentralityNPA(fCentNPA);
1258     esdCent->SetCentralityZNA(fCentZNA);
1259     esdCent->SetCentralityV0MvsFMD(fCentV0MvsFMD);
1260     esdCent->SetCentralityTKLvsV0M(fCentTKLvsV0M);
1261     esdCent->SetCentralityZEMvsZDC(fCentZEMvsZDC);
1262   }
1263
1264   // filling QA histograms
1265   if (fFillHistos) {    
1266     if ((fMB) && (abs(zvtx)<10))        fHOutMultCL1vsTKL->Fill(spdCorr,nTracklets);
1267
1268     if (fCVHN)   fHOutCentV0MCVHN->Fill(fCentV0M);
1269     if (fCVLN)   fHOutCentV0MCVLN->Fill(fCentV0M);
1270     if (fCCENT)  fHOutCentV0MCCENT->Fill(fCentV0M);
1271     if (fCSEMI)  fHOutCentV0MCSEMI->Fill(fCentV0M);
1272     if (fMSL) fHOutCentV0MMSL->Fill(fCentV0M);
1273     if (fMSH) fHOutCentV0MMSH->Fill(fCentV0M);
1274     if (fMUL) fHOutCentV0MMUL->Fill(fCentV0M);
1275     if (fMLL) fHOutCentV0MMLL->Fill(fCentV0M);
1276     if (fEJE) fHOutCentV0MEJE->Fill(fCentV0M);
1277     if (fEGA) fHOutCentV0MEGA->Fill(fCentV0M);
1278     if (fPHS) fHOutCentV0MPHS->Fill(fCentV0M);
1279
1280     if (((AliInputEventHandler*)(AliAnalysisManager::GetAnalysisManager()->GetInputEventHandler()))->IsEventSelected() & AliVEvent::kMB) { // fill the QA histograms only for MB events!
1281       fHOutQuality->Fill(fQuality);
1282       fHOutVertex->Fill(zvtx);
1283       fHOutVertexT0->Fill(zvtxT0);
1284       
1285       if (fQuality==0) {  
1286         fHOutCentV0M->Fill(fCentV0M);
1287         fHOutCentV0A->Fill(fCentV0A);
1288         fHOutCentV0C->Fill(fCentV0C);
1289         
1290         if (fCVHNbit)  fHOutCentV0MCVHNinMB->Fill(fCentV0M);
1291         if (fCVLNbit)  fHOutCentV0MCVLNinMB->Fill(fCentV0M);
1292         if (fCCENTbit) fHOutCentV0MCCENTinMB->Fill(fCentV0M);
1293         if (fCSEMIbit) fHOutCentV0MCSEMIinMB->Fill(fCentV0M);
1294         if (fMSL) fHOutCentV0MMSLinMB->Fill(fCentV0M);
1295         if (fMSH) fHOutCentV0MMSHinMB->Fill(fCentV0M);
1296         if (fMUL) fHOutCentV0MMULinMB->Fill(fCentV0M);
1297         if (fMLL) fHOutCentV0MMLLinMB->Fill(fCentV0M);
1298         if (fEJE) fHOutCentV0MEJEinMB->Fill(fCentV0M);
1299         if (fEGA) fHOutCentV0MEGAinMB->Fill(fCentV0M);
1300         if (fPHS) fHOutCentV0MPHSinMB->Fill(fCentV0M);
1301         
1302         fHOutCentFMD->Fill(fCentFMD);
1303         fHOutCentTRK->Fill(fCentTRK);
1304         fHOutCentTKL->Fill(fCentTKL);
1305         fHOutCentCL0->Fill(fCentCL0);
1306         fHOutCentCL1->Fill(fCentCL1);
1307         fHOutCentCND->Fill(fCentCND);
1308         fHOutCentNPA->Fill(fCentNPA);
1309         fHOutCentZNA->Fill(fCentZNA);
1310         fHOutCentV0MvsFMD->Fill(fCentV0MvsFMD);
1311         fHOutCentTKLvsV0M->Fill(fCentTKLvsV0M);
1312         fHOutCentZEMvsZDC->Fill(fCentZEMvsZDC);
1313         fHOutCentV0MvsCentCL1->Fill(fCentV0M,fCentCL1);
1314         fHOutCentV0MvsCentTRK->Fill(fCentV0M,fCentTRK);
1315         fHOutCentTRKvsCentCL1->Fill(fCentTRK,fCentCL1);
1316         fHOutCentV0MvsCentZDC->Fill(fCentV0M,fCentZEMvsZDC);
1317         fHOutCentV0AvsCentV0C->Fill(fCentV0A,fCentV0C);
1318         fHOutCentV0AvsCentTRK->Fill(fCentV0A,fCentTRK);
1319         fHOutCentV0AvsCentCND->Fill(fCentV0A,fCentCND);
1320         fHOutCentV0AvsCentCL1->Fill(fCentV0A,fCentCL1);
1321         fHOutCentV0CvsCentTRK->Fill(fCentV0C,fCentTRK);
1322         fHOutCentV0CvsCentCND->Fill(fCentV0C,fCentCND);
1323         fHOutCentV0CvsCentCL1->Fill(fCentV0C,fCentCL1);
1324         fHOutCentNPAvsCentV0A->Fill(fCentNPA,fCentV0A);
1325         fHOutCentNPAvsCentV0C->Fill(fCentNPA,fCentV0C);
1326         fHOutCentNPAvsCentTRK->Fill(fCentNPA,fCentTRK);
1327         fHOutCentNPAvsCentCND->Fill(fCentNPA,fCentCND);
1328         fHOutCentNPAvsCentCL1->Fill(fCentNPA,fCentCL1);
1329         fHOutCentZNAvsCentV0A->Fill(fCentZNA,fCentV0A);
1330         fHOutCentZNAvsCentV0C->Fill(fCentZNA,fCentV0C);
1331         fHOutCentZNAvsCentTRK->Fill(fCentZNA,fCentTRK);
1332         fHOutCentZNAvsCentCND->Fill(fCentZNA,fCentCND);
1333         fHOutCentZNAvsCentCL1->Fill(fCentZNA,fCentCL1);
1334
1335         fHOutMultV0AC->Fill(multV0A,multV0C);
1336         fHOutMultV0M->Fill(multV0ACorr+multV0CCorr);
1337         fHOutMultV0A->Fill(multV0ACorr);
1338         fHOutMultV0C->Fill(multV0CCorr);
1339         fHOutMultV0Mnc->Fill(multV0A+multV0C);
1340         fHOutMultV0Anc->Fill(multV0A);
1341         fHOutMultV0Cnc->Fill(multV0C);
1342         fHOutMultV0O->Fill(multV0AOnline+multV0COnline);
1343         fHOutMultV0Cells->Fill(nV0A,nV0C); 
1344         fHOutMultFMD->Fill(multFMDA+multFMDC);
1345         fHOutMultTRK->Fill(nTracks);
1346         fHOutMultTKL->Fill(nTracklets);
1347         fHOutMultCL0->Fill(nClusters[0]);
1348         fHOutMultCL1->Fill(spdCorr);
1349         fHOutMultCND->Fill(multCND);
1350         fHOutMultNPA->Fill(Npart);
1351         fHOutMultZNA->Fill(znaTower);
1352
1353         fHOutMultV0MvsZDN->Fill(v0Corr,(zncEnergy+znaEnergy));
1354         fHOutMultZEMvsZDN->Fill((zem1Energy+zem2Energy),(zncEnergy+znaEnergy));
1355         fHOutMultV0MvsZDC->Fill(v0Corr,(zncEnergy+znaEnergy+zpcEnergy+zpaEnergy));
1356         fHOutMultZEMvsZDC->Fill((zem1Energy+zem2Energy),(zncEnergy+znaEnergy+zpcEnergy+zpaEnergy));
1357         fHOutMultZEMvsZDCw->Fill((zem1Energy+zem2Energy),(zncEnergy+znaEnergy+zpcEnergy+zpaEnergy),fCentV0M);
1358         fHOutMultV0MvsCL1->Fill(v0Corr,spdCorr);
1359         fHOutMultV0MvsTRK->Fill(v0Corr,nTracks);
1360         fHOutMultTRKvsCL1->Fill(nTracks,spdCorr);
1361         fHOutMultV0MvsV0O->Fill(v0Corr,(multV0AOnline+multV0COnline));
1362         fHOutMultV0OvsCL1->Fill((multV0AOnline+multV0COnline),spdCorr);
1363         fHOutMultV0OvsTRK->Fill((multV0AOnline+multV0COnline),nTracks);
1364       } else if (fQuality%2 == 0) {
1365         fHOutCentV0Mqual1->Fill(fCentV0M);
1366         fHOutCentTRKqual1->Fill(fCentTRK);
1367         fHOutCentCL1qual1->Fill(fCentCL1);
1368         fHOutMultV0MvsCL1qual1->Fill(v0Corr,spdCorr);
1369         fHOutMultV0MvsTRKqual1->Fill(v0Corr,nTracks);
1370         fHOutMultTRKvsCL1qual1->Fill(nTracks,spdCorr);
1371       } else {
1372         fHOutCentV0Mqual2->Fill(fCentV0M);
1373         fHOutCentTRKqual2->Fill(fCentTRK);
1374         fHOutCentCL1qual2->Fill(fCentCL1);
1375         fHOutMultV0MvsCL1qual2->Fill(v0Corr,spdCorr);
1376         fHOutMultV0MvsTRKqual2->Fill(v0Corr,nTracks);
1377         fHOutMultTRKvsCL1qual2->Fill(nTracks,spdCorr);
1378       }
1379     }
1380     PostData(1, fOutputList); 
1381   }
1382 }
1383 //________________________________________________________________________
1384 void AliCentralitySelectionTask::Terminate(Option_t */*option*/)
1385 {
1386   // Terminate analysis
1387 }
1388 //________________________________________________________________________
1389 Int_t AliCentralitySelectionTask::SetupRun(const AliESDEvent* const esd)
1390 {
1391   // Setup files for run
1392
1393   if (!esd)
1394     return -1;
1395
1396   // check if something to be done
1397   if (fCurrentRun == esd->GetRunNumber())
1398     return 0;
1399   else
1400     fCurrentRun = esd->GetRunNumber();
1401
1402   TString fileName =(Form("%s/COMMON/CENTRALITY/data/centrality.root", AliAnalysisManager::GetOADBPath()));
1403   AliInfo(Form("Setup Centrality Selection for run %d with file %s\n",fCurrentRun,fileName.Data()));
1404
1405   AliOADBContainer *con = new AliOADBContainer("OADB");
1406   con->InitFromFile(fileName,"Centrality");
1407
1408   AliOADBCentrality*  centOADB = 0;
1409   centOADB = (AliOADBCentrality*)(con->GetObject(fCurrentRun));
1410   if (!centOADB) {
1411     AliWarning(Form("Centrality OADB does not exist for run %d, using Default \n",fCurrentRun ));
1412     centOADB  = (AliOADBCentrality*)(con->GetDefaultObject("oadbDefault"));
1413   }
1414
1415   // modes
1416   fUseScaling   = centOADB->UseScaling();
1417   fUseCleaning  = centOADB->UseCleaning();
1418
1419   // cuts
1420   fZVCut        = centOADB->ZVCut();
1421   fOutliersCut  = centOADB->OutliersCut(); 
1422
1423   // centrality histos
1424   fHtempV0M       = centOADB->V0hist(); 
1425   fHtempV0A       = centOADB->V0Ahist(); 
1426   fHtempV0C       = centOADB->V0Chist(); 
1427   fHtempTRK       = centOADB->TPChist();
1428   fHtempCL1       = centOADB->SPDhist();
1429   fHtempCND       = centOADB->CNDhist();
1430   fHtempNPA       = centOADB->NPAhist();
1431   fHtempZNA       = centOADB->ZNAhist();
1432   fHtempZEMvsZDC  = centOADB->ZEMvsZDChist();
1433   
1434   TString path = gSystem->ExpandPathName(fileName.Data());
1435   if (!fHtempV0M) AliWarning(Form("Calibration for V0M does not exist in %s", path.Data()));
1436   if (!fHtempV0A) AliWarning(Form("Calibration for V0A does not exist in %s", path.Data()));
1437   if (!fHtempV0C) AliWarning(Form("Calibration for V0C does not exist in %s", path.Data()));
1438   if (!fHtempTRK) AliWarning(Form("Calibration for TRK does not exist in %s", path.Data()));
1439   if (!fHtempCL1) AliWarning(Form("Calibration for CL1 does not exist in %s", path.Data()));
1440   if (!fHtempCND) AliWarning(Form("Calibration for CND does not exist in %s", path.Data()));
1441   if (!fHtempNPA) AliWarning(Form("Calibration for NPA does not exist in %s", path.Data()));
1442   if (!fHtempZNA) AliWarning(Form("Calibration for ZNA does not exist in %s", path.Data()));
1443   if (!fHtempZEMvsZDC) AliWarning(Form("Calibration for ZEMvsZDC does not exist in %s", path.Data()));
1444
1445   // scale factors
1446   fV0MScaleFactor    = centOADB->V0MScaleFactor();
1447   fSPDScaleFactor    = centOADB->SPDScaleFactor();
1448   fTPCScaleFactor    = centOADB->TPCScaleFactor();
1449   fV0MScaleFactorMC  = centOADB->V0MScaleFactorMC();
1450
1451   // outliers parameters
1452   fV0MSPDOutlierPar0 = centOADB->V0MSPDOutlierPar0();      
1453   fV0MSPDOutlierPar1 = centOADB->V0MSPDOutlierPar1();     
1454   fV0MTPCOutlierPar0 = centOADB->V0MTPCOutlierPar0();      
1455   fV0MTPCOutlierPar1 = centOADB->V0MTPCOutlierPar1();     
1456                                                    
1457   fV0MSPDSigmaOutlierPar0 = centOADB->V0MSPDSigmaOutlierPar0(); 
1458   fV0MSPDSigmaOutlierPar1 = centOADB->V0MSPDSigmaOutlierPar1(); 
1459   fV0MSPDSigmaOutlierPar2 = centOADB->V0MSPDSigmaOutlierPar2();
1460   fV0MTPCSigmaOutlierPar0 = centOADB->V0MTPCSigmaOutlierPar0(); 
1461   fV0MTPCSigmaOutlierPar1 = centOADB->V0MTPCSigmaOutlierPar1(); 
1462   fV0MTPCSigmaOutlierPar2 = centOADB->V0MTPCSigmaOutlierPar2(); 
1463                             
1464   fV0MZDCOutlierPar0 =      centOADB->V0MZDCOutlierPar0();      
1465   fV0MZDCOutlierPar1 =      centOADB->V0MZDCOutlierPar1();      
1466   fV0MZDCEcalOutlierPar0 =  centOADB->V0MZDCEcalOutlierPar0();  
1467   fV0MZDCEcalOutlierPar1 =  centOADB->V0MZDCEcalOutlierPar1();  
1468
1469
1470
1471   return 0;
1472 }
1473
1474
1475
1476 //________________________________________________________________________
1477 Bool_t AliCentralitySelectionTask::IsOutlierV0MSPD(Float_t spd, Float_t v0, Int_t cent) const
1478 {
1479   // Clean outliers
1480   Float_t val = fV0MSPDOutlierPar0 +  fV0MSPDOutlierPar1 * v0;
1481   Float_t spdSigma = fV0MSPDSigmaOutlierPar0 + fV0MSPDSigmaOutlierPar1*cent + fV0MSPDSigmaOutlierPar2*cent*cent;
1482   if ( TMath::Abs(spd-val) > fOutliersCut*spdSigma ) 
1483     return kTRUE;
1484   else 
1485     return kFALSE;
1486 }
1487
1488 //________________________________________________________________________
1489 Bool_t AliCentralitySelectionTask::IsOutlierV0MTPC(Int_t tracks, Float_t v0, Int_t cent) const
1490 {
1491   // Clean outliers
1492   Float_t val = fV0MTPCOutlierPar0 +  fV0MTPCOutlierPar1 * v0;
1493   Float_t tpcSigma = fV0MTPCSigmaOutlierPar0 + fV0MTPCSigmaOutlierPar1*cent + fV0MTPCSigmaOutlierPar2*cent*cent;
1494   if ( TMath::Abs(tracks-val) > fOutliersCut*tpcSigma ) 
1495     return kTRUE;
1496   else 
1497     return kFALSE;
1498 }
1499
1500 //________________________________________________________________________
1501 Bool_t AliCentralitySelectionTask::IsOutlierV0MZDC(Float_t zdc, Float_t v0) const
1502 {
1503   // Clean outliers
1504   Float_t val = fV0MZDCOutlierPar0 + fV0MZDCOutlierPar1 * v0;
1505   if (zdc >  val) 
1506     return kTRUE;
1507   else 
1508   return kFALSE;
1509 }
1510
1511 //________________________________________________________________________
1512 Bool_t AliCentralitySelectionTask::IsOutlierV0MZDCECal(Float_t zdc, Float_t v0) const
1513 {
1514   // Clean outliers
1515   Float_t val = fV0MZDCEcalOutlierPar0 + fV0MZDCEcalOutlierPar1 * v0;
1516   if (zdc >  val) 
1517     return kTRUE;
1518   else 
1519     return kFALSE;
1520 }
1521
1522