diff --git a/benchmarks/tracking/central_electrons.sh b/benchmarks/tracking/central_electrons.sh
index d3c3cace5be61bd9cfb6c682579357029992fb64..4f7eab126839bbcd555a81004628b6fd5edb2ffd 100644
--- a/benchmarks/tracking/central_electrons.sh
+++ b/benchmarks/tracking/central_electrons.sh
@@ -1,12 +1,17 @@
 #!/bin/bash
 
 function print_the_help {
-  echo "USAGE: ${0} [--rec-only]  "
+  echo "USAGE: ${0} [-f,--force] [-b,--background] [--rec-only] [--ana-only] "
   echo "  OPTIONS: "
+  echo "            --force      force rerun even if files exist"
+  echo "            --background include overlay of background"
   echo "            --rec-only   only run gaudi reconstruction part"
+  echo "            --ana-only   only run analysis part"
   exit 
 }
 
+BKG=
+FORCE=
 REC_ONLY=
 ANALYSIS_ONLY=
 
@@ -20,6 +25,14 @@ do
       shift # past argument
       print_the_help
       ;;
+    -b|--background)
+      shift # past argument
+      BKG=1
+      ;;
+    -f|--force)
+      shift # past argument
+      FORCE=1
+      ;;
     --rec-only)
       REC_ONLY=1
       shift # past value
@@ -61,6 +74,7 @@ export JUGGLER_FILE_NAME_TAG="central_electrons"
 export JUGGLER_GEN_FILE="${JUGGLER_FILE_NAME_TAG}.hepmc"
 
 export JUGGLER_SIM_FILE="sim_${JUGGLER_FILE_NAME_TAG}.edm4hep.root"
+export JUGGLER_BKG_FILE="bkg_${JUGGLER_FILE_NAME_TAG}.edm4hep.root"
 export JUGGLER_REC_FILE="rec_${JUGGLER_FILE_NAME_TAG}.root"
 
 echo "JUGGLER_N_EVENTS = ${JUGGLER_N_EVENTS}"
@@ -71,34 +85,64 @@ if [[ -z "${REC_ONLY}" && -z "${ANALYSIS_ONLY}" ]] ;
 then
 
   ## generate the input events
-  root -b -q "benchmarks/tracking/scripts/gen_central_electrons.cxx+(${JUGGLER_N_EVENTS}, \"${JUGGLER_FILE_NAME_TAG}.hepmc\")"
-  if [[ "$?" -ne "0" ]] ; then
-    echo "ERROR running script"
-    exit 1
+  if [ -n "${FORCE}" -o ! -f ${JUGGLER_FILE_NAME_TAG}.hepmc ] ; then
+    root -b -q "benchmarks/tracking/scripts/gen_central_electrons.cxx+(${JUGGLER_N_EVENTS}, \"${JUGGLER_FILE_NAME_TAG}.hepmc\")"
+    if [[ "$?" -ne "0" ]] ; then
+      echo "ERROR running script"
+      exit 1
+    fi
   fi
 
-  echo "Running geant4 simulation"
   ## run geant4 simulations
-  ddsim --runType batch \
-    --part.minimalKineticEnergy 1000*GeV  \
-    --filter.tracker edep0 \
-    -v WARNING \
-    --numberOfEvents ${JUGGLER_N_EVENTS} \
-    --compactFile ${DETECTOR_PATH}/${DETECTOR_CONFIG}.xml \
-    --inputFiles  ${JUGGLER_FILE_NAME_TAG}.hepmc \
-    --outputFile  ${JUGGLER_SIM_FILE}
-  if [[ "$?" -ne "0" ]] ; then
-    echo "ERROR running script"
-    exit 1
+  if [ -n "${FORCE}" -o ! -f ${JUGGLER_SIM_FILE} ] ; then
+    echo "Running geant4 simulation"
+    ddsim --runType batch \
+      --part.minimalKineticEnergy 1000*GeV  \
+      --filter.tracker edep0 \
+      -v WARNING \
+      --numberOfEvents ${JUGGLER_N_EVENTS} \
+      --compactFile ${DETECTOR_PATH}/${DETECTOR_CONFIG}.xml \
+      --inputFiles  ${JUGGLER_FILE_NAME_TAG}.hepmc \
+      --outputFile  ${JUGGLER_SIM_FILE}
+    if [[ "$?" -ne "0" ]] ; then
+      echo "ERROR running script"
+      exit 1
+    fi
+    rootls -t ${JUGGLER_SIM_FILE}
   fi
-fi
 
-rootls -t ${JUGGLER_SIM_FILE}
+  ## run geant4 simulations (background)
+  if [[ -n "${BKG}" ]] ; then
+    if [ -n "${FORCE}" -o ! -f ${JUGGLER_BKG_FILE} ] ; then
+      echo "Running geant4 simulation (background)"
+      ddsim --runType batch \
+        --part.minimalKineticEnergy 1000*GeV  \
+        --filter.tracker edep0 \
+        -v WARNING \
+        --numberOfEvents ${JUGGLER_N_EVENTS} \
+        --compactFile ${DETECTOR_PATH}/${DETECTOR_CONFIG}.xml \
+        --enableGun --gun.particle e- \
+        --gun.momentumMin 0.0*GeV --gun.momentumMax 10.0*GeV \
+        --gun.thetaMin 45.0*deg --gun.thetaMax 135.0*deg \
+        --gun.distribution "cos(theta)" \
+        --outputFile  ${JUGGLER_BKG_FILE}
+      if [[ "$?" -ne "0" ]] ; then
+        echo "ERROR running script"
+        exit 1
+      fi
+      rootls -t ${JUGGLER_BKG_FILE}
+    fi
+  fi
+fi
 
 if [[ -z "${ANALYSIS_ONLY}" ]] ;
 then
   # Need to figure out how to pass file name to juggler from the commandline
-  gaudirun.py benchmarks/tracking/options/track_reconstruction.py
+  if [[ -n "${BKG}" ]] ; then
+    gaudirun.py benchmarks/tracking/options/track_reconstruction_with_overlay.py
+  else
+    gaudirun.py benchmarks/tracking/options/track_reconstruction.py
+  fi
   if [[ "$?" -ne "0" ]] ; then
     echo "ERROR running juggler"
     exit 1
diff --git a/benchmarks/tracking/options/track_reconstruction_with_overlay.py b/benchmarks/tracking/options/track_reconstruction_with_overlay.py
new file mode 100644
index 0000000000000000000000000000000000000000..716360345eda4dd5feb996d00d1e8798ab089ebd
--- /dev/null
+++ b/benchmarks/tracking/options/track_reconstruction_with_overlay.py
@@ -0,0 +1,219 @@
+from Gaudi.Configuration import *
+
+from Configurables import ApplicationMgr, EICDataSvc, PodioOutput, GeoSvc
+from GaudiKernel import SystemOfUnits as units
+
+detector_name = "topside"
+if "DETECTOR_CONFIG" in os.environ :
+    detector_name = str(os.environ["DETECTOR_CONFIG"])
+
+detector_path = ""
+if "DETECTOR_PATH" in os.environ :
+    detector_path = str(os.environ["DETECTOR_PATH"])
+
+detector_version = 'default'
+if "DETECTOR_VERSION" in os.environ:
+    env_version = str(os.environ["DETECTOR_VERSION"])
+
+# todo add checks
+input_sim_file  = str(os.environ["JUGGLER_SIM_FILE"])
+input_bkg_file  = str(os.environ["JUGGLER_BKG_FILE"])
+output_rec_file = str(os.environ["JUGGLER_REC_FILE"])
+n_events = str(os.environ["JUGGLER_N_EVENTS"])
+
+## note: old version of material map is called material-maps.XXX, new version is materials-map.XXX
+##       these names are somewhat inconsistent, and should probably all be renamed to 'material-map.XXX' 
+##       FIXME
+geo_service = GeoSvc("GeoSvc",
+                     detectors=["{}/{}.xml".format(detector_path,detector_name)],
+                     materials="calibrations/materials-map.cbor",
+                     OutputLevel=WARNING)
+
+podioevent   = EICDataSvc("EventDataSvc", inputs=[input_sim_file], OutputLevel=WARNING)
+
+from Configurables import PodioInput
+
+from Configurables import Jug__Digi__SimTrackerHitsCollector as SimTrackerHitsCollector
+from Configurables import Jug__Digi__SiliconTrackerDigi as TrackerDigi
+
+from Configurables import Jug__Reco__TrackerHitReconstruction as TrackerHitReconstruction
+from Configurables import Jug__Reco__TrackingHitsCollector2 as TrackingHitsCollector
+
+from Configurables import Jug__Reco__TrackerSourceLinker as TrackerSourceLinker
+from Configurables import Jug__Reco__TrackParamTruthInit as TrackParamTruthInit
+from Configurables import Jug__Reco__TrackParamClusterInit as TrackParamClusterInit
+from Configurables import Jug__Reco__TrackParamVertexClusterInit as TrackParamVertexClusterInit
+
+from Configurables import Jug__Reco__TrackFindingAlgorithm as TrackFindingAlgorithm
+from Configurables import Jug__Reco__ParticlesFromTrackFit as ParticlesFromTrackFit
+from Configurables import Jug__Reco__TrackProjector as TrackProjector
+
+from Configurables import Jug__Reco__SimpleClustering as SimpleClustering
+
+algorithms = [ ]
+
+tracker_endcap_collections = [
+    'TrackerEndcapHits1',
+    'TrackerEndcapHits2',
+    'TrackerEndcapHits3',
+    'TrackerEndcapHits4',
+    'TrackerEndcapHits5',
+    'TrackerEndcapHits6'
+]
+tracker_barrel_collections = [
+    'TrackerBarrelHits'
+]
+vertex_barrel_collections = [
+    'VertexBarrelHits'
+]
+gem_endcap_collections = [
+    'GEMTrackerEndcapHits1',
+    'GEMTrackerEndcapHits2',
+    'GEMTrackerEndcapHits3'
+]
+vertex_endcap_collections = [
+    'VertexEndcapHits'
+]
+mpgd_barrel_collections = [
+    'MPGDTrackerBarrelHits1',
+    'MPGDTrackerBarrelHits2'
+]
+
+input_collections = ['MCParticles'] + tracker_endcap_collections + tracker_barrel_collections + vertex_barrel_collections + gem_endcap_collections
+
+input_collections += mpgd_barrel_collections
+
+podioinput = PodioInput("PodioReader",
+                        collections=input_collections)
+algorithms.append( podioinput )
+
+# use the pileuptool to specify the number of pileup
+from Configurables import ConstPileUp
+pileuptool = ConstPileUp("MyPileupTool", numPileUpEvents=10)
+
+# edm data from simulation: hits and positioned hits
+from Configurables import PileupHitMergeTool_edm4hep__SimTrackerHitCollection_ as PileupTrackHitMergeTool
+trackhitsmergetool = PileupTrackHitMergeTool("MyTrackHitMergeTool")
+# branchnames for the pileup
+trackhitsmergetool.pileupHitsBranch = "hits"
+# branchnames for the signal
+trackhitsmergetool.signalHits = "hits"
+# branchnames for the output
+trackhitsmergetool.mergedHits = "overlaidTrackHits"
+
+# algorithm for the overlay
+from Configurables import PileupOverlayAlg
+overlay = PileupOverlayAlg()
+overlay.pileupFilenames = [input_bkg_file]
+overlay.randomizePileup = False
+overlay.mergeTools = ["PileupTrackHitMergeTool/MyTrackHitMergeTool"]
+overlay.PileUpTool = pileuptool
+overlay.OutputLevel = DEBUG
+
+trk_b_coll = SimTrackerHitsCollector("trk_b_coll",
+        inputSimTrackerHits = tracker_barrel_collections,
+        outputSimTrackerHits = "TrackerBarrelAllHits")
+algorithms.append( trk_b_coll )
+trk_b_digi = TrackerDigi("trk_b_digi",
+        inputHitCollection = trk_b_coll.outputSimTrackerHits,
+        outputHitCollection = "TrackerBarrelRawHits",
+        timeResolution=8)
+algorithms.append( trk_b_digi )
+
+vtx_b_coll = SimTrackerHitsCollector("vtx_b_coll",
+        inputSimTrackerHits = vertex_barrel_collections,
+        outputSimTrackerHits = "VertexBarrelAllHits")
+algorithms.append( vtx_b_coll )
+vtx_b_digi = TrackerDigi("vtx_b_digi",
+        inputHitCollection = vtx_b_coll.outputSimTrackerHits,
+        outputHitCollection = "VertexBarrelRawHits",
+        timeResolution=8)
+algorithms.append( vtx_b_digi )
+
+mm_b_coll = SimTrackerHitsCollector("mm_b_coll",
+        inputSimTrackerHits = mpgd_barrel_collections,
+        outputSimTrackerHits = "MPGDTrackerBarrelAllHits")
+algorithms.append( mm_b_coll )
+mm_b_digi = TrackerDigi("mm_b_digi",
+        inputHitCollection = mm_b_coll.outputSimTrackerHits,
+        outputHitCollection = "MPGDTrackerBarrelRawHits",
+        timeResolution=8)
+algorithms.append( mm_b_digi )
+
+# Tracker and vertex reconstruction
+trk_b_reco = TrackerHitReconstruction("trk_b_reco",
+        inputHitCollection = trk_b_digi.outputHitCollection,
+        outputHitCollection="TrackerBarrelRecHits")
+algorithms.append( trk_b_reco )
+
+vtx_b_reco = TrackerHitReconstruction("vtx_b_reco",
+        inputHitCollection = vtx_b_digi.outputHitCollection,
+        outputHitCollection="VertexBarrelRecHits")
+algorithms.append( vtx_b_reco )
+
+mm_b_reco = TrackerHitReconstruction("mm_b_reco",
+        inputHitCollection = mm_b_digi.outputHitCollection,
+        outputHitCollection="MPGDTrackerBarrelRecHits")
+algorithms.append( mm_b_reco )
+
+input_tracking_hits = [
+    str(trk_b_reco.outputHitCollection),
+    str(vtx_b_reco.outputHitCollection),
+]
+input_tracking_hits.append(str(mm_b_reco.outputHitCollection))
+
+trk_hit_col = TrackingHitsCollector("trk_hit_col",
+        inputTrackingHits=input_tracking_hits,
+        trackingHits="trackingHits")
+algorithms.append( trk_hit_col )
+
+# Hit Source linker 
+sourcelinker = TrackerSourceLinker("sourcelinker",
+        inputHitCollection=trk_hit_col.trackingHits,
+        outputSourceLinks="TrackSourceLinks",
+        outputMeasurements="TrackMeasurements")
+algorithms.append( sourcelinker )
+
+## Track param init
+truth_trk_init = TrackParamTruthInit("truth_trk_init",
+        inputMCParticles="MCParticles",
+        outputInitialTrackParameters="InitTrackParams")
+algorithms.append( truth_trk_init )
+
+# Tracking algorithms
+trk_find_alg = TrackFindingAlgorithm("trk_find_alg",
+        inputSourceLinks = sourcelinker.outputSourceLinks,
+        inputMeasurements = sourcelinker.outputMeasurements,
+        inputInitialTrackParameters = "InitTrackParams",
+        outputTrajectories = "trajectories")
+algorithms.append( trk_find_alg )
+
+parts_from_fit = ParticlesFromTrackFit("parts_from_fit",
+        inputTrajectories = "trajectories",
+        outputParticles = "ReconstructedParticles",
+        outputTrackParameters = "outputTrackParameters")
+algorithms.append( parts_from_fit )
+
+trk_proj = TrackProjector("trajs_from_fit",
+        inputTrajectories = trk_find_alg.outputTrajectories,
+        outputTrackSegments = "outputTrackSegments")
+algorithms.append( trk_proj )
+
+out = PodioOutput("out", filename=output_rec_file)
+out.outputCommands = ["keep *", 
+        "drop BarrelTrackSourceLinks", 
+        "drop InitTrackParams",
+        "drop trajectories",
+        "drop outputSourceLinks",
+        "drop outputInitialTrackParameters",
+        "keep MCParticles"
+        ]
+algorithms.append(out)
+
+ApplicationMgr(
+    TopAlg = algorithms,
+    EvtSel = 'NONE',
+    EvtMax   = n_events,
+    ExtSvc = [podioevent,geo_service],
+    OutputLevel=WARNING
+ )