]> git.uio.no Git - u/mrichter/AliRoot.git/blobdiff - JETAN/DEV/AliAnalysisTaskJetCluster.cxx
adding requirement on tzero or vzero if needed
[u/mrichter/AliRoot.git] / JETAN / DEV / AliAnalysisTaskJetCluster.cxx
index 28e721eefa1d27191f97b641e7f7c9860dd6991c..e83e0f7133df850432fa33cacfe22308a20bb448 100644 (file)
@@ -66,6 +66,7 @@
 // get info on how fastjet was configured
 #include "fastjet/config.h"
 
+using std::vector;
 
 ClassImp(AliAnalysisTaskJetCluster)
 
@@ -102,7 +103,10 @@ AliAnalysisTaskJetCluster::AliAnalysisTaskJetCluster():
   fUseAODMCInput(kFALSE),
   fUseBackgroundCalc(kFALSE),
   fEventSelection(kFALSE),     
+  fRequireVZEROAC(kFALSE),     
+  fRequireTZEROvtx(kFALSE),
   fFilterMask(0),
+  fFilterMaskBestPt(0),
   fFilterType(0),
   fJetTypes(kJet),
   fTrackTypeRec(kTrackUndef),
@@ -122,6 +126,7 @@ AliAnalysisTaskJetCluster::AliAnalysisTaskJetCluster():
   fVtxR2Cut(1),
   fCentCutUp(0),
   fCentCutLo(0),
+  fStoreRhoLeadingTrackCorr(kFALSE),
   fNonStdBranch(""),
   fBackgroundBranch(""),
   fNonStdFile(""),
@@ -134,8 +139,15 @@ AliAnalysisTaskJetCluster::AliAnalysisTaskJetCluster():
   fhEffH1(0x0),
   fhEffH2(0x0),
   fhEffH3(0x0),
-  fUseTrMomentumSmearing(kFALSE),
+  fUseTrPtResolutionSmearing(kFALSE),
   fUseDiceEfficiency(kFALSE),
+  fDiceEfficiencyMinPt(-1.),
+  fUseTrPtResolutionFromOADB(kFALSE),
+  fUseTrEfficiencyFromOADB(kFALSE),
+  fPathTrPtResolution(""),
+  fPathTrEfficiency(""),
+  fChangeEfficiencyFraction(0.),
+  fEfficiencyFixed(1.),
   fRparam(1.0), 
   fAlgorithm(fastjet::kt_algorithm),
   fStrategy(fastjet::Best),
@@ -204,6 +216,30 @@ AliAnalysisTaskJetCluster::AliAnalysisTaskJetCluster():
   fh2PtNchNRan(0x0),
   fh2TracksLeadingJetPhiPtRan(0x0),
   fh2TracksLeadingJetPhiPtWRan(0x0),
+  fh3CentvsRhoLeadingTrackPt(0x0),
+  fh3CentvsSigmaLeadingTrackPt(0x0),
+  fh3MultvsRhoLeadingTrackPt(0x0),
+  fh3MultvsSigmaLeadingTrackPt(0x0),
+  fh3CentvsRhoLeadingTrackPtQ1(0x0),
+  fh3CentvsRhoLeadingTrackPtQ2(0x0),
+  fh3CentvsRhoLeadingTrackPtQ3(0x0),
+  fh3CentvsRhoLeadingTrackPtQ4(0x0),
+  fh3CentvsSigmaLeadingTrackPtQ1(0x0),
+  fh3CentvsSigmaLeadingTrackPtQ2(0x0),
+  fh3CentvsSigmaLeadingTrackPtQ3(0x0),
+  fh3CentvsSigmaLeadingTrackPtQ4(0x0),
+  fh3MultvsRhoLeadingTrackPtQ1(0x0),
+  fh3MultvsRhoLeadingTrackPtQ2(0x0),
+  fh3MultvsRhoLeadingTrackPtQ3(0x0),
+  fh3MultvsRhoLeadingTrackPtQ4(0x0),
+  fh3MultvsSigmaLeadingTrackPtQ1(0x0),
+  fh3MultvsSigmaLeadingTrackPtQ2(0x0),
+  fh3MultvsSigmaLeadingTrackPtQ3(0x0),
+  fh3MultvsSigmaLeadingTrackPtQ4(0x0),
+  fh3CentvsDeltaRhoLeadingTrackPtQ1(0x0),
+  fh3CentvsDeltaRhoLeadingTrackPtQ2(0x0),
+  fh3CentvsDeltaRhoLeadingTrackPtQ3(0x0),
+  fh3CentvsDeltaRhoLeadingTrackPtQ4(0x0),
   fh2PtGenPtSmeared(0x0),
   fp1Efficiency(0x0),
   fp1PtResolution(0x0),
@@ -233,8 +269,11 @@ AliAnalysisTaskJetCluster::AliAnalysisTaskJetCluster(const char* name):
   fUseAODTrackInput(kFALSE),
   fUseAODMCInput(kFALSE),
   fUseBackgroundCalc(kFALSE),
-  fEventSelection(kFALSE),                                                     
+  fEventSelection(kFALSE),
+  fRequireVZEROAC(kFALSE),     
+  fRequireTZEROvtx(kFALSE), 
   fFilterMask(0),
+  fFilterMaskBestPt(0),
   fFilterType(0),
   fJetTypes(kJet),
   fTrackTypeRec(kTrackUndef),
@@ -254,6 +293,7 @@ AliAnalysisTaskJetCluster::AliAnalysisTaskJetCluster(const char* name):
   fVtxR2Cut(1),
   fCentCutUp(0),
   fCentCutLo(0),
+  fStoreRhoLeadingTrackCorr(kFALSE),
   fNonStdBranch(""),
   fBackgroundBranch(""),
   fNonStdFile(""),
@@ -266,8 +306,15 @@ AliAnalysisTaskJetCluster::AliAnalysisTaskJetCluster(const char* name):
   fhEffH1(0x0),
   fhEffH2(0x0),
   fhEffH3(0x0),
-  fUseTrMomentumSmearing(kFALSE),
+  fUseTrPtResolutionSmearing(kFALSE),
   fUseDiceEfficiency(kFALSE),
+  fDiceEfficiencyMinPt(-1.),
+  fUseTrPtResolutionFromOADB(kFALSE),
+  fUseTrEfficiencyFromOADB(kFALSE),
+  fPathTrPtResolution(""),
+  fPathTrEfficiency(""),
+  fChangeEfficiencyFraction(0.),
+  fEfficiencyFixed(1.),
   fRparam(1.0), 
   fAlgorithm(fastjet::kt_algorithm),
   fStrategy(fastjet::Best),
@@ -336,6 +383,30 @@ AliAnalysisTaskJetCluster::AliAnalysisTaskJetCluster(const char* name):
   fh2PtNchNRan(0x0),
   fh2TracksLeadingJetPhiPtRan(0x0),
   fh2TracksLeadingJetPhiPtWRan(0x0),
+  fh3CentvsRhoLeadingTrackPt(0x0),
+  fh3CentvsSigmaLeadingTrackPt(0x0),
+  fh3MultvsRhoLeadingTrackPt(0x0),
+  fh3MultvsSigmaLeadingTrackPt(0x0),
+  fh3CentvsRhoLeadingTrackPtQ1(0x0),
+  fh3CentvsRhoLeadingTrackPtQ2(0x0),
+  fh3CentvsRhoLeadingTrackPtQ3(0x0),
+  fh3CentvsRhoLeadingTrackPtQ4(0x0),
+  fh3CentvsSigmaLeadingTrackPtQ1(0x0),
+  fh3CentvsSigmaLeadingTrackPtQ2(0x0),
+  fh3CentvsSigmaLeadingTrackPtQ3(0x0),
+  fh3CentvsSigmaLeadingTrackPtQ4(0x0),
+  fh3MultvsRhoLeadingTrackPtQ1(0x0),
+  fh3MultvsRhoLeadingTrackPtQ2(0x0),
+  fh3MultvsRhoLeadingTrackPtQ3(0x0),
+  fh3MultvsRhoLeadingTrackPtQ4(0x0),
+  fh3MultvsSigmaLeadingTrackPtQ1(0x0),
+  fh3MultvsSigmaLeadingTrackPtQ2(0x0),
+  fh3MultvsSigmaLeadingTrackPtQ3(0x0),
+  fh3MultvsSigmaLeadingTrackPtQ4(0x0),
+  fh3CentvsDeltaRhoLeadingTrackPtQ1(0x0),
+  fh3CentvsDeltaRhoLeadingTrackPtQ2(0x0),
+  fh3CentvsDeltaRhoLeadingTrackPtQ3(0x0),
+  fh3CentvsDeltaRhoLeadingTrackPtQ4(0x0),
   fh2PtGenPtSmeared(0x0),
   fp1Efficiency(0x0),
   fp1PtResolution(0x0),
@@ -402,8 +473,12 @@ void AliAnalysisTaskJetCluster::UserCreateOutputObjects()
       if(fJetTypes&kJetRan){
        fTCAJetsOutRan = new TClonesArray("AliAODJet", 0);
        fTCAJetsOutRan->SetName(Form("%s_%s",fNonStdBranch.Data(),"random"));
-       if(fUseDiceEfficiency || fUseTrMomentumSmearing)
-         fTCAJetsOutRan->SetName(Form("%s_%sDetector%d%d",fNonStdBranch.Data(),"random",fUseTrMomentumSmearing,fUseDiceEfficiency));
+       if(fUseDiceEfficiency || fUseTrPtResolutionSmearing) {
+         if(  fEfficiencyFixed < 1.)
+           fTCAJetsOutRan->SetName(Form("%s_%sDetector%d%dEffFixed%d",fNonStdBranch.Data(),"random",fUseTrPtResolutionSmearing,fUseDiceEfficiency,(int)(fEfficiencyFixed*100.)));
+         else
+           fTCAJetsOutRan->SetName(Form("%s_%sDetector%d%dFr%d",fNonStdBranch.Data(),"random",fUseTrPtResolutionSmearing,fUseDiceEfficiency,(int)(fChangeEfficiencyFraction*100.)));
+       }
        AddAODBranch("TClonesArray",&fTCAJetsOutRan,fNonStdFile.Data());
       }
 
@@ -411,17 +486,23 @@ void AliAnalysisTaskJetCluster::UserCreateOutputObjects()
        if(!AODEvent()->FindListObject(Form("%s_%s",AliAODJetEventBackground::StdBranchName(),fNonStdBranch.Data()))){
          fAODJetBackgroundOut = new AliAODJetEventBackground();
          fAODJetBackgroundOut->SetName(Form("%s_%s",AliAODJetEventBackground::StdBranchName(),fNonStdBranch.Data()));
-         if(fUseDiceEfficiency || fUseTrMomentumSmearing)
-           fAODJetBackgroundOut->SetName(Form("%s_%sDetector%d%d",AliAODJetEventBackground::StdBranchName(),fNonStdBranch.Data(),fUseTrMomentumSmearing,fUseDiceEfficiency));
-
+         if(fUseDiceEfficiency || fUseTrPtResolutionSmearing) {
+           if(  fEfficiencyFixed < 1.)
+             fAODJetBackgroundOut->SetName(Form("%s_%sDetector%d%dEffFixed%d",AliAODJetEventBackground::StdBranchName(),fNonStdBranch.Data(),fUseTrPtResolutionSmearing,fUseDiceEfficiency,(int)(fEfficiencyFixed*100.)));
+           else
+             fAODJetBackgroundOut->SetName(Form("%s_%sDetector%d%dFr%d",AliAODJetEventBackground::StdBranchName(),fNonStdBranch.Data(),fUseTrPtResolutionSmearing,fUseDiceEfficiency,(int)(fChangeEfficiencyFraction*100.)));
+         }
          AddAODBranch("AliAODJetEventBackground",&fAODJetBackgroundOut,fNonStdFile.Data());  
        }
       }
       // create the branch for the random cones with the same R 
       TString cName = Form("%sRandomConeSkip%02d",fNonStdBranch.Data(),fNSkipLeadingCone);
-      if(fUseDiceEfficiency || fUseTrMomentumSmearing)
-       cName = Form("%sDetector%d%d_RandomConeSkip%02d",fNonStdBranch.Data(),fUseTrMomentumSmearing,fUseDiceEfficiency,fNSkipLeadingCone);
-
+      if(fUseDiceEfficiency || fUseTrPtResolutionSmearing) {
+       if(  fEfficiencyFixed < 1.)
+         cName = Form("%sDetector%d%dEffFixed%d_RandomConeSkip%02d",fNonStdBranch.Data(),fUseTrPtResolutionSmearing,fUseDiceEfficiency,(int)(fEfficiencyFixed*100.),fNSkipLeadingCone);
+       else
+         cName = Form("%sDetector%d%dFr%d_RandomConeSkip%02d",fNonStdBranch.Data(),fUseTrPtResolutionSmearing,fUseDiceEfficiency,(int)(fChangeEfficiencyFraction*100.),fNSkipLeadingCone);
+      }
       if(fNRandomCones>0){
        if(fJetTypes&kRC){
          if(!AODEvent()->FindListObject(cName.Data())){
@@ -451,8 +532,6 @@ void AliAnalysisTaskJetCluster::UserCreateOutputObjects()
       }
     }
 
-  //  FitMomentumResolution();
-
 
   if(!fHistList)fHistList = new TList();
   fHistList->SetOwner();
@@ -598,6 +677,71 @@ void AliAnalysisTaskJetCluster::UserCreateOutputObjects()
   fh2TracksLeadingJetPhiPtWRan = new TH2F("fh2TracksLeadingJetPhiPtWRan","leading p_T vs delta phi to leading jet;#Delta#phi;p_{T} (GeV/c)",
                                       nBinPhi,binLimitsPhi,nBinPt,binLimitsPt);
 
+  if(fStoreRhoLeadingTrackCorr) {
+    fh3CentvsRhoLeadingTrackPt = new TH3F("fh3CentvsRhoLeadingTrackPt","centrality vs background density full event; centrality; #rho", 50,0.,100.,500,0.,250.,100,0.,100.);
+    fh3CentvsSigmaLeadingTrackPt = new TH3F("fh3CentvsSigmaLeadingTrackPt","centrality vs sigma full event; centrality; #sigma(#rho)", 50,0.,100.,50,0.,50.,100,0.,100.);
+    fh3MultvsRhoLeadingTrackPt = new TH3F("fh3MultvsRhoLeadingTrackPt","multiplicity vs background density full event; multiplicity; #rho", 100,0.,5000.,500,0.,250.,100,0.,100.);
+    fh3MultvsSigmaLeadingTrackPt = new TH3F("fh3MultvsSigmaLeadingTrackPt","multiplicity vs sigma full event; multiplicity; #sigma(#rho)", 100,0.,5000.,50,0.,50.,100,0.,100.);
+
+
+    fh3CentvsRhoLeadingTrackPtQ1 = new TH3F("fh3CentvsRhoLeadingTrackPtQ1","centrality vs background density Q1; centrality; #rho; leading p_{t}^{track}", 50,0.,100.,500,0.,250.,100,0.,100.);
+    fh3CentvsRhoLeadingTrackPtQ2 = new TH3F("fh3CentvsRhoLeadingTrackPtQ2","centrality vs background density Q2; centrality; #rho; leading p_{t}^{track}", 50,0.,100.,500,0.,250.,100,0.,100.);
+    fh3CentvsRhoLeadingTrackPtQ3 = new TH3F("fh3CentvsRhoLeadingTrackPtQ3","centrality vs background density Q3; centrality; #rho; leading p_{t}^{track}", 50,0.,100.,500,0.,250.,100,0.,100.);
+    fh3CentvsRhoLeadingTrackPtQ4 = new TH3F("fh3CentvsRhoLeadingTrackPtQ4","centrality vs background density Q4; centrality; #rho; leading p_{t}^{track}", 50,0.,100.,500,0.,250.,100,0.,100.);
+
+    fh3CentvsSigmaLeadingTrackPtQ1 = new TH3F("fh3CentvsSigmaLeadingTrackPtQ1","centrality vs background density Q1; centrality; #sigma(#rho); leading p_{t}^{track}", 50,0.,100.,50,0.,50.,100,0.,100.);
+    fh3CentvsSigmaLeadingTrackPtQ2 = new TH3F("fh3CentvsSigmaLeadingTrackPtQ2","centrality vs background density Q2; centrality; #sigma(#rho); leading p_{t}^{track}", 50,0.,100.,50,0.,50.,100,0.,100.);
+    fh3CentvsSigmaLeadingTrackPtQ3 = new TH3F("fh3CentvsSigmaLeadingTrackPtQ3","centrality vs background density Q3; centrality; #sigma(#rho); leading p_{t}^{track}", 50,0.,100.,50,0.,50.,100,0.,100.);
+    fh3CentvsSigmaLeadingTrackPtQ4 = new TH3F("fh3CentvsSigmaLeadingTrackPtQ4","centrality vs background density Q4; centrality; #sigma(#rho); leading p_{t}^{track}", 50,0.,100.,50,0.,50.,100,0.,100.);
+
+    fh3MultvsRhoLeadingTrackPtQ1 = new TH3F("fh3MultvsRhoLeadingTrackPtQ1","multiplicity vs background density Q1; multiplicity; #rho; leading p_{t}^{track}", 100,0.,5000.,500,0.,250.,100,0.,100.);
+    fh3MultvsRhoLeadingTrackPtQ2 = new TH3F("fh3MultvsRhoLeadingTrackPtQ2","multiplicity vs background density Q2; multiplicity; #rho; leading p_{t}^{track}", 100,0.,5000.,500,0.,250.,100,0.,100.);
+    fh3MultvsRhoLeadingTrackPtQ3 = new TH3F("fh3MultvsRhoLeadingTrackPtQ3","multiplicity vs background density Q3; multiplicity; #rho; leading p_{t}^{track}", 100,0.,5000.,500,0.,250.,100,0.,100.);
+    fh3MultvsRhoLeadingTrackPtQ4 = new TH3F("fh3MultvsRhoLeadingTrackPtQ4","multiplicity vs background density Q4; multiplicity; #rho; leading p_{t}^{track}", 100,0.,5000.,500,0.,250.,100,0.,100.);
+
+    fh3MultvsSigmaLeadingTrackPtQ1 = new TH3F("fh3MultvsSigmaLeadingTrackPtQ1","multiplicity vs background density Q1; multiplicity; #sigma(#rho); leading p_{t}^{track}", 100,0.,5000.,50,0.,50.,100,0.,100.);
+    fh3MultvsSigmaLeadingTrackPtQ2 = new TH3F("fh3MultvsSigmaLeadingTrackPtQ2","multiplicity vs background density Q2; multiplicity; #sigma(#rho); leading p_{t}^{track}", 100,0.,5000.,50,0.,50.,100,0.,100.);
+    fh3MultvsSigmaLeadingTrackPtQ3 = new TH3F("fh3MultvsSigmaLeadingTrackPtQ3","multiplicity vs background density Q3; multiplicity; #sigma(#rho); leading p_{t}^{track}", 100,0.,5000.,50,0.,50.,100,0.,100.);
+    fh3MultvsSigmaLeadingTrackPtQ4 = new TH3F("fh3MultvsSigmaLeadingTrackPtQ4","multiplicity vs background density Q4; multiplicity; #sigma(#rho); leading p_{t}^{track}", 100,0.,5000.,50,0.,50.,100,0.,100.);
+
+
+    fh3CentvsDeltaRhoLeadingTrackPtQ1 = new TH3F("fh3CentvsDeltaRhoLeadingTrackPtQ1","centrality vs background density Q1; centrality; #rho_{quadrant}-#rho_{full event}; leading p_{t}^{track}", 50,0.,100.,200,-10.,10.,100,0.,100.);
+    fh3CentvsDeltaRhoLeadingTrackPtQ2 = new TH3F("fh3CentvsDeltaRhoLeadingTrackPtQ2","centrality vs background density Q2; centrality; #rho_{quadrant}-#rho_{full event}; leading p_{t}^{track}", 50,0.,100.,200,-10.,10.,100,0.,100.);
+    fh3CentvsDeltaRhoLeadingTrackPtQ3 = new TH3F("fh3CentvsDeltaRhoLeadingTrackPtQ3","centrality vs background density Q3; centrality; #rho_{quadrant}-#rho_{full event}; leading p_{t}^{track}", 50,0.,100.,200,-10.,10.,100,0.,100.);
+    fh3CentvsDeltaRhoLeadingTrackPtQ4 = new TH3F("fh3CentvsDeltaRhoLeadingTrackPtQ4","centrality vs background density Q4; centrality; #rho_{quadrant}-#rho_{full event}; leading p_{t}^{track}", 50,0.,100.,200,-10.,10.,100,0.,100.);
+
+    fHistList->Add(fh3CentvsRhoLeadingTrackPt);
+    fHistList->Add(fh3CentvsSigmaLeadingTrackPt);
+    fHistList->Add(fh3MultvsRhoLeadingTrackPt);
+    fHistList->Add(fh3MultvsSigmaLeadingTrackPt);
+
+    fHistList->Add(fh3CentvsRhoLeadingTrackPtQ1);
+    fHistList->Add(fh3CentvsRhoLeadingTrackPtQ2);
+    fHistList->Add(fh3CentvsRhoLeadingTrackPtQ3);
+    fHistList->Add(fh3CentvsRhoLeadingTrackPtQ4);
+
+    fHistList->Add(fh3CentvsSigmaLeadingTrackPtQ1);
+    fHistList->Add(fh3CentvsSigmaLeadingTrackPtQ2);
+    fHistList->Add(fh3CentvsSigmaLeadingTrackPtQ3);
+    fHistList->Add(fh3CentvsSigmaLeadingTrackPtQ4);
+
+    fHistList->Add(fh3MultvsRhoLeadingTrackPtQ1);
+    fHistList->Add(fh3MultvsRhoLeadingTrackPtQ2);
+    fHistList->Add(fh3MultvsRhoLeadingTrackPtQ3);
+    fHistList->Add(fh3MultvsRhoLeadingTrackPtQ4);
+
+    fHistList->Add(fh3MultvsSigmaLeadingTrackPtQ1);
+    fHistList->Add(fh3MultvsSigmaLeadingTrackPtQ2);
+    fHistList->Add(fh3MultvsSigmaLeadingTrackPtQ3);
+    fHistList->Add(fh3MultvsSigmaLeadingTrackPtQ4);
+
+    fHistList->Add(fh3CentvsDeltaRhoLeadingTrackPtQ1);
+    fHistList->Add(fh3CentvsDeltaRhoLeadingTrackPtQ2);
+    fHistList->Add(fh3CentvsDeltaRhoLeadingTrackPtQ3);
+    fHistList->Add(fh3CentvsDeltaRhoLeadingTrackPtQ4);
+
+  }
+  
   //Detector level effects histos
   fh2PtGenPtSmeared = new TH2F("fh2PtGenPtSmeared","fh2PtGenPtSmeared",nBinPt,binLimitsPt,nBinPt,binLimitsPt);
 
@@ -703,7 +847,7 @@ void AliAnalysisTaskJetCluster::UserCreateOutputObjects()
   TH1::AddDirectory(oldStatus);
 }
 
-void AliAnalysisTaskJetCluster::Init()
+void AliAnalysisTaskJetCluster::LocalInit()
 {
   //
   // Initialization
@@ -711,6 +855,10 @@ void AliAnalysisTaskJetCluster::Init()
 
   if (fDebug > 1) printf("AnalysisTaskJetCluster::Init() \n");
 
+  if(fUseTrPtResolutionFromOADB) LoadTrPtResolutionRootFileFromOADB();
+  if(fUseTrEfficiencyFromOADB)   LoadTrEfficiencyRootFileFromOADB();
+
+
   FitMomentumResolution();
 
 }
@@ -757,9 +905,10 @@ void AliAnalysisTaskJetCluster::UserExec(Option_t */*option*/)
   }
 
   //Check if information is provided detector level effects
-  if(!fMomResH1 || !fMomResH2 || !fMomResH3) fUseTrMomentumSmearing = kFALSE;
-  if(!fhEffH1 || !fhEffH2 || !fhEffH3)       fUseDiceEfficiency = kFALSE;
-  
+  if(!fMomResH1 || !fMomResH2 || !fMomResH3) fUseTrPtResolutionSmearing = kFALSE;
+  if(!fhEffH1 || !fhEffH2 || !fhEffH3 )      fUseDiceEfficiency = kFALSE;
+  if(  fEfficiencyFixed < 1. )               fUseDiceEfficiency = kTRUE;  
+
   Bool_t selectEvent =  false;
   Bool_t physicsSelection = true;// handled by the framework(fInputHandler->IsEventSelected()&AliVEvent::kMB)==AliVEvent::kMB;
 
@@ -801,7 +950,28 @@ void AliAnalysisTaskJetCluster::UserExec(Option_t */*option*/)
       selectEvent = true;
     }
   }
+
+
+  Bool_t T0 = false;
+  Bool_t V0 = false;
+  const AliAODVZERO  *vzero = fAOD->GetVZEROData();
+  if(vzero){
+    if((vzero->GetTriggerChargeA()>0)&&(vzero->GetTriggerChargeC()>0)){
+      V0 = true;
+    }
+  }
   
+  const AliAODTZERO  *tzero = fAOD->GetTZEROData();
+  if(tzero){
+    if(TMath::Abs(tzero->GetT0VertexRaw())<100){
+      T0 = true;
+    }
+  }
+  
+  if(fRequireVZEROAC&&fRequireTZEROvtx)selectEvent = selectEvent&&V0&&T0;
+  else if(fRequireTZEROvtx)selectEvent = selectEvent&&T0;
+  else if(fRequireVZEROAC)selectEvent = selectEvent&&V0;
+
 
   if(!selectEvent){
     PostData(1, fHistList);
@@ -844,7 +1014,7 @@ void AliAnalysisTaskJetCluster::UserExec(Option_t */*option*/)
     // Carefull energy is not well determined in real data, should not matter for p_T scheme?
     // we take total momentum here
 
-    if((!fUseTrMomentumSmearing) && (!fUseDiceEfficiency)) {
+    if((!fUseTrPtResolutionSmearing) && (!fUseDiceEfficiency)) {
       //Add particles to fastjet in case we are not running toy model
       fastjet::PseudoJet jInp(vp->Px(),vp->Py(),vp->Pz(),vp->P());
       jInp.set_user_index(i);
@@ -854,67 +1024,76 @@ void AliAnalysisTaskJetCluster::UserExec(Option_t */*option*/)
 
       // Dice to decide if particle is kept or not - toy  model for efficiency
       //
-      Double_t rnd = fRandom->Uniform(1.);
-      Double_t pT = vp->Pt();
+      Double_t sumEff = 0.;
+      Double_t pT = 0.;
       Double_t eff[3] = {0.};
-      Double_t pTtmp = pT;
-      if(pT>10.) pTtmp = 10.;
-      Double_t eff1 = fhEffH1->GetBinContent(fhEffH1->FindBin(pTtmp));
-      Double_t eff2 = fhEffH2->GetBinContent(fhEffH2->FindBin(pTtmp));
-      Double_t eff3 = fhEffH3->GetBinContent(fhEffH3->FindBin(pTtmp));
       Int_t cat[3] = {0};
-      //Sort efficiencies from large to small
-      if(eff1>eff2 && eff1>eff3) { 
-       eff[0] = eff1; 
-       cat[0] = 1;
-       if(eff2>eff3) {
-         eff[1] = eff2;
-         eff[2] = eff3; 
-         cat[1] = 2;
-         cat[2] = 3;
-       } else {
-         eff[1] = eff3;
-         eff[2] = eff2; 
-         cat[1] = 3;
-         cat[2] = 2;
+
+      Double_t rnd = fRandom->Uniform(1.);
+      if(  fEfficiencyFixed<1. ) {
+       sumEff = fEfficiencyFixed;
+      } else {
+
+       pT = vp->Pt();
+       Double_t pTtmp = pT;
+       if(pT>10.) pTtmp = 10.;
+       Double_t eff1 = fhEffH1->GetBinContent(fhEffH1->FindBin(pTtmp));
+       Double_t eff2 = fhEffH2->GetBinContent(fhEffH2->FindBin(pTtmp));
+       Double_t eff3 = fhEffH3->GetBinContent(fhEffH3->FindBin(pTtmp));
+
+       //Sort efficiencies from large to small
+       if(eff1>eff2 && eff1>eff3) { 
+         eff[0] = eff1; 
+         cat[0] = 1;
+         if(eff2>eff3) {
+           eff[1] = eff2;
+           eff[2] = eff3; 
+           cat[1] = 2;
+           cat[2] = 3;
+         } else {
+           eff[1] = eff3;
+           eff[2] = eff2; 
+           cat[1] = 3;
+           cat[2] = 2;
+         }
        }
-      }
-      else if(eff2>eff1 && eff2>eff3) {
-       eff[0] = eff2;
-       cat[0] = 2;
-       if(eff1>eff3) {
-         eff[1] = eff1;
-         eff[2] = eff3; 
-         cat[1] = 1;
-         cat[2] = 3;
-       } else {
-         eff[1] = eff3;
-         eff[2] = eff1; 
-         cat[1] = 3;
-         cat[2] = 1;
+       else if(eff2>eff1 && eff2>eff3) {
+         eff[0] = eff2;
+         cat[0] = 2;
+         if(eff1>eff3) {
+           eff[1] = eff1;
+           eff[2] = eff3; 
+           cat[1] = 1;
+           cat[2] = 3;
+         } else {
+           eff[1] = eff3;
+           eff[2] = eff1; 
+           cat[1] = 3;
+           cat[2] = 1;
+         }
        }
-      }
-      else if(eff3>eff1 && eff3>eff2) {
-       eff[0] = eff3;
-       cat[0] = 3;
-       if(eff1>eff2) {
-         eff[1] = eff1;
-         eff[2] = eff2; 
-         cat[1] = 1;
-         cat[2] = 2;
-       } else {
-         eff[1] = eff2;
-         eff[2] = eff1; 
-         cat[1] = 2;
-         cat[2] = 1;
+       else if(eff3>eff1 && eff3>eff2) {
+         eff[0] = eff3;
+         cat[0] = 3;
+         if(eff1>eff2) {
+           eff[1] = eff1;
+           eff[2] = eff2; 
+           cat[1] = 1;
+           cat[2] = 2;
+         } else {
+           eff[1] = eff2;
+           eff[2] = eff1; 
+           cat[1] = 2;
+           cat[2] = 1;
+         }
        }
+       
+       sumEff = eff[0]+eff[1]+eff[2];
       }
-      
-      Double_t sumEff = eff[0]+eff[1]+eff[2];
       fp1Efficiency->Fill(vp->Pt(),sumEff);
-      if(rnd>sumEff) continue;
+      if(rnd>sumEff && pT > fDiceEfficiencyMinPt) continue;
 
-      if(fUseTrMomentumSmearing) {
+      if(fUseTrPtResolutionSmearing) {
        //Smear momentum of generated particle
        Double_t smear = 1.;
        //Select hybrid track category
@@ -976,11 +1155,11 @@ void AliAnalysisTaskJetCluster::UserExec(Option_t */*option*/)
     if(fTCAJetsOut){
       if(i == 0){
        fRef->Delete(); // make sure to delete before placement new...
-       if((!fUseTrMomentumSmearing) && (!fUseDiceEfficiency)) {
+       if((!fUseTrPtResolutionSmearing) && (!fUseDiceEfficiency)) {
          new(fRef) TRefArray(TProcessID::GetProcessWithUID(vp)); //TRefArray does not work with toy model ...
        } 
       }
-      if((!fUseTrMomentumSmearing) && (!fUseDiceEfficiency)) fRef->Add(vp);  //TRefArray does not work with toy model ...
+      if((!fUseTrPtResolutionSmearing) && (!fUseDiceEfficiency)) fRef->Add(vp);  //TRefArray does not work with toy model ...
     }
   }// recparticles
 
@@ -1119,18 +1298,40 @@ void AliAnalysisTaskJetCluster::UserExec(Option_t */*option*/)
       fh2NConstPt->Fill(tmpPt,constituents.size());
       // loop over constiutents and fill spectrum
    
+      AliVParticle *partLead = 0;
+      Float_t ptLead = -1;
+
       for(unsigned int ic = 0; ic < constituents.size();ic++){
        AliVParticle *part = (AliVParticle*)recParticles.At(constituents[ic].user_index());
        if(!part) continue;
-       
        fh1PtJetConstRec->Fill(part->Pt());
        if(aodOutJet){
-         if((!fUseTrMomentumSmearing) && (!fUseDiceEfficiency)) aodOutJet->AddTrack(fRef->At(constituents[ic].user_index()));
-         if(part->Pt()>fMaxTrackPtInJet)aodOutJet->SetTrigger(AliAODJet::kHighTrackPtTriggered);
+         if((!fUseTrPtResolutionSmearing) && (!fUseDiceEfficiency)) aodOutJet->AddTrack(fRef->At(constituents[ic].user_index()));
+         if(part->Pt()>fMaxTrackPtInJet){
+           aodOutJet->SetTrigger(AliAODJet::kHighTrackPtTriggered);
+         }
+       }
+       if(part->Pt()>ptLead){
+         ptLead = part->Pt();
+         partLead = part;
        }
        if(j==0)fh1PtJetConstLeadingRec->Fill(part->Pt());
       }
-      
+
+      AliAODTrack *aodT = 0;
+      if(partLead){
+       if(aodOutJet){
+         //set pT of leading constituent of jet
+         aodOutJet->SetPtLeading(partLead->Pt());
+         aodT = dynamic_cast<AliAODTrack*>(partLead);
+         if(aodT){
+           if(aodT->TestFilterBit(fFilterMaskBestPt)){
+             aodOutJet->SetTrigger(AliAODJet::kHighTrackPtBest);
+           }
+         }
+       }
+      }
+    
      // correlation
      Float_t tmpPhi =  tmpRec.Phi();
      Float_t tmpEta =  tmpRec.Eta();
@@ -1198,6 +1399,7 @@ void AliAnalysisTaskJetCluster::UserExec(Option_t */*option*/)
        }
        if(skip)continue;
        tmpRecC.SetBgEnergy(0,0); // this is use as temporary storage of the summed p_T below
+       tmpRecC.SetPtLeading(-1.);
        if(fTCARandomConesOut)new ((*fTCARandomConesOut)[nCone++]) AliAODJet(tmpRecC);
        if(fTCARandomConesOutRan)new ((*fTCARandomConesOutRan)[nConeRan++]) AliAODJet(tmpRecC);
      }// loop over random cones creation
@@ -1215,10 +1417,11 @@ void AliAnalysisTaskJetCluster::UserExec(Option_t */*option*/)
           if(jC&&jC->DeltaR(vp)<fRparam){
             if(vp->Pt()>fMaxTrackPtInJet)jC->SetTrigger(AliAODJet::kHighTrackPtTriggered);
             jC->SetBgEnergy(jC->ChargedBgEnergy()+vp->Pt(),0);
+            if(vp->Pt() > jC->GetPtLeading()) jC->SetPtLeading(vp->Pt());
           }
         }  
        }// add up energy in cone
-      
+
        // the randomized input changes eta and phi, but keeps the p_T
        if(i>=fNSkipLeadingRan){// eventually skip the leading particles
         Double_t pTR = vp->Pt();
@@ -1238,6 +1441,7 @@ void AliAnalysisTaskJetCluster::UserExec(Option_t */*option*/)
             if(jC&&jC->DeltaR(&vTmpRanR)<fRparam){
               if(vTmpRanR.Pt()>fMaxTrackPtInJet)jC->SetTrigger(AliAODJet::kHighTrackPtTriggered);
               jC->SetBgEnergy(jC->ChargedBgEnergy()+vTmpRanR.Pt(),0);
+              if(vTmpRanR.Pt() > jC->GetPtLeading()) jC->SetPtLeading(vTmpRanR.Pt());
             }
           }  
         }
@@ -1247,7 +1451,7 @@ void AliAnalysisTaskJetCluster::UserExec(Option_t */*option*/)
      Float_t jetArea = fRparam*fRparam*TMath::Pi();
      if(fTCARandomConesOut){
        for(int ir = 0;ir < fTCARandomConesOut->GetEntriesFast();ir++){
-        // rescale the momntum vectors for the random cones
+        // rescale the momentum vectors for the random cones
         
         AliAODJet *rC = (AliAODJet*)fTCARandomConesOut->At(ir);
         if(rC){
@@ -1291,14 +1495,11 @@ void AliAnalysisTaskJetCluster::UserExec(Option_t */*option*/)
    }// if(fNRandomCones
   
    //background estimates:all bckg jets(0) & wo the 2 hardest(1)
-
-
-
 
    if(fAODJetBackgroundOut){
      vector<fastjet::PseudoJet> jets2=sortedJets;
      if(jets2.size()>2) jets2.erase(jets2.begin(),jets2.begin()+2); 
+
      Double_t bkg1=0;
      Double_t sigma1=0.;
      Double_t meanarea1=0.;
@@ -1319,6 +1520,100 @@ void AliAnalysisTaskJetCluster::UserExec(Option_t */*option*/)
 
    }
   }
+
+  if(fStoreRhoLeadingTrackCorr) {
+    vector<fastjet::PseudoJet> jets3=sortedJets;
+    if(jets3.size()>2) jets3.erase(jets3.begin(),jets3.begin()+2); 
+
+    Double_t bkg2=0;
+    Double_t sigma2=0.;
+    Double_t meanarea2=0.;
+
+    clustSeq.get_median_rho_and_sigma(jets3, range, false, bkg2, sigma2, meanarea2, true);
+    fAODJetBackgroundOut->SetBackground(1,bkg2,sigma2,meanarea2);
+
+    //Get leading track in event
+    AliVParticle *leading = (AliVParticle*)recParticles.At(0);
+
+    fh3CentvsRhoLeadingTrackPt->Fill(cent,bkg2,leading->Pt());
+    fh3CentvsSigmaLeadingTrackPt->Fill(cent,sigma2,leading->Pt());
+    fh3MultvsRhoLeadingTrackPt->Fill(nCh,bkg2,leading->Pt());
+    fh3MultvsSigmaLeadingTrackPt->Fill(nCh,sigma2,leading->Pt());
+     
+
+    //Calculate rho with 4-vector area method for different quadrants with respect to the leading track in the event
+    //exclude 2 leading jets from event
+    //Quadrant 1: |phi_leadingTrack - phi_bkgCluster| < pi/2./2. - R  (Near side to leading track)
+    //Quadrant 2: pi/2 - (pi/2./2. - R) < |phi_leadingTrack - phi_bkgCluster| < pi/2 + (pi/2./2. - R) (Orthogonal to leading track)
+    //Quadrant 3: pi - (pi/2/2 - R) < |phi_leadingTrack - phi_bkgCluster| < pi + (pi/2./2. - R)  (Away side to leading track)
+    //Quadrant 4: 3/2*pi - (pi/2./2. - R) < |phi_leadingTrack - phi_bkgCluster| < 3/2*pi + (pi/2./2. - R)   (Orthogonal to leading track)
+
+    //Assuming R=0.2 for background clusters
+
+    Double_t bkg2Q[4]      = {0.};
+    Double_t sigma2Q[4]    = {0.};
+    Double_t meanarea2Q[4] = {0.};
+
+    //Get phi of leading track in event
+    Float_t phiLeadingTrack = leading->Phi();
+    Float_t phiStep = (TMath::Pi()/2./2. - 0.2);
+
+    //Quadrant1 - near side
+    phiMin = phiLeadingTrack - phiStep;
+    phiMax = phiLeadingTrack + phiStep;
+    fastjet::RangeDefinition rangeQ1(rapMin,rapMax, phiMin, phiMax);
+    clustSeq.get_median_rho_and_sigma(jets3, rangeQ1, false, bkg2Q[0], sigma2Q[0], meanarea2Q[0], true);
+
+    //Quadrant2 - orthogonal
+    Double_t phiQ2 = phiLeadingTrack + TMath::Pi()/2.;
+    if(phiQ2>TMath::TwoPi()) phiQ2 = phiQ2 - TMath::TwoPi();
+    phiMin = phiQ2 - phiStep;
+    phiMax = phiQ2 + phiStep;
+    fastjet::RangeDefinition rangeQ2(rapMin,rapMax, phiMin, phiMax);
+    clustSeq.get_median_rho_and_sigma(jets3, rangeQ2, false, bkg2Q[1], sigma2Q[1], meanarea2Q[1], true);
+
+    //Quadrant3 - away side
+    Double_t phiQ3 = phiLeadingTrack + TMath::Pi();
+    if(phiQ3>TMath::TwoPi()) phiQ3 = phiQ3 - TMath::TwoPi();
+    phiMin = phiQ3 - phiStep;
+    phiMax = phiQ3 + phiStep;
+    fastjet::RangeDefinition rangeQ3(rapMin,rapMax, phiMin, phiMax);
+    clustSeq.get_median_rho_and_sigma(jets3, rangeQ3, false, bkg2Q[2], sigma2Q[2], meanarea2Q[2], true);
+
+    //Quadrant4 - othogonal
+    Double_t phiQ4 = phiLeadingTrack + 3./2.*TMath::Pi();
+    if(phiQ4>TMath::TwoPi()) phiQ4 = phiQ4 - TMath::TwoPi();
+    phiMin = phiQ4 - phiStep;
+    phiMax = phiQ4 + phiStep;
+    fastjet::RangeDefinition rangeQ4(rapMin,rapMax, phiMin, phiMax);
+    clustSeq.get_median_rho_and_sigma(jets3, rangeQ4, false, bkg2Q[3], sigma2Q[3], meanarea2Q[3], true);
+
+    fh3CentvsRhoLeadingTrackPtQ1->Fill(cent,bkg2Q[0],leading->Pt());
+    fh3CentvsRhoLeadingTrackPtQ2->Fill(cent,bkg2Q[1],leading->Pt());
+    fh3CentvsRhoLeadingTrackPtQ3->Fill(cent,bkg2Q[2],leading->Pt());
+    fh3CentvsRhoLeadingTrackPtQ4->Fill(cent,bkg2Q[3],leading->Pt());
+
+    fh3CentvsSigmaLeadingTrackPtQ1->Fill(cent,sigma2Q[0],leading->Pt());
+    fh3CentvsSigmaLeadingTrackPtQ2->Fill(cent,sigma2Q[1],leading->Pt());
+    fh3CentvsSigmaLeadingTrackPtQ3->Fill(cent,sigma2Q[2],leading->Pt());
+    fh3CentvsSigmaLeadingTrackPtQ4->Fill(cent,sigma2Q[3],leading->Pt());
+
+    fh3MultvsRhoLeadingTrackPtQ1->Fill(nCh,bkg2Q[0],leading->Pt());
+    fh3MultvsRhoLeadingTrackPtQ2->Fill(nCh,bkg2Q[1],leading->Pt());
+    fh3MultvsRhoLeadingTrackPtQ3->Fill(nCh,bkg2Q[2],leading->Pt());
+    fh3MultvsRhoLeadingTrackPtQ4->Fill(nCh,bkg2Q[3],leading->Pt());
+
+    fh3MultvsSigmaLeadingTrackPtQ1->Fill(nCh,sigma2Q[0],leading->Pt());
+    fh3MultvsSigmaLeadingTrackPtQ2->Fill(nCh,sigma2Q[1],leading->Pt());
+    fh3MultvsSigmaLeadingTrackPtQ3->Fill(nCh,sigma2Q[2],leading->Pt());
+    fh3MultvsSigmaLeadingTrackPtQ4->Fill(nCh,sigma2Q[3],leading->Pt());
+
+    fh3CentvsDeltaRhoLeadingTrackPtQ1->Fill(cent,bkg2Q[0]-bkg2,leading->Pt());
+    fh3CentvsDeltaRhoLeadingTrackPtQ2->Fill(cent,bkg2Q[1]-bkg2,leading->Pt());
+    fh3CentvsDeltaRhoLeadingTrackPtQ3->Fill(cent,bkg2Q[2]-bkg2,leading->Pt());
+    fh3CentvsDeltaRhoLeadingTrackPtQ4->Fill(cent,bkg2Q[3]-bkg2,leading->Pt());
+
+  }
    
 
   
@@ -1560,6 +1855,7 @@ Int_t  AliAnalysisTaskJetCluster::GetListOfTracks(TList *list,Int_t type){
        if(fDebug>2)Printf("%s:%d No AOD",(char*)__FILE__,__LINE__);
        return iCount;
       }
+
       for(int it = 0;it < aod->GetNumberOfTracks();++it){
        AliAODTrack *tr = aod->GetTrack(it);
        Bool_t bGood = false;
@@ -1664,15 +1960,64 @@ Int_t  AliAnalysisTaskJetCluster::GetListOfTracks(TList *list,Int_t type){
   return iCount;
 }
 
+void AliAnalysisTaskJetCluster::LoadTrPtResolutionRootFileFromOADB() {
+
+  TFile *f = TFile::Open(fPathTrPtResolution.Data());
+
+  if(!f)return;
+
+  TProfile *fProfPtPtSigma1PtGlobSt     = (TProfile*)f->Get("fProfPtPtSigma1PtGlobSt");
+  TProfile *fProfPtPtSigma1PtGlobCnoITS = (TProfile*)f->Get("fProfPtPtSigma1PtGlobCnoITS");
+  TProfile *fProfPtPtSigma1PtGlobCnoSPD = (TProfile*)f->Get("fProfPtPtSigma1PtGlobCnoSPD");
+
+  SetSmearResolution(kTRUE);
+  SetMomentumResolutionHybrid(fProfPtPtSigma1PtGlobSt,fProfPtPtSigma1PtGlobCnoITS,fProfPtPtSigma1PtGlobCnoSPD);
+  
+
+}
+
+void AliAnalysisTaskJetCluster::LoadTrEfficiencyRootFileFromOADB() {
+
+  TFile *f = TFile::Open(fPathTrEfficiency.Data());
+  if(!f)return;
+
+  TH1D *hEffPosGlobSt = (TH1D*)f->Get("hEffPosGlobSt");
+  TH1D *hEffPosGlobCnoITS = (TH1D*)f->Get("hEffPosGlobCnoITS");
+  TH1D *hEffPosGlobCnoSPD = (TH1D*)f->Get("hEffPosGlobCnoSPD");
+
+  SetDiceEfficiency(kTRUE);
+
+  if(fChangeEfficiencyFraction>0.) {
+
+    TH1D *hEffPosGlobStMin = (TH1D*)hEffPosGlobSt->Clone("hEffPosGlobStMin");
+    
+    for(int bin=1; bin<=hEffPosGlobSt->GetNbinsX(); bin++) {
+      Double_t content = hEffPosGlobSt->GetBinContent(bin);
+      hEffPosGlobStMin->SetBinContent(bin,content-fChangeEfficiencyFraction);
+    }
+
+    SetEfficiencyHybrid(hEffPosGlobStMin,hEffPosGlobCnoITS,hEffPosGlobCnoSPD);
+
+  }
+  else
+    SetEfficiencyHybrid(hEffPosGlobSt,hEffPosGlobCnoITS,hEffPosGlobCnoSPD);
+
+}
+
 void AliAnalysisTaskJetCluster::SetMomentumResolutionHybrid(TProfile *p1, TProfile *p2, TProfile *p3) {
 
   //
   // set mom res profiles
   //
 
-  fMomResH1 = (TProfile*)p1->Clone("fMomResH1");
-  fMomResH2 = (TProfile*)p2->Clone("fMomResH2");
-  fMomResH3 = (TProfile*)p3->Clone("fMomResH3");
+  if(fMomResH1) delete fMomResH1;
+  if(fMomResH2) delete fMomResH2;
+  if(fMomResH3) delete fMomResH3;
+
+  fMomResH1 = new TProfile(*p1);//(TProfile*)p1->Clone("fMomResH1");
+  fMomResH2 = new TProfile(*p2);//(TProfile*)p2->Clone("fMomResH2");
+  fMomResH3 = new TProfile(*p3);//(TProfile*)p3->Clone("fMomResH3");
+
 }
 
 void AliAnalysisTaskJetCluster:: SetEfficiencyHybrid(TH1 *h1, TH1 *h2, TH1 *h3) {