From 030b69fb506b7dc7ff74cd013cde850dc0de9a35 Mon Sep 17 00:00:00 2001
From: Jean-Baptiste De Vivie De Regie <devivie@lal.in2p3.fr>
Date: Fri, 9 Dec 2016 12:46:20 +0100
Subject: [PATCH] writing very fwd density (IsolationAlgs-00-01-32)

	* Adding the very forward density (might be removed again later)
	* tagging as IsolationAlgs-00-01-32

06-12-2016 Jean-Baptiste de Vivie
	* bug fix (wrong deletion of a pointer)
	* tagging as IsolationAlgs-00-01-31

26-10-2016 Jean-Baptiste de Vivie
	* Modification to split the AODFix into one component per collection (so 4 : el, ph, fwdel, mu)
	* tagging as IsolationAlgs-00-01-30

26-10-2016 Jean-Baptiste de Vivie
	* Minor modif for IsolationCorrections config that detects automatically if mc or data
	* tagging as IsolationAlgs-00-01-29

26-10-2016 Jean-Baptiste de Vivie
	* Changes for AODFix and fix for CMakeList
	* tagging as IsolationAlgs-00-01-28

21-10-2016 Jean-Baptiste de Vivie
...
(Long ChangeLog diff - truncated)
---
 .../RecoAlgs/IsolationAlgs/CMakeLists.txt     |   5 +-
 .../RecoAlgs/IsolationAlgs/cmt/requirements   |   2 +
 .../IsolationAlgs/python/IsoAODFixGetter.py   | 229 ++++++++++++
 .../IsolationAlgs/python/IsoGetter.py         |  76 ++--
 .../IsoEventShapeOutputItemList_jobOptions.py |   2 +
 .../IsolationAlgs/src/IsolationBuilder.cxx    | 350 ++++++++++++++++--
 .../IsolationAlgs/src/IsolationBuilder.h      |  23 +-
 7 files changed, 634 insertions(+), 53 deletions(-)
 create mode 100644 Reconstruction/RecoAlgs/IsolationAlgs/python/IsoAODFixGetter.py

diff --git a/Reconstruction/RecoAlgs/IsolationAlgs/CMakeLists.txt b/Reconstruction/RecoAlgs/IsolationAlgs/CMakeLists.txt
index dcaf1d2a4b2..2b8b9a4efe1 100644
--- a/Reconstruction/RecoAlgs/IsolationAlgs/CMakeLists.txt
+++ b/Reconstruction/RecoAlgs/IsolationAlgs/CMakeLists.txt
@@ -15,7 +15,8 @@ atlas_depends_on_subdirs( PUBLIC
                           Event/xAOD/xAODEgamma
                           Event/xAOD/xAODMuon
                           Event/xAOD/xAODPrimitives
-                          Reconstruction/RecoTools/RecoToolInterfaces )
+                          Reconstruction/RecoTools/RecoToolInterfaces
+			  PhysicsAnalysis/ElectronPhotonID/IsolationCorrections )
 
 # External dependencies:
 find_package( Boost COMPONENTS filesystem thread system )
@@ -25,7 +26,7 @@ atlas_add_component( IsolationAlgs
                      src/*.cxx
                      src/components/*.cxx
                      INCLUDE_DIRS ${Boost_INCLUDE_DIRS}
-                     LINK_LIBRARIES ${Boost_LIBRARIES} GaudiKernel CaloEvent AthContainers AthenaBaseComps xAODEgamma xAODMuon xAODPrimitives RecoToolInterfaces )
+                     LINK_LIBRARIES ${Boost_LIBRARIES} GaudiKernel CaloEvent AthContainers AthenaBaseComps xAODEgamma xAODMuon xAODPrimitives RecoToolInterfaces IsolationCorrectionsLib)
 
 # Install files from the package:
 atlas_install_python_modules( python/*.py )
diff --git a/Reconstruction/RecoAlgs/IsolationAlgs/cmt/requirements b/Reconstruction/RecoAlgs/IsolationAlgs/cmt/requirements
index 0a1721b5641..c1a1d5e445d 100644
--- a/Reconstruction/RecoAlgs/IsolationAlgs/cmt/requirements
+++ b/Reconstruction/RecoAlgs/IsolationAlgs/cmt/requirements
@@ -16,6 +16,8 @@ use xAODEgamma         xAODEgamma-*         Event/xAOD
 use xAODPrimitives     xAODPrimitives-*     Event/xAOD
 use xAODMuon           xAODMuon-*           Event/xAOD
 use CaloEvent          CaloEvent-*          Calorimeter
+use IsolationCorrections IsolationCorrections-*   PhysicsAnalysis/ElectronPhotonID   
+
 ##
 
 branches src src/components doc python share
diff --git a/Reconstruction/RecoAlgs/IsolationAlgs/python/IsoAODFixGetter.py b/Reconstruction/RecoAlgs/IsolationAlgs/python/IsoAODFixGetter.py
new file mode 100644
index 00000000000..2fbb1563e95
--- /dev/null
+++ b/Reconstruction/RecoAlgs/IsolationAlgs/python/IsoAODFixGetter.py
@@ -0,0 +1,229 @@
+# Copyright (C) 2002-2017 CERN for the benefit of the ATLAS collaboration
+
+from AthenaCommon.Constants import INFO
+from AthenaCommon.AppMgr import ToolSvc
+from AthenaCommon.Logging import logging
+from RecExConfig.RecFlags import rec
+
+from egammaRec.Factories import ToolFactory, AlgFactory, getPropertyValue
+
+#---------------------------------------
+
+# tool to collect topo clusters in cone
+from ParticlesInConeTools.ParticlesInConeToolsConf import xAOD__CaloClustersInConeTool
+CaloClustersInConeTool = ToolFactory(xAOD__CaloClustersInConeTool,
+                                     CaloClusterLocation = "CaloCalTopoClusters")
+
+# configuration for ED computation
+def configureEDCorrection(tool):
+  """Configure tools and algorithm for energy density correction 
+     (only if doEnergyDensityCorrection = True)"""
+  # Return if doEnergyDensityCorrection is false
+  if not getPropertyValue(tool, 'doEnergyDensityCorrection'):
+    return
+  # Set OutputLevel to INFO or higher if tool has it too
+  OutputLevel = min(getPropertyValue(tool, 'OutputLevel'), INFO)
+  try:
+    from AthenaCommon.AppMgr import ToolSvc
+    from AthenaCommon.AlgSequence import AlgSequence
+    from EventShapeTools.EventDensityConfig import configEventDensityTool, EventDensityAlg
+    from JetRec.JetRecStandard import jtm
+    topSequence = AlgSequence()
+    if not hasattr(topSequence,'EDtpIsoVeryForwardAlg'):
+      tfcc = configEventDensityTool("EDtpIsoVeryForwardTool", jtm.emget,
+                                    radius          = 0.5,
+                                    AbsRapidityMin  = 2.5,
+                                    AbsRapidityMax  = 4.5,
+                                    OutputContainer = "TopoClusterIsoVeryForwardEventShape",
+                                    OutputLevel     = OutputLevel
+                                    )
+      ToolSvc += tfcc
+      topSequence += EventDensityAlg("EDtpIsoVeryForwardAlg", EventDensityTool = tfcc)
+
+  except Exception:
+    print '\nERROR: could not get handle to ED'
+    raise
+
+from CaloIdentifier import SUBCALO 
+from IsolationTool.IsolationToolConf import xAOD__CaloIsolationTool, xAOD__TrackIsolationTool
+CaloIsolationTool = ToolFactory(xAOD__CaloIsolationTool,name = "CaloIsolationTool",
+                                postInit                        = [configureEDCorrection],
+                                CaloFillRectangularClusterTool  = None,
+                                ClustersInConeTool              = CaloClustersInConeTool,
+                                PFlowObjectsInConeTool          = None,
+                                IsoLeakCorrectionTool           = None,
+                                saveOnlyRequestedCorrections    = True,
+                                EMCaloNums                      = [],
+                                HadCaloNums                     = [],
+                                UseEMScale                      = True,
+                                addCaloExtensionDecoration      = False,
+                                OutputLevel                     = 3)
+
+TrackIsolationTool = ToolFactory(xAOD__TrackIsolationTool, name = 'TrackIsolationTool')
+from AthenaCommon import CfgMgr
+tit = CfgMgr.xAOD__TrackIsolationTool('TrackIsolationTool')
+tit.TrackSelectionTool.maxZ0SinTheta = 3
+tit.TrackSelectionTool.minPt         = 1000
+tit.TrackSelectionTool.CutLevel      = "Loose"
+
+import ROOT, cppyy
+# Need to be sure base dict is loaded first.
+cppyy.loadDictionary('xAODCoreRflxDict')
+cppyy.loadDictionary('xAODPrimitivesDict')
+isoPar = ROOT.xAOD.Iso
+
+# The types to be computed
+IsoTypesEg =  [
+  [ isoPar.ptcone40, 
+    isoPar.ptcone30,
+    isoPar.ptcone20 ]
+  ]
+IsoTypesFe =  [
+  [ isoPar.topoetcone20, 
+    isoPar.topoetcone30,
+    isoPar.topoetcone40 ]
+  ]
+IsoTypesMu =  [
+  [ isoPar.topoetcone20, 
+    isoPar.topoetcone30,
+    isoPar.topoetcone40 ]
+  ]
+# And the corrections
+IsoCorEg = [
+  [ isoPar.coreTrackPtr ] 
+  ]
+#
+IsoCorFe = [
+  [ isoPar.coreCone, isoPar.pileupCorrection ] 
+  ]
+#
+IsoCorMu = [
+  [ isoPar.coreCone, isoPar.pileupCorrection ] 
+  ]
+
+# 
+
+from IsolationCorrections.IsolationCorrectionsConf import CP__IsolationCorrectionTool as ict
+leakTool = ict(name     = "LeakageCorrection",
+               CorrFile = "IsolationCorrections/isolation_ptcorrections_rel20_2.root") # in principle the default
+ToolSvc += leakTool
+
+from IsolationAlgs.IsolationAlgsConf import IsolationBuilder
+isoAODFixBuilderElectron = AlgFactory(IsolationBuilder,
+                                      name                  = "IsolationBuilderElectron",
+                                      ElectronCollectionContainerName    = "Electrons",
+                                      PhotonCollectionContainerName      = "",
+                                      MuonCollectionContainerName        = "",
+                                      FwdElectronCollectionContainerName = "",
+                                      CaloCellIsolationTool = None,
+                                      CaloTopoIsolationTool = None,
+                                      PFlowIsolationTool    = None,
+                                      TrackIsolationTool    = None, 
+                                      FeIsoTypes            = [[]] ,
+                                      FeCorTypes            = IsoCorFe,
+                                      EgIsoTypes            = [[]],
+                                      EgCorTypes            = IsoCorEg,
+                                      MuIsoTypes            = [[]] ,
+                                      MuCorTypes            = IsoCorMu,
+                                      IsAODFix              = True,
+                                      LeakageTool           = leakTool,
+                                      IsolateEl             = False,
+                                      OutputLevel           = 3)
+
+isoAODFixBuilderPhoton = AlgFactory(IsolationBuilder,
+                                    name                  = "IsolationBuilderPhoton",
+                                    ElectronCollectionContainerName    = "",
+                                    PhotonCollectionContainerName      = "Photons",
+                                    MuonCollectionContainerName        = "",
+                                    FwdElectronCollectionContainerName = "",
+                                    CaloCellIsolationTool = None,
+                                    CaloTopoIsolationTool = None,
+                                    PFlowIsolationTool    = None,
+                                    TrackIsolationTool    = TrackIsolationTool, 
+                                    FeIsoTypes            = [[]] ,
+                                    FeCorTypes            = IsoCorFe,
+                                    EgIsoTypes            = IsoTypesEg,
+                                    EgCorTypes            = IsoCorEg,
+                                    MuIsoTypes            = [[]] ,
+                                    MuCorTypes            = IsoCorMu,
+                                    IsAODFix              = True,
+                                    LeakageTool           = leakTool,
+                                    IsolateEl             = False,
+                                    OutputLevel           = 3)
+
+isoAODFixBuilderMuon = AlgFactory(IsolationBuilder,
+                                  name                  = "IsolationBuilderMuon",
+                                  ElectronCollectionContainerName    = "",
+                                  PhotonCollectionContainerName      = "",
+                                  MuonCollectionContainerName        = "Muons",
+                                  FwdElectronCollectionContainerName = "",
+                                  CaloCellIsolationTool = None,
+                                  CaloTopoIsolationTool = CaloIsolationTool,
+                                  PFlowIsolationTool    = None,
+                                  TrackIsolationTool    = None, 
+                                  FeIsoTypes            = [[]] ,
+                                  FeCorTypes            = IsoCorFe,
+                                  EgIsoTypes            = [[]],
+                                  EgCorTypes            = IsoCorEg,
+                                  MuIsoTypes            = IsoTypesMu,
+                                  MuCorTypes            = IsoCorMu,
+                                  CustomConfigurationNameMu = 'Core0p05',
+                                  IsAODFix              = True,
+                                  LeakageTool           = None,
+                                  IsolateEl             = False,
+                                  OutputLevel           = 3)
+
+
+isoAODFixBuilderFwdElectron = AlgFactory(IsolationBuilder,
+                                         name                  = "IsolationBuilderForwardElectron",
+                                         ElectronCollectionContainerName    = "",
+                                         PhotonCollectionContainerName      = "",
+                                         MuonCollectionContainerName        = "",
+                                         FwdElectronCollectionContainerName = "ForwardElectrons",
+                                         CaloCellIsolationTool = None,
+                                         CaloTopoIsolationTool = CaloIsolationTool,
+                                         PFlowIsolationTool    = None,
+                                         TrackIsolationTool    = None, 
+                                         FeIsoTypes            = IsoTypesFe,
+                                         FeCorTypes            = IsoCorFe,
+                                         EgIsoTypes            = [[]],
+                                         EgCorTypes            = IsoCorEg,
+                                         MuIsoTypes            = [[]],
+                                         MuCorTypes            = IsoCorMu,
+                                         IsAODFix              = True,
+                                         LeakageTool           = None,
+                                         IsolateEl             = False,
+                                         OutputLevel           = 3)
+
+
+from RecExConfig.Configured import Configured
+class isoAODFixGetter ( Configured ) :
+    
+  def __init__(self, intype = ""):
+    super( isoAODFixGetter, self ).__init__()
+    self.type = intype
+
+    mlog = logging.getLogger ('isoAODFixGetter.py::configure:')
+    mlog.info('entering')
+    
+    # configure iso here:
+    try:
+      if self.type == "Electrons":
+        self._isoBuilderHandle = isoAODFixBuilderElectron()
+      elif self.type == "Photons":
+        self._isoBuilderHandle = isoAODFixBuilderPhoton()
+      elif self.type == "ForwardElectrons":
+        self._isoBuilderHandle = isoAODFixBuilderFwdElectron()  
+      elif self.type == "Muons":
+        self._isoBuilderHandle = isoAODFixBuilderMuon()  
+      else:
+        mlog.error("wrong object for IsolationBuilder")
+    except Exception:
+      mlog.error("could not get handle to IsolationBuilder")
+      import traceback
+      print traceback.format_exc()
+      return False
+         
+  def isoBuilderHandle(self):
+    return self._BuilderHandle
+
diff --git a/Reconstruction/RecoAlgs/IsolationAlgs/python/IsoGetter.py b/Reconstruction/RecoAlgs/IsolationAlgs/python/IsoGetter.py
index 20845da18a1..4ea998107f7 100644
--- a/Reconstruction/RecoAlgs/IsolationAlgs/python/IsoGetter.py
+++ b/Reconstruction/RecoAlgs/IsolationAlgs/python/IsoGetter.py
@@ -15,8 +15,7 @@ from AthenaCommon.GlobalFlags import globalflags
 isMC = not globalflags.DataSource()=='data'
 from IsolationCorrections.IsolationCorrectionsConf import CP__IsolationCorrectionTool as ICT
 IsoCorrectionTool = ToolFactory(ICT,
-                                name = "NewLeakageCorrTool",
-                                IsMC = isMC)
+                                name = "NewLeakageCorrTool")
 
 doPFlow = False
 PFlowObjectsInConeTool = None
@@ -73,7 +72,7 @@ def configureEDCorrection(tool):
   try:
     from AthenaCommon.AppMgr import ToolSvc
     from AthenaCommon.AlgSequence import AlgSequence
-    from EventShapeTools.EventDensityConfig import configEventDensityTool, EventDensityAlg
+    from EventShapeTools.EventDensityConfig import configEventDensityTool, EventDensityAthAlg
     from JetRec.JetRecStandard import jtm
     topSequence = AlgSequence()
     if not hasattr(topSequence,'EDtpIsoCentralAlg'):
@@ -85,7 +84,7 @@ def configureEDCorrection(tool):
                                     OutputLevel     = OutputLevel
                                     )
       ToolSvc += tccc
-      topSequence += EventDensityAlg("EDtpIsoCentralAlg", EventDensityTool = tccc)
+      topSequence += EventDensityAthAlg("EDtpIsoCentralAlg", EventDensityTool = tccc)
 
     if not hasattr(topSequence,'EDtpIsoForwardAlg'):
       tfcc = configEventDensityTool("EDtpIsoForwardTool", jtm.emget,
@@ -96,7 +95,18 @@ def configureEDCorrection(tool):
                                     OutputLevel     = OutputLevel
                                     )
       ToolSvc += tfcc
-      topSequence += EventDensityAlg("EDtpIsoForwardAlg", EventDensityTool = tfcc)
+      topSequence += EventDensityAthAlg("EDtpIsoForwardAlg", EventDensityTool = tfcc)
+
+    if not hasattr(topSequence,'EDtpIsoVeryForwardAlg'):
+      tvfcc = configEventDensityTool("EDtpIsoVeryForwardTool", jtm.emget,
+                                     radius          = 0.5,
+                                     AbsRapidityMin  = 2.5,
+                                     AbsRapidityMax  = 4.5,
+                                     OutputContainer = "TopoClusterIsoVeryForwardEventShape",
+                                     OutputLevel     = OutputLevel
+                                   )
+      ToolSvc += tvfcc
+      topSequence += EventDensityAthAlg("EDtpIsoVeryForwardAlg", EventDensityTool = tvfcc)
 
     if doPFlow:
       if not hasattr(topSequence,'EDpfIsoCentralAlg'):
@@ -108,7 +118,7 @@ def configureEDCorrection(tool):
                                       OutputLevel     = OutputLevel
                                       )
         ToolSvc += tcpf
-        topSequence += EventDensityAlg("EDpfIsoCentralAlg", EventDensityTool = tcpf)
+        topSequence += EventDensityAthAlg("EDpfIsoCentralAlg", EventDensityTool = tcpf)
 
       if not hasattr(topSequence,'EDpfIsoForwardAlg'):
         tfpf = configEventDensityTool("EDpfIsoForwardTool", jtm.empflowget,
@@ -119,7 +129,7 @@ def configureEDCorrection(tool):
                                       OutputLevel     = OutputLevel
                                       )
         ToolSvc += tfpf
-        topSequence += EventDensityAlg("EDpfIsoForwardAlg", EventDensityTool = tfpf)
+        topSequence += EventDensityAthAlg("EDpfIsoForwardAlg", EventDensityTool = tfpf)
 
       ## Try a neutral density
       if not hasattr(topSequence,'EDnpfIsoCentralAlg'):
@@ -131,7 +141,7 @@ def configureEDCorrection(tool):
                                        OutputLevel     = OutputLevel
                                        )
         ToolSvc += tcnpf
-        topSequence += EventDensityAlg("EDnpfIsoCentralAlg", EventDensityTool = tcnpf)
+        topSequence += EventDensityAthAlg("EDnpfIsoCentralAlg", EventDensityTool = tcnpf)
 
       if not hasattr(topSequence,'EDnpfIsoForwardAlg'):
         tfnpf = configEventDensityTool("EDnpfIsoForwardTool", jtm.emnpflowget,
@@ -142,7 +152,7 @@ def configureEDCorrection(tool):
                                        OutputLevel     = OutputLevel
                                        )
         ToolSvc += tfnpf
-        topSequence += EventDensityAlg("EDnpfIsoForwardAlg", EventDensityTool = tfnpf)
+        topSequence += EventDensityAthAlg("EDnpfIsoForwardAlg", EventDensityTool = tfnpf)
 
   except Exception:
     print '\nERROR: could not get handle to ED'
@@ -186,6 +196,12 @@ IsoTypes =  [
     isoPar.ptcone30,
     isoPar.ptcone20 ]
   ]
+IsoTypesFe =  [
+  [ isoPar.topoetcone20, 
+    isoPar.topoetcone30,
+    isoPar.topoetcone40 ]
+  ]
+
 # The Default corrections : have to follow exactly the order of the iso var.
 IsoCorEg = [
   [ isoPar.core57cells, isoPar.ptCorrection ],
@@ -198,6 +214,9 @@ IsoCorMu = [
   [ isoPar.coreCone, isoPar.pileupCorrection ],
   [ isoPar.coreTrackPtr ] #still hard-coded
   ]
+IsoCorFe = [
+  [ isoPar.coreCone, isoPar.pileupCorrection ] 
+  ]
 
 if doPFlow:
   IsoTypes.append(  
@@ -214,31 +233,34 @@ isoBuilder = AlgFactory(IsolationBuilder,
                         CaloTopoIsolationTool = CaloIsolationTool,
                         PFlowIsolationTool    = CaloIsolationTool,
                         TrackIsolationTool    = TrackIsolationTool, 
+                        FeIsoTypes            = [[]] if not rec.doEgamma() else IsoTypesFe,
+                        FeCorTypes            = IsoCorFe,
 			EgIsoTypes            = [[]] if not rec.doEgamma() else IsoTypes,
                         EgCorTypes            = IsoCorEg,
 			MuIsoTypes            = [[]] if not rec.doMuon() else IsoTypes,
                         MuCorTypes            = IsoCorMu,
+                        LeakageTool           = None,
                         OutputLevel           = 3)
 
 from RecExConfig.Configured import Configured
 class isoGetter ( Configured ) :
  
-    def configure(self):
-        mlog = logging.getLogger ('isoGetter.py::configure:')
-        mlog.info('entering')        
-        
-         # configure iso here:
-        try:
-            self._isoBuilderHandle = isoBuilder()
-        except Exception:
-            mlog.error("could not get handle to IsolationBuilder")
-            import traceback
-            print traceback.format_exc()
-            return False
-         
-        #print self._isoBuilderHandle
-        return True
- 
-    def isoBuilderHandle(self):
-        return self._BuilderHandle
+  def configure(self):
+    mlog = logging.getLogger ('isoGetter.py::configure:')
+    mlog.info('entering')        
+    
+    # configure iso here:
+    try:
+      self._isoBuilderHandle = isoBuilder()
+    except Exception:
+      mlog.error("could not get handle to IsolationBuilder")
+      import traceback
+      print traceback.format_exc()
+      return False
+    
+    #print self._isoBuilderHandle
+    return True
+
+  def isoBuilderHandle(self):
+    return self._BuilderHandle
 
diff --git a/Reconstruction/RecoAlgs/IsolationAlgs/share/IsoEventShapeOutputItemList_jobOptions.py b/Reconstruction/RecoAlgs/IsolationAlgs/share/IsoEventShapeOutputItemList_jobOptions.py
index 0bdaa4fe38b..2c9a7712249 100644
--- a/Reconstruction/RecoAlgs/IsolationAlgs/share/IsoEventShapeOutputItemList_jobOptions.py
+++ b/Reconstruction/RecoAlgs/IsolationAlgs/share/IsoEventShapeOutputItemList_jobOptions.py
@@ -10,6 +10,8 @@ if rec.doWriteAOD() or rec.doWriteESD():
         IsoAODESList += ["xAOD::EventShapeAuxInfo#TopoClusterIsoCentralEventShapeAux."]
         IsoAODESList += ["xAOD::EventShape#TopoClusterIsoForwardEventShape"]
         IsoAODESList += ["xAOD::EventShapeAuxInfo#TopoClusterIsoForwardEventShapeAux."]
+        IsoAODESList += ["xAOD::EventShape#TopoClusterIsoVeryForwardEventShape"]
+        IsoAODESList += ["xAOD::EventShapeAuxInfo#TopoClusterIsoVeryForwardEventShapeAux."]
         if recAlgs.doEFlow():
             IsoAODESList += ["xAOD::EventShape#ParticleFlowIsoCentralEventShape"]
             IsoAODESList += ["xAOD::EventShapeAuxInfo#ParticleFlowIsoCentralEventShapeAux."]
diff --git a/Reconstruction/RecoAlgs/IsolationAlgs/src/IsolationBuilder.cxx b/Reconstruction/RecoAlgs/IsolationAlgs/src/IsolationBuilder.cxx
index 16230240a46..b667191ab89 100644
--- a/Reconstruction/RecoAlgs/IsolationAlgs/src/IsolationBuilder.cxx
+++ b/Reconstruction/RecoAlgs/IsolationAlgs/src/IsolationBuilder.cxx
@@ -23,10 +23,13 @@
 #include "xAODEgamma/EgammaxAODHelpers.h"
 #include "xAODEgamma/EgammaContainer.h"
 #include "xAODEgamma/ElectronContainer.h"
+#include "xAODEgamma/ElectronAuxContainer.h"
 #include "xAODEgamma/Electron.h"
 #include "xAODEgamma/PhotonContainer.h"
+#include "xAODEgamma/PhotonAuxContainer.h"
 #include "xAODEgamma/Photon.h"
 #include "xAODMuon/MuonContainer.h"
+#include "xAODMuon/MuonAuxContainer.h"
 #include "xAODMuon/Muon.h"
 #include "CaloEvent/CaloCellContainer.h"
 
@@ -44,20 +47,33 @@ IsolationBuilder::IsolationBuilder( const std::string& name,
   //
   // Property declaration
   // 
-  declareProperty("ElectronCollectionContainerName", m_ElectronContainerName = "Electrons");
-  declareProperty("PhotonCollectionContainerName",   m_PhotonContainerName   = "Photons");
-  declareProperty("MuonCollectionContainerName",     m_MuonContainerName     = "Muons");
-  declareProperty("CellCollectionName",              m_cellsName             ="AllCalo", "Name of container which contain calo cells");
+  declareProperty("ElectronCollectionContainerName",    m_ElectronContainerName    = "Electrons");
+  declareProperty("PhotonCollectionContainerName",      m_PhotonContainerName      = "Photons");
+  declareProperty("MuonCollectionContainerName",        m_MuonContainerName        = "Muons");
+  declareProperty("FwdElectronCollectionContainerName", m_FwdElectronContainerName = "ForwardElectrons");
+  declareProperty("CellCollectionName",              m_cellsName             = "AllCalo", "Name of container which contain calo cells");
   declareProperty("CaloCellIsolationTool",           m_cellIsolationTool,    "Handle of the calo cell IsolationTool");
   declareProperty("CaloTopoIsolationTool",           m_topoIsolationTool,    "Handle of the calo topo IsolationTool");
   declareProperty("PFlowIsolationTool",              m_pflowIsolationTool,   "Handle of the pflow IsolationTool");
   declareProperty("TrackIsolationTool",              m_trackIsolationTool,   "Handle of the track IsolationTool");
   declareProperty("useBremAssoc",                    m_useBremAssoc          = true, "use track to track assoc after brem");
+  declareProperty("FeIsoTypes",                      m_feisoInts, "The isolation types to do for forward electron: vector of vector of enum type Iso::IsolationType, stored as float");
+  declareProperty("FeCorTypes",                      m_fecorInts, "The correction types to do for forward electron iso: vector of vector of enum type Iso::IsolationCalo/TrackCorrection, stored as float");
   declareProperty("EgIsoTypes",                      m_egisoInts, "The isolation types to do for egamma: vector of vector of enum type Iso::IsolationType, stored as float");
   declareProperty("EgCorTypes",                      m_egcorInts, "The correction types to do for egamma iso: vector of vector of enum type Iso::IsolationCalo/TrackCorrection, stored as float");
   declareProperty("MuIsoTypes",                      m_muisoInts, "The isolation types to do for Muons : vector of vector of enum type Iso::IsolationType, stored as float");
   declareProperty("MuCorTypes",                      m_mucorInts, "The correction types to do for Muon iso: vector of vector of enum type Iso::IsolationCalo/TrackCorrection, stored as float");
-  declareProperty("CustomConfgurationName",          m_customConfig          = "", "use a custom configuration");
+  declareProperty("CustomConfigurationName",         m_customConfig          = "", "use a custom configuration");
+  declareProperty("CustomConfigurationNameEl",       m_customConfigEl        = "", "use a custom configuration for electron");
+  declareProperty("CustomConfigurationNamePh",       m_customConfigPh        = "", "use a custom configuration for photon");
+  declareProperty("CustomConfigurationNameFwd",      m_customConfigFwd       = "", "use a custom configuration for forward electron");
+  declareProperty("CustomConfigurationNameMu",       m_customConfigMu        = "", "use a custom configuration for muon");
+  declareProperty("IsAODFix",                        m_isAODFix              = false);
+  declareProperty("LeakageTool",                     m_leakTool,                   "Handle of the leakage Tool");
+  declareProperty("IsolateEl",                       m_isolateEl             = true, "since egIsoTypes is common for el and ph, a new flag to control individually : electron");
+  declareProperty("IsolatePh",                       m_isolatePh             = true, "since egIsoTypes is common for el and ph, a new flag to control individually : photon");
+  declareProperty("AllTrackRemoval",                 m_allTrackRemoval       = true);
+  
 }
 
 IsolationBuilder::~IsolationBuilder() {}
@@ -69,18 +85,18 @@ StatusCode IsolationBuilder::initialize()
   std::set<xAOD::Iso::IsolationFlavour> runIsoType;
 
   // Build helpers for all required isolations
-  // For egamma
+  // For egamma (central)
   unsigned int nI = m_egisoInts.size();
   if (nI > 0 && m_egisoInts[0].size() > 0) 
     ATH_MSG_INFO("Will built egamma isolation types : ");
   else 
     nI = 0;
-  for (unsigned int ii = 0; ii < nI; ii++) {
+  for (unsigned int ii = 0; ii < nI; ii++) { // the various types : should be 4 by default (etcone, topoetcone, ptcone, neflowisol)
     std::vector<xAOD::Iso::IsolationType> isoTypes;
     std::vector<SG::AuxElement::Decorator<float>*> Deco;
     xAOD::Iso::IsolationFlavour isoFlav =
       xAOD::Iso::numIsolationFlavours;
-    for (unsigned int jj = 0; jj < m_egisoInts[ii].size(); jj++) {
+    for (unsigned int jj = 0; jj < m_egisoInts[ii].size(); jj++) { // the size of the cone : should be 3 by default
       int egIso = int(m_egisoInts[ii][jj]+0.5);
       xAOD::Iso::IsolationType isoType = static_cast<xAOD::Iso::IsolationType>(egIso);
       isoTypes.push_back(isoType);
@@ -126,25 +142,64 @@ StatusCode IsolationBuilder::initialize()
     if (runIsoType.find(isoFlav) == runIsoType.end()) runIsoType.insert(isoFlav);
   }
 
+  // For forward electrons
+  nI = m_feisoInts.size();
+  if (nI > 0 && m_feisoInts[0].size() > 0) 
+    ATH_MSG_INFO("Will built forward electron isolation types : ");
+  else 
+    nI = 0;
+  for (unsigned int ii = 0; ii < nI; ii++) { // here, should be at most three (no ptcone)
+    std::vector<xAOD::Iso::IsolationType> isoTypes;
+    std::vector<SG::AuxElement::Decorator<float>*> Deco;
+    xAOD::Iso::IsolationFlavour isoFlav =
+      xAOD::Iso::numIsolationFlavours;
+    for (unsigned int jj = 0; jj < m_feisoInts[ii].size(); jj++) {
+      int egIso = int(m_feisoInts[ii][jj]+0.5);
+      xAOD::Iso::IsolationType isoType = static_cast<xAOD::Iso::IsolationType>(egIso);
+      isoTypes.push_back(isoType);
+      isoFlav = xAOD::Iso::isolationFlavour(isoType);
+      std::string isoName = xAOD::Iso::toString(isoType);
+      if (m_customConfig != "") {
+	isoName += "_"; isoName += m_customConfig;
+	Deco.push_back(new SG::AuxElement::Decorator<float>(isoName));
+      }
+      ATH_MSG_INFO(isoName);
+    }
+    if (isoFlav == xAOD::Iso::etcone || isoFlav == xAOD::Iso::topoetcone || isoFlav == xAOD::Iso::neflowisol) {
+      CaloIsoHelp cisoH;
+      cisoH.isoTypes = isoTypes;
+      if (m_customConfig != "") cisoH.isoDeco  = Deco;
+      xAOD::CaloCorrection corrlist;
+      for (unsigned int jj = 0; jj < m_fecorInts[ii].size(); jj++) {
+        int egCor = int(m_fecorInts[ii][jj]+0.5);
+	corrlist.calobitset.set(static_cast<unsigned int>(egCor));
+      }
+      cisoH.CorrList = corrlist;
+      m_feCaloIso.insert(std::make_pair(isoFlav,cisoH));
+    } else 
+      ATH_MSG_WARNING("Isolation flavour " << xAOD::Iso::toString(isoFlav) << " does not exist ! Check your inputs");
+    if (runIsoType.find(isoFlav) == runIsoType.end()) runIsoType.insert(isoFlav);
+  }
+
   // For muons
   nI = m_muisoInts.size();
   if (nI > 0 && m_muisoInts[0].size() > 0) 
     ATH_MSG_INFO("Will built muon isolation types : ");
   else 
     nI = 0;
-  for (unsigned int ii = 0; ii < nI; ii++) {
+  for (unsigned int ii = 0; ii < nI; ii++) {   // the flavor (cell, topo, eflow, track
     std::vector<xAOD::Iso::IsolationType> isoTypes;
     std::vector<SG::AuxElement::Decorator<float>*> Deco;
     xAOD::Iso::IsolationFlavour isoFlav =
       xAOD::Iso::numIsolationFlavours;
-    for (unsigned int jj = 0; jj < m_muisoInts[ii].size(); jj++) {
+    for (unsigned int jj = 0; jj < m_muisoInts[ii].size(); jj++) { // the cone size
       int muIso = int(m_muisoInts[ii][jj]+0.5);
       xAOD::Iso::IsolationType isoType = static_cast<xAOD::Iso::IsolationType>(muIso);
       isoTypes.push_back(isoType);
       isoFlav = xAOD::Iso::isolationFlavour(isoType);
       std::string isoName = xAOD::Iso::toString(isoType);
-      if (m_customConfig != "") {
-	isoName += "_"; isoName += m_customConfig;
+      if (m_customConfigMu != "") {
+	isoName += "_"; isoName += m_customConfigMu;
 	Deco.push_back(new SG::AuxElement::Decorator<float>(isoName));
       }
       ATH_MSG_INFO(isoName);
@@ -152,8 +207,8 @@ StatusCode IsolationBuilder::initialize()
 	isoName = "ptvarcone";
 	int is = 100*xAOD::Iso::coneSize(isoType);
 	isoName += std::to_string(is);
-	if (m_customConfig != "") {
-	  isoName += "_"; isoName += m_customConfig;
+	if (m_customConfigMu != "") {
+	  isoName += "_"; isoName += m_customConfigMu;
 	  Deco.push_back(new SG::AuxElement::Decorator<float>(isoName));
 	}
 	ATH_MSG_INFO(isoName);
@@ -162,18 +217,29 @@ StatusCode IsolationBuilder::initialize()
     if (isoFlav == xAOD::Iso::etcone || isoFlav == xAOD::Iso::topoetcone || isoFlav == xAOD::Iso::neflowisol) {
       CaloIsoHelp cisoH;
       cisoH.isoTypes = isoTypes;
-      if (m_customConfig != "") cisoH.isoDeco  = Deco;
+      cisoH.coreCorisoDeco = nullptr;
+      if (m_customConfigMu != "") cisoH.isoDeco  = Deco;
       xAOD::CaloCorrection corrlist;
       for (unsigned int jj = 0; jj < m_mucorInts[ii].size(); jj++) {
 	int muCor = int(m_mucorInts[ii][jj]+0.5);
 	corrlist.calobitset.set(static_cast<unsigned int>(muCor));
+	xAOD::Iso::IsolationCaloCorrection isoCor = static_cast<xAOD::Iso::IsolationCaloCorrection>(muCor);
+	if (m_customConfigMu != "" && (isoCor == xAOD::Iso::coreCone || isoCor == xAOD::Iso::coreMuon) ) {
+	  std::string isoCorName = "";
+	  if (isoFlav == xAOD::Iso::topoetcone || isoFlav == xAOD::Iso::neflowisol)
+	    isoCorName = xAOD::Iso::toString(isoFlav);
+	  isoCorName += xAOD::Iso::toString(isoCor);
+	  isoCorName += xAOD::Iso::toString(xAOD::Iso::coreEnergy); isoCorName += "Correction"; // hard coded since we never store the core area in fact
+	  isoCorName += "_"; isoCorName += m_customConfigMu;
+	  cisoH.coreCorisoDeco = new SG::AuxElement::Decorator<float>(isoCorName);
+      }
       }
       cisoH.CorrList = corrlist;
       m_muCaloIso.insert(std::make_pair(isoFlav,cisoH));
     } else if (isoFlav == xAOD::Iso::ptcone) {
       TrackIsoHelp tisoH;
       tisoH.isoTypes = isoTypes;
-      if (m_customConfig != "") tisoH.isoDeco  = Deco;
+      if (m_customConfigMu != "") tisoH.isoDeco  = Deco;
       xAOD::TrackCorrection corrlist;
       corrlist.trackbitset.set(static_cast<unsigned int>(xAOD::Iso::coreTrackPtr));
       tisoH.CorrList = corrlist;
@@ -196,6 +262,14 @@ StatusCode IsolationBuilder::initialize()
   if (!m_trackIsolationTool.empty() && runIsoType.find(xAOD::Iso::IsolationFlavour::ptcone) != runIsoType.end())
     ATH_CHECK(m_trackIsolationTool.retrieve());
 
+
+  // Also can apply leakage correction if AODFix :
+  // be carefull ! There is a leakage in topoetcone, etcone, ... : do not run this !!!!!! (useless)
+  if (m_isAODFix && !m_leakTool.empty()) {
+    ATH_MSG_INFO("Will run leakage corrections for photons and electrons");
+    ATH_CHECK(m_leakTool.retrieve());
+  }
+
   return StatusCode::SUCCESS;
 }
 
@@ -218,6 +292,8 @@ StatusCode IsolationBuilder::finalize()
   }
   itc = m_muCaloIso.begin(), itcE = m_muCaloIso.end();
   for (; itc != itcE; itc++) {
+    if (itc->second.coreCorisoDeco != nullptr)      // it is only there for coreCone or coreMuon and if decoration
+      delete itc->second.coreCorisoDeco;            // only implemented for muon calo iso for the time being...
     std::vector<SG::AuxElement::Decorator<float>*> vec = itc->second.isoDeco;
     for (unsigned int ii = 0; ii < vec.size(); ii++)
       delete vec[ii];
@@ -244,25 +320,99 @@ StatusCode IsolationBuilder::execute()
     }
   }
 
+  // If AODFix, first deep copy
+  if (m_isAODFix) {
+    if (m_ElectronContainerName.size()) {
+      if (!evtStore()->tryRetrieve<xAOD::ElectronContainer>(m_ElectronContainerName)) {
+	if( deepCopy<xAOD::ElectronContainer,xAOD::ElectronAuxContainer>(m_ElectronContainerName).isFailure()) {
+	  ATH_MSG_FATAL( "Couldn't deep copy electrons" );
+	  return StatusCode::FAILURE;
+	}
+      }
+    }
+    if (m_FwdElectronContainerName.size()) {
+      if (!evtStore()->tryRetrieve<xAOD::ElectronContainer>(m_FwdElectronContainerName)) {
+	if( deepCopy<xAOD::ElectronContainer,xAOD::ElectronAuxContainer>(m_FwdElectronContainerName).isFailure()) {
+	  ATH_MSG_FATAL( "Couldn't deep copy forward electrons" );
+	  return StatusCode::FAILURE;
+	}
+      }
+    }
+    if (m_PhotonContainerName.size()) {
+      if (!evtStore()->tryRetrieve<xAOD::PhotonContainer>(m_PhotonContainerName)) {
+	if( deepCopy<xAOD::PhotonContainer,xAOD::PhotonAuxContainer>(m_PhotonContainerName).isFailure()) {
+	  ATH_MSG_FATAL( "Couldn't deep copy photons" );
+	  return StatusCode::FAILURE;
+	}
+      }
+    }
+    if (m_MuonContainerName.size()) {
+      if (!evtStore()->tryRetrieve<xAOD::MuonContainer>(m_MuonContainerName)) {
+	if( deepCopy<xAOD::MuonContainer,xAOD::MuonAuxContainer>(m_MuonContainerName).isFailure()) {
+	  ATH_MSG_FATAL( "Couldn't deep copy muons" );
+	  return StatusCode::FAILURE;
+	}
+      }
+    }
+  }
+
   // Compute isolations
+  /*
   if (m_customConfig == "") { 
     // Here, default configurations with standard names can be are computed
     if (m_egCaloIso.size() != 0 || m_egTrackIso.size() != 0) {
-      CHECK(IsolateEgamma("electron"));
-      CHECK(IsolateEgamma("photon"));
+      if (m_ElectronContainerName.size()) CHECK(IsolateEgamma("electron"));
+      if (m_PhotonContainerName.size())   CHECK(IsolateEgamma("photon"));
     }
+    if (m_feCaloIso.size() != 0)
+      if (m_FwdElectronContainerName.size()) CHECK(IsolateEgamma("fwdelectron"));
+      
     if (m_muCaloIso.size() != 0 || m_muTrackIso.size() != 0)
-      CHECK(IsolateMuon());
+      if (m_MuonContainerName.size()) CHECK(IsolateMuon());
   } else {
     // Here, custom configurations with custom names can be computed
     if (m_egCaloIso.size() != 0 || m_egTrackIso.size() != 0) {
-      CHECK(DecorateEgamma("electron"));
-      CHECK(DecorateEgamma("photon"));
+      if (m_ElectronContainerName.size()) CHECK(DecorateEgamma("electron"));
+      if (m_PhotonContainerName.size())   CHECK(DecorateEgamma("photon"));
     }
+    if (m_feCaloIso.size() != 0)
+      if (m_FwdElectronContainerName.size()) CHECK(DecorateEgamma("fwdelectron"));
+
     if (m_muCaloIso.size() != 0 || m_muTrackIso.size() != 0)
+      if (m_MuonContainerName.size()) CHECK(DecorateMuon());
+  }
+  */
+  if (m_egCaloIso.size() != 0 || m_egTrackIso.size() != 0) {
+    if (m_isolateEl && m_ElectronContainerName.size()) {
+      if (m_customConfigEl == "")
+	CHECK(IsolateEgamma("electron"));
+      else
+	CHECK(DecorateEgamma("electron"));
+    }
+    if (m_isolatePh && m_PhotonContainerName.size()) {
+      if (m_customConfigPh == "")
+	CHECK(IsolateEgamma("photon"));
+      else
+	CHECK(DecorateEgamma("photon"));
+    }
+  }
+  if (m_feCaloIso.size() != 0 && m_FwdElectronContainerName.size()) {
+    if (m_customConfigFwd == "")
+      CHECK(IsolateEgamma("fwdelectron"));
+    else
+      CHECK(DecorateEgamma("fwdelectron"));
+  }
+  if ( (m_muCaloIso.size() != 0 || m_muTrackIso.size() != 0) && m_MuonContainerName.size()) {
+    if (m_customConfigMu == "")
+      CHECK(IsolateMuon());
+    else
       CHECK(DecorateMuon());
   }
 
+  if (m_isAODFix && !m_leakTool.empty())
+    CHECK(runLeakage());
+
+  
   return StatusCode::SUCCESS;
 }
 
@@ -289,6 +439,16 @@ StatusCode IsolationBuilder::IsolateEgamma(std::string egType) {
       ATH_MSG_DEBUG("PhotonContainer " << m_PhotonContainerName << " not available");
       return StatusCode::SUCCESS;
     }
+  } else if (egType == "fwdelectron") {
+    if (evtStore()->contains<xAOD::ElectronContainer>(m_FwdElectronContainerName)) {
+      if (evtStore()->retrieve(egC,m_FwdElectronContainerName).isFailure()) {
+	ATH_MSG_FATAL("Cannot retrieve forward electron container " << m_FwdElectronContainerName);
+	return StatusCode::FAILURE;
+      }
+    } else {
+      ATH_MSG_DEBUG("Forward ElectronContainer " << m_FwdElectronContainerName << " not available");
+      return StatusCode::SUCCESS;
+    }
   } else {
     ATH_MSG_WARNING("Unknown egamma type " << egType);
     return StatusCode::SUCCESS;
@@ -302,6 +462,10 @@ StatusCode IsolationBuilder::IsolateEgamma(std::string egType) {
     // Calo Isolation types
     //
     std::map<xAOD::Iso::IsolationFlavour,CaloIsoHelp>::iterator itc = m_egCaloIso.begin(), itcE = m_egCaloIso.end();
+    if (egType == "fwdelectron") {
+      itc  = m_feCaloIso.begin();
+      itcE = m_feCaloIso.end();
+    }
     for (; itc != itcE; itc++) {
       CaloIsoHelp isoH = itc->second;
       xAOD::Iso::IsolationFlavour flav = itc->first;
@@ -324,13 +488,30 @@ StatusCode IsolationBuilder::IsolateEgamma(std::string egType) {
       } else
 	ATH_MSG_WARNING("Call to CaloIsolationTool failed for flavour " << xAOD::Iso::toString(flav));
     }
+    if (egType == "fwdelectron")
+      return StatusCode::SUCCESS;
     // 
     // Track Isolation types
     //
     std::map<xAOD::Iso::IsolationFlavour,TrackIsoHelp>::iterator itt = m_egTrackIso.begin(), ittE = m_egTrackIso.end();
     for (; itt != ittE; itt++) {
       TrackIsoHelp isoH = itt->second;
-      const std::set<const xAOD::TrackParticle*> tracksToExclude = xAOD::EgammaHelpers::getTrackParticles(eg, m_useBremAssoc);
+      //const std::set<const xAOD::TrackParticle*> tracksToExclude = xAOD::EgammaHelpers::getTrackParticles(eg, m_useBremAssoc);
+      std::set<const xAOD::TrackParticle*> tracksToExclude;
+      if (eg->type() == xAOD::Type::Electron)
+	tracksToExclude = xAOD::EgammaHelpers::getTrackParticles(eg, m_useBremAssoc);
+      else {
+	if (m_allTrackRemoval) //New (from ??/??/16) : now this gives all tracks
+	  tracksToExclude = xAOD::EgammaHelpers::getTrackParticles(eg, m_useBremAssoc);
+	else { // this is just to be able to have the 2015+2016 default case (only tracks from first vertex)
+	  const xAOD::Photon *gam = dynamic_cast<const xAOD::Photon *>(eg);
+	  if (gam->nVertices() > 0) {
+	    const xAOD::Vertex *phvtx = gam->vertex(0);
+	    for (unsigned int itk = 0; itk < phvtx->nTrackParticles(); itk++)
+	      tracksToExclude.insert(m_useBremAssoc ? xAOD::EgammaHelpers::getOriginalTrackParticleFromGSF(phvtx->trackParticle(itk)) : phvtx->trackParticle(itk));
+	  }
+	}
+      }
       xAOD::Vertex *vx = 0;
       bool bsc = m_trackIsolationTool->decorateParticle(*eg, isoH.isoTypes, isoH.CorrList, vx, &tracksToExclude);
       if (bsc) {
@@ -450,6 +631,16 @@ StatusCode IsolationBuilder::DecorateEgamma(std::string egType) {
     } else {
       ATH_MSG_DEBUG("PhotonContainer " << m_PhotonContainerName << " not available");
       return StatusCode::SUCCESS;
+    } 
+  } else if (egType == "fwdelectron") {
+    if (evtStore()->contains<xAOD::ElectronContainer>(m_FwdElectronContainerName)) {
+      if (evtStore()->retrieve(egC,m_FwdElectronContainerName).isFailure()) {
+	ATH_MSG_FATAL("Cannot retrieve forward electron container " << m_FwdElectronContainerName);
+	return StatusCode::FAILURE;
+      }
+    } else {
+      ATH_MSG_DEBUG("Forward ElectronContainer " << m_FwdElectronContainerName << " not available");
+      return StatusCode::SUCCESS;
     }
   } else {
     ATH_MSG_WARNING("Unknown egamma type " << egType);
@@ -464,6 +655,10 @@ StatusCode IsolationBuilder::DecorateEgamma(std::string egType) {
     // Calo Isolation types
     //
     std::map<xAOD::Iso::IsolationFlavour,CaloIsoHelp>::iterator itc = m_egCaloIso.begin(), itcE = m_egCaloIso.end();
+    if (egType == "fwdelectron") {
+      itc  = m_feCaloIso.begin();
+      itcE = m_feCaloIso.end();
+    }
     for (; itc != itcE; itc++) {
       xAOD::CaloIsolation CaloIsoResult;
       CaloIsoHelp isoH = itc->second;
@@ -484,6 +679,8 @@ StatusCode IsolationBuilder::DecorateEgamma(std::string egType) {
       } else
 	ATH_MSG_WARNING("Call to CaloIsolationTool failed for custom flavour " << xAOD::Iso::toString(flav));
     }
+    if (egType == "fwdelectron")
+      return StatusCode::SUCCESS;
     // 
     // Track Isolation types
     //
@@ -548,6 +745,23 @@ StatusCode IsolationBuilder::DecorateMuon() {
 	  ATH_MSG_DEBUG("custom Iso " << xAOD::Iso::toString(isoH.isoTypes[i]) << " = " << iso/1e3);
 	  (*isoH.isoDeco[i])(*mu) = iso;
 	}
+	// not that nice. I expect a single core correction (e.g. for topoetcone, not coreMuon and coreCone together...)
+	xAOD::Iso::IsolationCaloCorrection icc = xAOD::Iso::coreCone;
+	bool ison = false;
+	if (CaloIsoResult.corrlist.calobitset.test(static_cast<unsigned int>(xAOD::Iso::coreMuon))) {
+	  icc  = xAOD::Iso::coreMuon;
+	  ison = true;
+	} else if (CaloIsoResult.corrlist.calobitset.test(static_cast<unsigned int>(xAOD::Iso::coreCone)))
+	  ison = true;
+	if (ison) {
+	  if (CaloIsoResult.coreCorrections.find(icc) != CaloIsoResult.coreCorrections.end()) {
+	    if (CaloIsoResult.coreCorrections[icc].find(xAOD::Iso::coreEnergy) != CaloIsoResult.coreCorrections[icc].end())
+	      (*isoH.coreCorisoDeco)(*mu) = CaloIsoResult.coreCorrections[icc][xAOD::Iso::coreEnergy];
+	    else
+	      ATH_MSG_WARNING("Cannot find the core energy correction for custom flavour " << xAOD::Iso::toString(flav));
+	  } else
+	    ATH_MSG_WARNING("Cannot find the core correction for custom flavour " << xAOD::Iso::toString(flav));
+	}
       } else
 	ATH_MSG_WARNING("Call to CaloIsolationTool failed for custom flavour " << xAOD::Iso::toString(flav));
     }
@@ -573,4 +787,94 @@ StatusCode IsolationBuilder::DecorateMuon() {
   }
 
   return StatusCode::SUCCESS;
-} 
+}
+
+StatusCode IsolationBuilder::runLeakage() {
+  
+  // Retrieve data
+  if (m_PhotonContainerName.size()) {
+    xAOD::PhotonContainer* photons = evtStore()->retrieve< xAOD::PhotonContainer >(m_PhotonContainerName);
+    if( !photons ) {
+      ATH_MSG_ERROR("Couldn't retrieve photon container with key: " << m_PhotonContainerName);
+      return StatusCode::FAILURE;
+    }
+    for (auto ph : *photons)           
+      m_leakTool->applyCorrection(*ph);
+  }
+    
+  if (m_ElectronContainerName.size()) {
+    xAOD::ElectronContainer* electrons = evtStore()->retrieve< xAOD::ElectronContainer >(m_ElectronContainerName);
+    if( !electrons ) {
+      ATH_MSG_ERROR("Couldn't retrieve electron container with key: " << m_ElectronContainerName);
+      return StatusCode::FAILURE;
+    }
+    for (auto el : *electrons)
+      m_leakTool->applyCorrection(*el);
+  }
+
+  return StatusCode::SUCCESS;
+}
+
+
+template< class CONTAINER, class AUXSTORE > StatusCode IsolationBuilder::deepCopy( const std::string& key ) const {
+    
+  // Let the user know what's happening:
+  ATH_MSG_VERBOSE( "Running deepCopy on container: " << key );
+  
+  // Decide which implementation to call:
+  if( evtStore()->template contains< AUXSTORE >( key + "Aux." ) ) {
+    if( deepCopyImp< CONTAINER, AUXSTORE >( key ).isFailure() ) {
+      ATH_MSG_FATAL( "Couldn't call deepCopyImp with concrete "
+		     "auxiliary store" );
+      return StatusCode::FAILURE;
+    }
+  } else if( evtStore()->template contains< xAOD::AuxContainerBase >( key +
+								      "Aux." ) ) {
+    if( deepCopyImp< CONTAINER,
+	xAOD::AuxContainerBase >( key ).isFailure() ) {
+      ATH_MSG_FATAL( "Couldn't call deepCopyImp with generic "
+		     "auxiliary store" );
+      return StatusCode::FAILURE;
+    }
+  } else {
+    ATH_MSG_FATAL( "Couldn't discover auxiliary store type for container \""
+		   << key << "\"" );
+    return StatusCode::FAILURE;
+  }
+  
+  // Return gracefully:
+  return StatusCode::SUCCESS;
+}
+  
+template< class CONTAINER, class AUXSTORE > StatusCode IsolationBuilder::deepCopyImp( const std::string& key ) const {
+  
+  // Retrieve the const container:
+  ATH_MSG_VERBOSE( "Will retrieve " << key);
+  const CONTAINER* c = 0;
+  ATH_CHECK( evtStore()->retrieve( c, key ) );
+  
+  // Create the new container:
+  ATH_MSG_VERBOSE( "Will create new containers" );
+  CONTAINER* copy = new CONTAINER();
+  AUXSTORE* copyAux = new AUXSTORE();
+  copy->setStore( copyAux );
+  
+  // Do the deep copy:
+  ATH_MSG_VERBOSE( "Will copy the object" );
+  for( auto oldObj : *c ) {
+    ATH_MSG_VERBOSE( "Now working on object " << oldObj);
+    auto newObj = new typename CONTAINER::base_value_type();
+    copy->push_back( newObj );
+    *newObj = *oldObj;
+  }
+  
+  // Do the overwrite:
+  ATH_MSG_VERBOSE( "Will overwrite the container" );
+  ATH_CHECK( evtStore()->overwrite( copy, key, true, false ) );
+  ATH_MSG_VERBOSE( "Will overwrite the aux container" );
+  ATH_CHECK( evtStore()->overwrite( copyAux, key + "Aux.", true, false ) );
+  ATH_MSG_VERBOSE( "Done" );
+  
+  // Return gracefully:
+  return StatusCode::SUCCESS;
+}
diff --git a/Reconstruction/RecoAlgs/IsolationAlgs/src/IsolationBuilder.h b/Reconstruction/RecoAlgs/IsolationAlgs/src/IsolationBuilder.h
index 981943672a7..7c603578c54 100644
--- a/Reconstruction/RecoAlgs/IsolationAlgs/src/IsolationBuilder.h
+++ b/Reconstruction/RecoAlgs/IsolationAlgs/src/IsolationBuilder.h
@@ -22,6 +22,8 @@
 #include "xAODPrimitives/IsolationFlavour.h"
 #include "RecoToolInterfaces/IsolationCommon.h"
 
+#include "IsolationCorrections/IIsolationCorrectionTool.h"
+
 namespace xAOD {
   class INeutralEFlowIsolationTool;
   class ICaloTopoClusterIsolationTool;
@@ -70,6 +72,7 @@ class IsolationBuilder
 
   /// Containers
   std::string m_ElectronContainerName;
+  std::string m_FwdElectronContainerName;
   std::string m_PhotonContainerName;
   std::string m_MuonContainerName;
 
@@ -92,28 +95,46 @@ class IsolationBuilder
   /** @brief Isolation types (for the alg. properties, only vector<vector<double>> available */
   std::vector<std::vector<double> > m_egisoInts, m_muisoInts;
   std::vector<std::vector<double> > m_egcorInts, m_mucorInts;
+  std::vector<std::vector<double> > m_fecorInts, m_feisoInts;
 
   struct CaloIsoHelp {
     std::vector<SG::AuxElement::Decorator<float>*> isoDeco;
+    std::vector<SG::AuxElement::Decorator<float>*> nonecoreCorisoDeco;
+    SG::AuxElement::Decorator<float>* coreCorisoDeco;
     std::vector<xAOD::Iso::IsolationType> isoTypes;
     xAOD::CaloCorrection CorrList;
   };
-  std::map<xAOD::Iso::IsolationFlavour,CaloIsoHelp> m_egCaloIso, m_muCaloIso;
+  std::map<xAOD::Iso::IsolationFlavour,CaloIsoHelp> m_egCaloIso, m_feCaloIso, m_muCaloIso;
 
   struct TrackIsoHelp {
     std::vector<SG::AuxElement::Decorator<float>*> isoDeco;
+    SG::AuxElement::Decorator<float>* coreCorisoDeco;
     std::vector<xAOD::Iso::IsolationType> isoTypes;
     xAOD::TrackCorrection CorrList;
   };
   std::map<xAOD::Iso::IsolationFlavour,TrackIsoHelp> m_egTrackIso, m_muTrackIso;
 
+  // individual flags to run or not computation on electron or photon
+  bool m_isolateEl, m_isolatePh;
+  
   // Compute the isolation variables
   StatusCode IsolateEgamma(std::string egType);
   StatusCode IsolateMuon();
 
   std::string m_customConfig;
+  std::string m_customConfigEl, m_customConfigPh, m_customConfigFwd, m_customConfigMu; // for the time being, only mu vs eg, no separation in eg
   StatusCode DecorateEgamma(std::string egType);
   StatusCode DecorateMuon();
+
+  // For an AODFix
+  bool m_isAODFix;
+  bool m_allTrackRemoval;
+  ToolHandle<CP::IIsolationCorrectionTool> m_leakTool;
+  StatusCode runLeakage();
+  // From Attila, for a deep copy
+  template< class CONTAINER, class AUXSTORE > StatusCode deepCopy( const std::string& key ) const;
+  template< class CONTAINER, class AUXSTORE > StatusCode deepCopyImp( const std::string& key ) const;
+
   
 }; 
 
-- 
GitLab