diff --git a/MuonSpectrometer/MuonDetDescr/MuonReadoutGeometry/MuonReadoutGeometry/GlobalUtilities.h b/MuonSpectrometer/MuonDetDescr/MuonReadoutGeometry/MuonReadoutGeometry/GlobalUtilities.h
index d1e3d31915b371bce88c8ced6bd58544ca325327..a34813ffff3cca750726f3dc17352dafeeced350 100644
--- a/MuonSpectrometer/MuonDetDescr/MuonReadoutGeometry/MuonReadoutGeometry/GlobalUtilities.h
+++ b/MuonSpectrometer/MuonDetDescr/MuonReadoutGeometry/MuonReadoutGeometry/GlobalUtilities.h
@@ -48,6 +48,18 @@ namespace MuonGM {
         const AmgVector(N) AminusB = posA - posB;
         return (AminusB.dot(dirB) - AminusB.dot(dirA) * dirDots) / divisor;
     }
+    /// Intersects a line parametrized as A + lambda * B with the (N-1) dimensional
+    /// hyperplane that's given in the Jordan normal form <P, N> - C = 0 
+    template <int N>
+    std::optional<double> intersect(const AmgVector(N)& pos, const AmgVector(N)& dir, 
+                                    const AmgVector(N)& planeNorm, const double offset) {
+        ///  <P, N> - C = 0
+        /// --> <A + lambda *B , N> - C = 0
+        /// --> lambda = (C - <A,N> ) / <N, B>
+        const double normDot = planeNorm.dot(dir); 
+        if (std::abs(normDot) < std::numeric_limits<double>::epsilon()) return std::nullopt;
+        return (offset - pos.dot(planeNorm)) / normDot;
+    }
     /// Returns the position of the most left bit which is set to 1
     template <typename T> constexpr int maxBit(const T &number) {
         for (int bit = sizeof(number) * 8 - 1; bit >= 0; --bit) {
diff --git a/MuonSpectrometer/MuonPhaseII/MuonDetDescr/MuonGeoModelTestR4/python/testGeoModel.py b/MuonSpectrometer/MuonPhaseII/MuonDetDescr/MuonGeoModelTestR4/python/testGeoModel.py
index df0a80a4d24a8c8bb8925c79143fe31cb0d65b0c..81143a7d84bc1e19e30dfc54cf3befe4f1623bfa 100644
--- a/MuonSpectrometer/MuonPhaseII/MuonDetDescr/MuonGeoModelTestR4/python/testGeoModel.py
+++ b/MuonSpectrometer/MuonPhaseII/MuonDetDescr/MuonGeoModelTestR4/python/testGeoModel.py
@@ -141,7 +141,6 @@ def executeTest(cfg, num_events = 1):
     DetDescCnvSvc.IdDictFromRDB = False
     DetDescCnvSvc.MuonIDFileName="IdDictParser/IdDictMuonSpectrometer_R.10.00.xml"
     DetDescCnvSvc.MuonIDFileName="IdDictParser/IdDictMuonSpectrometer_R.09.03.xml"
-    cfg.merge(setupHistSvcCfg(flags, out_file = args.outRootFile))
     
     cfg.printConfig(withDetails=True, summariseProps=True)
     if not cfg.run(num_events).isSuccess(): exit(1)
diff --git a/MuonSpectrometer/MuonPhaseII/MuonDetDescr/MuonGeoModelTestR4/python/testSensitiveDetectors.py b/MuonSpectrometer/MuonPhaseII/MuonDetDescr/MuonGeoModelTestR4/python/testSensitiveDetectors.py
index 98ff8f8e993d9abef7792da8ab4925e9e629038a..6364ef9e541dfa12091dd31be8624a52bd7fe944 100644
--- a/MuonSpectrometer/MuonPhaseII/MuonDetDescr/MuonGeoModelTestR4/python/testSensitiveDetectors.py
+++ b/MuonSpectrometer/MuonPhaseII/MuonDetDescr/MuonGeoModelTestR4/python/testSensitiveDetectors.py
@@ -35,11 +35,11 @@ if __name__=="__main__":
     from G4AtlasAlg.G4AtlasAlgConfig import G4AtlasAlgCfg
     cfg.merge(G4AtlasAlgCfg(flags))
     ### Keep the Volume debugger commented for the moment
-    #from G4DebuggingTools.PostIncludes import VolumeDebuggerAtlas
-    #cfg.merge(VolumeDebuggerAtlas(flags, name="G4UA::UserActionSvc", 
-    #                                     PrintGeometry = True,
-    #                                     TargetVolume="BMS7_MDT04_1_3_2"
-    #                                    ))
+    from G4DebuggingTools.PostIncludes import VolumeDebuggerAtlas
+    cfg.merge(VolumeDebuggerAtlas(flags, name="G4UA::UserActionSvc", 
+                                         PrintGeometry = True,
+                                         TargetVolume="BIS7_RPC26_7_0_1_1_1"
+                                        ))
     
     ## xAOD TruthParticle conversion
     from xAODEventInfoCnv.xAODEventInfoCnvConfig import EventInfoCnvAlgCfg
diff --git a/MuonSpectrometer/MuonPhaseII/MuonDetDescr/MuonReadoutGeometryR4/MuonReadoutGeometryR4/RpcReadoutElement.h b/MuonSpectrometer/MuonPhaseII/MuonDetDescr/MuonReadoutGeometryR4/MuonReadoutGeometryR4/RpcReadoutElement.h
index 3b6a9fb328ee70b3d8829369628bc8e26c487ea7..91aacdc158352966e7c4d459ca1e011c85830bd7 100644
--- a/MuonSpectrometer/MuonPhaseII/MuonDetDescr/MuonReadoutGeometryR4/MuonReadoutGeometryR4/RpcReadoutElement.h
+++ b/MuonSpectrometer/MuonPhaseII/MuonDetDescr/MuonReadoutGeometryR4/MuonReadoutGeometryR4/RpcReadoutElement.h
@@ -90,6 +90,8 @@ class RpcReadoutElement : public MuonReadoutElement {
     double stripEtaLength() const;
     /// Returns the length of a phi strip
     double stripPhiLength() const;
+    /// Returns the thickness of a RPC gasgap
+    double gasGapThickness() const;
 
 
 
@@ -117,6 +119,9 @@ class RpcReadoutElement : public MuonReadoutElement {
         static bool measuresPhi(const IdentifierHash& measHash);
 
         Amg::Transform3D fromGapToChamOrigin(const IdentifierHash& layerHash) const;
+        /// Returns the local strip position w.r.t. to the chamber origin
+        Amg::Vector3D chamberStripPos(const IdentifierHash& measHash) const;
+
 
 
         parameterBook m_pars{};
@@ -134,6 +139,9 @@ class RpcReadoutElement : public MuonReadoutElement {
         const unsigned int m_hashShiftDbl{m_pars.hasPhiStrips ? 1u :0u};
         const unsigned int m_hashShiftGap{m_hashShiftDbl + (nPhiPanels() <= m_doubletPhi ? 0u : 1u)};
         const unsigned int m_hashShiftStr{m_hashShiftGap + MuonGM::maxBit(nGasGaps()) + 1};
+
+        /// Distance between 2 gas gaps (Radial direction)
+        double m_gasThickness{0.};
         
  };
 std::ostream& operator<<(std::ostream& ostr, const RpcReadoutElement::parameterBook& pars);
diff --git a/MuonSpectrometer/MuonPhaseII/MuonDetDescr/MuonReadoutGeometryR4/MuonReadoutGeometryR4/RpcReadoutElement.icc b/MuonSpectrometer/MuonPhaseII/MuonDetDescr/MuonReadoutGeometryR4/MuonReadoutGeometryR4/RpcReadoutElement.icc
index 3c2c60ea31c9b2a2fea7bdb8af8f2e806236448d..42e72c1c963a7a9a522487d417aa5fafa355dc0b 100644
--- a/MuonSpectrometer/MuonPhaseII/MuonDetDescr/MuonReadoutGeometryR4/MuonReadoutGeometryR4/RpcReadoutElement.icc
+++ b/MuonSpectrometer/MuonPhaseII/MuonDetDescr/MuonReadoutGeometryR4/MuonReadoutGeometryR4/RpcReadoutElement.icc
@@ -9,6 +9,7 @@ namespace MuonGMR4 {
     constexpr unsigned int minOne = -1;
 
 inline double RpcReadoutElement::thickness() const { return 2.* m_pars.halfThickness; }
+inline double  RpcReadoutElement::gasGapThickness() const {return m_gasThickness; }
 inline int RpcReadoutElement::doubletZ() const{ return m_doubletZ; }
 inline int RpcReadoutElement::doubletR() const{ return m_doubletR; }
 inline int RpcReadoutElement::doubletPhi() const{ return m_doubletPhi; }
diff --git a/MuonSpectrometer/MuonPhaseII/MuonDetDescr/MuonReadoutGeometryR4/src/RpcReadoutElement.cxx b/MuonSpectrometer/MuonPhaseII/MuonDetDescr/MuonReadoutGeometryR4/src/RpcReadoutElement.cxx
index 0cc561e23534cc8de2e2ba80077f40090905b3f9..9488fac1761566d7813390644ae4a81272c84057 100644
--- a/MuonSpectrometer/MuonPhaseII/MuonDetDescr/MuonReadoutGeometryR4/src/RpcReadoutElement.cxx
+++ b/MuonSpectrometer/MuonPhaseII/MuonDetDescr/MuonReadoutGeometryR4/src/RpcReadoutElement.cxx
@@ -45,6 +45,8 @@ StatusCode RpcReadoutElement::initElement() {
                                     return toStation(store) * fromGapToChamOrigin(hash); 
                                  }));
     }
+    m_gasThickness = (chamberStripPos(createHash(1, 2, 1, false)) - 
+                      chamberStripPos(createHash(1, 1, 1, false))).mag();
     return StatusCode::SUCCESS;
 }
 
@@ -66,6 +68,16 @@ Amg::Vector3D RpcReadoutElement::stripPosition(const ActsGeometryContext& ctx, c
                  <<" is out of range. Maximum range "<<m_pars.layers.size());
    return Amg::Vector3D::Zero();
 }
+Amg::Vector3D RpcReadoutElement::chamberStripPos(const IdentifierHash& measHash) const {
+   const IdentifierHash lHash = layerHash(measHash);
+   unsigned int layIdx = static_cast<unsigned int>(lHash);
+   if (layIdx < m_pars.layers.size()) {
+      return  m_pars.layers[layIdx].stripPosition(stripNumber(measHash));
+   }
+   ATH_MSG_WARNING(__FILE__<<":"<<__LINE__<<" The layer hash "<<layIdx
+                 <<" is out of range. Maximum range "<<m_pars.layers.size());
+   return Amg::Vector3D::Zero();
+}
     
 
 
diff --git a/MuonSpectrometer/MuonPhaseII/MuonG4/MuonSensitiveDetectorsR4/python/SensitiveDetectorsCfg.py b/MuonSpectrometer/MuonPhaseII/MuonG4/MuonSensitiveDetectorsR4/python/SensitiveDetectorsCfg.py
index aeb73e6967b6322afa60bd362fc70934c3ecbb07..104a6cb02b66083213159e09b8ce0420ba6e027b 100644
--- a/MuonSpectrometer/MuonPhaseII/MuonG4/MuonSensitiveDetectorsR4/python/SensitiveDetectorsCfg.py
+++ b/MuonSpectrometer/MuonPhaseII/MuonG4/MuonSensitiveDetectorsR4/python/SensitiveDetectorsCfg.py
@@ -1,6 +1,7 @@
 # Copyright (C) 2002-2023 CERN for the benefit of the ATLAS collaboration
 from AthenaConfiguration.ComponentAccumulator import ComponentAccumulator
 from AthenaConfiguration.ComponentFactory import CompFactory
+
 def MdtSensitiveDetectorToolCfg(flags, name = "MdtSensitiveDetector", **kwargs):
     result = ComponentAccumulator()
     kwargs.setdefault("OutputCollectionNames", [ "xMdtSimHits"])
@@ -9,11 +10,23 @@ def MdtSensitiveDetectorToolCfg(flags, name = "MdtSensitiveDetector", **kwargs):
     result.setPrivateTools(the_tool)
     return result
 
+def RpcSensitiveDetectorToolCfg(flags, name = "RpcSensitiveDetector", **kwargs):
+    result = ComponentAccumulator()
+    kwargs.setdefault("OutputCollectionNames", [ "xRpcSimHits"])
+    kwargs.setdefault("LogicalVolumeNames", ["MuonR4::RpcGasGap"])
+    the_tool = CompFactory.MuonG4R4.RpcSensitiveDetectorTool(name, **kwargs)
+    result.setPrivateTools(the_tool)
+    return result
+
 def SetupSensitiveDetectorsCfg(flags):
     result = ComponentAccumulator()
     tools = []
+
     if flags.Detector.EnableMDT:
         tools += [result.popToolsAndMerge(MdtSensitiveDetectorToolCfg(flags))]
-    result.setPrivateTools(tools)
     
+    if flags.Detector.EnableRPC:
+        tools += [result.popToolsAndMerge(RpcSensitiveDetectorToolCfg(flags))]
+
+    result.setPrivateTools(tools)
     return result
\ No newline at end of file
diff --git a/MuonSpectrometer/MuonPhaseII/MuonG4/MuonSensitiveDetectorsR4/src/RpcSensitiveDetector.cxx b/MuonSpectrometer/MuonPhaseII/MuonG4/MuonSensitiveDetectorsR4/src/RpcSensitiveDetector.cxx
new file mode 100755
index 0000000000000000000000000000000000000000..3f7358129084670441aad1fd0f31eded79c3f73f
--- /dev/null
+++ b/MuonSpectrometer/MuonPhaseII/MuonG4/MuonSensitiveDetectorsR4/src/RpcSensitiveDetector.cxx
@@ -0,0 +1,203 @@
+/*
+  Copyright (C) 2002-2023 CERN for the benefit of the ATLAS collaboration
+*/
+
+#include "RpcSensitiveDetector.h"
+
+#include "G4ThreeVector.hh"
+#include "G4Trd.hh"
+#include "G4Geantino.hh"
+#include "G4ChargedGeantino.hh"
+
+#include "MCTruth/TrackHelper.h"
+#include <sstream>
+
+#include "GeoPrimitives/CLHEPtoEigenConverter.h"
+#include "xAODMuonSimHit/MuonSimHitAuxContainer.h"
+#include "MuonReadoutGeometryR4/StringUtils.h"
+#include "GaudiKernel/SystemOfUnits.h"
+
+inline std::ostream& operator<<(std::ostream& ostr, const G4StepPoint& step) {
+      ostr<<"position: "<<Amg::toString(Amg::Hep3VectorToEigen(step.GetPosition()),2)<<", ";
+      ostr<<"momentum: "<<Amg::toString(Amg::Hep3VectorToEigen(step.GetMomentum()),2)<<", ";
+      ostr<<"time: "<<step.GetGlobalTime()<<", ";
+      ostr<<"mass: "<<step.GetMass()<<", ";
+      ostr<<"kinetic energy: "<<step.GetKineticEnergy()<<", ";
+      ostr<<"charge: "<<step.GetCharge();
+      return ostr;
+}
+
+
+inline Amg::Transform3D getTransform(const G4VTouchable* history, unsigned int level) {
+   return Amg::Translation3D{Amg::Hep3VectorToEigen(history->GetTranslation(level))}*
+            Amg::CLHEPRotationToEigen(*history->GetRotation(level)).inverse();
+}
+
+using namespace MuonGMR4;
+
+
+// construction/destruction
+namespace MuonG4R4 {
+
+RpcSensitiveDetector::RpcSensitiveDetector(const std::string& name, const std::string& output_key,
+                         const MuonGMR4::MuonDetectorManager* detMgr):
+    G4VSensitiveDetector{name},
+    AthMessaging{name},
+    m_writeHandle{output_key},
+    m_detMgr{detMgr} {}
+
+void RpcSensitiveDetector::Initialize(G4HCofThisEvent*) {
+  if (m_writeHandle.isValid()) {
+      ATH_MSG_VERBOSE("Simulation hit container "<<m_writeHandle.fullKey()<<" is already written");
+      return;
+  }
+  if (!m_writeHandle.recordNonConst(std::make_unique<xAOD::MuonSimHitContainer>(),
+                                    std::make_unique<xAOD::MuonSimHitAuxContainer>()).isSuccess()) {
+      ATH_MSG_FATAL(__FILE__<<":"<<__LINE__<<" Failed to record "<<m_writeHandle.fullKey());
+      throw std::runtime_error("Container saving is impossible");
+  }
+  ATH_MSG_DEBUG("Output container "<<m_writeHandle.fullKey()<<" has been successfully created");
+}
+
+G4bool RpcSensitiveDetector::ProcessHits(G4Step* aStep,G4TouchableHistory*) {
+
+  G4Track* track = aStep->GetTrack();
+
+  /** RPCs sensitive to charged particle only */
+  if (track->GetDefinition()->GetPDGCharge() == 0.0) {
+    if (track->GetDefinition()!=G4Geantino::GeantinoDefinition()) {
+      return true;
+    }
+  }
+
+  G4Track* currentTrack = aStep->GetTrack();
+
+  // MDTs sensitive to charged particle only
+  if (currentTrack->GetDefinition()->GetPDGCharge() == 0.0) {
+    if (currentTrack->GetDefinition()!= G4Geantino::GeantinoDefinition()) return true;
+    else if (currentTrack->GetDefinition()==G4ChargedGeantino::ChargedGeantinoDefinition()) return true;
+  }
+
+  const G4TouchableHistory* touchHist = static_cast<const G4TouchableHistory*>(aStep->GetPreStepPoint()->GetTouchable());
+  const MuonGMR4::RpcReadoutElement* readOutEle = getReadoutElement(touchHist);
+  if (!readOutEle) return true;
+
+
+  const G4StepPoint* preStep = aStep->GetPreStepPoint();
+  const G4StepPoint* postStep = aStep->GetPostStepPoint();
+  G4VPhysicalVolume*  physVolPostStep = postStep->GetPhysicalVolume();
+  if (!physVolPostStep)  {
+    ATH_MSG_VERBOSE("No physical volume stored "<<(*postStep));
+    return true;
+  }
+  
+  const Amg::Transform3D globalToLocal = getTransform(currentTrack->GetTouchable(), 0).inverse();
+  ATH_MSG_VERBOSE(" Track is inside volume "
+                 <<currentTrack->GetTouchable()->GetHistory()->GetTopVolume()->GetName()
+                 <<" transformation: "<<to_string(globalToLocal));
+  // transform pre and post step positions to local positions
+  const Amg::Vector3D globalVertex1{Amg::Hep3VectorToEigen(preStep->GetPosition())};
+  const Amg::Vector3D globalVertex2{Amg::Hep3VectorToEigen(postStep->GetPosition())};
+
+  const Amg::Vector3D localVertex1{globalToLocal * globalVertex1};
+  const Amg::Vector3D localVertex2{globalToLocal * globalVertex2};
+
+  const Amg::Vector3D localDir = (localVertex2 - localVertex1).unit();
+  ATH_MSG_VERBOSE("Entry / exit point in "<<m_detMgr->idHelperSvc()->toStringDetEl(readOutEle->identify())
+               <<" "<<Amg::toString(localVertex1, 2)<<" / "<<Amg::toString(localVertex2, 2));
+  
+  
+  // The middle of the gas gap is at X= 0
+  std::optional<double> travelDist = MuonGM::intersect<3>(localVertex1, localDir, Amg::Vector3D::UnitX(), 0.);
+  if (!travelDist) return true;
+  const Amg::Vector3D locGapCross = localVertex1 + (*travelDist) * localDir;
+  ATH_MSG_VERBOSE("Propagation to the gas gap center: "<<Amg::toString(locGapCross, 2));
+  const Amg::Vector3D gapCenterCross = globalToLocal.inverse() * locGapCross;
+
+  const Identifier etaHitID = getIdentifier(readOutEle, 
+                                            gapCenterCross, false);
+  if (!etaHitID.is_valid()) {
+      ATH_MSG_VERBOSE("No valid hit found");
+      return true;
+  }
+  const double globalTime = preStep->GetGlobalTime() + (*travelDist) / preStep->GetVelocity();
+  const Amg::Transform3D& gapTrans{readOutEle->globalToLocalTrans(m_gctx, etaHitID)};
+  const Amg::Vector3D dir = gapTrans*(globalVertex2 - globalVertex1).unit();
+  const Amg::Vector3D locDir = gapTrans * gapCenterCross;
+  
+  xAOD::MuonSimHit* hit = new xAOD::MuonSimHit();
+  m_writeHandle->push_back(hit);  
+  
+  TrackHelper trHelp(aStep->GetTrack());
+  hit->setIdentifier(etaHitID); 
+  hit->setLocalPosition(xAOD::toStorage(locDir));  
+  hit->setLocalDirection(xAOD::toStorage(dir));
+  hit->setStepLength(aStep->GetStepLength());
+  hit->setGlobalTime(globalTime);
+  hit->setPdgId(currentTrack->GetDefinition()->GetPDGEncoding());
+  hit->setEnergyDeposit(aStep->GetTotalEnergyDeposit());
+  hit->setKineticEnergy(preStep->GetKineticEnergy());
+  hit->setGenParticleLink(trHelp.GetParticleLink());
+  return true;
+}
+
+Identifier RpcSensitiveDetector::getIdentifier(const MuonGMR4::RpcReadoutElement* readOutEle, 
+                                               const Amg::Vector3D& hitAtGapPlane, bool phiGap) const {
+  const RpcIdHelper& idHelper{m_detMgr->idHelperSvc()->rpcIdHelper()};
+
+  const Identifier firstChan = idHelper.channelID(readOutEle->identify(),
+                                                  readOutEle->doubletZ(),
+                                                  readOutEle->doubletPhi(), 1, phiGap, 1);
+   
+  /// The 10 mm are introduced to match the old & new readout geometry. Revert the offset
+  /// to ease the determination of the gasgap
+  const Amg::Vector3D locHitPos{Amg::Translation3D{0.,0.,10.} *
+                                 readOutEle->globalToLocalTrans(m_gctx, firstChan) * 
+                                 hitAtGapPlane};   
+  const double gapHalfWidth = readOutEle->stripEtaLength() / 2;
+  const double gapHalfLength = readOutEle->stripPhiLength()/ 2;
+  ATH_MSG_VERBOSE("Detector element: "<<m_detMgr->idHelperSvc()->toStringDetEl(firstChan)
+                <<" locPos: "<<Amg::toString(locHitPos, 2)
+                <<" gap thickness "<<readOutEle->gasGapThickness()
+                <<" gap width: "<<gapHalfWidth
+                <<" gap length: "<<gapHalfLength);
+  const int doubletPhi = locHitPos.x() < - gapHalfWidth ? readOutEle->doubletPhiMax() :
+                                                          readOutEle->doubletPhi();
+  const int gasGap = (std::abs(locHitPos.z()) /  readOutEle->gasGapThickness()) + 1;
+
+  return idHelper.channelID(readOutEle->identify(),
+                            readOutEle->doubletZ(),
+                            doubletPhi, gasGap, phiGap, 1);
+}
+const MuonGMR4::RpcReadoutElement* RpcSensitiveDetector::getReadoutElement(const G4TouchableHistory* touchHist) const {
+    /// The fourth volume is the envelope volume of the rpc gas gap
+   const std::string stationVolume = touchHist->GetVolume(3)->GetName();
+   
+   const std::vector<std::string> volumeTokens = tokenize(stationVolume, "_");
+   ATH_MSG_VERBOSE("Name of the station volume is "<<stationVolume);
+   if (volumeTokens.size() != 7) {
+      ATH_MSG_FATAL(__FILE__<<":"<<__LINE__<<" Cannot deduce the station name from "<<stationVolume);
+      throw std::runtime_error("Invalid station Identifier");
+   }
+   /// Find the Detector element from the Identifier
+    ///       <STATIONETA>_(<STATIONPHI>-1)_<DOUBLETR>_<DOUBLETPHI>_<DOUBLETZ>
+   const std::string stName = volumeTokens[0].substr(0,3);
+   const int stationEta = atoi(volumeTokens[2]);
+   const int stationPhi = atoi(volumeTokens[3]) + 1;
+   const int doubletR = atoi(volumeTokens[4]);
+   const int doubletPhi = atoi(volumeTokens[5]);
+   const int doubletZ = atoi(volumeTokens[6]);
+   const RpcIdHelper& idHelper{m_detMgr->idHelperSvc()->rpcIdHelper()};
+
+   const Identifier detElId = idHelper.padID(idHelper.stationNameIndex(stName), 
+                                             stationEta, stationPhi, doubletR, doubletZ, doubletPhi);
+   const RpcReadoutElement* readOutElem = m_detMgr->getRpcReadoutElement(detElId);
+   if (!readOutElem) {
+      ATH_MSG_FATAL(__FILE__<<":"<<__LINE__<<" Failed to retrieve a valid detector element from "
+                    <<m_detMgr->idHelperSvc()->toStringDetEl(detElId)<<" "<<stationVolume);
+      /// Keep the failure for the moment commented because there're few ID issues
+      /// throw std::runtime_error("Invalid detector Element");
+   }
+   return readOutElem;
+}
+}
\ No newline at end of file
diff --git a/MuonSpectrometer/MuonPhaseII/MuonG4/MuonSensitiveDetectorsR4/src/RpcSensitiveDetector.h b/MuonSpectrometer/MuonPhaseII/MuonG4/MuonSensitiveDetectorsR4/src/RpcSensitiveDetector.h
new file mode 100755
index 0000000000000000000000000000000000000000..9d23794c10649811a959604e6d4672b5f1b21287
--- /dev/null
+++ b/MuonSpectrometer/MuonPhaseII/MuonG4/MuonSensitiveDetectorsR4/src/RpcSensitiveDetector.h
@@ -0,0 +1,124 @@
+#ifndef MUONG4R4_RPCSensitiveDetector_H
+#define MUONG4R4_RPCSensitiveDetector_H
+/*
+  Copyright (C) 2002-2023 CERN for the benefit of the ATLAS collaboration
+*/
+
+/** @class RPCSensitiveDetector
+    @section RPCSensitiveDetector Class methods and properties
+
+The method RPCSensitiveDetector::ProcessHits is executed by the G4 kernel each
+time a particle crosses one of the RPC gas gaps.
+Navigating with the touchableHistory method GetHistoryDepth()
+through the hierarchy of volumes crossed by the particle,
+the Sensitive Detector determines the correct set of Simulation Identifiers
+to associate to each hit. The RPC SimIDs are 32-bit unsigned integers, built 
+using the MuonSimEvent/RpcHitIdHelper class
+which inherits from the MuonHitIdHelper base class. 
+
+We describe here how each field of the identifier is determined.
+
+1) stationName, stationEta, stationPhi: when a volume is found in the hierarchy
+   whose name contains the substring "station", the stationName is extracted
+   from the volume's name. stationEta and stationPhi values are calculated
+   starting from the volume's copy number.
+   
+2) doubletR:  when a volume is found in the hierarchy whose name contains
+   the substring "rpccomponent", its copy number is used to calulate the
+   doubletR of the hit.
+   
+3) gasGap: the substring "layer" is searched through the names of the volumes
+   in the hierarchy. When a the volume is found, its copy number is used,
+   together with the sign of the copy number of the "rpccomponent" to decide
+   what is the correct gasGap value.
+   
+4) doubletPhi: when the volume with the substring "gas volume" in its name is
+   found, its copy number is enough to determine the doubletPhi of the hit
+   if the hit is registered in a standard chamber. For special chambers (i.e.
+   chambers with only one gas gap per layer readout by two strip panels per
+   direction, or chambers in the ribs) some a special attribution is done
+   using also the stationName.
+   
+5) doubletZ: "gazGap" is the required substring in the volume's name. doubletZ
+   attribution is fairly complicated. For standard chambers it just uses the
+   sign of the copy number of the "rpccomponent", but severa
+   arrangements (based on the stationName and the technology name) are needed
+   for special chambers.
+   
+6) the doubletPhi calculated as described above needs one further correction
+   for some chambers to take into account chamber rotation before placement
+   in the spectrometer. This is done just before creating the hit.
+
+    @section Some notes:
+
+1) presently, chamber efficiency is assumed to be 100% at the level of the
+   Sensitive Detector, i.e. two hits (one in eta and one in phi) are created
+   each time the sensitive detector is created.
+
+2) the hits created as described above contain information about their local
+   position *in the gas gap* and the simulation identifier of *the strip panel*
+
+3) the strip number is not assigned by the sensitive detector, and must be
+   calculated by the digitization algorithm.
+
+4) points 2 and 3 are due to the necessity of not introducing any dependency
+   on the geometry description in the Sensitive Detector.
+
+5) the present version of the Sensitive Detector produces hits which are
+   different depending on whether hand coded G4 geometry or GeoModel is used
+   for the geometrical description of the setup. When hand coded geometry
+   is used, both GeoModel and the old MuonDetDescr can be used for digitization
+   (using a proper tag of RPC_Digitization), while if GeoModel is used in the
+   simulation, *it must be used* also in the digitization.
+
+6) for each hit, the time of flight (the G4 globalTime), is recorded and
+   associated to the hit.
+
+7) the RPCHit object contains: the SimID, the globalTime, the hit local position
+   and the track barcode.
+
+
+*/
+
+#include "G4VSensitiveDetector.hh"
+
+#include <StoreGate/WriteHandle.h>
+#include <AthenaBaseComps/AthMessaging.h>
+#include <MuonReadoutGeometryR4/MuonDetectorManager.h>
+#include <MuonReadoutGeometryR4/RpcReadoutElement.h>
+#include <xAODMuonSimHit/MuonSimHitContainer.h>
+#include <AthenaBaseComps/AthMessaging.h>
+
+namespace MuonG4R4 {
+
+class RpcSensitiveDetector : public G4VSensitiveDetector, public AthMessaging {
+
+
+public:
+     /** construction/destruction */
+    RpcSensitiveDetector(const std::string& name, const std::string& output_key,
+                         const MuonGMR4::MuonDetectorManager* detMgr);
+    ~RpcSensitiveDetector()=default;
+    
+    /** member functions */
+    void   Initialize(G4HCofThisEvent* HCE) override final;
+    G4bool ProcessHits(G4Step* aStep, G4TouchableHistory* ROhist) override final;
+
+    
+private:
+    /// Retrieves the matching readout element to a G4 hit
+    const MuonGMR4::RpcReadoutElement* getReadoutElement(const G4TouchableHistory* touchHist) const;
+    Identifier getIdentifier(const MuonGMR4::RpcReadoutElement* readOutEle, const Amg::Vector3D& hitAtGapPlane, bool phiGap) const;
+    /* For the moment use write handles because the sensitive detectors are 
+     *  managed by a service which must not have a data dependency
+    */
+    SG::WriteHandle<xAOD::MuonSimHitContainer> m_writeHandle;
+    /// Pointer to the acts Geometry context
+    const ActsGeometryContext m_gctx{};
+    /// Pointer to the underlying detector manager
+    const MuonGMR4::MuonDetectorManager* m_detMgr{nullptr};
+   
+};
+
+}
+#endif
diff --git a/MuonSpectrometer/MuonPhaseII/MuonG4/MuonSensitiveDetectorsR4/src/RpcSensitiveDetectorTool.cxx b/MuonSpectrometer/MuonPhaseII/MuonG4/MuonSensitiveDetectorsR4/src/RpcSensitiveDetectorTool.cxx
new file mode 100755
index 0000000000000000000000000000000000000000..9835ad79234899bef3d9dc8015e848a4d4a3a6f5
--- /dev/null
+++ b/MuonSpectrometer/MuonPhaseII/MuonG4/MuonSensitiveDetectorsR4/src/RpcSensitiveDetectorTool.cxx
@@ -0,0 +1,20 @@
+/*
+  Copyright (C) 2002-2023 CERN for the benefit of the ATLAS collaboration
+*/
+#include "RpcSensitiveDetectorTool.h"
+#include "RpcSensitiveDetector.h"
+
+namespace MuonG4R4 {
+
+RpcSensitiveDetectorTool::RpcSensitiveDetectorTool(const std::string& type, const std::string& name, const IInterface* parent)
+  : SensitiveDetectorBase( type , name , parent ) {}
+
+StatusCode RpcSensitiveDetectorTool::initialize() {
+    ATH_CHECK(SensitiveDetectorBase::initialize());
+    ATH_CHECK(detStore()->retrieve(m_detMgr));
+    return StatusCode::SUCCESS;
+}
+G4VSensitiveDetector* RpcSensitiveDetectorTool::makeSD() const {
+  return new RpcSensitiveDetector(name(), m_outputCollectionNames[0], m_detMgr);
+}
+}
diff --git a/MuonSpectrometer/MuonPhaseII/MuonG4/MuonSensitiveDetectorsR4/src/RpcSensitiveDetectorTool.h b/MuonSpectrometer/MuonPhaseII/MuonG4/MuonSensitiveDetectorsR4/src/RpcSensitiveDetectorTool.h
new file mode 100755
index 0000000000000000000000000000000000000000..c204185e267ae749443cabdd1597bfcd802151a7
--- /dev/null
+++ b/MuonSpectrometer/MuonPhaseII/MuonG4/MuonSensitiveDetectorsR4/src/RpcSensitiveDetectorTool.h
@@ -0,0 +1,30 @@
+/*
+  Copyright (C) 2002-2023 CERN for the benefit of the ATLAS collaboration
+*/
+#ifndef MUONG4R4_RpcSENSITIVEDETECTORTOOL_H
+#define MUONG4R4_RpcSENSITIVEDETECTORTOOL_H
+
+#include <GeoPrimitives/GeoPrimitives.h>
+
+#include <G4AtlasTools/SensitiveDetectorBase.h>
+#include <MuonReadoutGeometryR4/MuonDetectorManager.h>
+#include <xAODMuonSimHit/MuonSimHitContainer.h>
+
+namespace MuonG4R4 {
+
+class RpcSensitiveDetectorTool : public SensitiveDetectorBase {
+
+public:
+    RpcSensitiveDetectorTool(const std::string& type, const std::string& name, const IInterface *parent);
+    ~RpcSensitiveDetectorTool()=default;
+
+    StatusCode initialize() override final;
+protected:
+    G4VSensitiveDetector* makeSD() const override final;
+private:
+    const MuonGMR4::MuonDetectorManager* m_detMgr{nullptr};
+   
+};
+}
+
+#endif
diff --git a/MuonSpectrometer/MuonPhaseII/MuonG4/MuonSensitiveDetectorsR4/src/components/MuonSensitiveDetectorsR4_entries.cxx b/MuonSpectrometer/MuonPhaseII/MuonG4/MuonSensitiveDetectorsR4/src/components/MuonSensitiveDetectorsR4_entries.cxx
index 91568059cb6802df2677f6ddbd0036460a46f16e..5b303eb0832b30778222da1c8257b69ab51c9d3f 100644
--- a/MuonSpectrometer/MuonPhaseII/MuonG4/MuonSensitiveDetectorsR4/src/components/MuonSensitiveDetectorsR4_entries.cxx
+++ b/MuonSpectrometer/MuonPhaseII/MuonG4/MuonSensitiveDetectorsR4/src/components/MuonSensitiveDetectorsR4_entries.cxx
@@ -2,4 +2,7 @@
   Copyright (C) 2002-2023 CERN for the benefit of the ATLAS collaboration
 */
 #include "../MdtSensitiveDetectorTool.h"
+#include "../RpcSensitiveDetectorTool.h"
+
 DECLARE_COMPONENT( MuonG4R4::MdtSensitiveDetectorTool )
+DECLARE_COMPONENT( MuonG4R4::RpcSensitiveDetectorTool )