revamped jet chemistry task (O. Busch)
authorkleinb <kleinb@f7af4fe6-9843-0410-8265-dc069ae4e863>
Thu, 23 Jun 2011 19:01:27 +0000 (19:01 +0000)
committerkleinb <kleinb@f7af4fe6-9843-0410-8265-dc069ae4e863>
Thu, 23 Jun 2011 19:01:27 +0000 (19:01 +0000)
PWG4/JetTasks/AliAnalysisTaskFragmentationFunction.h
PWG4/JetTasks/AliAnalysisTaskJetChem.cxx
PWG4/JetTasks/AliAnalysisTaskJetChem.h
PWG4/PWG4JetTasksLinkDef.h
PWG4/macros/AddTaskJetChem.C

index ac4277c..7f528bf 100644 (file)
@@ -484,7 +484,7 @@ class AliAnalysisTaskFragmentationFunction : public AliAnalysisTaskSE {
        kBckgOut2JStat, kBckgOut3JStat, kBckgOutAJStat,  kBckgASide, kBckgASideWindow, kBckgPerpWindow};
 
  
- private:
+ protected:
   
   Int_t   GetListOfTracks(TList* list, Int_t type);
   Int_t          GetListOfJets(TList* list, Int_t type);
index 403fd4d..618a6a0 100644 (file)
@@ -1,3 +1,11 @@
+/*************************************************************************
+ *                                                                       *
+ *                                                                       *
+ *      Task for Jet Chemistry Analysis in PWG4 Jet Task Force Train     *
+ *                                                                       *
+ *                                                                       *
+ *************************************************************************/
+
 /**************************************************************************
  * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
  *                                                                        *
  * provided "as is" without express or implied warranty.                  *
  **************************************************************************/
 
-/* $Id:$ */
-
+/* $Id: */
 
-#include <TChain.h>
-#include <TFile.h>
-#include <TList.h>
-#include <TH1F.h>
-#include <TH2F.h>
-#include <TProfile.h>
-#include <TVector3.h>
-#include <TLorentzVector.h>
-#include <TMath.h>
-#include <TTree.h>
-#include <TDatabasePDG.h> 
+#include "TH2F.h"
+#include "TH3F.h"
+#include "TProfile.h"
+#include "THnSparse.h"
 
-#include "AliAnalysisTaskJetChem.h"
+#include "AliAnalysisHelperJetTasks.h"
 #include "AliAnalysisManager.h"
-#include "AliMCEventHandler.h"
-#include "AliAODInputHandler.h"
-#include "AliAODHandler.h"
-#include "AliAODJet.h"
-#include "AliAODTrack.h"
-#include "AliAODMCParticle.h"
-
+#include "AliAODHandler.h" 
+#include "AliAODInputHandler.h" 
+#include "AliESDEvent.h"
 #include "AliGenPythiaEventHeader.h"
-#include "AliAnalysisHelperJetTasks.h"
-#include "AliESDtrack.h"
-#include "AliExternalTrackParam.h"
-
-//#include <iostream> // TEST!!!
+#include "AliGenHijingEventHeader.h"
 
-
-//
-// Analysis class for jet chemistry studies
-// based on AliAnalysisTaskUE by Arian Abrahantes Quintana and Enesto Lopez
-// contact: Oliver Busch, o.busch@gsi.de
-
-// This class needs as input AOD with track and Jets
-// the output is a list of histograms
-//
-// AOD can be either connected to the InputEventHandler  
-// for a chain of AOD files 
-// or 
-// to the OutputEventHandler
-// for a chain of ESD files, so this case class should be 
-// in the train after the Jet finder
-//
+#include "AliAODEvent.h"
+#include "AliAODJet.h"
+#include "AliAODv0.h"
+#include "AliAODTrack.h"
 
 
+#include "AliPID.h" 
+#include "AliExternalTrackParam.h"
 
+#include "AliAnalysisTaskJetChem.h"
 
 ClassImp(AliAnalysisTaskJetChem)
 
-////////////////////////////////////////////////////////////////////////
-
-
-//____________________________________________________________________
-AliAnalysisTaskJetChem:: AliAnalysisTaskJetChem(const char* name): AliAnalysisTaskSE(name),
-fDebug(kFALSE),
-fDeltaAOD(kFALSE),
-fDeltaAODBranch(""),
-fAODBranch("jets"),
-fDeltaAODBranchMC(""),
-fAODBranchMC(""),
-fArrayJetsAOD(0x0), 
-fArrayJetsMC(0x0),
-fAOD(0x0),            
-fAODjets(0x0),
-fListOfHistos(0x0),
-fJetsOnFly(kFALSE),
-fUseLOConeJets(kFALSE),
-fUseLOConeMCJets(kFALSE),
-fUsePythiaJets(kFALSE),
-fConeRadius(0.),
-fTrackPtCutJF(0),
-fFilterBitJF(0x01),
-fRequireITSRefitJF(kFALSE),
-fRejectK0TracksJF(kFALSE),
-fJetPtCut(5.0),
-fJetEtaCut(0.5),
-fFilterBit(0x10),
-fTrackPtCut(0.),
-fTrackEtaCut(0.9),
-fUseOnFlyV0s(kFALSE),
-fCutnSigdEdx(0.), 
-fUseAODMCTracksForUE(kFALSE),
-fAreaReg(1.0), 
-fAvgTrials(1),
-fhPrimVertexNCont(0x0),
-fhPrimVertexRho(0x0),
-fhPrimVertexZ(0x0),
-fhNJets(0x0),
-fhNJetsMC(0x0),
-fhLeadingEta(0x0),
-fhLeadingNTracksVsEta(0x0),
-fhLeadingPtVsEta(0x0),
-fhLeadingPhi(0x0), 
-fhLeadingPt(0x0),
-fhLeadingPtDiffr(0x0),
-fhLeadingEtaMC(0x0), 
-fhLeadingPhiMC(0x0), 
-fhLeadingPtMC(0x0),
-fhLeadingPtMCDiffr(0x0),
-fhPhiEtaTracksNoCut(0x0),
-fhPtTracksNoCut(0x0),
-fhPhiEtaTracks(0x0),
-fhPtTracks(0x0),
-fhTrackMult(0x0),
-fhEtaMCTracks(0x0),            
-fhPhiMCTracks(0x0),            
-fhPtMCTracks(0x0),
-fhnTracksVsPtLeading(0x0),
-fhdNdEtaPhiDist(0x0),        
-fhRegionSumPtMaxVsEt(0x0),
-fhRegionMultMaxVsEt(0x0),     
-fhRegionSumPtMinVsEt(0x0),
-fhRegionMultMinVsEt(0x0),     
-fhNV0s(0x0),
-fhV0onFly(0x0),
-fhV0DCADaughters(0x0),
-fhV0Radius(0x0),
-fhV0DCAToVertex(0x0),
-fhV0DCAToVertexK0(0x0),
-fhV0InvMassK0(0x0),
-fhV0PtVsInvMassK0(0x0),
-fhV0InvMassK0JetEvt(0x0),
-fhV0InvMassLambda(0x0),
-fhV0InvMassAntiLambda(0x0),
-fhV0InvMassLambdaJetEvt(0x0),
-fhV0InvMassAntiLambdaJetEvt(0x0),
-fhdROpanK0VsPt(0x0),
-fhdPhiJetV0(0x0),
-fhdPhiJetK0(0x0),
-fhdRJetK0(0x0),
-fhdNdptV0(0x0),
-fhdNdptK0(0x0),
-fhPtVsEtaK0(0x0),
-fhV0InvMassK0DCA(0x0),
-fhV0InvMassK0DCAdEdx(0x0),
-fhV0PtVsInvMassK0DCAdEdx(0x0),
-fhV0InvMassK0DCAPID(0x0),
-fhV0InvMassLambdaDCAdEdx(0x0), 
-fhV0InvMassAntiLambdaDCAdEdx(0x0),
-fhdNdptK0DCA(0x0),
-fhdNdptK0DCAdEdx(0x0),
-fhV0InvMassK0Min(0x0),
-fhV0InvMassLambdaMin(0x0),
-fhV0InvMassAntiLambdaMin(0x0),
-fhV0InvMassK0Max(0x0),
-fhV0InvMassLambdaMax(0x0),
-fhV0InvMassAntiLambdaMax(0x0),
-fhV0InvMassK0Jet(0x0),
-fhV0InvMassLambdaJet(0x0),
-fhV0InvMassAntiLambdaJet(0x0),
-fhV0InvMassK0Lambda(0x0), 
-fhdNdptK0JetEvt(0x0),
-fhdNdzK0(0x0),  
-fhdNdzK05to10(0x0),  
-fhdNdzK010to20(0x0), 
-fhdNdzK020to30(0x0), 
-fhdNdzK030to40(0x0), 
-fhdNdzK040to60(0x0), 
-fhdNdxiK0(0x0),
-fhdNdzLambda(0x0),
-fhdNdzAntiLambda(0x0),
-fhdNdzK0Max(0x0),     
-fhdNdxiK0Max(0x0),    
-fhdNdzLambdaMax(0x0), 
-fhdNdxiLambdaMax(0x0),
-fhdNdptK0Max(0x0), 
-fhdNdptLambdaMax(0x0),
-fhdNdzK0Min(0x0),     
-fhdNdxiK0Min(0x0),    
-fhdNdzLambdaMin(0x0), 
-fhdNdxiLambdaMin(0x0),
-fhdNdptK0Min(0x0),
-fhdNdptLambdaMin(0x0),
-fhdNdzK0Jet(0x0),     
-fhdNdxiK0Jet(0x0),    
-fhdNdzLambdaJet(0x0), 
-fhdNdxiLambdaJet(0x0),
-fhdNdptK0Jet(0x0),
-fhdNdptLambdaJet(0x0),
-fhdEdxVsMomV0(0x0),
-fhdEdxVsMomV0pidEdx(0x0),
-fhdEdxVsMomV0piPID(0x0),
-fhdPhiJetK0MC(0x0),
-fhdRJetK0MC(0x0),
-fhdRV0MC(0x0),
-fhdNdptchPiMCMax(0x0),    
-fhdNdptK0MCMax(0x0),      
-fhdNdptchKMCMax(0x0),     
-fhdNdptpMCMax(0x0),       
-fhdNdptpBarMCMax(0x0),    
-fhdNdptLambdaMCMax(0x0),  
-fhdNdptLambdaBarMCMax(0x0),
-fhdNdptchPiMCMin(0x0),    
-fhdNdptK0MCMin(0x0),      
-fhdNdptchKMCMin(0x0),     
-fhdNdptpMCMin(0x0),       
-fhdNdptpBarMCMin(0x0),    
-fhdNdptLambdaMCMin(0x0),  
-fhdNdptLambdaBarMCMin(0x0),
-fhdNdptOmegaMCMin(0x0),   
-fhdNdptOmegaBarMCMin(0x0),
-fhdNdptchPiMCJet(0x0),    
-fhdNdptK0MCJet(0x0),      
-fhdNdptchKMCJet(0x0),     
-fhdNdptpMCJet(0x0),       
-fhdNdptpBarMCJet(0x0),    
-fhdNdptLambdaMCJet(0x0),  
-fhdNdptLambdaBarMCJet(0x0),
-fhPIDMC(0x0),
-fhPIDMC_quarkEv(0x0),
-fhPIDMC_gluonEv(0x0),
-fhPIDMCAll(0x0),
-fhPIDMCMin(0x0),
-fhPIDMCJet(0x0),
-fhPIDMCMotherK0(0x0),
-fhPIDMCGrandMotherK0(0x0),
-fhPIDMCMotherChK(0x0),
-fhPIDMCMotherK0Trans(0x0),
-fhPIDMCGrandMotherK0Trans(0x0),
-fhPIDMCMotherChKTrans(0x0),
-fhdNdptgammaMC(0x0),
-fhdNdptchPiMC(0x0),
-fhdNdptpi0MC(0x0),
-fhdNdptK0MC(0x0),
-fhdNdptchKMC(0x0),
-fhdNdptpMC(0x0),
-fhdNdptpBarMC(0x0),
-fhdNdptLambdaMC(0x0),
-fhdNdptLambdaBarMC(0x0),
-fhdNdptOmegaMC(0x0),
-fhdNdptOmegaBarMC(0x0),
-fhdNdxiMC(0x0),
-fhdNdxiK0MC(0x0),
-fhdNdxiK0MCJet(0x0),
-fhdNdzK0MC(0x0), 
-fhdNdzK0MCJet(0x0), 
-fhdNdptK0MCJetEvt(0x0),
-fhnJetsAODvsMC(0x0),
-fhLeadingPtAODvsMC(0x0),
-fhLeadingEtaAODvsMC(0x0),
-fhLeadingPhiAODvsMC(0x0),
-fhnTracksLeadingAODvsMC(0x0),
-fhLeadingdRAODMC(0x0),
-fhLeadingPtAODvsMCdRcut(0x0),
-fhdnTracksVsdPtLeadingAODMC(0x0),
-fhnTracksJetVsPtAOD(0x0),
-fhnTracksJetVsPtAODquarkEv(0x0),
-fhRadiusJetVsPtAOD(0x0), 
-fhnTracksJetVsPtMC(0x0),
-fhnTracksJetVsPtMCquarkEv(0x0),
-fhRadiusJetVsPtMC(0x0),
-fhnTracksJetVsPtMCK0(0x0),
-fhnTracksJetVsPtMCK0quarkEv(0x0),
-fhRadiusJetVsPtMCK0(0x0),
-fhnTracksJetVsPtAODK0(0x0),
-fhnTracksJetVsPtAODK0quarkEv(0x0),
-fhRadiusJetVsPtAODK0(0x0),
-fhnTracksJetVsPtAODpKch(0x0),
-fhRadiusJetVsPtAODpKch(0x0),
-fhPythiaProcess(0x0), 
-fhPythiaProcessK0(0x0),
-fhPythiaProcessKch(0x0),
-fhPythiaProcessp(0x0),
-fhPythiaProcesspbar(0x0),
-fhdNdzJets5to10(0x0),   
-fhdNdzJets10to20(0x0),  
-fhdNdzJets20to30(0x0),  
-fhdNdzJets30to40(0x0),  
-fhdNdzJets40to60(0x0),  
-fhdNdxiJets5to10(0x0),  
-fhdNdxiJets10to20(0x0),
-fhdNdxiJets20to30(0x0), 
-fhdNdxiJets30to40(0x0), 
-fhdNdxiJets40to60(0x0), 
-fhdNdptTracksJetPt5to10(0x0),
-fhdNdptTracksJetPt10to20(0x0),
-fhdNdptTracksJetPt20to30(0x0),
-fhdNdptTracksJetPt30to40(0x0),
-fhdNdptTracksJetPt40to60(0x0),
-fh1Xsec(0x0),
-fh1Trials(0x0),
-fpdgdb(0x0){
-  // Default constructor 
-
-  fAreaReg = 2*TMath::Pi()/6.0 * 2*fTrackEtaCut;
-  fpdgdb = TDatabasePDG::Instance(); 
-
-  // Output slot #1 writes into a TList container, 0 reserved for std AOD output
-  DefineOutput(1, TList::Class());
-}
-
-//______________________________________________________________
-Bool_t AliAnalysisTaskJetChem::UserNotify()
+//____________________________________________________________________________
+AliAnalysisTaskJetChem::AliAnalysisTaskJetChem()
+   : AliAnalysisTaskFragmentationFunction()
+   ,fK0Type(0)
+   ,fFilterMaskK0(0)
+   ,fListK0s(0)
+   ,fV0QAK0(0)
+   ,fFFHistosRecCutsK0Evt(0)      
+   ,fFFHistosIMK0AllEvt(0)        
+   ,fFFHistosIMK0Jet(0)           
+   ,fFFHistosIMK0Cone(0)
+   ,fFFHistosPhiCorrIMK0(0)
+   ,fFFIMNBinsJetPt(0)    
+   ,fFFIMJetPtMin(0) 
+   ,fFFIMJetPtMax(0)
+   ,fFFIMNBinsInvM(0) 
+   ,fFFIMInvMMin(0)   
+   ,fFFIMInvMMax(0)   
+   ,fFFIMNBinsPt(0)      
+   ,fFFIMPtMin(0)        
+   ,fFFIMPtMax(0)        
+   ,fFFIMNBinsXi(0)      
+   ,fFFIMXiMin(0)        
+   ,fFFIMXiMax(0)        
+   ,fFFIMNBinsZ(0)       
+   ,fFFIMZMin(0)         
+   ,fFFIMZMax(0)
+   ,fPhiCorrIMNBinsPt(0)
+   ,fPhiCorrIMPtMin(0)
+   ,fPhiCorrIMPtMax(0)
+   ,fPhiCorrIMNBinsPhi(0)
+   ,fPhiCorrIMPhiMin(0)
+   ,fPhiCorrIMPhiMax(0)
+   ,fPhiCorrIMNBinsInvM(0)
+   ,fPhiCorrIMInvMMin(0)
+   ,fPhiCorrIMInvMMax(0)
+   ,fh1K0Mult(0)
+   ,fh1dPhiJetK0(0)
 {
-  //
-  // read the cross sections
-  // and number of trials from pyxsec.root
-  // 
-
-  fAvgTrials = 1;
-  TTree *tree = AliAnalysisManager::GetAnalysisManager()->GetTree();
-  Float_t xsection = 0;
-  Float_t ftrials  = 1;
-  if(tree){
-    TFile *curfile = tree->GetCurrentFile();
-    if (!curfile) {
-      Error("Notify","No current file");
-      return kFALSE;
-    }
-    if(!fh1Xsec||!fh1Trials){
-      Printf("%s%d No Histogram fh1Xsec",(char*)__FILE__,__LINE__);
-      return kFALSE;
-    }
-    AliAnalysisHelperJetTasks::PythiaInfoFromFile(curfile->GetName(),xsection,ftrials);
-    fh1Xsec->Fill("<#sigma>",xsection);
-    // construct average trials 
-    Float_t nEntries = (Float_t)tree->GetTree()->GetEntries();
-    if(ftrials>=nEntries)fAvgTrials = ftrials/nEntries; 
-  }  
-  return kTRUE;
+   // default constructor
 }
 
-//____________________________________________________________________
-void  AliAnalysisTaskJetChem::UserCreateOutputObjects()
+//__________________________________________________________________________________________
+AliAnalysisTaskJetChem::AliAnalysisTaskJetChem(const char *name) 
+  : AliAnalysisTaskFragmentationFunction(name)
+  ,fK0Type(0)  
+  ,fFilterMaskK0(0)
+  ,fListK0s(0)
+  ,fV0QAK0(0)
+  ,fFFHistosRecCutsK0Evt(0)      
+  ,fFFHistosIMK0AllEvt(0)        
+  ,fFFHistosIMK0Jet(0)           
+  ,fFFHistosIMK0Cone(0)
+  ,fFFHistosPhiCorrIMK0(0)
+  ,fFFIMNBinsJetPt(0)    
+  ,fFFIMJetPtMin(0) 
+  ,fFFIMJetPtMax(0)
+  ,fFFIMNBinsInvM(0) 
+  ,fFFIMInvMMin(0)   
+  ,fFFIMInvMMax(0)   
+  ,fFFIMNBinsPt(0)      
+  ,fFFIMPtMin(0)        
+  ,fFFIMPtMax(0)        
+  ,fFFIMNBinsXi(0)      
+  ,fFFIMXiMin(0)        
+  ,fFFIMXiMax(0)        
+  ,fFFIMNBinsZ(0)       
+  ,fFFIMZMin(0)         
+  ,fFFIMZMax(0)
+  ,fPhiCorrIMNBinsPt(0)
+  ,fPhiCorrIMPtMin(0)
+  ,fPhiCorrIMPtMax(0)
+  ,fPhiCorrIMNBinsPhi(0)
+  ,fPhiCorrIMPhiMin(0)
+  ,fPhiCorrIMPhiMax(0)
+  ,fPhiCorrIMNBinsInvM(0)
+  ,fPhiCorrIMInvMMin(0)
+  ,fPhiCorrIMInvMMax(0)
+  ,fh1K0Mult(0)
+  ,fh1dPhiJetK0(0)     
 {
-  // Create the output container
-  //
-  AliInfo("UserCreateOutPutObjects()");
-  //
-  //  Histograms
-
-  OpenFile(1);
-  CreateHistos();
-  PostData(1, fListOfHistos); // PostData at least once
-
+  // constructor
+  
+  DefineOutput(1,TList::Class());  
 }
 
-//____________________________________________________________________
-void  AliAnalysisTaskJetChem::UserExec(Option_t */*option*/)
+//__________________________________________________________________________________________________________________________
+AliAnalysisTaskJetChem::AliAnalysisTaskJetChem(const  AliAnalysisTaskJetChem &copy)
+  : AliAnalysisTaskFragmentationFunction()
+  ,fK0Type(copy.fK0Type)
+  ,fFilterMaskK0(copy.fFilterMaskK0)
+  ,fListK0s(copy.fListK0s)
+  ,fV0QAK0(copy.fV0QAK0)
+  ,fFFHistosRecCutsK0Evt(copy.fFFHistosRecCutsK0Evt)      
+  ,fFFHistosIMK0AllEvt(copy.fFFHistosIMK0AllEvt)        
+  ,fFFHistosIMK0Jet(copy.fFFHistosIMK0Jet)           
+  ,fFFHistosIMK0Cone(copy.fFFHistosIMK0Cone)          
+  ,fFFHistosPhiCorrIMK0(copy.fFFHistosPhiCorrIMK0)          
+  ,fFFIMNBinsJetPt(copy.fFFIMNBinsJetPt)   
+  ,fFFIMJetPtMin(copy.fFFIMJetPtMin)     
+  ,fFFIMJetPtMax(copy.fFFIMJetPtMax)     
+  ,fFFIMNBinsInvM(copy.fFFIMNBinsInvM)  
+  ,fFFIMInvMMin(copy.fFFIMInvMMin)    
+  ,fFFIMInvMMax(copy.fFFIMInvMMax)    
+  ,fFFIMNBinsPt(copy.fFFIMNBinsPt)      
+  ,fFFIMPtMin(copy.fFFIMPtMin)        
+  ,fFFIMPtMax(copy.fFFIMPtMax)        
+  ,fFFIMNBinsXi(copy.fFFIMNBinsXi)      
+  ,fFFIMXiMin(copy.fFFIMXiMin)        
+  ,fFFIMXiMax(copy.fFFIMXiMax)        
+  ,fFFIMNBinsZ(copy.fFFIMNBinsZ)       
+  ,fFFIMZMin(copy.fFFIMZMin)         
+  ,fFFIMZMax(copy.fFFIMZMax) 
+  ,fPhiCorrIMNBinsPt(copy.fPhiCorrIMNBinsPt)
+  ,fPhiCorrIMPtMin(copy.fPhiCorrIMPtMin)
+  ,fPhiCorrIMPtMax(copy.fPhiCorrIMPtMax)
+  ,fPhiCorrIMNBinsPhi(copy.fPhiCorrIMNBinsPhi)
+  ,fPhiCorrIMPhiMin(copy.fPhiCorrIMPhiMin)
+  ,fPhiCorrIMPhiMax(copy.fPhiCorrIMPhiMax)
+  ,fPhiCorrIMNBinsInvM(copy.fPhiCorrIMNBinsInvM)
+  ,fPhiCorrIMInvMMin(copy.fPhiCorrIMInvMMin)
+  ,fPhiCorrIMInvMMax(copy.fPhiCorrIMInvMMax)
+  ,fh1K0Mult(copy.fh1K0Mult)
+  ,fh1dPhiJetK0(copy.fh1dPhiJetK0)
 {
-  // Execute analysis for current event
-  //
-  if (fDebug > 3) AliInfo( " Processing event..." );
-
-  TObject* handler = AliAnalysisManager::GetAnalysisManager()->GetInputEventHandler();
+  // copy constructor
   
-  if( handler && handler->InheritsFrom("AliAODInputHandler") ) {
-    fAOD  =  ((AliAODInputHandler*)handler)->GetEvent();
-    if (fDebug > 1) AliInfo("  ==== Tracks from AliAODInputHandler");
-    // Case when jets are reconstructed on the fly from AOD tracks
-    // (the Jet Finder is using the AliJetAODReader) of InputEventHandler
-    // and put in the OutputEventHandler AOD. Useful whe you want to reconstruct jets with
-    // different parameters to default ones stored in the AOD or to use a different algorithm
-    if( fJetsOnFly ) {
-      handler = AliAnalysisManager::GetAnalysisManager()->GetOutputEventHandler();
-      if( handler && handler->InheritsFrom("AliAODHandler") ) {
-        fAODjets = ((AliAODHandler*)handler)->GetAOD();
-        if (fDebug > 1) AliInfo("  ==== Jets from AliAODHandler");
-      }
-    } else {
-      fAODjets = fAOD;
-      if (fDebug > 1) AliInfo("  ==== Jets from AliAODInputHandler");
-    }
-  } else {
-    handler = AliAnalysisManager::GetAnalysisManager()->GetOutputEventHandler();
-    if( handler && handler->InheritsFrom("AliAODHandler") ) {
-      fAOD  = ((AliAODHandler*)handler)->GetAOD();
-      fAODjets = fAOD;
-      if (fDebug > 1) AliInfo("  ==== Tracks and Jets from AliAODHandler");
-    } else {
-      AliFatal("I can't get any AOD Event Handler");
-      return;
-    }
-  }
+}
 
-  // -------------
-
-  // fetch the pythia header info and get the trials
-  AliMCEventHandler* mcHandler = dynamic_cast<AliMCEventHandler*> (AliAnalysisManager::GetAnalysisManager()->GetMCtruthEventHandler());
-  Float_t nTrials = 1;
-  if (mcHandler) {  
-    AliMCEvent* mcEvent = mcHandler->MCEvent();
-    if (mcEvent) {
-      AliGenPythiaEventHeader*  pythiaGenHeader = AliAnalysisHelperJetTasks::GetPythiaEventHeader(mcEvent);
-      if(pythiaGenHeader){
-       nTrials = pythiaGenHeader->Trials();
-      }
-    }
+// _________________________________________________________________________________________________________________________________
+AliAnalysisTaskJetChem& AliAnalysisTaskJetChem::operator=(const AliAnalysisTaskJetChem& o)
+{
+  // assignment
+  
+  if(this!=&o){
+    AliAnalysisTaskFragmentationFunction::operator=(o);
+
+    fK0Type                         = o.fK0Type;
+    fFilterMaskK0                   = o.fFilterMaskK0;
+    fListK0s                        = o.fListK0s;
+    fV0QAK0                         = o.fV0QAK0;
+    fFFHistosRecCutsK0Evt           = o.fFFHistosRecCutsK0Evt;      
+    fFFHistosIMK0AllEvt             = o.fFFHistosIMK0AllEvt;        
+    fFFHistosIMK0Jet                = o.fFFHistosIMK0Jet;           
+    fFFHistosIMK0Cone               = o.fFFHistosIMK0Cone;          
+    fFFHistosPhiCorrIMK0            = o.fFFHistosPhiCorrIMK0;
+    fFFIMNBinsJetPt                 = o.fFFIMNBinsJetPt;    
+    fFFIMJetPtMin                   = o.fFFIMJetPtMin; 
+    fFFIMJetPtMax                   = o.fFFIMJetPtMax;
+    fFFIMNBinsPt                    = o.fFFIMNBinsPt;      
+    fFFIMPtMin                      = o.fFFIMPtMin;        
+    fFFIMPtMax                      = o.fFFIMPtMax;        
+    fFFIMNBinsXi                    = o.fFFIMNBinsXi;      
+    fFFIMXiMin                      = o.fFFIMXiMin;        
+    fFFIMXiMax                      = o.fFFIMXiMax;        
+    fFFIMNBinsZ                     = o.fFFIMNBinsZ;       
+    fFFIMZMin                       = o.fFFIMZMin;         
+    fFFIMZMax                       = o.fFFIMZMax;         
+    fPhiCorrIMNBinsPt               = o.fPhiCorrIMNBinsPt;
+    fPhiCorrIMPtMin                 = o.fPhiCorrIMPtMin;
+    fPhiCorrIMPtMax                 = o.fPhiCorrIMPtMax;
+    fPhiCorrIMNBinsPhi              = o.fPhiCorrIMNBinsPhi;
+    fPhiCorrIMPhiMin                = o.fPhiCorrIMPhiMin;
+    fPhiCorrIMPhiMax                = o.fPhiCorrIMPhiMax;
+    fPhiCorrIMNBinsInvM             = o.fPhiCorrIMNBinsInvM;
+    fPhiCorrIMInvMMin               = o.fPhiCorrIMInvMMin;
+    fPhiCorrIMInvMMax               = o.fPhiCorrIMInvMMax;
+    fh1K0Mult                       = o.fh1K0Mult;
+    fh1dPhiJetK0                    = o.fh1dPhiJetK0;
   }
-
-  fh1Trials->Fill("#sum{ntrials}",fAvgTrials); 
     
-  AnalyseEvent();
-  
-  // Post the data
-  PostData(1, fListOfHistos);
+  return *this;
 }
 
-//____________________________________________________________________
-void  AliAnalysisTaskJetChem::AnalyseEvent()
+//_______________________________________________
+AliAnalysisTaskJetChem::~AliAnalysisTaskJetChem()
 {
-  // check Trigger 
-
-  // TString firedTriggerClasses = fAOD->GetFiredTriggerClasses(); 
-  // AliInfo(Form("firedTriggerClasses %s",firedTriggerClasses.Data()));
-  // if(firedTriggerClasses.Length() > 0 && !firedTriggerClasses.Contains("CINT1B")) return;  
+  // destructor  
 
-  AliAODVertex *primVertex = fAOD->GetPrimaryVertex();  // is this from tracks or from SPD ? SPD should only be used if no track vertex
-  if(!primVertex){
-    AliInfo("no prim Vertex found - skip event");
-    fhPrimVertexNCont->Fill(-1);
-    return; 
-  }
-    
-  TString primVtxName(primVertex->GetName());
-  
-  if(primVtxName.CompareTo("TPCVertex",TString::kIgnoreCase) == 1){
-
-    AliInfo("found TPC prim Vertex  - skip event");
-    fhPrimVertexNCont->Fill(-1);
-    return; 
-  }
-    
-  Int_t nVertContributors = primVertex->GetNContributors();
-  fhPrimVertexNCont->Fill(nVertContributors);
-  
-  //cout<<" prim vertex name "<<primVertex->GetName()<<" nCont "<<nVertContributors<<endl;
-  
-  if(nVertContributors<1){ // eventually check if not SPD vertex ??? 
-    AliInfo("prim vertex no contributors - skip event");
-    return;
-  }
-    
-  Double_t vertX = primVertex->GetX();
-  Double_t vertY = primVertex->GetY();
-  Double_t vertZ = primVertex->GetZ();
-  
-  Double_t vertRho = TMath::Sqrt(vertX*vertX+vertY*vertY);
-  
-  fhPrimVertexRho->Fill(vertRho);
-  fhPrimVertexZ->Fill(vertZ);
-  
-  if(TMath::Abs(vertZ)>10){
-    AliInfo(Form("prim vertex z=%f - skip event",vertZ));
-    return; 
-  }
-  
-  
-  Int_t pythiaPID = GetPythiaProcessID();
-  
-  // ------------------------------------------------
-  // Find Leading Jets 1,2,3 
-  // (could be skipped if Jets are sort by Pt...)
-
-  //Int_t    index1 = -1;
-  //Int_t nTracksLeading = 0;
-
-  Int_t nJetsAOD = 0;
-  AliAODJet* leadingJetAOD = NULL; // non-zero if any jet in acc and leading jet survives pt cut 
-  Int_t indexLeadingAOD   = -1;
-  Double_t ptLeadingAOD   = 0.; 
-  Int_t indexMaxRegionAOD = 0; // initialize with '0', '+-1' transverse regions 
-
-  Int_t nJetsMC = 0;
-  AliAODJet* leadingJetMC = NULL; // non-zero if leading jet survives acc, pt cut 
-  Int_t indexLeadingMC   = -1;
-  Double_t ptLeadingMC   = 0.; 
-  Int_t indexMaxRegionMC = 0; // initialize with '0', '+-1' transverse regions 
-
-  Double_t   ptLeadingAODAllEta  = 0; 
-  AliAODJet* leadingJetAODAllEta = NULL; 
-  
-  if(fUseLOConeJets)   fArrayJetsAOD  = FindChargedParticleJets(); 
-  else{ // Use jets on AOD/DeltaAOD
-
-    if(fDeltaAOD){
-      if (fDebug > 1) AliInfo(" ==== Jets From  Delta-AODs !");
-      if (fDebug > 1) AliInfo(Form(" ====  Reading Branch: %s  ",fDeltaAODBranch.Data()));
-      fArrayJetsAOD = (TClonesArray*)fAODjets->GetList()->FindObject(fDeltaAODBranch.Data());
-      if (!fArrayJetsAOD){
-       AliFatal(" No jet-array! ");
-       return;
-      }
-    }
-    else{
-      if (fDebug > 1) AliInfo(" ==== AOD jets: Read Standard-AODs  !");
-      if (fDebug > 1) AliInfo(Form(" ====  Reading Branch: %s  ",fAODBranch.Data()));
-      
-      fArrayJetsAOD = (TClonesArray*)fAODjets->FindListObject(fAODBranch.Data());
-    }
-  }
   
+  if(fTracksRecCuts) delete fTracksRecCuts;
+  if(fJetsRecCuts) delete fJetsRecCuts;
+  if(fListK0s) delete fListK0s;
+  if(fCommonHistList) delete fCommonHistList;
 
-  if(fUseLOConeMCJets)    fArrayJetsMC = FindChargedParticleJetsMC(); 
-  else if(fUsePythiaJets) fArrayJetsMC = GetPythiaJets();
-  else if( fAODBranchMC.Length() || fDeltaAODBranchMC.Length() ){
-
-    if(fDeltaAOD){
-      if (fDebug > 1) AliInfo(" ==== MC Jets From  Delta-AODs !");
-      if (fDebug > 1) AliInfo(Form(" ====  Reading Branch: %s  ",fDeltaAODBranchMC.Data()));
-      fArrayJetsMC = (TClonesArray*)fAODjets->GetList()->FindObject(fDeltaAODBranchMC.Data());
-      if (!fArrayJetsMC){
-       AliFatal(" No MC jet-array! ");
-       return;
-      }
-    }
-    else{
-      if (fDebug > 1) AliInfo(" ==== MC jets: Read Standard-AODs  !");
-      if (fDebug > 1) AliInfo(Form(" ====  Reading Branch: %s  ",fAODBranchMC.Data()));
-      
-      fArrayJetsMC = (TClonesArray*)fAODjets->FindListObject(fAODBranchMC.Data());
-    }
-  }
-
-  if(fArrayJetsAOD) nJetsAOD = fArrayJetsAOD->GetEntries();
-  if(fArrayJetsMC)  nJetsMC  = fArrayJetsMC->GetEntries();
-
-  fhNJets->Fill(nJetsAOD);
-  fhNJetsMC->Fill(nJetsMC);
-  fhnJetsAODvsMC->Fill(nJetsMC,nJetsAOD);
+}
 
-  if(fDebug>1) AliInfo(Form("AOD %d jets",nJetsAOD)); 
+//________________________________________________________________________________________________________________________________
+AliAnalysisTaskJetChem::AliFragFuncHistosInvMass::AliFragFuncHistosInvMass(const char* name, 
+                                                                          Int_t nJetPt, Float_t jetPtMin, Float_t jetPtMax,  
+                                                                          Int_t nInvMass, Float_t invMassMin, Float_t invMassMax,
+                                                                          Int_t nPt, Float_t ptMin, Float_t ptMax,
+                                                                          Int_t nXi, Float_t xiMin, Float_t xiMax,
+                                                                          Int_t nZ , Float_t zMin , Float_t zMax )
+  : TObject()
+  ,fNBinsJetPt(nJetPt)
+  ,fJetPtMin(jetPtMin)
+  ,fJetPtMax(jetPtMax)
+  ,fNBinsInvMass(nInvMass)
+  ,fInvMassMin(invMassMin)  
+  ,fInvMassMax(invMassMax)
+  ,fNBinsPt(nPt) 
+  ,fPtMin(ptMin)   
+  ,fPtMax(ptMax)   
+  ,fNBinsXi(nXi) 
+  ,fXiMin(xiMin)   
+  ,fXiMax(xiMax)   
+  ,fNBinsZ(nZ)  
+  ,fZMin(zMin)    
+  ,fZMax(zMax)    
+  ,fh3TrackPt(0)
+  ,fh3Xi(0)
+  ,fh3Z(0)
+  ,fh1JetPt(0)
+  ,fNameFF(name)
+{
+  // default constructor
 
+}
 
-  // for xcheck: find leading jet in large eta range
-  for(Int_t i=0; i<nJetsAOD; i++){
+//______________________________________________________________________________________________________________
+AliAnalysisTaskJetChem::AliFragFuncHistosInvMass::AliFragFuncHistosInvMass(const AliFragFuncHistosInvMass& copy)
+  : TObject()
+  ,fNBinsJetPt(copy.fNBinsJetPt)
+  ,fJetPtMin(copy.fJetPtMin)
+  ,fJetPtMax(copy.fJetPtMax)
+  ,fNBinsInvMass(copy.fNBinsInvMass)
+  ,fInvMassMin(copy.fInvMassMin)  
+  ,fInvMassMax(copy.fInvMassMax)
+  ,fNBinsPt(copy.fNBinsPt) 
+  ,fPtMin(copy.fPtMin)   
+  ,fPtMax(copy.fPtMax)   
+  ,fNBinsXi(copy.fNBinsXi) 
+  ,fXiMin(copy.fXiMin)   
+  ,fXiMax(copy.fXiMax)   
+  ,fNBinsZ(copy.fNBinsZ)  
+  ,fZMin(copy.fZMin)    
+  ,fZMax(copy.fZMax)    
+  ,fh3TrackPt(copy.fh3TrackPt)
+  ,fh3Xi(copy.fh3Xi)
+  ,fh3Z(copy.fh3Z)
+  ,fh1JetPt(copy.fh1JetPt)
+  ,fNameFF(copy.fNameFF)
+{
+  // copy constructor
+}
 
-    AliAODJet* jet  = (AliAODJet*) fArrayJetsAOD->At(i);
-    Double_t jetPt  = jet->Pt();
-    
-    if(jetPt > ptLeadingAODAllEta){ 
-      ptLeadingAODAllEta  = jetPt;
-      leadingJetAODAllEta = jet;
-    }
+//______________________________________________________________________________________________________________________________________________________________________
+AliAnalysisTaskJetChem::AliFragFuncHistosInvMass& AliAnalysisTaskJetChem::AliFragFuncHistosInvMass::operator=(const AliAnalysisTaskJetChem::AliFragFuncHistosInvMass& o)
+{
+  // assignment
+  
+  if(this!=&o){
+    TObject::operator=(o);
+    fNBinsJetPt   = o.fNBinsJetPt;
+    fJetPtMin     = o.fJetPtMin;
+    fJetPtMax     = o.fJetPtMax;
+    fNBinsInvMass = o.fNBinsInvMass;
+    fInvMassMin   = o.fInvMassMin;  
+    fInvMassMax   = o.fInvMassMax;
+    fNBinsPt      = o.fNBinsPt; 
+    fPtMin        = o.fPtMin;   
+    fPtMax        = o.fPtMax;   
+    fNBinsXi      = o.fNBinsXi; 
+    fXiMin        = o.fXiMin;   
+    fXiMax        = o.fXiMax;   
+    fNBinsZ       = o.fNBinsZ;  
+    fZMin         = o.fZMin;    
+    fZMax         = o.fZMax;    
+    fh3TrackPt    = o.fh3TrackPt;
+    fh3Xi         = o.fh3Xi;
+    fh3Z          = o.fh3Z;
+    fh1JetPt      = o.fh1JetPt;
+    fNameFF       = o.fNameFF;
   }
-
-  if(leadingJetAODAllEta){
-    //cout<<" trackRefs entries "<<leadingJetAODAllEta->GetRefTracks()->GetEntriesFast()<<endl;
-    fhLeadingNTracksVsEta->Fill(leadingJetAODAllEta->Eta(),leadingJetAODAllEta->GetRefTracks()->GetEntriesFast());
-    fhLeadingPtVsEta->Fill(leadingJetAODAllEta->Eta(),leadingJetAODAllEta->Pt()); 
     
-  }
-
-  // find leading jet AOD
-  for(Int_t i=0; i<nJetsAOD; ++i){
-
-    AliAODJet* jet  = (AliAODJet*) fArrayJetsAOD->At(i); 
-    Double_t jetPt  = jet->Pt();
-    Double_t jetEta = jet->Eta();
-
-    if((jetPt > ptLeadingAOD) && (TMath::Abs(jetEta)<fJetEtaCut)){ 
-      ptLeadingAOD    = jetPt; 
-      indexLeadingAOD = i;
-    }
-  }
-
-  // find leading jet MC
-  for(Int_t i=0; i<nJetsMC; ++i){
+  return *this;
+}
 
-    AliAODJet* jet  = (AliAODJet*) fArrayJetsMC->At(i); 
-    Double_t jetPt  = jet->Pt();
-    Double_t jetEta = jet->Eta();
+//___________________________________________________________________________
+AliAnalysisTaskJetChem::AliFragFuncHistosInvMass::~AliFragFuncHistosInvMass()
+{ 
+  // destructor 
 
-    if((jetPt > ptLeadingMC) && (TMath::Abs(jetEta)<fJetEtaCut)){
-      ptLeadingMC    = jetPt; 
-      indexLeadingMC = i;
-    }
-  }
+  if(fh1JetPt)   delete fh1JetPt;
+  if(fh3TrackPt) delete fh3TrackPt;
+  if(fh3Xi)      delete fh3Xi;
+  if(fh3Z)       delete fh3Z;
+}
 
-  //cout<<" ptLeadingAOD "<<ptLeadingAOD<<" MC "<<ptLeadingMC<<endl;
+//_________________________________________________________________
+void AliAnalysisTaskJetChem::AliFragFuncHistosInvMass::DefineHistos()
+{
+  // book FF histos
 
-  if(indexLeadingAOD>=0){ // event with jet  
+  fh1JetPt   = new TH1F(Form("fh1FFJetPtIM%s", fNameFF.Data()),"",fNBinsJetPt,fJetPtMin,fJetPtMax);
+  fh3TrackPt = new TH3F(Form("fh3FFTrackPtIM%s",fNameFF.Data()),"",fNBinsJetPt, fJetPtMin, fJetPtMax, fNBinsInvMass, fInvMassMin, fInvMassMax, fNBinsPt, fPtMin, fPtMax);
+  fh3Xi      = new TH3F(Form("fh3FFXiIM%s", fNameFF.Data()),"", fNBinsJetPt, fJetPtMin, fJetPtMax, fNBinsInvMass, fInvMassMin, fInvMassMax, fNBinsXi, fXiMin, fXiMax);
+  fh3Z       = new TH3F(Form("fh3FFZIM%s", fNameFF.Data()),"", fNBinsJetPt, fJetPtMin, fJetPtMax, fNBinsInvMass, fInvMassMin, fInvMassMax, fNBinsZ, fZMin, fZMax);
 
-    Double_t etaLeadingAOD  = ((AliAODJet*) fArrayJetsAOD->At(indexLeadingAOD))->Eta();
-    Double_t phiLeadingAOD  = ((AliAODJet*) fArrayJetsAOD->At(indexLeadingAOD))->Phi();
-    Int_t nTracksLeadingAOD = ((AliAODJet*) fArrayJetsAOD->At(indexLeadingAOD))->GetRefTracks()->GetEntriesFast();
-    
-    if(fDebug>1) AliInfo(Form("\n Pt Leading AOD Jet = %6.1f eta=%5.3f, nTracks %d",ptLeadingAOD,etaLeadingAOD,nTracksLeadingAOD));
+  AliAnalysisTaskFragmentationFunction::SetProperties(fh1JetPt, "p_{t} (GeV/c)", "entries"); 
+  AliAnalysisTaskJetChem::SetProperties(fh3TrackPt,"jet p_{t} (GeV/c)","inv Mass (GeV/c^2)","p_{t} (GeV/c)");
+  AliAnalysisTaskJetChem::SetProperties(fh3Xi,"jet p_{t} (GeV/c)","inv Mass (GeV/c^2)","#xi");
+  AliAnalysisTaskJetChem::SetProperties(fh3Z,"jet p_{t} (GeV/c)","inv Mass (GeV/c^2)","z");
+}
 
-    fhLeadingEta->Fill(etaLeadingAOD);
+//________________________________________________________________________________________________________________________________
+void AliAnalysisTaskJetChem::AliFragFuncHistosInvMass::FillFF(Float_t trackPt, Float_t invM, Float_t jetPt, Bool_t incrementJetPt)
+{
+  // fill FF
  
-    if(TMath::Abs(etaLeadingAOD)<fJetEtaCut){ // leading jet eta cut
-
-      fhnTracksVsPtLeading->Fill(ptLeadingAOD,nTracksLeadingAOD);
-      fhLeadingPt->Fill(ptLeadingAOD);
-      if(IsDiffractiveEvent(pythiaPID))  fhLeadingPtDiffr->Fill(ptLeadingAOD);
-
-      if(ptLeadingAOD>fJetPtCut){ // leading jet pt cut
-       
-       leadingJetAOD = (AliAODJet*) fArrayJetsAOD->At(indexLeadingAOD);
-       
-       fhLeadingPhi->Fill(phiLeadingAOD);
-       
-       // ----------------------------------------------
-       // Find max and min regions
-       Double_t sumPtRegionPosit  = 0.;
-       Double_t sumPtRegionNegat  = 0.;
-       Int_t    nTrackRegionPosit = 0;
-       Int_t    nTrackRegionNegat = 0;
-      
-       Int_t nTracks = fAOD->GetNTracks();
-         
-       for (Int_t ipart=0; ipart<nTracks; ++ipart) {
-           
-         AliAODTrack* part = fAOD->GetTrack(ipart);
-         if ( !part->TestFilterBit(fFilterBit) ) continue; // track cut selection
-         if (!part->IsPrimaryCandidate()) continue; // reject whatever is not linked to collision point
-         // PID Selection: Reject everything but hadrons
-         Bool_t isHadron = part->GetMostProbablePID()==AliAODTrack::kPion || 
-           part->GetMostProbablePID()==AliAODTrack::kKaon || 
-           part->GetMostProbablePID()==AliAODTrack::kProton;
-         if (!isHadron ) continue;
-         if (!part->Charge() ) continue; //Only charged
-         if (part->Pt() < fTrackPtCut ) continue;
-         if(TMath::Abs(part->Eta()) > fTrackEtaCut ) continue;
-         
-         TVector3 partVect(part->Px(), part->Py(), part->Pz());
-           
-         Int_t region = IsTrackInsideRegion(leadingJetAOD,&partVect );  
-
-         if (region > 0) {
-           sumPtRegionPosit += part->Pt();
-           nTrackRegionPosit++;
-         }
-         if (region < 0) {
-           sumPtRegionNegat += part->Pt();
-           nTrackRegionNegat++;
-         }
-       } // tracks loop 
-
-       // fill sumPt and mult density for transverse regions
-       
-       if( sumPtRegionPosit > sumPtRegionNegat ) {
-         FillSumPtRegion( ptLeadingAOD, sumPtRegionPosit/fAreaReg, sumPtRegionNegat/fAreaReg );
-       } 
-       else {
-         FillSumPtRegion( ptLeadingAOD, sumPtRegionNegat/fAreaReg, sumPtRegionPosit/fAreaReg );
-       }
-       if (nTrackRegionPosit > nTrackRegionNegat ) {
-         FillMultRegion( ptLeadingAOD, nTrackRegionPosit/fAreaReg, nTrackRegionNegat/fAreaReg);
-       } 
-       else {
-         FillMultRegion( ptLeadingAOD, nTrackRegionNegat/fAreaReg, nTrackRegionPosit/fAreaReg);
-       }
-        
-       indexMaxRegionAOD = (sumPtRegionPosit > sumPtRegionNegat) ? 1 : -1;
-       
-      } // leading jet pt cut
-    } // leading jet eta cut  
-  } // jet event
-  
-
-  if(indexLeadingMC>=0){ // event with MC jet  
-
-    Double_t etaLeadingMC  = ((AliAODJet*) fArrayJetsMC->At(indexLeadingMC))->Eta();
-    Double_t phiLeadingMC  = ((AliAODJet*) fArrayJetsMC->At(indexLeadingMC))->Phi();
-    Int_t nTracksLeadingMC = ((AliAODJet*) fArrayJetsMC->At(indexLeadingMC))->GetRefTracks()->GetEntriesFast();
-    
-    if(fDebug>1) AliInfo(Form("\n Pt Leading MC Jet = %6.1f eta=%5.3f, nTracks %d",ptLeadingMC,etaLeadingMC,nTracksLeadingMC));
-
-    fhLeadingEtaMC->Fill(etaLeadingMC);
+  if(incrementJetPt) fh1JetPt->Fill(jetPt);    
+  fh3TrackPt->Fill(jetPt,invM,trackPt);
   
-    if(TMath::Abs(etaLeadingMC)<fJetEtaCut){ // leading jet eta cut
-
-      fhLeadingPtMC->Fill(ptLeadingMC);
-      if(IsDiffractiveEvent(pythiaPID)) fhLeadingPtMCDiffr->Fill(ptLeadingMC);
-
-      if(ptLeadingMC>fJetPtCut){ // leading jet pt cut
-       
-       leadingJetMC = (AliAODJet*) fArrayJetsMC->At(indexLeadingMC);
-   
-       fhLeadingPhiMC->Fill(phiLeadingMC); // -pi to pi
-
-       // ----------------------------------------------
-       // Find max and min regions
-       Double_t sumPtRegionPosit  = 0;
-       Double_t sumPtRegionNegat  = 0;
-       Int_t    nTrackRegionPosit = 0;
-       Int_t    nTrackRegionNegat = 0;
-
-       TClonesArray* farray = (TClonesArray*)fAOD->FindListObject("mcparticles");
-       
-       Int_t ntrks = farray->GetEntries();
-       if (fDebug>1) AliInfo(Form("In UE MC analysis tracks %d \n",ntrks));
-         
-       for(Int_t i =0 ; i < ntrks; i++){   
-           
-         AliAODMCParticle* mctrk = (AliAODMCParticle*)farray->At(i);
-         //Cuts
-         if (!(mctrk->IsPhysicalPrimary())) continue;
-         //if (!(mctrk->IsPrimary())) continue;
-         
-         if (mctrk->Charge() == 0 || mctrk->Charge()==-99) continue;
-           
-         if (mctrk->Pt() < fTrackPtCut ) continue;
-         if( TMath::Abs(mctrk->Eta()) > fTrackEtaCut ) continue;
-         
-         Bool_t isHadron = TMath::Abs(mctrk->GetPdgCode())==211 ||
-           TMath::Abs(mctrk->GetPdgCode())==2212 ||
-           TMath::Abs(mctrk->GetPdgCode())==321;
-         
-         if (!isHadron) continue;
-           
-         TVector3 partVect(mctrk->Px(), mctrk->Py(), mctrk->Pz());
-         
-         Int_t region = IsTrackInsideRegion(leadingJetMC,&partVect );  
-         
-         if (region > 0) {
-           sumPtRegionPosit += mctrk->Pt();
-           nTrackRegionPosit++;
-         }
-         if (region < 0) {
-           sumPtRegionNegat += mctrk->Pt();
-           nTrackRegionNegat++;
-         }
-       } // AliAODMCParticle loop
-       
-       indexMaxRegionMC = (sumPtRegionPosit > sumPtRegionNegat) ? 1 : -1;
-       
-      } // leading jet pt
-    } // leading jet eta 
-  } // found jet
-
-  Bool_t foundK0AOD   = kFALSE; // event with K0 ? 
-  Bool_t foundK0MC    = kFALSE; // event with K0 ? 
-
-  CheckV0s(leadingJetAOD,indexMaxRegionAOD,foundK0AOD); // here leadingJetAOD/MC nonzero if jet passes eta & pt cut
-  CheckMCParticles(leadingJetMC,indexMaxRegionMC,foundK0MC);
-  CompLeadingJets(leadingJetAOD,leadingJetMC,pythiaPID,foundK0AOD,foundK0MC);
-  FillReferencePlotsTracks();
-  FillReferenceFF(leadingJetAOD);
-
+  Double_t z = 0.;
+  if(jetPt>0) z = trackPt / jetPt;
+  Double_t xi = 0;
+  if(z>0) xi = TMath::Log(1/z);
   
-
-  if(fUseLOConeJets && fArrayJetsAOD){
-    fArrayJetsAOD->Delete(); // no 'Clear': AliAODjet contains TMomentum and TRefArray
-    delete fArrayJetsAOD;
-  }
-  if(fUseLOConeMCJets && fArrayJetsMC){
-    fArrayJetsMC->Delete(); // no 'Clear': AliAODjet contains TMomentum and TRefArray
-    delete fArrayJetsMC;
-  }
+  fh3Xi->Fill(jetPt,invM,xi);
+  fh3Z->Fill(jetPt,invM,z);
 }
 
-// __________________________________________________________________
-
-
-Double_t AliAnalysisTaskJetChem::GetJetRadius(const AliAODJet* jet, const Double_t energyFrac){
-
-  // calc jet radius containing fraction energyFrac of full jet pt  
-  
-  const Int_t kArraySize = 1000;
-  const Int_t kInitVal   = -999;
-
-  Int_t nTracks = jet->GetRefTracks()->GetEntriesFast();
-
-  if(nTracks>kArraySize){
-    AliError(Form("nTracks in jet %d exceeds max array size",nTracks));
-    return -1;
-  }
-
+//___________________________________________________________________________________
+void AliAnalysisTaskJetChem::AliFragFuncHistosInvMass::AddToOutput(TList* list) const
+{
+  // add histos to list
 
-  Double_t deltaR[kArraySize] = {kInitVal};
-  Double_t pt[kArraySize] = {kInitVal};
-  Int_t index[kArraySize] = {kInitVal};
-  for(int i=0; i<kArraySize; i++) index[i] = kInitVal;
+  list->Add(fh1JetPt);
   
+  list->Add(fh3TrackPt);
+  list->Add(fh3Xi);
+  list->Add(fh3Z);
+}
 
-  Double_t ptTot = 0;
-
-  for(int i=0; i<nTracks; i++){
-
-    AliAODTrack* track = (AliAODTrack*) jet->GetRefTracks()->At(i);
-
-    TLorentzVector *mom4Jet   = jet->MomentumVector();
-    TVector3 mom3Jet = mom4Jet->Vect();
-
-    Double_t trackMom[3];
-    track->PxPyPz(trackMom);
-    TVector3 mom3Track(trackMom);
-
-    Double_t dR = mom3Jet.DeltaR(mom3Track);
-
-    deltaR[i] = dR;
-    pt[i]     = track->Pt();
-
-    ptTot += pt[i];
-  }
-
-  //cout<<" ptTot "<<ptTot<<" jetPt "<<jet->Pt()<<endl; // Xcheck 
-
-  TMath::Sort(nTracks,deltaR,index,kFALSE); // sort in decreasing order 
-
-  Double_t ptSum = 0;
-
-  for(int i=0; i<nTracks; i++){
-
-    Int_t ind = index[i];
-    
-    Double_t ptTrack = pt[ind];
-    Double_t dR      = deltaR[ind];
-
-    ptSum += ptTrack;
-
-    if(ptSum >= ptTot*energyFrac) return dR;
-  }
-
-  return -1;
+// ---
+
+
+//_______________________________________________________________________________________________________
+AliAnalysisTaskJetChem::AliFragFuncHistosPhiCorrInvMass::AliFragFuncHistosPhiCorrInvMass(const char* name,
+                                                                                        Int_t nPt,   Float_t ptMin,   Float_t ptMax,
+                                                                                        Int_t nPhi,  Float_t phiMin,  Float_t phiMax,
+                                                                                        Int_t nInvMass,  Float_t invMassMin,  Float_t invMassMax)
+  : TObject()
+  ,fNBinsPt(nPt)
+  ,fPtMin(ptMin)
+  ,fPtMax(ptMax)
+  ,fNBinsPhi(nPhi)
+  ,fPhiMin(phiMin)
+  ,fPhiMax(phiMax)
+  ,fNBinsInvMass(nInvMass)
+  ,fInvMassMin(invMassMin)  
+  ,fInvMassMax(invMassMax)
+  ,fh3PhiCorr(0) 
+  ,fNamePhiCorr(name) 
+{
+  // default constructor
+}
 
+//____________________________________________________________________________________________________________________________________
+AliAnalysisTaskJetChem::AliFragFuncHistosPhiCorrInvMass::AliFragFuncHistosPhiCorrInvMass(const AliFragFuncHistosPhiCorrInvMass& copy)
+  : TObject()
+  ,fNBinsPt(copy.fNBinsPt)
+  ,fPtMin(copy.fPtMin)
+  ,fPtMax(copy.fPtMax)
+  ,fNBinsPhi(copy.fNBinsPhi)
+  ,fPhiMin(copy.fPhiMin)
+  ,fPhiMax(copy.fPhiMax)
+  ,fNBinsInvMass(copy.fNBinsInvMass)
+  ,fInvMassMin(copy.fInvMassMin)  
+  ,fInvMassMax(copy.fInvMassMax)
+  ,fh3PhiCorr(copy.fh3PhiCorr)
+  ,fNamePhiCorr(copy.fNamePhiCorr)
+{
+  // copy constructor
 }
 
-//____________________________________________________________________
-void AliAnalysisTaskJetChem::FillSumPtRegion( Double_t leadingE, Double_t ptMax, Double_t ptMin  )
+//________________________________________________________________________________________________________________________________________________________________________
+AliAnalysisTaskJetChem::AliFragFuncHistosPhiCorrInvMass& AliAnalysisTaskJetChem::AliFragFuncHistosPhiCorrInvMass::operator=(const AliAnalysisTaskJetChem::AliFragFuncHistosPhiCorrInvMass& o)
 {
-  // Fill sumPt of control regions
+  // assignment
+  
+  if(this!=&o){
+    TObject::operator=(o);
+    fNBinsPt      = o.fNBinsPt;
+    fPtMin        = o.fPtMin;
+    fPtMax        = o.fPtMax;
+    fNBinsPhi     = o.fNBinsPhi;
+    fPhiMin       = o.fPhiMin;
+    fPhiMax       = o.fPhiMax;
+    fNBinsInvMass = o.fNBinsInvMass;
+    fInvMassMin   = o.fInvMassMin;  
+    fInvMassMax   = o.fInvMassMax;
+    
+    fh3PhiCorr    = o.fh3PhiCorr;
+    fNamePhiCorr  = o.fNamePhiCorr;
+  }
   
-  fhRegionSumPtMaxVsEt->Fill( leadingE, ptMax );
-  fhRegionSumPtMinVsEt->Fill( leadingE, ptMin );
+  return *this;
 }
 
-//____________________________________________________________________
-void AliAnalysisTaskJetChem::FillMultRegion(Double_t leadingE, Double_t nTrackPtmax, Double_t nTrackPtmin)
+//_________________________________________________________________________________________
+AliAnalysisTaskJetChem::AliFragFuncHistosPhiCorrInvMass::~AliFragFuncHistosPhiCorrInvMass()
 {
-  // Fill Nch multiplicity of control regions
+  // destructor 
   
-  fhRegionMultMaxVsEt->Fill( leadingE, nTrackPtmax );
-  fhRegionMultMinVsEt->Fill( leadingE, nTrackPtmin );
-
+  if(fh3PhiCorr) delete fh3PhiCorr;
 }
 
-//____________________________________________________________________
-Int_t AliAnalysisTaskJetChem::IsTrackInsideRegion(const AliAODJet* aodjetVect,const TVector3 *partVect) 
-{  
-  // return de region in delta phi
-  // -1 negative delta phi 
-  //  1 positive delta phi
-  //  0 outside region
-
-  TLorentzVector* jetVectLorentz = aodjetVect->MomentumVector(); 
-  TVector3 jetVect = jetVectLorentz->Vect();
+//__________________________________________________________________________
+void AliAnalysisTaskJetChem::AliFragFuncHistosPhiCorrInvMass::DefineHistos()
+{
+  // book jet QA histos
 
-  static const Double_t k60rad  = 60.*TMath::Pi()/180.;
-  static const Double_t k120rad = 120.*TMath::Pi()/180.;
+  fh3PhiCorr  = new TH3F(Form("fh3PhiCorrIM%s", fNamePhiCorr.Data()), 
+                        Form("%s: p_{t} - #phi - m_{inv} distribution",fNamePhiCorr.Data()), 
+                        fNBinsPt, fPtMin, fPtMax, 
+                        fNBinsPhi, fPhiMin, fPhiMax,
+                        fNBinsInvMass, fInvMassMin, fInvMassMax);
   
-  Int_t region = 0;
-  if( TMath::Abs(partVect->Eta()) > fTrackEtaCut ) return 0;
-  // transverse regions
-  if (jetVect.DeltaPhi(*partVect) < -k60rad && jetVect.DeltaPhi(*partVect) > -k120rad ) region = -1;
-  if (jetVect.DeltaPhi(*partVect) > k60rad && jetVect.DeltaPhi(*partVect) < k120rad ) region = 1;
-    
-  return region;
+  AliAnalysisTaskJetChem::SetProperties(fh3PhiCorr, "p_{t} (GeV/c)", "#phi", "m_{inv} (GeV/c^2)"); 
 }
 
-//____________________________________________________________________
-TClonesArray* AliAnalysisTaskJetChem::FindChargedParticleJetsMC()
+//___________________________________________________________________________________________________________
+void AliAnalysisTaskJetChem::AliFragFuncHistosPhiCorrInvMass::FillPhiCorr(Float_t pt, Float_t phi, Float_t invM)
 {
-  // CDF jet finder: 
-  // loop over pt-ordered list of tracks, combine tracks within jet cone, recalc cone axis after each step
-  // based on implementation by Arian Abrahantes Quintana and Enesto Lopez
-  TClonesArray* farray = (TClonesArray*)fAOD->FindListObject("mcparticles");
-  if(!farray){
-    AliInfo("no mcparticles branch"); 
-    return 0;
-  }
-
-  Int_t nTracks = farray->GetEntries();
-
-  if( !nTracks ) return 0;
-  TObjArray tracks(nTracks);
+  // fill jet QA histos 
 
-  for (Int_t ipart=0; ipart<nTracks; ++ipart) {
+  fh3PhiCorr->Fill(pt, phi, invM);
+}
 
-    AliAODMCParticle* part = (AliAODMCParticle*)farray->At(ipart);
+//______________________________________________________________________________________________
+void AliAnalysisTaskJetChem::AliFragFuncHistosPhiCorrInvMass::AddToOutput(TList* list) const 
+{
+  // add histos to list
 
-    if(!part->IsPhysicalPrimary()) continue;
+  list->Add(fh3PhiCorr);
+}
 
-    // exclude neutrinos
-    Int_t   pdg  = TMath::Abs(part->GetPdgCode());
-    if((pdg == 12 || pdg == 14 || pdg == 16)) return kFALSE;
+//____________________________________________________
+void AliAnalysisTaskJetChem::UserCreateOutputObjects()
+{
+  // create output objects
 
-    if( !part->Charge() ) continue; // comment / uncomment here
-    fhEtaMCTracks->Fill(part->Eta());
+  if(fDebug > 1) Printf("AliAnalysisTaskFragmentationFunction::UserCreateOutputObjects()");
+  // create list of tracks and jets 
+  
+  fTracksRecCuts = new TList();
+  fTracksRecCuts->SetOwner(kFALSE);  
 
-    if(TMath::Abs(part->Eta()) > fTrackEtaCut) continue; 
-    fhPtMCTracks->Fill(part->Pt());
+  fJetsRecCuts = new TList();
+  fJetsRecCuts->SetOwner(kFALSE);
 
-    if( part->Pt() < fTrackPtCutJF ) continue;
-    fhPhiMCTracks->Fill(part->Phi());
+  fListK0s = new TList(); 
+  fListK0s->SetOwner(kFALSE);
 
-    tracks.AddLast(part);
-  
-  }
+  //
+  // Create histograms / output container
+  //
 
-  QSortTracks( tracks, 0, tracks.GetEntriesFast() );
+  OpenFile(1);
+  fCommonHistList = new TList();
   
-  nTracks = tracks.GetEntriesFast();
-
-  if( !nTracks ) return 0;
-  TObjArray *jets = new TObjArray(nTracks);
-  TIter itrack(&tracks);
-  while( nTracks ) {
-    // 2- Start with the highest pT particle ...
-    Float_t px,py,pz,pt; 
-
-    AliAODMCParticle* track = (AliAODMCParticle*)itrack.Next();
-
-    if( !track ) continue;
-    px = track->Px();
-    py = track->Py();
-    pz = track->Pz();
-    pt = track->Pt(); 
-    jets->AddLast( new AliAODJet(px,py,pz,pt) ); // Use the energy member to store Pt
-    AliAODJet* jet = (AliAODJet*)jets->Last();
-    jet->AddTrack(track);
-    tracks.Remove(track);
-   
+  Bool_t oldStatus = TH1::AddDirectoryStatus();
+  TH1::AddDirectory(kFALSE);
+       
 
-    TVector2 jetVect2(jet->Px(),jet->Py());
+  // histograms inherited from AliAnalysisTaskFragmentationFunction
 
-    // 3- Go to the next highest pT particle not already included...
-    AliAODMCParticle* track1;
-    while ( (track1  = (AliAODMCParticle*)(itrack.Next())) ) {
-      TVector2 trackVect2(track1->Px(),track1->Py());
-      Double_t dphi = trackVect2.DeltaPhi(jetVect2);
+  fh1EvtSelection            = new TH1F("fh1EvtSelection", "Event Selection", 6, -0.5, 5.5);
+  fh1EvtSelection->GetXaxis()->SetBinLabel(1,"ACCEPTED");
+  fh1EvtSelection->GetXaxis()->SetBinLabel(2,"event selection: rejected");
+  fh1EvtSelection->GetXaxis()->SetBinLabel(3,"event class: rejected");
+  fh1EvtSelection->GetXaxis()->SetBinLabel(4,"vertex Ncontr: rejected");
+  fh1EvtSelection->GetXaxis()->SetBinLabel(5,"vertex z: rejected");
+  fh1EvtSelection->GetXaxis()->SetBinLabel(6,"vertex type: rejected");
 
-      Double_t r = TMath::Sqrt( (jet->Eta()-track1->Eta())*(jet->Eta()-track1->Eta()) +
-                               dphi*dphi );
+  fh1EvtCent                = new TH1F("fh1EvtCent","centrality",100,0.,100.);
+  fh1VertexNContributors     = new TH1F("fh1VertexNContributors", "Vertex N contributors", 11,-.5, 10.5);
+  fh1VertexZ                 = new TH1F("fh1VertexZ", "Vertex z distribution", 30, -15., 15.);
+  fh1Xsec                    = new TProfile("fh1Xsec","xsec from pyxsec.root",1,0,1);
+  fh1Xsec->GetXaxis()->SetBinLabel(1,"<#sigma>");
+  fh1Trials                  = new TH1F("fh1Trials","trials from pyxsec.root",1,0,1);
+  fh1Trials->GetXaxis()->SetBinLabel(1,"#sum{ntrials}");
+  fh1PtHard                  = new TH1F("fh1PtHard","PYTHIA Pt hard;p_{T,hard}",350,-.5,349.5);
+  fh1PtHardTrials            = new TH1F("fh1PtHardTrials","PYTHIA Pt hard weight with trials;p_{T,hard}",350,-.5,349.5);
 
-      if( r < fConeRadius ) {
-        Double_t fPt   = jet->E()+track1->Pt();  // Scalar sum of Pt
-        // recalculating the centroid
-        Double_t eta = jet->Eta()*jet->E()/fPt + track1->Eta()*track1->Pt()/fPt;
-       jetVect2.SetMagPhi(jet->E()/fPt,jetVect2.Phi());
-       trackVect2.SetMagPhi(track1->Pt()/fPt,trackVect2.Phi());
-       
-       TVector2 sumVect2 = jetVect2+trackVect2;
-       Double_t phi = sumVect2.Phi();
+  fh1nRecJetsCuts            = new TH1F("fh1nRecJetsCuts","reconstructed jets per event",10,-0.5,9.5);
 
-        //jet->SetPtEtaPhiE( 1., eta, phi, fPt );
-       ((TLorentzVector*) jet->MomentumVector())->SetPtEtaPhiE(1,eta,phi,fPt);
 
-       jet->AddTrack(track1);
-        tracks.Remove(track1);
-      }
-    }
-    
-    tracks.Compress();
-
-    nTracks = tracks.GetEntries();
-    //   4- Continue until all particles are in a jet.
-    itrack.Reset();
-  } // end while nTracks
-  
-  // Convert to AODjets....
-  Int_t njets = jets->GetEntriesFast();
+  // histograms jetChem proper 
 
-  TClonesArray* aodjets = new TClonesArray("AliAODJet",njets);
-  aodjets->SetOwner(kTRUE);
+  fh1EvtMult                = new TH1F("fh1EvtMult","multiplicity",1200,0.,12000.);
+  fh1K0Mult                 = new TH1F("fh1K0Mult","K0 multiplicity",500,0.,500.);
+  fh1dPhiJetK0               = new TH1F("fh1dPhiJetK0","",640,-1,5.4);
 
-  Int_t count = 0;
-  for(Int_t ijet=0; ijet<njets; ++ijet) {
-    AliAODJet* jet = (AliAODJet*)jets->At(ijet);
 
-    if (jet->E() < fJetPtCut) continue;
-    Float_t px, py,pz,en; // convert to 4-vector
-    px = jet->E() * TMath::Cos(jet->Phi());  // Pt * cos(phi)
-    py = jet->E() * TMath::Sin(jet->Phi());  // Pt * sin(phi)
-    pz = jet->E() / TMath::Tan(2.0 * TMath::ATan(TMath::Exp(-jet->Eta())));
-    en = TMath::Sqrt(px * px + py * py + pz * pz);
-    jet->SetPxPyPzE(px,py,pz,en);
 
-
-    TClonesArray &tmpaodjets = *aodjets;
-    new(tmpaodjets[count++]) AliAODJet(*jet);
-    //aodjets->AddLast( new AliAODJet(*jet));    
-    //aodjets->AddLast( new AliAODJet(px, py, pz, en) );
-  }
-  // Order jets according to their pT .
-  QSortTracks( *aodjets, 0, aodjets->GetEntriesFast() );
+  fFFHistosRecCuts          = new AliFragFuncHistos("RecCuts", fFFNBinsJetPt, fFFJetPtMin, fFFJetPtMax, 
+                                                    fFFNBinsPt, fFFPtMin, fFFPtMax, 
+                                                    fFFNBinsXi, fFFXiMin, fFFXiMax,  
+                                                    fFFNBinsZ , fFFZMin , fFFZMax);
   
-  // debug
-  //if (fDebug>3) AliInfo(Form(" %d Charged jets found\n",njets));
+  fV0QAK0                   = new AliFragFuncQATrackHistos("V0QAK0",fQATrackNBinsPt, fQATrackPtMin, fQATrackPtMax, 
+                                                           fQATrackNBinsEta, fQATrackEtaMin, fQATrackEtaMax,
+                                                           fQATrackNBinsPhi, fQATrackPhiMin, fQATrackPhiMax, 
+                                                           fQATrackHighPtThreshold);
   
-  jets->Delete(); // OB - should I cleanup or leave it to garbage collection ? 
-  delete jets;
-
-  return aodjets;
-}
-// ___________________________________________________________________
-
-Bool_t AliAnalysisTaskJetChem::IsTrackFromK0(const Int_t indexTrack){
-
-  // check wether track with index is from K0 (= V0 with proper inv mass)
-
-  TClonesArray* tracks = fAOD->GetTracks();
-
-  for(int i=0; i<fAOD->GetNumberOfV0s(); i++){ // loop over V0s
-    
-    AliAODv0* v0 = fAOD->GetV0(i);
-    
-    Bool_t isOnFly = v0->GetOnFlyStatus();
-
-    if( (fUseOnFlyV0s && !isOnFly) || (!fUseOnFlyV0s && isOnFly) ) continue;
+  fFFHistosRecCutsK0Evt      = new AliFragFuncHistos("RecCutsK0Evt", fFFNBinsJetPt, fFFJetPtMin, fFFJetPtMax, 
+                                                    fFFNBinsPt, fFFPtMin, fFFPtMax, 
+                                                    fFFNBinsXi, fFFXiMin, fFFXiMax,  
+                                                    fFFNBinsZ , fFFZMin , fFFZMax);
+  
+  
+  fFFHistosIMK0AllEvt        = new AliFragFuncHistosInvMass("K0AllEvt", fFFIMNBinsJetPt, fFFIMJetPtMin, fFFIMJetPtMax, 
+                                                           fFFIMNBinsInvM,fFFIMInvMMin,fFFIMInvMMax,
+                                                           fFFIMNBinsPt, fFFIMPtMin, fFFIMPtMax, 
+                                                           fFFIMNBinsXi, fFFIMXiMin, fFFIMXiMax,  
+                                                           fFFIMNBinsZ , fFFIMZMin , fFFIMZMax);
+  
+  fFFHistosIMK0Jet           = new AliFragFuncHistosInvMass("K0Jet", fFFIMNBinsJetPt, fFFIMJetPtMin, fFFIMJetPtMax, 
+                                                           fFFIMNBinsInvM,fFFIMInvMMin,fFFIMInvMMax,
+                                                           fFFIMNBinsPt, fFFIMPtMin, fFFIMPtMax, 
+                                                           fFFIMNBinsXi, fFFIMXiMin, fFFIMXiMax,  
+                                                           fFFIMNBinsZ , fFFIMZMin , fFFIMZMax);
+  
+  
+  fFFHistosIMK0Cone          = new AliFragFuncHistosInvMass("K0Cone", fFFIMNBinsJetPt, fFFIMJetPtMin, fFFIMJetPtMax, 
+                                                           fFFIMNBinsInvM,fFFIMInvMMin,fFFIMInvMMax,
+                                                           fFFIMNBinsPt, fFFIMPtMin, fFFIMPtMax, 
+                                                           fFFIMNBinsXi, fFFIMXiMin, fFFIMXiMax,  
+                                                           fFFIMNBinsZ , fFFIMZMin , fFFIMZMax);
+  
+  fFFHistosPhiCorrIMK0       = new AliFragFuncHistosPhiCorrInvMass("K0",fPhiCorrIMNBinsPt, fPhiCorrIMPtMin, fPhiCorrIMPtMax, 
+                                                                  fPhiCorrIMNBinsPhi, fPhiCorrIMPhiMin, fPhiCorrIMPhiMax,  
+                                                                  fPhiCorrIMNBinsInvM , fPhiCorrIMInvMMin , fPhiCorrIMInvMMax);
+  
+  
+  fV0QAK0->DefineHistos();
+  fFFHistosRecCuts->DefineHistos();
+  fFFHistosRecCutsK0Evt->DefineHistos();
+  fFFHistosIMK0AllEvt->DefineHistos();
+  fFFHistosIMK0Jet->DefineHistos();
+  fFFHistosIMK0Cone->DefineHistos();
+  fFFHistosPhiCorrIMK0->DefineHistos();
 
-    Double_t massK0 = v0->MassK0Short();
+  const Int_t saveLevel = 5;
+  if(saveLevel>0){
     
-    // FIXME: here should go cuts (at least V0 quality)
-
-    if(IsK0InvMass(massK0)){
-      
-      AliAODTrack *trackPos = (AliAODTrack *) (v0->GetSecondaryVtx()->GetDaughter(0));
-      AliAODTrack *trackNeg = (AliAODTrack *) (v0->GetSecondaryVtx()->GetDaughter(1));   
-
-      Int_t indexPos = tracks->IndexOf(trackPos);
-      Int_t indexNeg = tracks->IndexOf(trackNeg);
-     
-
-      if(indexPos == indexTrack){ return kTRUE;}
-      if(indexNeg == indexTrack){ return kTRUE;}
-
+    fCommonHistList->Add(fh1EvtSelection);
+    fCommonHistList->Add(fh1EvtCent);
+    fCommonHistList->Add(fh1VertexNContributors);
+    fCommonHistList->Add(fh1VertexZ);
+    fCommonHistList->Add(fh1Xsec);
+    fCommonHistList->Add(fh1Trials);
+    fCommonHistList->Add(fh1PtHard);
+    fCommonHistList->Add(fh1PtHardTrials);
+    fCommonHistList->Add(fh1nRecJetsCuts);
+    fCommonHistList->Add(fh1EvtMult);
+    fCommonHistList->Add(fh1K0Mult);
+    fCommonHistList->Add(fh1dPhiJetK0);
+
+    fV0QAK0->AddToOutput(fCommonHistList);
+    fFFHistosRecCuts->AddToOutput(fCommonHistList);
+    fFFHistosRecCutsK0Evt->AddToOutput(fCommonHistList);
+
+    fFFHistosIMK0AllEvt->AddToOutput(fCommonHistList);
+    fFFHistosIMK0Jet->AddToOutput(fCommonHistList);
+    fFFHistosIMK0Cone->AddToOutput(fCommonHistList);
+    fFFHistosPhiCorrIMK0->AddToOutput(fCommonHistList);
+  }
+  
+  // =========== Switch on Sumw2 for all histos ===========
+  for (Int_t i=0; i<fCommonHistList->GetEntries(); ++i){
+    TH1 *h1 = dynamic_cast<TH1*>(fCommonHistList->At(i));
+    if (h1) h1->Sumw2();
+    else{
+      THnSparse *hnSparse = dynamic_cast<THnSparse*>(fCommonHistList->At(i));
+      if(hnSparse) hnSparse->Sumw2();
     }
   }
-
-  return kFALSE;
-
+  
+  TH1::AddDirectory(oldStatus);
 }
 
-//____________________________________________________________________
-TClonesArray*  AliAnalysisTaskJetChem::FindChargedParticleJets()
+//_______________________________________________
+void AliAnalysisTaskJetChem::UserExec(Option_t *) 
 {
+  // Main loop
+  // Called for each event
 
-  // CDF jet finder: 
-  // loop over pt-ordered list of tracks, combine tracks within jet cone, recalc cone axis after each step 
-  // based on implementation by Arian Abrahantes Quintana and Enesto Lopez
-  // ref: PHYSICAL REVIEW D 65 092002, CDF Collaboration
-
-  //  1 - Order all charged particles according to their pT .
-  Int_t nTracks = fAOD->GetNTracks();
-
-  if( !nTracks ) return 0;
-  TObjArray tracks(nTracks);
-
-  // standard cuts + ITS refit = stdrd PWG4 cut (compose them for productions with old ESDFilter task) 
-
-  for (Int_t ipart=0; ipart<nTracks; ++ipart) {
-    AliAODTrack* part = fAOD->GetTrack( ipart );
+  if(fDebug > 1) Printf("Analysis event #%5d", (Int_t) fEntry);
 
-    UInt_t status = part->GetStatus();
-
-    if( !part->TestFilterBit(fFilterBitJF) ) continue; // track cut selection
-    if(fRequireITSRefitJF &&  ((status&AliESDtrack::kITSrefit)==0)) continue;  
-
-    if(!part->Charge() ) continue;
-
-    if(TMath::Abs(part->Eta()) > fTrackEtaCut) continue; 
-    
-    if(fRejectK0TracksJF && IsTrackFromK0(ipart)) continue; 
-
-    if( part->Pt() < fTrackPtCutJF ) continue;
-
-    tracks.AddLast(part);
+  AliInputEventHandler* inputHandler = (AliInputEventHandler*)
+    ((AliAnalysisManager::GetAnalysisManager())->GetInputEventHandler());
+  if(inputHandler->IsEventSelected() & AliVEvent::kMB){
+    if(fDebug > 1)  Printf(" Trigger Selection: event ACCEPTED ... ");
+    fh1EvtSelection->Fill(1.);
+  } else {
+    fh1EvtSelection->Fill(0.);
+    if(inputHandler->InheritsFrom("AliESDInputHandler") && fUsePhysicsSelection){ // PhysicsSelection only with ESD input
+      if (fDebug > 1 ) Printf(" Trigger Selection: event REJECTED ... ");
+      PostData(1, fCommonHistList);
+      return;
+    }
   }
 
-  QSortTracks( tracks, 0, tracks.GetEntriesFast() );
+  fESD = dynamic_cast<AliESDEvent*>(InputEvent());
+  if(!fESD){
+    if(fDebug>3) Printf("%s:%d ESDEvent not found in the input", (char*)__FILE__,__LINE__);
+  }
   
-  nTracks = tracks.GetEntriesFast();
+  fMCEvent = MCEvent();
+  if(!fMCEvent){
+    if(fDebug>3) Printf("%s:%d MCEvent not found in the input", (char*)__FILE__,__LINE__);
+  }
   
-  if( !nTracks ) return 0;
-  TObjArray *jets = new TObjArray(nTracks);
-  TIter itrack(&tracks);
-  while( nTracks ) {
-    // 2- Start with the highest pT particle ...
-    Float_t px,py,pz,pt; 
-    AliAODTrack* track = (AliAODTrack*)itrack.Next();
-    if( !track ) continue;
-    px = track->Px();
-    py = track->Py();
-    pz = track->Pz();
-    pt = track->Pt(); 
-    //jets->AddLast( new TLorentzVector(px, py, pz, pt) ); // Use the energy member to store Pt
-    jets->AddLast( new AliAODJet(px,py,pz,pt) ); // Use the energy member to store Pt
-   //TLorentzVector* jet = (TLorentzVector*)jets->Last();
-    AliAODJet* jet = (AliAODJet*)jets->Last();
-    jet->AddTrack(track);
-    tracks.Remove( track );
-
-    TVector2 jetVect2(jet->Px(),jet->Py()); // OB
-
-    // 3- Go to the next highest pT particle not already included...
-    AliAODTrack* track1;
-    while ( (track1  = (AliAODTrack*)(itrack.Next())) ) {
-      //Double_t dphi = TVector2::Phi_mpi_pi(jet->Phi()-track1->Phi()); // OB remove
-      TVector2 trackVect2(track1->Px(),track1->Py());
-      Double_t dphi = trackVect2.DeltaPhi(jetVect2);
-
-      Double_t r = TMath::Sqrt( (jet->Eta()-track1->Eta())*(jet->Eta()-track1->Eta()) +
-                               dphi*dphi );
-
-      if( r < fConeRadius ) {
-        Double_t fPt   = jet->E()+track1->Pt();  // Scalar sum of Pt
-        // recalculating the centroid
-        Double_t eta = jet->Eta()*jet->E()/fPt + track1->Eta()*track1->Pt()/fPt;
-
-        //Double_t phi = jet->Phi()*jet->E()/fPt + track1->Phi()*track1->Pt()/fPt; // OB - remove
-       // OB - recalc phi via weighted 2vectors
-       jetVect2.SetMagPhi(jet->E()/fPt,jetVect2.Phi());
-       trackVect2.SetMagPhi(track1->Pt()/fPt,trackVect2.Phi());
-       
-       TVector2 sumVect2 = jetVect2+trackVect2;
-       Double_t phi = sumVect2.Phi();
-
-        //jet->SetPtEtaPhiE( 1., eta, phi, fPt );
-       ((TLorentzVector*) jet->MomentumVector())->SetPtEtaPhiE(1,eta,phi,fPt);
-
-       jet->AddTrack(track1);
-        tracks.Remove(track1);
-      }
+  // get AOD event from input/ouput
+  TObject* handler = AliAnalysisManager::GetAnalysisManager()->GetInputEventHandler();
+  if( handler && handler->InheritsFrom("AliAODInputHandler") ) {
+    fAOD  =  ((AliAODInputHandler*)handler)->GetEvent();
+    if(fUseAODInputJets) fAODJets = fAOD;
+    if (fDebug > 1)  Printf("%s:%d AOD event from input", (char*)__FILE__,__LINE__);
+  }
+  else {
+    handler = AliAnalysisManager::GetAnalysisManager()->GetOutputEventHandler();
+    if( handler && handler->InheritsFrom("AliAODHandler") ) {
+      fAOD = ((AliAODHandler*)handler)->GetAOD();
+      fAODJets = fAOD;
+      if (fDebug > 1)  Printf("%s:%d AOD event from output", (char*)__FILE__,__LINE__);
     }
-    
-    tracks.Compress();
-
-    nTracks = tracks.GetEntries();
-    //   4- Continue until all particles are in a jet.
-    itrack.Reset();
-  } // end while nTracks
-  
-  // Convert to AODjets....
-  Int_t njets = jets->GetEntriesFast();
-
-  TClonesArray* aodjets = new TClonesArray("AliAODJet",njets);
-  aodjets->SetOwner(kTRUE);
-  
-  Int_t count = 0;
-
-  for(Int_t ijet=0; ijet<njets; ++ijet) {
-    //TLorentzVector* jet = (TLorentzVector*)jets->At(ijet);
-    AliAODJet* jet = (AliAODJet*)jets->At(ijet);
-
-    if (jet->E() < fJetPtCut) continue;
-    Float_t px, py,pz,en; // convert to 4-vector
-    px = jet->E() * TMath::Cos(jet->Phi());  // Pt * cos(phi)
-    py = jet->E() * TMath::Sin(jet->Phi());  // Pt * sin(phi)
-    pz = jet->E() / TMath::Tan(2.0 * TMath::ATan(TMath::Exp(-jet->Eta())));
-    en = TMath::Sqrt(px * px + py * py + pz * pz);
-    jet->SetPxPyPzE(px,py,pz,en);
-
-    TClonesArray &tmpaodjets = *aodjets;
-    new(tmpaodjets[count++]) AliAODJet(*jet);
-    //aodjets->AddLast( new AliAODJet(*jet));    
-    //aodjets->AddLast( new AliAODJet(px, py, pz, en) );
   }
-  // Order jets according to their pT .
-  //QSortTracks( *aodjets, 0, aodjets->GetEntriesFast() ); // not for TClonesArray
-  
-  if (fDebug>3) AliInfo(Form(" %d Charged jets found - after cuts %d \n",njets,count));
   
-  jets->Delete(); // OB: cleanup
-  delete jets;
-
-
-  return aodjets;
-}
-
-//____________________________________________________________________
-
-TClonesArray*  AliAnalysisTaskJetChem::GetPythiaJets(){
-
-  // return Pythia jets in jet acc
-
-  // note: present verion of AliMCEventHandler expects "AliAOD.root" (not "AliAODs.root"), otherwise galice.root will not be found in proper dir
-  AliMCEventHandler* mcHandler = dynamic_cast<AliMCEventHandler*> (AliAnalysisManager::GetAnalysisManager()->GetMCtruthEventHandler());
-  if (!mcHandler) {
-    Printf("ERROR: Could not retrieve MC event handler");
-    return 0;
+  if(!fAODJets && !fUseAODInputJets){ // case we have AOD in input & output and want jets from output
+    TObject* outHandler = AliAnalysisManager::GetAnalysisManager()->GetOutputEventHandler();
+    if( outHandler && outHandler->InheritsFrom("AliAODHandler") ) {
+      fAODJets = ((AliAODHandler*)outHandler)->GetAOD();
+      if (fDebug > 1)  Printf("%s:%d jets from output AOD", (char*)__FILE__,__LINE__);
+    }
   }
   
-  AliMCEvent* mcEvent = mcHandler->MCEvent();
-  if (!mcEvent) {
-    Printf("ERROR: Could not retrieve MC event");
-    return 0;
+  if(fNonStdFile.Length()!=0){
+    // case we have an AOD extension - fetch the jets from the extended output
+    
+    AliAODHandler *aodH = dynamic_cast<AliAODHandler*>(AliAnalysisManager::GetAnalysisManager()->GetOutputEventHandler());
+    fAODExtension = (aodH?aodH->GetExtension(fNonStdFile.Data()):0);    
+    if(!fAODExtension){
+      if(fDebug>1)Printf("AODExtension not found for %s",fNonStdFile.Data());
+    }
   }
   
-  
-  AliGenPythiaEventHeader*  pythiaGenHeader = AliAnalysisHelperJetTasks::GetPythiaEventHeader(mcEvent);
-  if(!pythiaGenHeader){
-    Printf("ERROR: Could not retrieve pythiaGenHeader");
-    return 0;
+  if(!fAOD){
+    Printf("%s:%d AODEvent not found", (char*)__FILE__,__LINE__);
+    return;
+  }
+  if(!fAODJets){
+    Printf("%s:%d AODEvent with jet branch not found", (char*)__FILE__,__LINE__);
+    return;
   }
   
-  //Get Jets from MC header
-  Int_t nPythiaGenJets = pythiaGenHeader->NTriggerJets();
-  
-  TClonesArray* aodjets = new TClonesArray("AliAODJet",nPythiaGenJets);
-  aodjets->SetOwner(kTRUE);
-
-  int count = 0;
-  for(int ip = 0; ip<nPythiaGenJets; ++ip){
-    Float_t p[4];
-    pythiaGenHeader->TriggerJet(ip,p);
-    TVector3 tempVect(p[0],p[1],p[2]);
-    if ( TMath::Abs(tempVect.Eta())>fJetEtaCut ) continue;
-
-    TClonesArray &tmpaodjets = *aodjets;
-    new(tmpaodjets[count++]) AliAODJet(p[0],p[1],p[2],p[3]);
-    //aodjets->AddLast(new AliAODJet(p[0],p[1],p[2],p[3]));
-   }
-  
-  return aodjets;
-}
-
 
-//____________________________________________________________________
-void  AliAnalysisTaskJetChem::QSortTracks(TObjArray &a, Int_t first, Int_t last)
-{
-  // Sort array of TObjArray of tracks by Pt using a quicksort algorithm.
-  
-  static TObject *tmp;
-  static int i;           // "static" to save stack space
-  int j;
+  // event selection  *****************************************
   
-  while (last - first > 1) {
-    i = first;
-    j = last;
-    for (;;) {
-      while (++i < last && ((AliVParticle*)a[i])->Pt() > ((AliVParticle*)a[first])->Pt() )
-        ;
-      while (--j > first && ((AliVParticle*)a[j])->Pt() < ((AliVParticle*)a[first])->Pt() )
-        ;
-      if (i >= j)
-        break;
-      
-      tmp  = a[i];
-      a[i] = a[j];
-      a[j] = tmp;
+  Double_t centPercent = -1;
+  if(fEventClass>0){
+    Int_t cl = 0;
+    if(handler->InheritsFrom("AliAODInputHandler")){ 
+      // since it is not supported by the helper task define own classes
+      centPercent = fAOD->GetHeader()->GetCentrality();
+      cl = 1;
+      if(centPercent>10) cl = 2;
+      if(centPercent>30) cl = 3;
+      if(centPercent>50) cl = 4;
     }
-    if (j == first) {
-      ++first;
-      continue;
+    else {
+      cl = AliAnalysisHelperJetTasks::EventClass();
+      centPercent = fESD->GetCentrality()->GetCentralityPercentile("V0M"); // OB added
     }
-    tmp = a[first];
-    a[first] = a[j];
-    a[j] = tmp;
-    if (j - first < last - (j + 1)) {
-      QSortTracks(a, first, j);
-      first = j + 1;   // QSortTracks(j + 1, last);
-    } else {
-      QSortTracks(a, j + 1, last);
-      last = j;        // QSortTracks(first, j);
+    
+    if(cl!=fEventClass){
+      // event not in selected event class, reject event
+      if (fDebug > 1) Printf("%s:%d event not in selected event class: event REJECTED ...",(char*)__FILE__,__LINE__);
+      fh1EvtSelection->Fill(2.);
+      PostData(1, fCommonHistList);
+      return;
     }
   }
-}
-   
-//------------------------------------------------------------------
-
-TH1F* AliAnalysisTaskJetChem::CreatePIDhisto(const char* name){
   
-  // create histogram
+  // *** vertex cut ***
+  AliAODVertex* primVtx = fAOD->GetPrimaryVertex();
+  Int_t nTracksPrim = primVtx->GetNContributors();
+  fh1VertexNContributors->Fill(nTracksPrim);
   
-  TH1F* result = new TH1F(name,"",60,0,60);
-  result->SetOption("E");
-
-  // bin equal Geant ID
-
-  result->GetXaxis()->SetBinLabel(result->GetXaxis()->FindBin(1),"photon");
-  result->GetXaxis()->SetBinLabel(result->GetXaxis()->FindBin(2),"e+");
-  result->GetXaxis()->SetBinLabel(result->GetXaxis()->FindBin(3),"e-");
-  result->GetXaxis()->SetBinLabel(result->GetXaxis()->FindBin(4),"e-neutrino"); 
-  result->GetXaxis()->SetBinLabel(result->GetXaxis()->FindBin(5),"mu+");
-  result->GetXaxis()->SetBinLabel(result->GetXaxis()->FindBin(6),"mu-");
-  result->GetXaxis()->SetBinLabel(result->GetXaxis()->FindBin(7),"pi0");
-  result->GetXaxis()->SetBinLabel(result->GetXaxis()->FindBin(8),"pi+");
-  result->GetXaxis()->SetBinLabel(result->GetXaxis()->FindBin(9),"pi-");
-  result->GetXaxis()->SetBinLabel(result->GetXaxis()->FindBin(10),"K long");
-  result->GetXaxis()->SetBinLabel(result->GetXaxis()->FindBin(11),"K+");
-  result->GetXaxis()->SetBinLabel(result->GetXaxis()->FindBin(12),"K-");
-  result->GetXaxis()->SetBinLabel(result->GetXaxis()->FindBin(13),"n");
-  result->GetXaxis()->SetBinLabel(result->GetXaxis()->FindBin(14),"p");
-  result->GetXaxis()->SetBinLabel(result->GetXaxis()->FindBin(15),"anti-proton");
-  result->GetXaxis()->SetBinLabel(result->GetXaxis()->FindBin(16),"K short");
-  result->GetXaxis()->SetBinLabel(result->GetXaxis()->FindBin(17),"eta");
-  result->GetXaxis()->SetBinLabel(result->GetXaxis()->FindBin(18),"Lambda");
-  result->GetXaxis()->SetBinLabel(result->GetXaxis()->FindBin(19),"Sigma+");
-  result->GetXaxis()->SetBinLabel(result->GetXaxis()->FindBin(20),"Sigma0");
-  result->GetXaxis()->SetBinLabel(result->GetXaxis()->FindBin(21),"Sigma-");
-  result->GetXaxis()->SetBinLabel(result->GetXaxis()->FindBin(22),"Xi0");
-  result->GetXaxis()->SetBinLabel(result->GetXaxis()->FindBin(23),"Xi-");
-  result->GetXaxis()->SetBinLabel(result->GetXaxis()->FindBin(24),"Omega-");
-  result->GetXaxis()->SetBinLabel(result->GetXaxis()->FindBin(25),"anti-neutron");
-  result->GetXaxis()->SetBinLabel(result->GetXaxis()->FindBin(26),"anti-Lambda");
-  result->GetXaxis()->SetBinLabel(result->GetXaxis()->FindBin(27),"anti-Sigma-");
-  result->GetXaxis()->SetBinLabel(result->GetXaxis()->FindBin(28),"anti-Sigma0");
-  result->GetXaxis()->SetBinLabel(result->GetXaxis()->FindBin(29),"anti-Sigma+"); 
-  result->GetXaxis()->SetBinLabel(result->GetXaxis()->FindBin(30),"anti-Xi0");
-  result->GetXaxis()->SetBinLabel(result->GetXaxis()->FindBin(31),"anti-Xi+");
-  result->GetXaxis()->SetBinLabel(result->GetXaxis()->FindBin(32),"anti-Omega+");
-  result->GetXaxis()->SetBinLabel(result->GetXaxis()->FindBin(33),"tau+");
-  result->GetXaxis()->SetBinLabel(result->GetXaxis()->FindBin(34),"tau-");
-  result->GetXaxis()->SetBinLabel(result->GetXaxis()->FindBin(35),"D+");
-  result->GetXaxis()->SetBinLabel(result->GetXaxis()->FindBin(36),"D-");
-  result->GetXaxis()->SetBinLabel(result->GetXaxis()->FindBin(37),"D0");
-  result->GetXaxis()->SetBinLabel(result->GetXaxis()->FindBin(38),"anti-D0");
-  result->GetXaxis()->SetBinLabel(result->GetXaxis()->FindBin(39),"Ds+");
-  result->GetXaxis()->SetBinLabel(result->GetXaxis()->FindBin(40),"anti Ds-");
-  result->GetXaxis()->SetBinLabel(result->GetXaxis()->FindBin(41),"Lamba_c+");
-  result->GetXaxis()->SetBinLabel(result->GetXaxis()->FindBin(42),"W+");
-  result->GetXaxis()->SetBinLabel(result->GetXaxis()->FindBin(43),"W-");
-  result->GetXaxis()->SetBinLabel(result->GetXaxis()->FindBin(44),"Z0");
-  result->GetXaxis()->SetBinLabel(result->GetXaxis()->FindBin(45),"d");
-  result->GetXaxis()->SetBinLabel(result->GetXaxis()->FindBin(46),"t");
-  result->GetXaxis()->SetBinLabel(result->GetXaxis()->FindBin(47),"alpha");
-  result->GetXaxis()->SetBinLabel(result->GetXaxis()->FindBin(48),"G_nu");
-  result->GetXaxis()->SetBinLabel(result->GetXaxis()->FindBin(kPDGpm311Bin),"K0/#bar{K0}");
-  result->GetXaxis()->SetBinLabel(result->GetXaxis()->FindBin(kPDG333Bin),"phi");
-  result->GetXaxis()->SetBinLabel(result->GetXaxis()->FindBin(kPDGpm313Bin),"K*(892)0");
-  result->GetXaxis()->SetBinLabel(result->GetXaxis()->FindBin(kPDGp323Bin),"K*(892)+");
-  result->GetXaxis()->SetBinLabel(result->GetXaxis()->FindBin(kPDGm323Bin),"K*(892)-");
-  result->GetXaxis()->SetBinLabel(result->GetXaxis()->FindBin(kPDGNeutrinoBin),"nu");
-  result->GetXaxis()->SetBinLabel(result->GetXaxis()->FindBin(kPDGCharmedBaryonBin),"charmed baryon");
-  result->GetXaxis()->SetBinLabel(result->GetXaxis()->FindBin(kPDGQuarkBin),"q/#bar{q}");
-  result->GetXaxis()->SetBinLabel(result->GetXaxis()->FindBin(kPDGDiQuarkBin),"q #bar{q}");
-  result->GetXaxis()->LabelsOption("v"); // "u" ?
+  if (fDebug > 1) Printf("%s:%d primary vertex selection: %d", (char*)__FILE__,__LINE__,nTracksPrim);
+  if(!nTracksPrim){
+    if (fDebug > 1) Printf("%s:%d primary vertex selection: event REJECTED...",(char*)__FILE__,__LINE__); 
+    fh1EvtSelection->Fill(3.);
+    PostData(1, fCommonHistList);
+    return;
+  }
   
-  return result;
-}
-
-//------------------------------------------------------------------
-
-TH1F*  AliAnalysisTaskJetChem::CreatePythiaIDhisto(const char* name){
+  fh1VertexZ->Fill(primVtx->GetZ());
   
-  // create histogram
+  if(TMath::Abs(primVtx->GetZ())>fMaxVertexZ){
+    if (fDebug > 1) Printf("%s:%d primary vertex z = %f: event REJECTED...",(char*)__FILE__,__LINE__,primVtx->GetZ()); 
+    fh1EvtSelection->Fill(4.);
+    PostData(1, fCommonHistList);
+    return; 
+  }
   
-  TH1F* result = new TH1F(name,"",22,0,22);
-  result->SetOption("E");
-
-  result->GetXaxis()->SetBinLabel(kPythiaPIDP11Bin,"qq #rightarrow qq");             // ISUB = 11 
-  result->GetXaxis()->SetBinLabel(kPythiaPIDP12Bin,"q#bar{q} #rightarrow q#bar{q}"); // ISUB = 12
-  result->GetXaxis()->SetBinLabel(kPythiaPIDP13Bin,"q#bar{q} #rightarrow gg");       // ISUB = 13
-  result->GetXaxis()->SetBinLabel(kPythiaPIDP28Bin,"qg #rightarrow qg ");            // ISUB = 28 
-  result->GetXaxis()->SetBinLabel(kPythiaPIDP53Bin,"gg #rightarrow q#bar{q}");       // ISUB = 53
-  result->GetXaxis()->SetBinLabel(kPythiaPIDP68Bin,"gg #rightarrow gg");             // ISUB = 68
-  result->GetXaxis()->SetBinLabel(kPythiaPIDP92Bin,"SD");                            // ISUB = 92
-  result->GetXaxis()->SetBinLabel(kPythiaPIDP93Bin,"SD");                            // ISUB = 93
-  result->GetXaxis()->SetBinLabel(kPythiaPIDP94Bin,"DD");                            // ISUB = 94
-  result->GetXaxis()->SetBinLabel(kPythiaPIDP95Bin,"low pt (MPI)");                  // ISUB = 95
-  result->GetXaxis()->SetBinLabel(kPythiaPIDPOtherBin,"other");                      // ISUB = XX
-
-  result->GetXaxis()->LabelsOption("v"); // "u" ?
+  TString primVtxName(primVtx->GetName());
   
-  return result;
-}
-
-
-//------------------------------------------------------------------
-
-void AliAnalysisTaskJetChem::FillPythiaIDhisto(TH1F* h, const Int_t PID){
-
-  // fill Pyhtia PID histogram
-
-  Int_t bin = -1; 
-  if(PID == 11)      bin = kPythiaPIDP11Bin;
-  else if(PID == 12) bin = kPythiaPIDP12Bin;
-  else if(PID == 13) bin = kPythiaPIDP13Bin;
-  else if(PID == 28) bin = kPythiaPIDP28Bin;
-  else if(PID == 53) bin = kPythiaPIDP53Bin;
-  else if(PID == 68) bin = kPythiaPIDP68Bin;
-  else if(PID == 92) bin = kPythiaPIDP92Bin;
-  else if(PID == 93) bin = kPythiaPIDP93Bin;
-  else if(PID == 94) bin = kPythiaPIDP94Bin;
-  else if(PID == 95) bin = kPythiaPIDP95Bin;
-  else{ 
-    if(PID != -1) AliInfo(Form("unknown PID %d",PID));
+  if(primVtxName.CompareTo("TPCVertex",TString::kIgnoreCase) == 1){
+    if (fDebug > 1) Printf("%s:%d primary vertex selection: TPC vertex, event REJECTED...",(char*)__FILE__,__LINE__);
+    fh1EvtSelection->Fill(5.);
+    PostData(1, fCommonHistList);
     return;
   }
   
-  h->Fill(h->GetBinCenter(bin),1);
-}
-
-
-//____________________________________________________________________
+  if (fDebug > 1) Printf("%s:%d event ACCEPTED ...",(char*)__FILE__,__LINE__); 
+  fh1EvtSelection->Fill(0.);
+  fh1EvtCent->Fill(centPercent);
+    
 
+  //___ get MC information __________________________________________________________________
 
-void  AliAnalysisTaskJetChem::FillPIDhisto(TH1F* hist, Int_t pdg, Float_t weight){
+  Double_t ptHard = 0.;
+  Double_t nTrials = 1; // trials for MC trigger weight for real data
   
-  // convert pdg code to Geant ID and fill corresponding bin 
+  if(fMCEvent){
+     AliGenEventHeader* genHeader = fMCEvent->GenEventHeader();
+     AliGenPythiaEventHeader*  pythiaGenHeader = dynamic_cast<AliGenPythiaEventHeader*>(genHeader);
+     AliGenHijingEventHeader*  hijingGenHeader = 0x0;
 
-  Int_t fGID = fpdgdb->ConvertPdgToGeant3(pdg);
+     if(pythiaGenHeader){
+        if(fDebug>3) Printf("%s:%d pythiaGenHeader found", (char*)__FILE__,__LINE__);
+        nTrials = pythiaGenHeader->Trials();
+        ptHard  = pythiaGenHeader->GetPtHard();
 
-  //cout<<" pdg "<<pdg<<" fGID "<<fGID<<endl;
+        fh1PtHard->Fill(ptHard);
+        fh1PtHardTrials->Fill(ptHard,nTrials);
 
-  if(TMath::Abs(pdg) == 311) fGID = kPDGpm311Bin; 
-  if(pdg == 333)             fGID = kPDG333Bin; 
-  if(TMath::Abs(pdg) == 313) fGID = kPDGpm313Bin; 
-  if(pdg == 323)             fGID = kPDGp323Bin; 
-  if(pdg == -323)            fGID = kPDGm323Bin; 
 
-  if(TMath::Abs(pdg)==12 || TMath::Abs(pdg)==14 || TMath::Abs(pdg)==16) fGID = kPDGNeutrinoBin;
+     } else { // no pythia, hijing?
 
-  if(TMath::Abs(pdg)==4122) fGID = kPDGCharmedBaryonBin;
+        if(fDebug>3) Printf("%s:%d no pythiaGenHeader found", (char*)__FILE__,__LINE__);
 
-  if(1<=TMath::Abs(pdg) && TMath::Abs(pdg)<=6) fGID = kPDGQuarkBin; 
+         hijingGenHeader = dynamic_cast<AliGenHijingEventHeader*>(genHeader);
+         if(!hijingGenHeader){
+            Printf("%s:%d no pythiaGenHeader or hjingGenHeader found", (char*)__FILE__,__LINE__);
+         } else {
+           if(fDebug>3) Printf("%s:%d hijingGenHeader found", (char*)__FILE__,__LINE__);
+        }
+     }
 
-  if(TMath::Abs(pdg)==1103 || TMath::Abs(pdg)==2101 || TMath::Abs(pdg)==2103 || TMath::Abs(pdg)==2203 || 
-     TMath::Abs(pdg)==3101 || TMath::Abs(pdg)==3103 || TMath::Abs(pdg)==3201 || TMath::Abs(pdg)==3203 || 
-     TMath::Abs(pdg)==3303 || TMath::Abs(pdg)==4101 || TMath::Abs(pdg)==4103 || TMath::Abs(pdg)==4201 || 
-     TMath::Abs(pdg)==4203 || TMath::Abs(pdg)==4301 || TMath::Abs(pdg)==4303 || TMath::Abs(pdg)==4403 || 
-     TMath::Abs(pdg)==5101 || TMath::Abs(pdg)==5103 || TMath::Abs(pdg)==5201 || TMath::Abs(pdg)==5203 || 
-     TMath::Abs(pdg)==5301 || TMath::Abs(pdg)==5303 || TMath::Abs(pdg)==5401 || TMath::Abs(pdg)==5403 || 
-     TMath::Abs(pdg)==5503)  fGID = kPDGDiQuarkBin;
-    
-
-  hist->Fill(fGID,weight);
-
-  if(fGID>hist->GetBinLowEdge(hist->GetNbinsX()+1) || 
-     fGID<hist->GetBinLowEdge(1)){
-    
-    AliError(Form("fGID %d for pdg %d exceeding histo limits ",fGID,pdg));
-  }
-
-  if(fGID == 0){
-    AliInfo(Form("fGID 0 for pdg %d ",pdg));
+     fh1Trials->Fill("#sum{ntrials}",fAvgTrials);
   }
 
-}
 
-//____________________________________________________________________
-void  AliAnalysisTaskJetChem::CreateHistos(){
+  //____ fetch jets ______________________________________________________________
 
-  // make histos
-  
-  fListOfHistos = new TList();
-  fListOfHistos->SetOwner(kTRUE);  
-
-  fhPrimVertexNCont = new TH1F("hPrimVertexNCont","",52,-2,50);
-  fhPrimVertexNCont->SetXTitle("");
-  fhPrimVertexNCont->Sumw2();
-  fListOfHistos->Add( fhPrimVertexNCont ); 
-
-  fhPrimVertexRho = new TH1F("hPrimVertexRho","",100,0,1);
-  fhPrimVertexRho->SetXTitle("");
-  fhPrimVertexRho->Sumw2();
-  fListOfHistos->Add( fhPrimVertexRho ); 
-
-  fhPrimVertexZ = new TH1F("hPrimVertexZ","",40,-20,20);
-  fhPrimVertexZ->SetXTitle("");
-  fhPrimVertexZ->Sumw2();
-  fListOfHistos->Add( fhPrimVertexZ ); 
-
-  fhNJets = new TH1F("hNJets","",10, 0, 10);
-  fhNJets->SetXTitle("# of jets");
-  fhNJets->Sumw2();
-  fListOfHistos->Add( fhNJets ); 
-
-  fhNJetsMC = new TH1F("hNJetsMC","",10,0,10);
-  fhNJetsMC->SetXTitle("# of jets");
-  fhNJetsMC->Sumw2();
-  fListOfHistos->Add( fhNJetsMC ); 
-
-  fhLeadingEta = new TH1F("hLeadingEta","Leading Jet eta",12,-0.6,0.6);
-  fhLeadingEta->SetXTitle("eta");
-  fhLeadingEta->SetYTitle("dN/deta");
-  fhLeadingEta->Sumw2();
-  fListOfHistos->Add(fhLeadingEta); 
-
-  fhLeadingNTracksVsEta = new TH2F("hLeadingNTracksVsEta","",20,-1.0,1.0,20,0,20);
-  fhLeadingNTracksVsEta->SetXTitle("eta");
-  fhLeadingNTracksVsEta->SetYTitle("# of tracks");
-  fhLeadingNTracksVsEta->Sumw2();
-  fListOfHistos->Add(fhLeadingNTracksVsEta); 
-
-  fhLeadingPtVsEta = new TH2F("hLeadingPtVsEta","",20,-1.0,1.0,50,0,50);
-  fhLeadingPtVsEta->SetXTitle("eta");
-  fhLeadingPtVsEta->SetYTitle("# of tracks");
-  fhLeadingPtVsEta->Sumw2();
-  fListOfHistos->Add(fhLeadingPtVsEta); 
-
-
-  fhLeadingPhi = new TH1F("hLeadingPhi","Leading Jet phi",63,0,6.3);
-  fhLeadingPhi->SetXTitle("phi");
-  fhLeadingPhi->SetYTitle("dN/dphi");
-  fhLeadingPhi->Sumw2();
-  fListOfHistos->Add(fhLeadingPhi);
-
-
-  fhLeadingPt  = new TH1F("hLeadingPt","leading Jet p_{T}",50,0,50);
-  fhLeadingPt->SetXTitle("p_{T} (GeV/c)");
-  fhLeadingPt->SetYTitle("dN/dp_{T} (1/GeV)");
-  fhLeadingPt->Sumw2();
-  fListOfHistos->Add( fhLeadingPt ); 
-
-  fhLeadingPtDiffr  = new TH1F("hLeadingPtDiffr","leading Jet p_{T}",50,0,50);
-  fhLeadingPtDiffr->SetXTitle("P_{T} (GeV/c)");
-  fhLeadingPtDiffr->SetYTitle("dN/dp_{T} (1/GeV)");
-  fhLeadingPtDiffr->Sumw2();
-  fListOfHistos->Add( fhLeadingPtDiffr ); 
-
-  fhLeadingEtaMC = new TH1F("hLeadingEtaMC","Leading Jet eta",12,-0.6,0.6);
-  fhLeadingEtaMC->SetXTitle("eta");
-  fhLeadingEtaMC->SetYTitle("dN/deta");
-  fhLeadingEtaMC->Sumw2();
-  fListOfHistos->Add(fhLeadingEtaMC); 
-
-  fhLeadingPhiMC = new TH1F("hLeadingPhiMC","Leading Jet phi",63,0,6.3);
-  fhLeadingPhiMC->SetXTitle("phi");
-  fhLeadingPhiMC->SetYTitle("dN/dphi");
-  fhLeadingPhiMC->Sumw2();
-  fListOfHistos->Add(fhLeadingPhiMC);
-
-  fhLeadingPtMC  = new TH1F("hLeadingPtMC","Leading Jet p_{T}",50,0,50);
-  fhLeadingPtMC->SetXTitle("p_{T} (GeV/c)");
-  fhLeadingPtMC->SetYTitle("dN/dp_{T} (1/GeV)");
-  fhLeadingPtMC->Sumw2();
-  fListOfHistos->Add( fhLeadingPtMC ); 
-
-  fhLeadingPtMCDiffr  = new TH1F("hLeadingPtMCDiffr","Leading Jet p_{T}",50,0,50);
-  fhLeadingPtMCDiffr->SetXTitle("p_{T} (GeV/c)");
-  fhLeadingPtMCDiffr->SetYTitle("dN/dp_{T} (1/GeV)");
-  fhLeadingPtMCDiffr->Sumw2();
-  fListOfHistos->Add( fhLeadingPtMCDiffr ); 
-
-  fhPhiEtaTracksNoCut = new TH2F("hPhiEtaTracksNoCut","phi vs eta tracks",20,-1.0,1.0,63,0,6.3);
-  fhPhiEtaTracksNoCut->SetXTitle("eta");
-  fhPhiEtaTracksNoCut->SetYTitle("phi");
-  fhPhiEtaTracksNoCut->Sumw2();
-  fListOfHistos->Add(fhPhiEtaTracksNoCut);
-
-  fhPtTracksNoCut = new TH1F("hPtTracksNoCut","p_{T} tracks",150,0,150);
-  fhPtTracksNoCut->SetXTitle("p_{T} (GeV)");
-  fhPtTracksNoCut->SetYTitle("dN/dp_{T} (1/GeV)");
-  fhPtTracksNoCut->Sumw2();
-  fListOfHistos->Add(fhPtTracksNoCut);
-
-  fhPhiEtaTracks = new TH2F("hPhiEtaTracks","phi vs eta tracks",20,-1.0,1.0,63,0,6.3);
-  fhPhiEtaTracks->SetXTitle("eta");
-  fhPhiEtaTracks->SetYTitle("phi");
-  fhPhiEtaTracks->Sumw2();
-  fListOfHistos->Add(fhPhiEtaTracks);
-
-  fhPtTracks = new TH1F("hPtTracks","p_{T} tracks",150,0,150);
-  fhPtTracks->SetXTitle("P{T} (GeV)");
-  fhPtTracks->SetYTitle("dN/dp_{T} (1/GeV)");
-  fhPtTracks->Sumw2();
-  fListOfHistos->Add(fhPtTracks);
-
-  fhTrackMult = new TH1F("hTrackMult","",150,0,150);
-  fhTrackMult->SetXTitle("n Tracks");
-  fhTrackMult->SetYTitle("counts");
-  fhTrackMult->Sumw2();
-  fListOfHistos->Add(fhTrackMult);
-
-  fhEtaMCTracks = new TH1F("hEtaMCTracks","eta tracks",30,-1.5,1.5);
-  fhEtaMCTracks->SetXTitle("eta");
-  fhEtaMCTracks->SetYTitle("dN/deta");
-  fhEtaMCTracks->Sumw2();
-  fListOfHistos->Add(fhEtaMCTracks);
-
-  fhPhiMCTracks = new TH1F("hPhiMCTracks","phi tracks",63,0,6.3);
-  fhPhiMCTracks->SetXTitle("phi");
-  fhPhiMCTracks->SetYTitle("dN/dphi");
-  fhPhiMCTracks->Sumw2();
-  fListOfHistos->Add(fhPhiMCTracks);
-
-  fhPtMCTracks = new TH1F("hPtMCTracks","p_{T} tracks",50,0,50);
-  fhPtMCTracks->SetXTitle("p_{T} (GeV)");
-  fhPtMCTracks->SetYTitle("dN/dp_{T} (1/GeV)");
-  fhPtMCTracks->Sumw2();
-  fListOfHistos->Add(fhPtMCTracks);
-
-  fhnTracksVsPtLeading  = new TH2F("hnTracksVsPtLeading","",50,0.,50.,20,-0.5,19.5);
-  fhnTracksVsPtLeading->SetXTitle("p_{T} (GeV/c)");
-  fhnTracksVsPtLeading->SetYTitle("n tracks");
-  fhnTracksVsPtLeading->Sumw2();
-  fListOfHistos->Add( fhnTracksVsPtLeading );
-
-  fhdNdEtaPhiDist  = new TH1F("hdNdEtaPhiDist","Charge particle density |#eta|< 1 vs #Delta#phi",  120, 0.,   2.*TMath::Pi());
-  fhdNdEtaPhiDist->SetXTitle("#Delta#phi");
-  fhdNdEtaPhiDist->SetYTitle("dN{ch}/d#etad#phi");
-  fhdNdEtaPhiDist->Sumw2();
-  fListOfHistos->Add( fhdNdEtaPhiDist );
-
-  fhRegionSumPtMaxVsEt = new TH1F("hRegionSumPtMaxVsEt","P{T}^{90, max} vs Leading Jet P{T}",50,0.,50.);
-  fhRegionSumPtMaxVsEt->SetXTitle("P{T} (GeV/c)");
-  fhRegionSumPtMaxVsEt->Sumw2();
-  fListOfHistos->Add( fhRegionSumPtMaxVsEt ); 
-
-  fhRegionMultMaxVsEt = new TH1F("hRegionMultMaxVsEt","N{ch}^{90, max} vs Leading Jet P{T}",50,0.,50.);
-  fhRegionMultMaxVsEt->SetXTitle("E (GeV hRegionAveSumPtVsEt/c)");
-  fhRegionMultMaxVsEt->Sumw2();
-  fListOfHistos->Add( fhRegionMultMaxVsEt );
-
-  fhRegionSumPtMinVsEt = new TH1F("hRegionSumPtMinVsEt","P{T}^{90, min} vs Leading Jet P{T}",50,0.,50.);
-  fhRegionSumPtMinVsEt->SetXTitle("P{T} (GeV/c)");
-  fhRegionSumPtMinVsEt->Sumw2();
-  fListOfHistos->Add( fhRegionSumPtMinVsEt );
-    
-  fhRegionMultMinVsEt = new TH1F("hRegionMultMinVsEt","N{ch}^{90, min} vs Leading Jet P{T}",50,0.,50.);
-  fhRegionMultMinVsEt->SetXTitle("E (GeV/c)");
-  fhRegionMultMinVsEt->Sumw2();
-  fListOfHistos->Add( fhRegionMultMinVsEt );
-  
-  // V0s 
-
-  fhNV0s = new TH1F("hNV0s","n V0s",50,0,50);
-  fhNV0s->SetXTitle("n V0s");
-  fhNV0s->Sumw2();
-  fListOfHistos->Add(fhNV0s);
-
-  fhV0onFly = new TH1F("hV0onFly","on-the-fly V0",5,0,5);
-  fhV0onFly->SetXTitle("is on-the-fly V0");
-  fhV0onFly->Sumw2();
-  fListOfHistos->Add(fhV0onFly);
-
-  fhV0DCADaughters = new TH1F("hV0DCADaughters","V0 DCA daughters",200,0,2.0);
-  fhV0DCADaughters->SetXTitle("V0 DCA daughters");
-  fhV0DCADaughters->Sumw2();
-  fListOfHistos->Add(fhV0DCADaughters); 
-
-  fhV0Radius = new TH1F("hV0Radius","V0 radius",2500,0,250);
-  fhV0Radius->SetXTitle("V0 radius");
-  fhV0Radius->Sumw2();
-  fListOfHistos->Add(fhV0Radius); 
-
-  fhV0DCAToVertex = new TH1F("hV0DCAToVertex","",100,0,10);
-  fhV0DCAToVertex->SetXTitle("V0 DCA (cm)");
-  fhV0DCAToVertex->Sumw2();
-  fListOfHistos->Add(fhV0DCAToVertex);
-
-  fhV0DCAToVertexK0 = new TH1F("hV0DCAToVertexK0","",100,0,10); 
-  fhV0DCAToVertexK0->SetXTitle("V0 DCA (cm)");
-  fhV0DCAToVertexK0->Sumw2();
-  fListOfHistos->Add(fhV0DCAToVertexK0);
-
-  fhV0InvMassK0 = new TH1F("hV0InvMassK0","",2000,0,2);
-  fhV0InvMassK0->SetXTitle("inv mass (GeV)");
-  fhV0InvMassK0->Sumw2();
-  fListOfHistos->Add(fhV0InvMassK0); 
-
-  fhV0PtVsInvMassK0 = new TH2F("hV0PtVsInvMassK0","",100,0.45,0.55,100,0.,50.);
-  fhV0PtVsInvMassK0->SetXTitle("inv mass (GeV)");
-  fhV0PtVsInvMassK0->SetYTitle("p_{T} (GeV)");
-  fhV0PtVsInvMassK0->Sumw2();
-  fListOfHistos->Add(fhV0PtVsInvMassK0); 
-
-  fhV0InvMassK0JetEvt = new TH1F("hV0InvMassK0JetEvt","",2000,0,2);
-  fhV0InvMassK0JetEvt->SetXTitle("inv mass (GeV)");
-  fhV0InvMassK0JetEvt->Sumw2();
-  fListOfHistos->Add(fhV0InvMassK0JetEvt); 
-
-  fhV0InvMassLambda = new TH1F("hV0InvMassLambda","",2000,0,2);
-  fhV0InvMassLambda->SetXTitle("inv mass (GeV)");
-  fhV0InvMassLambda->Sumw2();
-  fListOfHistos->Add(fhV0InvMassLambda);
-
-  fhV0InvMassAntiLambda = new TH1F("hV0InvMassAntiLambda","",2000,0,2);
-  fhV0InvMassAntiLambda->SetXTitle("inv mass (GeV)");
-  fhV0InvMassAntiLambda->Sumw2();
-  fListOfHistos->Add(fhV0InvMassAntiLambda);
-
-  fhV0InvMassLambdaJetEvt = new TH1F("hV0InvMassLambdaJetEvt","",2000,0,2);
-  fhV0InvMassLambdaJetEvt->SetXTitle("inv mass (GeV)");
-  fhV0InvMassLambdaJetEvt->Sumw2();
-  fListOfHistos->Add(fhV0InvMassLambdaJetEvt);
-
-  fhV0InvMassAntiLambdaJetEvt = new TH1F("hV0InvMassAntiLambdaJetEvt","",2000,0,2);
-  fhV0InvMassAntiLambdaJetEvt->SetXTitle("inv mass (GeV)");
-  fhV0InvMassAntiLambdaJetEvt->Sumw2();
-  fListOfHistos->Add(fhV0InvMassAntiLambdaJetEvt);
-
-  fhdROpanK0VsPt = new TH2F("hdROpanK0VsPt","V0 dR vs pt",100,0,10,100,0,10);
-  fhdROpanK0VsPt->SetXTitle("opening angle R (rad)");
-  fhdROpanK0VsPt->Sumw2();
-  fListOfHistos->Add(fhdROpanK0VsPt);
-
-  fhdPhiJetV0 = new TH1F("hdPhiJetV0","",640,-3.2,3.2);
-  fhdPhiJetV0->SetXTitle("#Delta #phi V0 - jet");
-  fhdPhiJetV0->Sumw2();
-  fListOfHistos->Add(fhdPhiJetV0); 
-
-  fhdPhiJetK0 = new TH1F("hdPhiJetK0","",640,-3.2,3.2);
-  fhdPhiJetK0->SetXTitle("#Delta #phi V0 - jet");
-  fhdPhiJetK0->Sumw2();
-  fListOfHistos->Add(fhdPhiJetK0); 
-
-  fhdRJetK0 = new TH1F("hdRJetK0","dN/dR K0-jet",500,0,5);
-  fhdRJetK0->SetXTitle("#Delta R K0 - jet");
-  fhdRJetK0->SetYTitle("1/N_{jet} dN/dR");
-  fhdRJetK0->Sumw2();
-  fListOfHistos->Add(fhdRJetK0); 
-
-  fhdNdptV0 = new TH1F("hdNdptV0","dN/dpt V0",100,0,10);
-  fhdNdptV0->SetXTitle("p_{T} (GeV/c)");
-  fhdNdptV0->SetYTitle("1/N_{jet} dN/dpt (1/GeV)");
-  fhdNdptV0->Sumw2();
-  fListOfHistos->Add(fhdNdptV0);
-
-  fhdNdptK0 = new TH1F("hdNdptK0","",100,0,10);
-  fhdNdptK0->SetXTitle("p_{T} (GeV/c)");
-  fhdNdptK0->SetYTitle("1/N_{jet} dN/dpt (1/GeV)");
-  fhdNdptK0->Sumw2();
-  fListOfHistos->Add(fhdNdptK0); 
-
-  fhPtVsEtaK0 = new TH2F("hPtVsEtaK0","",20,-1,1,100,0,10);
-  fhPtVsEtaK0->SetXTitle("#eta");  
-  fhPtVsEtaK0->SetYTitle("p_{T} (GeV/c)");
-  fhPtVsEtaK0->Sumw2();
-  fListOfHistos->Add(fhPtVsEtaK0); 
-
-  fhV0InvMassK0DCA = new TH1F("hV0InvMassK0DCA","",2000,0,2);
-  fhV0InvMassK0DCA->SetXTitle("inv mass (GeV)");
-  fhV0InvMassK0DCA->Sumw2();
-  fListOfHistos->Add(fhV0InvMassK0DCA); 
-
-  fhV0InvMassK0DCAdEdx = new TH1F("hV0InvMassK0DCAdEdx","",2000,0,2);
-  fhV0InvMassK0DCAdEdx->SetXTitle("inv mass (GeV)");
-  fhV0InvMassK0DCAdEdx->Sumw2();
-  fListOfHistos->Add(fhV0InvMassK0DCAdEdx); 
-
-  fhV0PtVsInvMassK0DCAdEdx = new TH2F("hV0PtVsInvMassK0DCAdEdx","",100,0.45,0.55,100,0.,50.);
-  fhV0PtVsInvMassK0DCAdEdx->SetXTitle("inv mass (GeV)");
-  fhV0PtVsInvMassK0DCAdEdx->SetYTitle("p_{T} (GeV)");
-  fhV0PtVsInvMassK0DCAdEdx->Sumw2();
-  fListOfHistos->Add(fhV0PtVsInvMassK0DCAdEdx); 
-
-  fhV0InvMassK0DCAPID = new TH1F("hV0InvMassK0DCAPID","",2000,0,2);
-  fhV0InvMassK0DCAPID->SetXTitle("inv mass (GeV)");
-  fhV0InvMassK0DCAPID->Sumw2();
-  fListOfHistos->Add(fhV0InvMassK0DCAPID); 
-
-  fhV0InvMassLambdaDCAdEdx = new TH1F("hV0InvMassLambdaDCAdEdx","",2000,0,2);
-  fhV0InvMassLambdaDCAdEdx->SetXTitle("inv mass (GeV)");
-  fhV0InvMassLambdaDCAdEdx->Sumw2();
-  fListOfHistos->Add(fhV0InvMassLambdaDCAdEdx); 
-
-  fhV0InvMassAntiLambdaDCAdEdx = new TH1F("hV0InvMassAntiLambdaDCAdEdx","",2000,0,2);
-  fhV0InvMassAntiLambdaDCAdEdx->SetXTitle("inv mass (GeV)");
-  fhV0InvMassAntiLambdaDCAdEdx->Sumw2();
-  fListOfHistos->Add(fhV0InvMassAntiLambdaDCAdEdx); 
-
-  fhdNdptK0DCA = new TH1F("hdNdptK0DCA","",100,0,10);
-  fhdNdptK0DCA->SetXTitle("p_{T} (GeV/c)");
-  fhdNdptK0DCA->SetYTitle("1/N_{jet} dN/dpt (1/GeV)");
-  fhdNdptK0DCA->Sumw2();
-  fListOfHistos->Add(fhdNdptK0DCA); 
-
-  fhdNdptK0DCAdEdx = new TH1F("hdNdptK0DCAdEdx","",100,0,10);
-  fhdNdptK0DCAdEdx->SetXTitle("p_{T} (GeV/c)");
-  fhdNdptK0DCAdEdx->SetYTitle("1/N_{jet} dN/dpt (1/GeV)");
-  fhdNdptK0DCAdEdx->Sumw2();
-  fListOfHistos->Add(fhdNdptK0DCAdEdx); 
-
-  fhV0InvMassK0Min = new TH1F("hV0InvMassK0Min","",2000,0,2);
-  fhV0InvMassK0Min->SetXTitle("inv mass (GeV)");
-  fhV0InvMassK0Min->Sumw2();
-  fListOfHistos->Add(fhV0InvMassK0Min); 
-
-  fhV0InvMassLambdaMin = new TH1F("hV0InvMassLambdaMin","",2000,0,2);
-  fhV0InvMassLambdaMin->SetXTitle("inv mass (GeV)");
-  fhV0InvMassLambdaMin->Sumw2();
-  fListOfHistos->Add(fhV0InvMassLambdaMin); 
-
-  fhV0InvMassAntiLambdaMin = new TH1F("hV0InvMassAntiLambdaMin","",2000,0,2);
-  fhV0InvMassAntiLambdaMin->SetXTitle("inv mass (GeV)");
-  fhV0InvMassAntiLambdaMin->Sumw2();
-  fListOfHistos->Add(fhV0InvMassAntiLambdaMin);
-
-  fhV0InvMassK0Max = new TH1F("hV0InvMassK0Max","",2000,0,2);
-  fhV0InvMassK0Max->SetXTitle("inv mass (GeV)");
-  fhV0InvMassK0Max->Sumw2();
-  fListOfHistos->Add(fhV0InvMassK0Max); 
-
-  fhV0InvMassLambdaMax = new TH1F("hV0InvMassLambdaMax","",2000,0,2);
-  fhV0InvMassLambdaMax->SetXTitle("inv mass (GeV)");
-  fhV0InvMassLambdaMax->Sumw2();
-  fListOfHistos->Add(fhV0InvMassLambdaMax); 
-
-  fhV0InvMassAntiLambdaMax = new TH1F("hV0InvMassAntiLambdaMax","",2000,0,2);
-  fhV0InvMassAntiLambdaMax->SetXTitle("inv mass (GeV)");
-  fhV0InvMassAntiLambdaMax->Sumw2();
-  fListOfHistos->Add(fhV0InvMassAntiLambdaMax); 
-
-  fhV0InvMassK0Jet = new TH1F("hV0InvMassK0Jet","",2000,0,2);
-  fhV0InvMassK0Jet->SetXTitle("inv mass (GeV)");
-  fhV0InvMassK0Jet->Sumw2();
-  fListOfHistos->Add(fhV0InvMassK0Jet); 
-
-  fhV0InvMassLambdaJet = new TH1F("hV0InvMassLambdaJet","",2000,0,2);
-  fhV0InvMassLambdaJet->SetXTitle("inv mass (GeV)");
-  fhV0InvMassLambdaJet->Sumw2();
-  fListOfHistos->Add(fhV0InvMassLambdaJet); 
-
-  fhV0InvMassAntiLambdaJet = new TH1F("hV0InvMassAntiLambdaJet","",2000,0,2);
-  fhV0InvMassAntiLambdaJet->SetXTitle("inv mass (GeV)");
-  fhV0InvMassAntiLambdaJet->Sumw2();
-  fListOfHistos->Add(fhV0InvMassAntiLambdaJet); 
-
-  fhV0InvMassK0Lambda = new TH1F("hV0InvMassK0Lambda","",2000,0,2);
-  fhV0InvMassK0Lambda->SetXTitle("inv mass (GeV)");
-  fhV0InvMassK0Lambda->Sumw2();
-  fListOfHistos->Add(fhV0InvMassK0Lambda); 
-
-  fhdNdptK0JetEvt = new TH1F("hdNdptK0JetEvt","",100,0,10);
-  fhdNdptK0JetEvt->SetXTitle("p_{T} (GeV/c)");
-  fhdNdptK0JetEvt->SetYTitle("1/N_{jet} dN/dpt (1/GeV)");
-  fhdNdptK0JetEvt->Sumw2();
-  fListOfHistos->Add(fhdNdptK0JetEvt); 
-  
-  fhdNdzK0 = new TH1F("hdNdzK0","",150,0,1.5);
-  fhdNdzK0->SetXTitle("z");
-  fhdNdzK0->SetYTitle("1/N_{jet} dN/dz");
-  fhdNdzK0->Sumw2();
-  fListOfHistos->Add(fhdNdzK0);
-
-  fhdNdzK05to10 = new TH1F("hdNdzK05to10","",150,0,1.5);
-  fhdNdzK05to10->SetXTitle("z");
-  fhdNdzK05to10->SetYTitle("1/N_{jet} dN/dz");
-  fhdNdzK05to10->Sumw2();
-  fListOfHistos->Add(fhdNdzK05to10);
-
-  fhdNdzK010to20 = new TH1F("hdNdzK010to20","",150,0,1.5);
-  fhdNdzK010to20->SetXTitle("z");
-  fhdNdzK010to20->SetYTitle("1/N_{jet} dN/dz");
-  fhdNdzK010to20->Sumw2();
-  fListOfHistos->Add(fhdNdzK010to20);
-
-  fhdNdzK020to30 = new TH1F("hdNdzK020to30","",150,0,1.5);
-  fhdNdzK020to30->SetXTitle("z");
-  fhdNdzK020to30->SetYTitle("1/N_{jet} dN/dz");
-  fhdNdzK020to30->Sumw2();
-  fListOfHistos->Add(fhdNdzK020to30);
-
-  fhdNdzK030to40 = new TH1F("hdNdzK030to40","",150,0,1.5);
-  fhdNdzK030to40->SetXTitle("z");
-  fhdNdzK030to40->SetYTitle("1/N_{jet} dN/dz");
-  fhdNdzK030to40->Sumw2();
-  fListOfHistos->Add(fhdNdzK030to40);
-
-  fhdNdzK040to60 = new TH1F("hdNdzK040to60","",150,0,1.5);
-  fhdNdzK040to60->SetXTitle("z");
-  fhdNdzK040to60->SetYTitle("1/N_{jet} dN/dz");
-  fhdNdzK040to60->Sumw2();
-  fListOfHistos->Add(fhdNdzK040to60);
-
-  fhdNdxiK0 = new TH1F("hdNdxiK0","",100,0,10);
-  fhdNdxiK0->SetXTitle("xi");
-  fhdNdxiK0->SetYTitle("1/N_{jet} dN/dxi");
-  fhdNdxiK0->Sumw2();
-  fListOfHistos->Add(fhdNdxiK0);
-  
-  fhdNdzLambda = new TH1F("hdNdzLambda","",150,0,1.5);
-  fhdNdzLambda->SetXTitle("z");
-  fhdNdzLambda->SetYTitle("1/N_{jet} dN/dz");
-  fhdNdzLambda->Sumw2();
-  fListOfHistos->Add(fhdNdzLambda);
-
-  fhdNdzAntiLambda = new TH1F("hdNdzAntiLambda","",150,0,1.5);
-  fhdNdzAntiLambda->SetXTitle("z");
-  fhdNdzAntiLambda->SetYTitle("1/N_{jet} dN/dz");
-  fhdNdzAntiLambda->Sumw2();
-  fListOfHistos->Add(fhdNdzAntiLambda);
-
-  fhdNdzK0Max = new TH1F("hdNdzK0Max","",150,0,1.5);
-  fhdNdzK0Max->SetXTitle("z");
-  fhdNdzK0Max->SetYTitle("1/N_{jet} dN/dz");
-  fhdNdzK0Max->Sumw2();
-  fListOfHistos->Add(fhdNdzK0Max);
-
-  fhdNdxiK0Max = new TH1F("hdNdxiK0Max","",100,0,10);
-  fhdNdxiK0Max->SetXTitle("xi");
-  fhdNdxiK0Max->SetYTitle("1/N_{jet} dN/dxi");
-  fhdNdxiK0Max->Sumw2();
-  fListOfHistos->Add(fhdNdxiK0Max);
-
-  fhdNdzLambdaMax = new TH1F("hdNdzLambdaMax","",150,0,1.5);
-  fhdNdzLambdaMax->SetXTitle("z");
-  fhdNdzLambdaMax->SetYTitle("1/N_{jet} dN/dz");
-  fhdNdzLambdaMax->Sumw2();
-  fListOfHistos->Add(fhdNdzLambdaMax); 
-
-  fhdNdxiLambdaMax = new TH1F("hdNdxiLambdaMax","",700,0,7);
-  fhdNdxiLambdaMax->SetXTitle("xi");
-  fhdNdxiLambdaMax->SetYTitle("1/N_{jet} dN/dxi");
-  fhdNdxiLambdaMax->Sumw2();
-  fListOfHistos->Add(fhdNdxiLambdaMax); 
-
-  fhdNdptK0Max = new TH1F("hdNdptK0Max","",100,0,10);
-  fhdNdptK0Max->SetXTitle("p_{T} (GeV/c)");
-  fhdNdptK0Max->SetYTitle("1/N_{jet} dN/dpt (1/GeV)");
-  fhdNdptK0Max->Sumw2();
-  fListOfHistos->Add(fhdNdptK0Max); 
-
-  fhdNdptLambdaMax = new TH1F("hdNdptLambdaMax","",100,0,10);
-  fhdNdptLambdaMax->SetXTitle("p_{T} (GeV/c)");
-  fhdNdptLambdaMax->SetYTitle("1/N_{jet} dN/dpt (1/GeV)");
-  fhdNdptLambdaMax->Sumw2();
-  fListOfHistos->Add(fhdNdptLambdaMax); 
-
-  fhdNdzK0Min = new TH1F("hdNdzK0Min","",150,0,1.5);
-  fhdNdzK0Min->SetXTitle("z");
-  fhdNdzK0Min->SetYTitle("1/N_{jet} dN/dz");
-  fhdNdzK0Min->Sumw2();
-  fListOfHistos->Add(fhdNdzK0Min); 
-
-  fhdNdxiK0Min = new TH1F("hdNdxiK0Min","",100,0,10);
-  fhdNdxiK0Min->SetXTitle("xi");
-  fhdNdxiK0Min->SetYTitle("1/N_{jet} dN/dxi");
-  fhdNdxiK0Min->Sumw2();
-  fListOfHistos->Add(fhdNdxiK0Min); 
-
-  fhdNdzLambdaMin = new TH1F("hdNdzLambdaMin","",150,0,1.5);
-  fhdNdzLambdaMin->SetXTitle("z");
-  fhdNdzLambdaMin->SetYTitle("1/N_{jet} dN/dz");
-  fhdNdzLambdaMin->Sumw2();
-  fListOfHistos->Add(fhdNdzLambdaMin); 
-
-  fhdNdxiLambdaMin = new TH1F("hdNdxiLambdaMin","",700,0,7);
-  fhdNdxiLambdaMin->SetXTitle("xi");
-  fhdNdxiLambdaMin->SetYTitle("1/N_{jet} dN/dxi");
-  fhdNdxiLambdaMin->Sumw2();
-  fListOfHistos->Add(fhdNdxiLambdaMin); 
-
-  fhdNdptK0Min = new TH1F("hdNdptK0Min","",100,0,10);
-  fhdNdptK0Min->SetXTitle("p_{T} (GeV/c)");
-  fhdNdptK0Min->SetYTitle("1/N_{jet} dN/dpt (1/GeV)");
-  fhdNdptK0Min->Sumw2();
-  fListOfHistos->Add(fhdNdptK0Min); 
-
-  fhdNdptLambdaMin = new TH1F("hdNdptLambdaMin","",100,0,10);
-  fhdNdptLambdaMin->SetXTitle("p_{T} (GeV/c)");
-  fhdNdptLambdaMin->SetYTitle("1/N_{jet} dN/dpt (1/GeV)");
-  fhdNdptLambdaMin->Sumw2();
-  fListOfHistos->Add(fhdNdptLambdaMin); 
-
-  fhdNdzK0Jet = new TH1F("hdNdzK0Jet","",150,0,1.5);
-  fhdNdzK0Jet->SetXTitle("z");
-  fhdNdzK0Jet->SetYTitle("1/N_{jet} dN/dz");
-  fhdNdzK0Jet->Sumw2();
-  fListOfHistos->Add(fhdNdzK0Jet);
-
-  fhdNdxiK0Jet = new TH1F("hdNdxiK0Jet","",100,0,10);
-  fhdNdxiK0Jet->SetXTitle("xi");
-  fhdNdxiK0Jet->SetYTitle("1/N_{jet} dN/dxi");
-  fhdNdxiK0Jet->Sumw2();
-  fListOfHistos->Add(fhdNdxiK0Jet); 
-
-  fhdNdzLambdaJet = new TH1F("hdNdzLambdaJet","",150,0,1.5);
-  fhdNdzLambdaJet->SetXTitle("z");
-  fhdNdzLambdaJet->SetYTitle("1/N_{jet} dN/dz");
-  fhdNdzLambdaJet->Sumw2();
-  fListOfHistos->Add(fhdNdzLambdaJet); 
-
-  fhdNdxiLambdaJet = new TH1F("hdNdxiLambdaJet","",700,0,7);
-  fhdNdxiLambdaJet->SetXTitle("xi");
-  fhdNdxiLambdaJet->SetYTitle("1/N_{jet} dN/dxi");
-  fhdNdxiLambdaJet->Sumw2();
-  fListOfHistos->Add(fhdNdxiLambdaJet); 
-
-  fhdNdptK0Jet = new TH1F("hdNdptK0Jet","",100,0,10);
-  fhdNdptK0Jet->SetXTitle("p_{T} (GeV/c)");
-  fhdNdptK0Jet->SetYTitle("1/N_{jet} dN/dpt (1/GeV)");
-  fhdNdptK0Jet->Sumw2();
-  fListOfHistos->Add(fhdNdptK0Jet); 
-
-  fhdNdptLambdaJet = new TH1F("hdNdptLambdaJet","",100,0,10);
-  fhdNdptLambdaJet->SetXTitle("p_{T} (GeV/c)");
-  fhdNdptLambdaJet->SetYTitle("1/N_{jet} dN/dpt (1/GeV)");
-  fhdNdptLambdaJet->Sumw2();
-  fListOfHistos->Add(fhdNdptLambdaJet); 
-
-  fhdEdxVsMomV0 = new TH2F("hdEdxVsMomV0","",200,-10,10,200,0,2000);
-  fhdEdxVsMomV0->SetXTitle("mom (GeV/c)");
-  fhdEdxVsMomV0->SetYTitle("dE/dx");
-  fhdEdxVsMomV0->Sumw2();
-  fListOfHistos->Add(fhdEdxVsMomV0); 
-
-  fhdEdxVsMomV0pidEdx = new TH2F("hdEdxVsMomV0pidEdx","",200,-10,10,200,0,2000);
-  fhdEdxVsMomV0pidEdx->SetXTitle("mom (GeV/c)");
-  fhdEdxVsMomV0pidEdx->SetYTitle("dE/dx");
-  fhdEdxVsMomV0pidEdx->Sumw2();
-  fListOfHistos->Add(fhdEdxVsMomV0pidEdx); 
-
-  fhdEdxVsMomV0piPID = new TH2F("hdEdxVsMomV0piPID","",200,-10,10,200,0,2000);
-  fhdEdxVsMomV0piPID->SetXTitle("mom (GeV/c)");
-  fhdEdxVsMomV0piPID->SetYTitle("dE/dx");
-  fhdEdxVsMomV0piPID->Sumw2();
-  fListOfHistos->Add(fhdEdxVsMomV0piPID); 
-
-  fhdPhiJetK0MC = new TH1F("hdPhiJetK0MC","",640,-3.2,3.2);
-  fhdPhiJetK0MC->SetXTitle("#Delta #phi K0 - jet");
-  fhdPhiJetK0MC->SetYTitle("1/N_{jet} dN/dphi");
-  fhdPhiJetK0MC->Sumw2();
-  fListOfHistos->Add(fhdPhiJetK0MC); 
-
-  fhdRJetK0MC   = new TH1F("hdRJetK0MC","dN/R K0-jet",500,0,5);
-  fhdRJetK0MC->SetXTitle("#Delta R K0 - jet");
-  fhdRJetK0MC->SetYTitle("1/N_{jet} dN/dR");
-  fhdRJetK0MC->Sumw2();
-  fListOfHistos->Add(fhdRJetK0MC); 
-
-  fhdRV0MC =  new TH1F("hdRV0MC","",500,0.,1.);
-  fhdRV0MC->SetXTitle("#Delta R");
-  fhdRV0MC->SetYTitle("");
-  fhdRV0MC->Sumw2();
-  fListOfHistos->Add(fhdRV0MC); 
-
-  fhdNdptchPiMCMax = new TH1F("hdNdptchPiMCMax","",100,0,10);
-  fhdNdptchPiMCMax->SetXTitle("p_{T} (GeV/c)");
-  fhdNdptchPiMCMax->SetYTitle("1/N_{jet} dN/dpt (1/GeV)");
-  fhdNdptchPiMCMax->Sumw2();
-  fListOfHistos->Add(fhdNdptchPiMCMax); 
-
-  fhdNdptK0MCMax = new TH1F("hdNdptK0MCMax","",100,0,10);
-  fhdNdptK0MCMax->SetXTitle("p_{T} (GeV/c)");
-  fhdNdptK0MCMax->SetYTitle("1/N_{jet} dN/dpt (1/GeV)");
-  fhdNdptK0MCMax->Sumw2();
-  fListOfHistos->Add(fhdNdptK0MCMax); 
-
-  fhdNdptchKMCMax = new TH1F("hdNdptchKMCMax","",100,0,10);
-  fhdNdptchKMCMax->SetXTitle("p_{T} (GeV/c)");
-  fhdNdptchKMCMax->SetYTitle("1/N_{jet} dN/dpt (1/GeV)");
-  fhdNdptchKMCMax->Sumw2();
-  fListOfHistos->Add(fhdNdptchKMCMax); 
-
-  fhdNdptpMCMax = new TH1F("hdNdptpMCMax","",100,0,10);
-  fhdNdptpMCMax->SetXTitle("p_{T} (GeV/c)");
-  fhdNdptpMCMax->SetYTitle("1/N_{jet} dN/dpt (1/GeV)");
-  fhdNdptpMCMax->Sumw2();
-  fListOfHistos->Add(fhdNdptpMCMax); 
-
-  fhdNdptpBarMCMax = new TH1F("hdNdptpBarMCMax","",100,0,10);
-  fhdNdptpBarMCMax->SetXTitle("p_{T} (GeV/c)");
-  fhdNdptpBarMCMax->SetYTitle("1/N_{jet} dN/dpt (1/GeV)");
-  fhdNdptpBarMCMax->Sumw2();
-  fListOfHistos->Add(fhdNdptpBarMCMax); 
-
-  fhdNdptLambdaMCMax = new TH1F("hdNdptLambdaMCMax","",100,0,10);
-  fhdNdptLambdaMCMax->SetXTitle("p_{T} (GeV/c)");
-  fhdNdptLambdaMCMax->SetYTitle("1/N_{jet} dN/dpt (1/GeV)");
-  fhdNdptLambdaMCMax->Sumw2();
-  fListOfHistos->Add(fhdNdptLambdaMCMax); 
-
-  fhdNdptLambdaBarMCMax = new TH1F("hdNdptLambdaBarMCMax","",100,0,10);
-  fhdNdptLambdaBarMCMax->SetXTitle("p_{T} (GeV/c)");
-  fhdNdptLambdaBarMCMax->SetYTitle("1/N_{jet} dN/dpt (1/GeV)");
-  fhdNdptLambdaBarMCMax->Sumw2();
-  fListOfHistos->Add(fhdNdptLambdaBarMCMax); 
-
-
-  fhdNdptchPiMCMin = new TH1F("hdNdptchPiMCMin","",100,0,10);
-  fhdNdptchPiMCMin->SetXTitle("p_{T} (GeV/c)");
-  fhdNdptchPiMCMin->SetYTitle("1/N_{jet} dN/dpt (1/GeV)");
-  fhdNdptchPiMCMin->Sumw2();
-  fListOfHistos->Add(fhdNdptchPiMCMin); 
-  fhdNdptK0MCMin = new TH1F("hdNdptK0MCMin","",100,0,10);
-  fhdNdptK0MCMin->SetXTitle("p_{T} (GeV/c)");
-  fhdNdptK0MCMin->SetYTitle("1/N_{jet} dN/dpt (1/GeV)");
-  fhdNdptK0MCMin->Sumw2();
-  fListOfHistos->Add(fhdNdptK0MCMin); 
-
-  fhdNdptchKMCMin = new TH1F("hdNdptchKMCMin","",100,0,10);
-  fhdNdptchKMCMin->SetXTitle("p_{T} (GeV/c)");
-  fhdNdptchKMCMin->SetYTitle("1/N_{jet} dN/dpt (1/GeV)");
-  fhdNdptchKMCMin->Sumw2();
-  fListOfHistos->Add(fhdNdptchKMCMin); 
-
-  fhdNdptpMCMin = new TH1F("hdNdptpMCMin","",100,0,10);
-  fhdNdptpMCMin->SetXTitle("p_{T} (GeV/c)");
-  fhdNdptpMCMin->SetYTitle("1/N_{jet} dN/dpt (1/GeV)");
-  fhdNdptpMCMin->Sumw2();
-  fListOfHistos->Add(fhdNdptpMCMin); 
-
-  fhdNdptpBarMCMin = new TH1F("hdNdptpBarMCMin","",100,0,10);
-  fhdNdptpBarMCMin->SetXTitle("p_{T} (GeV/c)");
-  fhdNdptpBarMCMin->SetYTitle("1/N_{jet} dN/dpt (1/GeV)");
-  fhdNdptpBarMCMin->Sumw2();
-  fListOfHistos->Add(fhdNdptpBarMCMin); 
-
-  fhdNdptLambdaMCMin = new TH1F("hdNdptLambdaMCMin","",100,0,10);
-  fhdNdptLambdaMCMin->SetXTitle("p_{T} (GeV/c)");
-  fhdNdptLambdaMCMin->SetYTitle("1/N_{jet} dN/dpt (1/GeV)");
-  fhdNdptLambdaMCMin->Sumw2();
-  fListOfHistos->Add(fhdNdptLambdaMCMin); 
-
-  fhdNdptLambdaBarMCMin = new TH1F("hdNdptLambdaBarMCMin","",100,0,10);
-  fhdNdptLambdaBarMCMin->SetXTitle("p_{T} (GeV/c)");
-  fhdNdptLambdaBarMCMin->SetYTitle("1/N_{jet} dN/dpt (1/GeV)");
-  fhdNdptLambdaBarMCMin->Sumw2();
-  fListOfHistos->Add(fhdNdptLambdaBarMCMin); 
-
-  fhdNdptOmegaMCMin = new TH1F("hdNdptOmegaMCMin","",100,0,10);;
-  fhdNdptOmegaMCMin->SetXTitle("p_{T} (GeV/c)");
-  fhdNdptOmegaMCMin->SetYTitle("1/N_{event} dN/dpt (1/GeV)");
-  fhdNdptOmegaMCMin->Sumw2();
-  fListOfHistos->Add(fhdNdptOmegaMCMin); 
-
-  fhdNdptOmegaBarMCMin = new TH1F("hdNdptOmegaBarMCMin","",100,0,10);;
-  fhdNdptOmegaBarMCMin->SetXTitle("p_{T} (GeV/c)");
-  fhdNdptOmegaBarMCMin->SetYTitle("1/N_{event} dN/dpt (1/GeV)");
-  fhdNdptOmegaBarMCMin->Sumw2();
-  fListOfHistos->Add(fhdNdptOmegaBarMCMin); 
-
-  fhdNdptchPiMCJet = new TH1F("hdNdptchPiMCJet","",100,0,10);
-  fhdNdptchPiMCJet->SetXTitle("p_{T} (GeV/c)");
-  fhdNdptchPiMCJet->SetYTitle("1/N_{jet} dN/dpt (1/GeV)");
-  fhdNdptchPiMCJet->Sumw2();
-  fListOfHistos->Add(fhdNdptchPiMCJet); 
-
-  fhdNdptK0MCJet = new TH1F("hdNdptK0MCJet","",100,0,10);
-  fhdNdptK0MCJet->SetXTitle("p_{T} (GeV/c)");
-  fhdNdptK0MCJet->SetYTitle("1/N_{jet} dN/dpt (1/GeV)");
-  fhdNdptK0MCJet->Sumw2();
-  fListOfHistos->Add(fhdNdptK0MCJet); 
-
-  fhdNdptchKMCJet = new TH1F("hdNdptchKMCJet","",100,0,10);
-  fhdNdptchKMCJet->SetXTitle("p_{T} (GeV/c)");
-  fhdNdptchKMCJet->SetYTitle("1/N_{jet} dN/dpt (1/GeV)");
-  fhdNdptchKMCJet->Sumw2();
-  fListOfHistos->Add(fhdNdptchKMCJet);
-
-  fhdNdptpMCJet = new TH1F("hdNdptpMCJet","",100,0,10);
-  fhdNdptpMCJet->SetXTitle("p_{T} (GeV/c)");
-  fhdNdptpMCJet->SetYTitle("1/N_{jet} dN/dpt (1/GeV)");
-  fhdNdptpMCJet->Sumw2();
-  fListOfHistos->Add(fhdNdptpMCJet); 
-
-  fhdNdptpBarMCJet = new TH1F("hdNdptpBarMCJet","",100,0,10);
-  fhdNdptpBarMCJet->SetXTitle("p_{T} (GeV/c)");
-  fhdNdptpBarMCJet->SetYTitle("1/N_{jet} dN/dpt (1/GeV)");
-  fhdNdptpBarMCJet->Sumw2();
-  fListOfHistos->Add(fhdNdptpBarMCJet);
-
-  fhdNdptLambdaMCJet = new TH1F("hdNdptLambdaMCJet","",100,0,10);
-  fhdNdptLambdaMCJet->SetXTitle("p_{T} (GeV/c)");
-  fhdNdptLambdaMCJet->SetYTitle("1/N_{jet} dN/dpt (1/GeV)");
-  fhdNdptLambdaMCJet->Sumw2();
-  fListOfHistos->Add(fhdNdptLambdaMCJet); 
-
-  fhdNdptLambdaBarMCJet = new TH1F("hdNdptLambdaBarMCJet","",100,0,10);
-  fhdNdptLambdaBarMCJet->SetXTitle("p_{T} (GeV/c)");
-  fhdNdptLambdaBarMCJet->SetYTitle("1/N_{jet} dN/dpt (1/GeV)");
-  fhdNdptLambdaBarMCJet->Sumw2();
-  fListOfHistos->Add(fhdNdptLambdaBarMCJet); 
-
-  fhPIDMC = CreatePIDhisto("hPIDMC");
-  fhPIDMC->Sumw2();
-  fListOfHistos->Add(fhPIDMC); 
-
-  fhPIDMC_quarkEv = CreatePIDhisto("hPIDMC_quarkEv");
-  fhPIDMC_quarkEv->Sumw2();
-  fListOfHistos->Add(fhPIDMC_quarkEv); 
-
-  fhPIDMC_gluonEv = CreatePIDhisto("hPIDMC_gluonEv");
-  fhPIDMC_gluonEv->Sumw2();
-  fListOfHistos->Add(fhPIDMC_gluonEv); 
-
-  fhPIDMCAll = CreatePIDhisto("hPIDMCAll");
-  fhPIDMCAll->Sumw2();
-  fListOfHistos->Add(fhPIDMCAll); 
-
-  fhPIDMCMin = CreatePIDhisto("hPIDMCMin");
-  fhPIDMCMin->Sumw2();
-  fListOfHistos->Add(fhPIDMCMin); 
-
-  fhPIDMCJet = CreatePIDhisto("hPIDMCJet");
-  fhPIDMCJet->Sumw2();
-  fListOfHistos->Add(fhPIDMCJet); 
-
-  fhPIDMCMotherK0 = CreatePIDhisto("hPIDMCMotherK0"); 
-  fhPIDMCMotherK0->Sumw2();
-  fListOfHistos->Add(fhPIDMCMotherK0); 
-
-  fhPIDMCGrandMotherK0 = CreatePIDhisto("hPIDMCGrandMotherK0"); 
-  fhPIDMCGrandMotherK0->Sumw2();
-  fListOfHistos->Add(fhPIDMCGrandMotherK0); 
-
-  fhPIDMCMotherChK = CreatePIDhisto("hPIDMCMotherChK"); 
-  fhPIDMCMotherChK->Sumw2();
-  fListOfHistos->Add(fhPIDMCMotherChK); 
-
-  fhPIDMCMotherK0Trans = CreatePIDhisto("hPIDMCMotherK0Trans"); 
-  fhPIDMCMotherK0Trans->Sumw2();
-  fListOfHistos->Add(fhPIDMCMotherK0Trans); 
-
-  fhPIDMCGrandMotherK0Trans = CreatePIDhisto("hPIDMCGrandMotherK0Trans"); 
-  fhPIDMCGrandMotherK0Trans->Sumw2();
-  fListOfHistos->Add(fhPIDMCGrandMotherK0Trans); 
-
-  fhPIDMCMotherChKTrans = CreatePIDhisto("hPIDMCMotherChKTrans"); 
-  fhPIDMCMotherChKTrans->Sumw2();
-  fListOfHistos->Add(fhPIDMCMotherChKTrans); 
-
-  fhdNdptgammaMC = new TH1F("hdNdptgammaMC","",100,0,10);
-  fhdNdptgammaMC->SetXTitle("p_{T} (GeV/c)");
-  fhdNdptgammaMC->SetYTitle("1/N_{event} dN/dpt (1/GeV)");
-  fhdNdptgammaMC->Sumw2();
-  fListOfHistos->Add(fhdNdptgammaMC); 
-
-  fhdNdptchPiMC  = new TH1F("hdNdptchPiMC","",100,0,10);;
-  fhdNdptchPiMC->SetXTitle("p_{T} (GeV/c)");
-  fhdNdptchPiMC->SetYTitle("1/N_{event} dN/dpt (1/GeV)");
-  fhdNdptchPiMC->Sumw2();
-  fListOfHistos->Add(fhdNdptchPiMC); 
-
-  fhdNdptpi0MC    = new TH1F("hdNdptpi0MC","",100,0,10);;
-  fhdNdptpi0MC->SetXTitle("p_{T} (GeV/c)");
-  fhdNdptpi0MC->SetYTitle("1/N_{event} dN/dpt (1/GeV)");
-  fhdNdptpi0MC->Sumw2();
-  fListOfHistos->Add(fhdNdptpi0MC); 
-
-  fhdNdptK0MC = new TH1F("hdNdptK0MC","",100,0,10);;
-  fhdNdptK0MC->SetXTitle("p_{T} (GeV/c)");
-  fhdNdptK0MC->SetYTitle("1/N_{event} dN/dpt (1/GeV)");
-  fhdNdptK0MC->Sumw2();
-  fListOfHistos->Add(fhdNdptK0MC); 
-
-  fhdNdptchKMC    = new TH1F("hdNdptchKMC","",100,0,10);;
-  fhdNdptchKMC->SetXTitle("p_{T} (GeV/c)");
-  fhdNdptchKMC->SetYTitle("1/N_{event} dN/dpt (1/GeV)");
-  fhdNdptchKMC->Sumw2();
-  fListOfHistos->Add(fhdNdptchKMC); 
-
-  fhdNdptpMC    = new TH1F("hdNdptpMC","",100,0,10);;
-  fhdNdptpMC->SetXTitle("p_{T} (GeV/c)");
-  fhdNdptpMC->SetYTitle("1/N_{event} dN/dpt (1/GeV)");
-  fhdNdptpMC->Sumw2();
-  fListOfHistos->Add(fhdNdptpMC); 
-
-  fhdNdptpBarMC    = new TH1F("hdNdptpBarMC","",100,0,10);;
-  fhdNdptpBarMC->SetXTitle("p_{T} (GeV/c)");
-  fhdNdptpBarMC->SetYTitle("1/N_{event} dN/dpt (1/GeV)");
-  fhdNdptpBarMC->Sumw2();
-  fListOfHistos->Add(fhdNdptpBarMC); 
-
-  fhdNdptLambdaMC    = new TH1F("hdNdptLambdaMC","",100,0,10);;
-  fhdNdptLambdaMC->SetXTitle("p_{T} (GeV/c)");
-  fhdNdptLambdaMC->SetYTitle("1/N_{event} dN/dpt (1/GeV)");
-  fhdNdptLambdaMC->Sumw2();
-  fListOfHistos->Add(fhdNdptLambdaMC); 
-
-  fhdNdptLambdaBarMC    = new TH1F("hdNdptLambdaBarMC","",100,0,10);;
-  fhdNdptLambdaBarMC->SetXTitle("p_{T} (GeV/c)");
-  fhdNdptLambdaBarMC->SetYTitle("1/N_{event} dN/dpt (1/GeV)");
-  fhdNdptLambdaBarMC->Sumw2();
-  fListOfHistos->Add(fhdNdptLambdaBarMC); 
-
-  fhdNdptOmegaMC = new TH1F("hdNdptOmegaMC","",100,0,10);;
-  fhdNdptOmegaMC->SetXTitle("p_{T} (GeV/c)");
-  fhdNdptOmegaMC->SetYTitle("1/N_{event} dN/dpt (1/GeV)");
-  fhdNdptOmegaMC->Sumw2();
-  fListOfHistos->Add(fhdNdptOmegaMC); 
-
-  fhdNdptOmegaBarMC = new TH1F("hdNdptOmegaBarMC","",100,0,10);;
-  fhdNdptOmegaBarMC->SetXTitle("p_{T} (GeV/c)");
-  fhdNdptOmegaBarMC->SetYTitle("1/N_{event} dN/dpt (1/GeV)");
-  fhdNdptOmegaBarMC->Sumw2();
-  fListOfHistos->Add(fhdNdptOmegaBarMC); 
-
-  fhdNdxiMC = new TH1F("hdNdxiMC","",100,0,10); 
-  fhdNdxiMC->SetXTitle("#xi");
-  fhdNdxiMC->Sumw2();
-  fListOfHistos->Add(fhdNdxiMC);            
-
-  fhdNdxiK0MC = new TH1F("hdNdxiK0MC","",100,0,10); 
-  fhdNdxiK0MC->SetXTitle("#xi");
-  fhdNdxiK0MC->Sumw2();
-  fListOfHistos->Add(fhdNdxiK0MC);            
-
-  fhdNdxiK0MCJet = new TH1F("hdNdxiK0MCJet","",100,0,10); 
-  fhdNdxiK0MCJet->SetXTitle("#xi");
-  fhdNdxiK0MCJet->Sumw2();
-  fListOfHistos->Add(fhdNdxiK0MCJet);
-
-  fhdNdzK0MC = new TH1F("hdNdzK0MC","",150,0,1.5);
-  fhdNdzK0MC->SetXTitle("z");
-  fhdNdzK0MC->SetYTitle("1/N_{jet} dN/dz");
-  fhdNdzK0MC->Sumw2();
-  fListOfHistos->Add(fhdNdzK0MC);
-
-  fhdNdzK0MCJet = new TH1F("hdNdzK0MCJet","",150,0,1.5);
-  fhdNdzK0MCJet->SetXTitle("z");
-  fhdNdzK0MCJet->SetYTitle("1/N_{jet} dN/dz");
-  fhdNdzK0MCJet->Sumw2();
-  fListOfHistos->Add(fhdNdzK0MCJet);
-
-  fhdNdptK0MCJetEvt = new TH1F("hdNdptK0MCJetEvt","",100,0,10);;
-  fhdNdptK0MCJetEvt->SetXTitle("p_{T} (GeV/c)");
-  fhdNdptK0MCJetEvt->SetYTitle("1/N_{event} dN/dpt (1/GeV)");
-  fhdNdptK0MCJetEvt->Sumw2();
-  fListOfHistos->Add(fhdNdptK0MCJetEvt);
-
-  fhnJetsAODvsMC = new TH2F("hnJetsAODvsMC","",20,0,20,20,0.,20.); 
-  fhnJetsAODvsMC->SetXTitle("n jets MC");
-  fhnJetsAODvsMC->SetXTitle("n jets AOD");
-  fhnJetsAODvsMC->Sumw2();
-  fListOfHistos->Add(fhnJetsAODvsMC);
-
-  fhLeadingPtAODvsMC = new TH2F("hLeadingPtAODvsMC","",20,0,20,20,0.,20.); 
-  fhLeadingPtAODvsMC->SetXTitle("p_{T} MC (GeV/c)");
-  fhLeadingPtAODvsMC->SetYTitle("p_{T} AOD (GeV/c)");
-  fhLeadingPtAODvsMC->Sumw2();
-  fListOfHistos->Add(fhLeadingPtAODvsMC);
-
-  fhLeadingEtaAODvsMC = new TH2F("hLeadingEtaAODvsMC","",20,-1,1,20,-1.,1.); 
-  fhLeadingEtaAODvsMC->SetXTitle("#eta MC");
-  fhLeadingEtaAODvsMC->SetYTitle("#eta AOD");
-  fhLeadingEtaAODvsMC->Sumw2();
-  fListOfHistos->Add(fhLeadingEtaAODvsMC);
-
-  fhLeadingPhiAODvsMC = new TH2F("hLeadingPhiAODvsMC","",63,0,6.3,63,0.,6.3); 
-  fhLeadingPhiAODvsMC->SetXTitle("#phi MC");
-  fhLeadingPhiAODvsMC->SetYTitle("#phi AOD");
-  fhLeadingPhiAODvsMC->Sumw2();
-  fListOfHistos->Add(fhLeadingPhiAODvsMC);
-
-
-  fhnTracksLeadingAODvsMC = new TH2F("hnTracksLeadingAODvsMC","",20,0.,20,20,-0.5,19.5);
-  fhnTracksLeadingAODvsMC->SetXTitle("nTracks MC");
-  fhnTracksLeadingAODvsMC->SetYTitle("nTracks AOD");
-  fhnTracksLeadingAODvsMC->Sumw2();
-  fListOfHistos->Add(fhnTracksLeadingAODvsMC);
-
-  fhLeadingdRAODMC = new TH1F("hLeadingdRAODMC","",40,0.,4);
-  fhLeadingdRAODMC->SetXTitle("#Delta R");
-  fhLeadingdRAODMC->Sumw2();
-  fListOfHistos->Add(fhLeadingdRAODMC);
-
-
-  fhLeadingPtAODvsMCdRcut = new TH2F("hLeadingPtAODvsMCdRcut","",40,0,20,40,0.,20.); 
-  fhLeadingPtAODvsMCdRcut->SetXTitle("p_{T} MC (GeV/c)");
-  fhLeadingPtAODvsMCdRcut->SetYTitle("p_{T} AOD (GeV/c)");
-  fhLeadingPtAODvsMCdRcut->Sumw2();
-  fListOfHistos->Add(fhLeadingPtAODvsMCdRcut);
-
-
-  fhdnTracksVsdPtLeadingAODMC = new TH2F("hdnTracksVsdPtLeadingAODMC","",40,-10.,10,20,-10.5,9.5);
-  fhdnTracksVsdPtLeadingAODMC->SetXTitle("#Delta Pt AOD-MC (GeV/c)");
-  fhdnTracksVsdPtLeadingAODMC->SetYTitle("#Delta N AOD-MC");
-  fhdnTracksVsdPtLeadingAODMC->Sumw2();
-  fListOfHistos->Add(fhdnTracksVsdPtLeadingAODMC);
-
-  fhnTracksJetVsPtAOD = new TH2F("hnTracksJetVsPtAOD","",50,0,50,20,-0.5,19.5);
-  fhnTracksJetVsPtAOD->SetXTitle("p_{T} (GeV/c)");
-  fhnTracksJetVsPtAOD->SetYTitle("nTracks");
-  fhnTracksJetVsPtAOD->Sumw2();
-  fListOfHistos->Add(fhnTracksJetVsPtAOD);
-
-  fhnTracksJetVsPtAODquarkEv = new TH2F("hnTracksJetVsPtAODquarkEv","",50,0,50,20,-0.5,19.5);
-  fhnTracksJetVsPtAODquarkEv->SetXTitle("p_{T} (GeV/c)");
-  fhnTracksJetVsPtAODquarkEv->SetYTitle("nTracks");
-  fhnTracksJetVsPtAODquarkEv->Sumw2();
-  fListOfHistos->Add(fhnTracksJetVsPtAODquarkEv);
-
-  fhRadiusJetVsPtAOD = new TH2F("hRadiusJetVsPtAOD","",50,0,50,10,0.,1);
-  fhRadiusJetVsPtAOD->SetXTitle("p_{T} (GeV/c)");
-  fhRadiusJetVsPtAOD->SetYTitle("radius");
-  fhRadiusJetVsPtAOD->Sumw2();
-  fListOfHistos->Add(fhRadiusJetVsPtAOD);
-
-  fhnTracksJetVsPtMC = new TH2F("hnTracksJetVsPtMC","",50,0,50,20,-0.5,19.5);
-  fhnTracksJetVsPtMC->SetXTitle("p_{T} (GeV/c)");
-  fhnTracksJetVsPtMC->SetYTitle("nTracks");
-  fhnTracksJetVsPtMC->Sumw2();
-  fListOfHistos->Add(fhnTracksJetVsPtMC);
-
-  fhnTracksJetVsPtMCquarkEv = new TH2F("hnTracksJetVsPtMCquarkEv","",50,0,50,20,-0.5,19.5);
-  fhnTracksJetVsPtMCquarkEv->SetXTitle("p_{T} (GeV/c)");
-  fhnTracksJetVsPtMCquarkEv->SetYTitle("nTracks");
-  fhnTracksJetVsPtMCquarkEv->Sumw2();
-  fListOfHistos->Add(fhnTracksJetVsPtMCquarkEv);
-
-  fhRadiusJetVsPtMC = new TH2F("hRadiusJetVsPtMC","",50,0,50,10,0.,1);
-  fhRadiusJetVsPtMC->SetXTitle("p_{T} (GeV/c)");
-  fhRadiusJetVsPtMC->SetYTitle("radius");
-  fhRadiusJetVsPtMC->Sumw2();
-  fListOfHistos->Add(fhRadiusJetVsPtMC);
-
-  fhnTracksJetVsPtMCK0 = new TH2F("hnTracksJetVsPtMCK0","",50,0,50,20,-0.5,19.5);
-  fhnTracksJetVsPtMCK0->SetXTitle("p_{T} (GeV/c)");
-  fhnTracksJetVsPtMCK0->SetYTitle("nTracks");
-  fhnTracksJetVsPtMCK0->Sumw2();
-  fListOfHistos->Add(fhnTracksJetVsPtMCK0);
-
-
-  fhnTracksJetVsPtMCK0quarkEv = new TH2F("hnTracksJetVsPtMCK0quarkEv","",50,0,50,20,-0.5,19.5);
-  fhnTracksJetVsPtMCK0quarkEv->SetXTitle("p_{T} (GeV/c)");
-  fhnTracksJetVsPtMCK0quarkEv->SetYTitle("nTracks");
-  fhnTracksJetVsPtMCK0quarkEv->Sumw2();
-  fListOfHistos->Add(fhnTracksJetVsPtMCK0quarkEv);
-
-  fhRadiusJetVsPtMCK0 = new TH2F("hRadiusJetVsPtMCK0","",50,0,50,10,0.,1);
-  fhRadiusJetVsPtMCK0->SetXTitle("p_{T} (GeV/c)");
-  fhRadiusJetVsPtMCK0->SetYTitle("radius");
-  fhRadiusJetVsPtMCK0->Sumw2();
-  fListOfHistos->Add(fhRadiusJetVsPtMCK0);
-
-  fhnTracksJetVsPtAODK0 = new TH2F("hnTracksJetVsPtAODK0","",50,0,50,20,-0.5,19.5);
-  fhnTracksJetVsPtAODK0->SetXTitle("p_{T} (GeV/c)");
-  fhnTracksJetVsPtAODK0->SetYTitle("nTracks AODK0");
-  fhnTracksJetVsPtAODK0->Sumw2();
-  fListOfHistos->Add(fhnTracksJetVsPtAODK0);
-
-  fhnTracksJetVsPtAODK0quarkEv = new TH2F("hnTracksJetVsPtAODK0quarkEv","",50,0,50,20,-0.5,19.5);
-  fhnTracksJetVsPtAODK0quarkEv->SetXTitle("p_{T} (GeV/c)");
-  fhnTracksJetVsPtAODK0quarkEv->SetYTitle("nTracks AODK0quarkEv");
-  fhnTracksJetVsPtAODK0quarkEv->Sumw2();
-  fListOfHistos->Add(fhnTracksJetVsPtAODK0quarkEv);
-
-  fhRadiusJetVsPtAODK0 = new TH2F("hRadiusJetVsPtAODK0","",50,0,50,10,0.,1);
-  fhRadiusJetVsPtAODK0->SetXTitle("p_{T} (GeV/c)");
-  fhRadiusJetVsPtAODK0->SetYTitle("radius");
-  fhRadiusJetVsPtAODK0->Sumw2();
-  fListOfHistos->Add(fhRadiusJetVsPtAODK0);
-
-  fhnTracksJetVsPtAODpKch = new TH2F("hnTracksJetVsPtAODpKch","",50,0,50,20,-0.5,19.5);
-  fhnTracksJetVsPtAODpKch->SetXTitle("p_{T} (GeV/c)");
-  fhnTracksJetVsPtAODpKch->SetYTitle("nTracks AODpKch");
-  fhnTracksJetVsPtAODpKch->Sumw2();
-  fListOfHistos->Add(fhnTracksJetVsPtAODpKch);
-
-  fhRadiusJetVsPtAODpKch = new TH2F("hRadiusJetVsPtAODpKch","",50,0,50,20,-0.5,19.5);
-  fhRadiusJetVsPtAODpKch->SetXTitle("p_{T} (GeV/c)");
-  fhRadiusJetVsPtAODpKch->SetYTitle("Radius AODpKch");
-  fhRadiusJetVsPtAODpKch->Sumw2();
-  fListOfHistos->Add(fhRadiusJetVsPtAODpKch);
-
-  fhPythiaProcess      = CreatePythiaIDhisto("hPythiaProcess");
-  fListOfHistos->Add(fhPythiaProcess);  
-
-  fhPythiaProcessK0   = CreatePythiaIDhisto("hPythiaProcessK0");
-  fListOfHistos->Add(fhPythiaProcessK0);  
-
-  fhPythiaProcessKch   = CreatePythiaIDhisto("hPythiaProcessKch");
-  fListOfHistos->Add(fhPythiaProcessKch);  
-
-  fhPythiaProcessp    = CreatePythiaIDhisto("hPythiaProcessp");
-  fListOfHistos->Add(fhPythiaProcessp);  
-  
-  fhPythiaProcesspbar = CreatePythiaIDhisto("hPythiaProcesspbar");
-  fListOfHistos->Add(fhPythiaProcesspbar);  
-  fhdNdzJets5to10 = new TH1F("hdNdzJets5to10","",25,0,1.25);
-  fhdNdzJets5to10->SetXTitle("z");
-  fhdNdzJets5to10->SetYTitle("dN/dz");
-  fhdNdzJets5to10->Sumw2();
-  fListOfHistos->Add(fhdNdzJets5to10);
-
-  fhdNdzJets10to20 = new TH1F("hdNdzJets10to20","",25,0,1.25);
-  fhdNdzJets10to20->SetXTitle("z");
-  fhdNdzJets10to20->SetYTitle("dN/dz");
-  fhdNdzJets10to20->Sumw2();
-  fListOfHistos->Add(fhdNdzJets10to20);
-
-  fhdNdzJets20to30 = new TH1F("hdNdzJets20to30","",25,0,1.25);
-  fhdNdzJets20to30->SetXTitle("z");
-  fhdNdzJets20to30->SetYTitle("dN/dz");
-  fhdNdzJets20to30->Sumw2();
-  fListOfHistos->Add(fhdNdzJets20to30);
-
-  fhdNdzJets30to40 = new TH1F("hdNdzJets30to40","",25,0,1.25);
-  fhdNdzJets30to40->SetXTitle("z");
-  fhdNdzJets30to40->SetYTitle("dN/dz");
-  fhdNdzJets30to40->Sumw2();
-  fListOfHistos->Add(fhdNdzJets30to40);
-
-  fhdNdzJets40to60 = new TH1F("hdNdzJets40to60","",25,0,1.25);
-  fhdNdzJets40to60->SetXTitle("z");
-  fhdNdzJets40to60->SetYTitle("dN/dz");
-  fhdNdzJets40to60->Sumw2();
-  fListOfHistos->Add(fhdNdzJets40to60);
+  Int_t nJCuts = GetListOfJets(fJetsRecCuts, kJetsRecAcceptance);
+  Int_t nRecJetsCuts = 0;
+  if(nJCuts>=0) nRecJetsCuts = fJetsRecCuts->GetEntries();
+  if(fDebug>2)Printf("%s:%d Selected Rec jets after cuts: %d %d",(char*)__FILE__,__LINE__,nJCuts,nRecJetsCuts);
+  if(nRecJetsCuts != nJCuts) Printf("%s:%d Mismatch selected Rec jets after cuts: %d %d",(char*)__FILE__,__LINE__,nJCuts,nRecJetsCuts);
+  fh1nRecJetsCuts->Fill(nRecJetsCuts);
   
 
-  fhdNdxiJets5to10 = new TH1F("hdNdxiJets5to10","",70,0,7);
-  fhdNdxiJets5to10->SetXTitle("z");
-  fhdNdxiJets5to10->SetYTitle("dN/dz");
-  fhdNdxiJets5to10->Sumw2();
-  fListOfHistos->Add(fhdNdxiJets5to10);
-
-  fhdNdxiJets10to20 = new TH1F("hdNdxiJets10to20","",70,0,7);
-  fhdNdxiJets10to20->SetXTitle("z");
-  fhdNdxiJets10to20->SetYTitle("dN/dz");
-  fhdNdxiJets10to20->Sumw2();
-  fListOfHistos->Add(fhdNdxiJets10to20);
-
-  fhdNdxiJets20to30 = new TH1F("hdNdxiJets20to30","",70,0,7);
-  fhdNdxiJets20to30->SetXTitle("z");
-  fhdNdxiJets20to30->SetYTitle("dN/dz");
-  fhdNdxiJets20to30->Sumw2();
-  fListOfHistos->Add(fhdNdxiJets20to30);
-
-  fhdNdxiJets30to40 = new TH1F("hdNdxiJets30to40","",70,0,7);
-  fhdNdxiJets30to40->SetXTitle("z");
-  fhdNdxiJets30to40->SetYTitle("dN/dz");
-  fhdNdxiJets30to40->Sumw2();
-  fListOfHistos->Add(fhdNdxiJets30to40);
-
-  fhdNdxiJets40to60 = new TH1F("hdNdxiJets40to60","",70,0,7);
-  fhdNdxiJets40to60->SetXTitle("z");
-  fhdNdxiJets40to60->SetYTitle("dN/dz");
-  fhdNdxiJets40to60->Sumw2();
-  fListOfHistos->Add(fhdNdxiJets40to60);
-
-  fhdNdptTracksJetPt5to10 = new TH1F("hdNdptTracksJetPt5to10","",250,0,25);
-  fhdNdptTracksJetPt5to10->SetXTitle("p_{T} (GeV)");
-  fhdNdptTracksJetPt5to10->SetYTitle("dN/dp_{T} 1/GeV");
-  fhdNdptTracksJetPt5to10->Sumw2();
-  fListOfHistos->Add(fhdNdptTracksJetPt5to10);
-
-  fhdNdptTracksJetPt10to20 = new TH1F("hdNdptTracksJetPt10to20","",25,0,25);
-  fhdNdptTracksJetPt10to20->SetXTitle("p_{T} (GeV)");
-  fhdNdptTracksJetPt10to20->SetYTitle("dN/dp_{T} 1/GeV");
-  fhdNdptTracksJetPt10to20->Sumw2();
-  fListOfHistos->Add(fhdNdptTracksJetPt10to20);
-
-  fhdNdptTracksJetPt20to30 = new TH1F("hdNdptTracksJetPt20to30","",25,0,25);
-  fhdNdptTracksJetPt20to30->SetXTitle("p_{T} (GeV)");
-  fhdNdptTracksJetPt20to30->SetYTitle("dN/dp_{T} 1/GeV");
-  fhdNdptTracksJetPt20to30->Sumw2();
-  fListOfHistos->Add(fhdNdptTracksJetPt20to30);
-
-  fhdNdptTracksJetPt30to40 = new TH1F("hdNdptTracksJetPt30to40","",25,0,25);
-  fhdNdptTracksJetPt30to40->SetXTitle("p_{T} (GeV)");
-  fhdNdptTracksJetPt30to40->SetYTitle("dN/dp_{T} 1/GeV");
-  fhdNdptTracksJetPt30to40->Sumw2();
-  fListOfHistos->Add(fhdNdptTracksJetPt30to40);
-
-  fhdNdptTracksJetPt40to60 = new TH1F("hdNdptTracksJetPt40to60","",25,0,25);
-  fhdNdptTracksJetPt40to60->SetXTitle("p_{T} (GeV)");
-  fhdNdptTracksJetPt40to60->SetYTitle("dN/dp_{T} 1/GeV");
-  fhdNdptTracksJetPt40to60->Sumw2();
-  fListOfHistos->Add(fhdNdptTracksJetPt40to60);
-
-
-  fh1Xsec = new TProfile("h1Xsec","xsec from pyxsec.root",1,0,1); 
-  fh1Xsec->SetXTitle("<#sigma>");
-  fh1Xsec->Sumw2();
-  fListOfHistos->Add( fh1Xsec ); 
-  
-  fh1Trials = new TH1F("h1Trials","trials from pyxsec.root",1,0,1);
-  fh1Trials->SetXTitle("#sum{ntrials}");
-  fh1Trials->Sumw2();
-  fListOfHistos->Add( fh1Trials ); 
+  //____ fetch particles __________________________________________________________
  
-//   fSettingsTree   = new TTree("JetChemAnalysisSettings","Analysis Settings");
-//   fSettingsTree->Branch("fUseLOConeJets",&fUseLOConeJets,"UseLOConeJets/O");
-//   fSettingsTree->Branch("fUseLOConeMCJets",&fUseLOConeMCJets,"UseLOConeMCJets/O");
-//   fSettingsTree->Branch("fUsePythiaJets",&fUsePythiaJets,"UsePythiaJets/O");
-//   fSettingsTree->Branch("fConeRadius", &fConeRadius,"Rad/D");
-//   fSettingsTree->Branch("fTrackPtCutJF", &fTrackPtCutJF,"TrackPt/D");
-//   fSettingsTree->Branch("fFilterBitJF",&fFilterBitJF,"FilterBitJF/i");
-//   fSettingsTree->Branch("fRequireITSRefitJF",&fRequireITSRefitJF,"RequireITSRefitJF/O");
-//   fSettingsTree->Branch("fRejectK0TracksJF",&fRejectK0TracksJF,"RejectK0TracksJF/O");
-//   fSettingsTree->Branch("fJetPtCut",&fJetPtCut,"JetPtCut/D");
-//   fSettingsTree->Branch("fJetEtaCut",&fJetEtaCut,"JetEtaCut/D");
-//   fSettingsTree->Branch("fFilterBit",&fFilterBit,"FilterBit/i");
-//   fSettingsTree->Branch("fTrackPtCut",&fTrackPtCut,"TrackPtCut/D");
-//   fSettingsTree->Branch("fTrackEtaCut",&fTrackEtaCut,"TrackEtaCut/D");
-//   fSettingsTree->Branch("fUseOnFlyV0s",&fUseOnFlyV0s,"UseOnFlyV0s/O");
-//   fSettingsTree->Branch("fCutnSigdEdx",&fCutnSigdEdx,"CutnSigdEdx/D");
-//   fSettingsTree->Branch("fUseAODMCTracksForUE",&fUseAODMCTracksForUE,"UseAODMCTracksForUE/O");
-//   fSettingsTree->Branch("fAreaReg",&fAreaReg,"AreaReg/D");
-//   fSettingsTree->Branch("fAvgTrials",&fAvgTrials,"AvgTrials/D");
-  // fListOfHistos->Add(fSettingsTree);
-  
-}
+  Int_t nTCuts = GetListOfTracks(fTracksRecCuts, kTrackAODCuts);
+  if(fDebug>2)Printf("%s:%d Selected Rec tracks after cuts: %d %d",(char*)__FILE__,__LINE__,nTCuts,fTracksRecCuts->GetEntries());
+  if(fTracksRecCuts->GetEntries() != nTCuts) 
+    Printf("%s:%d Mismatch selected Rec tracks after cuts: %d %d",(char*)__FILE__,__LINE__,nTCuts,fTracksRecCuts->GetEntries());
+  fh1EvtMult->Fill(fTracksRecCuts->GetEntries());
 
-//____________________________________________________________________
-void  AliAnalysisTaskJetChem::Terminate(Option_t */*option*/){
+  Int_t nK0s = GetListOfK0s(fListK0s,fK0Type);
+  if(fDebug>2)Printf("%s:%d Selected Rec tracks after cuts: %d %d",(char*)__FILE__,__LINE__,nK0s,fListK0s->GetEntries());
+  if(nK0s != fListK0s->GetEntries()) Printf("%s:%d Mismatch selected K0s: %d %d",(char*)__FILE__,__LINE__,nK0s,fListK0s->GetEntries());
+  fh1K0Mult->Fill(fListK0s->GetEntries());
 
-  // Terminate analysis
-  //
-  
-  WriteSettings();
+   // ___ V0 QA + K0 pt spectra all events _______________________________________________
   
-  // Update pointers reading them from the output slot
-  fListOfHistos = dynamic_cast<TList*> (GetOutputData(1));
-  if( !fListOfHistos ) {
-    AliError("Histogram List is not available");
-    return;
-  }
-
-  fhLeadingPt          = (TH1F*)fListOfHistos->FindObject("hLeadingPt");
-  fhRegionSumPtMaxVsEt = (TH1F*)fListOfHistos->FindObject("hRegionSumPtMaxVsEt");
-  fhRegionSumPtMinVsEt = (TH1F*)fListOfHistos->FindObject("hRegionSumPtMinVsEt");
-  fhRegionMultMaxVsEt  = (TH1F*)fListOfHistos->FindObject("hRegionMultMaxVsEt");
-  fhRegionMultMinVsEt  = (TH1F*)fListOfHistos->FindObject("hRegionMultMinVsEt");
-    
-   
-  fhRegionSumPtMaxVsEt->Divide(fhLeadingPt);
-  fhRegionSumPtMinVsEt->Divide(fhLeadingPt);
-  fhRegionMultMaxVsEt->Divide(fhLeadingPt);
-  fhRegionMultMinVsEt->Divide(fhLeadingPt);
-
+  if(fListK0s->GetEntries()>0){
     
-  //Get Normalization
-  fh1Xsec                = (TProfile*) fListOfHistos->FindObject("fh1Xsec");
-  fh1Trials              = (TH1F*)fListOfHistos->FindObject("fh1Trials");
-  
-  //std::cout<<" fh1Xsec "<<fh1Xsec<<"  fh1Trials "<<fh1Trials<<std::endl;
-
-  if(fh1Xsec && fh1Trials){
+    for(Int_t it=0; it<fListK0s->GetSize(); ++it){ // loop all K0s 
+      
+      AliAODv0* v0 = dynamic_cast<AliAODv0*>(fListK0s->At(it));
+      
+      Float_t trackPt       = v0->Pt();
+      Bool_t incrementJetPt = (it==0) ? kTRUE : kFALSE;
+      Double_t invM         = v0->MassK0Short();
+      Double_t jetPt        = fFFIMJetPtMin; // assign pro forma jet energy
 
-    Double_t xsec = fh1Xsec->GetBinContent(1);
-    Double_t ntrials = fh1Trials->GetBinContent(1);
-    Double_t normFactor = xsec/ntrials;
-    Printf("xSec %f nTrials %f Norm %f \n",xsec,ntrials,normFactor);
-    fhLeadingPt->Scale(normFactor);
+      fV0QAK0->FillTrackQA(v0->Eta(), TVector2::Phi_0_2pi(v0->Phi()), v0->Pt()); 
+      fFFHistosIMK0AllEvt->FillFF(trackPt, invM, jetPt, incrementJetPt);
+    }
   }
 
-  
-  if(fDebug > 1) AliInfo("End analysis");
-}
 
-// ----------------------------------------------------------------------------
+  //____ fill FF histos  __________________________________________________________
 
-void  AliAnalysisTaskJetChem::WriteSettings(){ 
+  for(Int_t ij=0; ij<nRecJetsCuts; ++ij){
 
-  // fill settings tree
-  // not on GRID !
+    AliAODJet* jet = dynamic_cast<AliAODJet*>(fJetsRecCuts->At(ij));
+    
+    if(ij==0){ // leading jet
+      
+      TList* jettracklist = new TList();
+      Double_t sumPt      = 0.;
+      
+      if(GetFFRadius()<=0){
+       GetJetTracksTrackrefs(jettracklist, jet);
+      } else {
+       GetJetTracksPointing(fTracksRecCuts, jettracklist, jet, GetFFRadius(), sumPt);
+      }
+      
+      for(Int_t it=0; it<jettracklist->GetSize(); ++it){
 
-  //fSettingsTree->Fill(); 
+       AliVParticle* trackVP = dynamic_cast<AliVParticle*>(jettracklist->At(it));
+       
+       Float_t jetPt   = jet->Pt();
+       Float_t trackPt = trackVP->Pt();
+       
+       Bool_t incrementJetPt = (it==0) ? kTRUE : kFALSE;
+       
+       fFFHistosRecCuts->FillFF(trackPt, jetPt, incrementJetPt);
+       if(nK0s>0) fFFHistosRecCutsK0Evt->FillFF(trackPt, jetPt, incrementJetPt);
+      }
+      
+      delete jettracklist;
+    }
+  }
 
-}
 
-// ---------------------------------------------------------------------------
+  for(Int_t ij=0; ij<nRecJetsCuts; ++ij){ // jet loop
 
-Bool_t AliAnalysisTaskJetChem::IsK0InvMass(const Double_t mass) const {
+    AliAODJet* jet = dynamic_cast<AliAODJet*>(fJetsRecCuts->At(ij));
+    
+    if(ij==0){ // leading jet
+      
+      //fQAJetHistosRecCutsLeading->FillJetQA( jet->Eta(), TVector2::Phi_0_2pi(jet->Phi()), jet->Pt() );
+      
+      Float_t jetPt   = jet->Pt();
 
-  // K0 mass ?  Use STAR params for bin counting & mass fit 
+      for(Int_t it=0; it<fListK0s->GetSize(); ++it){ // loop all K0s 
 
-  const Double_t massK0    = 0.497;        // from fits
-  const Double_t sigmaK0   = 0.0046;
-  const Double_t nSigmaSignal = 3.5; // STAR parameters for bin counting
+       AliAODv0* v0 = dynamic_cast<AliAODv0*>(fListK0s->At(it));
 
-  if((massK0-nSigmaSignal*sigmaK0)<=mass &&  
-     mass<=(massK0 + nSigmaSignal*sigmaK0)) return kTRUE;
+       Double_t v0Mom[3];
+       v0->PxPyPz(v0Mom);
+       TVector3 v0MomVect(v0Mom);
 
-  return kFALSE;
-}
+       Double_t dPhiJetK0 = (jet->MomentumVector()->Vect()).DeltaPhi(v0MomVect);
+       
+       Float_t trackPt       = v0->Pt();
+       Bool_t incrementJetPt = (it==0) ? kTRUE : kFALSE;
+       Double_t invM         = v0->MassK0Short();
+       
+       fFFHistosIMK0Jet->FillFF(trackPt, invM, jetPt, incrementJetPt);
+       fFFHistosPhiCorrIMK0->FillPhiCorr(trackPt,TVector2::Phi_0_2pi(dPhiJetK0),invM);
 
-// ---------------------------------------------------------------------------
+       if(dPhiJetK0<fh1dPhiJetK0->GetXaxis()->GetXmin()) dPhiJetK0 += 2*TMath::Pi();
+       fh1dPhiJetK0->Fill(dPhiJetK0);
 
-Bool_t AliAnalysisTaskJetChem::IsLambdaInvMass(const Double_t mass) const{
+      }
+      if(fListK0s->GetSize() == 0){ // no K0: increment jet pt spectrum 
 
-  // Lambda mass ?  
+       Bool_t incrementJetPt = kTRUE;
+       fFFHistosIMK0Jet->FillFF(-1, -1, jetPt, incrementJetPt);
+      }
+      
 
+      TList* jetConeK0list = new TList();
+      Double_t sumPt   = 0.;
+      
+      GetJetTracksPointing(fListK0s, jetConeK0list, jet, GetFFRadius(), sumPt);
+      if(fDebug>2)Printf("%s:%d nK0s total: %d, in jet cone: %d,FFRadius %f ",(char*)__FILE__,__LINE__,nK0s,jetConeK0list->GetEntries(),GetFFRadius());
+      
+      for(Int_t it=0; it<jetConeK0list->GetSize(); ++it){ // loop K0s in jet cone
 
-  if(1.1<mass && mass<1.13) return kTRUE; // FIXME - adjust range from fit
-  return kFALSE;
-}
+       AliAODv0* v0 = dynamic_cast<AliAODv0*>(jetConeK0list->At(it));
 
-// ---------------------------------------------------------------------------
+       Double_t invM           = v0->MassK0Short();
+       Float_t  trackPt        = v0->Pt();
+       Bool_t   incrementJetPt = (it==0) ? kTRUE : kFALSE;
+       
+       //std::cout<<" trackPt "<<trackPt<<" invM "<<invM<<std::endl;
+       fFFHistosIMK0Cone->FillFF(trackPt, invM, jetPt, incrementJetPt);
+      }
+     if(jetConeK0list->GetSize() == 0){ // no K0: increment jet pt spectrum 
 
-Bool_t AliAnalysisTaskJetChem::IsAcceptedDCAK0(/*const Double_t dca*/) const{
+       Bool_t incrementJetPt = kTRUE;
+       fFFHistosIMK0Cone->FillFF(-1, -1, jetPt, incrementJetPt);
+      }
 
-  // DCA cut 
+      delete jetConeK0list;
+    }
+  }
+  
+  fTracksRecCuts->Clear();
+  fJetsRecCuts->Clear();
+  fListK0s->Clear();
 
-  return kTRUE; // FIXME - adjust cut
+  //Post output data.
+  PostData(1, fCommonHistList);    
 }
 
-
-// ---------------------------------------------------------------------------
-
-Bool_t AliAnalysisTaskJetChem::IsAcceptedDCALambda(/*const Double_t dca*/) const {
-
-  // DCA cut
-
-  return kTRUE; // FIXME - adjust cut
+// ____________________________________________________________________________________________
+void AliAnalysisTaskJetChem::SetProperties(TH3F* h,const char* x, const char* y, const char* z)
+{
+  //Set properties of histos (x,y and z title)
+
+  h->SetXTitle(x);
+  h->SetYTitle(y);
+  h->SetZTitle(z);
+  h->GetXaxis()->SetTitleColor(1);
+  h->GetYaxis()->SetTitleColor(1);
+  h->GetZaxis()->SetTitleColor(1);
 }
 
-// ---------------------------------------------------------------------------
-
+// ____________________________________________________________________________________________________________________________________________
 Bool_t AliAnalysisTaskJetChem::IsAccepteddEdx(const Double_t mom,const Double_t signal, AliPID::EParticleType n, const Double_t cutnSig) const{
 
   // apply TPC dE/dx cut similar as in AliTPCpidESD 
@@ -2730,56 +1010,59 @@ Bool_t AliAnalysisTaskJetChem::IsAccepteddEdx(const Double_t mom,const Double_t
   return kTRUE;
 }
 
-// ----------------------------------------------------------------------------
+//___________________________________________________________________
+Bool_t AliAnalysisTaskJetChem::IsK0InvMass(const Double_t mass) const
+{
+  // K0 mass ? Use FF histo limits
+  
+  if(fFFIMInvMMin <= mass && mass <= fFFIMInvMMax) return kTRUE;
 
-void AliAnalysisTaskJetChem::CheckV0s(AliAODJet* jetVect, Int_t maxPtRegionIndex,Bool_t& foundK0){
+  return kFALSE;
+}
 
-  // loop over AOD V0s, fill masses etc
+//_____________________________________________________________________________________
+Int_t AliAnalysisTaskJetChem::GetListOfK0s(TList *list, const Int_t type)
+{
+  // fill list of V0s selected according to type
 
-  Int_t nV0 = fAOD->GetNumberOfV0s();
-  fhNV0s->Fill(nV0);  
+  if(!list){
+    if(fDebug>1) Printf("%s:%d no input list", (char*)__FILE__,__LINE__);
+    return -1;
+  }
+
+  if(type==kTrackUndef) return 0;
 
   for(int i=0; i<fAOD->GetNumberOfV0s(); i++){ // loop over V0s
     
     AliAODv0* v0 = fAOD->GetV0(i);
 
     Bool_t isOnFly = v0->GetOnFlyStatus();
-    isOnFly ? fhV0onFly->Fill(1) : fhV0onFly->Fill(0);
     
-    if((fUseOnFlyV0s && isOnFly)  || (!fUseOnFlyV0s && !isOnFly) ){ // 'offline' V0s :  using vertex tracks, on-fly: take decay vertex into account during mom fit 
-      
-      Double_t massK0         = v0->MassK0Short();
-      Double_t massLambda     = v0->MassLambda();
-      Double_t massAntiLambda = v0->MassAntiLambda();
+    if(!isOnFly &&  (type == kOnFly || type == kOnFlyPID || type == kOnFlydEdx || type == kOnFlyPrim)) continue; 
+    if( isOnFly &&  (type == kOffl  || type == kOfflPID  || type == kOffldEdx  || type == kOfflPrim))  continue; 
 
-      fhV0InvMassK0->Fill(massK0);
-      fhV0InvMassLambda->Fill(massLambda);
-      fhV0InvMassAntiLambda->Fill(massAntiLambda);
+    Double_t massK0 = v0->MassK0Short();
 
-      if(!(IsK0InvMass(massK0))) continue; // moved invMass cut for HI - otherwise too slow
+    if(!(IsK0InvMass(massK0))) continue; // moved invMass cut for HI - otherwise too slow
 
-      //Double_t radiusV0       = v0->RadiusV0(); // slow
+    if(type == kOnFlyPID || type == kOfflPID){
+      
+      AliAODTrack *trackPos = (AliAODTrack *) (v0->GetSecondaryVtx()->GetDaughter(0)); // slow 
+      AliAODTrack *trackNeg = (AliAODTrack *) (v0->GetSecondaryVtx()->GetDaughter(1)); // slow  
+     
+      // AOD pid - cuts strongly into signal
 
-      Double_t etaV0          = v0->Eta();
-      Double_t fDCAV0ToVertex = v0->DcaV0ToPrimVertex();
-      Double_t fDCADaughters  = v0->DcaV0Daughters();
+      AliAODTrack::AODTrkPID_t mpPIDNeg = trackNeg->GetMostProbablePID();
+      AliAODTrack::AODTrkPID_t mpPIDPos = trackPos->GetMostProbablePID();
+       
+      if(!( (mpPIDNeg == AliAODTrack::kPion) && (mpPIDPos == AliAODTrack::kPion) ) ) continue;
+    }
+   
+    if(type == kOnFlydEdx || type == kOffldEdx){
 
-      Double_t v0Mom[3];
-      v0->PxPyPz(v0Mom);
-      TVector3 v0MomVect(v0Mom);
-      
       AliAODTrack *trackPos = (AliAODTrack *) (v0->GetSecondaryVtx()->GetDaughter(0)); // slow 
       AliAODTrack *trackNeg = (AliAODTrack *) (v0->GetSecondaryVtx()->GetDaughter(1)); // slow  
 
-      Double_t posMom[3], negMom[3];
-      trackPos->PxPyPz(posMom);
-      trackNeg->PxPyPz(negMom);
-      
-      TVector3 posMomVec(posMom);
-      TVector3 negMomVec(negMom);
-      
-      Double_t dROpanV0 = posMomVec.DeltaR(negMomVec);
-      
       AliAODPid*  aodPidPos = trackPos->GetDetPid();
       AliAODPid*  aodPidNeg = trackNeg->GetDetPid();
 
@@ -2788,1094 +1071,29 @@ void AliAnalysisTaskJetChem::CheckV0s(AliAODJet* jetVect, Int_t maxPtRegionIndex
 
       Double_t momPos  = trackPos->P();
       Double_t momNeg  = trackNeg->P();
-           
-      fhdEdxVsMomV0->Fill(momPos,dEdxPos);
-      fhdEdxVsMomV0->Fill(momNeg,dEdxNeg);
-      Double_t pV0       = TMath::Sqrt(v0->Ptot2V0());
-      Double_t ptV0      = TMath::Sqrt(v0->Pt2V0());
-     
-    
-      fhV0DCADaughters->Fill(fDCADaughters);
-      //fhV0Radius->Fill(radiusV0);
-      fhV0DCAToVertex->Fill(fDCAV0ToVertex);
-      fhdNdptV0->Fill(ptV0);
-
-      fhV0PtVsInvMassK0->Fill(massK0,ptV0); 
-
-
-      // K0 signal before cuts 
-
-      if(IsK0InvMass(massK0)){
-       
-       fhdNdptK0->Fill(ptV0);
-       
-       fhPtVsEtaK0->Fill(etaV0,ptV0);
-
-       Double_t dRV0MC = AssociateV0MC(&v0MomVect,310); // K0
-       if(dRV0MC < 0) dRV0MC = 0.99;
-       fhdRV0MC->Fill(dRV0MC);
-      }
-      if(IsAcceptedDCAK0(/*fDCAV0ToVertex*/)){
-       
-       fhV0InvMassK0DCA->Fill(massK0);
-       if(IsK0InvMass(massK0)) fhdNdptK0DCA->Fill(ptV0);
-
-       
-       if(IsAccepteddEdx(momPos,dEdxPos,AliPID::kPion,fCutnSigdEdx)) fhdEdxVsMomV0pidEdx->Fill(momPos,dEdxPos);
-       if(IsAccepteddEdx(momNeg,dEdxNeg,AliPID::kPion,fCutnSigdEdx)) fhdEdxVsMomV0pidEdx->Fill(-1*momNeg,dEdxNeg);
-
-       if(IsAccepteddEdx(momPos,dEdxPos,AliPID::kPion,fCutnSigdEdx) && 
-          IsAccepteddEdx(momNeg,dEdxNeg,AliPID::kPion,fCutnSigdEdx)){ 
-         
-         fhV0InvMassK0DCAdEdx->Fill(massK0);
-         fhV0PtVsInvMassK0DCAdEdx->Fill(massK0,ptV0); 
-
-         if(IsK0InvMass(massK0)){
-           fhdNdptK0DCAdEdx->Fill(ptV0);
-           fhV0DCAToVertexK0->Fill(fDCAV0ToVertex);
-           fhdROpanK0VsPt->Fill(ptV0,dROpanV0);
-           foundK0 = kTRUE;
-         }
-       }
-       
-
-       // AOD pid - cuts strongly into signal
-
-       AliAODTrack::AODTrkPID_t mpPIDNeg = trackNeg->GetMostProbablePID();
-       AliAODTrack::AODTrkPID_t mpPIDPos = trackPos->GetMostProbablePID();
-       
-       if(mpPIDNeg == AliAODTrack::kPion) fhdEdxVsMomV0piPID->Fill(momPos,dEdxPos);
-       if(mpPIDPos == AliAODTrack::kPion) fhdEdxVsMomV0piPID->Fill(-1*momNeg,dEdxNeg);
-
-
-       if( (mpPIDNeg == AliAODTrack::kPion) && (mpPIDPos == AliAODTrack::kPion) ){
-
-         fhV0InvMassK0DCAPID->Fill(massK0);
-       }       
-      }
-      
-
-      if(IsLambdaInvMass(massLambda) || (IsLambdaInvMass(massAntiLambda))){
-       
-       fhV0InvMassK0Lambda->Fill(massK0);      
-      }
-
-
-      // Lambda 
-
-      if(IsAcceptedDCALambda(/*fDCAV0ToVertex*/) && 
-        IsAccepteddEdx(momPos,dEdxPos,AliPID::kProton,fCutnSigdEdx) && 
-        IsAccepteddEdx(momNeg,dEdxNeg,AliPID::kPion,fCutnSigdEdx)){
-    
-       fhV0InvMassLambdaDCAdEdx->Fill(massLambda);
-      }
-      
-      if(IsAcceptedDCALambda(/*fDCAV0ToVertex*/) && 
-        IsAccepteddEdx(momNeg,dEdxNeg,AliPID::kProton,fCutnSigdEdx) && 
-        IsAccepteddEdx(momPos,dEdxPos,AliPID::kPion,fCutnSigdEdx)){
-    
-       fhV0InvMassAntiLambdaDCAdEdx->Fill(massAntiLambda);
-      }
-
-
-      // jet events 
-
-      if(jetVect){
-       Int_t regionV0vect = IsTrackInsideRegion(jetVect,&v0MomVect);
-
-       // calc xi 
-       Double_t jetpt     = jetVect->Pt();
-       Double_t dPhiJetV0 = (jetVect->MomentumVector()->Vect()).DeltaPhi(v0MomVect);
-       Double_t dRJetV0   = (jetVect->MomentumVector()->Vect()).DeltaR(v0MomVect);             
-
-       Double_t z         = pV0/jetpt;
-       Double_t xi        = TMath::Log(1/z);
-
-       fhdPhiJetV0->Fill(dPhiJetV0);
-
-       // K0
-
-       if(IsAcceptedDCAK0(/*fDCAV0ToVertex*/) && 
-          IsAccepteddEdx(momPos,dEdxPos,AliPID::kPion,fCutnSigdEdx) && 
-          IsAccepteddEdx(momNeg,dEdxNeg,AliPID::kPion,fCutnSigdEdx)){
-         
-         fhV0InvMassK0JetEvt->Fill(massK0);
-
-         if(IsK0InvMass(massK0)){
-           fhdNdptK0JetEvt->Fill(ptV0);
-           fhdNdzK0->Fill(z); 
-           fhdNdxiK0->Fill(xi);
-
-           fhdPhiJetK0->Fill(dPhiJetV0);
-           fhdRJetK0->Fill(dRJetV0);
-
-           if(5<jetpt && jetpt<=10)       fhdNdzK05to10->Fill(z);
-           else if(10<jetpt && jetpt<=20) fhdNdzK010to20->Fill(z);
-           else if(20<jetpt && jetpt<=30) fhdNdzK020to30->Fill(z);
-           else if(30<jetpt && jetpt<=40) fhdNdzK030to40->Fill(z);
-           else if(40<jetpt && jetpt<=60) fhdNdzK040to60->Fill(z);
-         }
-       }
-       
-
-       // Lambda
-       
-       if(IsAcceptedDCALambda(/*fDCAV0ToVertex*/) && 
-          IsAccepteddEdx(momPos,dEdxPos,AliPID::kProton,fCutnSigdEdx) && 
-          IsAccepteddEdx(momNeg,dEdxNeg,AliPID::kPion,fCutnSigdEdx)){
-         
-         fhV0InvMassLambdaJetEvt->Fill(massLambda);
-
-         if(IsLambdaInvMass(massLambda)){
-           fhdNdzLambda->Fill(z); 
-         }
-       }
-       
-       if(IsAcceptedDCALambda(/*fDCAV0ToVertex*/) && 
-          IsAccepteddEdx(momNeg,dEdxNeg,AliPID::kProton,fCutnSigdEdx) && 
-          IsAccepteddEdx(momPos,dEdxPos,AliPID::kPion,fCutnSigdEdx)){
-         
-         fhV0InvMassAntiLambdaJetEvt->Fill(massAntiLambda);
-
-         if(IsLambdaInvMass(massAntiLambda)){
-           fhdNdzAntiLambda->Fill(z); 
-         }
-       }
-       
-       // fill histos max region 
-       
-       if(regionV0vect != 0 && regionV0vect == maxPtRegionIndex){ // max region
-       
-         // K0
-
-         if(IsAcceptedDCAK0(/*fDCAV0ToVertex*/) && 
-            IsAccepteddEdx(momPos,dEdxPos,AliPID::kPion,fCutnSigdEdx) && 
-            IsAccepteddEdx(momNeg,dEdxNeg,AliPID::kPion,fCutnSigdEdx)){ // K0 cuts
-           
-           fhV0InvMassK0Max->Fill(massK0);
-           
-           if(IsK0InvMass(massK0)){
-             fhdNdzK0Max->Fill(z);
-             fhdNdxiK0Max->Fill(xi);
-             fhdNdptK0Max->Fill(ptV0);
-           }
-         }
-
-         // Lambda
-         
-         if(IsAcceptedDCALambda(/*fDCAV0ToVertex*/) && 
-            IsAccepteddEdx(momPos,dEdxPos,AliPID::kProton,fCutnSigdEdx) && 
-            IsAccepteddEdx(momNeg,dEdxNeg,AliPID::kPion,fCutnSigdEdx)){
-           
-           fhV0InvMassLambdaMax->Fill(massLambda);
-         }
-
-         if(IsAcceptedDCALambda(/*fDCAV0ToVertex*/) && 
-            IsAccepteddEdx(momNeg,dEdxNeg,AliPID::kProton,fCutnSigdEdx) && 
-            IsAccepteddEdx(momPos,dEdxPos,AliPID::kPion,fCutnSigdEdx)){
-           
-           fhV0InvMassAntiLambdaMax->Fill(massAntiLambda);
-         }
-         
-         if(IsAcceptedDCALambda(/*fDCAV0ToVertex*/) && 
-            ((IsAccepteddEdx(momPos,dEdxPos,AliPID::kProton,fCutnSigdEdx) && 
-              IsAccepteddEdx(momNeg,dEdxNeg,AliPID::kPion,fCutnSigdEdx)) ||
-             (IsAccepteddEdx(momNeg,dEdxNeg,AliPID::kProton,fCutnSigdEdx) && 
-              IsAccepteddEdx(momPos,dEdxPos,AliPID::kPion,fCutnSigdEdx)))){
-              
-           if(IsLambdaInvMass(massLambda) || (IsLambdaInvMass(massAntiLambda))){
-             fhdNdzLambdaMax->Fill(z);
-             fhdNdxiLambdaMax->Fill(xi);
-             fhdNdptLambdaMax->Fill(ptV0);
-           }
-         }
-       }
-       
-       // fill histos min region
-
-       if(regionV0vect != 0 && regionV0vect != maxPtRegionIndex){ // min region 
-       
-         // K0
-  
-         if(IsAcceptedDCAK0(/*fDCAV0ToVertex*/) && 
-            IsAccepteddEdx(momPos,dEdxPos,AliPID::kPion,fCutnSigdEdx) && 
-            IsAccepteddEdx(momNeg,dEdxNeg,AliPID::kPion,fCutnSigdEdx)){
-           
-           fhV0InvMassK0Min->Fill(massK0);
-           
-           if(IsK0InvMass(massK0)){
-             fhdNdzK0Min->Fill(z);
-             fhdNdxiK0Min->Fill(xi);
-             fhdNdptK0Min->Fill(ptV0);
-           }
-         }
-       
-
-         // Lambda
-         
-         
-         if(IsAcceptedDCALambda(/*fDCAV0ToVertex*/) && 
-            IsAccepteddEdx(momPos,dEdxPos,AliPID::kProton,fCutnSigdEdx) && 
-            IsAccepteddEdx(momNeg,dEdxNeg,AliPID::kPion,fCutnSigdEdx)){
-           
-           fhV0InvMassLambdaMin->Fill(massLambda);
-         }
-
-         if(IsAcceptedDCALambda(/*fDCAV0ToVertex*/) && 
-            IsAccepteddEdx(momNeg,dEdxNeg,AliPID::kProton,fCutnSigdEdx) && 
-            IsAccepteddEdx(momPos,dEdxPos,AliPID::kPion,fCutnSigdEdx)){
-           
-           fhV0InvMassAntiLambdaMin->Fill(massAntiLambda);
-         }
-         
-         if(IsAcceptedDCALambda(/*fDCAV0ToVertex*/) && 
-            ((IsAccepteddEdx(momPos,dEdxPos,AliPID::kProton,fCutnSigdEdx) && 
-              IsAccepteddEdx(momNeg,dEdxNeg,AliPID::kPion,fCutnSigdEdx)) ||
-             (IsAccepteddEdx(momNeg,dEdxNeg,AliPID::kProton,fCutnSigdEdx) && 
-              IsAccepteddEdx(momPos,dEdxPos,AliPID::kPion,fCutnSigdEdx)))){
-           
-           if(IsLambdaInvMass(massLambda) || (IsLambdaInvMass(massAntiLambda))){
-             fhdNdzLambdaMin->Fill(z);
-             fhdNdxiLambdaMin->Fill(xi);
-             fhdNdptLambdaMin->Fill(ptV0);
-           }
-         }
-       }
-       
-
-       // jet region 
-
-       if(regionV0vect == 0){ // jet region
-       
-         //Double_t dRJetV0 = jetVect->DeltaR(v0MomVect); 
-         
-         if(dRJetV0 <= fConeRadius){ 
-         
-           // fill histos jet region  
-
-           // K0
-
-           if(IsAcceptedDCAK0(/*fDCAV0ToVertex*/) && 
-              IsAccepteddEdx(momPos,dEdxPos,AliPID::kPion,fCutnSigdEdx) && 
-              IsAccepteddEdx(momNeg,dEdxNeg,AliPID::kPion,fCutnSigdEdx)){ // K0 cuts
-             
-             fhV0InvMassK0Jet->Fill(massK0);
-           
-             if(IsK0InvMass(massK0)){
-               fhdNdzK0Jet->Fill(z);
-               fhdNdxiK0Jet->Fill(xi);
-               fhdNdptK0Jet->Fill(ptV0);
-             }
-           }
-
-
-           // Lambda 
-         
-           if(IsAcceptedDCALambda(/*fDCAV0ToVertex*/) && 
-              IsAccepteddEdx(momPos,dEdxPos,AliPID::kProton,fCutnSigdEdx) && 
-              IsAccepteddEdx(momNeg,dEdxNeg,AliPID::kPion,fCutnSigdEdx)){
-           
-             fhV0InvMassLambdaJet->Fill(massLambda);
-           }     
-           if(IsAcceptedDCALambda(/*fDCAV0ToVertex*/) && 
-              IsAccepteddEdx(momNeg,dEdxNeg,AliPID::kProton,fCutnSigdEdx) && 
-              IsAccepteddEdx(momPos,dEdxPos,AliPID::kPion,fCutnSigdEdx)){
-             
-             fhV0InvMassAntiLambdaJet->Fill(massAntiLambda);
-           }
-           
-           if(IsAcceptedDCALambda(/*fDCAV0ToVertex*/) && 
-              ((IsAccepteddEdx(momPos,dEdxPos,AliPID::kProton,fCutnSigdEdx) && 
-                IsAccepteddEdx(momNeg,dEdxNeg,AliPID::kPion,fCutnSigdEdx)) ||
-               (IsAccepteddEdx(momNeg,dEdxNeg,AliPID::kProton,fCutnSigdEdx) && 
-                IsAccepteddEdx(momPos,dEdxPos,AliPID::kPion,fCutnSigdEdx)))){
-             
-             if(IsLambdaInvMass(massLambda) || (IsLambdaInvMass(massAntiLambda))){
-               fhdNdzLambdaJet->Fill(z);
-               fhdNdxiLambdaJet->Fill(xi);
-               fhdNdptLambdaJet->Fill(ptV0);
-             }
-           }
-         }
-       }
-      }
-    }
-  }
-}
-
-// ----------------------------------------------------------------------------
-
-void AliAnalysisTaskJetChem::CheckMCParticles(AliAODJet* jetVect,Int_t maxPtRegionIndex,Bool_t& isK0Event){
-
-  // histos for generated MC  
-
-  TClonesArray* farray = (TClonesArray*)fAOD->FindListObject("mcparticles");
-       
-  if(!farray){
-    AliInfo("no mcparticles branch"); 
-    return;
-  }
-
-  Int_t ntrks = farray->GetEntries();
-  if (fDebug>1) AliInfo(Form("check MC particles, tracks %d \n",ntrks));
-
-  Int_t pythiaPID = GetPythiaProcessID();
-
-  isK0Event          = kFALSE;
-  Bool_t ispEvent    = kFALSE;
-  Bool_t ispBarEvent = kFALSE;
-  Bool_t isKchEvent  = kFALSE;
-
-  Bool_t isQuarkHardScatteringEvent = IsQuarkHardScatteringEvent(pythiaPID);
-  Bool_t isGluonHardScatteringEvent = IsGluonHardScatteringEvent(pythiaPID);
-
-  for(Int_t i =0 ; i<ntrks; i++){ // mc tracks loop
-
-    AliAODMCParticle* mctrk = (AliAODMCParticle*)farray->At(i);
-    
-    Double_t trackPt = mctrk->Pt(); 
-
-    //Cuts
-    if (!(mctrk->IsPhysicalPrimary())) continue;
-    if ((trackPt < fTrackPtCut) || (TMath::Abs(mctrk->Eta()) > fTrackEtaCut)) continue; 
-
-    // OB PID histo
-    Int_t pdg = mctrk->GetPdgCode();
-    
-
-    FillPIDhisto(fhPIDMC,pdg);
-    if(isQuarkHardScatteringEvent) FillPIDhisto(fhPIDMC_quarkEv,pdg);
-    if(isGluonHardScatteringEvent) FillPIDhisto(fhPIDMC_gluonEv,pdg);
-
-    if(pdg == 111)              fhdNdptpi0MC->Fill(trackPt);
-    if(pdg == 22)               fhdNdptgammaMC->Fill(trackPt);
-    if(TMath::Abs(pdg) == 211)  fhdNdptchPiMC->Fill(trackPt);
-    if(pdg == 310)              fhdNdptK0MC->Fill(trackPt);
-    if(TMath::Abs(pdg) == 321)  fhdNdptchKMC->Fill(trackPt);
-    if(pdg == 2212)             fhdNdptpMC->Fill(trackPt);
-    if(pdg == -2212)            fhdNdptpBarMC->Fill(trackPt);
-    if(pdg == 3122)             fhdNdptLambdaMC->Fill(trackPt);
-    if(pdg == -3122)            fhdNdptLambdaBarMC->Fill(trackPt);
-    if(pdg == 3332)             fhdNdptOmegaMC->Fill(trackPt);
-    if(pdg == -3332)            fhdNdptOmegaBarMC->Fill(trackPt);
-    
-    if(pdg == 310)             isK0Event   = kTRUE;
-    if(TMath::Abs(pdg) == 321) isKchEvent  = kTRUE;
-    if(pdg == 2212)            ispEvent    = kTRUE;
-    if(pdg == -2212)           ispBarEvent = kTRUE;
-    
-    
-    Int_t pdgMotherChK = -1;
-    Int_t pdgMotherK0  = -1;
-    Int_t pdgGrandMotherK0 = -1;
-
-    if(TMath::Abs(pdg) == 321){ // chK
-      Int_t labelMother = mctrk->GetMother();
-      AliAODMCParticle* mctrkMother    = NULL;
-      
-      for(Int_t k=0 ; k<ntrks; k++){
-      
-       AliAODMCParticle* mctrk2 = (AliAODMCParticle*)farray->At(k);
-       if(mctrk2->GetLabel() == labelMother) mctrkMother = mctrk2;
-      }
-      
-      if(mctrkMother) pdgMotherChK = mctrkMother->GetPdgCode();
-      FillPIDhisto(fhPIDMCMotherChK,pdgMotherChK);
-
-      //printf("pdgMotherChK %d \n",pdgMotherChK);
-    }
-    
-    if(pdg == 310){ // K0
-      Int_t labelMother = mctrk->GetMother();
-      AliAODMCParticle* mctrkMother = NULL;
-      
-      for(Int_t k=0 ; k<ntrks; k++){
-       AliAODMCParticle* mctrk2 = (AliAODMCParticle*)farray->At(k);
-       if(mctrk2->GetLabel() == labelMother) mctrkMother = mctrk2;
-      }
-      
-      if(mctrkMother) pdgMotherK0 = mctrkMother->GetPdgCode();
-      FillPIDhisto(fhPIDMCMotherK0,pdgMotherK0);
-
-      Int_t labelGrandMother = -1; 
-      if(mctrkMother) mctrkMother->GetMother();
-      AliAODMCParticle* mctrkGrandMother = NULL;
-      
-      for(Int_t k=0 ; k<ntrks; k++){
-       AliAODMCParticle* mctrk2 = (AliAODMCParticle*)farray->At(k);
-       if(mctrk2->GetLabel() == labelGrandMother) mctrkGrandMother = mctrk2;
-      }
-      
-      if(mctrkGrandMother) pdgGrandMotherK0 = mctrkGrandMother->GetPdgCode();
-      FillPIDhisto(fhPIDMCGrandMotherK0,pdgGrandMotherK0);
-    }
-    
-    if(jetVect){ // jet event 
-
-      FillPIDhisto(fhPIDMCAll,pdg);
-    
-      TVector3 partVect(mctrk->Px(), mctrk->Py(), mctrk->Pz());
-      
-      Int_t region = IsTrackInsideRegion(jetVect, &partVect );  
-      
-      //Double_t deltaPhi = jetVect->DeltaPhi(partVect); //+k270rad;
-      Double_t deltaPhi = (jetVect->MomentumVector()->Vect()).DeltaPhi(partVect);
-      if( deltaPhi > 2.*TMath::Pi() )  deltaPhi-= 2.*TMath::Pi();
-      //Double_t deltaR = jetVect->DeltaR(partVect); //+k270rad;
-      Double_t deltaR = (jetVect->MomentumVector()->Vect()).DeltaR(partVect);
-
-      // calc xi 
-      Double_t jetpt = jetVect->Pt(); 
-      //Double_t pV0   = partVect.Mag(); 
-      //Double_t ptV0  = partVect.Pt();
-   
-
-      Double_t z  = trackPt / jetpt;
-      Double_t xi = TMath::Log(1/z);
-
-      if(!(mctrk->Charge() == 0 || mctrk->Charge()==-99)) fhdNdxiMC->Fill(xi);
-
-      if(pdg == 310){ // K0
-       fhdPhiJetK0MC->Fill(deltaPhi);
-       fhdRJetK0MC->Fill(deltaR);
-       fhdNdxiK0MC->Fill(xi);
-       fhdNdzK0MC->Fill(z);
-       fhdNdptK0MCJetEvt->Fill(trackPt);
-      }
-
-      
-      if(region != 0 && region == maxPtRegionIndex){ // max region
-       
-       if(TMath::Abs(pdg) == 211)  fhdNdptchPiMCMax->Fill(trackPt);
-       if(pdg == 310)              fhdNdptK0MCMax->Fill(trackPt);
-       if(TMath::Abs(pdg) == 321)  fhdNdptchKMCMax->Fill(trackPt);
-       if(pdg == 2212)             fhdNdptpMCMax->Fill(trackPt);
-       if(pdg == -2212)            fhdNdptpBarMCMax->Fill(trackPt);
-       if(pdg == 3122)             fhdNdptLambdaMCMax->Fill(trackPt);
-       if(pdg == -3122)            fhdNdptLambdaBarMCMax->Fill(trackPt);
-      }
-      
-      if(region != 0 && region != maxPtRegionIndex){ // min region
-       
-       FillPIDhisto(fhPIDMCMin,pdg);
-       
-       if(TMath::Abs(pdg) == 211)  fhdNdptchPiMCMin->Fill(trackPt);
-       if(pdg == 310)              fhdNdptK0MCMin->Fill(trackPt);
-       if(TMath::Abs(pdg) == 321)  fhdNdptchKMCMin->Fill(trackPt);
-       if(pdg == 2212)             fhdNdptpMCMin->Fill(trackPt);
-       if(pdg == -2212)            fhdNdptpBarMCMin->Fill(trackPt);
-       if(pdg == 3122)             fhdNdptLambdaMCMin->Fill(trackPt);
-       if(pdg == -3122)            fhdNdptLambdaBarMCMin->Fill(trackPt);
-
-       if(pdg == 3332)  fhdNdptOmegaMCMin->Fill(trackPt);
-       if(pdg == -3332) fhdNdptOmegaBarMCMin->Fill(trackPt);
-      }
-
-      // trans region
-      if(region != 0){
-       if(TMath::Abs(pdg) == 321) FillPIDhisto(fhPIDMCMotherChKTrans,pdgMotherChK);
-       if(pdg == 310){
-         FillPIDhisto(fhPIDMCMotherK0Trans,pdgMotherK0);
-         FillPIDhisto(fhPIDMCGrandMotherK0Trans,pdgGrandMotherK0);
-       }
-      }
-   
-      if(region == 0){ //  jet region ?
-       
-       FillPIDhisto(fhPIDMCJet,pdg);
-       
-       //Double_t dRJetV0 = jetVect->DeltaR(partVect); 
-       Double_t dRJetV0 =  (jetVect->MomentumVector()->Vect()).DeltaR(partVect);
-       
-       if(dRJetV0 <= fConeRadius){ 
-         
-         if(pdg == 310){ // K0
-           
-           fhdNdptK0MCJet->Fill(trackPt);
-           fhdNdxiK0MCJet->Fill(xi);
-           fhdNdzK0MCJet->Fill(z);
-         }
-       
-         if(TMath::Abs(pdg) == 211)  fhdNdptchPiMCJet->Fill(trackPt);
-         if(TMath::Abs(pdg) == 321)  fhdNdptchKMCJet->Fill(trackPt);
-         if(pdg == 2212)             fhdNdptpMCJet->Fill(trackPt);
-         if(pdg == -2212)            fhdNdptpBarMCJet->Fill(trackPt);
-         if(pdg == 3122)             fhdNdptLambdaMCJet->Fill(trackPt);
-         if(pdg == -3122)            fhdNdptLambdaBarMCJet->Fill(trackPt);  
-       }
-      }
-    }
-  }
-
-
-  FillPythiaIDhisto(fhPythiaProcess,pythiaPID);
-  if(isK0Event)   FillPythiaIDhisto(fhPythiaProcessK0,pythiaPID);
-  if(isKchEvent)  FillPythiaIDhisto(fhPythiaProcessKch,pythiaPID);
-  if(ispEvent)    FillPythiaIDhisto(fhPythiaProcessp,pythiaPID);
-  if(ispBarEvent) FillPythiaIDhisto(fhPythiaProcesspbar,pythiaPID);
-
-}
-       
-// ----------------------------------------------------------------------------
-
-Double_t AliAnalysisTaskJetChem::AssociateV0MC(const TVector3* V0Mom,const Int_t pdg){
-
-  // find closest MC gen. particle for V0 vector
-
-  TClonesArray* farray = (TClonesArray*)fAOD->FindListObject("mcparticles");
-       
-  if(!farray){
-    AliInfo("no mcparticles branch"); 
-    return -1;
-  }
-
-  Double_t dRmin = -1;
-
-  Int_t ntrks = farray->GetEntries();
-
-  for(Int_t i =0 ; i<ntrks; i++){
-
-    AliAODMCParticle* mctrk = (AliAODMCParticle*)farray->At(i);
-    
-    //Cuts
-    if (!(mctrk->IsPhysicalPrimary())) continue;
-
-    Int_t pdgtrk = mctrk->GetPdgCode();
-    
-    if(pdgtrk != pdg) continue;
-
-    TVector3 partVect(mctrk->Px(), mctrk->Py(), mctrk->Pz());
-
-    Double_t dR = V0Mom->DeltaR(partVect);
-
-    if(dRmin<0) dRmin = dR; // initialize
-    
-    if(dR < dRmin) dRmin = dR;
-      
-  }
-  
-  return dRmin;
-}
-
-// ----------------------------------------------------------------------------
-
-void AliAnalysisTaskJetChem::CompLeadingJets(AliAODJet* jetLeadingAOD,AliAODJet* jetLeadingMC,const Int_t pythiaPID,
-                                            const Bool_t foundK0AOD,const Bool_t foundK0MC){
-
-  // leading jet properties
-
-  Double_t ptLeadingAOD      = -1;
-  Double_t etaLeadingAOD     = -1;
-  Double_t phiLeadingAOD     = -1;
-  Int_t    nTracksLeadingAOD = -1;
-
-  Double_t ptLeadingMC       = -1; 
-  Double_t etaLeadingMC      = -1;
-  Double_t phiLeadingMC      = -1;
-  Int_t    nTracksLeadingMC  = -1;
-
-  if(jetLeadingAOD){
-    if(jetLeadingAOD->Pt()>fJetPtCut && TMath::Abs(jetLeadingAOD->Eta())<fJetEtaCut){
-
-      ptLeadingAOD      = jetLeadingAOD->Pt();
-      etaLeadingAOD     = jetLeadingAOD->Eta();
-      phiLeadingAOD     = jetLeadingAOD->Phi();
-      nTracksLeadingAOD = jetLeadingAOD->GetRefTracks()->GetEntriesFast();
-
-      Double_t radiusAOD = GetJetRadius(jetLeadingAOD,0.8);
-      fhnTracksJetVsPtAOD->Fill(ptLeadingAOD,nTracksLeadingAOD);
-      fhRadiusJetVsPtAOD->Fill(ptLeadingAOD,radiusAOD);
-      if(IsQuarkHardScatteringEvent(pythiaPID)) fhnTracksJetVsPtAODquarkEv->Fill(ptLeadingAOD,nTracksLeadingAOD);
-      
-      if(foundK0AOD){
-       
-       fhnTracksJetVsPtAODK0->Fill(ptLeadingAOD,nTracksLeadingAOD);
-       if(IsQuarkHardScatteringEvent(pythiaPID)) fhnTracksJetVsPtAODK0quarkEv->Fill(ptLeadingAOD,nTracksLeadingAOD);
-       fhRadiusJetVsPtAODK0->Fill(ptLeadingAOD,radiusAOD);
-      }
-      
-
-      // check if p/Kch in jet 
-
-      Bool_t foundpKch = kFALSE; 
-      Int_t nTracksJet = jetLeadingAOD->GetRefTracks()->GetEntriesFast();
-
-      for(int i=0; i<nTracksJet; i++){
 
-       AliAODTrack* track = (AliAODTrack*) jetLeadingAOD->GetRefTracks()->At(i);
-
-       Double_t mom  = track->P();
-       
-       AliAODPid* aodPid = track->GetDetPid();
-       Double_t   dEdx   = aodPid->GetTPCsignal();
-
-       if(IsAccepteddEdx(mom,dEdx,AliPID::kKaon,fCutnSigdEdx) ||  
-          IsAccepteddEdx(mom,dEdx,AliPID::kProton,fCutnSigdEdx)){
-
-         foundpKch = kTRUE; 
-       }
-      } // track loop
+      Int_t cutnSigdEdx = 2;
+      if(! (IsAccepteddEdx(momPos,dEdxPos,AliPID::kPion,cutnSigdEdx)) && (IsAccepteddEdx(momNeg,dEdxNeg,AliPID::kPion,cutnSigdEdx)) ) continue;
       
-      if(foundpKch){
-       fhnTracksJetVsPtAODpKch->Fill(ptLeadingAOD,nTracksLeadingAOD);
-       fhRadiusJetVsPtAODpKch->Fill(ptLeadingAOD,radiusAOD);
-      }
-    }
-  }
-
-
-  if(jetLeadingMC){
-    if(jetLeadingMC->Pt()>fJetPtCut && TMath::Abs(jetLeadingMC->Eta())<fJetEtaCut){
-
-      ptLeadingMC      = jetLeadingMC->Pt();
-      etaLeadingMC     = jetLeadingMC->Eta();
-      phiLeadingMC     = jetLeadingMC->Phi();
-      nTracksLeadingMC = jetLeadingMC->GetRefTracks()->GetEntriesFast();    
-      
-      Double_t radiusMC  = GetJetRadius(jetLeadingMC,0.8);
-      fhnTracksJetVsPtMC->Fill(ptLeadingMC,nTracksLeadingMC);
-      fhRadiusJetVsPtMC->Fill(ptLeadingMC,radiusMC);
-      if(IsQuarkHardScatteringEvent(pythiaPID)) fhnTracksJetVsPtMCquarkEv->Fill(ptLeadingMC,nTracksLeadingMC);
-      
-      if(foundK0MC){
-       
-       fhnTracksJetVsPtMCK0->Fill(ptLeadingMC,nTracksLeadingMC);
-       if(IsQuarkHardScatteringEvent(pythiaPID)) fhnTracksJetVsPtMCK0quarkEv->Fill(ptLeadingMC,nTracksLeadingMC);
-       fhRadiusJetVsPtMCK0->Fill(ptLeadingMC,radiusMC);
-      }
-    }
-  }
-  
-  if(jetLeadingAOD && jetLeadingMC){
-    
-    //std::cout<<" comp: leading jetPt AOD "<<ptLeadingAOD<<" MC "<<ptLeadingMC<<std::endl;
-    //if(jetLeadingAOD && jetLeadingMC) 
-    //std::cout<<" leading jet eta AOD "<<jetLeadingAOD->Eta()<<" MC "<<jetLeadingMC->Eta()<<std::endl;
-
+    }   
     
-    if(jetLeadingMC->Pt()>fJetPtCut && jetLeadingAOD->Pt()>fJetPtCut && 
-       TMath::Abs(jetLeadingMC->Eta())<fJetEtaCut && TMath::Abs(jetLeadingAOD->Eta())<fJetEtaCut){
+    if(type == kOnFlyPrim || type == kOfflPrim){
       
-      fhLeadingPtAODvsMC->Fill(ptLeadingMC,ptLeadingAOD);
-      fhLeadingEtaAODvsMC->Fill(etaLeadingMC,etaLeadingAOD);
-      fhLeadingPhiAODvsMC->Fill(phiLeadingMC,phiLeadingAOD);
-      fhnTracksLeadingAODvsMC->Fill(nTracksLeadingMC,nTracksLeadingAOD);
-      
-      TLorentzVector *mom4MC  = jetLeadingMC->MomentumVector();
-      TLorentzVector *mom4AOD = jetLeadingAOD->MomentumVector();
-      
-      Double_t dR = mom4MC->DeltaR(*mom4AOD);
-      fhLeadingdRAODMC->Fill(dR);
-      
-      if(dR<0.4) fhLeadingPtAODvsMCdRcut->Fill(ptLeadingMC,ptLeadingAOD);
+      AliAODTrack *trackPos = (AliAODTrack *) (v0->GetSecondaryVtx()->GetDaughter(0)); // slow 
+      AliAODTrack *trackNeg = (AliAODTrack *) (v0->GetSecondaryVtx()->GetDaughter(1)); // slow  
       
-      Double_t dPt = ptLeadingAOD - ptLeadingMC;
-      Double_t dnTracks = nTracksLeadingAOD - nTracksLeadingMC;
+      //std::cout<<"  filer map trackPos "<<trackPos->GetFilterMap()<<" trackNeg "<<trackNeg->GetFilterMap()<<std::endl;
       
-      fhdnTracksVsdPtLeadingAODMC->Fill(dPt,dnTracks);  
+      if((fFilterMaskK0>0) && !(trackPos->TestFilterBit(fFilterMaskK0)))   continue;
+      if((fFilterMaskK0>0) && !(trackNeg->TestFilterBit(fFilterMaskK0)))   continue;
     }
-  }
-  
-}
-
-// ---------------------------------------------------------------------------
-
-// void AliAnalysisTaskJetChem::CheckK0(){
-
-//   TClonesArray* farray = (TClonesArray*)fAOD->FindListObject("mcparticles");
-  
-//   if(!farray){
-//     AliInfo("no mcparticles branch"); 
-//     return;
-//   }
-  
-  
-//   Int_t ntrks = farray->GetEntries();
-    
-//   for(Int_t i=0; i<ntrks; i++){ // trk loop
-    
-//     AliAODMCParticle* mctrk = (AliAODMCParticle*)farray->At(i);
-    
-//     Int_t pdg = mctrk->GetPdgCode();
-    
-//     if(pdg != 310) continue; // K0
-
-//     cout<<" MC K0: physPrim "<<mctrk->IsPhysicalPrimary()<<endl;
-
-//     Int_t labelMother = mctrk->GetMother();
-
-//     cout<<" found K00, label mother "<<labelMother<<endl;
-
-//     AliAODMCParticle* mctrkMother = NULL;
-//     Int_t pdgMother = -1;
-    
-//     for(Int_t k=0 ; k<ntrks; k++){
-      
-//       mctrkMother = (AliAODMCParticle*)farray->At(k);
-//       if(mctrkMother->GetLabel() == labelMother) break;
-//     }
     
-//     pdgMother = mctrkMother->GetPdgCode();
-//     cout<<" K0 mother pdg "<<pdgMother<<" GID "<<fpdgdb->ConvertPdgToGeant3(pdgMother)<<" isPrimary "<<mctrkMother->IsPrimary()<<endl;
-//     //cout<<" mother name "<<mctrkMother->GetName()<<endl;
-    
-//   }
-// }
-
-// ---------------------------------------------------------------------------
-
-// void AliAnalysisTaskJetChem::CheckK0Stack(){
-  
-//   AliMCEventHandler* mcHandler = dynamic_cast<AliMCEventHandler*> (AliAnalysisManager::GetAnalysisManager()->GetMCtruthEventHandler());
-//   if (!mcHandler) {
-//     Printf("ERROR: Could not retrieve MC event handler");
-//     return;
-//   }
-  
-//   AliMCEvent* mcEvent = mcHandler->MCEvent();
-//   if (!mcEvent) {
-//     Printf("ERROR: Could not retrieve MC event");
-//     return;
-//   }
-  
-//    AliStack* mcStack = mcEvent->Stack();//Load Stack
-//    AliGenPythiaEventHeader*  pythiaGenHeader = AliAnalysisHelperJetTasks::GetPythiaEventHeader(mcEvent);
-
-//    Int_t nTracksMC = mcStack->GetNtrack();
-
-//    Bool_t foundK0 = kFALSE;
-
-//    for (Int_t iTracks = 0; iTracks < nTracksMC; iTracks++) {
-//      //Cuts
-//      if(!(mcStack->IsPhysicalPrimary(iTracks))) continue;
-     
-//      TParticle* mctrk = mcStack->Particle(iTracks);
-        
-//      Int_t pdg = mctrk->GetPdgCode();
-
-//     if ((mctrk->Pt() < 2*fTrackPtCut) || (TMath::Abs(mctrk->Eta()) > fTrackEtaCut )) continue; 
-
-//      if(pdg == 310) foundK0 = kTRUE;
-
-//      if(pdg != 310) continue; // K0 short
-
-//      cout<<" check K0 "<<endl;
-
-//      Int_t indexMother = -999;
-//      TParticle* mctrkMother = mctrk;
-//      TParticle* mctrkFirstGenMother = NULL; // last mother which is not the primary proton
-
-//      Int_t nAncestors = 0;
-
-//      while(indexMother != -1){
-
-//        indexMother = mctrkMother->GetFirstMother();
-     
-//        if(indexMother != -1){
-//      mctrkFirstGenMother = mctrkMother;
-//      mctrkMother  =  mcStack->Particle(indexMother);
-//      nAncestors++;
-//        }
-
-//        cout<<" nAncestors "<<nAncestors<<" pdg mother "<<mctrkMother->GetPdgCode()<<" name "<<mctrkMother->GetName()<<endl;
-//      }
-//      cout<<" pdg firstGenMother "<<mctrkFirstGenMother->GetPdgCode()<<" name "<<mctrkFirstGenMother->GetName()<<endl;
-
-
-
-//      cout<<" pythiaGenHeader "<<pythiaGenHeader<<endl;
-//      cout<<" Pythia Process type "<<pythiaGenHeader->ProcessType()<<" ptHard "<<pythiaGenHeader->GetPtHard()<<endl;
-
-
-//      fhPythiaProcess->Fill(pythiaGenHeader->ProcessType());
-//      if(foundK0) fhPythiaProcess_K0->Fill(pythiaGenHeader->ProcessType());
-
-//      //Int_t indexGrandMother      = mctrkMother->GetFirstMother();
-//      //cout<<" indexGrandMother "<<indexGrandMother<<endl;
-     
-//      //if(indexGrandMother>-1){
-//      //  TParticle* mctrkGrandMother = mcStack->Particle(indexGrandMother);
-//      //  cout<<" pdg grandMother "<<mctrkGrandMother->GetPdgCode()<<" name "<<mctrkGrandMother->GetName()<<endl;
-//      // }
-
-//    }
-// }
-
-// ---------------------------------------------------------------------------------
-
-Bool_t  AliAnalysisTaskJetChem::IsQuarkHardScatteringEvent(const Int_t PID){
-
-  // Pythia Manual sec. 8.2.1 : 
-  // if Pythia PID = 92,93,94   event is diffractive 
-  // if Pythia PID = 13, 28, 68 hard scattering products are gg, qg 
-  // if Pythia PID = 11, 12, 53 hard scattering products are qq,  q\bar{q}  
-
-  if(PID == 92 || PID == 93 || PID == 94 || PID==95) return kFALSE;
-  else if(PID == 13 || PID == 28 || PID == 68) return kFALSE;
-  else if(PID == 11 || PID == 12 || PID == 53) return kTRUE;
-  else{
-    AliInfo(Form("unknown Pythia PID %d",PID));
-  }
-
-  return kFALSE;
-}
-
-
-// ---------------------------------------------------------------------------------
-
-Bool_t  AliAnalysisTaskJetChem::IsGluonHardScatteringEvent(const Int_t PID){
-
-  // Pythia Manual sec. 8.2.1 : 
-  // if Pythia PID = 92,93,94   event is diffractive 
-  // if Pythia PID = 95: low pt event (MPI) 
-  // if Pythia PID = 13, 28, 68 hard scattering products are gg, qg 
-  // if Pythia PID = 11, 12, 53 hard scattering products are qq,  q\bar{q}  
-
-
-  if(PID == 92 || PID == 93 || PID == 94 || PID == 95) return kFALSE;
-  else if(PID == 13 || PID == 68) return kTRUE;
-  else if(PID == 28) return kFALSE; // mixed gq final state
-  else if(PID == 11 || PID == 12 || PID == 53) return kFALSE;
-  else{
-    AliInfo(Form("unknown Pythia PID %d",PID));
-  }
-
-  return kFALSE;
-}
-
-// ---------------------------------------------------------------------------------
-
-Bool_t  AliAnalysisTaskJetChem::IsDiffractiveEvent(const Int_t PID){
-
-  // Pythia Manual sec. 8.2.1 : 
-  // if Pythia PID = 92,93,94   event is diffractive 
-  // if Pythia PID = 13, 28, 68 hard scattering products are gg, qg 
-  // if Pythia PID = 11, 12, 53 hard scattering products are qq,  q\bar{q}  
-
-  if(PID == -1) return kFALSE;
-
-  if(PID == 13 || PID == 28 || PID == 68) return kFALSE;
-  else if(PID == 11 || PID == 12 || PID == 53) return kFALSE;
-  else if(PID == 92 || PID == 93 || PID == 94 || PID==95) return kTRUE;
-  else{
-    AliInfo(Form("unknown Pythia PID %d",PID));
+    list->Add(v0);
   }
-
-  return kFALSE;
-
-}
-
-
-// ----------------------------------------------------------------------------------
-
-Int_t AliAnalysisTaskJetChem::GetPythiaProcessID(){
-
-  // Pythia PID for this event 
-  
-  AliMCEventHandler* mcHandler = dynamic_cast<AliMCEventHandler*> (AliAnalysisManager::GetAnalysisManager()->GetMCtruthEventHandler());
-  if (!mcHandler) {
-    //Printf("ERROR: Could not retrieve MC event handler");
-    return -1;
-  }
-  
-  AliMCEvent* mcEvent = mcHandler->MCEvent();
-  if (!mcEvent) {
-    AliInfo("could not retrieve MC event");
-    return -1;
-  }
-  
-  //AliGenPythiaEventHeader*  pythiaGenHeader = AliAnalysisHelperJetTasks::GetPythiaEventHeader(mcEvent); // not OK for HiJing
-  
-  AliGenEventHeader* genHeader = fMCEvent->GenEventHeader();
-  AliGenPythiaEventHeader*  pythiaGenHeader = dynamic_cast<AliGenPythiaEventHeader*>(genHeader);
-
-  if (!pythiaGenHeader) {
-    AliInfo("Could not retrieve pythiaEventHeader");
-    return -1;
-  }
-
-  Int_t  pid = pythiaGenHeader->ProcessType(); 
-
-  return pid;
-
-}
-
-// ----------------------------------------------------------------------------------
-
-void AliAnalysisTaskJetChem::GetJetTracksResum(TList* list, AliAODJet* jet, const Double_t radius){
-  
-
-  if(!jet) return; // no jet in acc in event
-
-  // list of AOD tracks in jet cone, using cone axis and distance axis-track (and not trackrefs)
-
-  Int_t nTracks = fAOD->GetNTracks();
-
-  if(!nTracks) return;
-
-  Double_t jetMom[3];
-  jet->PxPyPz(jetMom);
-  TVector3 jet3mom(jetMom);
-    
-  
-  // standard cuts + ITS refit = stdrd PWG4 cut (compose them for productions with old ESDFilter task) 
-
-  for (Int_t itrack=0; itrack<nTracks; itrack++) {
-
-    AliAODTrack* track = fAOD->GetTrack(itrack);
-    
-    UInt_t status = track->GetStatus();
-
-    if(!track->TestFilterBit(fFilterBitJF)) continue; // track cut selection
-    if(fRequireITSRefitJF &&  ((status&AliESDtrack::kITSrefit)==0)) continue;  
-
-    if(TMath::Abs(track->Eta()) > fTrackEtaCut) continue; 
-    if( track->Pt() < fTrackPtCutJF ) continue;
-
-    Double_t trackMom[3];
-    track->PxPyPz(trackMom);
-    TVector3 track3mom(trackMom);
-
-    Double_t dR = jet3mom.DeltaR(track3mom);
-
-    if(dR<radius){
-
-      list->Add(track);
-      
-    }
-  }
-}
-
-// ----------------------------------------------------------------------------------
-
-void AliAnalysisTaskJetChem::GetJetTracksTrackrefs(TList* list, AliAODJet* jet){
-  
-  if(!jet) return; // no jet in acc in event
-
-  // list of AOD tracks in jet cone, using trackrefs
   
-  Int_t nTracks = jet->GetRefTracks()->GetEntriesFast();
-
-  if(!nTracks) return;
-  
-  // standard cuts + ITS refit = stdrd PWG4 cut (compose them for productions with old ESDFilter task) 
-
-  for (Int_t itrack=0; itrack<nTracks; itrack++) {
-
-    AliAODTrack* track = (AliAODTrack*) jet->GetRefTracks()->At(itrack);
-
-    UInt_t status = track->GetStatus();
-
-    if(!track->TestFilterBit(fFilterBitJF)) continue; // track cut selection
-    if(fRequireITSRefitJF &&  ((status&AliESDtrack::kITSrefit)==0)) continue;  
-
-    if(TMath::Abs(track->Eta()) > fTrackEtaCut) continue; 
-    if( track->Pt() < fTrackPtCutJF ) continue;
-        
-    list->Add(track);
-  }
-
-  //cout<<" trackrefs Size "<<nTracks<<" acc track list size "<<list->GetEntries()<<endl;
-
-}
-
-
-// ----------------------------------------------------------------------------------
-
-void AliAnalysisTaskJetChem::FillReferenceFF(AliAODJet* jet){
-  
-
-  if(!jet) return;
-
-  TList* jetTracks = new TList(); // FIXME - avoid new/delete
-  //GetJetTracksResum(jetTracks,jet,0.7); 
-  GetJetTracksTrackrefs(jetTracks,jet); 
-
-  Double_t jetpt = jet->Pt();
-
-  TIter next(jetTracks);
-  while(AliAODTrack* track = static_cast<AliAODTrack*>(next())){
-    
-    Double_t trackpt = track->Pt();
-    Double_t z       = trackpt/jetpt;
-    Double_t xi      = TMath::Log(1/z);
-    
-    //cout<<" trackpt "<<trackpt<<" jetpt "<<jetpt<<" z "<<z<<" xi "<<xi<<endl;
-
-    if(5<jetpt && jetpt<=10){
-      fhdNdptTracksJetPt5to10->Fill(trackpt);
-      fhdNdzJets5to10->Fill(z);
-      fhdNdxiJets5to10->Fill(xi);
-    }
-    else if(10<jetpt && jetpt<=20){
-      fhdNdptTracksJetPt10to20->Fill(trackpt);
-      fhdNdzJets10to20->Fill(z);
-      fhdNdxiJets10to20->Fill(xi);
-    }
-    else if(20<jetpt && jetpt<=30){
-      fhdNdptTracksJetPt20to30->Fill(trackpt);
-      fhdNdzJets20to30->Fill(z);
-      fhdNdxiJets20to30->Fill(xi);
-    }
-    else if(30<jetpt && jetpt<=40){
-      fhdNdptTracksJetPt30to40->Fill(trackpt);
-      fhdNdzJets30to40->Fill(z);
-      fhdNdxiJets30to40->Fill(xi);
-    }
-    else if(40<jetpt && jetpt<=60){
-      fhdNdptTracksJetPt40to60->Fill(trackpt);
-      fhdNdzJets40to60->Fill(z);
-      fhdNdxiJets40to60->Fill(xi);
-    }
-  }
-
-  
-  delete jetTracks;
+  Int_t nK0s = list->GetSize();
   
+  return nK0s;
 }
 
-// -------------------------------------------------------------
-
-void AliAnalysisTaskJetChem::FillReferencePlotsTracks(){
-
-  // eta/phi & pt tracks before/after cuts
-  // track multiplicity / evt
-
-  Int_t nTracks     = fAOD->GetNTracks();
-  Int_t countTracks = 0;
-
-  // standard cuts + ITS refit = stdrd PWG4 cut (compose them for productions with old ESDFilter task) 
-
-  for(Int_t itrack=0; itrack<nTracks; itrack++){
-    
-    AliAODTrack* track = fAOD->GetTrack(itrack);
-    
-    Double_t trackPt  = track->Pt();
-    Double_t trackPhi = track->Phi();
-    Double_t trackEta = track->Eta();
-
-    fhPhiEtaTracksNoCut->Fill(trackEta,trackPhi);
-    fhPtTracksNoCut->Fill(trackPt);
-    
-    UInt_t status = track->GetStatus();
-    
-    if(!track->TestFilterBit(fFilterBitJF)) continue; // track cut selection
-    if(fRequireITSRefitJF &&  ((status&AliESDtrack::kITSrefit)==0)) continue;  
-
-    fhPhiEtaTracks->Fill(trackEta,trackPhi);
-    fhPtTracks->Fill(trackPt);
-
-    if(TMath::Abs(track->Eta()) > fTrackEtaCut) continue; 
-    if( track->Pt() < fTrackPtCutJF ) continue;
-
-    countTracks++;
-    
-  }
-
-  fhTrackMult->Fill(countTracks);
-}
 
index 0a7ed9d..0394a82 100644 (file)
@@ -1,7 +1,12 @@
-// Analysis task for jet chemistry analysis
-// Oliver Busch, o.busch@gsi.de
-// UE task & CDF jet finder based on UE task by Arian Abrahantes Quintana and Enesto Lopez
-
+/*************************************************************************
+ *                                                                       *
+ * Task for Jet Chemistry Analysis in PWG4 Jet Task Force Train          *
+ *                                                                       *
+ *                                                                       *
+ * contact: Oliver Busch                                                 *
+ * o.busch@gsi.de                                                        *
+ *                                                                       *
+ *************************************************************************/
 
 #ifndef ALIANALYSISTASKJETCHEM_H
 #define ALIANALYSISTASKJETCHEM_H
 /* Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
  * See cxx source for full Copyright notice                               */
 
-#include "AliAnalysisTaskSE.h" 
-#include "AliPID.h"
-
-//class AliAnalysisTaskSE; 
-class AliESDEvent;
-class AliAODEvent;
-class TH1F;
-class TH2F;
-class TH1I;
-class TProfile;
-class TVector3;
-class TTree;
-class AliAODJet;
-class AliAODTrack;
-class AliPID;
-class TDatabasePDG;
-
-class  AliAnalysisTaskJetChem : public AliAnalysisTaskSE
-  {
-  public:
-    AliAnalysisTaskJetChem(const char* name="AliAnalysisTaskJetChem");
-    virtual           ~AliAnalysisTaskJetChem() {;}
-    
-    // Implementation of interface methods
-    virtual     Bool_t UserNotify();
-    virtual     void   UserCreateOutputObjects();
-    virtual     void   UserExec(Option_t *option);
-    virtual     void   Terminate(Option_t *);
-    
-    //  Setters
-    virtual     void   SetDebugLevel( Int_t level )  { fDebug = level; }
-    // Read deltaAODs
-    void   ReadDeltaAOD()                       { fDeltaAOD = kTRUE; }
-    void   SelectDeltaAODBranch(const char* val){ fDeltaAODBranch = val;   }
-    void   SelectAODBranch(const char* val)     { fAODBranch = val;   }
-    void   SelectDeltaAODBranchMC(const char* val){ fDeltaAODBranchMC = val;   }
-    void   SelectAODBranchMC(const char* val)     { fAODBranchMC = val;   }
-
-    void   SetJetsOnFly( Bool_t val )           { fJetsOnFly = val;  }
-
-    // use internal jet finder 
-
-    void   SetUseLOConeJets( )              { fUseLOConeJets = kTRUE; }
-    void   SetUseLOConeMCJets()             { fUseLOConeMCJets = kTRUE;  fUsePythiaJets = kFALSE;}
-    void   SetUsePythiaJets()               { fUseLOConeMCJets = kFALSE; fUsePythiaJets = kTRUE; }
-    void   SetConeRadius( Double_t val )    { fConeRadius = val; }
-    void   SetTrackPtCutJF( Double_t val )  { fTrackPtCutJF = val; }
-
-    void   SetFilterBitJF( UInt_t val )     { fFilterBitJF = val; }
-    void   SetRequireITSRefitJF()           { fRequireITSRefitJF = kTRUE; }
-    void   SetRejectK0TracksJF()            { fRejectK0TracksJF  = kTRUE; }
-
-    // Jet cuts
-    void   SetJetPtCut( Double_t val )    { fJetPtCut = val;  }
-    void   SetJetEtaCut( Double_t val )   { fJetEtaCut = val; }
-
-    // track cuts
-    void   SetFilterBit( UInt_t val )     { fFilterBit = val;  }
-    void   SetTrackPtCut( Double_t val )  { fTrackPtCut = val; }
-    void   SetTrackEtaCut( Double_t val ) { fTrackEtaCut = val; }
-
-    // K0 cuts 
-
-    void SetUseOnFlyV0s()                 { fUseOnFlyV0s = kTRUE; }
-    void SetCutnSigdEdx( Double_t val)    { fCutnSigdEdx = val; }
+/* $Id$ */
 
+#include "AliAnalysisTaskFragmentationFunction.h"
 
-    // Setters for MC
-    void  SetUseAODMCTracksForUE(){ fUseAODMCTracksForUE = kTRUE;}
+class AliAnalysisTaskJetChem : public AliAnalysisTaskFragmentationFunction {
 
+ public:
+  
+  //----------------------------------------
+  class AliFragFuncHistosInvMass : public TObject
+  {
     
-  private:
-    AliAnalysisTaskJetChem(const  AliAnalysisTaskJetChem &det);
-    AliAnalysisTaskJetChem&   operator=(const  AliAnalysisTaskJetChem &det);
-    
-    void   AnalyseEvent();
-    Int_t  IsTrackInsideRegion(const AliAODJet* aodjetVect, const TVector3 *partVect);
-
-    void   FillPIDhisto(TH1F* hist,Int_t pdg,Float_t weight=1);
-    TH1F*  CreatePIDhisto(const char* name);
-    TH1F*  CreatePythiaIDhisto(const char* name);
-    void   FillPythiaIDhisto(TH1F* h, const Int_t PID);
-    void   CreateHistos();
-    void   FillSumPtRegion( Double_t leadingE, Double_t ptMax, Double_t ptMin );
-    void   FillMultRegion(Double_t leadingE, Double_t nTrackPtmax, Double_t nTrackPtmin);
-    TClonesArray* FindChargedParticleJets();
-    TClonesArray* FindChargedParticleJetsMC();
-    TClonesArray* GetPythiaJets();
-    void   QSortTracks(TObjArray &a, Int_t first, Int_t last);
-    void   WriteSettings();
+    public:
     
-    Bool_t IsK0InvMass(const Double_t mass) const;
-    Bool_t IsLambdaInvMass(const Double_t mass) const;
-    Bool_t IsAcceptedDCAK0(/*const Double_t dca*/) const;
-    Bool_t IsAcceptedDCALambda(/*const Double_t dca*/) const;
-    Bool_t IsAccepteddEdx(const Double_t mom, const Double_t dEdx, AliPID::EParticleType n, const Double_t cutnSig) const;
-    void   CheckV0s(AliAODJet* jetVect, Int_t maxPtRegionIndex, Bool_t& foundK0);
-    void   CheckMCParticles(AliAODJet* jetVect, Int_t maxPtRegionIndex,Bool_t& isK0event);
-    Double_t AssociateV0MC(const TVector3* V0Mom,const Int_t pdg);
-    void CompLeadingJets(AliAODJet* jetLeadingAOD,AliAODJet* jetLeadingMC,const Int_t pythiaPID,
-                        const Bool_t foundK0AID,const Bool_t foundK0MC);
-    Double_t GetJetRadius(const AliAODJet* jet, const Double_t energyFrac);
-    Bool_t IsTrackFromK0(const Int_t indexTrackCheck);
-
-    Bool_t IsQuarkHardScatteringEvent(const Int_t PID);
-    Bool_t IsGluonHardScatteringEvent(const Int_t PID);
-    Bool_t IsDiffractiveEvent(const Int_t PID);
-    Int_t  GetPythiaProcessID();
-
-    void GetJetTracksResum(TList* list, AliAODJet* jet, const Double_t radius);
-    void GetJetTracksTrackrefs(TList* list, AliAODJet* jet);
-    void FillReferenceFF(AliAODJet* jet);
-    void FillReferencePlotsTracks(); 
-
-    Int_t   fDebug;                 //  Debug flag
-
-    Bool_t      fDeltaAOD;          //  Read jets from delta AOD 
-    TString     fDeltaAODBranch;    //  Jet branch name from delta AOD
-    TString     fAODBranch;         //  Jet branch name from standard AOD
-    TString     fDeltaAODBranchMC;  //  MC Jet branch name from delta AOD
-    TString     fAODBranchMC;       //  MC Jet branch name from standard AOD
-
-    TClonesArray*  fArrayJetsAOD;   //  Array of AOD Jets from delta AOD
-    TClonesArray*  fArrayJetsMC;    //  Array of MC  Jets from delta AOD
+    AliFragFuncHistosInvMass(const char* name = "FFIMhistos", 
+                            Int_t nJetPt = 0, Float_t jetPtMin = 0, Float_t jetPtMax = 0,
+                            Int_t nInvMass = 0, Float_t invMassMin=0, Float_t invMassMax=0,  
+                            Int_t nPt = 0, Float_t ptMin = 0, Float_t ptMax = 0,
+                            Int_t nXi = 0, Float_t xiMin = 0, Float_t xiMax = 0,
+                            Int_t nZ  = 0, Float_t zMin  = 0, Float_t zMax  = 0);
+    AliFragFuncHistosInvMass(const AliFragFuncHistosInvMass& copy);
+    AliFragFuncHistosInvMass& operator=(const AliFragFuncHistosInvMass &o);
+    virtual ~AliFragFuncHistosInvMass();
     
+    virtual void DefineHistos();
+    virtual void FillFF(Float_t trackPt, Float_t invM, Float_t jetPt,Bool_t incrementJetPt);
+    virtual void AddToOutput(TList* list) const;
 
-    AliAODEvent*  fAOD;             //! AOD Event 
-    AliAODEvent*  fAODjets;         //! AOD Event for reconstructed on the fly (see ConnectInputData()
-    TList*  fListOfHistos;          //  Output list of histograms
-    Bool_t   fJetsOnFly;            // if jets are reconstructed on the fly from AOD tracks 
+  private:
 
-    Bool_t   fUseLOConeJets;        // Use LO cone finder instead of jets from AOD
-    Bool_t   fUseLOConeMCJets;      // Use LO cone finder on aodmcparticles
-    Bool_t   fUsePythiaJets;        // use pythia internal jet finder output
+    Int_t   fNBinsJetPt;    // FF histos bins
+    Float_t fJetPtMin;      // FF histos limits
+    Float_t fJetPtMax;      // FF histos limits
+    Int_t   fNBinsInvMass;  // FF histos bins
+    Float_t fInvMassMin;    // FF histos limits
+    Float_t fInvMassMax;    // FF histos limits
+    Int_t   fNBinsPt;       // FF histos bins
+    Float_t fPtMin;         // FF histos limits
+    Float_t fPtMax;         // FF histos limits
+    Int_t   fNBinsXi;       // FF histos bins
+    Float_t fXiMin;         // FF histos limits
+    Float_t fXiMax;         // FF histos limits
+    Int_t   fNBinsZ;        // FF histos bins
+    Float_t fZMin;          // FF histos limits
+    Float_t fZMax;          // FF histos limits
+  
+    TH3F*   fh3TrackPt;     //! FF: track transverse momentum 
+    TH3F*   fh3Xi;          //! FF: xi 
+    TH3F*   fh3Z;           //! FF: z  
+    TH1F*   fh1JetPt;       //! jet pt 
 
-    Double_t fConeRadius;          // if selected Cone-like region type, set Radius (=0.7 default)
-    Double_t fTrackPtCutJF;        // track lower pt for JF     
-    UInt_t   fFilterBitJF;         // track filter for JetFinder 
-    Bool_t   fRequireITSRefitJF;   // additional ITS refit requirement in JF 
-    Bool_t   fRejectK0TracksJF;    // exclude tracks from K0 decay in JF 
+    TString fNameFF;        // histo names prefix
+    
+    ClassDef(AliFragFuncHistosInvMass, 1);
+  };
+  
+
+ //----------------------------------------
+  class AliFragFuncHistosPhiCorrInvMass : public TObject
+  {
+                                  
+    public:
+    
+    AliFragFuncHistosPhiCorrInvMass(const char* name = "FFPhiCorrIMhistos", 
+                                   Int_t nPt = 0, Float_t ptMin = 0, Float_t ptMax = 0,
+                                   Int_t nPhi = 0, Float_t phiMin = 0, Float_t phiMax = 0,
+                                   Int_t nInvMass = 0, Float_t invMassMin=0, Float_t invMassMax=0);
 
+    AliFragFuncHistosPhiCorrInvMass(const AliFragFuncHistosPhiCorrInvMass& copy);
+    AliFragFuncHistosPhiCorrInvMass& operator=(const AliFragFuncHistosPhiCorrInvMass &o);
+    virtual ~AliFragFuncHistosPhiCorrInvMass();
+    
+    virtual void DefineHistos();
+    virtual void FillPhiCorr(Float_t pt, Float_t phi, Float_t invM);
+    virtual void AddToOutput(TList* list) const;
 
-    // Jet cuts 
-    Double_t   fJetPtCut;          // Min Pt for charged Particle Jet
-    Double_t   fJetEtaCut;         // |jet1 eta| < fJet1EtaCut   (fAnaType = 1,2,3)
+  private:
 
-    // track cuts
-    UInt_t   fFilterBit;           // Select tracks from an specific track cut (default 0xFF all track selected)
-    Double_t fTrackPtCut;          // Pt cut of tracks in the regions
-    Double_t fTrackEtaCut;         // Eta cut on tracks in the regions (fRegionType=1)
+    Int_t   fNBinsPt;       // FF histos bins
+    Float_t fPtMin;         // FF histos limits
+    Float_t fPtMax;         // FF histos limits
 
-    // K0 cuts 
-    Bool_t fUseOnFlyV0s;           // on-the-fly V0s versus vertex track 'offline' V0s
-    Double_t fCutnSigdEdx;         // TPC dEdx cut 
+    Int_t   fNBinsPhi;      // FF histos bins
+    Float_t fPhiMin;        // FF histos limits
+    Float_t fPhiMax;        // FF histos limits
+    
+    Int_t   fNBinsInvMass;  // FF histos bins
+    Float_t fInvMassMin;    // FF histos limits
+    Float_t fInvMassMax;    // FF histos limits
+  
+    TH3F*   fh3PhiCorr;     //! FF: phi correlation histo 
 
-    // For MC
-    Bool_t fUseAODMCTracksForUE;   // use aodmcparticles branch to determine max/min transverse region
-    Double_t fAreaReg;             // Area of the region To be used as normalization factor when filling histograms
-    Double_t fAvgTrials;           // average trials used to fill the fh1Triasl histogram in case we do not have trials on a event by event basis
+    TString fNamePhiCorr;   // histo names prefix
     
+    ClassDef(AliFragFuncHistosPhiCorrInvMass, 1);
+  };
+  
+  //----------------------------------------
 
-    // Histograms    ( are owned by fListOfHistos TList )
-    TH1F*  fhPrimVertexNCont;        //!
-    TH1F*  fhPrimVertexRho;          //!
-    TH1F*  fhPrimVertexZ;            //!
-    TH1F*  fhNJets;                  //!
-    TH1F*  fhNJetsMC;                //!
-    TH1F*  fhLeadingEta;             //!
-    TH2F*  fhLeadingNTracksVsEta;    //!
-    TH2F*  fhLeadingPtVsEta;         //!
-    TH1F*  fhLeadingPhi;             //!  
-    TH1F*  fhLeadingPt;              //!
-    TH1F*  fhLeadingPtDiffr;         //!
-    TH1F*  fhLeadingEtaMC;           //!
-    TH1F*  fhLeadingPhiMC;           //!  
-    TH1F*  fhLeadingPtMC;            //!
-    TH1F*  fhLeadingPtMCDiffr;       //!
-    TH2F*  fhPhiEtaTracksNoCut;      //!
-    TH1F*  fhPtTracksNoCut;          //!
-    TH2F*  fhPhiEtaTracks;           //!
-    TH1F*  fhPtTracks;               //!
-    TH1F*  fhTrackMult;              //!
-
-    TH1F*  fhEtaMCTracks;            //!
-    TH1F*  fhPhiMCTracks;            //!
-    TH1F*  fhPtMCTracks;             //!
-    TH2F*  fhnTracksVsPtLeading;     //!
+  AliAnalysisTaskJetChem(); 
+  AliAnalysisTaskJetChem(const char *name);
+  AliAnalysisTaskJetChem(const  AliAnalysisTaskJetChem &copy);
+  AliAnalysisTaskJetChem& operator=(const  AliAnalysisTaskJetChem &o);
+  virtual ~AliAnalysisTaskJetChem();
+  
+  virtual void   UserCreateOutputObjects();
+  virtual void   UserExec(Option_t *option);
+
+  static  void   SetProperties(TH3F* h,const char* x, const char* y,const char* z);
+
+  Bool_t IsAccepteddEdx(const Double_t mom,const Double_t signal, AliPID::EParticleType n, const Double_t cutnSig) const;
+  Bool_t IsK0InvMass(const Double_t mass) const; 
+  Int_t  GetListOfK0s(TList *list, const Int_t type);
+  virtual void SetK0Type(Int_t i){ fK0Type = i; }
+  virtual void SetFilterMaskK0(UInt_t i) {fFilterMaskK0 = i;}
+
+
+  void   SetFFInvMassHistoBins(Int_t nJetPt = 19, Float_t jetPtMin = 5, Float_t jetPtMax = 100, 
+                              Int_t nInvM = 50, Float_t invMMin = 0.450,  Float_t invMMax = 0.550,
+                              Int_t nPt = 20, Float_t ptMin = 0., Float_t ptMax = 20., 
+                              Int_t nXi = 35, Float_t xiMin = 0., Float_t xiMax = 7.,
+                              Int_t nZ = 11,  Float_t zMin = 0.,  Float_t zMax = 1.1)
+  { fFFIMNBinsJetPt = nJetPt; fFFIMJetPtMin = jetPtMin; fFFIMJetPtMax = jetPtMax; 
+    fFFIMNBinsInvM = nInvM; fFFIMInvMMin = invMMin; fFFIMInvMMax = invMMax; fFFIMNBinsPt = nPt; fFFIMPtMin = ptMin; fFFIMPtMax = ptMax; 
+    fFFIMNBinsXi = nXi; fFFIMXiMin = xiMin; fFFIMXiMax = xiMax; fFFIMNBinsZ  = nZ;  fFFIMZMin  = zMin;  fFFIMZMax  = zMax; }
+
+  void   SetPhiCorrInvMassHistoBins(Int_t nPt = 40, Float_t ptMin = 0., Float_t ptMax = 20., 
+                                   Int_t nPhi = 20, Float_t phiMin = 0., Float_t phiMax = 2*TMath::Pi(),
+                                   Int_t nInvM = 50, Float_t invMMin = 0.450,  Float_t invMMax = 0.550)
+                                   
+  { fPhiCorrIMNBinsPt = nPt; fPhiCorrIMPtMin = ptMin; fPhiCorrIMPtMax = ptMax;
+    fPhiCorrIMNBinsPhi = nPhi; fPhiCorrIMPhiMin = phiMin; fPhiCorrIMPhiMax = phiMax;
+    fPhiCorrIMNBinsInvM = nInvM; fPhiCorrIMInvMMin = invMMin; fPhiCorrIMInvMMax = invMMax;
+  }
+  
+  
+  // consts
+  enum { kTrackUndef =0, kOnFly, kOnFlyPID, kOnFlydEdx, kOnFlyPrim, kOffl, kOfflPID, kOffldEdx, kOfflPrim };  
+  
+ private:
+  
+  Int_t fK0Type;                                           //! K0 cuts
+  UInt_t fFilterMaskK0;                                    //! K0 legs cuts
+  TList* fListK0s;                                         //! K0 list 
+
+  AliFragFuncQATrackHistos*  fV0QAK0;                      //! track QA: V0s in K0 inv mass range
+  AliFragFuncHistos*         fFFHistosRecCutsK0Evt;        //! inclusive FF for K0 evt
+  AliFragFuncHistosInvMass*  fFFHistosIMK0AllEvt;          //! K0 pt spec for all events
+  AliFragFuncHistosInvMass*  fFFHistosIMK0Jet;             //! K0 FF all dPhi   
+  AliFragFuncHistosInvMass*  fFFHistosIMK0Cone;            //! K0 FF jet cone   
+  AliFragFuncHistosPhiCorrInvMass*  fFFHistosPhiCorrIMK0;  //! K0 corelation to jet axis 
+
+  // histogram bins  
+
+  Int_t   fFFIMNBinsJetPt;    // FF histos bins
+  Float_t fFFIMJetPtMin;      // FF histos limits
+  Float_t fFFIMJetPtMax;      // FF histos limits
+
+  Int_t   fFFIMNBinsInvM;     // FF histos bins
+  Float_t fFFIMInvMMin;       // FF histos bins
+  Float_t fFFIMInvMMax;       // FF histos bins
+
+  Int_t   fFFIMNBinsPt;       // FF histos bins
+  Float_t fFFIMPtMin;         // FF histos limits
+  Float_t fFFIMPtMax;         // FF histos limits
+
+  Int_t   fFFIMNBinsXi;       // FF histos bins
+  Float_t fFFIMXiMin;         // FF histos limits
+  Float_t fFFIMXiMax;         // FF histos limits
+
+  Int_t   fFFIMNBinsZ;        // FF histos bins
+  Float_t fFFIMZMin;          // FF histos limits
+  Float_t fFFIMZMax;          // FF histos limits
+
+
+  Int_t fPhiCorrIMNBinsPt;    // FF histos bins
+  Float_t fPhiCorrIMPtMin;    // FF histos limits
+  Float_t fPhiCorrIMPtMax;    // FF histos limits
+
+  Int_t fPhiCorrIMNBinsPhi;   // FF histos bins
+  Float_t fPhiCorrIMPhiMin;   // FF histos limits
+  Float_t fPhiCorrIMPhiMax;   // FF histos limits
+               
+  Int_t fPhiCorrIMNBinsInvM;  // FF histos bins
+  Float_t fPhiCorrIMInvMMin;  // FF histos limits
+  Float_t fPhiCorrIMInvMMax;  // FF histos limits
   
-    TH1F*  fhdNdEtaPhiDist;          //!
-   
-    TH1F*  fhRegionSumPtMaxVsEt;     //!
-    TH1F*  fhRegionMultMaxVsEt;      //!
-    TH1F*  fhRegionSumPtMinVsEt;     //!
-    TH1F*  fhRegionMultMinVsEt;      //!
-
-
-    TH1F* fhNV0s;                     //!
-    TH1F* fhV0onFly;                  //!
-    TH1F* fhV0DCADaughters;           //!
-    TH1F* fhV0Radius;                 //!
-    TH1F* fhV0DCAToVertex;            //!
-    TH1F* fhV0DCAToVertexK0;          //!
-    
-    TH1F* fhV0InvMassK0;               //!
-    TH2F* fhV0PtVsInvMassK0;           //!
-    TH1F* fhV0InvMassK0JetEvt;         //!       
-    TH1F* fhV0InvMassLambda;           //!
-    TH1F* fhV0InvMassAntiLambda;       //!
-    TH1F* fhV0InvMassLambdaJetEvt;     //!
-    TH1F* fhV0InvMassAntiLambdaJetEvt; //!
-
-    TH2F* fhdROpanK0VsPt;             //!
-    TH1F* fhdPhiJetV0;                //!
-    TH1F* fhdPhiJetK0;                //!
-    TH1F* fhdRJetK0;                  //!
-
-    TH1F* fhdNdptV0;                  //!
-    TH1F* fhdNdptK0;                  //!
-    TH2F* fhPtVsEtaK0;                //!
-
-    TH1F* fhV0InvMassK0DCA;              //!
-    TH1F* fhV0InvMassK0DCAdEdx;          //!
-    TH2F* fhV0PtVsInvMassK0DCAdEdx;      //!
-
-    TH1F* fhV0InvMassK0DCAPID;           //!
-    TH1F* fhV0InvMassLambdaDCAdEdx;      //!
-    TH1F* fhV0InvMassAntiLambdaDCAdEdx;  //!
-    TH1F* fhdNdptK0DCA;                  //!
-    TH1F* fhdNdptK0DCAdEdx;              //!
-
-    TH1F* fhV0InvMassK0Min;          //!
-    TH1F* fhV0InvMassLambdaMin;      //!
-    TH1F* fhV0InvMassAntiLambdaMin;  //!
-
-    TH1F* fhV0InvMassK0Max;          //!
-    TH1F* fhV0InvMassLambdaMax;      //!
-    TH1F* fhV0InvMassAntiLambdaMax;  //!
-
-    TH1F* fhV0InvMassK0Jet;          //!
-    TH1F* fhV0InvMassLambdaJet;      //!
-    TH1F* fhV0InvMassAntiLambdaJet;  //!
-
-    TH1F* fhV0InvMassK0Lambda;       //!
-
-    TH1F* fhdNdptK0JetEvt;           //!
-
-    TH1F* fhdNdzK0;                   //!
-    TH1F* fhdNdzK05to10;              //!
-    TH1F* fhdNdzK010to20;             //!
-    TH1F* fhdNdzK020to30;             //!
-    TH1F* fhdNdzK030to40;             //!
-    TH1F* fhdNdzK040to60;             //!
-
-    TH1F* fhdNdxiK0;                 //!
-
-    TH1F* fhdNdzLambda;              //!
-    TH1F* fhdNdzAntiLambda;          //!
-
-    TH1F* fhdNdzK0Max;               //!
-    TH1F* fhdNdxiK0Max;              //!
-    TH1F* fhdNdzLambdaMax;           //!
-    TH1F* fhdNdxiLambdaMax;          //!
-    
-    TH1F* fhdNdptK0Max;              //!
-    TH1F* fhdNdptLambdaMax;          //!
-    
-    TH1F* fhdNdzK0Min;               //!
-    TH1F* fhdNdxiK0Min;              //!
-    TH1F* fhdNdzLambdaMin;           //!
-    TH1F* fhdNdxiLambdaMin;          //!
-    
-    TH1F* fhdNdptK0Min;              //!
-    TH1F* fhdNdptLambdaMin;          //!
 
-    TH1F* fhdNdzK0Jet;               //!
-    TH1F* fhdNdxiK0Jet;              //!
-    TH1F* fhdNdzLambdaJet;           //!
-    TH1F* fhdNdxiLambdaJet;          //!
 
-    TH1F* fhdNdptK0Jet;              //!
-    TH1F* fhdNdptLambdaJet;          //!
-    
-    TH2F* fhdEdxVsMomV0;             //!
-    TH2F* fhdEdxVsMomV0pidEdx;       //! 
-    TH2F* fhdEdxVsMomV0piPID;        //! 
-
-    TH1F* fhdPhiJetK0MC;             //!
-    TH1F* fhdRJetK0MC;               //!
-    TH1F* fhdRV0MC;                  //!
-
-    TH1F* fhdNdptchPiMCMax;               //!
-    TH1F* fhdNdptK0MCMax;                 //!
-    TH1F* fhdNdptchKMCMax;                //!
-    TH1F* fhdNdptpMCMax;                  //!
-    TH1F* fhdNdptpBarMCMax;               //!
-    TH1F* fhdNdptLambdaMCMax;             //!
-    TH1F* fhdNdptLambdaBarMCMax;          //!
-
-    TH1F* fhdNdptchPiMCMin;               //!
-    TH1F* fhdNdptK0MCMin;                 //!
-    TH1F* fhdNdptchKMCMin;                //!
-    TH1F* fhdNdptpMCMin;                  //!
-    TH1F* fhdNdptpBarMCMin;               //!
-    TH1F* fhdNdptLambdaMCMin;             //!
-    TH1F* fhdNdptLambdaBarMCMin;          //!
-    TH1F* fhdNdptOmegaMCMin;              //!
-    TH1F* fhdNdptOmegaBarMCMin;           //!
-
-    TH1F* fhdNdptchPiMCJet;               //!
-    TH1F* fhdNdptK0MCJet;                 //!
-    TH1F* fhdNdptchKMCJet;                //!
-    TH1F* fhdNdptpMCJet;                  //!
-    TH1F* fhdNdptpBarMCJet;               //!
-    TH1F* fhdNdptLambdaMCJet;             //!
-    TH1F* fhdNdptLambdaBarMCJet;          //!
-
-
-    // kine tree 
-    TH1F* fhPIDMC;                    //!
-    TH1F* fhPIDMC_quarkEv;            //!
-    TH1F* fhPIDMC_gluonEv;            //!
-    TH1F* fhPIDMCAll;                 //!
-    TH1F* fhPIDMCMin;                 //!
-    TH1F* fhPIDMCJet;                 //!
-    TH1F* fhPIDMCMotherK0;            //!
-    TH1F* fhPIDMCGrandMotherK0;       //!
-    TH1F* fhPIDMCMotherChK;           //!
-    TH1F* fhPIDMCMotherK0Trans;       //!
-    TH1F* fhPIDMCGrandMotherK0Trans;  //!
-    TH1F* fhPIDMCMotherChKTrans;      //!
-
-    TH1F* fhdNdptgammaMC;             //!
-    TH1F* fhdNdptchPiMC;              //!
-    TH1F* fhdNdptpi0MC;               //!
-    TH1F* fhdNdptK0MC;                //!
-    TH1F* fhdNdptchKMC;               //!
-    TH1F* fhdNdptpMC;                 //!
-    TH1F* fhdNdptpBarMC;              //!
-    TH1F* fhdNdptLambdaMC;            //!
-    TH1F* fhdNdptLambdaBarMC;         //!
-    TH1F* fhdNdptOmegaMC;             //!
-    TH1F* fhdNdptOmegaBarMC;          //!
-
-    TH1F* fhdNdxiMC;                 //!
-    TH1F* fhdNdxiK0MC;               //!
-    TH1F* fhdNdxiK0MCJet;            //!
-
-    TH1F* fhdNdzK0MC;                //!
-    TH1F* fhdNdzK0MCJet;             //!
-    TH1F* fhdNdptK0MCJetEvt;         //!
-      
-    TH2F* fhnJetsAODvsMC;             //!
-    TH2F* fhLeadingPtAODvsMC;         //!
-    TH2F* fhLeadingEtaAODvsMC;        //!
-    TH2F* fhLeadingPhiAODvsMC;        //!
-    TH2F* fhnTracksLeadingAODvsMC;    //!
-    
-    TH1F* fhLeadingdRAODMC;            //!
-    TH2F* fhLeadingPtAODvsMCdRcut;     //!
-    TH2F* fhdnTracksVsdPtLeadingAODMC; //!
-
-    TH2F* fhnTracksJetVsPtAOD;         //!
-    TH2F* fhnTracksJetVsPtAODquarkEv;  //!
-    TH2F* fhRadiusJetVsPtAOD;          //!
-    TH2F* fhnTracksJetVsPtMC;          //!
-    TH2F* fhnTracksJetVsPtMCquarkEv;   //!
-    TH2F* fhRadiusJetVsPtMC;           //!
-
-    TH2F* fhnTracksJetVsPtMCK0;         //!
-    TH2F* fhnTracksJetVsPtMCK0quarkEv;  //!
-    TH2F* fhRadiusJetVsPtMCK0;          //!
-
-    TH2F* fhnTracksJetVsPtAODK0;         //!
-    TH2F* fhnTracksJetVsPtAODK0quarkEv;  //!
-    TH2F* fhRadiusJetVsPtAODK0;          //!
-    TH2F* fhnTracksJetVsPtAODpKch;       //!
-    TH2F* fhRadiusJetVsPtAODpKch;        //!
-
-    TH1F* fhPythiaProcess;           //!
-    TH1F* fhPythiaProcessK0;         //!
-    TH1F* fhPythiaProcessKch;        //!
-    TH1F* fhPythiaProcessp;          //!
-    TH1F* fhPythiaProcesspbar;       //!
-
-    TH1F* fhdNdzJets5to10;  //!
-    TH1F* fhdNdzJets10to20; //!
-    TH1F* fhdNdzJets20to30; //!
-    TH1F* fhdNdzJets30to40; //!
-    TH1F* fhdNdzJets40to60; //!
-
-    TH1F* fhdNdxiJets5to10;  //!
-    TH1F* fhdNdxiJets10to20; //!
-    TH1F* fhdNdxiJets20to30; //!
-    TH1F* fhdNdxiJets30to40; //!
-    TH1F* fhdNdxiJets40to60; //!
-
-    TH1F* fhdNdptTracksJetPt5to10;  //!
-    TH1F* fhdNdptTracksJetPt10to20; //!
-    TH1F* fhdNdptTracksJetPt20to30; //!
-    TH1F* fhdNdptTracksJetPt30to40; //!
-    TH1F* fhdNdptTracksJetPt40to60; //!
-
-
-    TProfile* fh1Xsec;               //!
-    TH1F*  fh1Trials;                //!
-    
-    //TTree* fSettingsTree;            //! Fast Settings saving
-    
-    TDatabasePDG* fpdgdb;            //!
+  // Histograms
 
-    enum PythiaPIDHistoBin{kPythiaPIDP11Bin=1, kPythiaPIDP12Bin=3, kPythiaPIDP13Bin=5, kPythiaPIDP28Bin=7, 
-                          kPythiaPIDP53Bin=9, kPythiaPIDP68Bin=11, kPythiaPIDP92Bin=13, kPythiaPIDP93Bin=15, 
-                          kPythiaPIDP94Bin=17,kPythiaPIDP95Bin=19, kPythiaPIDPOtherBin=21}; 
+  TH1F* fh1K0Mult;                  //!
+  TH1F* fh1dPhiJetK0;               //!
 
 
-    enum PIDHistoBin{kPDGpm311Bin=48,kPDG333Bin=49,kPDGpm313Bin=50,kPDGp323Bin=51,kPDGm323Bin=52,
-                    kPDGNeutrinoBin=53,kPDGCharmedBaryonBin=54,kPDGQuarkBin=55,kPDGDiQuarkBin=56};
-   
-    ClassDef( AliAnalysisTaskJetChem, 2); // Analysis task for jet chemistry analysis 
-  };
+  ClassDef(AliAnalysisTaskJetChem, 3);
+};
 
 #endif
-
-    
index ee38c37..594407e 100644 (file)
@@ -21,7 +21,6 @@
 #pragma link C++ class AliPWG4HighPtTrackQA+;
 #pragma link C++ class AliAnalysisTaskPWG4PidDetEx+;
 #pragma link C++ class AliJetSpectrumUnfolding+;
-#pragma link C++ class AliAnalysisTaskJetChem+;
 #pragma link C++ class AliAnalysisTaskFragmentationFunction+;
 #pragma link C++ class AliAnalysisTaskFragmentationFunction::AliFragFuncHistos+;
 #pragma link C++ class AliAnalysisTaskFragmentationFunction::AliFragFuncDiJetHistos+;
@@ -29,6 +28,9 @@
 #pragma link C++ class AliAnalysisTaskFragmentationFunction::AliFragFuncIntraJetHistos+;
 #pragma link C++ class AliAnalysisTaskFragmentationFunction::AliFragFuncQATrackHistos+;
 #pragma link C++ class AliAnalysisTaskFragmentationFunction::AliFragFuncQAJetHistos+;
+#pragma link C++ class AliAnalysisTaskJetChem+;
+#pragma link C++ class AliAnalysisTaskJetChem::AliFragFuncHistosInvMass+;
+#pragma link C++ class AliAnalysisTaskJetChem::AliFragFuncHistosPhiCorrInvMass+;
 #pragma link C++ class AliAnalysisTaskMinijet+;
 #pragma link C++ class AliUEHist+;
 #pragma link C++ class AliUEHistograms+;
index d12ac39..b9772df 100644 (file)
@@ -1,92 +1,92 @@
-AliAnalysisTask *AddTaskJetChem(UInt_t filterMask){
-
-
+AliAnalysisTaskJetChem *AddTaskJetChem(){
+        
+  AliAnalysisTaskJetChem *ff=0;
+   
+  ff = AddTask("clustersAOD_ANTIKT04_B1_Filter00256_Cut00150_Skip00",0,AliAnalysisTaskJetChem::kOnFly);
+    
+  return ff;
+}
 
-  //cout<<" OB : add JetTasks inlcude path ! "<<endl;
-  //gSystem->AddIncludePath("-I$ALICE_ROOT/PWG4/JetTasks");
+// _______________________________________________________________________________________
 
-  //get the current analysis manager
+AliAnalysisTaskJetChem *AddTask(const char* recJetsBranch,Int_t eventClass, Int_t K0type)
+{
+  // Creates a jet chem task,
+  // configures it and adds it to the analysis manager.
+  
+   
+  // Get the pointer to the existing analysis manager via the static access method.
+  //==============================================================================
   AliAnalysisManager *mgr = AliAnalysisManager::GetAnalysisManager();
   if (!mgr) {
-    Error("AddTask_obusch_jets", "No analysis manager found.");
-    return 0;
+    ::Error("AddTaskJetChem", "No analysis manager to connect to.");
+    return NULL;
   }
-
-  // physics event selection task - not for AOD analysis
-
-  //gROOT->LoadMacro("$ALICE_ROOT/ANALYSIS/macros/AddTaskPhysicsSelection.C");
-  //AliPhysicsSelectionTask* physSelTask = AddTaskPhysicsSelection(); 
-  //AliPhysicsSelection* physSel = physSelTask->GetPhysicsSelection();
-  //physSel->AddBackgroundIdentification(new AliBackgroundSelection()); 
-  //physSel->SetAnalyzeMC();
-
-  /*
-  //============= Set Task Name ===================
-  TString taskName=("AliAnalysisTaskJetChem.cxx+");
-  //===============================================
-  //            Load the task
-  gROOT->LoadMacro(taskName.Data());
-  if (gProof){
-    TString taskSO=gSystem->pwd();
-    taskSO+="/";
-    taskSO+=taskName(0,taskName.First('.'))+"_cxx.so";
-    gProof->Exec(Form("gSystem->Load(\"%s\")",taskSO.Data()),kTRUE);
+  
+  // Check the analysis type using the event handlers connected to the analysis manager.
+  //==============================================================================
+  if (!mgr->GetInputEventHandler()) {
+    ::Error("AddTaskJetChem", "This task requires an input event handler");
+    return NULL;
   }
-  */
-  //========= Add task to the ANALYSIS manager =====
-  AliAnalysisTaskJetChem *task = new AliAnalysisTaskJetChem;
   
-  // configure task 
-
-
-  //task->SetUseLOConeJets();
-  //task->SetUseLOConeMCJets();
+  TString type = mgr->GetInputEventHandler()->GetDataType(); // can be "ESD" or "AOD"
+  Printf("Data Type: %s", type.Data());
   
-  //task->SetUsePythiaJets();
-  task->SetConeRadius(0.4);
-  task->SetTrackPtCutJF(0.150); //
-  task->SetFilterBitJF(0x01);   // official PWG4 high pt filter bit 0x10, but not all AliAnalysisTaskESDFilter had configured this ESDTrackCut
-  task->SetRequireITSRefitJF(); // 0x01 + ITS refit = 0x10 
-  task->SetRejectK0TracksJF();  // uncomment for K0 analysis running jet finder in task - modifies jet spectrum
-
-  task->SetJetPtCut(2.0);
-  //task->SetJetPtCut(0.150); // lower pt cut: for plot of diffractive contribution to jet spectrum (goes up to 2 GeV ...) 
-  task->SetJetEtaCut(0.5);
-
-  task->SetFilterBit((UInt_t) 0X01); // std AOD track cuts
-  task->SetTrackPtCut(0.150); 
-  task->SetTrackEtaCut(0.9);
-
-  task->SetUseOnFlyV0s(); 
-  task->SetCutnSigdEdx(2); 
-  task->ReadDeltaAOD(); // uncomment for DeltaAODs
-  task->SelectDeltaAODBranch(Form("jetsAOD_FASTJET04_Filter%05d",filterMask)); 
-  //task->SelectAODBranch("jetsAOD_FASTKT04");
-  //task->SelectAODBranch("jetsAOD_UA107");
-  //task->SelectAODBranch("jets");
-
-  task->SelectCollisionCandidates(); // either here or in userExec of task - but not for AODs ...
-
-  //  AliLog::SetGlobalLogLevel(AliLog::kInfo);    // kInfo // kDebug // kFatal
-  //  task->SetDebugLevel(10);
-  //  mgr->SetDebugLevel(10);
-
-  mgr->AddTask(task);
-
-
-  //================================================
-  //              data containers
-  //================================================
-  //            find input container
-  //below the trunk version
-  AliAnalysisDataContainer *cinput  = mgr->GetCommonInputContainer();
-  AliAnalysisDataContainer *coutput = mgr->CreateContainer("PWG4_JetChem", 
-                                                          TList::Class(), AliAnalysisManager::kOutputContainer,
-                                                          Form("%s:PWG4_JetChem",AliAnalysisManager::GetCommonFileName()));
+  
+  // Create the task and configure it.
+  //===========================================================================
+  AliAnalysisTaskJetChem *task = new AliAnalysisTaskJetChem("bla");
 
-  mgr->ConnectInput(task,0,cinput );
-  mgr->ConnectOutput(task,1,coutput);
-  return task;
+  Int_t debug = -1; // debug level
+  if(debug>=0) task->SetDebugLevel(debug);
+  
+  TString branchRecJets(recJetsBranch);
+  if(!branchRecJets.Contains("noRecJets")) task->SetBranchRecJets(branchRecJets);
+  
+  task->SetFilterMask(256); // for h+/h-
+  task->SetEventClass(eventClass);
+  
+  task->SetK0Type(K0type);
+
+  if(K0type == AliAnalysisTaskJetChem::kOnFlyPrim || AliAnalysisTaskJetChem::kOfflPrim) task->SetFilterMaskK0(256);
+
+  task->SetFFRadius(); 
+
+  task->SetTrackCuts(0.15, -0.75, 0.75, 0., 2*TMath::Pi());
+  task->SetJetCuts(5., -0.35, 0.35, 0., 2*TMath::Pi());
+
+  // Define histo bins
+   task->SetFFHistoBins();
+   task->SetQATrackHistoBins();
+   task->SetFFInvMassHistoBins();
+   task->SetPhiCorrInvMassHistoBins();
+   
+   mgr->AddTask(task);
+
+   // Create ONLY the output containers for the data produced by the task.
+   // Get and connect other common input/output containers via the manager as below
+   //==============================================================================
+
+   TString strK0type;
+   if(K0type ==  AliAnalysisTaskJetChem::kOnFly)     strK0type = "OnFly";
+   if(K0type ==  AliAnalysisTaskJetChem::kOnFlyPID)  strK0type = "OnFlyPID";
+   if(K0type ==  AliAnalysisTaskJetChem::kOnFlydEdx) strK0type = "OnFlydEdx";
+   if(K0type ==  AliAnalysisTaskJetChem::kOnFlyPrim) strK0type = "OnFlyPrim";
+   if(K0type ==  AliAnalysisTaskJetChem::kOffl)      strK0type = "Offl";
+   if(K0type ==  AliAnalysisTaskJetChem::kOfflPID)   strK0type = "OfflPID";
+   if(K0type ==  AliAnalysisTaskJetChem::kOffldEdx)  strK0type = "OffldEdx";
+   if(K0type ==  AliAnalysisTaskJetChem::kOfflPrim)  strK0type = "OfflPrim";
+   
+   TString listName(Form("PWG4_JetChem_%s_%s_cl%d",branchRecJets.Data(),strK0type.Data(),eventClass));
+                   
+   AliAnalysisDataContainer *coutput_FragFunc = mgr->CreateContainer(listName, 
+                                                                    TList::Class(),
+                                                                    AliAnalysisManager::kOutputContainer,
+                                                                    Form("%s:PWG4_JetChem",AliAnalysisManager::GetCommonFileName()));
+   
+   mgr->ConnectInput  (task, 0, mgr->GetCommonInputContainer());
+   mgr->ConnectOutput (task, 1, coutput_FragFunc);
+   
+   return task;
 }