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 "AliPhysicsSelectionTask.h"
  72 #include "AliPhysicsSelection.h"
  73 #include "AliBackgroundSelection.h"
  74 #include "AliESDUtils.h"
  75
  76 ClassImp(AliCentralitySelectionTask)
  77
  78
  79 //________________________________________________________________________
  80 AliCentralitySelectionTask::AliCentralitySelectionTask():
  81 AliAnalysisTaskSE(),
  82   fAnalysisInput("ESD"),
  83   fIsMCInput(kFALSE),
  84   fCurrentRun(-1),
  85   fUseScaling(0),
  86   fUseCleaning(0),
  87   fFillHistos(0),
  88   fV0MScaleFactor(0),
  89   fSPDScaleFactor(0),
  90   fTPCScaleFactor(0),
  91   fV0MScaleFactorMC(0),
  92   fV0MSPDOutlierPar0(0),
  93   fV0MSPDOutlierPar1(0),
  94   fV0MTPCOutlierPar0(0),
  95   fV0MTPCOutlierPar1(0),
  96   fV0MSPDSigmaOutlierPar0(0),
  97   fV0MSPDSigmaOutlierPar1(0),
  98   fV0MSPDSigmaOutlierPar2(0),
  99   fV0MTPCSigmaOutlierPar0(0),
 100   fV0MTPCSigmaOutlierPar1(0),
 101   fV0MTPCSigmaOutlierPar2(0),
 102   fV0MZDCOutlierPar0(0),
 103   fV0MZDCOutlierPar1(0),
 104   fV0MZDCEcalOutlierPar0(0),
 105   fV0MZDCEcalOutlierPar1(0),
 106   fTrackCuts(0),
 107   fZVCut(10),
 108   fOutliersCut(5),
 109   fQuality(999),
 110   fIsSelected(0),
 111   fMSL(0),
 112   fMSH(0),
 113   fMUL(0),
 114   fMLL(0),
 115   fEJE(0),
 116   fEGA(0),
 117   fPHS(0),
 118   fCVHN(0),
 119   fCVLN(0),
 120   fCVHNbit(0),
 121   fCVLNbit(0),
 122   fCCENT(0),
 123   fCSEMI(0),
 124   fCCENTbit(0),
 125   fCSEMIbit(0),
 126   fCentV0M(0),
 127   fCentFMD(0),
 128   fCentTRK(0),
 129   fCentTKL(0),
 130   fCentCL0(0),
 131   fCentCL1(0),
 132   fCentV0MvsFMD(0),
 133   fCentTKLvsV0M(0),
 134   fCentZEMvsZDC(0),
 135   fHtempV0M(0),
 136   fHtempFMD(0),
 137   fHtempTRK(0),
 138   fHtempTKL(0),
 139   fHtempCL0(0),
 140   fHtempCL1(0),
 141   fHtempV0MvsFMD(0),
 142   fHtempTKLvsV0M(0),
 143   fHtempZEMvsZDC(0),
 144   fOutputList(0),
 145   fHOutCentV0M     (0),
 146   fHOutCentV0MCVHN(0),
 147   fHOutCentV0MCVLN(0),
 148   fHOutCentV0MCVHNinMB(0),
 149   fHOutCentV0MCVLNinMB(0),
 150   fHOutCentV0MCCENT(0),
 151   fHOutCentV0MCSEMI(0),
 152   fHOutCentV0MCCENTinMB(0),
 153   fHOutCentV0MCSEMIinMB(0),
 154   fHOutCentV0MMSL(0),
 155   fHOutCentV0MMSH(0),
 156   fHOutCentV0MMUL(0),
 157   fHOutCentV0MMLL(0),
 158   fHOutCentV0MEJE(0),
 159   fHOutCentV0MEGA(0),
 160   fHOutCentV0MPHS(0),
 161   fHOutCentV0MMSLinMB(0),
 162   fHOutCentV0MMSHinMB(0),
 163   fHOutCentV0MMULinMB(0),
 164   fHOutCentV0MMLLinMB(0),
 165   fHOutCentV0MEJEinMB(0),
 166   fHOutCentV0MEGAinMB(0),
 167   fHOutCentV0MPHSinMB(0),
 168   fHOutCentFMD     (0),
 169   fHOutCentTRK     (0),
 170   fHOutCentTKL     (0),
 171   fHOutCentCL0     (0),
 172   fHOutCentCL1     (0),
 173   fHOutCentV0MvsFMD(0),
 174   fHOutCentTKLvsV0M(0),
 175   fHOutCentZEMvsZDC(0),
 176   fHOutCentV0MvsCentCL1(0),
 177   fHOutCentV0MvsCentTRK(0),
 178   fHOutCentTRKvsCentCL1(0),
 179   fHOutCentV0MvsCentZDC(0),
 180   fHOutMultV0M(0),
 181   fHOutMultV0O(0),
 182   fHOutMultFMD(0),
 183   fHOutMultTRK(0),
 184   fHOutMultTKL(0),
 185   fHOutMultCL0(0),
 186   fHOutMultCL1(0),
 187   fHOutMultV0MvsZDN(0),
 188   fHOutMultZEMvsZDN(0),
 189   fHOutMultV0MvsZDC(0),
 190   fHOutMultZEMvsZDC(0),
 191   fHOutMultZEMvsZDCw(0),
 192   fHOutMultV0MvsCL1(0),
 193   fHOutMultV0MvsTRK(0),
 194   fHOutMultTRKvsCL1(0),
 195   fHOutMultV0MvsV0O(0),
 196   fHOutMultV0OvsCL1(0),
 197   fHOutMultV0OvsTRK(0),
 198   fHOutCentV0Mqual1(0),
 199   fHOutCentTRKqual1(0),
 200   fHOutCentCL1qual1(0),
 201   fHOutMultV0MvsCL1qual1(0),
 202   fHOutMultV0MvsTRKqual1(0),
 203   fHOutMultTRKvsCL1qual1(0),
 204   fHOutCentV0Mqual2(0),
 205   fHOutCentTRKqual2(0),
 206   fHOutCentCL1qual2(0),
 207   fHOutMultV0MvsCL1qual2(0),
 208   fHOutMultV0MvsTRKqual2(0),
 209   fHOutMultTRKvsCL1qual2(0),
 210   fHOutQuality(0),
 211   fHOutVertex(0),
 212   fHOutVertexT0(0)
 213 {
 214   // Default constructor
 215   AliInfo("Centrality Selection enabled.");
 216
 217   fUseScaling=kTRUE;
 218   fUseCleaning=kTRUE;
 219   fFillHistos=kFALSE;
 220   fBranchNames="ESD:AliESDRun.,AliESDHeader.,AliESDZDC.,AliESDFMD.,AliESDVZERO.,AliESDTZERO."
 221     ",SPDVertex.,TPCVertex.,PrimaryVertex.,AliMultiplicity.,Tracks ";
 222 }
 223
 224 //________________________________________________________________________
 225 AliCentralitySelectionTask::AliCentralitySelectionTask(const char *name):
 226   AliAnalysisTaskSE(name),
 227   fAnalysisInput("ESD"),
 228   fIsMCInput(kFALSE),
 229   fCurrentRun(-1),
 230   fUseScaling(0),
 231   fUseCleaning(0),
 232   fFillHistos(0),
 233   fV0MScaleFactor(0),
 234   fSPDScaleFactor(0),
 235   fTPCScaleFactor(0),
 236   fV0MScaleFactorMC(0),
 237   fV0MSPDOutlierPar0(0),
 238   fV0MSPDOutlierPar1(0),
 239   fV0MTPCOutlierPar0(0),
 240   fV0MTPCOutlierPar1(0),
 241   fV0MSPDSigmaOutlierPar0(0),
 242   fV0MSPDSigmaOutlierPar1(0),
 243   fV0MSPDSigmaOutlierPar2(0),
 244   fV0MTPCSigmaOutlierPar0(0),
 245   fV0MTPCSigmaOutlierPar1(0),
 246   fV0MTPCSigmaOutlierPar2(0),
 247   fV0MZDCOutlierPar0(0),
 248   fV0MZDCOutlierPar1(0),
 249   fV0MZDCEcalOutlierPar0(0),
 250   fV0MZDCEcalOutlierPar1(0),
 251   fTrackCuts(0),
 252   fZVCut(10),
 253   fOutliersCut(5),
 254   fQuality(999),
 255   fIsSelected(0),
 256   fMSL(0),
 257   fMSH(0),
 258   fMUL(0),
 259   fMLL(0),
 260   fEJE(0),
 261   fEGA(0),
 262   fPHS(0),
 263   fCVHN(0),
 264   fCVLN(0),
 265   fCVHNbit(0),
 266   fCVLNbit(0),
 267   fCCENT(0),
 268   fCSEMI(0),
 269   fCCENTbit(0),
 270   fCSEMIbit(0),
 271   fCentV0M(0),
 272   fCentFMD(0),
 273   fCentTRK(0),
 274   fCentTKL(0),
 275   fCentCL0(0),
 276   fCentCL1(0),
 277   fCentV0MvsFMD(0),
 278   fCentTKLvsV0M(0),
 279   fCentZEMvsZDC(0),
 280   fHtempV0M(0),
 281   fHtempFMD(0),
 282   fHtempTRK(0),
 283   fHtempTKL(0),
 284   fHtempCL0(0),
 285   fHtempCL1(0),
 286   fHtempV0MvsFMD(0),
 287   fHtempTKLvsV0M(0),
 288   fHtempZEMvsZDC(0),
 289   fOutputList(0),
 290   fHOutCentV0M     (0),
 291   fHOutCentV0MCVHN(0),
 292   fHOutCentV0MCVLN(0),
 293   fHOutCentV0MCVHNinMB(0),
 294   fHOutCentV0MCVLNinMB(0),
 295   fHOutCentV0MCCENT(0),
 296   fHOutCentV0MCSEMI(0),
 297   fHOutCentV0MCCENTinMB(0),
 298   fHOutCentV0MCSEMIinMB(0),
 299   fHOutCentV0MMSL(0),
 300   fHOutCentV0MMSH(0),
 301   fHOutCentV0MMUL(0),
 302   fHOutCentV0MMLL(0),
 303   fHOutCentV0MEJE(0),
 304   fHOutCentV0MEGA(0),
 305   fHOutCentV0MPHS(0),
 306   fHOutCentV0MMSLinMB(0),
 307   fHOutCentV0MMSHinMB(0),
 308   fHOutCentV0MMULinMB(0),
 309   fHOutCentV0MMLLinMB(0),
 310   fHOutCentV0MEJEinMB(0),
 311   fHOutCentV0MEGAinMB(0),
 312   fHOutCentV0MPHSinMB(0),
 313   fHOutCentFMD     (0),
 314   fHOutCentTRK     (0),
 315   fHOutCentTKL     (0),
 316   fHOutCentCL0     (0),
 317   fHOutCentCL1     (0),
 318   fHOutCentV0MvsFMD(0),
 319   fHOutCentTKLvsV0M(0),
 320   fHOutCentZEMvsZDC(0),
 321   fHOutCentV0MvsCentCL1(0),
 322   fHOutCentV0MvsCentTRK(0),
 323   fHOutCentTRKvsCentCL1(0),
 324   fHOutCentV0MvsCentZDC(0),
 325   fHOutMultV0M(0),
 326   fHOutMultV0O(0),
 327   fHOutMultFMD(0),
 328   fHOutMultTRK(0),
 329   fHOutMultTKL(0),
 330   fHOutMultCL0(0),
 331   fHOutMultCL1(0),
 332   fHOutMultV0MvsZDN(0),
 333   fHOutMultZEMvsZDN(0),
 334   fHOutMultV0MvsZDC(0),
 335   fHOutMultZEMvsZDC(0),
 336   fHOutMultZEMvsZDCw(0),
 337   fHOutMultV0MvsCL1(0),
 338   fHOutMultV0MvsTRK(0),
 339   fHOutMultTRKvsCL1(0),
 340   fHOutMultV0MvsV0O(0),
 341   fHOutMultV0OvsCL1(0),
 342   fHOutMultV0OvsTRK(0),
 343   fHOutCentV0Mqual1(0),
 344   fHOutCentTRKqual1(0),
 345   fHOutCentCL1qual1(0),
 346   fHOutMultV0MvsCL1qual1(0),
 347   fHOutMultV0MvsTRKqual1(0),
 348   fHOutMultTRKvsCL1qual1(0),
 349   fHOutCentV0Mqual2(0),
 350   fHOutCentTRKqual2(0),
 351   fHOutCentCL1qual2(0),
 352   fHOutMultV0MvsCL1qual2(0),
 353   fHOutMultV0MvsTRKqual2(0),
 354   fHOutMultTRKvsCL1qual2(0),
 355   fHOutQuality(0),
 356   fHOutVertex(0),
 357   fHOutVertexT0(0)
 358 {
 359   // Default constructor
 360   AliInfo("Centrality Selection enabled.");
 361   //DefineOutput(1, TList::Class());
 362   fUseScaling=kTRUE;
 363   fUseCleaning=kTRUE;
 364   fFillHistos=kFALSE;
 365   fBranchNames="ESD:AliESDRun.,AliESDHeader.,AliESDZDC.,AliESDFMD.,AliESDVZERO.,AliESDTZERO."
 366     ",SPDVertex.,TPCVertex.,PrimaryVertex.,AliMultiplicity.,Tracks ";
 367 }
 368
 369 //________________________________________________________________________
 370 AliCentralitySelectionTask& AliCentralitySelectionTask::operator=(const AliCentralitySelectionTask& c)
 371 {
 372   // Assignment operator
 373   if (this!=&c) {
 374     AliAnalysisTaskSE::operator=(c);
 375   }
 376   return *this;
 377 }
 378
 379 //________________________________________________________________________
 380 AliCentralitySelectionTask::AliCentralitySelectionTask(const AliCentralitySelectionTask& ana):
 381   AliAnalysisTaskSE(ana),
 382   fAnalysisInput(ana.fAnalysisInput),
 383   fIsMCInput(ana.fIsMCInput),
 384   fCurrentRun(ana.fCurrentRun),
 385   fUseScaling(ana.fUseScaling),
 386   fUseCleaning(ana.fUseCleaning),
 387   fFillHistos(ana.fFillHistos),
 388   fV0MScaleFactor(ana.fV0MScaleFactor),
 389   fSPDScaleFactor(ana.fSPDScaleFactor),
 390   fTPCScaleFactor(ana.fTPCScaleFactor),
 391   fV0MScaleFactorMC(ana.fV0MScaleFactorMC),
 392   fV0MSPDOutlierPar0(ana.fV0MSPDOutlierPar0),
 393   fV0MSPDOutlierPar1(ana.fV0MSPDOutlierPar1),
 394   fV0MTPCOutlierPar0(ana.fV0MTPCOutlierPar0),
 395   fV0MTPCOutlierPar1(ana.fV0MTPCOutlierPar1),
 396   fV0MSPDSigmaOutlierPar0(ana.fV0MSPDSigmaOutlierPar0),
 397   fV0MSPDSigmaOutlierPar1(ana.fV0MSPDSigmaOutlierPar1),
 398   fV0MSPDSigmaOutlierPar2(ana.fV0MSPDSigmaOutlierPar2),
 399   fV0MTPCSigmaOutlierPar0(ana.fV0MTPCSigmaOutlierPar0),
 400   fV0MTPCSigmaOutlierPar1(ana.fV0MTPCSigmaOutlierPar1),
 401   fV0MTPCSigmaOutlierPar2(ana.fV0MTPCSigmaOutlierPar2),
 402   fV0MZDCOutlierPar0(ana.fV0MZDCOutlierPar0),
 403   fV0MZDCOutlierPar1(ana.fV0MZDCOutlierPar1),
 404   fV0MZDCEcalOutlierPar0(ana.fV0MZDCEcalOutlierPar0),
 405   fV0MZDCEcalOutlierPar1(ana.fV0MZDCEcalOutlierPar1),
 406   fTrackCuts(ana.fTrackCuts),
 407   fZVCut(ana.fZVCut),
 408   fOutliersCut(ana.fOutliersCut),
 409   fQuality(ana.fQuality),
 410   fIsSelected(ana.fIsSelected),
 411   fMSL(ana.fMSL),
 412   fMSH(ana.fMSH),
 413   fMUL(ana.fMUL),
 414   fMLL(ana.fMLL),
 415   fEJE(ana.fEJE),
 416   fEGA(ana.fEGA),
 417   fPHS(ana.fPHS),
 418   fCVHN(ana.fCVHN),
 419   fCVLN(ana.fCVLN),
 420   fCVHNbit(ana.fCVHNbit),
 421   fCVLNbit(ana.fCVLNbit),
 422   fCCENT(ana.fCCENT),
 423   fCSEMI(ana.fCSEMI),
 424   fCCENTbit(ana.fCCENTbit),
 425   fCSEMIbit(ana.fCSEMIbit),
 426   fCentV0M(ana.fCentV0M),
 427   fCentFMD(ana.fCentFMD),
 428   fCentTRK(ana.fCentTRK),
 429   fCentTKL(ana.fCentTKL),
 430   fCentCL0(ana.fCentCL0),
 431   fCentCL1(ana.fCentCL1),
 432   fCentV0MvsFMD(ana.fCentV0MvsFMD),
 433   fCentTKLvsV0M(ana.fCentTKLvsV0M),
 434   fCentZEMvsZDC(ana.fCentZEMvsZDC),
 435   fHtempV0M(ana.fHtempV0M),
 436   fHtempFMD(ana.fHtempFMD),
 437   fHtempTRK(ana.fHtempTRK),
 438   fHtempTKL(ana.fHtempTKL),
 439   fHtempCL0(ana.fHtempCL0),
 440   fHtempCL1(ana.fHtempCL1),
 441   fHtempV0MvsFMD(ana.fHtempV0MvsFMD),
 442   fHtempTKLvsV0M(ana.fHtempTKLvsV0M),
 443   fHtempZEMvsZDC(ana.fHtempZEMvsZDC),
 444   fOutputList(ana.fOutputList),
 445   fHOutCentV0M     (ana.fHOutCentV0M     ),
 446   fHOutCentV0MCVHN(ana.fHOutCentV0MCVHN),
 447   fHOutCentV0MCVLN(ana.fHOutCentV0MCVLN),
 448   fHOutCentV0MCVHNinMB(ana.fHOutCentV0MCVHNinMB),
 449   fHOutCentV0MCVLNinMB(ana.fHOutCentV0MCVLNinMB),
 450   fHOutCentV0MCCENT(ana.fHOutCentV0MCCENT),
 451   fHOutCentV0MCSEMI(ana.fHOutCentV0MCSEMI),
 452   fHOutCentV0MCCENTinMB(ana.fHOutCentV0MCCENTinMB),
 453   fHOutCentV0MCSEMIinMB(ana.fHOutCentV0MCSEMIinMB),
 454   fHOutCentV0MMSL(ana.fHOutCentV0MMSL),
 455   fHOutCentV0MMSH(ana.fHOutCentV0MMSH),
 456   fHOutCentV0MMUL(ana.fHOutCentV0MMUL),
 457   fHOutCentV0MMLL(ana.fHOutCentV0MMLL),
 458   fHOutCentV0MEJE(ana.fHOutCentV0MEJE),
 459   fHOutCentV0MEGA(ana.fHOutCentV0MEGA),
 460   fHOutCentV0MPHS(ana.fHOutCentV0MPHS),
 461   fHOutCentV0MMSLinMB(ana.fHOutCentV0MMSLinMB),
 462   fHOutCentV0MMSHinMB(ana.fHOutCentV0MMSHinMB),
 463   fHOutCentV0MMULinMB(ana.fHOutCentV0MMULinMB),
 464   fHOutCentV0MMLLinMB(ana.fHOutCentV0MMLLinMB),
 465   fHOutCentV0MEJEinMB(ana.fHOutCentV0MEJEinMB),
 466   fHOutCentV0MEGAinMB(ana.fHOutCentV0MEGAinMB),
 467   fHOutCentV0MPHSinMB(ana.fHOutCentV0MPHSinMB),
 468   fHOutCentFMD     (ana.fHOutCentFMD     ),
 469   fHOutCentTRK     (ana.fHOutCentTRK     ),
 470   fHOutCentTKL     (ana.fHOutCentTKL     ),
 471   fHOutCentCL0     (ana.fHOutCentCL0     ),
 472   fHOutCentCL1     (ana.fHOutCentCL1     ),
 473   fHOutCentV0MvsFMD(ana.fHOutCentV0MvsFMD),
 474   fHOutCentTKLvsV0M(ana.fHOutCentTKLvsV0M),
 475   fHOutCentZEMvsZDC(ana.fHOutCentZEMvsZDC),
 476   fHOutCentV0MvsCentCL1(ana.fHOutCentV0MvsCentCL1),
 477   fHOutCentV0MvsCentTRK(ana.fHOutCentV0MvsCentTRK),
 478   fHOutCentTRKvsCentCL1(ana.fHOutCentTRKvsCentCL1),
 479   fHOutCentV0MvsCentZDC(ana.fHOutCentV0MvsCentZDC),
 480   fHOutMultV0M(ana.fHOutMultV0M),
 481   fHOutMultV0O(ana.fHOutMultV0O),
 482   fHOutMultFMD(ana.fHOutMultFMD),
 483   fHOutMultTRK(ana.fHOutMultTRK),
 484   fHOutMultTKL(ana.fHOutMultTKL),
 485   fHOutMultCL0(ana.fHOutMultCL0),
 486   fHOutMultCL1(ana.fHOutMultCL1),
 487   fHOutMultV0MvsZDN(ana.fHOutMultV0MvsZDN),
 488   fHOutMultZEMvsZDN(ana.fHOutMultZEMvsZDN),
 489   fHOutMultV0MvsZDC(ana.fHOutMultV0MvsZDC),
 490   fHOutMultZEMvsZDC(ana.fHOutMultZEMvsZDC),
 491   fHOutMultZEMvsZDCw(ana.fHOutMultZEMvsZDCw),
 492   fHOutMultV0MvsCL1(ana.fHOutMultV0MvsCL1),
 493   fHOutMultV0MvsTRK(ana.fHOutMultV0MvsTRK),
 494   fHOutMultTRKvsCL1(ana.fHOutMultTRKvsCL1),
 495   fHOutMultV0MvsV0O(ana.fHOutMultV0MvsV0O),
 496   fHOutMultV0OvsCL1(ana.fHOutMultV0OvsCL1),
 497   fHOutMultV0OvsTRK(ana.fHOutMultV0OvsTRK),
 498   fHOutCentV0Mqual1(ana.fHOutCentV0Mqual1),
 499   fHOutCentTRKqual1(ana.fHOutCentTRKqual1),
 500   fHOutCentCL1qual1(ana.fHOutCentCL1qual1),
 501   fHOutMultV0MvsCL1qual1(ana.fHOutMultV0MvsCL1qual1),
 502   fHOutMultV0MvsTRKqual1(ana.fHOutMultV0MvsTRKqual1),
 503   fHOutMultTRKvsCL1qual1(ana.fHOutMultTRKvsCL1qual1),
 504   fHOutCentV0Mqual2(ana.fHOutCentV0Mqual2),
 505   fHOutCentTRKqual2(ana.fHOutCentTRKqual2),
 506   fHOutCentCL1qual2(ana.fHOutCentCL1qual2),
 507   fHOutMultV0MvsCL1qual2(ana.fHOutMultV0MvsCL1qual2),
 508   fHOutMultV0MvsTRKqual2(ana.fHOutMultV0MvsTRKqual2),
 509   fHOutMultTRKvsCL1qual2(ana.fHOutMultTRKvsCL1qual2),
 510   fHOutQuality(ana.fHOutQuality),
 511   fHOutVertex(ana.fHOutVertex),
 512   fHOutVertexT0(ana.fHOutVertexT0)
 513 {
 514   // Copy Constructor
 515
 516 }
 517
 518 //________________________________________________________________________
 519 AliCentralitySelectionTask::~AliCentralitySelectionTask()
 520 {
 521   // Destructor
 522   if (fOutputList && !AliAnalysisManager::GetAnalysisManager()->IsProofMode()) delete fOutputList;
 523   if (fTrackCuts) delete fTrackCuts;
 524 }
 525
 526 //________________________________________________________________________
 527 void AliCentralitySelectionTask::UserCreateOutputObjects()
 528 {
 529   // Create the output containers
 530   if(fDebug>1) printf("AnalysisCentralitySelectionTask::UserCreateOutputObjects() \n");
 531   AliLog::SetClassDebugLevel("AliCentralitySelectionTask", AliLog::kInfo);
 532
 533   if (fFillHistos) {
 534     fOutputList = new TList();
 535     fOutputList->SetOwner();
 536     fHOutCentV0M     = new TH1F("fHOutCentV0M","fHOutCentV0M; Centrality V0",505,0,101);
 537     fHOutCentV0MCVHN= new TH1F("fHOutCentV0M_CVHN","fHOutCentV0M_CVHN; Centrality V0",505,0,101);
 538     fHOutCentV0MCVLN= new TH1F("fHOutCentV0M_CVLN","fHOutCentV0M_CVLN; Centrality V0",505,0,101);
 539     fHOutCentV0MCVHNinMB= new TH1F("fHOutCentV0M_CVHNinMB","fHOutCentV0M_CVHN; Centrality V0",505,0,101);
 540     fHOutCentV0MCVLNinMB= new TH1F("fHOutCentV0M_CVLNinMB","fHOutCentV0M_CVLN; Centrality V0",505,0,101);
 541     fHOutCentV0MCCENT= new TH1F("fHOutCentV0M_CCENT","fHOutCentV0M_CCENT; Centrality V0",505,0,101);
 542     fHOutCentV0MCSEMI= new TH1F("fHOutCentV0M_CSEMI","fHOutCentV0M_CSEMI; Centrality V0",505,0,101);
 543     fHOutCentV0MCCENTinMB= new TH1F("fHOutCentV0M_CCENTinMB","fHOutCentV0M_CCENT; Centrality V0",505,0,101);
 544     fHOutCentV0MCSEMIinMB= new TH1F("fHOutCentV0M_CSEMIinMB","fHOutCentV0M_CSEMI; Centrality V0",505,0,101);
 545     fHOutCentV0MMSL= new TH1F("fHOutCentV0M_MSL","fHOutCentV0M_MSL; Centrality V0",505,0,101);
 546     fHOutCentV0MMSH= new TH1F("fHOutCentV0M_MSH","fHOutCentV0M_MSH; Centrality V0",505,0,101);
 547     fHOutCentV0MMUL= new TH1F("fHOutCentV0M_MUL","fHOutCentV0M_MUL; Centrality V0",505,0,101);
 548     fHOutCentV0MMLL= new TH1F("fHOutCentV0M_MLL","fHOutCentV0M_MLL; Centrality V0",505,0,101);
 549     fHOutCentV0MEJE= new TH1F("fHOutCentV0M_EJE","fHOutCentV0M_EJE; Centrality V0",505,0,101);
 550     fHOutCentV0MEGA= new TH1F("fHOutCentV0M_EGA","fHOutCentV0M_EGA; Centrality V0",505,0,101);
 551     fHOutCentV0MPHS= new TH1F("fHOutCentV0M_PHS","fHOutCentV0M_PHS; Centrality V0",505,0,101);
 552     fHOutCentV0MMSLinMB= new TH1F("fHOutCentV0M_MSLinMB","fHOutCentV0M_MSLinMB; Centrality V0",505,0,101);
 553     fHOutCentV0MMSHinMB= new TH1F("fHOutCentV0M_MSHinMB","fHOutCentV0M_MSHinMB; Centrality V0",505,0,101);
 554     fHOutCentV0MMULinMB= new TH1F("fHOutCentV0M_MULinMB","fHOutCentV0M_MULinMB; Centrality V0",505,0,101);
 555     fHOutCentV0MMLLinMB= new TH1F("fHOutCentV0M_MLLinMB","fHOutCentV0M_MLLinMB; Centrality V0",505,0,101);
 556     fHOutCentV0MEJEinMB= new TH1F("fHOutCentV0M_EJEinMB","fHOutCentV0M_EJEinMB; Centrality V0",505,0,101);
 557     fHOutCentV0MEGAinMB= new TH1F("fHOutCentV0M_EGAinMB","fHOutCentV0M_EGAinMB; Centrality V0",505,0,101);
 558     fHOutCentV0MPHSinMB= new TH1F("fHOutCentV0M_PHSinMB","fHOutCentV0M_PHSinMB; Centrality V0",505,0,101);
 559     fHOutCentFMD     = new TH1F("fHOutCentFMD","fHOutCentFMD; Centrality FMD",505,0,101);
 560     fHOutCentTRK     = new TH1F("fHOutCentTRK","fHOutCentTRK; Centrality TPC",505,0,101);
 561     fHOutCentTKL     = new TH1F("fHOutCentTKL","fHOutCentTKL; Centrality tracklets",505,0,101);
 562     fHOutCentCL0     = new TH1F("fHOutCentCL0","fHOutCentCL0; Centrality SPD inner",505,0,101);
 563     fHOutCentCL1     = new TH1F("fHOutCentCL1","fHOutCentCL1; Centrality SPD outer",505,0,101);
 564     fHOutCentV0MvsFMD= new TH1F("fHOutCentV0MvsFMD","fHOutCentV0MvsFMD; Centrality V0 vs FMD",505,0,101);
 565     fHOutCentTKLvsV0M= new TH1F("fHOutCentTKLvsV0M","fHOutCentTKLvsV0M; Centrality tracklets vs V0",505,0,101);
 566     fHOutCentZEMvsZDC= new TH1F("fHOutCentZEMvsZDC","fHOutCentZEMvsZDC; Centrality ZEM vs ZDC",505,0,101);
 567     fHOutCentV0MvsCentCL1= new TH2F("fHOutCentV0MvsCentCL1","fHOutCentV0MvsCentCL1; Cent V0; Cent SPD",505,0,101,505,0,101);
 568     fHOutCentV0MvsCentTRK= new TH2F("fHOutCentV0MvsCentTRK","fHOutCentV0MvsCentTRK; Cent V0; Cent TPC",505,0,101,505,0,101);
 569     fHOutCentTRKvsCentCL1= new TH2F("fHOutCentTRKvsCentCL1","fHOutCentTRKvsCentCL1; Cent TPC; Cent SPD",505,0,101,505,0,101);
 570     fHOutCentV0MvsCentZDC= new TH2F("fHOutCentV0MvsCentZDC","fHOutCentV0MvsCentZDC; Cent V0; Cent ZDC",505,0,101,505,0,101);
 571     fHOutMultV0M = new TH1F("fHOutMultV0M","fHOutMultV0M; Multiplicity V0",25000,0,30000);
 572     fHOutMultV0O = new TH1F("fHOutMultV0O","fHOutMultV0O; Multiplicity V0",40000,0,40000);
 573     fHOutMultFMD = new TH1F("fHOutMultFMD","fHOutMultFMD; Multiplicity FMD",24000,0,24000);
 574     fHOutMultTRK = new TH1F("fHOutMultTRK","fHOutMultTRK; Multiplicity TPC",4000,0,4000);
 575     fHOutMultTKL = new TH1F("fHOutMultTKL","fHOutMultTKL; Multiplicity tracklets",5000,0,5000);
 576     fHOutMultCL0 = new TH1F("fHOutMultCL0","fHOutMultCL0; Multiplicity SPD inner",7000,0,7000);
 577     fHOutMultCL1 = new TH1F("fHOutMultCL1","fHOutMultCL1; Multiplicity SPD outer",7000,0,7000);
 578     fHOutMultV0MvsZDN = new TH2F("fHOutMultV0MvsZDN","fHOutMultV0MvsZDN; Multiplicity V0; Energy ZDC-N",500,0,30000,500,0,180000);
 579     fHOutMultZEMvsZDN = new TH2F("fHOutMultZEMvsZDN","fHOutMultZEMvsZDN; Energy ZEM; Energy ZDC-N",500,0,2500,500,0,180000);
 580     fHOutMultV0MvsZDC = new TH2F("fHOutMultV0MvsZDC","fHOutMultV0MvsZDC; Multiplicity V0; Energy ZDC",500,0,30000,500,0,200000);
 581     fHOutMultZEMvsZDC = new TH2F("fHOutMultZEMvsZDC","fHOutMultZEMvsZDC; Energy ZEM; Energy ZDC",500,0,2500,500,0,200000);
 582     fHOutMultZEMvsZDCw = new TH2F("fHOutMultZEMvsZDCw","fHOutMultZEMvsZDCw; Energy ZEM; Energy ZDC (weigthed with V0 percentile)",500,0,2500,500,0,200000);
 583     fHOutMultV0MvsCL1 = new TH2F("fHOutMultV0MvsCL1","fHOutMultV0MvsCL1; Multiplicity V0; Multiplicity SPD outer",2500,0,30000,700,0,7000);
 584     fHOutMultV0MvsTRK = new TH2F("fHOutMultV0MvsTRK","fHOutMultV0MvsTRK; Multiplicity V0; Multiplicity TPC",2500,0,30000,400,0,4000);
 585     fHOutMultTRKvsCL1 = new TH2F("fHOutMultTRKvsCL1","fHOutMultTRKvsCL1; Multiplicity TPC; Multiplicity SPD outer",400,0,4000,700,0,7000);
 586     fHOutMultV0MvsV0O = new TH2F("fHOutMultV0MvsV0O","fHOutMultV0MvsV0O; Multiplicity V0; Multiplicity V0 Online",500,0,30000,500,0,40000);
 587     fHOutMultV0OvsCL1 = new TH2F("fHOutMultV0OvsCL1","fHOutMultV0OvsCL1; Multiplicity V0; Multiplicity SPD outer",500,0,40000,700,0,7000);
 588     fHOutMultV0OvsTRK = new TH2F("fHOutMultV0OvsTRK","fHOutMultV0OvsTRK; Multiplicity V0; Multiplicity TPC",500,0,40000,400,0,4000);
 589     fHOutMultV0MvsV0O = new TH2F("fHOutMultV0MvsV0O","fHOutMultV0MvsV0O; Multiplicity V0; Multiplicity V0 Online",500,0,30000,500,0,30000);
 590     fHOutMultV0OvsCL1 = new TH2F("fHOutMultV0OvsCL1","fHOutMultV0OvsCL1; Multiplicity V0; Multiplicity SPD outer",2500,0,30000,700,0,7000);
 591     fHOutMultV0OvsTRK = new TH2F("fHOutMultV0OvsTRK","fHOutMultV0OvsTRK; Multiplicity V0; Multiplicity TPC",2500,0,30000,400,0,4000);
 592
 593     fHOutCentV0Mqual1 = new TH1F("fHOutCentV0M_qual1","fHOutCentV0M_qual1; Centrality V0",505,0,101);
 594     fHOutCentTRKqual1 = new TH1F("fHOutCentTRK_qual1","fHOutCentTRK_qual1; Centrality TPC",505,0,101);
 595     fHOutCentCL1qual1 = new TH1F("fHOutCentCL1_qual1","fHOutCentCL1_qual1; Centrality SPD outer",505,0,101);
 596     fHOutMultV0MvsCL1qual1 = new TH2F("fHOutMultV0MvsCL1_qual1","fHOutMultV0MvsCL1_qual1; Multiplicity V0; Multiplicity SPD outer",2500,0,25000,700,0,7000);
 597     fHOutMultV0MvsTRKqual1 = new TH2F("fHOutMultV0MvsTRK_qual1","fHOutMultV0MvsTRK_qual1; Multiplicity V0; Multiplicity TPC",2500,0,25000,400,0,4000);
 598     fHOutMultTRKvsCL1qual1 = new TH2F("fHOutMultTRKvsCL1_qual1","fHOutMultTRKvsCL1_qual1; Multiplicity TPC; Multiplicity SPD outer",400,0,4000,700,0,7000);
 599
 600     fHOutCentV0Mqual2 = new TH1F("fHOutCentV0M_qual2","fHOutCentV0M_qual2; Centrality V0",505,0,101);
 601     fHOutCentTRKqual2 = new TH1F("fHOutCentTRK_qual2","fHOutCentTRK_qual2; Centrality TPC",505,0,101);
 602     fHOutCentCL1qual2 = new TH1F("fHOutCentCL1_qual2","fHOutCentCL1_qual2; Centrality SPD outer",505,0,101);
 603     fHOutMultV0MvsCL1qual2 = new TH2F("fHOutMultV0MvsCL1_qual2","fHOutMultV0MvsCL1_qual2; Multiplicity V0; Multiplicity SPD outer",2500,0,25000,700,0,7000);
 604     fHOutMultV0MvsTRKqual2 = new TH2F("fHOutMultV0MvsTRK_qual2","fHOutMultV0MvsTRK_qual2; Multiplicity V0; Multiplicity TPC",2500,0,25000,400,0,4000);
 605     fHOutMultTRKvsCL1qual2 = new TH2F("fHOutMultTRKvsCL1_qual2","fHOutMultTRKvsCL1_qual2; Multiplicity TPC; Multiplicity SPD outer",400,0,4000,700,0,7000);
 606
 607     fHOutQuality = new TH1F("fHOutQuality", "fHOutQuality", 100,-0.5,99.5);
 608     fHOutVertex  = new TH1F("fHOutVertex", "fHOutVertex", 100,-20,20);
 609     fHOutVertexT0  = new TH1F("fHOutVertexT0", "fHOutVertexT0", 100,-20,20);
 610
 611     fOutputList->Add(  fHOutCentV0M     );
 612     fOutputList->Add(  fHOutCentV0MCVHN);
 613     fOutputList->Add(  fHOutCentV0MCVLN);
 614     fOutputList->Add(  fHOutCentV0MCVHNinMB);
 615     fOutputList->Add(  fHOutCentV0MCVLNinMB);
 616     fOutputList->Add(  fHOutCentV0MCCENT);
 617     fOutputList->Add(  fHOutCentV0MCSEMI);
 618     fOutputList->Add(  fHOutCentV0MCCENTinMB);
 619     fOutputList->Add(  fHOutCentV0MCSEMIinMB);
 620     fOutputList->Add(  fHOutCentV0MMSL    );
 621     fOutputList->Add(  fHOutCentV0MMSH    );
 622     fOutputList->Add(  fHOutCentV0MMUL    );
 623     fOutputList->Add(  fHOutCentV0MMLL    );
 624     fOutputList->Add(  fHOutCentV0MEJE    );
 625     fOutputList->Add(  fHOutCentV0MEGA    );
 626     fOutputList->Add(  fHOutCentV0MPHS   );
 627     fOutputList->Add(  fHOutCentV0MMSLinMB);
 628     fOutputList->Add(  fHOutCentV0MMSHinMB);
 629     fOutputList->Add(  fHOutCentV0MMULinMB);
 630     fOutputList->Add(  fHOutCentV0MMLLinMB);
 631     fOutputList->Add(  fHOutCentV0MEJEinMB);
 632     fOutputList->Add(  fHOutCentV0MEGAinMB);
 633     fOutputList->Add(  fHOutCentV0MPHSinMB);
 634     fOutputList->Add(  fHOutCentFMD     );
 635     fOutputList->Add(  fHOutCentTRK     );
 636     fOutputList->Add(  fHOutCentTKL     );
 637     fOutputList->Add(  fHOutCentCL0     );
 638     fOutputList->Add(  fHOutCentCL1     );
 639     fOutputList->Add(  fHOutCentV0MvsFMD);
 640     fOutputList->Add(  fHOutCentTKLvsV0M);
 641     fOutputList->Add(  fHOutCentZEMvsZDC);
 642     fOutputList->Add(  fHOutCentV0MvsCentCL1);
 643     fOutputList->Add(  fHOutCentV0MvsCentTRK);
 644     fOutputList->Add(  fHOutCentTRKvsCentCL1);
 645     fOutputList->Add(  fHOutCentV0MvsCentZDC);
 646     fOutputList->Add(  fHOutMultV0M);
 647     fOutputList->Add(  fHOutMultV0O);
 648     fOutputList->Add(  fHOutMultFMD);
 649     fOutputList->Add(  fHOutMultTRK);
 650     fOutputList->Add(  fHOutMultTKL);
 651     fOutputList->Add(  fHOutMultCL0);
 652     fOutputList->Add(  fHOutMultCL1);
 653     fOutputList->Add(  fHOutMultV0MvsZDN);
 654     fOutputList->Add(  fHOutMultZEMvsZDN);
 655     fOutputList->Add(  fHOutMultV0MvsZDC);
 656     fOutputList->Add(  fHOutMultZEMvsZDC);
 657     fOutputList->Add(  fHOutMultZEMvsZDCw);
 658     fOutputList->Add(  fHOutMultV0MvsCL1);
 659     fOutputList->Add(  fHOutMultV0MvsTRK);
 660     fOutputList->Add(  fHOutMultTRKvsCL1);
 661     fOutputList->Add(  fHOutMultV0MvsV0O);
 662     fOutputList->Add(  fHOutMultV0OvsCL1);
 663     fOutputList->Add(  fHOutMultV0OvsTRK);
 664     fOutputList->Add(  fHOutCentV0Mqual1 );
 665     fOutputList->Add(  fHOutCentTRKqual1 );
 666     fOutputList->Add(  fHOutCentCL1qual1 );
 667     fOutputList->Add(  fHOutMultV0MvsCL1qual1);
 668     fOutputList->Add(  fHOutMultV0MvsTRKqual1);
 669     fOutputList->Add(  fHOutMultTRKvsCL1qual1);
 670     fOutputList->Add(  fHOutCentV0Mqual2 );
 671     fOutputList->Add(  fHOutCentTRKqual2 );
 672     fOutputList->Add(  fHOutCentCL1qual2 );
 673     fOutputList->Add(  fHOutMultV0MvsCL1qual2);
 674     fOutputList->Add(  fHOutMultV0MvsTRKqual2);
 675     fOutputList->Add(  fHOutMultTRKvsCL1qual2);
 676     fOutputList->Add(  fHOutQuality );
 677     fOutputList->Add(  fHOutVertex );
 678     fOutputList->Add(  fHOutVertexT0 );
 679
 680     PostData(1, fOutputList);
 681   }
 682
 683   fTrackCuts = AliESDtrackCuts::GetStandardTPCOnlyTrackCuts();
 684 }
 685
 686 //________________________________________________________________________
 687 void AliCentralitySelectionTask::UserExec(Option_t */*option*/)
 688 {
 689   // Execute analysis for current event:
 690   if(fDebug>1) printf(" **** AliCentralitySelectionTask::UserExec() \n");
 691
 692   Float_t  zncEnergy = 0.;          //  ZNC Energy
 693   Float_t  zpcEnergy = 0.;          //  ZPC Energy
 694   Float_t  znaEnergy = 0.;          //  ZNA Energy
 695   Float_t  zpaEnergy = 0.;          //  ZPA Energy
 696   Float_t  zem1Energy = 0.;         //  ZEM1 Energy
 697   Float_t  zem2Energy = 0.;         //  ZEM2 Energy
 698   Bool_t   zdcEnergyCal = kFALSE;   // if zdc is calibrated (in pass2)
 699
 700   Int_t    nTracks = 0;             //  no. tracks
 701   Int_t    nTracklets = 0;          //  no. tracklets
 702   Int_t    nClusters[6] = {0};      //  no. clusters on 6 ITS layers
 703   Int_t    nChips[2];               //  no. chips on 2 SPD layers
 704   Float_t  spdCorr =0;              //  corrected spd2 multiplicity
 705
 706   Float_t  multV0A  = 0;            //  multiplicity from V0 reco side A
 707   Float_t  multV0C  = 0;            //  multiplicity from V0 reco side C
 708   Short_t  multV0AOnline  = 0;      //  multiplicity from V0 reco side A
 709   Short_t  multV0COnline  = 0;      //  multiplicity from V0 reco side C
 710   Float_t  v0Corr = 0;               // corrected V0 multiplicity (used for MC)
 711
 712   Float_t  multFMDA = 0;            //  multiplicity from FMD on detector A
 713   Float_t  multFMDC = 0;            //  multiplicity from FMD on detector C
 714
 715   Float_t zvtx =0;                  // z-vertex SPD
 716   Int_t zvtxNcont =0;               // contributors to z-vertex SPD
 717
 718   Float_t zvtxT0 =0;                // z-vertex T0
 719
 720
 721   AliCentrality *esdCent = 0;
 722
 723   if(fAnalysisInput.CompareTo("ESD")==0){
 724
 725     AliVEvent* event = InputEvent();
 726     AliESDEvent* esd = dynamic_cast<AliESDEvent*>(event);
 727     if (!esd) {
 728       AliError("No ESD Event");
 729       return;
 730     }
 731
 732     LoadBranches();
 733
 734     if (SetupRun(esd)<0) {
 735       AliError("Centrality File not available for this run");
 736       return;
 737     }
 738
 739     esdCent = esd->GetCentrality();
 740
 741     // ***** V0 info
 742     AliESDVZERO* esdV0 = esd->GetVZEROData();
 743     if (!esdV0)
 744       {
 745         AliError("AliESDVZERO not available");
 746         return;
 747       }
 748     multV0A=esdV0->GetMTotV0A();
 749     multV0C=esdV0->GetMTotV0C();
 750     v0Corr = multV0A+multV0C;
 751
 752     multV0AOnline=esdV0->GetTriggerChargeA();
 753     multV0COnline=esdV0->GetTriggerChargeC();
 754
 755     // ***** T0 info
 756     const AliESDTZERO* esdT0 = esd->GetESDTZERO();
 757     if (!esdT0)
 758       {
 759         AliError("AliESDTZERO not available");
 760         return;
 761       }
 762     Int_t trig=esdT0->GetT0Trig();
 763     Bool_t kT0BB = kFALSE;
 764     if(trig&1) kT0BB=kTRUE;
 765     zvtxT0=esdT0->GetT0zVertex();
 766
 767
 768     // ***** Trigger info
 769     fIsSelected = ((esdV0->GetV0ADecision()==1) && (esdV0->GetV0CDecision()==1));
 770     TString trigStr(esd->GetFiredTriggerClasses());
 771
 772     fCVHN=kFALSE; fCVLN=kFALSE; fCCENT=kFALSE; fCSEMI=kFALSE;
 773     fMSL=kFALSE;  fMSH=kFALSE;  fMUL=kFALSE;   fMLL=kFALSE;
 774     fEJE=kFALSE;  fEGA=kFALSE;  fPHS=kFALSE;
 775
 776     if ( (trigStr.Contains("-B-")) &&  (trigStr.Contains("CVHN")) && (fIsSelected))
 777       fCVHN=kTRUE;
 778     if ( (trigStr.Contains("-B-")) &&  (trigStr.Contains("CVLN")) && (fIsSelected))
 779       fCVLN=kTRUE;
 780     if ( (trigStr.Contains("-B-")) &&  (trigStr.Contains("CCENT")) && (fIsSelected))
 781       fCCENT=kTRUE;
 782     if ( (trigStr.Contains("-B-")) &&  (trigStr.Contains("CSEMI")) && (fIsSelected))
 783       fCSEMI=kTRUE;
 784
 785     if ( (trigStr.Contains("-B-")) &&  (trigStr.Contains("CPBI1MSL")) && (fIsSelected))
 786       fMSL=kTRUE;
 787     if ( (trigStr.Contains("-B-")) &&  (trigStr.Contains("CPBI1MSH")) && (fIsSelected))
 788       fMSH=kTRUE;
 789     if ( (trigStr.Contains("-B-")) &&  (trigStr.Contains("CPBI1MUL")) && (fIsSelected))
 790       fMUL=kTRUE;
 791     if ( (trigStr.Contains("-B-")) &&  (trigStr.Contains("CPBI1MLL")) && (fIsSelected))
 792       fMLL=kTRUE;
 793     if ( (trigStr.Contains("-B-")) &&  (trigStr.Contains("CPBI1EJE")) && (fIsSelected))
 794       fEJE=kTRUE;
 795     if ( (trigStr.Contains("-B-")) &&  (trigStr.Contains("CPBI1EGA")) && (fIsSelected))
 796       fEGA=kTRUE;
 797     if ( (trigStr.Contains("-B-")) &&  (trigStr.Contains("CPBI1PHS")) && (fIsSelected))
 798       fPHS=kTRUE;
 799
 800     fCVHNbit=kFALSE;    fCVLNbit=kFALSE;       fCCENTbit=kFALSE;    fCSEMIbit=kFALSE;
 801     if (esdV0->GetTriggerBits() & (1<<8))
 802       fCVHNbit=kTRUE;
 803     if (esdV0->GetTriggerBits() & (1<<6))
 804       fCVLNbit=kTRUE;
 805
 806     if (kT0BB && fCVHNbit)
 807       fCCENTbit=kTRUE;
 808     if (kT0BB && fCVLNbit)
 809       fCSEMIbit=kTRUE;
 810
 811
 812     // ***** Vertex Info
 813     const AliESDVertex* vtxESD = esd->GetPrimaryVertexSPD();
 814     zvtx        = vtxESD->GetZ();
 815     zvtxNcont   = vtxESD->GetNContributors();
 816
 817     // ***** CB info (tracklets, clusters, chips)
 818     //nTracks    = event->GetNumberOfTracks();
 819     nTracks    = fTrackCuts ? (Short_t)fTrackCuts->GetReferenceMultiplicity(esd,kTRUE):-1;
 820
 821     const AliMultiplicity *mult = esd->GetMultiplicity();
 822
 823     nTracklets = mult->GetNumberOfTracklets();
 824
 825     for(Int_t ilay=0; ilay<6; ilay++){
 826       nClusters[ilay] = mult->GetNumberOfITSClusters(ilay);
 827     }
 828
 829     for(Int_t ilay=0; ilay<2; ilay++){
 830       nChips[ilay] = mult->GetNumberOfFiredChips(ilay);
 831     }
 832
 833     spdCorr = AliESDUtils::GetCorrSPD2(nClusters[1],zvtx);
 834
 835     // ***** FMD info
 836     AliESDFMD *fmd = esd->GetFMDData();
 837     Float_t totalMultA = 0;
 838     Float_t totalMultC = 0;
 839     const Float_t fFMDLowCut = 0.4;
 840
 841     for(UShort_t det=1;det<=3;det++) {
 842       Int_t nRings = (det==1 ? 1 : 2);
 843       for (UShort_t ir = 0; ir < nRings; ir++) {
 844         Char_t   ring = (ir == 0 ? 'I' : 'O');
 845         UShort_t nsec = (ir == 0 ? 20  : 40);
 846         UShort_t nstr = (ir == 0 ? 512 : 256);
 847         for(UShort_t sec =0; sec < nsec;  sec++)  {
 848           for(UShort_t strip = 0; strip < nstr; strip++) {
 849
 850             Float_t fmdMult = fmd->Multiplicity(det,ring,sec,strip);
 851             if(fmdMult == 0 || fmdMult == AliESDFMD::kInvalidMult) continue;
 852
 853             Float_t nParticles=0;
 854
 855             if(fmdMult > fFMDLowCut) {
 856               nParticles = 1.;
 857             }
 858
 859             if (det<3) totalMultA = totalMultA + nParticles;
 860             else totalMultC = totalMultC + nParticles;
 861
 862           }
 863         }
 864       }
 865     }
 866     multFMDA = totalMultA;
 867     multFMDC = totalMultC;
 868
 869     // ***** ZDC info
 870     AliESDZDC *esdZDC = esd->GetESDZDC();
 871     zdcEnergyCal = esdZDC->AliESDZDC::TestBit(AliESDZDC::kEnergyCalibratedSignal);
 872     if (zdcEnergyCal) {
 873       zncEnergy = (Float_t) (esdZDC->GetZDCN1Energy());
 874       zpcEnergy = (Float_t) (esdZDC->GetZDCP1Energy());
 875       znaEnergy = (Float_t) (esdZDC->GetZDCN2Energy());
 876       zpaEnergy = (Float_t) (esdZDC->GetZDCP2Energy());
 877     } else {
 878       zncEnergy = (Float_t) (esdZDC->GetZDCN1Energy())/8.;
 879       zpcEnergy = (Float_t) (esdZDC->GetZDCP1Energy())/8.;
 880       znaEnergy = (Float_t) (esdZDC->GetZDCN2Energy())/8.;
 881       zpaEnergy = (Float_t) (esdZDC->GetZDCP2Energy())/8.;
 882     }
 883     zem1Energy = (Float_t) (esdZDC->GetZDCEMEnergy(0))/8.;
 884     zem2Energy = (Float_t) (esdZDC->GetZDCEMEnergy(1))/8.;
 885
 886   }
 887
 888   else if(fAnalysisInput.CompareTo("AOD")==0){
 889     //AliAODEvent *aod =  dynamic_cast<AliAODEvent*> (InputEvent());
 890     // to be implemented
 891     printf("  AOD analysis not yet implemented!!!\n\n");
 892     return;
 893   }
 894
 895   // ***** Scaling for MC
 896   if (fIsMCInput) {
 897     fUseScaling=kFALSE;
 898     v0Corr  = Short_t((multV0A+multV0C)  * fV0MScaleFactorMC);
 899   }
 900   // ***** Scaling for Data
 901   if (fUseScaling) {
 902     v0Corr  = Short_t(v0Corr /   fV0MScaleFactor);
 903     spdCorr = spdCorr / fSPDScaleFactor;
 904     nTracks = Int_t(nTracks /   fTPCScaleFactor);
 905   }
 906
 907   // ***** Centrality Selection
 908   if(fHtempV0M) fCentV0M = fHtempV0M->GetBinContent(fHtempV0M->FindBin((v0Corr)));
 909   if(fHtempFMD) fCentFMD = fHtempFMD->GetBinContent(fHtempFMD->FindBin((multFMDA+multFMDC)));
 910   if(fHtempTRK) fCentTRK = fHtempTRK->GetBinContent(fHtempTRK->FindBin(nTracks));
 911   if(fHtempTKL) fCentTKL = fHtempTKL->GetBinContent(fHtempTKL->FindBin(nTracklets));
 912   if(fHtempCL0) fCentCL0 = fHtempCL0->GetBinContent(fHtempCL0->FindBin(nClusters[0]));
 913   if(fHtempCL1) fCentCL1 = fHtempCL1->GetBinContent(fHtempCL1->FindBin(spdCorr));
 914
 915   if(fHtempV0MvsFMD) fCentV0MvsFMD = fHtempV0MvsFMD->GetBinContent(fHtempV0MvsFMD->FindBin((multV0A+multV0C)));
 916   if(fHtempTKLvsV0M) fCentTKLvsV0M = fHtempTKLvsV0M->GetBinContent(fHtempTKLvsV0M->FindBin(nTracklets));
 917   if(fHtempZEMvsZDC) fCentZEMvsZDC = fHtempZEMvsZDC->GetBinContent(fHtempZEMvsZDC->FindBin(zem1Energy+zem2Energy,zncEnergy+znaEnergy+zpcEnergy+zpaEnergy));
 918
 919
 920   // ***** Cleaning
 921   if (fUseCleaning) {
 922     fQuality=0;
 923
 924     // ***** vertex
 925     if (TMath::Abs(zvtx)>fZVCut || zvtxNcont<1) fQuality += 1;
 926
 927     // ***** outliers, skip in case of MC input
 928     if (!fIsMCInput) {
 929       // **** V0 vs SPD
 930       if (IsOutlierV0MSPD(spdCorr, v0Corr, int(fCentV0M))) fQuality  += 2;
 931       // ***** V0 vs TPC
 932       if (IsOutlierV0MTPC(nTracks, v0Corr, int(fCentV0M))) fQuality  += 4;
 933       // ***** V0 vs ZDC
 934       if (IsOutlierV0MZDC((zncEnergy+znaEnergy+zpcEnergy+zpaEnergy), v0Corr) &&
 935           (zdcEnergyCal==kFALSE) ) fQuality  += 8;
 936       if (IsOutlierV0MZDCECal((zncEnergy+znaEnergy+zpcEnergy+zpaEnergy), v0Corr) &&
 937           (zdcEnergyCal==kTRUE) ) fQuality  += 8;
 938     }
 939   } else {
 940     fQuality = 0;
 941   }
 942
 943
 944   if (esdCent) {
 945     esdCent->SetQuality(fQuality);
 946     esdCent->SetCentralityV0M(fCentV0M);
 947     esdCent->SetCentralityFMD(fCentFMD);
 948     esdCent->SetCentralityTRK(fCentTRK);
 949     esdCent->SetCentralityTKL(fCentTKL);
 950     esdCent->SetCentralityCL0(fCentCL0);
 951     esdCent->SetCentralityCL1(fCentCL1);
 952     esdCent->SetCentralityV0MvsFMD(fCentV0MvsFMD);
 953     esdCent->SetCentralityTKLvsV0M(fCentTKLvsV0M);
 954     esdCent->SetCentralityZEMvsZDC(fCentZEMvsZDC);
 955   }
 956
 957   // filling QA histograms
 958   if (fFillHistos) {
 959     if (fCVHN)   fHOutCentV0MCVHN->Fill(fCentV0M);
 960     if (fCVLN)   fHOutCentV0MCVLN->Fill(fCentV0M);
 961     if (fCCENT)  fHOutCentV0MCCENT->Fill(fCentV0M);
 962     if (fCSEMI)  fHOutCentV0MCSEMI->Fill(fCentV0M);
 963     if (fMSL) fHOutCentV0MMSL->Fill(fCentV0M);
 964     if (fMSH) fHOutCentV0MMSH->Fill(fCentV0M);
 965     if (fMUL) fHOutCentV0MMUL->Fill(fCentV0M);
 966     if (fMLL) fHOutCentV0MMLL->Fill(fCentV0M);
 967     if (fEJE) fHOutCentV0MEJE->Fill(fCentV0M);
 968     if (fEGA) fHOutCentV0MEGA->Fill(fCentV0M);
 969     if (fPHS) fHOutCentV0MPHS->Fill(fCentV0M);
 970
 971     if (((AliInputEventHandler*)(AliAnalysisManager::GetAnalysisManager()->GetInputEventHandler()))->IsEventSelected() & AliVEvent::kMB) { // fill the QA histograms only for MB events!
 972       fHOutQuality->Fill(fQuality);
 973       fHOutVertex->Fill(zvtx);
 974       fHOutVertexT0->Fill(zvtxT0);
 975
 976       if (fQuality==0) {
 977         fHOutCentV0M->Fill(fCentV0M);
 978
 979         if (fCVHNbit)  fHOutCentV0MCVHNinMB->Fill(fCentV0M);
 980         if (fCVLNbit)  fHOutCentV0MCVLNinMB->Fill(fCentV0M);
 981         if (fCCENTbit) fHOutCentV0MCCENTinMB->Fill(fCentV0M);
 982         if (fCSEMIbit) fHOutCentV0MCSEMIinMB->Fill(fCentV0M);
 983         if (fMSL) fHOutCentV0MMSLinMB->Fill(fCentV0M);
 984         if (fMSH) fHOutCentV0MMSHinMB->Fill(fCentV0M);
 985         if (fMUL) fHOutCentV0MMULinMB->Fill(fCentV0M);
 986         if (fMLL) fHOutCentV0MMLLinMB->Fill(fCentV0M);
 987         if (fEJE) fHOutCentV0MEJEinMB->Fill(fCentV0M);
 988         if (fEGA) fHOutCentV0MEGAinMB->Fill(fCentV0M);
 989         if (fPHS) fHOutCentV0MPHSinMB->Fill(fCentV0M);
 990
 991
 992         fHOutCentFMD->Fill(fCentFMD);
 993         fHOutCentTRK->Fill(fCentTRK);
 994         fHOutCentTKL->Fill(fCentTKL);
 995         fHOutCentCL0->Fill(fCentCL0);
 996         fHOutCentCL1->Fill(fCentCL1);
 997         fHOutCentV0MvsFMD->Fill(fCentV0MvsFMD);
 998         fHOutCentTKLvsV0M->Fill(fCentTKLvsV0M);
 999         fHOutCentZEMvsZDC->Fill(fCentZEMvsZDC);
1000         fHOutCentV0MvsCentCL1->Fill(fCentV0M,fCentCL1);
1001         fHOutCentV0MvsCentTRK->Fill(fCentV0M,fCentTRK);
1002         fHOutCentTRKvsCentCL1->Fill(fCentTRK,fCentCL1);
1003         fHOutCentV0MvsCentZDC->Fill(fCentV0M,fCentZEMvsZDC);
1004         fHOutMultV0M->Fill(multV0A+multV0C);
1005         fHOutMultV0O->Fill(multV0AOnline+multV0COnline);
1006         fHOutMultFMD->Fill(multFMDA+multFMDC);
1007         fHOutMultTRK->Fill(nTracks);
1008         fHOutMultTKL->Fill(nTracklets);
1009         fHOutMultCL0->Fill(nClusters[0]);
1010         fHOutMultCL1->Fill(spdCorr);
1011         fHOutMultV0MvsZDN->Fill(v0Corr,(zncEnergy+znaEnergy));
1012         fHOutMultZEMvsZDN->Fill((zem1Energy+zem2Energy),(zncEnergy+znaEnergy));
1013         fHOutMultV0MvsZDC->Fill(v0Corr,(zncEnergy+znaEnergy+zpcEnergy+zpaEnergy));
1014         fHOutMultZEMvsZDC->Fill((zem1Energy+zem2Energy),(zncEnergy+znaEnergy+zpcEnergy+zpaEnergy));
1015         fHOutMultZEMvsZDCw->Fill((zem1Energy+zem2Energy),(zncEnergy+znaEnergy+zpcEnergy+zpaEnergy),fCentV0M);
1016         fHOutMultV0MvsCL1->Fill(v0Corr,spdCorr);
1017         fHOutMultV0MvsTRK->Fill(v0Corr,nTracks);
1018         fHOutMultTRKvsCL1->Fill(nTracks,spdCorr);
1019         fHOutMultV0MvsV0O->Fill(v0Corr,(multV0AOnline+multV0COnline));
1020         fHOutMultV0OvsCL1->Fill((multV0AOnline+multV0COnline),spdCorr);
1021         fHOutMultV0OvsTRK->Fill((multV0AOnline+multV0COnline),nTracks);
1022       } else if (fQuality%2 == 0) {
1023         fHOutCentV0Mqual1->Fill(fCentV0M);
1024         fHOutCentTRKqual1->Fill(fCentTRK);
1025         fHOutCentCL1qual1->Fill(fCentCL1);
1026         fHOutMultV0MvsCL1qual1->Fill(v0Corr,spdCorr);
1027         fHOutMultV0MvsTRKqual1->Fill(v0Corr,nTracks);
1028         fHOutMultTRKvsCL1qual1->Fill(nTracks,spdCorr);
1029       } else {
1030         fHOutCentV0Mqual2->Fill(fCentV0M);
1031         fHOutCentTRKqual2->Fill(fCentTRK);
1032         fHOutCentCL1qual2->Fill(fCentCL1);
1033         fHOutMultV0MvsCL1qual2->Fill(v0Corr,spdCorr);
1034         fHOutMultV0MvsTRKqual2->Fill(v0Corr,nTracks);
1035         fHOutMultTRKvsCL1qual2->Fill(nTracks,spdCorr);
1036       }
1037     }
1038     PostData(1, fOutputList);
1039   }
1040 }
1041 //________________________________________________________________________
1042 void AliCentralitySelectionTask::Terminate(Option_t */*option*/)
1043 {
1044   // Terminate analysis
1045 }
1046 //________________________________________________________________________
1047 Int_t AliCentralitySelectionTask::SetupRun(const AliESDEvent* const esd)
1048 {
1049   // Setup files for run
1050
1051   if (!esd)
1052     return -1;
1053
1054   // check if something to be done
1055   if (fCurrentRun == esd->GetRunNumber())
1056     return 0;
1057   else
1058     fCurrentRun = esd->GetRunNumber();
1059
1060   TString fileName =(Form("%s/COMMON/CENTRALITY/data/centrality.root", AliAnalysisManager::GetOADBPath()));
1061   AliInfo(Form("Setup Centrality Selection for run %d with file %s\n",fCurrentRun,fileName.Data()));
1062
1063   AliOADBContainer *con = new AliOADBContainer("OADB");
1064   con->InitFromFile(fileName,"Centrality");
1065
1066   AliOADBCentrality*  centOADB = 0;
1067   centOADB = (AliOADBCentrality*)(con->GetObject(fCurrentRun));
1068   if (!centOADB) {
1069     AliWarning(Form("Centrality OADB does not exist for run %d, using Default \n",fCurrentRun ));
1070     centOADB  = (AliOADBCentrality*)(con->GetDefaultObject("oadbDefault"));
1071   }
1072
1073   // modes
1074   fUseScaling   = centOADB->UseScaling();
1075   fUseCleaning  = centOADB->UseCleaning();
1076
1077   // cuts
1078   fZVCut        = centOADB->ZVCut();
1079   fOutliersCut  = centOADB->OutliersCut();
1080
1081   // centrality histos
1082   fHtempV0M       = centOADB->V0hist();
1083   fHtempTRK       = centOADB->TPChist();
1084   fHtempCL1       = centOADB->SPDhist();
1085   fHtempZEMvsZDC  = centOADB->ZEMvsZDChist();
1086
1087   TString path = gSystem->ExpandPathName(fileName.Data());
1088   if (!fHtempV0M) AliWarning(Form("Calibration for V0M does not exist in %s", path.Data()));
1089   if (!fHtempTRK) AliWarning(Form("Calibration for TRK does not exist in %s", path.Data()));
1090   if (!fHtempCL1) AliWarning(Form("Calibration for CL1 does not exist in %s", path.Data()));
1091   if (!fHtempZEMvsZDC) AliWarning(Form("Calibration for ZEMvsZDC does not exist in %s", path.Data()));
1092
1093   // scale factors
1094   fV0MScaleFactor    = centOADB->V0MScaleFactor();
1095   fSPDScaleFactor    = centOADB->SPDScaleFactor();
1096   fTPCScaleFactor    = centOADB->TPCScaleFactor();
1097   fV0MScaleFactorMC  = centOADB->V0MScaleFactorMC();
1098
1099   // outliers parameters
1100   fV0MSPDOutlierPar0 = centOADB->V0MSPDOutlierPar0();
1101   fV0MSPDOutlierPar1 = centOADB->V0MSPDOutlierPar1();
1102   fV0MTPCOutlierPar0 = centOADB->V0MTPCOutlierPar0();
1103   fV0MTPCOutlierPar1 = centOADB->V0MTPCOutlierPar1();
1104
1105   fV0MSPDSigmaOutlierPar0 = centOADB->V0MSPDSigmaOutlierPar0();
1106   fV0MSPDSigmaOutlierPar1 = centOADB->V0MSPDSigmaOutlierPar1();
1107   fV0MSPDSigmaOutlierPar2 = centOADB->V0MSPDSigmaOutlierPar2();
1108   fV0MTPCSigmaOutlierPar0 = centOADB->V0MTPCSigmaOutlierPar0();
1109   fV0MTPCSigmaOutlierPar1 = centOADB->V0MTPCSigmaOutlierPar1();
1110   fV0MTPCSigmaOutlierPar2 = centOADB->V0MTPCSigmaOutlierPar2();
1111
1112   fV0MZDCOutlierPar0 =      centOADB->V0MZDCOutlierPar0();
1113   fV0MZDCOutlierPar1 =      centOADB->V0MZDCOutlierPar1();
1114   fV0MZDCEcalOutlierPar0 =  centOADB->V0MZDCEcalOutlierPar0();
1115   fV0MZDCEcalOutlierPar1 =  centOADB->V0MZDCEcalOutlierPar1();
1116
1117
1118
1119   return 0;
1120 }
1121
1122
1123
1124 //________________________________________________________________________
1125 Bool_t AliCentralitySelectionTask::IsOutlierV0MSPD(Float_t spd, Float_t v0, Int_t cent) const
1126 {
1127   // Clean outliers
1128   Float_t val = fV0MSPDOutlierPar0 +  fV0MSPDOutlierPar1 * v0;
1129   Float_t spdSigma = fV0MSPDSigmaOutlierPar0 + fV0MSPDSigmaOutlierPar1*cent + fV0MSPDSigmaOutlierPar2*cent*cent;
1130   if ( TMath::Abs(spd-val) > fOutliersCut*spdSigma )
1131     return kTRUE;
1132   else
1133     return kFALSE;
1134 }
1135
1136 //________________________________________________________________________
1137 Bool_t AliCentralitySelectionTask::IsOutlierV0MTPC(Int_t tracks, Float_t v0, Int_t cent) const
1138 {
1139   // Clean outliers
1140   Float_t val = fV0MTPCOutlierPar0 +  fV0MTPCOutlierPar1 * v0;
1141   Float_t tpcSigma = fV0MTPCSigmaOutlierPar0 + fV0MTPCSigmaOutlierPar1*cent + fV0MTPCSigmaOutlierPar2*cent*cent;
1142   if ( TMath::Abs(tracks-val) > fOutliersCut*tpcSigma )
1143     return kTRUE;
1144   else
1145     return kFALSE;
1146 }
1147
1148 //________________________________________________________________________
1149 Bool_t AliCentralitySelectionTask::IsOutlierV0MZDC(Float_t zdc, Float_t v0) const
1150 {
1151   // Clean outliers
1152   Float_t val = fV0MZDCOutlierPar0 + fV0MZDCOutlierPar1 * v0;
1153   if (zdc >  val)
1154     return kTRUE;
1155   else
1156   return kFALSE;
1157 }
1158
1159 //________________________________________________________________________
1160 Bool_t AliCentralitySelectionTask::IsOutlierV0MZDCECal(Float_t zdc, Float_t v0) const
1161 {
1162   // Clean outliers
1163   Float_t val = fV0MZDCEcalOutlierPar0 + fV0MZDCEcalOutlierPar1 * v0;
1164   if (zdc >  val)
1165     return kTRUE;
1166   else
1167     return kFALSE;
1168 }
1169
1170