]> git.uio.no Git - u/mrichter/AliRoot.git/commitdiff
add calibrated combined vzero EP for LHC10h data
authorrbertens <rbertens@cern.ch>
Tue, 29 Jul 2014 14:25:24 +0000 (16:25 +0200)
committerrbertens <rbertens@cern.ch>
Tue, 29 Jul 2014 14:41:51 +0000 (16:41 +0200)
multiple updates for jet v2 task (should have been multiple commits actually)
a) added calibrated combined vzerp EP for LHC10h data
calibration is forked from the procedure used in the flow package (to avoid dependency)
combination relies on weight (chi) evaluated from event plane resolution of vzeroa and vzeroc
b) added functions to evaluated q-vectors (component wise) for different sub-detectors
c) converted helper functions to static functions
d) added missing debug flags
e) removed deprecated code (event selection)
f) fixed infinite loop bug in QA routine which occured when expected runnumbers did not coindice with detected runnumbers
g) fixed devision by zero bug for events with 2 (or less) accepted tracks

PWGJE/EMCALJetTasks/UserTasks/AliAnalysisTaskJetV2.cxx
PWGJE/EMCALJetTasks/UserTasks/AliAnalysisTaskJetV2.h

index caf9d4aea055c0d35f4781f21a69bf12a11b84c2..c7a515ef57ea26614465860d3558f9324061d86a 100644 (file)
 #include <TH1F.h>
 #include <TH2F.h>
 #include <TProfile.h>
+#include <TFile.h>
 // aliroot includes
 #include <AliAnalysisTask.h>
 #include <AliAnalysisManager.h>
 #include <AliCentrality.h>
 #include <AliVVertex.h>
 #include <AliVTrack.h>
+#include <AliVVZERO.h>
 #include <AliESDEvent.h>
 #include <AliAODEvent.h>
 #include <AliAODTrack.h>
+#include <AliOADBContainer.h>
 // emcal jet framework includes
 #include <AliPicoTrack.h>
 #include <AliEmcalJet.h>
@@ -62,7 +65,8 @@ using namespace std;
 ClassImp(AliAnalysisTaskJetV2)
 
 AliAnalysisTaskJetV2::AliAnalysisTaskJetV2() : AliAnalysisTaskEmcalJet("AliAnalysisTaskJetV2", kTRUE), 
-    fDebug(0), fRunToyMC(kFALSE), fLocalInit(0), fAttachToEvent(kTRUE), fFillHistograms(kTRUE), fFillQAHistograms(kTRUE), fReduceBinsXByFactor(-1.), fReduceBinsYByFactor(-1.), fNoEventWeightsForQC(kTRUE), fCentralityClasses(0), fExpectedRuns(0), fExpectedSemiGoodRuns(0), fUserSuppliedV2(0), fUserSuppliedV3(0), fUserSuppliedR2(0), fUserSuppliedR3(0), fTracksCont(0), fClusterCont(0), fJetsCont(0), fLeadingJet(0), fNAcceptedTracks(0), fNAcceptedTracksQCn(0), fFitModulationType(kNoFit), fFitGoodnessTest(kChi2Poisson), fQCRecovery(kTryFit), fUsePtWeight(kTRUE), fUsePtWeightErrorPropagation(kTRUE), fDetectorType(kTPC), fAnalysisType( kCharged), fFitModulationOptions("QWLI"), fRunModeType(kGrid), fDataType(kESD), fCollisionType(kPbPb), fRandom(0), fRunNumber(-1), fMappedRunNumber(0), fInCentralitySelection(-1), fFitModulation(0), fFitControl(0), fMinPvalue(0.01), fMaxPvalue(1), fNameSmallRho(""), fCachedRho(0), fSoftTrackMinPt(0.15), fSoftTrackMaxPt(5.), fSemiGoodJetMinPhi(0.), fSemiGoodJetMaxPhi(4.), fSemiGoodTrackMinPhi(0.), fSemiGoodTrackMaxPhi(4.), fAbsVertexZ(10), fHistCentrality(0), fHistVertexz(0), fHistRunnumbersPhi(0), fHistRunnumbersEta(0), fHistPvalueCDFROOT(0), fHistPvalueCDFROOTCent(0), fHistChi2ROOTCent(0), fHistPChi2Root(0),  fHistPvalueCDF(0), fHistPvalueCDFCent(0), fHistChi2Cent(0), fHistPChi2(0), fHistKolmogorovTest(0), fHistKolmogorovTestCent(0), fHistPKolmogorov(0), fHistRhoStatusCent(0), fHistUndeterminedRunQA(0), fMinDisanceRCtoLJ(0), fMaxCones(-1), fExcludeLeadingJetsFromFit(1.), fRebinSwapHistoOnTheFly(kTRUE), fPercentageOfFits(10.), fUseV0EventPlaneFromHeader(kTRUE), /*fExplicitOutlierCut(-1),*/ fOutputList(0), fOutputListGood(0), fOutputListBad(0), fHistAnalysisSummary(0), fHistSwap(0), fProfV2(0), fProfV2Cumulant(0), fProfV3(0), fProfV3Cumulant(0), fHistPsiControl(0), fHistPsiSpread(0), fHistPsiVZEROA(0), fHistPsiVZEROC(0), fHistPsiVZERO(0), fHistPsiTPC(0), fHistPsiVZEROAV0M(0), fHistPsiVZEROCV0M(0), fHistPsiVZEROVV0M(0), fHistPsiTPCiV0M(0), fHistPsiVZEROATRK(0), fHistPsiVZEROCTRK(0), fHistPsiVZEROTRK(0), fHistPsiTPCTRK(0), fHistRhoVsMult(0), fHistRhoVsCent(0), fHistRhoAVsMult(0), fHistRhoAVsCent(0) {
+    fDebug(0), fRunToyMC(kFALSE), fLocalInit(0), fAttachToEvent(kTRUE), fFillHistograms(kTRUE), fFillQAHistograms(kTRUE), fReduceBinsXByFactor(-1.), fReduceBinsYByFactor(-1.), fNoEventWeightsForQC(kTRUE), fCentralityClasses(0), fExpectedRuns(0), fExpectedSemiGoodRuns(0), fUserSuppliedV2(0), fUserSuppliedV3(0), fUserSuppliedR2(0), fUserSuppliedR3(0), fTracksCont(0), fClusterCont(0), fJetsCont(0), fLeadingJet(0), fNAcceptedTracks(0), fNAcceptedTracksQCn(0), fFitModulationType(kNoFit), fFitGoodnessTest(kChi2Poisson), fQCRecovery(kTryFit), fUsePtWeight(kTRUE), fUsePtWeightErrorPropagation(kTRUE), fDetectorType(kTPC), fAnalysisType( kCharged), fFitModulationOptions("QWLI"), fRunModeType(kGrid), fDataType(kESD), fCollisionType(kPbPb), fRandom(0), fRunNumber(-1), fMappedRunNumber(0), fInCentralitySelection(-1), fFitModulation(0), fFitControl(0), fMinPvalue(0.01), fMaxPvalue(1), fNameSmallRho(""), fCachedRho(0), fSoftTrackMinPt(0.15), fSoftTrackMaxPt(5.), fSemiGoodJetMinPhi(0.), fSemiGoodJetMaxPhi(4.), fSemiGoodTrackMinPhi(0.), fSemiGoodTrackMaxPhi(4.), fAbsVertexZ(10), fHistCentrality(0), fHistVertexz(0), fHistRunnumbersPhi(0), fHistRunnumbersEta(0), fHistPvalueCDFROOT(0), fHistPvalueCDFROOTCent(0), fHistChi2ROOTCent(0), fHistPChi2Root(0),  fHistPvalueCDF(0), fHistPvalueCDFCent(0), fHistChi2Cent(0), fHistPChi2(0), fHistKolmogorovTest(0), fHistKolmogorovTestCent(0), fHistPKolmogorov(0), fHistRhoStatusCent(0), fHistUndeterminedRunQA(0), fMinDisanceRCtoLJ(0), fMaxCones(-1), fExcludeLeadingJetsFromFit(1.), fRebinSwapHistoOnTheFly(kTRUE), fPercentageOfFits(10.), fOutputList(0), fOutputListGood(0), fOutputListBad(0), fHistAnalysisSummary(0), fHistSwap(0), fProfV2(0), fProfV2Cumulant(0), fProfV3(0), fProfV3Cumulant(0), fHistPsiControl(0), fHistPsiSpread(0), fHistPsiVZEROA(0), fHistPsiVZEROC(0), fHistPsiVZERO(0), fHistPsiTPC(0), fHistPsiVZEROAV0M(0), fHistPsiVZEROCV0M(0), fHistPsiVZEROVV0M(0), fHistPsiTPCiV0M(0), fHistPsiVZEROATRK(0), fHistPsiVZEROCTRK(0), fHistPsiVZEROTRK(0), fHistPsiTPCTRK(0), fHistRhoVsMult(0), fHistRhoVsCent(0), fHistRhoAVsMult(0), fHistRhoAVsCent(0), fVZEROgainEqualization(0x0), fVZEROgainEqualizationPerRing(kFALSE), fChi2A(0x0), fChi2C(0x0), fChi3A(0x0), fChi3C(0x0), fOADB(0x0)
+{
     for(Int_t i(0); i < 10; i++) {
         fProfV2Resolution[i] = 0;
         fProfV3Resolution[i] = 0;
@@ -96,11 +100,27 @@ AliAnalysisTaskJetV2::AliAnalysisTaskJetV2() : AliAnalysisTaskEmcalJet("AliAnaly
         fHistJetPsi2Pt[i] = 0;
         fHistJetPsi2PtRho0[i] = 0;
    }
-    // default constructor
+   for(Int_t i(0); i < 9; i++) {
+       for(Int_t j(0); j < 2; j++) {
+           for(Int_t k(0); k < 2; k++) {
+               fMeanQ[i][j][k] = 0.; 
+               fWidthQ[i][j][k] = 0.;  
+               fMeanQv3[i][j][k] = 0.; 
+               fWidthQv3[i][j][k] = 0.;
+           }
+       }
+   }
+   for(Int_t i(0); i < 4; i++) {
+       fVZEROApol[i] = 0.;
+       fVZEROCpol[i] = 0.;
+   }
+   for(Int_t i(0); i < 8; i++) fUseVZERORing[i] = kTRUE;
+   // default constructor
 }
 //_____________________________________________________________________________
 AliAnalysisTaskJetV2::AliAnalysisTaskJetV2(const char* name, runModeType type) : AliAnalysisTaskEmcalJet(name, kTRUE),
-  fDebug(0), fRunToyMC(kFALSE), fLocalInit(0), fAttachToEvent(kTRUE), fFillHistograms(kTRUE), fFillQAHistograms(kTRUE), fReduceBinsXByFactor(-1.), fReduceBinsYByFactor(-1.), fNoEventWeightsForQC(kTRUE), fCentralityClasses(0), fExpectedRuns(0), fExpectedSemiGoodRuns(0), fUserSuppliedV2(0), fUserSuppliedV3(0), fUserSuppliedR2(0), fUserSuppliedR3(0), fTracksCont(0), fClusterCont(0), fJetsCont(0), fLeadingJet(0), fNAcceptedTracks(0), fNAcceptedTracksQCn(0), fFitModulationType(kNoFit), fFitGoodnessTest(kChi2Poisson), fQCRecovery(kTryFit), fUsePtWeight(kTRUE), fUsePtWeightErrorPropagation(kTRUE), fDetectorType(kTPC), fAnalysisType(kCharged), fFitModulationOptions("QWLI"), fRunModeType(type), fDataType(kESD), fCollisionType(kPbPb), fRandom(0), fRunNumber(-1), fMappedRunNumber(0), fInCentralitySelection(-1), fFitModulation(0), fFitControl(0), fMinPvalue(0.01), fMaxPvalue(1), fNameSmallRho(""), fCachedRho(0), fSoftTrackMinPt(0.15), fSoftTrackMaxPt(5.), fSemiGoodJetMinPhi(0.), fSemiGoodJetMaxPhi(4.), fSemiGoodTrackMinPhi(0.), fSemiGoodTrackMaxPhi(4.), fAbsVertexZ(10), fHistCentrality(0), fHistVertexz(0), fHistRunnumbersPhi(0), fHistRunnumbersEta(0), fHistPvalueCDFROOT(0), fHistPvalueCDFROOTCent(0), fHistChi2ROOTCent(0), fHistPChi2Root(0),  fHistPvalueCDF(0), fHistPvalueCDFCent(0), fHistChi2Cent(0), fHistPChi2(0), fHistKolmogorovTest(0), fHistKolmogorovTestCent(0), fHistPKolmogorov(0), fHistRhoStatusCent(0), fHistUndeterminedRunQA(0), fMinDisanceRCtoLJ(0), fMaxCones(-1), fExcludeLeadingJetsFromFit(1.), fRebinSwapHistoOnTheFly(kTRUE), fPercentageOfFits(10.), fUseV0EventPlaneFromHeader(kTRUE), /*fExplicitOutlierCut(-1),*/ fOutputList(0), fOutputListGood(0), fOutputListBad(0), fHistAnalysisSummary(0), fHistSwap(0), fProfV2(0), fProfV2Cumulant(0), fProfV3(0), fProfV3Cumulant(0), fHistPsiControl(0), fHistPsiSpread(0), fHistPsiVZEROA(0), fHistPsiVZEROC(0), fHistPsiVZERO(0), fHistPsiTPC(0), fHistPsiVZEROAV0M(0), fHistPsiVZEROCV0M(0), fHistPsiVZEROVV0M(0), fHistPsiTPCiV0M(0), fHistPsiVZEROATRK(0), fHistPsiVZEROCTRK(0), fHistPsiVZEROTRK(0), fHistPsiTPCTRK(0), fHistRhoVsMult(0), fHistRhoVsCent(0), fHistRhoAVsMult(0), fHistRhoAVsCent(0) {
+  fDebug(0), fRunToyMC(kFALSE), fLocalInit(0), fAttachToEvent(kTRUE), fFillHistograms(kTRUE), fFillQAHistograms(kTRUE), fReduceBinsXByFactor(-1.), fReduceBinsYByFactor(-1.), fNoEventWeightsForQC(kTRUE), fCentralityClasses(0), fExpectedRuns(0), fExpectedSemiGoodRuns(0), fUserSuppliedV2(0), fUserSuppliedV3(0), fUserSuppliedR2(0), fUserSuppliedR3(0), fTracksCont(0), fClusterCont(0), fJetsCont(0), fLeadingJet(0), fNAcceptedTracks(0), fNAcceptedTracksQCn(0), fFitModulationType(kNoFit), fFitGoodnessTest(kChi2Poisson), fQCRecovery(kTryFit), fUsePtWeight(kTRUE), fUsePtWeightErrorPropagation(kTRUE), fDetectorType(kTPC), fAnalysisType(kCharged), fFitModulationOptions("QWLI"), fRunModeType(type), fDataType(kESD), fCollisionType(kPbPb), fRandom(0), fRunNumber(-1), fMappedRunNumber(0), fInCentralitySelection(-1), fFitModulation(0), fFitControl(0), fMinPvalue(0.01), fMaxPvalue(1), fNameSmallRho(""), fCachedRho(0), fSoftTrackMinPt(0.15), fSoftTrackMaxPt(5.), fSemiGoodJetMinPhi(0.), fSemiGoodJetMaxPhi(4.), fSemiGoodTrackMinPhi(0.), fSemiGoodTrackMaxPhi(4.), fAbsVertexZ(10), fHistCentrality(0), fHistVertexz(0), fHistRunnumbersPhi(0), fHistRunnumbersEta(0), fHistPvalueCDFROOT(0), fHistPvalueCDFROOTCent(0), fHistChi2ROOTCent(0), fHistPChi2Root(0),  fHistPvalueCDF(0), fHistPvalueCDFCent(0), fHistChi2Cent(0), fHistPChi2(0), fHistKolmogorovTest(0), fHistKolmogorovTestCent(0), fHistPKolmogorov(0), fHistRhoStatusCent(0), fHistUndeterminedRunQA(0), fMinDisanceRCtoLJ(0), fMaxCones(-1), fExcludeLeadingJetsFromFit(1.), fRebinSwapHistoOnTheFly(kTRUE), fPercentageOfFits(10.), fOutputList(0), fOutputListGood(0), fOutputListBad(0), fHistAnalysisSummary(0), fHistSwap(0), fProfV2(0), fProfV2Cumulant(0), fProfV3(0), fProfV3Cumulant(0), fHistPsiControl(0), fHistPsiSpread(0), fHistPsiVZEROA(0), fHistPsiVZEROC(0), fHistPsiVZERO(0), fHistPsiTPC(0), fHistPsiVZEROAV0M(0), fHistPsiVZEROCV0M(0), fHistPsiVZEROVV0M(0), fHistPsiTPCiV0M(0), fHistPsiVZEROATRK(0), fHistPsiVZEROCTRK(0), fHistPsiVZEROTRK(0), fHistPsiTPCTRK(0), fHistRhoVsMult(0), fHistRhoVsCent(0), fHistRhoAVsMult(0), fHistRhoAVsCent(0), fVZEROgainEqualization(0x0), fVZEROgainEqualizationPerRing(kFALSE), fChi2A(0x0), fChi2C(0x0), fChi3A(0x0), fChi3C(0x0), fOADB(0x0)
+{
     for(Int_t i(0); i < 10; i++) {
         fProfV2Resolution[i] = 0;
         fProfV3Resolution[i] = 0;
@@ -133,7 +153,23 @@ AliAnalysisTaskJetV2::AliAnalysisTaskJetV2(const char* name, runModeType type) :
         fHistJetEtaRho[i] = 0;
         fHistJetPsi2Pt[i] = 0;
         fHistJetPsi2PtRho0[i] = 0;
-    }
+   }
+   for(Int_t i(0); i < 9; i++) {
+       for(Int_t j(0); j < 2; j++) {
+           for(Int_t k(0); k < 2; k++) {
+               fMeanQ[i][j][k] = 0.; 
+               fWidthQ[i][j][k] = 0.;  
+               fMeanQv3[i][j][k] = 0.; 
+               fWidthQv3[i][j][k] = 0.;
+           }
+       }
+   }
+   for(Int_t i(0); i < 4; i++) {
+       fVZEROApol[i] = 0.;
+       fVZEROCpol[i] = 0.;
+   }
+   for(Int_t i(0); i < 8; i++) fUseVZERORing[i] = kTRUE;
+
     // constructor
     DefineInput(0, TChain::Class());
     DefineOutput(1, TList::Class());
@@ -167,6 +203,14 @@ AliAnalysisTaskJetV2::~AliAnalysisTaskJetV2()
     if(fExpectedRuns)           {delete fExpectedRuns;          fExpectedRuns = 0x0;}
     if(fExpectedSemiGoodRuns)   {delete fExpectedSemiGoodRuns;  fExpectedSemiGoodRuns = 0x0;}
     if(fFitControl)             {delete fFitControl;            fFitControl = 0x0;}
+    if(fVZEROgainEqualization)  {delete fVZEROgainEqualization; fVZEROgainEqualization = 0x0;}
+    if(fChi2A)                  {delete fChi2A;                 fChi2A = 0x0;}
+    if(fChi2C)                  {delete fChi2C;                 fChi2C = 0x0;}
+    if(fChi3A)                  {delete fChi3A;                 fChi3A = 0x0;}
+    if(fChi3C)                  {delete fChi3C;                 fChi3C = 0x0;}
+    if(fOADB && !fOADB->IsZombie()) {
+        fOADB->Close();        fOADB = 0x0;
+    } else if (fOADB) fOADB = 0x0;
 }
 //_____________________________________________________________________________
 void AliAnalysisTaskJetV2::ExecOnce()
@@ -193,17 +237,22 @@ Bool_t AliAnalysisTaskJetV2::Notify()
     if(fRunNumber != InputEvent()->GetRunNumber()) {
         fRunNumber = InputEvent()->GetRunNumber();        // set the current run number
         if(fDebug > 0) printf("__FUNC__ %s > NEW RUNNUMBER DETECTED \n ", __func__);
-        // check if this is 10h or 11h data. for 10h we don't want to change the acceptance
+        // check if this is 10h or 11h data
         switch (fCollisionType) {
             case kPbPb10h : {
+                if(fDebug > 0) printf(" LHC10h data, assuming full acceptance, reading VZERO calibration DB \n ");
+                // for 10h data the vzero event plane calibration needs to be cached
+                ReadVZEROCalibration2010h(); 
+                // no need to change rho or acceptance for 10h, so we're done
                 return kTRUE;
             } break;
-            default : break;
+            default :  {
+                if(fDebug > 0) printf(" checking runnumber to adjust acceptance on the fly \n");           
+            } break;
         }
         // reset the cuts. should be a pointless operation except for the case where the run number changes
         // from semi-good back to good on one node, which is not a likely scenario (unless trains will
         // run as one masterjob)
-        AliAnalysisTaskEmcal::SetTrackPhiLimits(-10., 10.);
         switch (fAnalysisType) {
             case kCharged: {
                 AliAnalysisTaskEmcalJet::SetJetPhiLimits(-10., 10.);   
@@ -211,7 +260,9 @@ Bool_t AliAnalysisTaskJetV2::Notify()
             case kFull: {
                 AliAnalysisTaskEmcalJet::SetJetPhiLimits(1.405 + GetJetRadius(), 3.135 - GetJetRadius());
             } break;
-            default: break;
+            default: {
+                AliAnalysisTaskEmcal::SetTrackPhiLimits(-10., 10.);
+            } break;
         }
         if(fCachedRho) {                // if there's a cached rho, it's the default, so switch back
             if(fDebug > 0) printf("__FUNC__ %s > replacing rho with cached rho \n ", __func__);
@@ -573,7 +624,7 @@ void AliAnalysisTaskJetV2::UserCreateOutputObjects()
 //_____________________________________________________________________________
 Bool_t AliAnalysisTaskJetV2::Run()
 {
-    // user exec: execute once for each event
+    // called for each accepted event (call made from user exec of parent class)
     if(fDebug > 0) printf("__FILE__ = %s \n __LINE __ %i , __FUNC__ %s \n ", __FILE__, __LINE__, __func__);
     if(!fTracks||!fJets||!fRho) return kFALSE;
     if(!fLocalInit) fLocalInit = InitializeAnalysis();
@@ -707,41 +758,68 @@ Bool_t AliAnalysisTaskJetV2::Run()
     return kTRUE;
 }
 //_____________________________________________________________________________
+Double_t AliAnalysisTaskJetV2::CalculateEventPlaneChi(Double_t res)
+{
+    // return chi for given resolution FIXME needs citation and explanation (29072014)
+    Double_t chi(2.), delta(1.), con(0.626657);
+    for (Int_t i(0); i < 15; i++) {
+        chi = ((con*chi*TMath::Exp(-chi*chi/4.)*(TMath::BesselI0(chi*chi/4.)+TMath::BesselI1(chi*chi/4.))) < res) ? chi + delta : chi - delta;
+        delta = delta / 2.;
+    }
+    return chi;
+}
+//_____________________________________________________________________________
 void AliAnalysisTaskJetV2::CalculateEventPlaneVZERO(Double_t vzero[2][2]) const 
 {
-    // get the vzero event plane
+    // get the vzero event plane (a and c separately)
     if(fDebug > 0) printf("__FILE__ = %s \n __LINE __ %i , __FUNC__ %s \n ", __FILE__, __LINE__, __func__);
-    if(fUseV0EventPlaneFromHeader) {    // use the vzero from the header
-        Double_t a(0), b(0), c(0), d(0), e(0), f(0), g(0), h(0);
-        vzero[0][0] = InputEvent()->GetEventplane()->CalculateVZEROEventPlane(InputEvent(), 8, 2, a, b);
-        vzero[1][0] = InputEvent()->GetEventplane()->CalculateVZEROEventPlane(InputEvent(), 9, 2, c, d);
-        vzero[0][1] = InputEvent()->GetEventplane()->CalculateVZEROEventPlane(InputEvent(), 8, 3, e, f);
-        vzero[1][1] = InputEvent()->GetEventplane()->CalculateVZEROEventPlane(InputEvent(), 9, 3, g, h);
-        return;
+    switch (fCollisionType) {
+        case kPbPb10h : {
+            // for 10h data, get the calibrated q-vector from the database
+            Double_t QA2[] = {-999., -999.};
+            Double_t QA3[] = {-999., -999.};
+            Double_t QC2[] = {-999., -999.};
+            Double_t QC3[] = {-999., -999.};
+            CalculateQvectorVZERO(QA2, QA3, QC2, QC3);
+            vzero[0][0] = .5*TMath::ATan2(QA2[1], QA2[0]);
+            vzero[1][0] = .5*TMath::ATan2(QC2[1], QC2[0]);
+            vzero[0][1] = (1./3.)*TMath::ATan2(QA3[1], QA3[0]);
+            vzero[1][1] = (1./3.)*TMath::ATan2(QC3[1], QC3[0]);
+        } break;
+        default: {
+            // by default use the ep from the event header (make sure EP selection task is enabeled!)
+            Double_t a(0), b(0), c(0), d(0), e(0), f(0), g(0), h(0);
+            vzero[0][0] = InputEvent()->GetEventplane()->CalculateVZEROEventPlane(InputEvent(), 8, 2, a, b);
+            vzero[1][0] = InputEvent()->GetEventplane()->CalculateVZEROEventPlane(InputEvent(), 9, 2, c, d);
+            vzero[0][1] = InputEvent()->GetEventplane()->CalculateVZEROEventPlane(InputEvent(), 8, 3, e, f);
+            vzero[1][1] = InputEvent()->GetEventplane()->CalculateVZEROEventPlane(InputEvent(), 9, 3, g, h);
+            return;
+        }
     }
-    // grab the vzero event plane without recentering
+}
+//_____________________________________________________________________________
+void AliAnalysisTaskJetV2::CalculateEventPlaneCombinedVZERO(Double_t* comb) const
+{
+    // return the combined vzero event plane
     if(fDebug > 0) printf("__FILE__ = %s \n __LINE __ %i , __FUNC__ %s \n ", __FILE__, __LINE__, __func__);
-    Double_t qxa2(0), qya2(0), qxc2(0), qyc2(0);    // for psi2
-    Double_t qxa3(0), qya3(0), qxc3(0), qyc3(0);    // for psi3
-    for(Int_t iVZERO(0); iVZERO < 64; iVZERO++) {
-        Double_t phi(TMath::PiOver4()*(.5+iVZERO%8)), /* eta(0), */ weight(InputEvent()->GetVZEROEqMultiplicity(iVZERO));
-        if(iVZERO<32) {
-            qxa2 += weight*TMath::Cos(2.*phi);
-            qya2 += weight*TMath::Sin(2.*phi);
-            qxa3 += weight*TMath::Cos(3.*phi);
-            qya3 += weight*TMath::Sin(3.*phi);
-        }
-        else {
-            qxc2 += weight*TMath::Cos(2.*phi);
-            qyc2 += weight*TMath::Sin(2.*phi);
-            qxc3 += weight*TMath::Cos(3.*phi);
-            qyc3 += weight*TMath::Sin(3.*phi);
-       }
+    switch (fCollisionType) {
+        // for 10h data call calibration info
+        case kPbPb10h : {
+            // get the calibrated q-vectors
+            Double_t Q2[] = {-999., -999.};            
+            Double_t Q3[] = {-999., -999.};
+            // return if something isn't ok from the calibration side
+            CalculateQvectorCombinedVZERO(Q2, Q3);
+            comb[0] = .5*TMath::ATan2(Q2[1], Q2[0]);
+            comb[1] = (1./3.)*TMath::ATan2(Q3[1], Q3[0]);
+        } break;
+        default : {
+            // for all other types use calibrated event plane from the event header
+            Double_t a(0), b(0), c(0), d(0);
+            comb[0] = InputEvent()->GetEventplane()->CalculateVZEROEventPlane(InputEvent(), 10, 2, a, b);
+            comb[1] = InputEvent()->GetEventplane()->CalculateVZEROEventPlane(InputEvent(), 10, 3, c, d);
+        } break;
     }
-    vzero[0][0] = .5*TMath::ATan2(qya2, qxa2);
-    vzero[1][0] = .5*TMath::ATan2(qyc2, qxc2);
-    vzero[0][1] = (1./3.)*TMath::ATan2(qya3, qxa3);
-    vzero[1][1] = (1./3.)*TMath::ATan2(qyc3, qxc3);
 }
 //_____________________________________________________________________________
 void AliAnalysisTaskJetV2::CalculateEventPlaneTPC(Double_t* tpc)
@@ -769,21 +847,12 @@ void AliAnalysisTaskJetV2::CalculateEventPlaneTPC(Double_t* tpc)
    }
    tpc[0] = .5*TMath::ATan2(qy2, qx2);
    tpc[1] = (1./3.)*TMath::ATan2(qy3, qx3);
-} 
-//_____________________________________________________________________________
-void AliAnalysisTaskJetV2::CalculateEventPlaneCombinedVZERO(Double_t* comb) const
-{
-    // grab the combined vzero event plane
-    if(fDebug > 0) printf("__FILE__ = %s \n __LINE __ %i , __FUNC__ %s \n ", __FILE__, __LINE__, __func__);
-    Double_t a(0), b(0), c(0), d(0);
-    comb[0] = InputEvent()->GetEventplane()->CalculateVZEROEventPlane(InputEvent(), 10, 2, a, b);
-    comb[1] = InputEvent()->GetEventplane()->CalculateVZEROEventPlane(InputEvent(), 10, 3, c, d);
 }
 //_____________________________________________________________________________
 void AliAnalysisTaskJetV2::CalculateEventPlaneResolution(Double_t vzero[2][2], Double_t* vzeroComb, Double_t* tpc)
 {
     // fill the profiles for the resolution parameters
-    if(fDebug > 1) printf("__FILE__ = %s \n __LINE __ %i , __FUNC__ %s \n ", __FILE__, __LINE__, __func__);
+    if(fDebug > 0) printf("__FILE__ = %s \n __LINE __ %i , __FUNC__ %s \n ", __FILE__, __LINE__, __func__);
     fProfV2Resolution[fInCentralitySelection]->Fill(2., TMath::Cos(2.*(vzero[0][0] - vzero[1][0])));
     fProfV2Resolution[fInCentralitySelection]->Fill(3., TMath::Cos(2.*(vzero[1][0] - vzero[0][0])));
     fProfV2Resolution[fInCentralitySelection]->Fill(4., TMath::Cos(2.*(vzero[0][0] - tpc[0])));
@@ -831,6 +900,109 @@ void AliAnalysisTaskJetV2::CalculateEventPlaneResolution(Double_t vzero[2][2], D
    fProfV3Resolution[fInCentralitySelection]->Fill(10., TMath::Cos(3.*(tpca3 - tpcb3))); 
 }   
 //_____________________________________________________________________________
+void AliAnalysisTaskJetV2::CalculateQvectorVZERO(Double_t Qa2[2], Double_t Qc2[2], Double_t Qa3[2], Double_t Qc3[2]) const
+{
+    // return the calibrated 2nd and 3rd order q-vectors for vzeroa and vzeroc
+    // function takes arrays as arguments, which correspond to vzero info in the following way
+    // 
+    // Qa2[0] = Qx2 for vzero A         Qa2[1] = Qy2 for vzero A (etc)
+    
+    if(fDebug > 0) printf("__FILE__ = %s \n __LINE __ %i , __FUNC__ %s \n ", __FILE__, __LINE__, __func__);
+    // placeholders for geometric information
+    Double_t phi(-999.), weight(-999.); 
+    // reset placeholders for Q-vector components
+    Qa2[0] = 0.;    Qc2[0] = 0.;    Qa3[0] = 0.;    Qc3[0] = 0.;
+    Qa2[1] = 0.;    Qc2[1] = 0.;    Qa3[1] = 0.;    Qc3[1] = 0.;
+    
+    for(Int_t i(0); i < 64; i++) {
+        // loop over all scintillators, construct Q-vectors in the same loop
+        phi     = TMath::PiOver4()*(0.5+i%8);
+        weight  = 0.;
+        // note that disabled rings have already been excluded in ReadVZEROCalibration2010h
+        if(i<32) {    // v0c side
+            if(i < 8) weight = InputEvent()->GetVZEROData()->GetMultiplicity(i)*fVZEROCpol[0]/fVZEROgainEqualization->GetBinContent(1+i);
+            else if (i < 16 ) weight = InputEvent()->GetVZEROData()->GetMultiplicity(i)*fVZEROCpol[1]/fVZEROgainEqualization->GetBinContent(1+i);
+            else if (i < 24 ) weight = InputEvent()->GetVZEROData()->GetMultiplicity(i)*fVZEROCpol[2]/fVZEROgainEqualization->GetBinContent(1+i);
+            else if (i < 32 ) weight = InputEvent()->GetVZEROData()->GetMultiplicity(i)*fVZEROCpol[3]/fVZEROgainEqualization->GetBinContent(1+i);
+            // fill Q-vectors for v0c side
+            Qc2[0]+=weight*TMath::Cos(2.*phi);
+            Qc3[0]+=weight*TMath::Cos(3.*phi);
+            Qc2[1]+=weight*TMath::Sin(2.*phi);
+            Qc3[1]+=weight*TMath::Sin(3.*phi);
+        } else {       // v0a side
+            if( i < 40) weight = InputEvent()->GetVZEROData()->GetMultiplicity(i)*fVZEROApol[0]/fVZEROgainEqualization->GetBinContent(1+i);
+            else if ( i < 48 ) weight = InputEvent()->GetVZEROData()->GetMultiplicity(i)*fVZEROApol[1]/fVZEROgainEqualization->GetBinContent(1+i);
+            else if ( i < 56 ) weight = InputEvent()->GetVZEROData()->GetMultiplicity(i)*fVZEROApol[2]/fVZEROgainEqualization->GetBinContent(1+i);
+            else if ( i < 64 ) weight = InputEvent()->GetVZEROData()->GetMultiplicity(i)*fVZEROApol[3]/fVZEROgainEqualization->GetBinContent(1+i);
+            // fill Q-vectors for v0a side
+            Qa2[0]+=weight*TMath::Cos(2.*phi);
+            Qa3[0]+=weight*TMath::Cos(3.*phi);
+            Qa2[1]+=weight*TMath::Sin(2.*phi);
+            Qa3[1]+=weight*TMath::Sin(3.*phi);
+        }
+    }
+    // get the cache index and read the correction terms from the cache
+    Int_t VZEROcentralityBin(GetVZEROCentralityBin());
+    Double_t Qx2amean = fMeanQ[VZEROcentralityBin][1][0];
+    Double_t Qx2arms  = fWidthQ[VZEROcentralityBin][1][0];
+    Double_t Qy2amean = fMeanQ[VZEROcentralityBin][1][1];
+    Double_t Qy2arms  = fWidthQ[VZEROcentralityBin][1][1];
+
+    Double_t Qx2cmean = fMeanQ[VZEROcentralityBin][0][0];
+    Double_t Qx2crms  = fWidthQ[VZEROcentralityBin][0][0];
+    Double_t Qy2cmean = fMeanQ[VZEROcentralityBin][0][1];
+    Double_t Qy2crms  = fWidthQ[VZEROcentralityBin][0][1];     
+
+    Double_t Qx3amean = fMeanQv3[VZEROcentralityBin][1][0];
+    Double_t Qx3arms  = fWidthQv3[VZEROcentralityBin][1][0];
+    Double_t Qy3amean = fMeanQv3[VZEROcentralityBin][1][1];
+    Double_t Qy3arms  = fWidthQv3[VZEROcentralityBin][1][1];
+
+    Double_t Qx3cmean = fMeanQv3[VZEROcentralityBin][0][0];
+    Double_t Qx3crms  = fWidthQv3[VZEROcentralityBin][0][0];
+    Double_t Qy3cmean = fMeanQv3[VZEROcentralityBin][0][1];
+    Double_t Qy3crms  = fWidthQv3[VZEROcentralityBin][0][1];   
+
+    // update the weighted q-vectors with the re-centered values
+    Qa2[0] = (Qa2[0] - Qx2amean)/Qx2arms;
+    Qa2[1] = (Qa2[1] - Qy2amean)/Qy2arms;
+    Qc2[0] = (Qc2[0] - Qx2cmean)/Qx2crms;
+    Qc2[1] = (Qc2[1] - Qy2cmean)/Qy2crms;
+
+    Qa3[0] = (Qa3[0] - Qx3amean)/Qx3arms;
+    Qa3[1] = (Qa3[1] - Qy3amean)/Qy3arms;
+    Qc3[0] = (Qc3[0] - Qx3cmean)/Qx3crms;
+    Qc3[1] = (Qc3[0] - Qy3cmean)/Qy3crms;
+}
+//_____________________________________________________________________________
+void AliAnalysisTaskJetV2::CalculateQvectorCombinedVZERO(Double_t Q2[2], Double_t Q3[2]) const
+{
+    // calculate calibrated q-vector of the combined vzeroa, vzeroc system
+    // this is somewhat ugly as CalculateQvectorCombinedVZERO is called more than once per event
+    // but for now it will have to do ...
+    if(fDebug > 0) printf("__FILE__ = %s \n __LINE __ %i , __FUNC__ %s \n ", __FILE__, __LINE__, __func__);
+
+    // first step: retrieve the q-vectors component-wise per vzero detector
+    Double_t QA2[] = {-999., -999.};
+    Double_t QA3[] = {-999., -999.};
+    Double_t QC2[] = {-999., -999.};
+    Double_t QC3[] = {-999., -999.};
+    CalculateQvectorVZERO(QA2, QA3, QC2, QC3);
+
+    // get cache index and retrieve the chi weights for this centrality
+    Int_t VZEROcentralityBin(GetVZEROCentralityBin());
+    Double_t chi2A(fChi2A->At(VZEROcentralityBin));
+    Double_t chi2C(fChi2C->At(VZEROcentralityBin));
+    Double_t chi3A(fChi3A->At(VZEROcentralityBin));
+    Double_t chi3C(fChi3C->At(VZEROcentralityBin));
+
+    // combine the vzera and vzeroc signal
+    Q2[0] = chi2A*chi2A*QA2[0]+chi2C*chi2C*QC2[0];
+    Q2[1] = chi2A*chi2A*QA2[1]+chi2C*chi2C*QC2[1];
+    Q3[0] = chi3A*chi3A*QA3[0]+chi3C*chi3C*QC3[0];
+    Q3[1] = chi3A*chi3A*QC3[1]+chi3C*chi3C*QC3[1];
+}
+//_____________________________________________________________________________
 void AliAnalysisTaskJetV2::CalculateRandomCone(Float_t &pt, Float_t &eta, Float_t &phi, 
         AliParticleContainer* tracksCont, AliClusterContainer* clusterCont, AliEmcalJet* jet) const
 {
@@ -960,7 +1132,8 @@ void AliAnalysisTaskJetV2::QCnDiffentialFlowVectors(
         TClonesArray* pois, TArrayD* ptBins, Bool_t vpart, Double_t* repn, Double_t* impn, 
         Double_t *mp, Double_t *reqn, Double_t *imqn, Double_t* mq, Int_t n) 
 {
-    // get  unweighted differential flow vectors
+     if(fDebug > 0) printf("__FILE__ = %s \n __LINE __ %i , __FUNC__ %s \n ", __FILE__, __LINE__, __func__);
+   // get  unweighted differential flow vectors
     Int_t iPois(pois->GetEntriesFast());
     if(vpart) {
         for(Int_t i(0); i < iPois; i++) {
@@ -1185,7 +1358,6 @@ Bool_t AliAnalysisTaskJetV2::CorrectRho(Double_t psi2, Double_t psi3)
     Double_t lowBound(0.), upBound(TMath::TwoPi());     // bounds for fit
     if(GetParticleContainer()->GetParticlePhiMin() > lowBound) lowBound = GetParticleContainer()->GetParticlePhiMin();
     if(GetParticleContainer()->GetParticlePhiMax() < upBound) upBound = GetParticleContainer()->GetParticlePhiMax();
-
     fHistSwap->Reset(); // clear the histogram
     TH1F _tempSwap;     // on stack for quick access
     TH1F _tempSwapN;    // on stack for quick access, bookkeeping histogram
@@ -1199,26 +1371,26 @@ Bool_t AliAnalysisTaskJetV2::CorrectRho(Double_t psi2, Double_t psi3)
     // non poissonian error when using pt weights
     Double_t totalpts(0.), totalptsquares(0.), totalns(0.);
     for(Int_t i(0); i < iTracks; i++) {
-            AliVTrack* track = static_cast<AliVTrack*>(fTracks->At(i));
-            if(fExcludeLeadingJetsFromFit > 0 &&( (TMath::Abs(track->Eta() - excludeInEta) < GetJetContainer()->GetJetRadius()*fExcludeLeadingJetsFromFit ) || (TMath::Abs(track->Eta()) - GetJetContainer()->GetJetRadius() - GetJetContainer()->GetJetEtaMax() ) > 0 )) continue;
-            if(!PassesCuts(track) || track->Pt() > fSoftTrackMaxPt || track->Pt() < fSoftTrackMinPt) continue;
-            if(fUsePtWeight) {
-                _tempSwap.Fill(track->Phi(), track->Pt());
-                if(fUsePtWeightErrorPropagation) {
-                    totalpts += track->Pt();
-                   totalptsquares += track->Pt()*track->Pt();
-                   totalns += 1;
-                    _tempSwapN.Fill(track->Phi());
-                }
+        AliVTrack* track = static_cast<AliVTrack*>(fTracks->At(i));
+        if(fExcludeLeadingJetsFromFit > 0 &&( (TMath::Abs(track->Eta() - excludeInEta) < GetJetContainer()->GetJetRadius()*fExcludeLeadingJetsFromFit ) || (TMath::Abs(track->Eta()) - GetJetContainer()->GetJetRadius() - GetJetContainer()->GetJetEtaMax() ) > 0 )) continue;
+        if(!PassesCuts(track) || track->Pt() > fSoftTrackMaxPt || track->Pt() < fSoftTrackMinPt) continue;
+        if(fUsePtWeight) {
+            _tempSwap.Fill(track->Phi(), track->Pt());
+            if(fUsePtWeightErrorPropagation) {
+                totalpts += track->Pt();
+                totalptsquares += track->Pt()*track->Pt();
+               totalns += 1;
+                _tempSwapN.Fill(track->Phi());
             }
-            else _tempSwap.Fill(track->Phi());
+        }
+        else _tempSwap.Fill(track->Phi());
     }
     if(fUsePtWeight && fUsePtWeightErrorPropagation) {
         // in the case of pt weights overwrite the poissonian error estimate which is assigned by root by a more sophisticated appraoch
         // the assumption here is that the bin error will be dominated by the uncertainty in the mean pt in a bin and in the uncertainty
         // of the number of tracks in a bin, the first of which will be estimated from the sample standard deviation of all tracks in the 
         // event, for the latter use a poissonian estimate. the two contrubitions are assumed to be uncorrelated
-        if(totalns < 1) return kFALSE; // not one track passes the cuts
+        if(totalns < 2) return kFALSE; // not one track passes the cuts > 2 avoids possible division by 0 later on
         for(Int_t l = 0; l < _tempSwap.GetNbinsX(); l++) {
             if(_tempSwapN.GetBinContent(l+1) == 0) {
                 _tempSwap.SetBinContent(l+1,0);
@@ -1240,7 +1412,6 @@ Bool_t AliAnalysisTaskJetV2::CorrectRho(Double_t psi2, Double_t psi3)
             }
         }
     }
-
     fFitModulation->SetParameter(0, fLocalRho->GetVal());
     switch (fFitModulationType) {
         case kNoFit : { 
@@ -1299,7 +1470,7 @@ Bool_t AliAnalysisTaskJetV2::CorrectRho(Double_t psi2, Double_t psi3)
     // the quality of the fit is evaluated from 1 - the cdf of the chi square distribution
     // three methods are available, all with their drawbacks. all are stored, one is selected to do the cut
     Int_t NDF(_tempSwap.GetXaxis()->GetNbins()-freeParams);
-    if(NDF == 0) return kFALSE;
+    if(NDF == 0 || (float)NDF <= 0.) return kFALSE;
     Double_t CDF(1.-ChiSquareCDF(NDF, ChiSquare(_tempSwap, fFitModulation)));
     Double_t CDFROOT(1.-ChiSquareCDF(NDF, fFitModulation->GetChisquare()));
     Double_t CDFKolmogorov(KolmogorovTest(_tempSwap, fFitModulation));
@@ -1349,9 +1520,10 @@ Bool_t AliAnalysisTaskJetV2::CorrectRho(Double_t psi2, Double_t psi3)
             fHistRhoStatusCent->Fill(fCent, -1);
         }
     }
-    if(CDF >= fMinPvalue && CDF <= fMaxPvalue && ( fFitModulation->GetMinimum(0, TMath::TwoPi()) > 0)) {       // fit quality. not that although with limited acceptance the fit is performed on just
-                       // part of phase space, the requirement that energy desntiy is larger than zero is applied
-                       // to the FULL spectrum
+    if(CDF >= fMinPvalue && CDF <= fMaxPvalue && ( fFitModulation->GetMinimum(0, TMath::TwoPi()) > 0)) {       
+        // fit quality. not that although with limited acceptance the fit is performed on just
+        // part of phase space, the requirement that energy desntiy is larger than zero is applied
+        // to the FULL spectrum
         fHistRhoStatusCent->Fill(fCent, 0.);
         // for LOCAL didactic purposes, save the  best and the worst fits
         // this routine can produce a lot of output histograms (it's not memory 'safe') and will not work on GRID 
@@ -1469,9 +1641,6 @@ Bool_t AliAnalysisTaskJetV2::PassesCuts(AliVEvent* event)
             break; }
     } 
     if(fInCentralitySelection<0) return kFALSE;     // should be null op
-/*    if(fExplicitOutlierCut == 2010 || fExplicitOutlierCut == 2011) {
-       if(!PassesCuts(fExplicitOutlierCut)) return kFALSE;
-    }*/
     // see if input containers are filled
     if(fTracks->GetEntries() < 1) return kFALSE;
     if(fRho->GetVal() <= 0 ) return kFALSE;
@@ -1479,29 +1648,6 @@ Bool_t AliAnalysisTaskJetV2::PassesCuts(AliVEvent* event)
     return kTRUE;
 }
 //_____________________________________________________________________________
-/*Bool_t AliAnalysisTaskJetV2::PassesCuts(Int_t year) 
-{
-    // additional centrality cut based on relation between tpc and global multiplicity
-    if(fDebug > 0) printf("__FILE__ = %s \n __LINE __ %i , __FUNC__ %s \n ", __FILE__, __LINE__, __func__);
-    AliAODEvent* event(dynamic_cast<AliAODEvent*>(InputEvent()));
-    if(!event) return kFALSE;
-    Int_t multTPC(0), multGlob(0), nTracks(InputEvent()->GetNumberOfTracks());
-    for(Int_t iTracks = 0; iTracks < nTracks; iTracks++) { 
-        AliAODTrack* track = event->GetTrack(iTracks);
-        if(!track) continue;
-        if (!track || track->Pt() < .2 || track->Pt() > 5.0 || TMath::Abs(track->Eta()) > .8 || track->GetTPCNcls() < 70 || !track->GetDetPid() || track->GetDetPid()->GetTPCsignal() < 10.0)  continue;  // general quality cut
-        if (track->TestFilterBit(1) && track->Chi2perNDF() > 0.2) multTPC++;
-        if (!track->TestFilterBit(16) || track->Chi2perNDF() < 0.1) continue;
-        Double_t b[2] = {-99., -99.};
-        Double_t bCov[3] = {-99., -99., -99.};
-        AliAODTrack copy(*track);
-        if (copy.PropagateToDCA(event->GetPrimaryVertex(), event->GetMagneticField(), 100., b, bCov) && TMath::Abs(b[0]) < 0.3 && TMath::Abs(b[1]) < 0.3) multGlob++;
-    }
-    if(year == 2010 && multTPC > (-40.3+1.22*multGlob) && multTPC < (32.1+1.59*multGlob)) return kTRUE;
-    if(year == 2011  && multTPC > (-36.73 + 1.48*multGlob) && multTPC < (62.87 + 1.78*multGlob)) return kTRUE;
-    return kFALSE;
-}*/
-//_____________________________________________________________________________
 void AliAnalysisTaskJetV2::FillHistogramsAfterSubtraction(Double_t psi2, Double_t vzero[2][2], Double_t* vzeroComb, Double_t* tpc)
 {
     // fill histograms 
@@ -1653,7 +1799,8 @@ void AliAnalysisTaskJetV2::FillJetHistograms(Double_t psi2)
 void AliAnalysisTaskJetV2::FillQAHistograms(AliVTrack* vtrack) const
 {
     // fill qa histograms for pico tracks
-    if(!vtrack) return;
+     if(fDebug > 1) printf("__FILE__ = %s \n __LINE __ %i , __FUNC__ %s \n ", __FILE__, __LINE__, __func__);
+   if(!vtrack) return;
     AliPicoTrack* track = static_cast<AliPicoTrack*>(vtrack);
     fHistRunnumbersPhi->Fill(fMappedRunNumber, track->Phi());
     fHistRunnumbersEta->Fill(fMappedRunNumber, track->Eta());
@@ -1675,13 +1822,15 @@ void AliAnalysisTaskJetV2::FillQAHistograms(AliVTrack* vtrack) const
 void AliAnalysisTaskJetV2::FillQAHistograms(AliVEvent* vevent) 
 {
     // fill qa histograms for events
+    if(fDebug > 0) printf("__FILE__ = %s \n __LINE __ %i , __FUNC__ %s \n ", __FILE__, __LINE__, __func__);
     if(!vevent) return;
     fHistVertexz->Fill(vevent->GetPrimaryVertex()->GetZ());
     fHistCentrality->Fill(fCent);
     Int_t runNumber(InputEvent()->GetRunNumber());
-    for(fMappedRunNumber = 0; fExpectedRuns->GetSize()+1; fMappedRunNumber++) {
-        if(fExpectedRuns->At(fMappedRunNumber) == runNumber) break;
+    for(fMappedRunNumber = 0; fMappedRunNumber < fExpectedRuns->GetSize()+1; fMappedRunNumber++) {
+        if(fExpectedRuns->At(fMappedRunNumber) == runNumber) return;
     }
+    printf("\n > TASK %s CANNOT IDENTIFY RUN - CONFIGURATION COULD BE INCORRECT < \n", GetName());
 }
 //_____________________________________________________________________________
 void AliAnalysisTaskJetV2::FillAnalysisSummaryHistogram() const
@@ -1728,8 +1877,6 @@ void AliAnalysisTaskJetV2::FillAnalysisSummaryHistogram() const
     fHistAnalysisSummary->SetBinContent(41, (int)fRebinSwapHistoOnTheFly);
     fHistAnalysisSummary->GetXaxis()->SetBinLabel(42, "fUsePtWeight");
     fHistAnalysisSummary->SetBinContent(42, (int)fUsePtWeight);
-//    fHistAnalysisSummary->GetXaxis()->SetBinLabel(43, "fExplicitOutlierCut");
-//    fHistAnalysisSummary->SetBinContent(43, fExplicitOutlierCut);
     fHistAnalysisSummary->GetXaxis()->SetBinLabel(44, "fSoftTrackMinPt");
     fHistAnalysisSummary->SetBinContent(44, fSoftTrackMinPt);
     fHistAnalysisSummary->GetXaxis()->SetBinLabel(45, "fSoftTrackMaxPt");
@@ -1745,7 +1892,7 @@ void AliAnalysisTaskJetV2::Terminate(Option_t *)
     // terminate
     switch (fRunModeType) {
         case kLocal : {
-        printf("__FILE__ = %s \n __LINE __ %i , __FUNC__ %s \n ", __FILE__, __LINE__, __func__);
+        if(fDebug > 0) printf("__FILE__ = %s \n __LINE __ %i , __FUNC__ %s \n ", __FILE__, __LINE__, __func__);
         AliAnalysisTaskJetV2::Dump();
         for(Int_t i(0); i < fHistAnalysisSummary->GetXaxis()->GetNbins(); i++) printf( " > flag: %s \t content %.2f \n", fHistAnalysisSummary->GetXaxis()->GetBinLabel(1+i), fHistAnalysisSummary->GetBinContent(1+i));
         } break;
@@ -1756,6 +1903,7 @@ void AliAnalysisTaskJetV2::Terminate(Option_t *)
 void AliAnalysisTaskJetV2::SetModulationFit(TF1* fit) 
 {
     // set modulation fit
+    if(fDebug > 0) printf("__FILE__ = %s \n __LINE __ %i , __FUNC__ %s \n ", __FILE__, __LINE__, __func__);
     if (fFitModulation) delete fFitModulation;
     fFitModulation = fit; 
 }
@@ -1763,6 +1911,7 @@ void AliAnalysisTaskJetV2::SetModulationFit(TF1* fit)
 void AliAnalysisTaskJetV2::SetUseControlFit(Bool_t c)
 {
     // set control fit
+    if(fDebug > 0) printf("__FILE__ = %s \n __LINE __ %i , __FUNC__ %s \n ", __FILE__, __LINE__, __func__);
     if (fFitControl) delete fFitControl;
     if (c) {
         fFitControl = new TF1("controlFit", "pol0", 0, TMath::TwoPi());
@@ -1773,6 +1922,7 @@ TH1F* AliAnalysisTaskJetV2::GetResolutionFromOuptutFile(detectorType det, Int_t
 {
     // INTERFACE METHOD FOR OUTPUTFILE
     // get the detector resolution, user has ownership of the returned histogram
+    if(fDebug > 0) printf("__FILE__ = %s \n __LINE __ %i , __FUNC__ %s \n ", __FILE__, __LINE__, __func__);
     if(!fOutputList) {
         printf(" > Please add fOutputList first < \n");
         return 0x0;
@@ -1820,6 +1970,7 @@ TH1F* AliAnalysisTaskJetV2::CorrectForResolutionDiff(TH1F* v, detectorType det,
 {
     // INTERFACE METHOD FOR OUTPUT FILE
     // correct the supplied differential vn histogram v for detector resolution
+    if(fDebug > 0) printf("__FILE__ = %s \n __LINE __ %i , __FUNC__ %s \n ", __FILE__, __LINE__, __func__);
     TH1F* r(GetResolutionFromOuptutFile(det, h, cen));
     if(!r) {
         printf(" > Couldn't find resolution < \n");
@@ -1837,6 +1988,7 @@ TH1F* AliAnalysisTaskJetV2::CorrectForResolutionInt(TH1F* v, detectorType det, T
     // INTERFACE METHOD FOR OUTPUT FILE
     // correct the supplied intetrated vn histogram v for detector resolution
     // integrated vn must have the same centrality binning as the resolotion correction
+    if(fDebug > 0) printf("__FILE__ = %s \n __LINE __ %i , __FUNC__ %s \n ", __FILE__, __LINE__, __func__);
     TH1F* r(GetResolutionFromOuptutFile(det, h, cen));
     v->Divide(v, r);
     return v;
@@ -1845,6 +1997,7 @@ TH1F* AliAnalysisTaskJetV2::CorrectForResolutionInt(TH1F* v, detectorType det, T
 TH1F* AliAnalysisTaskJetV2::GetDifferentialQC(TProfile* refCumulants, TProfile* diffCumlants, TArrayD* ptBins, Int_t h)
 {
     // get differential QC
+    if(fDebug > 0) printf("__FILE__ = %s \n __LINE __ %i , __FUNC__ %s \n ", __FILE__, __LINE__, __func__);
     Double_t r(refCumulants->GetBinContent(h-1)); // v2 reference flow
     if(r > 0) r = TMath::Sqrt(r);
     TH1F* qc = new TH1F(Form("QC2v%i", h), Form("QC2v%i", h), ptBins->GetSize()-1, ptBins->GetArray());
@@ -1862,3 +2015,166 @@ TH1F* AliAnalysisTaskJetV2::GetDifferentialQC(TProfile* refCumulants, TProfile*
 }
 
 //_____________________________________________________________________________
+void AliAnalysisTaskJetV2::ReadVZEROCalibration2010h()
+{
+    // necessary for calibration of 10h vzero event plane. code copied from flow package 
+    // (duplicate, but i didn't want to introduce an ulgy dependency )
+    // this function is only called when the runnumber changes 
+    if(fDebug > 0) printf("__FILE__ = %s \n __LINE __ %i , __FUNC__ %s \n ", __FILE__, __LINE__, __func__);
+
+    // 1) check if the proper chi weights for merging vzero a and vzero c ep are present
+    // if not, use sane defaults. centrality binning is equal to that given in the fVZEROcentralityBin snippet
+    //
+    // chi values can be calculated using the static helper function 
+    // AliAnalysisTaskJetV2::CalculateEventPlaneChi(Double_t res) where res is the event plane
+    // resolution in a given centrality bin
+    //
+    // the resolutions that were used for these defaults are
+    // this might need a bit of updating as they were read 'by-eye' from a performance plot ..
+    // Double_t R2VZEROA[] = {.35, .40, .48, .50, .48, .45, .38, .26, .16};
+    // Double_t R2VZEROC[] = {.45, .60, .70, .73, .68, .60, .40, .36, .17};
+    // Double_t R3VZEROA[] = {.22, .23, .22, .19, .15, .12, .08, .00, .00};
+    // Double_t R3VZEROC[] = {.30, .30, .28, .25, .22, .17, .11, .00, .00};
+
+    Double_t chiC2[] = {0.771423, 1.10236, 1.38116, 1.48077, 1.31964, 1.10236, 0.674622, 0.600403, 0.273865};
+    Double_t chiA2[] = {0.582214, 0.674622, 0.832214, 0.873962, 0.832214, 0.771423, 0.637146, 0.424255, 0.257385};
+    Double_t chiC3[] = {0.493347, 0.493347, 0.458557, 0.407166, 0.356628, 0.273865, 0.176208, 6.10352e-05, 6.10352e-05};
+    Double_t chiA3[] = {0.356628, 0.373474, 0.356628, 0.306702, 0.24115, 0.192322, 0.127869, 6.10352e-05, 6.10352e-05};
+
+    if(!fChi2A) fChi2A = new TArrayD(9, chiA2);
+    if(!fChi2C) fChi2C = new TArrayD(9, chiC2);
+    if(!fChi3A) fChi3A = new TArrayD(9, chiA3);
+    if(!fChi3C) fChi3C = new TArrayD(9, chiC3);
+
+    // 2) open database file
+    fOADB = TFile::Open("$ALICE_ROOT/OADB/PWGCF/VZERO/VZEROcalibEP.root");
+    if(fOADB->IsZombie()){
+       printf("OADB file $ALICE_ROOT/OADB/PWGCF/VZERO/VZEROcalibEP.root cannot be opened, CALIBRATION FAILED !");
+       return;
+    }
+
+    AliOADBContainer *cont = (AliOADBContainer*) fOADB->Get("hMultV0BefCorr");
+    if(!cont){
+        // see if database is readable
+       printf("OADB object hMultV0BefCorr is not available in the file\n");
+       return; 
+    }
+    Int_t run(fRunNumber);
+    if(!(cont->GetObject(run))){
+        // if the run isn't recognized fall back to a default run
+       printf("OADB object hMultVZEROBefCorr is not available for run %i (used default run 137366)\n",run);
+       run = 137366;
+    }
+    // step 3) get the proper multiplicity weights from the vzero signal
+    fVZEROgainEqualization = ((TH2F*)cont->GetObject(run))->ProfileX();
+    if(!fVZEROgainEqualization) {
+        AliFatal(Form("%s: Fatal error, couldn't read fVZEROgainEqualization from OADB object < \n", GetName()));
+        return;
+    }
+
+    TF1* fpol0 = new TF1("fpol0","pol0"); 
+    if(fVZEROgainEqualizationPerRing) {
+        // do the calibration per ring
+        // start with the vzero c rings (segments 0 through 31)
+        fVZEROgainEqualization->Fit(fpol0, "", "", 0, 8);
+        (fUseVZERORing[0]) ? SetVZEROCpol(0, fpol0->GetParameter(0)) : SetVZEROCpol(0, 0.);
+        fVZEROgainEqualization->Fit(fpol0, "", "", 8, 16);
+        (fUseVZERORing[1]) ? SetVZEROCpol(1, fpol0->GetParameter(0)) : SetVZEROCpol(1, 0.);
+        fVZEROgainEqualization->Fit(fpol0, "", "", 16, 24);
+        (fUseVZERORing[2]) ? SetVZEROCpol(2, fpol0->GetParameter(0)) : SetVZEROCpol(2, 0.);
+        fVZEROgainEqualization->Fit(fpol0, "", "", 24, 32);
+        (fUseVZERORing[3]) ? SetVZEROCpol(3, fpol0->GetParameter(0)) : SetVZEROCpol(3, 0.);
+        // same thing for vero A
+        fVZEROgainEqualization->Fit(fpol0, "", "", 32, 40);
+        (fUseVZERORing[4]) ? SetVZEROApol(0, fpol0->GetParameter(0)) : SetVZEROApol(0, 0.);
+        fVZEROgainEqualization->Fit(fpol0, "", "", 40, 48);
+        (fUseVZERORing[5]) ? SetVZEROApol(1, fpol0->GetParameter(0)) : SetVZEROApol(1, 0.);
+        fVZEROgainEqualization->Fit(fpol0, "", "", 48, 56);
+        (fUseVZERORing[6]) ? SetVZEROApol(2, fpol0->GetParameter(0)) : SetVZEROApol(2, 0.);
+        fVZEROgainEqualization->Fit(fpol0, "", "", 56, 64);
+        (fUseVZERORing[7]) ? SetVZEROApol(3, fpol0->GetParameter(0)) : SetVZEROApol(3, 0.);
+    } else {
+        // do the calibration in one go. the calibration will still be 
+        // stored per ring, but each ring has the same weight now
+        // this should be the default for the analysis as the database is tuned to this configuration
+       fVZEROgainEqualization->Fit(fpol0,"","",0,31);
+       for(Int_t i(0); i < 4; i++) SetVZEROCpol(i, fpol0->GetParameter(0));
+       fVZEROgainEqualization->Fit(fpol0,"","",32,64);
+       for(Int_t i(0); i < 4; i++) SetVZEROApol(i, fpol0->GetParameter(0));
+    }
+
+    // step 4) extract the information to re-weight the q-vectors 
+    for(Int_t iside=0;iside<2;iside++){
+       for(Int_t icoord=0;icoord<2;icoord++){
+           for(Int_t i=0;i  < 9;i++){
+               char namecont[100];
+               if(iside==0 && icoord==0)
+                 snprintf(namecont,100,"hQxc2_%i",i);
+               else if(iside==1 && icoord==0)
+                 snprintf(namecont,100,"hQxa2_%i",i);
+               else if(iside==0 && icoord==1)
+                 snprintf(namecont,100,"hQyc2_%i",i);
+               else if(iside==1 && icoord==1)
+                 snprintf(namecont,100,"hQya2_%i",i);
+
+               cont = (AliOADBContainer*) fOADB->Get(namecont);
+               if(!cont){
+                   printf("OADB object %s is not available in the file\n",namecont);
+                   return;     
+               }
+       
+               if(!(cont->GetObject(run))){
+                   printf("OADB object %s is not available for run %i (used run 137366)\n",namecont,run);
+                   run = 137366;
+               }
+
+                // store info for all centralities to cache
+                fMeanQ[i][iside][icoord] = ((TH1F *) cont->GetObject(run))->GetMean();
+               fWidthQ[i][iside][icoord] = ((TH1F *) cont->GetObject(run))->GetRMS();
+
+               //for v3
+               if(iside==0 && icoord==0)
+                 snprintf(namecont,100,"hQxc3_%i",i);
+               else if(iside==1 && icoord==0)
+                 snprintf(namecont,100,"hQxa3_%i",i);
+               else if(iside==0 && icoord==1)
+                 snprintf(namecont,100,"hQyc3_%i",i);
+               else if(iside==1 && icoord==1)
+                 snprintf(namecont,100,"hQya3_%i",i);
+
+               cont = (AliOADBContainer*) fOADB->Get(namecont);
+               if(!cont){
+                   printf("OADB object %s is not available in the file\n",namecont);
+                   return;     
+               }
+               
+               if(!(cont->GetObject(run))){
+                   printf("OADB object %s is not available for run %i (used run 137366)\n",namecont,run);
+                   run = 137366;
+               }
+                // store info for all centralities to cache
+               fMeanQv3[i][iside][icoord] = ((TH1F *) cont->GetObject(run))->GetMean();
+               fWidthQv3[i][iside][icoord] = ((TH1F *) cont->GetObject(run))->GetRMS();
+           }
+       }
+    }
+    // cleanup. the opened file is closed in the destructor, otherwise fVZEROgainEqualization is no longer available
+    delete fpol0;
+}
+//_____________________________________________________________________________
+Int_t AliAnalysisTaskJetV2::GetVZEROCentralityBin() const
+{
+    // return cache index number corresponding to the event centrality
+    if(fDebug > 0) printf("__FILE__ = %s \n __LINE __ %i , __FUNC__ %s \n ", __FILE__, __LINE__, __func__);
+    Float_t v0Centr(InputEvent()->GetCentrality()->GetCentralityPercentile("V0M"));
+    if(v0Centr < 5) return 0;
+    else if(v0Centr < 10) return 1;
+    else if(v0Centr < 20) return  2;
+    else if(v0Centr < 30) return  3;
+    else if(v0Centr < 40) return  4;
+    else if(v0Centr < 50) return  5;
+    else if(v0Centr < 60) return  6;
+    else if(v0Centr < 70) return  7;
+    else return 8;
+}
+//_____________________________________________________________________________
index 59de265c17fa548205fe4de71713fcbca203be09..fc7f77a01117b60ac58f21d86b411bdc4f051657 100644 (file)
@@ -17,6 +17,7 @@
 #include <AliJetContainer.h>
 #include <AliParticleContainer.h>
 
+class TFile;
 class TF1;
 class THF1;
 class THF2;
@@ -48,11 +49,11 @@ class AliAnalysisTaskJetV2 : public AliAnalysisTaskEmcalJet {
         virtual Bool_t          Run();
         TH1F*                   BookTH1F(const char* name, const char* x, Int_t bins, Double_t min, Double_t max, Int_t c = -1, Bool_t append = kTRUE);
         TH2F*                   BookTH2F(const char* name, const char* x, const char* y, Int_t binsx, Double_t minx, Double_t maxx, Int_t binsy, Double_t miny, Double_t maxy, Int_t c = -1, Bool_t append = kTRUE);
-        /* inline */    Double_t PhaseShift(Double_t x) const {  
+        /* inline */    static Double_t PhaseShift(Double_t x) {  
             while (x>=TMath::TwoPi())x-=TMath::TwoPi();
             while (x<0.)x+=TMath::TwoPi();
             return x; }
-        /* inline */    Double_t PhaseShift(Double_t x, Double_t n) const {
+        /* inline */    static Double_t PhaseShift(Double_t x, Double_t n) {
             x = PhaseShift(x);
             if(TMath::Nint(n)==2) while (x>TMath::Pi()) x-=TMath::Pi();
             if(TMath::Nint(n)==3) {
@@ -60,13 +61,13 @@ class AliAnalysisTaskJetV2 : public AliAnalysisTaskEmcalJet {
                 if(x>TMath::TwoPi()/n) x = TMath::TwoPi()-(x+TMath::TwoPi()/n);
             }
             return x; }
-        /* inline */    Double_t ChiSquarePDF(Int_t ndf, Double_t x) const {
+        /* inline */    static Double_t ChiSquarePDF(Int_t ndf, Double_t x) {
             Double_t n(ndf/2.), denom(TMath::Power(2, n)*TMath::Gamma(n));
             if (denom!=0)  return ((1./denom)*TMath::Power(x, n-1)*TMath::Exp(-x/2.)); 
             return -999; }
         // note that the cdf of the chisquare distribution is the normalized lower incomplete gamma function
-        /* inline */    Double_t ChiSquareCDF(Int_t ndf, Double_t x) const { return TMath::Gamma(ndf/2., x/2.); }
-        /* inline */    Double_t ChiSquare(TH1& histo, TF1* func) const {
+        /* inline */    static Double_t ChiSquareCDF(Int_t ndf, Double_t x) { return TMath::Gamma(ndf/2., x/2.); }
+        /* inline */    static Double_t ChiSquare(TH1& histo, TF1* func) {
             // evaluate the chi2 using a poissonian error estimate on bins
             Double_t chi2(0.);
             for(Int_t i(0); i < histo.GetXaxis()->GetNbins(); i++) {
@@ -75,7 +76,7 @@ class AliAnalysisTaskJetV2 : public AliAnalysisTaskEmcalJet {
             }
            return chi2;
         }
-        /* inline*/ Double_t KolmogorovTest(TH1F& histo, TF1* func) const {
+        /* inline */ Double_t KolmogorovTest(TH1F& histo, TF1* func) const {
             // return the probability from a Kolmogorov test
             return .5;
             TH1F test(histo);       // stack copy of test statistic
@@ -83,7 +84,7 @@ class AliAnalysisTaskJetV2 : public AliAnalysisTaskEmcalJet {
             if(fFitGoodnessTest == kKolmogorovTOY) return histo.TH1::KolmogorovTest((&test), "X");
             return histo.TH1::KolmogorovTest((&test));
         }
+
         // setters - analysis setup
         void                    SetDebugMode(Int_t d)                           {fDebug = d;}
         void                    SetRunToyMC(Bool_t t)                           {fRunToyMC = t; }
@@ -123,9 +124,23 @@ class AliAnalysisTaskJetV2 : public AliAnalysisTaskEmcalJet {
         void                    SetExcludeLeadingJetsFromFit(Float_t n)         {fExcludeLeadingJetsFromFit = n; }
         void                    SetRebinSwapHistoOnTheFly(Bool_t r)             {fRebinSwapHistoOnTheFly = r; }
         void                    SetSaveThisPercentageOfFits(Float_t p)          {fPercentageOfFits = p; }
-        void                    SetUseV0EventPlaneFromHeader(Bool_t h)          {fUseV0EventPlaneFromHeader = h;}
-//        void                    SetExplicitOutlierCutForYear(Int_t y)           {fExplicitOutlierCut = y;}
-        // getters - these are used as well by AliAnalyisTaskJetFlow, so be careful when changing them
+        // setters specific to the vzero calibration for 10h data        
+        void                    SetVZEROApol(Int_t ring, Float_t f)             {fVZEROApol[ring]=f;}
+        void                    SetVZEROCpol(Int_t ring, Float_t f)             {fVZEROCpol[ring]=f;}
+        void                    SetVZEROgainEqualizationPerRing(Bool_t s)       {fVZEROgainEqualizationPerRing = s;}
+        void                    SetUseVZERORing(Int_t i, Bool_t u) {
+            // exclude vzero rings: 0 through 7 can be excluded by calling this setter multiple times
+            // 0 corresponds to segment ID 0 through 7, etc
+            fUseVZERORing[i] = u;
+            fVZEROgainEqualizationPerRing = kTRUE;       // must be true for this option
+        }
+
+        void                    SetChi2VZEROA(TArrayD* a)                       { fChi2A = a;}
+        void                    SetChi2VZEROC(TArrayD* a)                       { fChi2C = a;}
+        void                    SetChi3VZEROA(TArrayD* a)                       { fChi3A = a;}
+        void                    SetChi3VZEROC(TArrayD* a)                       { fChi3C = a;}
+
+        // getters 
         TString                 GetJetsName() const                             {return GetJetContainer()->GetArrayName(); }
         TString                 GetTracksName() const                           {return GetParticleContainer()->GetArrayName(); }
         TString                 GetLocalRhoName() const                         {return fLocalRhoName; }
@@ -152,16 +167,19 @@ class AliAnalysisTaskJetV2 : public AliAnalysisTaskEmcalJet {
             return leadingJet;
         }
         void                    ExecMe()                                {ExecOnce();}
-        AliAnalysisTaskJetV2* ReturnMe()                                {return this;}
+        AliAnalysisTaskJetV2*   ReturnMe()                              {return this;}
         // local cuts
         void                    SetSoftTrackMinMaxPt(Float_t min, Float_t max)          {fSoftTrackMinPt = min; fSoftTrackMaxPt = max;}
         void                    SetSemiGoodJetMinMaxPhi(Double_t a, Double_t b)         {fSemiGoodJetMinPhi = a; fSemiGoodJetMaxPhi = b;}
         void                    SetSemiGoodTrackMinMaxPhi(Double_t a, Double_t b)       {fSemiGoodTrackMinPhi = a; fSemiGoodTrackMaxPhi = b;}
         // numerical evaluations
+        static Double_t         CalculateEventPlaneChi(Double_t res);
         void                    CalculateEventPlaneVZERO(Double_t vzero[2][2]) const;
-        void                    CalculateEventPlaneTPC(Double_t* tpc);
         void                    CalculateEventPlaneCombinedVZERO(Double_t* comb) const;
+        void                    CalculateEventPlaneTPC(Double_t* tpc);
         void                    CalculateEventPlaneResolution(Double_t vzero[2][2], Double_t* vzeroComb, Double_t* tpc);
+        void                    CalculateQvectorVZERO(Double_t Qa2[2], Double_t Qc2[2], Double_t Qa3[2], Double_t Qc3[2]) const;
+        void                    CalculateQvectorCombinedVZERO(Double_t Q2[2], Double_t Q3[2]) const;
         void                    CalculateRandomCone(
                 Float_t &pt, 
                 Float_t &eta, 
@@ -172,7 +190,7 @@ class AliAnalysisTaskJetV2 : public AliAnalysisTaskEmcalJet {
                 ) const;
         Double_t                CalculateQC2(Int_t harm);
         Double_t                CalculateQC4(Int_t harm);
-        // helper calculations for the q-cumulant analysis, also used by AliAnalyisTaskJetFlow
+        // helper calculations for the q-cumulant analysis
         void                    QCnQnk(Int_t n, Int_t k, Double_t &reQ, Double_t &imQ);
         void                    QCnDiffentialFlowVectors(
             TClonesArray* pois, TArrayD* ptBins, Bool_t vpart, Double_t* repn, Double_t* impn, 
@@ -184,7 +202,7 @@ class AliAnalysisTaskJetV2 : public AliAnalysisTaskEmcalJet {
         Bool_t                  QCnRecovery(Double_t psi2, Double_t psi3);
         // analysis details
         Bool_t                  CorrectRho(Double_t psi2, Double_t psi3);
-        // event and track selection, also used by AliAnalyisTaskJetFlow
+        // event and track selection
         /* inline */    Bool_t PassesCuts(AliVParticle* track) const    { return AcceptTrack(track, 0); }
         /* inline */    Bool_t PassesCuts(AliEmcalJet* jet)             { return AcceptJet(jet, 0); }
         /* inline */    Bool_t PassesCuts(AliVCluster* clus) const      { return AcceptCluster(clus, 0); }
@@ -194,7 +212,6 @@ class AliAnalysisTaskJetV2 : public AliAnalysisTaskEmcalJet {
             return (jet && jet->Pt() > 1. && jet->Eta() > minEta && jet->Eta() < maxEta && jet->Phi() > minPhi && jet->Phi() < maxPhi && jet->Area() > .557*GetJetRadius()*GetJetRadius()*TMath::Pi());
         }
         Bool_t                  PassesCuts(AliVEvent* event);
-        /*Bool_t                  PassesCuts(Int_t year);*/
         Bool_t                  PassesCuts(const AliVCluster* track) const;
         // filling histograms
         void                    FillHistogramsAfterSubtraction(Double_t psi2, Double_t vzero[2][2], Double_t* vzeroComb, Double_t* tpc);
@@ -214,6 +231,8 @@ class AliAnalysisTaskJetV2 : public AliAnalysisTaskEmcalJet {
         TH1F*                   CorrectForResolutionDiff(TH1F* v, detectorType detector, TArrayD* cen, Int_t c, Int_t h = 2);
         TH1F*                   CorrectForResolutionInt(TH1F* v, detectorType detector, TArrayD* cen, Int_t h = 2);
         TH1F*                   GetDifferentialQC(TProfile* refCumulants, TProfile* diffCumlants, TArrayD* ptBins, Int_t h);
+        void                    ReadVZEROCalibration2010h();
+        Int_t                   GetVZEROCentralityBin() const;
     private:
         // analysis flags and settings
         Int_t                   fDebug;                 // debug level (0 none, 1 fcn calls, 2 verbose)
@@ -293,8 +312,6 @@ class AliAnalysisTaskJetV2 : public AliAnalysisTaskEmcalJet {
         Float_t                 fExcludeLeadingJetsFromFit;    // exclude n leading jets from fit
         Bool_t                  fRebinSwapHistoOnTheFly;       // rebin swap histo on the fly
         Float_t                 fPercentageOfFits;      // save this percentage of fits
-        Bool_t                  fUseV0EventPlaneFromHeader;    // use the vzero event plane from the header
-//        Int_t                   fExplicitOutlierCut;    // cut on correlation of tpc and global multiplicity
         // transient object pointers
         TList*                  fOutputList;            //! output list
         TList*                  fOutputListGood;        //! output list for local analysis
@@ -361,11 +378,27 @@ class AliAnalysisTaskJetV2 : public AliAnalysisTaskEmcalJet {
         // in plane, out of plane jet spectra
         TH2F*                   fHistJetPsi2Pt[10];             //! event plane dependence of jet pt
         TH2F*                   fHistJetPsi2PtRho0[10];         //! event plane dependence of jet pt vs rho_0
+        // vzero event plane calibration cache for 10h data
+        Float_t                 fMeanQ[9][2][2];                //! recentering
+        Float_t                 fWidthQ[9][2][2];               //! recentering
+        Float_t                 fMeanQv3[9][2][2];              //! recentering
+        Float_t                 fWidthQv3[9][2][2];             //! recentering
+        TH1*                    fVZEROgainEqualization;         //! equalization histo
+        Bool_t                  fVZEROgainEqualizationPerRing;  // per ring vzero gain calibration
+        Float_t                 fVZEROApol[4];                  //! calibration info per ring
+        Float_t                 fVZEROCpol[4];                  //! calibration info per ring
+        Bool_t                  fUseVZERORing[8];               // kTRUE means the ring is included
+        TArrayD*                fChi2A;                         // chi vs cent for vzero A ep_2
+        TArrayD*                fChi2C;                         // chi vs cent for vzero C ep_2
+        TArrayD*                fChi3A;                         // chi vs cent for vzero A ep_3
+        TArrayD*                fChi3C;                         // chi vs cent for vzero C ep_3
+        TFile*                  fOADB;                          //! fOADB
+
 
         AliAnalysisTaskJetV2(const AliAnalysisTaskJetV2&);                  // not implemented
         AliAnalysisTaskJetV2& operator=(const AliAnalysisTaskJetV2&);       // not implemented
 
-        ClassDef(AliAnalysisTaskJetV2, 2);
+        ClassDef(AliAnalysisTaskJetV2, 3);
 };
 
 #endif