add possibility to smear cluster energy, only for MC, some stetic changes in header...
authorgconesab <gconesab@f7af4fe6-9843-0410-8265-dc069ae4e863>
Wed, 9 Feb 2011 10:37:48 +0000 (10:37 +0000)
committergconesab <gconesab@f7af4fe6-9843-0410-8265-dc069ae4e863>
Wed, 9 Feb 2011 10:37:48 +0000 (10:37 +0000)
PWG4/PartCorrBase/AliCaloTrackReader.cxx
PWG4/PartCorrBase/AliCaloTrackReader.h

index 03fc257..a0bec6a 100755 (executable)
@@ -31,9 +31,9 @@
 
 // --- ROOT system ---
 #include "TFile.h"
+#include "TRandom3.h"
 
 //---- ANALYSIS system ----
-#include "AliCaloTrackReader.h"
 #include "AliMCEvent.h"
 #include "AliAODMCHeader.h"
 #include "AliGenPythiaEventHeader.h"
 #include "AliVParticle.h"
 #include "AliMixedEvent.h"
 #include "AliESDtrack.h"
-#include "AliEMCALRecoUtils.h"
 #include "AliESDtrackCuts.h"
 #include "AliTriggerAnalysis.h"
 
+//---- PartCorr/EMCAL ---
+#include "AliEMCALRecoUtils.h"
+#include "AliCaloTrackReader.h"
+
 ClassImp(AliCaloTrackReader)
   
   
@@ -59,7 +62,8 @@ ClassImp(AliCaloTrackReader)
     fEMCALCells(0x0), fPHOSCells(0x0),
     fInputEvent(0x0), fOutputEvent(0x0),fMC(0x0),
     fFillCTS(0),fFillEMCAL(0),fFillPHOS(0),
-    fFillEMCALCells(0),fFillPHOSCells(0),  fRemoveSuspiciousClusters(kFALSE),
+    fFillEMCALCells(0),fFillPHOSCells(0),  
+    fRemoveSuspiciousClusters(kFALSE), fSmearClusterEnergy(kFALSE), fRandom(),
 //    fSecondInputAODTree(0x0), fSecondInputAODEvent(0x0),
 //    fSecondInputFileName(""),fSecondInputFirstEvent(0), 
 //    fAODCTSNormalInputEntries(0), fAODEMCALNormalInputEntries(0), 
@@ -73,6 +77,7 @@ ClassImp(AliCaloTrackReader)
     fEMCALClustersListName(""),fZvtxCut(0.), 
     fDoEventSelection(kFALSE),   fDoV0ANDEventSelection(kFALSE),
     fTriggerAnalysis (new AliTriggerAnalysis), fCentralityClass("V0M"),fCentralityOpt(10)
+   
 {
   //Ctor
   
@@ -165,9 +170,9 @@ Bool_t AliCaloTrackReader::ComparePtHardAndJetPt(){
       //Compare jet pT and pt Hard
       //if(fDebug > 1) printf("AliMCAnalysisUtils:: %d pycell jet pT %f\n",ijet, jet->Pt());
       if(jet->Pt() > fPtHardAndJetPtFactor * ptHard) {
-       printf("AliMCAnalysisUtils::PythiaEventHeader: Njets: %d, pT Hard %2.2f, pycell jet pT %2.2f, rejection factor %1.1f\n",
-              nTriggerJets, ptHard, jet->Pt(), fPtHardAndJetPtFactor);
-       return kFALSE;
+        printf("AliMCAnalysisUtils::PythiaEventHeader: Njets: %d, pT Hard %2.2f, pycell jet pT %2.2f, rejection factor %1.1f\n",
+               nTriggerJets, ptHard, jet->Pt(), fPtHardAndJetPtFactor);
+        return kFALSE;
       }
     }
     if(jet) delete jet; 
@@ -335,6 +340,12 @@ void AliCaloTrackReader::InitParameters()
   //Centrality
   fCentralityBin[0]=fCentralityBin[1]=-1;
   
+  //Cluster smearing
+  fSmearClusterEnergy   = kFALSE;
+  fSmearClusterParam[0] = 0.07; // * sqrt E term
+  fSmearClusterParam[1] = 0.02; // * E term
+  fSmearClusterParam[2] = 0.00; // constant
+  
 }
 
 //________________________________________________________________
@@ -493,7 +504,7 @@ Bool_t AliCaloTrackReader::FillInputEvent(const Int_t iEntry, const char * curre
         }
       }//at least one cluster
       else {
-        printf("AliCaloTrackReader::FillInputEvent() - No clusters in event\n");
+        //printf("AliCaloTrackReader::FillInputEvent() - No clusters in event\n");
         //Remove events with  vertex (0,0,0), bad vertex reconstruction
         if(TMath::Abs(fVertex[0][0]) < 1.e-6 && TMath::Abs(fVertex[0][1]) < 1.e-6 && TMath::Abs(fVertex[0][2]) < 1.e-6) return kFALSE;
         
@@ -785,6 +796,7 @@ void AliCaloTrackReader::FillInputEMCALAlgorithm(AliVCluster * clus, const Int_t
 //    }
   }
   
+  
   TLorentzVector momentum ;
   
   clus->GetMomentum(momentum, fVertex[vindex]);      
@@ -828,6 +840,14 @@ void AliCaloTrackReader::FillInputEMCALAlgorithm(AliVCluster * clus, const Int_t
       //printf("Linearity Corrected Energy %f\n",clus->E());  
     }          
     
+    //In case of MC analysis, to match resolution/calibration in real data
+    if(fSmearClusterEnergy){
+      Float_t energy    = clus->E();
+      Float_t rdmEnergy = fRandom.Gaus(energy,fSmearClusterParam[0]*TMath::Sqrt(energy)+fSmearClusterParam[1]*energy+fSmearClusterParam[0]);
+      clus->SetE(rdmEnergy);
+      if(fDebug > 2) printf("\t Energy %f, smeared %f\n", energy, clus->E());
+    }
+    
     if (fMixedEvent) 
       clus->SetID(iclus) ; 
     
index 553f7b3..5b137f5 100755 (executable)
 
 // --- ROOT system ---
 #include "TObject.h" 
-class TObjArray ; 
 #include "TString.h"
+#include "TRandom3.h"
+class TObjArray ; 
 class TTree ;
 
 //--- ANALYSIS system ---
+#include "AliVCaloCells.h"
 class AliStack ; 
 class AliHeader ; 
 class AliGenEventHeader ; 
@@ -31,14 +33,15 @@ class AliVEvent;
 class AliAODEvent;
 class AliMCEvent;
 class AliMixedEvent;
-#include "AliVCaloCells.h"
-#include "AliFiducialCut.h"
 class AliAODMCHeader;
-#include "AliCalorimeterUtils.h"
 class AliESDtrackCuts;
 class AliCentrality;
 class AliTriggerAnalysis;
 
+// --- PartCorr
+#include "AliCalorimeterUtils.h"
+#include "AliFiducialCut.h"
+
 class AliCaloTrackReader : public TObject {
 
 public: 
@@ -49,113 +52,113 @@ private:
   AliCaloTrackReader & operator = (const AliCaloTrackReader & g) ;//cpy assignment
 
 public:
-  enum inputDataType {kESD, kAOD, kMC};
-    
-  //Select generated events, depending on comparison of pT hard and jets.
-  virtual Bool_t ComparePtHardAndJetPt() ;
-  virtual Bool_t IsPtHardAndJetPtComparisonSet() const {return  fComparePtHardAndJetPt ;}
-  virtual void SetPtHardAndJetPtComparison(Bool_t compare) { fComparePtHardAndJetPt = compare ;}       
-  virtual Float_t GetPtHardAndJetFactor() const {return  fPtHardAndJetPtFactor ;}
-  virtual void SetPtHardAndJetPtFactor(Float_t factor) { fPtHardAndJetPtFactor = factor ;}             
-       
-  virtual void InitParameters();
-  virtual void Print(const Option_t * opt) const;
+  
+  //--------------------------------
+  // General methods
+  //--------------------------------
+
+  virtual void    Init();
+  virtual void    InitParameters();
+  virtual void    Print(const Option_t * opt) const;
+  virtual void    ResetLists();
 
-  virtual Int_t GetDebug()         const { return fDebug ; }
-  virtual void  SetDebug(Int_t d)        { fDebug = d ; }
-  virtual Int_t GetDataType()      const { return fDataType ; }
-  virtual void  SetDataType(Int_t data ) { fDataType = data ; }
+  virtual Int_t   GetDebug()           const { return fDebug          ; }
+  virtual void    SetDebug(Int_t d)          { fDebug = d             ; }
+  
+  enum inputDataType {kESD, kAOD, kMC};
+  virtual Int_t   GetDataType()        const { return fDataType       ; }
+  virtual void    SetDataType(Int_t data )   { fDataType = data       ; }
 
-  virtual Int_t   GetEventNumber()     const {return fEventNumber ; }
+  virtual Int_t   GetEventNumber()     const {return fEventNumber     ; }
   virtual TString GetCurrentFileName() const {return fCurrentFileName ; }
        
-  //Minimum pt setters and getters 
-  virtual Float_t  GetEMCALPtMin() const { return fEMCALPtMin  ; }
-  virtual Float_t  GetPHOSPtMin()  const { return fPHOSPtMin  ; }
-  virtual Float_t  GetCTSPtMin()   const { return fCTSPtMin  ; }
+  TString         GetTaskName()        const {return fTaskName        ; }
+  void            SetTaskName(TString name)  {fTaskName = name        ; }
+    
 
-  virtual void SetEMCALPtMin(Float_t  pt) { fEMCALPtMin = pt ; }
-  virtual void SetPHOSPtMin(Float_t  pt)  { fPHOSPtMin = pt ; }
-  virtual void SetCTSPtMin(Float_t  pt)   { fCTSPtMin = pt ; }
+  //---------------------------------------
+  //Input/output event setters and getters
+  //---------------------------------------
+  virtual void    SetInputEvent(AliVEvent* const input) ;
+  virtual void    SetOutputEvent(AliAODEvent* const aod) {fOutputEvent = aod;}
+  virtual void    SetMC(AliMCEvent* const mc)            {fMC  = mc;}
+  virtual void    SetInputOutputMCEvent(AliVEvent* /*esd*/, AliAODEvent* /*aod*/, AliMCEvent* /*mc*/) {;}
+  
+  // Delta AODs
+  virtual TList * GetAODBranchList()      const { return fAODBranchList        ; }
+  void    SetDeltaAODFileName(TString name )    { fDeltaAODFileName = name     ; }
+  TString GetDeltaAODFileName()           const { return fDeltaAODFileName     ; }
+  void    SwitchOnWriteDeltaAOD()               { fWriteOutputDeltaAOD = kTRUE ; }
+  void    SwitchOffWriteDeltaAOD()              { fWriteOutputDeltaAOD = kFALSE; }
+  Bool_t  WriteDeltaAODToFile()           const { return fWriteOutputDeltaAOD  ; } 
+  
+  //------------------------------------------------------------
+  //Clusters/Tracks arrays filtering/filling methods and switchs 
+  //------------------------------------------------------------
+  
+  //Minimum pt setters and getters 
+  virtual Float_t  GetEMCALPtMin()      const { return fEMCALPtMin ; }
+  virtual Float_t  GetPHOSPtMin()       const { return fPHOSPtMin  ; }
+  virtual Float_t  GetCTSPtMin()        const { return fCTSPtMin   ; }
+  
+  virtual void     SetEMCALPtMin(Float_t  pt) { fEMCALPtMin = pt   ; }
+  virtual void     SetPHOSPtMin(Float_t  pt)  { fPHOSPtMin  = pt   ; }
+  virtual void     SetCTSPtMin(Float_t  pt)   { fCTSPtMin   = pt   ; }  
+  
+  // Fidutial cuts  
+  virtual AliFiducialCut * GetFiducialCut()   {if(!fFiducialCut) fFiducialCut = new AliFiducialCut(); 
+    return  fFiducialCut ;}
+  virtual void SetFiducialCut(AliFiducialCut * const fc) { fFiducialCut = fc ;}
+  virtual Bool_t IsFiducialCutOn()     const  { return fCheckFidCut ; }
+  virtual void SwitchOnFiducialCut()          { fCheckFidCut = kTRUE;  fFiducialCut = new AliFiducialCut();}
+  virtual void SwitchOffFiducialCut()         { fCheckFidCut = kFALSE;}
   
-  //Input setters and getters
+  // Cluster origin
   Bool_t IsEMCALCluster(AliVCluster *clus) const;
-  Bool_t IsPHOSCluster (AliVCluster *clus)  const;
-  void SwitchOnOldAODs()   {fOldAOD = kTRUE  ; }
-  void SwitchOffOldAODs()  {fOldAOD = kFALSE ; }
+  Bool_t IsPHOSCluster (AliVCluster *clus) const;
+  //Patch for cluster origin for Old AODs implementation
+  void   SwitchOnOldAODs()         { fOldAOD = kTRUE     ; }
+  void   SwitchOffOldAODs()        { fOldAOD = kFALSE    ; }
   
-  Bool_t IsCTSSwitchedOn()  const { return fFillCTS ; }
-  void SwitchOnCTS()    {fFillCTS = kTRUE  ; }
-  void SwitchOffCTS()   {fFillCTS = kFALSE ; }
+  // Cluster/track/cells switchs
+  Bool_t IsCTSSwitchedOn()   const { return fFillCTS     ; }
+  void   SwitchOnCTS()             { fFillCTS = kTRUE    ; }
+  void   SwitchOffCTS()            { fFillCTS = kFALSE   ; }
 
-  Bool_t IsEMCALSwitchedOn() const { return fFillEMCAL ; }
-  void SwitchOnEMCAL()  {fFillEMCAL = kTRUE  ; }
-  void SwitchOffEMCAL() {fFillEMCAL = kFALSE ; }
+  Bool_t IsEMCALSwitchedOn() const { return fFillEMCAL   ; }
+  void   SwitchOnEMCAL()           { fFillEMCAL = kTRUE  ; }
+  void   SwitchOffEMCAL()          { fFillEMCAL = kFALSE ; }
 
-  Bool_t IsPHOSSwitchedOn()  const { return fFillPHOS ; }
-  void SwitchOnPHOS()   {fFillPHOS = kTRUE  ; }
-  void SwitchOffPHOS()  {fFillPHOS = kFALSE ; }
+  Bool_t IsPHOSSwitchedOn()  const { return fFillPHOS    ; }
+  void   SwitchOnPHOS()            { fFillPHOS = kTRUE   ; }
+  void   SwitchOffPHOS()           { fFillPHOS = kFALSE  ; }
 
-  Bool_t IsEMCALCellsSwitchedOn() const { return fFillEMCALCells ; }
-  void SwitchOnEMCALCells()  {fFillEMCALCells = kTRUE  ; }
-  void SwitchOffEMCALCells() {fFillEMCALCells = kFALSE ; }
+  Bool_t IsEMCALCellsSwitchedOn() const { return fFillEMCALCells   ; }
+  void   SwitchOnEMCALCells()           { fFillEMCALCells = kTRUE  ; }
+  void   SwitchOffEMCALCells()          { fFillEMCALCells = kFALSE ; }
 
-  Bool_t IsPHOSCellsSwitchedOn()  const { return fFillPHOSCells ; }
-  void SwitchOnPHOSCells()   {fFillPHOSCells = kTRUE  ; }
-  void SwitchOffPHOSCells()  {fFillPHOSCells = kFALSE ; }
+  Bool_t IsPHOSCellsSwitchedOn()  const { return fFillPHOSCells   ; }
+  void   SwitchOnPHOSCells()            { fFillPHOSCells = kTRUE  ; }
+  void   SwitchOffPHOSCells()           { fFillPHOSCells = kFALSE ; }
 
+  // Filling/ filtering methods
   virtual Bool_t FillInputEvent(const Int_t iEntry, const char *currentFileName)  ;
-  virtual void FillInputCTS() ;
-  virtual void FillInputEMCAL() ;
-  virtual void FillInputEMCALAlgorithm(AliVCluster * clus, const Int_t iclus) ;
-  virtual void FillInputPHOS() ;
-  virtual void FillInputEMCALCells() ;
-  virtual void FillInputPHOSCells() ;
-  void SetEMCALClusterListName(TString &name) {fEMCALClustersListName = name;}
-
-  virtual TList * GetAODBranchList() const { return fAODBranchList ; }
-
-  virtual TObjArray* GetAODCTS()   const {return fAODCTS   ;}
-  virtual TObjArray* GetAODEMCAL() const {return fAODEMCAL ;}
-  virtual TObjArray* GetAODPHOS()  const {return fAODPHOS  ;}
-  virtual AliVCaloCells* GetEMCALCells()  const { return fEMCALCells ;}
-  virtual AliVCaloCells* GetPHOSCells()   const { return fPHOSCells  ;}
-
-  //Get MC  informatio
-  //Kinematics and galice.root available 
-  virtual AliStack*    GetStack()      const ;
-  virtual AliHeader*   GetHeader()     const ;
-  virtual AliGenEventHeader* GetGenEventHeader() const ;
-  //Filtered kinematics in AOD 
-  virtual TClonesArray*   GetAODMCParticles(Int_t input = 0) const ;
-  virtual AliAODMCHeader* GetAODMCHeader(Int_t input = 0)    const ;
-       
-  virtual AliVEvent*     GetInputEvent()  const {return fInputEvent;}
-  virtual AliAODEvent*   GetOutputEvent() const {return fOutputEvent;}
-  virtual AliMCEvent*    GetMC()          const {return fMC;}
-  virtual AliMixedEvent* GetMixedEvent()  const {return fMixedEvent;}
-  virtual Int_t          GetNMixedEvent() const {return fNMixedEvent ; } 
-
-  virtual Double_t       GetBField() const { return 0.;}
-       
-  virtual void Init();
-
-//  virtual void SetInputEvent(AliVEvent* const input)  {fInputEvent  = input;}
-  virtual void SetInputEvent(AliVEvent* const input) ;
-  virtual void SetOutputEvent(AliAODEvent* const aod) {fOutputEvent = aod;}
-  virtual void SetMC(AliMCEvent* const mc)            {fMC  = mc;}
+  virtual void   FillInputCTS() ;
+  virtual void   FillInputEMCAL() ;
+  virtual void   FillInputEMCALAlgorithm(AliVCluster * clus, const Int_t iclus) ;
+  virtual void   FillInputPHOS() ;
+  virtual void   FillInputEMCALCells() ;
+  virtual void   FillInputPHOSCells() ;
+  void           SetEMCALClusterListName(TString &name) {fEMCALClustersListName = name ; }
+
+  // Arrayes with clusters/track/cells access methods
+  //FIXME CHANGE NAMES, remove AOD, now they are VClusters, VTracks ...
+  virtual TObjArray*     GetAODCTS()      const { return fAODCTS     ; }
+  virtual TObjArray*     GetAODEMCAL()    const { return fAODEMCAL   ; }
+  virtual TObjArray*     GetAODPHOS()     const { return fAODPHOS    ; }
+  virtual AliVCaloCells* GetEMCALCells()  const { return fEMCALCells ; }
+  virtual AliVCaloCells* GetPHOSCells()   const { return fPHOSCells  ; }
 
-  virtual void ResetLists();
-
-  virtual AliFiducialCut * GetFiducialCut() {if(!fFiducialCut) fFiducialCut = new AliFiducialCut(); 
-         return  fFiducialCut ;}
-  virtual void SetFiducialCut(AliFiducialCut * const fc) { fFiducialCut = fc ;}
-  virtual Bool_t IsFiducialCutOn()       {return fCheckFidCut ; }
-  virtual void SwitchOnFiducialCut()     { fCheckFidCut = kTRUE;  fFiducialCut = new AliFiducialCut();}
-  virtual void SwitchOffFiducialCut()    { fCheckFidCut = kFALSE;}
-       
-  virtual void SetInputOutputMCEvent(AliVEvent* /*esd*/, AliAODEvent* /*aod*/, AliMCEvent* /*mc*/) {;}
-       
   //Methods for mixing with external input file (AOD)
   //virtual TTree* GetSecondInputAODTree() const {return  fSecondInputAODTree ; } 
   //virtual void SetSecondInputAODTree(TTree * tree) {fSecondInputAODTree = tree ;
@@ -176,101 +179,139 @@ public:
 //  Int_t GetAODPHOSNormalInputEntries()  {if(!fSecondInputAODTree) { fAODPHOSNormalInputEntries  = fAODPHOS->GetEntriesFast() ;}
 //                                                                              return fAODPHOSNormalInputEntries ; }
        
+   
+  //-------------------------------------
+  // Event/track selection methods
+  //-------------------------------------
+  
+  void    SetFiredTriggerClassName(TString name ) { fFiredTriggerClassName = name ; }
+  TString GetFiredTriggerClassName()        const { return fFiredTriggerClassName ; }
+  virtual TString GetFiredTriggerClasses()        { return ""                     ; } // look the ESD/AOD reader 
+  
+  void    SwitchOnEventSelection()                { fDoEventSelection = kTRUE  ; }
+  void    SwitchOffEventSelection()               { fDoEventSelection = kFALSE ; }
+  Bool_t  IsEventSelectionDone()            const { return fDoEventSelection   ; } 
+  
+  void    SwitchOnV0ANDSelection()                { fDoV0ANDEventSelection = kTRUE  ; }
+  void    SwitchOffV0ANDSelection()               { fDoV0ANDEventSelection = kFALSE ; }
+  Bool_t  IsV0ANDEventSelectionDone()       const { return fDoV0ANDEventSelection   ; } 
+  
   // Track selection
-  ULong_t GetTrackStatus() const   {return fTrackStatus ; }
-  void SetTrackStatus(ULong_t bit) { fTrackStatus = bit ; }            
+  ULong_t GetTrackStatus()                  const {return fTrackStatus      ; }
+  void    SetTrackStatus(ULong_t bit)             { fTrackStatus = bit      ; }                
        
-  AliESDtrackCuts* GetTrackCuts()          const  { return fESDtrackCuts    ; }
+  AliESDtrackCuts* GetTrackCuts()           const { return fESDtrackCuts    ; }
   void    SetTrackCuts(AliESDtrackCuts * cuts)    { fESDtrackCuts = cuts    ; }                  
-  Int_t   GetTrackMultiplicity()           const  { return fTrackMult       ; }
-  Float_t GetTrackMultiplicityEtaCut()     const  { return fTrackMultEtaCut ; }
+  Int_t   GetTrackMultiplicity()            const { return fTrackMult       ; }
+  Float_t GetTrackMultiplicityEtaCut()      const { return fTrackMultEtaCut ; }
   void    SetTrackMultiplicityEtaCut(Float_t eta) { fTrackMultEtaCut = eta  ; }                
   
-  //MC switchs
-  void SwitchOnStack()              { fReadStack          = kTRUE  ; }
-  void SwitchOffStack()             { fReadStack          = kFALSE ; }
-  void SwitchOnAODMCParticles()     { fReadAODMCParticles = kTRUE  ; }
-  void SwitchOffAODMCParticles()    { fReadAODMCParticles = kFALSE ; }
-  Bool_t ReadStack()          const { return fReadStack            ; }
-  Bool_t ReadAODMCParticles() const { return fReadAODMCParticles   ; }
-       
-  void SetDeltaAODFileName(TString name ) {fDeltaAODFileName = name ; }
-  TString GetDeltaAODFileName() const {return fDeltaAODFileName ; }
-
-  void SetFiredTriggerClassName(TString name ) {fFiredTriggerClassName = name ; }
-  TString GetFiredTriggerClassName() const {return fFiredTriggerClassName ; }
-  virtual TString GetFiredTriggerClasses() {return "";}
-
-  void AnalyzeOnlyLED()     {fAnaLED = kTRUE;}
-  void AnalyzeOnlyPhysics() {fAnaLED = kFALSE;}
-       
-  TString  GetTaskName() const {return fTaskName;}
-  void SetTaskName(TString name) {fTaskName = name;}
-
-  AliCalorimeterUtils * GetCaloUtils() const {return fCaloUtils ; }
-  void SetCaloUtils(AliCalorimeterUtils * caloutils) { fCaloUtils = caloutils ; }
+  // Calorimeter specific and patches
+  void    AnalyzeOnlyLED()                        {fAnaLED = kTRUE  ; }
+  void    AnalyzeOnlyPhysics()                    {fAnaLED = kFALSE ; }
   
+  void    SwitchOnCaloFilterPatch()               {fCaloFilterPatch = kTRUE ; fFillCTS = kFALSE    ; }
+  void    SwitchOffCaloFilterPatch()              {fCaloFilterPatch = kFALSE                       ; }
+  Bool_t  IsCaloFilterPatchOn()             const {if(fDataType == kAOD) { return fCaloFilterPatch ; } 
+                                                   else                  { return kFALSE           ; } }
+       
+  //-------------------------------
   //Vertex methods
-  virtual void GetVertex(Double_t v[3]) const ;
+  //-------------------------------
+  virtual void      GetVertex(Double_t v[3])        const ;
   virtual Double_t* GetVertex(const Int_t evtIndex) const {return fVertex[evtIndex];}
-  virtual void GetVertex(Double_t vertex[3], const Int_t evtIndex) const ;
-  virtual void FillVertexArray();
-  virtual Bool_t CheckForPrimaryVertex();
+  virtual void      GetVertex(Double_t vertex[3],   const Int_t evtIndex) const ;
+  virtual void      FillVertexArray();
+  virtual Bool_t    CheckForPrimaryVertex();
  // virtual void       GetSecondInputAODVertex(Double_t *) const {;}
-  virtual Float_t GetZvertexCut() const {return fZvtxCut ;} //cut on vertex position  
-  virtual void SetZvertexCut(Float_t zcut=10.){fZvtxCut=zcut ;} //cut on vertex position
+  virtual Float_t   GetZvertexCut()                 const {return fZvtxCut ;} //cut on vertex position  
+  virtual void      SetZvertexCut(Float_t zcut=10.)       {fZvtxCut=zcut   ;} //cut on vertex position
 
-  void SwitchOnWriteDeltaAOD()  {fWriteOutputDeltaAOD = kTRUE ;  }
-  void SwitchOffWriteDeltaAOD() {fWriteOutputDeltaAOD = kFALSE ; }
-  Bool_t WriteDeltaAODToFile() const {return fWriteOutputDeltaAOD ; } 
+  //------------------------
+  // Centrality
+  //------------------------
+  virtual AliCentrality* GetCentrality() const {return 0x0;} //Actual method to recover the pointer is in the ESD/AODReader
+  void    SetCentralityClass(TString name)     { fCentralityClass   = name       ;}
+  void    SetCentralityOpt(Int_t opt)          { fCentralityOpt     = opt        ;}
+  TString GetCentralityClass()           const { return fCentralityClass         ;}
+  Int_t   GetCentralityOpt()             const { return fCentralityOpt           ;}
+  Int_t   GetEventCentrality()           const ;
+  void    SetCentralityBin(Int_t min, Int_t max) //Set the centrality bin to select the event. If used, then need to get percentile
+                                               {fCentralityBin[0]=min; fCentralityBin[1]=max;  
+                                                if(min>=0 && max > 0) fCentralityOpt = 100 ; }
+  Float_t GetCentralityBin(Int_t i)     const  { if(i < 0 || i > 1) return 0 ; 
+                                                else return fCentralityBin[i]   ; }
   
-  void SwitchOnSuspiciousClustersRemoval()  {fRemoveSuspiciousClusters = kTRUE ;  }
-  void SwitchOffSuspiciousClustersRemoval() {fRemoveSuspiciousClusters = kFALSE ; }
-  Bool_t IsSuspiciousClustersRemovalOn() const {return fRemoveSuspiciousClusters ; } 
-
-  virtual void FillInputVZERO(){;}
-  Int_t GetV0Signal(Int_t i) const { return fV0ADC[i];}
-  Int_t GetV0Multiplicity(Int_t i)   const { return fV0Mul[i];}
-
-  void SwitchOnCaloFilterPatch()  {fCaloFilterPatch = kTRUE ; fFillCTS = kFALSE; }
-  void SwitchOffCaloFilterPatch() {fCaloFilterPatch = kFALSE ; }
-  Bool_t IsCaloFilterPatchOn()    {if(fDataType == kAOD) { return fCaloFilterPatch ; } 
-                                   else                  { return kFALSE ; } }
+  //-------------------------------------
+  // Other methods
+  //-------------------------------------
+  AliCalorimeterUtils * GetCaloUtils()            const { return fCaloUtils      ; }
+  void   SetCaloUtils(AliCalorimeterUtils * caloutils)  { fCaloUtils = caloutils ; }  
   
-  void SwitchOnEventSelection()  {fDoEventSelection = kTRUE  ; }
-  void SwitchOffEventSelection() {fDoEventSelection = kFALSE ; }
-  Bool_t IsEventSelectionDone()  { return fDoEventSelection  ; } 
+  void   SwitchOnSuspiciousClustersRemoval()            { fRemoveSuspiciousClusters = kTRUE ;  }
+  void   SwitchOffSuspiciousClustersRemoval()           { fRemoveSuspiciousClusters = kFALSE ; }
+  Bool_t IsSuspiciousClustersRemovalOn()          const { return fRemoveSuspiciousClusters   ; } 
   
-  void SwitchOnV0ANDSelection()  {fDoV0ANDEventSelection = kTRUE  ; }
-  void SwitchOffV0ANDSelection() {fDoV0ANDEventSelection = kFALSE ; }
-  Bool_t IsV0ANDEventSelectionDone() { return fDoV0ANDEventSelection  ; } 
-
-  //Centrality
-  virtual AliCentrality* GetCentrality() const {return 0x0;} //Actual method to recover the pointer is in the ESD/AODReader
-  void SetCentralityClass(TString name)    { fCentralityClass   = name       ;}
-  void SetCentralityOpt(Int_t opt)         { fCentralityOpt     = opt        ;}
-  TString GetCentralityClass()      const  { return fCentralityClass         ;}
-  Int_t   GetCentralityOpt()        const  { return fCentralityOpt           ;}
-  Int_t   GetEventCentrality()      const ;
-  void    SetCentralityBin(Int_t min, Int_t max) //Set the centrality bin to select the event. If used, then need to get percentile
-    {fCentralityBin[0]=min; fCentralityBin[1]=max;  if(min>=0 && max > 0) fCentralityOpt = 100; }
-  Float_t GetCentralityBin(Int_t i) const  { if(i < 0 || i > 1) return 0 ; else return fCentralityBin[i] ; }
+  //Use only for MC
+  void    SwitchOnClusterEnergySmearing()               { fSmearClusterEnergy = kTRUE    ; }
+  void    SwitchOffClusterEnergySmearing()              { fSmearClusterEnergy = kFALSE   ; }
+  Bool_t  IsClusterEnergySmeared()                const { return fSmearClusterEnergy     ; }   
+  void    SetSmearingParameters(Int_t i, Float_t param) { if(i < 3)fSmearClusterParam[i] = param  ; }
+  
+  virtual void FillInputVZERO()                         { ; } // done in the AOD/ESD reader
+  Int_t   GetV0Signal(Int_t i)                    const { return fV0ADC[i] ; }
+  Int_t   GetV0Multiplicity(Int_t i)              const { return fV0Mul[i] ; }
+  
+  virtual Double_t GetBField()                    const { return 0.        ; } // get in the AOD/ESD reader
+  
+  //------------------------------------------------
+  // MC analysis specific methods
+  //-------------------------------------------------
+  
+  //Kinematics and galice.root available 
+  virtual AliStack*          GetStack()          const ;
+  virtual AliHeader*         GetHeader()         const ;
+  virtual AliGenEventHeader* GetGenEventHeader() const ;
+  
+  //Filtered kinematics in AOD 
+  virtual TClonesArray*     GetAODMCParticles(Int_t input = 0) const ;
+  virtual AliAODMCHeader*   GetAODMCHeader(Int_t input = 0)    const ;
+       
+  virtual AliVEvent*        GetInputEvent()  const { return fInputEvent  ; }
+  virtual AliAODEvent*      GetOutputEvent() const { return fOutputEvent ; }
+  virtual AliMCEvent*       GetMC()          const { return fMC          ; }
+  virtual AliMixedEvent*    GetMixedEvent()  const { return fMixedEvent  ; }
+  virtual Int_t             GetNMixedEvent() const { return fNMixedEvent ; } 
+  
+  void   SwitchOnStack()                           { fReadStack          = kTRUE  ; }
+  void   SwitchOffStack()                          { fReadStack          = kFALSE ; }
+  void   SwitchOnAODMCParticles()                  { fReadAODMCParticles = kTRUE  ; }
+  void   SwitchOffAODMCParticles()                 { fReadAODMCParticles = kFALSE ; }
+  Bool_t ReadStack()                         const { return fReadStack            ; }
+  Bool_t ReadAODMCParticles()                const { return fReadAODMCParticles   ; }
+       
+  //Select generated events, depending on comparison of pT hard and jets.
+  virtual Bool_t  ComparePtHardAndJetPt() ;
+  virtual Bool_t  IsPtHardAndJetPtComparisonSet()       const { return  fComparePtHardAndJetPt   ; }
+  virtual void    SetPtHardAndJetPtComparison(Bool_t compare) { fComparePtHardAndJetPt = compare ; }   
+  virtual Float_t GetPtHardAndJetFactor()               const { return  fPtHardAndJetPtFactor    ; }
+  virtual void    SetPtHardAndJetPtFactor(Float_t factor)     { fPtHardAndJetPtFactor = factor   ; }           
   
-  //MC reader methods:
+  //MC reader methods, declared there to allow compilation, they are only used in the MC reader:
   
   virtual void AddNeutralParticlesArray(TArrayI & /*array*/)  { ; }  
   virtual void AddChargedParticlesArray(TArrayI & /*array*/)  { ; } 
   virtual void AddStatusArray(TArrayI & /*array*/)            { ; }
   
-  virtual void SwitchOnPi0Decay()         { ; } 
-  virtual void SwitchOffPi0Decay()        { ; } 
-  virtual void SwitchOnStatusSelection()  { ; }
-  virtual void SwitchOffStatusSelection() { ; }
-  virtual void SwitchOnOverlapCheck()     { ; }
-  virtual void SwitchOffOverlapCheck()    { ; }
+  virtual void SwitchOnPi0Decay()                             { ; } 
+  virtual void SwitchOffPi0Decay()                            { ; } 
+  virtual void SwitchOnStatusSelection()                      { ; }
+  virtual void SwitchOffStatusSelection()                     { ; }
+  virtual void SwitchOnOverlapCheck()                         { ; }
+  virtual void SwitchOffOverlapCheck()                        { ; }
   
-  virtual void SetEMCALOverlapAngle(Float_t /*angle*/)  { ; }
-  virtual void SetPHOSOverlapAngle(Float_t /*angle*/)   { ; }
+  virtual void SetEMCALOverlapAngle(Float_t /*angle*/)        { ; }
+  virtual void SetPHOSOverlapAngle(Float_t /*angle*/)         { ; }
 
   
  protected:
@@ -289,6 +330,7 @@ public:
   Float_t          fPHOSPtMin;      // pT Threshold on phos clusters
 
   TList          * fAODBranchList ; //! List with AOD branches created and needed in analysis  
+  //FIXME CHANGE NAMES, remove AOD, now they are VClusters, VTracks ...
   TObjArray      * fAODCTS ;        //! temporal referenced array with tracks
   TObjArray      * fAODEMCAL ;      //! temporal referenced array with EMCAL CaloClusters
   TObjArray      * fAODPHOS ;       //! temporal referenced array with PHOS CaloClusters
@@ -304,7 +346,10 @@ public:
   Bool_t           fFillPHOS;       // use data from PHOS
   Bool_t           fFillEMCALCells; // use data from EMCAL
   Bool_t           fFillPHOSCells;  // use data from PHOS
-  Bool_t           fRemoveSuspiciousClusters; //Remove high energy clusters with low number of cells
+  Bool_t           fRemoveSuspiciousClusters; // Remove high energy clusters with low number of cells
+  Bool_t           fSmearClusterEnergy;       // Smear cluster energy, to be done only for simulated data to match real data
+  Float_t          fSmearClusterParam[3];     // Smearing parameters
+  TRandom3         fRandom;                   // Random generator
   
 //  TTree *        fSecondInputAODTree;    // Tree with second input AOD, for mixing analysis. 
 //  AliAODEvent*   fSecondInputAODEvent;   //! pointer to second input AOD event.
@@ -353,7 +398,7 @@ public:
   Int_t            fCentralityOpt;      // Option for the returned value of the centrality, possible options 5, 10, 100
   Int_t            fCentralityBin[2];   // Minimum and maximum value of the centrality for the analysis
   
-  ClassDef(AliCaloTrackReader,24)
+  ClassDef(AliCaloTrackReader,25)
 } ;