]> git.uio.no Git - u/mrichter/AliRoot.git/blobdiff - ANALYSIS/AliCentralitySelectionTask.cxx
Changes in control histograms.
[u/mrichter/AliRoot.git] / ANALYSIS / AliCentralitySelectionTask.cxx
index 818a3e03469863302fea3a70ceffe6d449d4422d..4c15b1ff7c3cd80af7fd1c892820f2b452ac258b 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,7 @@
 #include "AliESDZDC.h"
 #include "AliESDFMD.h"
 #include "AliESDVZERO.h"
-#include "AliESDCentrality.h"
+#include "AliCentrality.h"
 #include "AliESDtrackCuts.h"
 #include "AliMultiplicity.h"
 #include "AliAODHandler.h"
@@ -61,7 +65,7 @@
 #include "AliPhysicsSelectionTask.h"
 #include "AliPhysicsSelection.h"
 #include "AliBackgroundSelection.h"
-#include "AliCentralitySelectionTask.h"
+#include "AliESDUtils.h"
 
 ClassImp(AliCentralitySelectionTask)
 
@@ -76,7 +80,13 @@ AliAnalysisTaskSE(),
   fFile2(0),
   fCurrentRun(-1),
   fRunNo(-1),
+  fLowRunN(0),
+  fHighRunN(0),
+  fUseScaling(0),
   fTrackCuts(0),
+  fZVCut(10),
+  fOutliersCut(5),
+  fQuality(0),
   fCentV0M(0),
   fCentFMD(0),
   fCentTRK(0),
@@ -104,11 +114,42 @@ AliAnalysisTaskSE(),
   fHOutCentCL1     (0),
   fHOutCentV0MvsFMD(0),
   fHOutCentTKLvsV0M(0),
-  fHOutCentZEMvsZDC(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
   AliInfo("Centrality Selection enabled.");
-  DefineOutput(1, TList::Class());
+  fLowRunN =136851;
+  fHighRunN=139517;
+
+  for (Int_t i=0; i<(fHighRunN-fLowRunN); i++) {
+    V0MScaleFactor[i]=0.0;
+    SPDScaleFactor[i]=0.0;
+    TPCScaleFactor[i]=0.0;
+    V0MScaleFactorMC[i]=0.0;
+  }
+  fUseScaling=kTRUE;
 }   
 
 //________________________________________________________________________
@@ -121,7 +162,13 @@ AliCentralitySelectionTask::AliCentralitySelectionTask(const char *name):
   fFile2(0),
   fCurrentRun(-1),
   fRunNo(-1),
+  fLowRunN(0),
+  fHighRunN(0),
+  fUseScaling(0),
   fTrackCuts(0),
+  fZVCut(10),
+  fOutliersCut(5),
+  fQuality(0),
   fCentV0M(0),
   fCentFMD(0),
   fCentTRK(0),
@@ -149,11 +196,43 @@ AliCentralitySelectionTask::AliCentralitySelectionTask(const char *name):
   fHOutCentCL1     (0),
   fHOutCentV0MvsFMD(0),
   fHOutCentTKLvsV0M(0),
-  fHOutCentZEMvsZDC(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
+    // Default constructor
+    fLowRunN =136851;
+    fHighRunN=139517;
+
   AliInfo("Centrality Selection enabled.");
   DefineOutput(1, TList::Class());
+  for (Int_t i=0; i<(fHighRunN-fLowRunN); i++) {
+    V0MScaleFactor[i]=0.0;
+    SPDScaleFactor[i]=0.0;
+    TPCScaleFactor[i]=0.0;
+    V0MScaleFactorMC[i]=0.0;
+  }
+  fUseScaling=kTRUE;
 }
 
 //________________________________________________________________________
@@ -176,7 +255,13 @@ AliCentralitySelectionTask::AliCentralitySelectionTask(const AliCentralitySelect
   fFile2(ana.fFile2),
   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),
@@ -204,7 +289,29 @@ AliCentralitySelectionTask::AliCentralitySelectionTask(const AliCentralitySelect
   fHOutCentCL1     (ana.fHOutCentCL1     ),
   fHOutCentV0MvsFMD(ana.fHOutCentV0MvsFMD),
   fHOutCentTKLvsV0M(ana.fHOutCentTKLvsV0M),
-  fHOutCentZEMvsZDC(ana.fHOutCentZEMvsZDC)
+  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  
 }
@@ -214,6 +321,7 @@ AliCentralitySelectionTask::~AliCentralitySelectionTask()
 {
   // Destructor  
   if (fOutputList && !AliAnalysisManager::GetAnalysisManager()->IsProofMode()) delete fOutputList;
+  if (fTrackCuts) delete fTrackCuts;
 }  
 
 //________________________________________________________________________
@@ -225,15 +333,41 @@ void AliCentralitySelectionTask::UserCreateOutputObjects()
 
   fOutputList = new TList();
   fOutputList->SetOwner();
-  fHOutCentV0M     = new TH1F("fHOutCentV0M","fHOutCentV0M; Centrality V0",101,-0.5,100.5);
-  fHOutCentFMD     = new TH1F("fHOutCentFMD","fHOutCentFMD; Centrality FMD",101,-0.5,100.5);
-  fHOutCentTRK     = new TH1F("fHOutCentTRK","fHOutCentTRK; Centrality TPC",101,-0.5,100.5);
-  fHOutCentTKL     = new TH1F("fHOutCentTKL","fHOutCentTKL; Centrality tracklets",101,-0.5,100.5);
-  fHOutCentCL0     = new TH1F("fHOutCentCL0","fHOutCentCL0; Centrality SPD inner",101,-0.5,100.5);
-  fHOutCentCL1     = new TH1F("fHOutCentCL1","fHOutCentCL1; Centrality SPD outer",101,-0.5,100.5);
-  fHOutCentV0MvsFMD= new TH1F("fHOutCentV0MvsFMD","fHOutCentV0MvsFMD; Centrality V0 vs FMD",101,-0.5,100.5);
-  fHOutCentTKLvsV0M= new TH1F("fHOutCentTKLvsV0M","fHOutCentTKLvsV0M; Centrality tracklets vs V0",101,-0.5,100.5);
-  fHOutCentZEMvsZDC= new TH1F("fHOutCentZEMvsZDC","fHOutCentZEMvsZDC; Centrality ZEM vs ZDC",101,-0.5,100.5);
+  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);
+
+  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     );
@@ -244,8 +378,36 @@ void AliCentralitySelectionTask::UserCreateOutputObjects()
   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();
 }
 
 //________________________________________________________________________
@@ -254,16 +416,16 @@ 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
   
   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
 
@@ -277,13 +439,17 @@ 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 (!esd) {
+       AliError("No ESD Event");
+       return;
+    }
+    
     if (fRunNo<=0) {
       if (SetupRun(esd)<0)
          AliFatal("Centrality File not available for this run");
@@ -296,9 +462,8 @@ void AliCentralitySelectionTask::UserExec(Option_t */*option*/)
     multV0A=esdV0->GetMTotV0A();
     multV0C=esdV0->GetMTotV0C();
 
-    float v0CorrR;
-    v0Corr = GetCorrV0(esd,v0CorrR,fRunNo);
-    v0Corr = (Short_t)v0Corr;
+    Float_t v0CorrR;
+    v0Corr = (Short_t)AliESDUtils::GetCorrV0(esd,v0CorrR);
     v0CorrResc = (Short_t)v0CorrR;
 
     // ***** Vertex Info
@@ -307,7 +472,6 @@ void AliCentralitySelectionTask::UserExec(Option_t */*option*/)
 
     // ***** CB info (tracklets, clusters, chips)
     //nTracks    = event->GetNumberOfTracks();     
-    fTrackCuts = AliESDtrackCuts::GetStandardTPCOnlyTrackCuts();
     nTracks    = fTrackCuts ? (Short_t)fTrackCuts->GetReferenceMultiplicity(esd,kTRUE):-1;
 
     const AliMultiplicity *mult = esd->GetMultiplicity();
@@ -322,7 +486,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();
@@ -364,8 +528,8 @@ void AliCentralitySelectionTask::UserExec(Option_t */*option*/)
     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));
+    zem1Energy = (Float_t) (esdZDC->GetZDCEMEnergy(0))/8.;
+    zem2Energy = (Float_t) (esdZDC->GetZDCEMEnergy(1))/8.;
     
   }   
   else if(fAnalysisInput.CompareTo("AOD")==0){
@@ -373,48 +537,105 @@ void AliCentralitySelectionTask::UserExec(Option_t */*option*/)
     // 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);
-
-  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);
+
+
+  // ***** Cleaning
+  fQuality=0;
+  fZVCut=10;
+  fOutliersCut=6;
+  
+  // ***** vertex
+  if (TMath::Abs(zvtx)>fZVCut) fQuality += 1;   
+
+  // ***** outliers
+  // **** V0 vs SPD
+  printf("AnalysisCentralitySelectionTask::centrality is %f\n",fCentV0M);
+
+  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)) 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); 
 }
@@ -424,6 +645,12 @@ void AliCentralitySelectionTask::ReadCentralityHistos(TString fCentfilename)
 {
   //  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"));
@@ -432,6 +659,15 @@ void AliCentralitySelectionTask::ReadCentralityHistos(TString fCentfilename)
   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();
 }  
 
 //________________________________________________________________________
@@ -439,11 +675,22 @@ void AliCentralitySelectionTask::ReadCentralityHistos2(TString fCentfilename2)
 {
   //  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  = (TH1F*) (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();
 }
 
 //________________________________________________________________________
@@ -471,149 +718,388 @@ 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("$ALICE_ROOT/ANALYSIS/macros/AliCentralityBy1D_%d.root", fRunNo));
-  TString fileName2(Form("$ALICE_ROOT/ANALYSIS/macros/AliCentralityByFunction_%d.root", fRunNo));
-  Bool_t isRunKnown = kFALSE;
+  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));
   
-  // Check if run exists
-  fFile  = TFile::Open(fileName);
-  if (fFile) {
-    isRunKnown = kTRUE;
-  }
-  // if not, take the last one  
-  if (!isRunKnown) {
-    //     fileName += (dynamic_cast<TObjString*>(fFileList->Last()))->GetString();
-    AliError(Form("Run %d not known to centrality selection!", fCurrentRun));
+  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;
+}
+
+//________________________________________________________________________
+Bool_t AliCentralitySelectionTask::IsOutlierV0MSPD(Float_t spd, Float_t v0, Int_t cent)
+{
+  TF1 *V0MSPDfun = new TF1("V0MSPDfun","-0.143789+ 0.288874*x",0,25000);
+  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};
+
+  if ( TMath::Abs(spd-V0MSPDfun->Eval(v0)) > fOutliersCut*SPDsigma[cent] ) 
+    return kTRUE;
+  else 
+    return kFALSE;
+}
+
+//________________________________________________________________________
+Bool_t AliCentralitySelectionTask::IsOutlierV0MTPC(Int_t tracks, Float_t v0, Int_t cent)
+{
+  TF1 *V0MTPCfun = new TF1("V0MTPCfun","-0.540691+0.128358*x",0,25000);
+  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 ( TMath::Abs(tracks-V0MTPCfun->Eval(v0)) > fOutliersCut*TPCsigma[cent] ) 
+    return kTRUE;
+  else 
+    return kFALSE;
+}
+
+//________________________________________________________________________
+Bool_t AliCentralitySelectionTask::IsOutlierV0MZDC(Float_t zdc, Float_t v0)
+{
+  TF1 *fun1 = new TF1("fun1","6350-0.26*x",0,25000);
+  TF1 *fun2 = new TF1("fun2","5580",0,25000);
+
+  if ( (zdc > fun1->Eval(v0)) || (zdc > fun2->Eval(v0)) )
+    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()));
-    ReadCentralityHistos(fileName.Data());
-    ReadCentralityHistos2(fileName2.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, int run) 
+
+//________________________________________________________________________  
+void AliCentralitySelectionTask::MyInitScaleFactor () 
 {
-  // correct V0 non-linearity, prepare a version rescaled to SPD2 corr
-  Double_t *par0;
-  Double_t *par1;
-  Double_t *par2;
-  
-  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 }; 
-  
+  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;
   
-  if (run==137161) {
-    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;
-    if (s<0) {
-      printf("FPE %d %.2f %.2f %.2e\n",i,esdV0->GetMultiplicity(i),b,(b*b-4.*par2[i]*esdV0->GetMultiplicity(i)*esdV0->GetMultiplicity(i)));
-      n = -b;
+  // 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;
+         }
+       }
+      }
     }
-    else {
-      n = (-b + TMath::Sqrt(s));
+  } // 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;
+         }
+       }
+      }
     }
-    multChCorr[i] = 2.*esdV0->GetMultiplicity(i)/n*par0[i];
-    multCorr += multChCorr[i];
-    multCorr2 += (multChCorr[i]/par0[i]/64.);
-  }
-  v0CorrResc =  multCorr2;
-  return multCorr;
+  } 
+
+  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;
+  
+  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;
+
 }
+