]> git.uio.no Git - u/mrichter/AliRoot.git/blobdiff - ANALYSIS/AliCentralitySelectionTask.cxx
changing default for cut on tracks which have shared TPC clusters to ***OFF*** (was...
[u/mrichter/AliRoot.git] / ANALYSIS / AliCentralitySelectionTask.cxx
index fdcce8b9d6c67c6599db3410d9fdc97b1a998109..c7dc6e20c2a4594691d53a84802b474f08168988 100644 (file)
 //   author: Alberica Toia
 //*****************************************************
 
+#include "AliCentralitySelectionTask.h"
+
 #include <TTree.h>
 #include <TList.h>
 #include <TH1F.h>
 #include <TH2F.h>
+#include <TF1.h>
 #include <TProfile.h>
 #include <TFile.h>
 #include <TObjString.h>
@@ -30,6 +33,7 @@
 #include <TCanvas.h>
 #include <TROOT.h>
 #include <TDirectory.h>
+#include <TSystem.h>
 #include <iostream>
 
 #include "AliAnalysisManager.h"
@@ -41,7 +45,8 @@
 #include "AliESDZDC.h"
 #include "AliESDFMD.h"
 #include "AliESDVZERO.h"
-#include "AliESDCentrality.h"
+#include "AliCentrality.h"
+#include "AliESDtrackCuts.h"
 #include "AliMultiplicity.h"
 #include "AliAODHandler.h"
 #include "AliAODEvent.h"
@@ -60,7 +65,7 @@
 #include "AliPhysicsSelectionTask.h"
 #include "AliPhysicsSelection.h"
 #include "AliBackgroundSelection.h"
-#include "AliCentralitySelectionTask.h"
+#include "AliESDUtils.h"
 
 ClassImp(AliCentralitySelectionTask)
 
@@ -73,12 +78,15 @@ AliAnalysisTaskSE(),
   fIsMCInput(kFALSE),
   fFile(0),
   fFile2(0),
-  fCentfilename(""),
-  fCentfilename2(""),
-  fFileList(new TList),
-  fFileList2(new TList),
   fCurrentRun(-1),
   fRunNo(-1),
+  fLowRunN(0),
+  fHighRunN(0),
+  fUseScaling(0),
+  fTrackCuts(0),
+  fZVCut(10),
+  fOutliersCut(5),
+  fQuality(0),
   fCentV0M(0),
   fCentFMD(0),
   fCentTRK(0),
@@ -96,13 +104,52 @@ AliAnalysisTaskSE(),
   fHtempCL1(0),
   fHtempV0MvsFMD(0),
   fHtempTKLvsV0M(0),
-  fHtempZEMvsZDC(0)
+  fHtempZEMvsZDC(0),
+  fOutputList(0),
+  fHOutCentV0M     (0),
+  fHOutCentFMD     (0),
+  fHOutCentTRK     (0),
+  fHOutCentTKL     (0),
+  fHOutCentCL0     (0),
+  fHOutCentCL1     (0),
+  fHOutCentV0MvsFMD(0),
+  fHOutCentTKLvsV0M(0),
+  fHOutCentZEMvsZDC(0),
+  fHOutMultV0M(0),
+  fHOutMultV0R(0),
+  fHOutMultFMD(0),
+  fHOutMultTRK(0),
+  fHOutMultTKL(0),
+  fHOutMultCL0(0),
+  fHOutMultCL1(0),
+  fHOutMultV0MvsZDN(0),
+  fHOutMultZEMvsZDN(0),
+  fHOutMultV0MvsZDC(0),
+  fHOutMultZEMvsZDC(0),
+  fHOutMultV0MvsCL1(0),
+  fHOutMultV0MvsTRK(0),
+  fHOutMultTRKvsCL1(0),
+  fHOutCentV0M_qual1(0),
+  fHOutCentTRK_qual1(0),
+  fHOutCentCL1_qual1(0),
+  fHOutCentV0M_qual2(0),
+  fHOutCentTRK_qual2(0),
+  fHOutCentCL1_qual2(0),
+  fHOutQuality(0),
+  fHOutVertex(0)
 {   
   // Default constructor
-  fFileList->SetOwner();
-  fFileList2->SetOwner();
-
   AliInfo("Centrality Selection enabled.");
+  fLowRunN =136851;
+  fHighRunN=139517;
+
+  for (Int_t i=0; i < 2667; i++) {
+    V0MScaleFactor[i]=0.0;
+    SPDScaleFactor[i]=0.0;
+    TPCScaleFactor[i]=0.0;
+    V0MScaleFactorMC[i]=0.0;
+  }
+  fUseScaling=kTRUE;
 }   
 
 //________________________________________________________________________
@@ -113,12 +160,15 @@ AliCentralitySelectionTask::AliCentralitySelectionTask(const char *name):
   fIsMCInput(kFALSE),
   fFile(0),
   fFile2(0),
-  fCentfilename(""),
-  fCentfilename2(""),
-  fFileList(new TList),
-  fFileList2(new TList),
   fCurrentRun(-1),
   fRunNo(-1),
+  fLowRunN(0),
+  fHighRunN(0),
+  fUseScaling(0),
+  fTrackCuts(0),
+  fZVCut(10),
+  fOutliersCut(5),
+  fQuality(0),
   fCentV0M(0),
   fCentFMD(0),
   fCentTRK(0),
@@ -136,13 +186,53 @@ AliCentralitySelectionTask::AliCentralitySelectionTask(const char *name):
   fHtempCL1(0),
   fHtempV0MvsFMD(0),
   fHtempTKLvsV0M(0),
-  fHtempZEMvsZDC(0)
+  fHtempZEMvsZDC(0),
+  fOutputList(0),
+  fHOutCentV0M     (0),
+  fHOutCentFMD     (0),
+  fHOutCentTRK     (0),
+  fHOutCentTKL     (0),
+  fHOutCentCL0     (0),
+  fHOutCentCL1     (0),
+  fHOutCentV0MvsFMD(0),
+  fHOutCentTKLvsV0M(0),
+  fHOutCentZEMvsZDC(0),
+  fHOutMultV0M(0),
+  fHOutMultV0R(0),
+  fHOutMultFMD(0),
+  fHOutMultTRK(0),
+  fHOutMultTKL(0),
+  fHOutMultCL0(0),
+  fHOutMultCL1(0),
+  fHOutMultV0MvsZDN(0),
+  fHOutMultZEMvsZDN(0),
+  fHOutMultV0MvsZDC(0),
+  fHOutMultZEMvsZDC(0),
+  fHOutMultV0MvsCL1(0),
+  fHOutMultV0MvsTRK(0),
+  fHOutMultTRKvsCL1(0),
+  fHOutCentV0M_qual1(0),
+  fHOutCentTRK_qual1(0),
+  fHOutCentCL1_qual1(0),
+  fHOutCentV0M_qual2(0),
+  fHOutCentTRK_qual2(0),
+  fHOutCentCL1_qual2(0),
+  fHOutQuality(0),
+  fHOutVertex(0)
 {
-  // Default constructor
-  fFileList->SetOwner();
-  fFileList2->SetOwner();
+    // Default constructor
+    fLowRunN =136851;
+    fHighRunN=139517;
 
   AliInfo("Centrality Selection enabled.");
+  DefineOutput(1, TList::Class());
+  for (Int_t i=0; i<2667; i++) {
+    V0MScaleFactor[i]=0.0;
+    SPDScaleFactor[i]=0.0;
+    TPCScaleFactor[i]=0.0;
+    V0MScaleFactorMC[i]=0.0;
+  }
+  fUseScaling=kTRUE;
 }
 
 //________________________________________________________________________
@@ -163,12 +253,15 @@ AliCentralitySelectionTask::AliCentralitySelectionTask(const AliCentralitySelect
   fIsMCInput(ana.fIsMCInput),
   fFile(ana.fFile),
   fFile2(ana.fFile2),
-  fCentfilename(ana.fCentfilename),
-  fCentfilename2(ana.fCentfilename2),
-  fFileList(ana.fFileList),
-  fFileList2(ana.fFileList2),
   fCurrentRun(ana.fCurrentRun),
   fRunNo(ana.fRunNo),
+  fLowRunN(ana.fLowRunN),
+  fHighRunN(ana.fHighRunN),
+  fUseScaling(ana.fUseScaling),
+  fTrackCuts(ana.fTrackCuts),
+  fZVCut(ana.fZVCut),
+  fOutliersCut(ana.fOutliersCut),
+  fQuality(ana.fQuality),
   fCentV0M(ana.fCentV0M),
   fCentFMD(ana.fCentFMD),
   fCentTRK(ana.fCentTRK),
@@ -186,27 +279,56 @@ AliCentralitySelectionTask::AliCentralitySelectionTask(const AliCentralitySelect
   fHtempCL1(ana.fHtempCL1),
   fHtempV0MvsFMD(ana.fHtempV0MvsFMD),
   fHtempTKLvsV0M(ana.fHtempTKLvsV0M),
-  fHtempZEMvsZDC(ana.fHtempZEMvsZDC)
+  fHtempZEMvsZDC(ana.fHtempZEMvsZDC),
+  fOutputList(ana.fOutputList),
+  fHOutCentV0M     (ana.fHOutCentV0M     ),
+  fHOutCentFMD     (ana.fHOutCentFMD     ),
+  fHOutCentTRK     (ana.fHOutCentTRK     ),
+  fHOutCentTKL     (ana.fHOutCentTKL     ),
+  fHOutCentCL0     (ana.fHOutCentCL0     ),
+  fHOutCentCL1     (ana.fHOutCentCL1     ),
+  fHOutCentV0MvsFMD(ana.fHOutCentV0MvsFMD),
+  fHOutCentTKLvsV0M(ana.fHOutCentTKLvsV0M),
+  fHOutCentZEMvsZDC(ana.fHOutCentZEMvsZDC),
+  fHOutMultV0M(ana.fHOutMultV0M),
+  fHOutMultV0R(ana.fHOutMultV0R),
+  fHOutMultFMD(ana.fHOutMultFMD),
+  fHOutMultTRK(ana.fHOutMultTRK),
+  fHOutMultTKL(ana.fHOutMultTKL),
+  fHOutMultCL0(ana.fHOutMultCL0),
+  fHOutMultCL1(ana.fHOutMultCL1),
+  fHOutMultV0MvsZDN(ana.fHOutMultV0MvsZDN),
+  fHOutMultZEMvsZDN(ana.fHOutMultZEMvsZDN),
+  fHOutMultV0MvsZDC(ana.fHOutMultV0MvsZDC),
+  fHOutMultZEMvsZDC(ana.fHOutMultZEMvsZDC),
+  fHOutMultV0MvsCL1(ana.fHOutMultV0MvsCL1),
+  fHOutMultV0MvsTRK(ana.fHOutMultV0MvsTRK),
+  fHOutMultTRKvsCL1(ana.fHOutMultTRKvsCL1),
+  fHOutCentV0M_qual1(ana.fHOutCentV0M_qual1),
+  fHOutCentTRK_qual1(ana.fHOutCentTRK_qual1),
+  fHOutCentCL1_qual1(ana.fHOutCentCL1_qual1),
+  fHOutCentV0M_qual2(ana.fHOutCentV0M_qual2),
+  fHOutCentTRK_qual2(ana.fHOutCentTRK_qual2),
+  fHOutCentCL1_qual2(ana.fHOutCentCL1_qual2),
+  fHOutQuality(ana.fHOutQuality),
+  fHOutVertex(ana.fHOutVertex)
 {
   // Copy Constructor  
+    for (Int_t i=0; i<2667; i++) {
+       V0MScaleFactor[i]=0.0;
+       SPDScaleFactor[i]=0.0;
+       TPCScaleFactor[i]=0.0;
+       V0MScaleFactorMC[i]=0.0;
+    }
+
 }
  
 //________________________________________________________________________
 AliCentralitySelectionTask::~AliCentralitySelectionTask()
 {
-  // Destructor
-  
-  if (fFileList) {
-    fFileList->Clear();
-    delete fFileList;
-  }
-  fFileList = NULL;
-
-  if (fFileList2) {
-    fFileList2->Clear();
-    delete fFileList2;
-  }
-  fFileList2 = NULL;
+  // Destructor  
+  if (fOutputList && !AliAnalysisManager::GetAnalysisManager()->IsProofMode()) delete fOutputList;
+  if (fTrackCuts) delete fTrackCuts;
 }  
 
 //________________________________________________________________________
@@ -214,12 +336,85 @@ void AliCentralitySelectionTask::UserCreateOutputObjects()
 {  
   // Create the output containers
   if(fDebug>1) printf("AnalysisCentralitySelectionTask::UserCreateOutputObjects() \n");
+  AliLog::SetClassDebugLevel("AliCentralitySelectionTask", AliLog::kInfo);
 
-  if (fFileList->GetEntries() < 1) {
-    AliError("No Inputfiles Added");
-  }
+  fOutputList = new TList();
+  fOutputList->SetOwner();
+  fHOutCentV0M     = new TH1F("fHOutCentV0M","fHOutCentV0M; Centrality V0",501,0,101);
+  fHOutCentFMD     = new TH1F("fHOutCentFMD","fHOutCentFMD; Centrality FMD",501,0,101);
+  fHOutCentTRK     = new TH1F("fHOutCentTRK","fHOutCentTRK; Centrality TPC",501,0,101);
+  fHOutCentTKL     = new TH1F("fHOutCentTKL","fHOutCentTKL; Centrality tracklets",501,0,101);
+  fHOutCentCL0     = new TH1F("fHOutCentCL0","fHOutCentCL0; Centrality SPD inner",501,0,101);
+  fHOutCentCL1     = new TH1F("fHOutCentCL1","fHOutCentCL1; Centrality SPD outer",501,0,101);
+  fHOutCentV0MvsFMD= new TH1F("fHOutCentV0MvsFMD","fHOutCentV0MvsFMD; Centrality V0 vs FMD",501,0,101);
+  fHOutCentTKLvsV0M= new TH1F("fHOutCentTKLvsV0M","fHOutCentTKLvsV0M; Centrality tracklets vs V0",501,0,101);
+  fHOutCentZEMvsZDC= new TH1F("fHOutCentZEMvsZDC","fHOutCentZEMvsZDC; Centrality ZEM vs ZDC",501,0,101);
 
-  AliLog::SetClassDebugLevel("AliCentralitySelectionTask", AliLog::kInfo);
+  fHOutMultV0M = new TH1F("fHOutMultV0M","fHOutMultV0M; Multiplicity V0",25000,0,25000);
+  fHOutMultV0R = new TH1F("fHOutMultV0R","fHOutMultV0R; Multiplicity V0",25000,0,25000);
+  fHOutMultFMD = new TH1F("fHOutMultFMD","fHOutMultFMD; Multiplicity FMD",24000,0,24000);
+  fHOutMultTRK = new TH1F("fHOutMultTRK","fHOutMultTRK; Multiplicity TPC",4000,0,4000);
+  fHOutMultTKL = new TH1F("fHOutMultTKL","fHOutMultTKL; Multiplicity tracklets",5000,0,5000);
+  fHOutMultCL0 = new TH1F("fHOutMultCL0","fHOutMultCL0; Multiplicity SPD inner",7000,0,7000);
+  fHOutMultCL1 = new TH1F("fHOutMultCL1","fHOutMultCL1; Multiplicity SPD outer",7000,0,7000);
+  fHOutMultV0MvsZDN = new TH2F("fHOutMultV0MvsZDN","fHOutMultV0MvsZDN; Multiplicity V0; Energy ZDC-N",500,0,25000,500,0,180000);
+  fHOutMultZEMvsZDN = new TH2F("fHOutMultZEMvsZDN","fHOutMultZEMvsZDN; Energy ZEM; Energy ZDC-N",500,0,2500,500,0,180000);
+  fHOutMultV0MvsZDC = new TH2F("fHOutMultV0MvsZDC","fHOutMultV0MvsZDC; Multiplicity V0; Energy ZDC",500,0,25000,500,0,200000);
+  fHOutMultZEMvsZDC = new TH2F("fHOutMultZEMvsZDC","fHOutMultZEMvsZDC; Energy ZEM; Energy ZDC",500,0,2500,500,0,200000);
+  fHOutMultV0MvsCL1 = new TH2F("fHOutMultV0MvsCL1","fHOutMultV0MvsCL1; Multiplicity V0; Multiplicity SPD outer",2500,0,25000,700,0,7000);
+  fHOutMultV0MvsTRK = new TH2F("fHOutMultV0MvsTRK","fHOutMultV0MvsTRK; Multiplicity V0; Multiplicity TPC",2500,0,25000,400,0,4000);
+  fHOutMultTRKvsCL1 = new TH2F("fHOutMultTRKvsCL1","fHOutMultTRKvsCL1; Multiplicity TPC; Multiplicity SPD outer",400,0,4000,700,0,7000);
+
+  fHOutCentV0M_qual1 = new TH1F("fHOutCentV0M_qual1","fHOutCentV0M_qual1; Centrality V0",501,0,100);
+  fHOutCentTRK_qual1 = new TH1F("fHOutCentTRK_qual1","fHOutCentTRK_qual1; Centrality TPC",501,0,100);
+  fHOutCentCL1_qual1 = new TH1F("fHOutCentCL1_qual1","fHOutCentCL1_qual1; Centrality SPD outer",501,0,100);
+
+  fHOutCentV0M_qual2 = new TH1F("fHOutCentV0M_qual2","fHOutCentV0M_qual2; Centrality V0",101,-0.1,100.1);
+  fHOutCentTRK_qual2 = new TH1F("fHOutCentTRK_qual2","fHOutCentTRK_qual2; Centrality TPC",101,-0.1,100.1);
+  fHOutCentCL1_qual2 = new TH1F("fHOutCentCL1_qual2","fHOutCentCL1_qual2; Centrality SPD outer",101,-0.1,100.1);
+
+  fHOutQuality = new TH1F("fHOutQuality", "fHOutQuality", 10,-0.5,9.5);
+  fHOutVertex  = new TH1F("fHOutVertex", "fHOutVertex", 100,-20,20);
+
+  fOutputList->Add(  fHOutCentV0M     );
+  fOutputList->Add(  fHOutCentFMD     );
+  fOutputList->Add(  fHOutCentTRK     );
+  fOutputList->Add(  fHOutCentTKL     );
+  fOutputList->Add(  fHOutCentCL0     );
+  fOutputList->Add(  fHOutCentCL1     );
+  fOutputList->Add(  fHOutCentV0MvsFMD);
+  fOutputList->Add(  fHOutCentTKLvsV0M);
+  fOutputList->Add(  fHOutCentZEMvsZDC);
+  fOutputList->Add(  fHOutMultV0M); 
+  fOutputList->Add(  fHOutMultV0R); 
+  fOutputList->Add(  fHOutMultFMD); 
+  fOutputList->Add(  fHOutMultTRK); 
+  fOutputList->Add(  fHOutMultTKL); 
+  fOutputList->Add(  fHOutMultCL0); 
+  fOutputList->Add(  fHOutMultCL1); 
+  fOutputList->Add(  fHOutMultV0MvsZDN);
+  fOutputList->Add(  fHOutMultZEMvsZDN);
+  fOutputList->Add(  fHOutMultV0MvsZDC);
+  fOutputList->Add(  fHOutMultZEMvsZDC);
+  fOutputList->Add(  fHOutMultV0MvsCL1);
+  fOutputList->Add(  fHOutMultV0MvsTRK);
+  fOutputList->Add(  fHOutMultTRKvsCL1);
+  fOutputList->Add(  fHOutCentV0M_qual1 );
+  fOutputList->Add(  fHOutCentTRK_qual1 );
+  fOutputList->Add(  fHOutCentCL1_qual1 );                   
+  fOutputList->Add(  fHOutCentV0M_qual2 );
+  fOutputList->Add(  fHOutCentTRK_qual2 );
+  fOutputList->Add(  fHOutCentCL1_qual2 );
+  fOutputList->Add(  fHOutQuality );
+  fOutputList->Add(  fHOutVertex );
+
+
+  fTrackCuts = AliESDtrackCuts::GetStandardTPCOnlyTrackCuts();
+
+  PostData(1, fOutputList); 
+
+  MyInitScaleFactor();
+  if (fIsMCInput) MyInitScaleFactorMC();
 }
 
 //________________________________________________________________________
@@ -228,16 +423,17 @@ void AliCentralitySelectionTask::UserExec(Option_t */*option*/)
   // Execute analysis for current event:
   if(fDebug>1) printf(" **** AliCentralitySelectionTask::UserExec() \n");
   
-  Float_t  zncEnergy;          //  ZNC Energy
-  Float_t  zpcEnergy;          //  ZPC Energy
-  Float_t  znaEnergy;          //  ZNA Energy
-  Float_t  zpaEnergy;          //  ZPA Energy
+  Float_t  zncEnergy = 0.;          //  ZNC Energy
+  Float_t  zpcEnergy = 0.;          //  ZPC Energy
+  Float_t  znaEnergy = 0.;          //  ZNA Energy
+  Float_t  zpaEnergy = 0.;          //  ZPA Energy
   Float_t  zem1Energy = 0.;         //  ZEM1 Energy
   Float_t  zem2Energy = 0.;         //  ZEM2 Energy
-  
+  Bool_t   zdcEnergyCal = kFALSE;   // if zdc is calibrated (in pass2)
+
   Int_t    nTracks = 0;             //  no. tracks
   Int_t    nTracklets = 0;          //  no. tracklets
-  Int_t    nClusters[6];            //  no. clusters on 6 ITS layers
+  Int_t    nClusters[6] = {0};      //  no. clusters on 6 ITS layers
   Int_t    nChips[2];               //  no. chips on 2 SPD layers
   Float_t  spdCorr =0;              //  corrected spd2 multiplicity
 
@@ -251,15 +447,21 @@ void AliCentralitySelectionTask::UserExec(Option_t */*option*/)
 
   Float_t zvtx =0;                  // z-vertex SPD
  
-  AliESDCentrality *esdCent = 0;
+  AliCentrality *esdCent = 0;
 
   if(fAnalysisInput.CompareTo("ESD")==0){
 
     AliVEvent* event = InputEvent();
     AliESDEvent* esd = dynamic_cast<AliESDEvent*>(event);
-
-    if (SetupRun(esd))   
-      return;
+    if (!esd) {
+       AliError("No ESD Event");
+       return;
+    }
+    
+    if (fRunNo<=0) {
+      if (SetupRun(esd)<0)
+         AliFatal("Centrality File not available for this run");
+    }
 
     esdCent = esd->GetCentrality();
 
@@ -268,8 +470,8 @@ void AliCentralitySelectionTask::UserExec(Option_t */*option*/)
     multV0A=esdV0->GetMTotV0A();
     multV0C=esdV0->GetMTotV0C();
 
-    float v0CorrR;
-    v0Corr = (Short_t) (GetCorrV0(esd,v0CorrR));
+    Float_t v0CorrR;
+    v0Corr = (Short_t)AliESDUtils::GetCorrV0(esd,v0CorrR);
     v0CorrResc = (Short_t)v0CorrR;
 
     // ***** Vertex Info
@@ -277,7 +479,8 @@ void AliCentralitySelectionTask::UserExec(Option_t */*option*/)
     zvtx        = vtxESD->GetZ(); 
 
     // ***** CB info (tracklets, clusters, chips)
-    nTracks    = event->GetNumberOfTracks();     
+    //nTracks    = event->GetNumberOfTracks();     
+    nTracks    = fTrackCuts ? (Short_t)fTrackCuts->GetReferenceMultiplicity(esd,kTRUE):-1;
 
     const AliMultiplicity *mult = esd->GetMultiplicity();
 
@@ -291,7 +494,7 @@ void AliCentralitySelectionTask::UserExec(Option_t */*option*/)
       nChips[ilay] = mult->GetNumberOfFiredChips(ilay);
     }
 
-    spdCorr = GetCorrSPD2(nClusters[1],zvtx);    
+    spdCorr = AliESDUtils::GetCorrSPD2(nClusters[1],zvtx);    
 
     // ***** FMD info
     AliESDFMD *fmd = esd->GetFMDData();
@@ -329,94 +532,179 @@ void AliCentralitySelectionTask::UserExec(Option_t */*option*/)
     
     // ***** ZDC info
     AliESDZDC *esdZDC = esd->GetESDZDC();
-    zncEnergy = (Float_t) (esdZDC->GetZDCN1Energy());
-    zpcEnergy = (Float_t) (esdZDC->GetZDCP1Energy());
-    znaEnergy = (Float_t) (esdZDC->GetZDCN2Energy());
-    zpaEnergy = (Float_t) (esdZDC->GetZDCP2Energy());
-    zem1Energy = (Float_t) (esdZDC->GetZDCEMEnergy(0));
-    zem2Energy = (Float_t) (esdZDC->GetZDCEMEnergy(1));
-    
+    zdcEnergyCal = esdZDC->AliESDZDC::TestBit(AliESDZDC::kEnergyCalibratedSignal);
+    if (zdcEnergyCal) {      
+      zncEnergy = (Float_t) (esdZDC->GetZDCN1Energy());
+      zpcEnergy = (Float_t) (esdZDC->GetZDCP1Energy());
+      znaEnergy = (Float_t) (esdZDC->GetZDCN2Energy());
+      zpaEnergy = (Float_t) (esdZDC->GetZDCP2Energy());
+    } else {
+      zncEnergy = (Float_t) (esdZDC->GetZDCN1Energy())/8.;
+      zpcEnergy = (Float_t) (esdZDC->GetZDCP1Energy())/8.;
+      znaEnergy = (Float_t) (esdZDC->GetZDCN2Energy())/8.;
+      zpaEnergy = (Float_t) (esdZDC->GetZDCP2Energy())/8.;
+    }
+    zem1Energy = (Float_t) (esdZDC->GetZDCEMEnergy(0))/8.;
+    zem2Energy = (Float_t) (esdZDC->GetZDCEMEnergy(1))/8.;
+  
   }   
   else if(fAnalysisInput.CompareTo("AOD")==0){
     //AliAODEvent *aod =  dynamic_cast<AliAODEvent*> (InputEvent());
     // to be implemented
     printf("  AOD analysis not yet implemented!!!\n\n");
     return;
+  } 
+
+
+  // ***** Scaling
+  if (fUseScaling) {
+    Float_t temp_scalefactorV0M = MyGetScaleFactor(fCurrentRun,0);
+    Float_t temp_scalefactorSPD = MyGetScaleFactor(fCurrentRun,1);
+    Float_t temp_scalefactorTPC = MyGetScaleFactor(fCurrentRun,2);
+    v0Corr  = Short_t(v0Corr / temp_scalefactorV0M);
+    spdCorr = spdCorr / temp_scalefactorSPD;
+    nTracks = Int_t(nTracks / temp_scalefactorTPC);
+  }
+  // ***** Scaling for MC
+  if (fIsMCInput) {
+    Float_t temp_scalefactorV0M = MyGetScaleFactorMC(fCurrentRun);
+    v0Corr  = Short_t((multV0A+multV0C)  * temp_scalefactorV0M);
   }
 
   // ***** Centrality Selection
-  if(fHtempV0M)  fCentV0M = fHtempV0M->GetBinContent(fHtempV0M->FindBin((v0Corr)));
-  ///  else     printf("  Centrality by V0 not available!!!\n\n");
+  if(fHtempV0M) fCentV0M = fHtempV0M->GetBinContent(fHtempV0M->FindBin((v0Corr)));
   if(fHtempFMD) fCentFMD = fHtempFMD->GetBinContent(fHtempFMD->FindBin((multFMDA+multFMDC)));
-  //  else     printf("  Centrality by FMD not available!!!\n\n");
   if(fHtempTRK) fCentTRK = fHtempTRK->GetBinContent(fHtempTRK->FindBin(nTracks));
-  //  else     printf("  Centrality by TRK not available!!!\n\n");
   if(fHtempTKL) fCentTKL = fHtempTKL->GetBinContent(fHtempTKL->FindBin(nTracklets));
-  //  else     printf("  Centrality by TKL not available!!!\n\n");
   if(fHtempCL0) fCentCL0 = fHtempCL0->GetBinContent(fHtempCL0->FindBin(nClusters[0]));
-  //  else     printf("  Centrality by CL0 not available!!!\n\n");
   if(fHtempCL1) fCentCL1 = fHtempCL1->GetBinContent(fHtempCL1->FindBin(spdCorr));
-  ///  else     printf("  Centrality by CL1 not available!!!\n\n");
   
   if(fHtempV0MvsFMD) fCentV0MvsFMD = fHtempV0MvsFMD->GetBinContent(fHtempV0MvsFMD->FindBin((multV0A+multV0C)));
-  //  else     printf("  Centrality by V0 vs FMD not available!!!\n\n");
   if(fHtempTKLvsV0M) fCentTKLvsV0M = fHtempTKLvsV0M->GetBinContent(fHtempTKLvsV0M->FindBin(nTracklets));
-  //  else     printf("  Centrality by V0 vs TKL not available!!!\n\n");
-  if(fHtempZEMvsZDC) fCentZEMvsZDC = fHtempZEMvsZDC->GetBinContent(fHtempZEMvsZDC->FindBin((zem1Energy+zem2Energy)/1000.));
-  //  else     printf("  Centrality by ZEM vs ZDC not available!!!\n\n");
-
-  esdCent->SetCentralityV0M(fCentV0M);
-  esdCent->SetCentralityFMD(fCentFMD);
-  esdCent->SetCentralityTRK(fCentTRK);
-  esdCent->SetCentralityTKL(fCentTKL);
-  esdCent->SetCentralityCL0(fCentCL0);
-  esdCent->SetCentralityCL1(fCentCL1);
-  esdCent->SetCentralityV0MvsFMD(fCentV0MvsFMD);
-  esdCent->SetCentralityTKLvsV0M(fCentTKLvsV0M);
-  esdCent->SetCentralityZEMvsZDC(fCentZEMvsZDC);
-}
+  if(fHtempZEMvsZDC) fCentZEMvsZDC = fHtempZEMvsZDC->GetBinContent(fHtempZEMvsZDC->FindBin(zem1Energy+zem2Energy,zncEnergy+znaEnergy+zpcEnergy+zpaEnergy));
 
-//________________________________________________________________________
-void AliCentralitySelectionTask::ReadCentralityHistos() 
-{
-//  Read centrality histograms
-    TDirectory *owd = gDirectory;
-    fFile  = TFile::Open(fCentfilename);
-    owd->cd();
-    fHtempV0M  = (TH1F*) (fFile->Get("hmultV0_percentile"));
-    fHtempFMD  = (TH1F*) (fFile->Get("hmultFMD_percentile"));
-    fHtempTRK  = (TH1F*) (fFile->Get("hNtracks_percentile"));
-    fHtempTKL  = (TH1F*) (fFile->Get("hNtracklets_percentile"));
-    fHtempCL0  = (TH1F*) (fFile->Get("hNclusters0_percentile"));
-    fHtempCL1  = (TH1F*) (fFile->Get("hNclusters1_percentile"));
-}  
+  // ***** Cleaning
+  fQuality=0;
+  fZVCut=10;
+  fOutliersCut=6;
+  
+  // ***** vertex
+  if (TMath::Abs(zvtx)>fZVCut) fQuality += 1;   
 
-//________________________________________________________________________
-void AliCentralitySelectionTask::ReadCentralityHistos2() 
-{
-//  Read centrality histograms
-    TDirectory *owd = gDirectory;
-    fFile2  = TFile::Open(fCentfilename2);
-    owd->cd();
-    fHtempV0MvsFMD =  (TH1F*) (fFile2->Get("hmultV0vsmultFMD_all_percentile"));
-    fHtempTKLvsV0M  = (TH1F*) (fFile2->Get("hNtrackletsvsmultV0_all_percentile"));
-    fHtempZEMvsZDC  = (TH1F*) (fFile2->Get("hEzemvsEzdc_all_percentile"));
+  // ***** outliers
+  // **** V0 vs SPD
+  if (IsOutlierV0MSPD(spdCorr, v0Corr, int(fCentV0M))) fQuality  += 2;
+  // ***** V0 vs TPC
+  if (IsOutlierV0MTPC(nTracks, v0Corr, int(fCentV0M))) fQuality  += 4;
+  // ***** V0 vs ZDC
+  if (IsOutlierV0MZDC((zncEnergy+znaEnergy+zpcEnergy+zpaEnergy), v0Corr) && 
+      (zdcEnergyCal==kFALSE)) fQuality  += 8;
+
+  if (esdCent) {
+      esdCent->SetQuality(fQuality);
+      esdCent->SetCentralityV0M(fCentV0M);
+      esdCent->SetCentralityFMD(fCentFMD);
+      esdCent->SetCentralityTRK(fCentTRK);
+      esdCent->SetCentralityTKL(fCentTKL);
+      esdCent->SetCentralityCL0(fCentCL0);
+      esdCent->SetCentralityCL1(fCentCL1);
+      esdCent->SetCentralityV0MvsFMD(fCentV0MvsFMD);
+      esdCent->SetCentralityTKLvsV0M(fCentTKLvsV0M);
+      esdCent->SetCentralityZEMvsZDC(fCentZEMvsZDC);
+  }
+
+  fHOutQuality->Fill(fQuality);
+  fHOutVertex->Fill(zvtx);
+
+  fHOutMultV0MvsZDN->Fill(v0Corr,(zncEnergy+znaEnergy));
+  fHOutMultZEMvsZDN->Fill((zem1Energy+zem2Energy),(zncEnergy+znaEnergy));
+  fHOutMultV0MvsZDC->Fill(v0Corr,(zncEnergy+znaEnergy+zpcEnergy+zpaEnergy));
+  fHOutMultZEMvsZDC->Fill((zem1Energy+zem2Energy),(zncEnergy+znaEnergy+zpcEnergy+zpaEnergy));
+  fHOutMultV0MvsCL1->Fill(v0Corr,spdCorr);
+  fHOutMultV0MvsTRK->Fill(v0Corr,nTracks);
+  
+  if (fQuality==0) {  
+    fHOutCentV0M->Fill(fCentV0M);
+    fHOutCentFMD->Fill(fCentFMD);
+    fHOutCentTRK->Fill(fCentTRK);
+    fHOutCentTKL->Fill(fCentTKL);
+    fHOutCentCL0->Fill(fCentCL0);
+    fHOutCentCL1->Fill(fCentCL1);
+    fHOutCentV0MvsFMD->Fill(fCentV0MvsFMD);
+    fHOutCentTKLvsV0M->Fill(fCentTKLvsV0M);
+    fHOutCentZEMvsZDC->Fill(fCentZEMvsZDC);
+    fHOutMultV0M->Fill(v0Corr);
+    fHOutMultV0R->Fill(multV0A+multV0C);
+    fHOutMultFMD->Fill((multFMDA+multFMDC));
+    fHOutMultTRK->Fill(nTracks);
+    fHOutMultTKL->Fill(nTracklets);
+    fHOutMultCL0->Fill(nClusters[0]);
+    fHOutMultCL1->Fill(spdCorr);
+    fHOutMultTRKvsCL1->Fill(nTracks,spdCorr);
+  } else if (fQuality ==1) {
+    fHOutCentV0M_qual1->Fill(fCentV0M);
+    fHOutCentTRK_qual1->Fill(fCentTRK);
+    fHOutCentCL1_qual1->Fill(fCentCL1);
+  } else {
+    fHOutCentV0M_qual2->Fill(fCentV0M);
+    fHOutCentTRK_qual2->Fill(fCentTRK);
+    fHOutCentCL1_qual2->Fill(fCentCL1);
+  }
+
+  PostData(1, fOutputList); 
 }
 
 //________________________________________________________________________
-void AliCentralitySelectionTask::SetPercentileFile(TString filename) 
+void AliCentralitySelectionTask::ReadCentralityHistos(TString fCentfilename) 
 {
-// Set the percentile file name
-  fCentfilename = filename;
-  ReadCentralityHistos();
-}
+  //  Read centrality histograms
+  TDirectory *owd = gDirectory;
+  // Check if the file is present
+  TString path = gSystem->ExpandPathName(fCentfilename.Data());
+  if (gSystem->AccessPathName(path)) {
+     AliError(Form("File %s does not exist", path.Data()));
+     return;
+  }
+  fFile  = TFile::Open(fCentfilename);
+  owd->cd();
+  fHtempV0M  = (TH1F*) (fFile->Get("hmultV0_percentile"));
+  fHtempFMD  = (TH1F*) (fFile->Get("hmultFMD_percentile"));
+  fHtempTRK  = (TH1F*) (fFile->Get("hNtracks_percentile"));
+  fHtempTKL  = (TH1F*) (fFile->Get("hNtracklets_percentile"));
+  fHtempCL0  = (TH1F*) (fFile->Get("hNclusters0_percentile"));
+  fHtempCL1  = (TH1F*) (fFile->Get("hNclusters1_percentile"));
+
+  if (!fHtempV0M) AliWarning(Form("Calibration for V0M does not exist in %s", path.Data()));
+  if (!fHtempFMD) AliWarning(Form("Calibration for FMD does not exist in %s", path.Data()));
+  if (!fHtempTRK) AliWarning(Form("Calibration for TRK does not exist in %s", path.Data()));
+  if (!fHtempTKL) AliWarning(Form("Calibration for TKL does not exist in %s", path.Data()));
+  if (!fHtempCL0) AliWarning(Form("Calibration for CL0 does not exist in %s", path.Data()));
+  if (!fHtempCL1) AliWarning(Form("Calibration for CL1 does not exist in %s", path.Data()));
+  
+  owd->cd();
+}  
 
 //________________________________________________________________________
-void AliCentralitySelectionTask::SetPercentileFile2(TString filename
+void AliCentralitySelectionTask::ReadCentralityHistos2(TString fCentfilename2
 {
-// Set the percentile file name
-  fCentfilename2 = filename;
-  ReadCentralityHistos2();
+  //  Read centrality histograms
+  TDirectory *owd = gDirectory;
+  TString path = gSystem->ExpandPathName(fCentfilename2.Data());
+  if (gSystem->AccessPathName(path)) {
+     AliError(Form("File %s does not exist", path.Data()));
+     return;
+  }   
+  fFile2  = TFile::Open(fCentfilename2);
+  owd->cd();
+  fHtempV0MvsFMD =  (TH1F*) (fFile2->Get("hmultV0vsmultFMD_all_percentile"));
+  fHtempTKLvsV0M  = (TH1F*) (fFile2->Get("hNtrackletsvsmultV0_all_percentile"));
+  fHtempZEMvsZDC  = (TH2F*) (fFile2->Get("hEzemvsEzdc_all_percentile"));
+
+  if (!fHtempV0MvsFMD) AliWarning(Form("Calibration for V0MvsFMD does not exist in %s", path.Data()));
+  if (!fHtempTKLvsV0M) AliWarning(Form("Calibration for TKLvsV0M does not exist in %s", path.Data()));
+  if (!fHtempZEMvsZDC) AliWarning(Form("Calibration for ZEMvsZDC does not exist in %s", path.Data()));
+  
+  owd->cd();
 }
 
 //________________________________________________________________________
@@ -444,150 +732,387 @@ Int_t AliCentralitySelectionTask::SetupRun(AliESDEvent* esd)
   
   AliInfo(Form("Setup Centrality Selection for run %d\n",fCurrentRun));
 
-  fRunNo = fCurrentRun;
-
   // CHANGE HERE FOR RUN RANGES
-  if ( fRunNo == 137162 ) fRunNo = 137161;
-  else if ( fRunNo == 137365 ) fRunNo = 137366;
-  else if ( fRunNo > 137366 ) fRunNo = 137366;
+  if ( fCurrentRun <= 137165 ) fRunNo = 137161;
+  else fRunNo = 137366;
   // CHANGE HERE FOR RUN RANGES
+  
+  TString fileName(Form("%s/COMMON/CENTRALITY/data/AliCentralityBy1D_%d.root", AliAnalysisManager::GetOADBPath(), fRunNo));
+  TString fileName2(Form("%s/COMMON/CENTRALITY/data/AliCentralityByFunction_%d.root", AliAnalysisManager::GetOADBPath(), fRunNo));
+  
+  AliInfo(Form("Centrality Selection for run %d is initialized with %s", fCurrentRun, fileName.Data()));
+  ReadCentralityHistos(fileName.Data());
+  ReadCentralityHistos2(fileName2.Data());
+  if (!fFile && !fFile2) {
+     AliFatal(Form("Run %d not known to centrality selection!", fCurrentRun));       
+     return -1;
+  }   
+  return 0;
+}
 
-  TString runName(Form("%d", fRunNo));
-  TString fileName("");
-  Bool_t isRunKnown = kFALSE;
+//________________________________________________________________________
+Bool_t AliCentralitySelectionTask::IsOutlierV0MSPD(Float_t spd, Float_t v0, Int_t cent)
+{
+  Float_t val= -0.143789 + 0.288874 * v0;
+  Float_t SPDsigma[100]={231.483, 189.446, 183.359, 179.923, 174.229, 170.309, 165.021, 
+                        160.84, 159.33, 154.453, 151.644, 148.337, 145.215, 142.353, 
+                        139.351, 136, 133.838, 129.885, 127.36, 125.032, 122.21, 120.3, 
+                        117.766, 114.77, 113.1, 110.268, 107.463, 105.293, 102.845, 
+                        100.835, 98.9632, 97.3287, 93.6887, 92.1066, 89.3224, 87.8382, 
+                        86.04, 83.6431, 81.9655, 80.0491, 77.8208, 76.4716, 74.2165, 
+                        72.2752, 70.4875, 68.9414, 66.8622, 65.022, 63.5134, 61.8228, 
+                        59.7166, 58.5008, 56.2789, 54.9615, 53.386, 51.2165, 49.4842, 
+                        48.259, 47.1129, 45.3115, 43.8486, 41.9207, 40.5754, 39.3872, 
+                        38.1897, 36.5401, 35.1283, 33.9702, 32.6429, 31.3612, 29.5876, 
+                        28.9319, 27.564, 26.0443, 25.2836, 23.9753, 22.8936, 21.5665, 
+                        20.7048, 19.8016, 18.7095, 18.1144, 17.2095, 16.602, 16.3233, 
+                        15.7185, 15.3006, 14.7432, 14.4174, 14.0805, 13.7638, 13.7638, 
+                        13.7638, 13.7638, 13.7638, 13.7638, 13.7638, 13.7638, 13.7638, 18.0803};
 
-  // Check if run is in fileList
-  // if not, take the last name in the list
-  for ( Int_t idx=0 ; idx < fFileList->GetEntries(); ++idx ) {
+  if ( TMath::Abs(spd-val) > fOutliersCut*SPDsigma[cent] ) 
+    return kTRUE;
+  else 
+    return kFALSE;
+}
 
-    TString str((dynamic_cast<TObjString*>(fFileList->At(idx)))->GetString());
-    if (str.Contains(runName)) {
-      fileName += str;
-      isRunKnown = kTRUE;
-      break;
-    }
-  }
+//________________________________________________________________________
+Bool_t AliCentralitySelectionTask::IsOutlierV0MTPC(Int_t tracks, Float_t v0, Int_t cent)
+{
+  Float_t val = -0.540691 + 0.128358 * v0;
+  Float_t TPCsigma[100]={106.439, 89.2834, 86.7568, 85.3641, 83.379, 81.6093, 79.3189, 
+                        78.0616, 77.2167, 75.0021, 73.9957, 72.0926, 71.0442, 69.8395, 
+                        68.1169, 66.6676, 66.0038, 64.2284, 63.3845, 61.7439, 60.642, 
+                        59.5383, 58.3696, 57.0227, 56.0619, 54.7108, 53.8382, 52.3398, 
+                        51.5297, 49.9488, 49.1126, 48.208, 46.8566, 45.7724, 44.7829, 
+                        43.8726, 42.7499, 41.9307, 40.6874, 39.9619, 38.5534, 38.0884, 
+                        36.6141, 35.7482, 34.8804, 34.1769, 33.1278, 32.3435, 31.4783, 
+                        30.2587, 29.4741, 28.8575, 27.9298, 26.9752, 26.1675, 25.1234, 
+                        24.4702, 23.6843, 22.9764, 21.8579, 21.2924, 20.3241, 19.8296, 
+                        19.2465, 18.4474, 17.7216, 16.8956, 16.342, 15.626, 15.0329, 
+                        14.3911, 13.9301, 13.254, 12.6745, 12.2436, 11.7776, 11.1795, 
+                        10.673, 10.27, 9.95646, 9.50939, 9.26162, 8.95315, 8.73439, 
+                        8.67375, 8.43029, 8.34818, 8.33484, 8.40709, 8.3974, 8.32814, 
+                        8.32814, 8.32814, 8.32814, 8.32814, 8.32814, 8.32814, 8.32814, 8.32814, 12.351};
 
-  if (!isRunKnown) {
-    if (fFileList->Last()) {
-      fileName += (dynamic_cast<TObjString*>(fFileList->Last()))->GetString();
-      AliError(Form("Run %d not known to centrality selection!", fCurrentRun));
-    }
+  if ( TMath::Abs(tracks-val) > fOutliersCut*TPCsigma[cent] ) 
+    return kTRUE;
+  else 
+    return kFALSE;
+}
+
+//________________________________________________________________________
+Bool_t AliCentralitySelectionTask::IsOutlierV0MZDC(Float_t zdc, Float_t v0)
+{
+  Float_t val1 = 6350. - 0.26 * v0;
+  Float_t val2 = 5580.;
+  if ((zdc >  val1) || (zdc > val2)) 
+    return kTRUE;
+  else 
+    return kFALSE;
+}
+
+//________________________________________________________________________  
+Float_t AliCentralitySelectionTask::MyGetScaleFactor(Int_t runnumber, Int_t flag) 
+{
+  if (! (runnumber >= fLowRunN && runnumber <=fHighRunN)) {
+    cout << "MyGetScaleFactor error in run number range " << runnumber << endl;
+    return 0.0;
   }
 
-  if (fileName.Contains(".root")) {
-    AliInfo(Form("Centrality Selection for run %d is initialized with %s", fCurrentRun, fileName.Data()));
-    SetPercentileFile(fileName.Data());
-    return 0;
+  Float_t scalefactor=0.0;
+  if (flag==0)
+    scalefactor = V0MScaleFactor[runnumber - fLowRunN]; // subtracting reference offset index
+  else if (flag==1)
+    scalefactor = SPDScaleFactor[runnumber - fLowRunN]; // subtracting reference offset index
+  else if (flag==2)
+    scalefactor = TPCScaleFactor[runnumber - fLowRunN]; // subtracting reference offset index
+
+  return scalefactor;
+
+}
+
+//________________________________________________________________________  
+Float_t AliCentralitySelectionTask::MyGetScaleFactorMC(Int_t runnumber) 
+{
+  if (! (runnumber >= fLowRunN && runnumber <=fHighRunN)) {
+    cout << "MyGetScaleFactor error in run number range " << runnumber << endl;
+    return 0.0;
   }
-  
-  return -1;
+
+  Float_t scalefactor= V0MScaleFactorMC[runnumber - fLowRunN]; // subtracting reference offset index
+  return scalefactor;
+
 }
-//________________________________________________________________________
-Float_t AliCentralitySelectionTask::GetCorrV0(const AliESDEvent* esd, float &v0CorrResc) 
+
+//________________________________________________________________________  
+void AliCentralitySelectionTask::MyInitScaleFactor () 
 {
-  // correct V0 non-linearity, prepare a version rescaled to SPD2 corr
-  Double_t *par0;
-  Double_t *par1;
-  Double_t *par2;
+  for (int i=0; i<(fHighRunN-fLowRunN); i++) V0MScaleFactor[i] = 0.0;
+  for (int i=0; i<(fHighRunN-fLowRunN); i++) SPDScaleFactor[i] = 0.0;
+  for (int i=0; i<(fHighRunN-fLowRunN); i++) TPCScaleFactor[i] = 0.0;
   
-  Double_t par0_137161[64] = { 6.71e-02 , 6.86e-02 , 7.06e-02 , 6.32e-02 , 
-                       5.91e-02 , 6.07e-02 , 5.78e-02 , 5.73e-02 , 5.91e-02 , 6.22e-02 , 
-                       5.90e-02 , 6.11e-02 , 5.55e-02 , 5.29e-02 , 5.19e-02 , 5.56e-02 , 
-                       6.25e-02 , 7.03e-02 , 5.64e-02 , 5.81e-02 , 4.57e-02 , 5.30e-02 , 
-                       5.13e-02 , 6.43e-02 , 6.27e-02 , 6.48e-02 , 6.07e-02 , 1.01e-01 , 
-                       6.68e-02 , 7.16e-02 , 6.36e-02 , 5.95e-02 , 2.52e-02 , 2.82e-02 , 
-                       2.56e-02 , 2.86e-02 , 2.82e-02 , 2.10e-02 , 2.13e-02 , 2.32e-02 , 
-                       2.75e-02 , 4.34e-02 , 3.78e-02 , 4.52e-02 , 4.11e-02 , 3.89e-02 , 
-                       4.10e-02 , 3.73e-02 , 4.51e-02 , 5.07e-02 , 5.42e-02 , 4.74e-02 , 
-                       4.33e-02 , 4.44e-02 , 4.64e-02 , 3.01e-02 , 6.38e-02 , 5.26e-02 , 
-                       4.99e-02 , 5.26e-02 , 5.47e-02 , 3.84e-02 , 5.00e-02 , 5.20e-02 };
-  Double_t par1_137161[64] = { -6.68e-05 , -7.78e-05 , -6.88e-05 , -5.92e-05 , 
-                       -2.43e-05 , -3.54e-05 , -2.91e-05 , -1.99e-05 , -1.40e-05 , -4.01e-05 , 
-                       -2.29e-05 , -3.68e-05 , -2.53e-05 , -2.44e-06 , -9.22e-06 , -1.51e-05 , 
-                       -2.80e-05 , -2.34e-05 , -1.72e-05 , -1.81e-05 , -1.29e-05 , -2.65e-05 , 
-                       -1.61e-05 , -2.86e-05 , -1.74e-05 , -4.23e-05 , -3.41e-05 , -1.05e-04 , 
-                       -2.76e-05 , -4.71e-05 , -3.06e-05 , -2.32e-05 , -1.55e-06 , 2.15e-05 , 
-                       1.40e-05 , 2.16e-05 , 1.21e-05 , 3.05e-06 , 1.67e-05 , -3.84e-06 , 
-                       3.09e-06 , 1.50e-05 , 3.47e-06 , 4.87e-06 , -3.71e-07 , -1.75e-06 , 
-                       -1.80e-06 , 9.99e-06 , -6.46e-06 , -4.91e-06 , 1.33e-05 , -2.52e-07 , 
-                       -3.85e-06 , 4.94e-06 , -2.48e-07 , -1.20e-05 , 2.07e-06 , 6.12e-06 , 
-                       -1.18e-06 , 4.54e-06 , -1.54e-05 , -1.25e-05 , 1.46e-06 , -6.67e-06 };
-  Double_t par2_137161[64] = { 1.29e-08 , 1.51e-08 , 1.43e-08 , 1.11e-08 , 
-                       5.04e-09 , 6.99e-09 , 5.58e-09 , 4.15e-09 , 4.00e-09 , 8.22e-09 , 
-                       4.97e-09 , 7.66e-09 , 4.91e-09 , 1.10e-09 , 2.64e-09 , 3.64e-09 , 
-                       5.76e-09 , 5.46e-09 , 3.38e-09 , 3.47e-09 , 2.43e-09 , 4.13e-09 , 
-                       2.80e-09 , 5.80e-09 , 3.86e-09 , 7.46e-09 , 5.98e-09 , 2.58e-08 , 
-                       5.50e-09 , 8.72e-09 , 5.23e-09 , 4.37e-09 , 2.33e-09 , -6.01e-10 , 
-                       3.99e-11 , -2.02e-10 , 7.67e-10 , 2.03e-09 , 1.17e-10 , 2.56e-09 , 
-                       1.16e-09 , -4.75e-10 , 1.28e-09 , 1.23e-09 , 1.62e-09 , 1.61e-09 , 
-                       1.93e-09 , 2.97e-10 , 2.21e-09 , 2.16e-09 , 5.22e-10 , 1.03e-09 , 
-                       1.56e-09 , 5.00e-10 , 1.01e-09 , 2.93e-09 , 1.05e-09 , 9.96e-11 , 
-                       1.21e-09 , 7.45e-10 , 3.07e-09 , 2.31e-09 , 6.70e-10 , 1.89e-09 };
-
-  Double_t par0_137366[64] = { 7.12e-02 , 7.34e-02 , 7.39e-02 , 6.54e-02 , 6.11e-02 , 6.31e-02 , 6.15e-02 , 
-                              6.00e-02 , 6.10e-02 , 6.49e-02 , 6.17e-02 , 6.33e-02 , 6.00e-02 , 5.48e-02 , 
-                              5.44e-02 , 5.81e-02 , 6.49e-02 , 7.07e-02 , 5.91e-02 , 6.18e-02 , 4.82e-02 , 
-                              5.67e-02 , 5.36e-02 , 6.60e-02 , 6.37e-02 , 6.78e-02 , 6.31e-02 , 1.04e-01 , 
-                              6.91e-02 , 7.32e-02 , 6.61e-02 , 6.16e-02 , 2.64e-02 , 2.81e-02 , 2.64e-02 , 
-                              2.85e-02 , 2.87e-02 , 2.18e-02 , 2.19e-02 , 2.43e-02 , 2.81e-02 , 4.37e-02 , 
-                              3.90e-02 , 4.66e-02 , 4.24e-02 , 4.09e-02 , 4.21e-02 , 3.88e-02 , 4.83e-02 , 
-                              5.23e-02 , 5.44e-02 , 4.85e-02 , 4.42e-02 , 4.58e-02 , 4.74e-02 , 3.14e-02 , 
-                              6.31e-02 , 5.30e-02 , 5.01e-02 , 5.33e-02 , 5.70e-02 , 3.95e-02 , 4.98e-02 , 5.31e-02 };
-  Double_t par1_137366[64] = { -6.99e-05 , -6.99e-05 , -6.94e-05 , -6.55e-05 , -3.55e-05 , -4.50e-05 , 
-                              -3.10e-05 , -2.81e-05 , -2.29e-05 , -3.89e-05 , -2.53e-05 , -4.25e-05 ,
-                              -1.87e-05 , -2.01e-05 , -1.53e-05 , -2.14e-05 , -2.86e-05 , -4.70e-05 ,
-                              -2.23e-05 , -3.30e-05 ,-9.74e-06 , -2.62e-05 , -1.76e-05 , -2.38e-05 , 
-                              -2.40e-05 , -3.43e-05 , -2.75e-05 , -6.86e-05 ,-2.35e-05 , -4.45e-05 , 
-                              -2.51e-05 , -2.20e-05 , -1.25e-16 , -2.04e-17 , -2.06e-17 , -3.74e-19 ,
-                              -1.18e-18 , -2.02e-15 , -3.78e-06 , -1.26e-06 , -2.71e-06 , -6.23e-17 , 
-                              -7.39e-08 , -1.76e-16 , -8.98e-06 , -4.10e-18 , -1.34e-05 , -1.06e-16 , 
-                              -3.34e-06 , -1.04e-05 , -5.28e-06 , -7.34e-06 , -1.05e-05 , -7.68e-06 ,
-                              -1.78e-05 , -1.19e-05 , -1.78e-05 , -1.34e-06 , -9.23e-06 , -3.34e-06 ,
-                              -8.02e-06 , -1.39e-05 , -1.38e-05 , -1.40e-05 };
-  Double_t par2_137366[64] = { 1.41e-08 , 1.47e-08 , 1.48e-08 , 1.24e-08 , 6.82e-09 , 8.73e-09 , 6.26e-09 , 
-                              5.53e-09 , 5.40e-09 , 7.93e-09 , 5.49e-09 , 8.77e-09 , 4.21e-09 , 3.93e-09 , 
-                              3.60e-09 , 4.67e-09 , 5.59e-09 , 8.81e-09 , 3.89e-09 , 6.19e-09 , 1.97e-09 , 
-                              4.38e-09 , 3.26e-09 , 5.00e-09 , 4.58e-09 , 6.39e-09 , 5.03e-09 , 1.30e-08 , 
-                              4.95e-09 , 8.26e-09 , 4.57e-09 , 4.10e-09 , 2.35e-09 , 2.30e-09 , 2.15e-09 , 
-                              2.27e-09 , 2.17e-09 , 2.27e-09 , 2.97e-09 , 2.25e-09 , 1.69e-09 , 1.44e-09 , 
-                              1.66e-09 , 1.75e-09 , 2.88e-09 , 1.82e-09 , 3.64e-09 , 1.80e-09 , 1.71e-09 , 
-                              2.66e-09 , 3.01e-09 , 1.95e-09 , 2.64e-09 , 2.42e-09 , 3.68e-09 , 2.66e-09 , 
-                              3.92e-09 , 1.18e-09 , 2.26e-09 , 1.57e-09 , 2.02e-09 , 2.71e-09 , 2.99e-09 , 3.04e-09 }; 
+  // scale factors determined from <V0 charge> on a run-by-run basis
+   V0MScaleFactor[310] = 0.956841;
+   V0MScaleFactor[311] = 0.958274;
+   V0MScaleFactor[514] = 1.0046;
+   V0MScaleFactor[515] = 0.983535;
+   V0MScaleFactor[579] = 0.988185;
+   V0MScaleFactor[580] = 0.983351;
+   V0MScaleFactor[581] = 0.989013;
+   V0MScaleFactor[583] = 0.990056;
+   V0MScaleFactor[588] = 0.974438;
+   V0MScaleFactor[589] = 0.981572;
+   V0MScaleFactor[590] = 0.989316;
+   V0MScaleFactor[592] = 0.98345;
+   V0MScaleFactor[688] = 0.993647;
+   V0MScaleFactor[690] = 0.994758;
+   V0MScaleFactor[693] = 0.989569;
+   V0MScaleFactor[698] = 0.993119;
+   V0MScaleFactor[744] = 0.989583;
+   V0MScaleFactor[757] = 0.990377;
+   V0MScaleFactor[787] = 0.990176;
+   V0MScaleFactor[788] = 0.98723;
+   V0MScaleFactor[834] = 1.00403;
+   V0MScaleFactor[835] = 0.994376;
+   V0MScaleFactor[840] = 0.99759;
+   V0MScaleFactor[842] = 1.01031;
+   V0MScaleFactor[900] = 0.996216;
+   V0MScaleFactor[901] = 0.994205;
+   V0MScaleFactor[992] = 0.998479;
+   V0MScaleFactor[997] = 1.00078;
+   V0MScaleFactor[1299] = 1.00515;
+   V0MScaleFactor[1303] = 1.00094;
+   V0MScaleFactor[1339] = 0.986596;
+   V0MScaleFactor[1346] = 0.972226;
+   V0MScaleFactor[1349] = 0.960358;
+   V0MScaleFactor[1350] = 0.970023;
+   V0MScaleFactor[1374] = 1.00575;
+   V0MScaleFactor[1545] = 1.00471;
+   V0MScaleFactor[1587] = 1.00611;
+   V0MScaleFactor[1588] = 1.00976;
+   V0MScaleFactor[1618] = 1.00771;
+   V0MScaleFactor[1682] = 1.01622;
+   V0MScaleFactor[1727] = 1.01305;
+   V0MScaleFactor[1728] = 1.00388;
+   V0MScaleFactor[1731] = 1.00673;
+   V0MScaleFactor[1732] = 1.00916;
+   V0MScaleFactor[1770] = 1.0092;
+   V0MScaleFactor[1773] = 1.00728;
+   V0MScaleFactor[1786] = 1.01655;
+   V0MScaleFactor[1787] = 1.00672;
+   V0MScaleFactor[1801] = 0.983339;
+   V0MScaleFactor[1802] = 1.00754;
+   V0MScaleFactor[1811] = 1.00608;
+   V0MScaleFactor[1815] = 1.01227;
+   V0MScaleFactor[1879] = 0.99907;
+   V0MScaleFactor[1881] = 0.995696;
+   V0MScaleFactor[2186] = 1.00559;
+   V0MScaleFactor[2187] = 1.00631;
+   V0MScaleFactor[2254] = 1.01512;
+   V0MScaleFactor[2256] = 0.998727;
+   V0MScaleFactor[2321] = 1.00701;
+   SPDScaleFactor[310] = 1.00211;
+   SPDScaleFactor[311] = 1.00067;
+   SPDScaleFactor[514] = 1.02268;
+   SPDScaleFactor[515] = 0.994902;
+   SPDScaleFactor[579] = 1.00215;
+   SPDScaleFactor[580] = 0.993421;
+   SPDScaleFactor[581] = 1.00129;
+   SPDScaleFactor[583] = 1.00242;
+   SPDScaleFactor[588] = 0.984762;
+   SPDScaleFactor[589] = 0.994355;
+   SPDScaleFactor[590] = 1.00073;
+   SPDScaleFactor[592] = 0.995889;
+   SPDScaleFactor[688] = 0.994532;
+   SPDScaleFactor[690] = 0.998307;
+   SPDScaleFactor[693] = 0.994052;
+   SPDScaleFactor[698] = 0.993224;
+   SPDScaleFactor[744] = 0.993279;
+   SPDScaleFactor[757] = 0.992494;
+   SPDScaleFactor[787] = 0.992678;
+   SPDScaleFactor[788] = 0.996563;
+   SPDScaleFactor[834] = 1.01116;
+   SPDScaleFactor[835] = 0.993108;
+   SPDScaleFactor[840] = 0.997574;
+   SPDScaleFactor[842] = 1.01829;
+   SPDScaleFactor[900] = 0.999438;
+   SPDScaleFactor[901] = 0.995849;
+   SPDScaleFactor[992] = 0.999227;
+   SPDScaleFactor[997] = 1.00575;
+   SPDScaleFactor[1299] = 0.99877;
+   SPDScaleFactor[1303] = 0.999682;
+   SPDScaleFactor[1339] = 0.978198;
+   SPDScaleFactor[1346] = 0.964178;
+   SPDScaleFactor[1349] = 0.959439;
+   SPDScaleFactor[1350] = 0.956945;
+   SPDScaleFactor[1374] = 0.994434;
+   SPDScaleFactor[1545] = 1.0016;
+   SPDScaleFactor[1587] = 1.00153;
+   SPDScaleFactor[1588] = 1.00698;
+   SPDScaleFactor[1618] = 1.00554;
+   SPDScaleFactor[1682] = 1.0123;
+   SPDScaleFactor[1727] = 1.011;
+   SPDScaleFactor[1728] = 1.00143;
+   SPDScaleFactor[1731] = 1.00486;
+   SPDScaleFactor[1732] = 1.00646;
+   SPDScaleFactor[1770] = 1.00515;
+   SPDScaleFactor[1773] = 1.00485;
+   SPDScaleFactor[1786] = 1.01215;
+   SPDScaleFactor[1787] = 1.00419;
+   SPDScaleFactor[1801] = 0.983327;
+   SPDScaleFactor[1802] = 1.00529;
+   SPDScaleFactor[1811] = 1.00367;
+   SPDScaleFactor[1815] = 1.01045;
+   SPDScaleFactor[1879] = 0.996374;
+   SPDScaleFactor[1881] = 0.988827;
+   SPDScaleFactor[2186] = 1.00354;
+   SPDScaleFactor[2187] = 1.00397;
+   SPDScaleFactor[2254] = 1.01138;
+   SPDScaleFactor[2256] = 0.996641;
+   SPDScaleFactor[2321] = 1.00357;
+   TPCScaleFactor[310] = 1.00434;
+   TPCScaleFactor[311] = 1.0056;
+   TPCScaleFactor[514] = 1.02185;
+   TPCScaleFactor[515] = 1.0036;
+   TPCScaleFactor[579] = 1.00607;
+   TPCScaleFactor[580] = 1.00183;
+   TPCScaleFactor[581] = 1.00693;
+   TPCScaleFactor[583] = 1.00746;
+   TPCScaleFactor[588] = 0.990524;
+   TPCScaleFactor[589] = 0.998582;
+   TPCScaleFactor[590] = 1.00618;
+   TPCScaleFactor[592] = 1.00088;
+   TPCScaleFactor[688] = 1.00598;
+   TPCScaleFactor[690] = 1.00658;
+   TPCScaleFactor[693] = 1.00144;
+   TPCScaleFactor[698] = 1.00279;
+   TPCScaleFactor[744] = 1.00122;
+   TPCScaleFactor[757] = 1.002;
+   TPCScaleFactor[787] = 0.997818;
+   TPCScaleFactor[788] = 0.994583;
+   TPCScaleFactor[834] = 1.01508;
+   TPCScaleFactor[835] = 1.00218;
+   TPCScaleFactor[840] = 1.00569;
+   TPCScaleFactor[842] = 1.01789;
+   TPCScaleFactor[900] = 1.00739;
+   TPCScaleFactor[901] = 1.00462;
+   TPCScaleFactor[992] = 1.00502;
+   TPCScaleFactor[997] = 1.00943;
+   TPCScaleFactor[1299] = 1.00438;
+   TPCScaleFactor[1303] = 0.996701;
+   TPCScaleFactor[1339] = 0.978641;
+   TPCScaleFactor[1346] = 0.968906;
+   TPCScaleFactor[1349] = 0.954311;
+   TPCScaleFactor[1350] = 0.958764;
+   TPCScaleFactor[1374] = 0.997899;
+   TPCScaleFactor[1545] = 0.992;
+   TPCScaleFactor[1587] = 0.992635;
+   TPCScaleFactor[1588] = 1.00207;
+   TPCScaleFactor[1618] = 1.00126;
+   TPCScaleFactor[1682] = 1.00324;
+   TPCScaleFactor[1727] = 1.00042;
+   TPCScaleFactor[1728] = 0.978881;
+   TPCScaleFactor[1731] = 0.999818;
+   TPCScaleFactor[1732] = 1.00109;
+   TPCScaleFactor[1770] = 0.99935;
+   TPCScaleFactor[1773] = 0.998531;
+   TPCScaleFactor[1786] = 0.999125;
+   TPCScaleFactor[1787] = 0.998479;
+   TPCScaleFactor[1801] = 0.9775;
+   TPCScaleFactor[1802] = 0.999095;
+   TPCScaleFactor[1811] = 0.998197;
+   TPCScaleFactor[1815] = 1.00413;
+   TPCScaleFactor[1879] = 0.990916;
+   TPCScaleFactor[1881] = 0.987241;
+   TPCScaleFactor[2186] = 1.00048;
+   TPCScaleFactor[2187] = 1.00057;
+   TPCScaleFactor[2254] = 1.00588;
+   TPCScaleFactor[2256] = 0.991503;
+   TPCScaleFactor[2321] = 1.00057;
+
+
+  // set all missing values to the value of the run before it ....
+  for (int i=0; i<(fHighRunN-fLowRunN); i++) {    
+    if (V0MScaleFactor[i] == 0.0) {     
+      if (i==0) {
+       V0MScaleFactor[i] = 1.00;
+      } else {
+       // search for last run number with non-zero value
+       for (int j=i-1;j>=0;j--) {
+         if (V0MScaleFactor[j] != 0.0) {
+           V0MScaleFactor[i] = V0MScaleFactor[j];
+           break;
+         }
+       }
+      }
+    }
+  } // end loop over checking all run-numbers 
+
+  for (int i=0; i<(fHighRunN-fLowRunN); i++) {    
+    if (SPDScaleFactor[i] == 0.0) {     
+      if (i==0) {
+       SPDScaleFactor[i] = 1.00;
+      } else {
+       for (int j=i-1;j>=0;j--) {
+         if (SPDScaleFactor[j] != 0.0) {
+           SPDScaleFactor[i] = SPDScaleFactor[j];
+           break;
+         }
+       }
+      }
+    }
+  } 
+
+  for (int i=0; i<(fHighRunN-fLowRunN); i++) {    
+    if (TPCScaleFactor[i] == 0.0) {     
+      if (i==0) {
+       TPCScaleFactor[i] = 1.00;
+      } else {
+       for (int j=i-1;j>=0;j--) {
+         if (TPCScaleFactor[j] != 0.0) {
+           TPCScaleFactor[i] = TPCScaleFactor[j];
+           break;
+         }
+       }
+      }
+    }
+  } 
   
+
+  //    for (int i=0; i<(fHighRunN-fLowRunN); i++) cout << "Scale Factor = " << V0MScaleFactor[i] << " for Run " << i+fLowRunN << endl;
   
-  if (esd->GetRunNumber() <= 137165) {
-    par0=par0_137161;
-    par1=par1_137161;
-    par2=par2_137161;
-  }  else  {
-    par0=par0_137366;
-    par1=par1_137366;
-    par2=par2_137366;
- }
-  //
-  Float_t multCorr = 0;
-  Float_t multCorr2 = 0;
-  Float_t multChCorr[64];
-  AliESDVZERO* esdV0 = esd->GetVZEROData();
-  for(Int_t i = 0; i < 64; ++i) {
-    Double_t b = (esdV0->GetMultiplicity(i)*par1[i]-par0[i]);
-    Double_t s = (b*b-4.*par2[i]*esdV0->GetMultiplicity(i)*esdV0->GetMultiplicity(i));
-    Double_t n = (s<0) ? -b : (-b + TMath::Sqrt(s));
-    multChCorr[i] = 2.*esdV0->GetMultiplicity(i)/n*par0[i];
-    multCorr += multChCorr[i];
-    multCorr2 += (multChCorr[i]/par0[i]/64.);
-  }
-  v0CorrResc =  multCorr2;
-  return multCorr;
+  return;
+
 }
 
-//____________________________________________________________________
-Float_t AliCentralitySelectionTask::GetCorrSPD2(Float_t spd2raw,Float_t zv)
+
+//________________________________________________________________________  
+void AliCentralitySelectionTask::MyInitScaleFactorMC() 
 {
-  // renormalize N spd2 clusters at given Zv to acceptance at Zv=0
-  const double pars[] = {8.10030e-01,-2.80364e-03,-7.19504e-04};
-  zv -= pars[0];
-  float corr = 1 + zv*(pars[1] + zv*pars[2]);
-  return corr>0 ? spd2raw/corr : -1;
+  for (int i=0; i<(fHighRunN-fLowRunN); i++) V0MScaleFactorMC[i] = 0.0;
+  // scale factors determined from <V0 charge> on a run-by-run basis
+  V0MScaleFactor[0] = 0.75108;
+  // set all missing values to the value of the run before it ....
+  for (int i=0; i<(fHighRunN-fLowRunN); i++) {    
+    if (V0MScaleFactorMC[i] == 0.0) {     
+      if (i==0) {
+       V0MScaleFactorMC[i] = 1.00;
+      } else {
+       // search for last run number with non-zero value
+       for (int j=i-1;j>=0;j--) {
+         if (V0MScaleFactorMC[j] != 0.0) {
+           V0MScaleFactorMC[i] = V0MScaleFactorMC[j];
+           break;
+         }
+       }
+      }
+    }
+  } // end loop over checking all run-numbers 
+
+
+  //    for (int i=0; i<(fHighRunN-fLowRunN); i++) cout << "Scale Factor = " << V0MScaleFactorMC[i] << " for Run " << i+fLowRunN << endl;
+  
+  return;
+
 }
+