diff --git a/Tracking/TrkDetDescr/TrkSurfaces/TrkSurfaces/AnnulusBounds.h b/Tracking/TrkDetDescr/TrkSurfaces/TrkSurfaces/AnnulusBounds.h
index 65e7b4ba24d32ddf115133ccdb66dcefc11ab00e..fb82d09e564b1aa0ccc917282ccf4027dc02378f 100644
--- a/Tracking/TrkDetDescr/TrkSurfaces/TrkSurfaces/AnnulusBounds.h
+++ b/Tracking/TrkDetDescr/TrkSurfaces/TrkSurfaces/AnnulusBounds.h
@@ -41,7 +41,7 @@ namespace Trk {
  @author Marcin.Wolter@cern.ch
  */
 
-class AnnulusBounds : public SurfaceBounds
+class AnnulusBounds final: public SurfaceBounds
 {
 
 public:
diff --git a/Tracking/TrkDetDescr/TrkSurfaces/TrkSurfaces/DiscTrapezoidalBounds.h b/Tracking/TrkDetDescr/TrkSurfaces/TrkSurfaces/DiscTrapezoidalBounds.h
index a5a858435077d52b68a08f2e76561229c1613914..dc8bab54f2e995bc08f143e0758c7f5bc2f97d61 100644
--- a/Tracking/TrkDetDescr/TrkSurfaces/TrkSurfaces/DiscTrapezoidalBounds.h
+++ b/Tracking/TrkDetDescr/TrkSurfaces/TrkSurfaces/DiscTrapezoidalBounds.h
@@ -38,7 +38,7 @@ namespace Trk {
    @author Noemi.Calace@cern.ch , Andreas.Salzburger@cern.ch
 */
 
-class DiscTrapezoidalBounds : public SurfaceBounds
+class DiscTrapezoidalBounds final : public SurfaceBounds
 {
 
 public:
diff --git a/Tracking/TrkDetDescr/TrkSurfaces/TrkSurfaces/NoBounds.h b/Tracking/TrkDetDescr/TrkSurfaces/TrkSurfaces/NoBounds.h
index 19d82c451c054118757379b4496fb9b9c8c2294b..a30e1d443b56c3c7a07426518460fd8bafd9574b 100644
--- a/Tracking/TrkDetDescr/TrkSurfaces/TrkSurfaces/NoBounds.h
+++ b/Tracking/TrkDetDescr/TrkSurfaces/TrkSurfaces/NoBounds.h
@@ -26,7 +26,7 @@ namespace Trk {
  @author Andreas.Salzburger@cern.ch
  */
 
-class NoBounds : public SurfaceBounds
+class NoBounds final : public SurfaceBounds
 {
 public:
   /**Default Constructor*/
@@ -36,7 +36,7 @@ public:
   ~NoBounds() = default;
 
   /**Equality operator */
-  virtual bool operator==(const SurfaceBounds& sbo) const override;
+  virtual bool operator==(const SurfaceBounds& sbo) const override final;
 
   /** Return SurfaceBounds for persistency */
   virtual SurfaceBounds::BoundsType type() const override { return SurfaceBounds::Other; }
@@ -67,71 +67,8 @@ public:
   /** Output Method for std::ostream */
   virtual std::ostream& dump(std::ostream& sl) const override final;
 
-private:
 };
-
-inline bool
-NoBounds::operator==(const SurfaceBounds&) const
-{
-  return true;
-}
-
-inline bool
-NoBounds::inside(const Amg::Vector2D&, double, double) const
-{
-  return true;
-}
-
-inline bool
-NoBounds::inside(const Amg::Vector2D&, const BoundaryCheck&) const
-{
-  return true;
-}
-
-inline bool
-NoBounds::insideLoc1(const Amg::Vector2D&, double) const
-{
-  return true;
-}
-
-inline bool
-NoBounds::insideLoc2(const Amg::Vector2D&, double) const
-{
-  return true;
-}
-
-inline double
-NoBounds::minDistance(const Amg::Vector2D&) const
-{
-  return 0.;
-}
-
-inline NoBounds*
-NoBounds::clone() const
-{
-  return new NoBounds();
-}
-
-inline double
-NoBounds::r() const
-{
-  return 0.;
-}
-
-inline MsgStream&
-NoBounds::dump(MsgStream& sl) const
-{
-  sl << "Trk::NoBounds ... boundless surface" << endmsg;
-  return sl;
-}
-
-inline std::ostream&
-NoBounds::dump(std::ostream& sl) const
-{
-  sl << "Trk::NoBounds ... boundless surface" << std::endl;
-  return sl;
-}
-
 } // end of namespace
 
+#include "TrkSurfaces/NoBounds.icc"
 #endif // TRKSURFACES_NOBOUNDS_H
diff --git a/Tracking/TrkDetDescr/TrkSurfaces/TrkSurfaces/NoBounds.icc b/Tracking/TrkDetDescr/TrkSurfaces/TrkSurfaces/NoBounds.icc
new file mode 100644
index 0000000000000000000000000000000000000000..406f6af93de586c6b7a5e4241178958d05c53fae
--- /dev/null
+++ b/Tracking/TrkDetDescr/TrkSurfaces/TrkSurfaces/NoBounds.icc
@@ -0,0 +1,69 @@
+/*
+  Copyright (C) 2002-2020 CERN for the benefit of the ATLAS collaboration
+*/
+
+namespace Trk {
+
+inline bool
+NoBounds::operator==(const SurfaceBounds&) const
+{
+  return true;
+}
+
+inline bool
+NoBounds::inside(const Amg::Vector2D&, double, double) const
+{
+  return true;
+}
+
+inline bool
+NoBounds::inside(const Amg::Vector2D&, const BoundaryCheck&) const
+{
+  return true;
+}
+
+inline bool
+NoBounds::insideLoc1(const Amg::Vector2D&, double) const
+{
+  return true;
+}
+
+inline bool
+NoBounds::insideLoc2(const Amg::Vector2D&, double) const
+{
+  return true;
+}
+
+inline double
+NoBounds::minDistance(const Amg::Vector2D&) const
+{
+  return 0.;
+}
+
+inline NoBounds*
+NoBounds::clone() const
+{
+  return new NoBounds();
+}
+
+inline double
+NoBounds::r() const
+{
+  return 0.;
+}
+
+inline MsgStream&
+NoBounds::dump(MsgStream& sl) const
+{
+  sl << "Trk::NoBounds ... boundless surface" << endmsg;
+  return sl;
+}
+
+inline std::ostream&
+NoBounds::dump(std::ostream& sl) const
+{
+  sl << "Trk::NoBounds ... boundless surface" << std::endl;
+  return sl;
+}
+
+} // end of namespace
diff --git a/Tracking/TrkDetDescr/TrkSurfaces/TrkSurfaces/RectangleBounds.h b/Tracking/TrkDetDescr/TrkSurfaces/TrkSurfaces/RectangleBounds.h
index 2165659c860c1adca4815e1f55d8e4a44abac0bc..3d13b6bd2c982a65c7d11bdcc2d9bec79d85b904 100644
--- a/Tracking/TrkDetDescr/TrkSurfaces/TrkSurfaces/RectangleBounds.h
+++ b/Tracking/TrkDetDescr/TrkSurfaces/TrkSurfaces/RectangleBounds.h
@@ -1,5 +1,5 @@
 /*
-  Copyright (C) 2002-2017 CERN for the benefit of the ATLAS collaboration
+  Copyright (C) 2002-2020 CERN for the benefit of the ATLAS collaboration
 */
 
 ///////////////////////////////////////////////////////////////////
@@ -34,7 +34,7 @@ namespace Trk {
 
  @author Andreas.Salzburger@cern.ch */
 
-class RectangleBounds : public SurfaceBounds
+class RectangleBounds final : public SurfaceBounds
 {
 
 public:
@@ -62,7 +62,7 @@ public:
   RectangleBounds& operator=(const RectangleBounds& recbo);
 
   /**Equality operator*/
-  virtual bool operator==(const SurfaceBounds& sbo) const override;
+  virtual bool operator==(const SurfaceBounds& sbo) const override final;
 
   /**Virtual constructor*/
   virtual RectangleBounds* clone() const override;
@@ -113,109 +113,6 @@ private:
   std::vector<TDD_real_t> m_boundValues;
 };
 
-inline RectangleBounds*
-RectangleBounds::clone() const
-{
-  return new RectangleBounds(*this);
-}
-
-inline bool
-RectangleBounds::inside(const Amg::Vector2D& locpo, double tol1, double tol2) const
-{
-  return ((fabs(locpo[locX]) < m_boundValues[RectangleBounds::bv_halfX] + tol1) &&
-          (fabs(locpo[locY]) < m_boundValues[RectangleBounds::bv_halfY] + tol2));
-}
-
-inline bool
-RectangleBounds::inside(const Amg::Vector2D& locpo, const BoundaryCheck& bchk) const
-{
-  if (bchk.bcType == 0)
-    return RectangleBounds::inside(locpo, bchk.toleranceLoc1, bchk.toleranceLoc2);
-
-  // a fast FALSE
-  double max_ell = bchk.lCovariance(0, 0) > bchk.lCovariance(1, 1) ? bchk.lCovariance(0, 0) : bchk.lCovariance(1, 1);
-  double limit = bchk.nSigmas * sqrt(max_ell);
-  if (!RectangleBounds::inside(locpo, limit, limit))
-    return false;
-  // a fast TRUE
-  double min_ell = bchk.lCovariance(0, 0) < bchk.lCovariance(1, 1) ? bchk.lCovariance(0, 0) : bchk.lCovariance(1, 1);
-  limit = bchk.nSigmas * sqrt(min_ell);
-  if (RectangleBounds::inside(locpo, limit, limit))
-    return true;
-
-  // compute KDOP and axes for surface polygon
-  std::vector<KDOP> elementKDOP(4);
-  std::vector<Amg::Vector2D> elementP(4);
-  float theta = (bchk.lCovariance(1, 0) != 0 && (bchk.lCovariance(1, 1) - bchk.lCovariance(0, 0)) != 0)
-                  ? .5 * bchk.FastArcTan(2 * bchk.lCovariance(1, 0) / (bchk.lCovariance(1, 1) - bchk.lCovariance(0, 0)))
-                  : 0.;
-  sincosCache scResult = bchk.FastSinCos(theta);
-  AmgMatrix(2, 2) rotMatrix;
-  rotMatrix << scResult.cosC, scResult.sinC, -scResult.sinC, scResult.cosC;
-  // ellipse is always at (0,0), surface is moved to ellipse position and then rotated
-  Amg::Vector2D p;
-  p << m_boundValues[RectangleBounds::bv_halfX], m_boundValues[RectangleBounds::bv_halfY];
-  elementP[0] = (rotMatrix * (p - locpo));
-  p << m_boundValues[RectangleBounds::bv_halfX], -m_boundValues[RectangleBounds::bv_halfY];
-  elementP[1] = (rotMatrix * (p - locpo));
-  p << -m_boundValues[RectangleBounds::bv_halfX], m_boundValues[RectangleBounds::bv_halfY];
-  elementP[2] = (rotMatrix * (p - locpo));
-  p << -m_boundValues[RectangleBounds::bv_halfX], -m_boundValues[RectangleBounds::bv_halfY];
-  elementP[3] = (rotMatrix * (p - locpo));
-  std::vector<Amg::Vector2D> axis = {
-    elementP[0] - elementP[1], elementP[0] - elementP[2], elementP[0] - elementP[3], elementP[1] - elementP[2]
-  };
-  bchk.ComputeKDOP(elementP, axis, elementKDOP);
-  // compute KDOP for error ellipse
-  std::vector<KDOP> errelipseKDOP(4);
-  bchk.ComputeKDOP(bchk.EllipseToPoly(3), axis, errelipseKDOP);
-  // check if KDOPs overlap and return result
-  return bchk.TestKDOPKDOP(elementKDOP, errelipseKDOP);
-}
-
-inline bool
-RectangleBounds::insideLoc1(const Amg::Vector2D& locpo, double tol1) const
-{
-  return (fabs(locpo[locX]) < m_boundValues[RectangleBounds::bv_halfX] + tol1);
-}
-
-inline bool
-RectangleBounds::insideLoc2(const Amg::Vector2D& locpo, double tol2) const
-{
-  return (fabs(locpo[locY]) < m_boundValues[RectangleBounds::bv_halfY] + tol2);
-}
-
-inline double
-RectangleBounds::halflengthPhi() const
-{
-  return this->halflengthX();
-}
-
-inline double
-RectangleBounds::halflengthEta() const
-{
-  return this->halflengthY();
-}
-
-inline double
-RectangleBounds::halflengthX() const
-{
-  return m_boundValues[RectangleBounds::bv_halfX];
-}
-
-inline double
-RectangleBounds::halflengthY() const
-{
-  return m_boundValues[RectangleBounds::bv_halfY];
-}
-
-inline double
-RectangleBounds::r() const
-{
-  return sqrt(m_boundValues[RectangleBounds::bv_halfX] * m_boundValues[RectangleBounds::bv_halfX] +
-              m_boundValues[RectangleBounds::bv_halfY] * m_boundValues[RectangleBounds::bv_halfY]);
-}
-
 } // end of namespace
-
+#include "TrkSurfaces/RectangleBounds.icc"
 #endif // TRKSURFACES_RECTANGLEBOUNDS_H
diff --git a/Tracking/TrkDetDescr/TrkSurfaces/TrkSurfaces/RectangleBounds.icc b/Tracking/TrkDetDescr/TrkSurfaces/TrkSurfaces/RectangleBounds.icc
new file mode 100644
index 0000000000000000000000000000000000000000..37e8827c59d8b9a4b6d0cea98a231749aab0ff34
--- /dev/null
+++ b/Tracking/TrkDetDescr/TrkSurfaces/TrkSurfaces/RectangleBounds.icc
@@ -0,0 +1,130 @@
+/*
+  Copyright (C) 2002-2020 CERN for the benefit of the ATLAS collaboration
+*/
+
+namespace Trk {
+inline RectangleBounds*
+RectangleBounds::clone() const
+{
+  return new RectangleBounds(*this);
+}
+
+inline bool
+RectangleBounds::inside(const Amg::Vector2D& locpo,
+                        double tol1,
+                        double tol2) const
+{
+  return (
+    (fabs(locpo[locX]) < m_boundValues[RectangleBounds::bv_halfX] + tol1) &&
+    (fabs(locpo[locY]) < m_boundValues[RectangleBounds::bv_halfY] + tol2));
+}
+
+inline bool
+RectangleBounds::inside(const Amg::Vector2D& locpo,
+                        const BoundaryCheck& bchk) const
+{
+  if (bchk.bcType == 0)
+    return RectangleBounds::inside(
+      locpo, bchk.toleranceLoc1, bchk.toleranceLoc2);
+
+  // a fast FALSE
+  double max_ell = bchk.lCovariance(0, 0) > bchk.lCovariance(1, 1)
+                     ? bchk.lCovariance(0, 0)
+                     : bchk.lCovariance(1, 1);
+  double limit = bchk.nSigmas * sqrt(max_ell);
+  if (!RectangleBounds::inside(locpo, limit, limit))
+    return false;
+  // a fast TRUE
+  double min_ell = bchk.lCovariance(0, 0) < bchk.lCovariance(1, 1)
+                     ? bchk.lCovariance(0, 0)
+                     : bchk.lCovariance(1, 1);
+  limit = bchk.nSigmas * sqrt(min_ell);
+  if (RectangleBounds::inside(locpo, limit, limit))
+    return true;
+
+  // compute KDOP and axes for surface polygon
+  std::vector<KDOP> elementKDOP(4);
+  std::vector<Amg::Vector2D> elementP(4);
+  float theta =
+    (bchk.lCovariance(1, 0) != 0 &&
+     (bchk.lCovariance(1, 1) - bchk.lCovariance(0, 0)) != 0)
+      ? .5 * bchk.FastArcTan(2 * bchk.lCovariance(1, 0) /
+                             (bchk.lCovariance(1, 1) - bchk.lCovariance(0, 0)))
+      : 0.;
+  sincosCache scResult = bchk.FastSinCos(theta);
+  AmgMatrix(2, 2) rotMatrix;
+  rotMatrix << scResult.cosC, scResult.sinC, -scResult.sinC, scResult.cosC;
+  // ellipse is always at (0,0), surface is moved to ellipse position and then
+  // rotated
+  Amg::Vector2D p;
+  p << m_boundValues[RectangleBounds::bv_halfX],
+    m_boundValues[RectangleBounds::bv_halfY];
+  elementP[0] = (rotMatrix * (p - locpo));
+  p << m_boundValues[RectangleBounds::bv_halfX],
+    -m_boundValues[RectangleBounds::bv_halfY];
+  elementP[1] = (rotMatrix * (p - locpo));
+  p << -m_boundValues[RectangleBounds::bv_halfX],
+    m_boundValues[RectangleBounds::bv_halfY];
+  elementP[2] = (rotMatrix * (p - locpo));
+  p << -m_boundValues[RectangleBounds::bv_halfX],
+    -m_boundValues[RectangleBounds::bv_halfY];
+  elementP[3] = (rotMatrix * (p - locpo));
+  std::vector<Amg::Vector2D> axis = { elementP[0] - elementP[1],
+                                      elementP[0] - elementP[2],
+                                      elementP[0] - elementP[3],
+                                      elementP[1] - elementP[2] };
+  bchk.ComputeKDOP(elementP, axis, elementKDOP);
+  // compute KDOP for error ellipse
+  std::vector<KDOP> errelipseKDOP(4);
+  bchk.ComputeKDOP(bchk.EllipseToPoly(3), axis, errelipseKDOP);
+  // check if KDOPs overlap and return result
+  return bchk.TestKDOPKDOP(elementKDOP, errelipseKDOP);
+}
+
+inline bool
+RectangleBounds::insideLoc1(const Amg::Vector2D& locpo, double tol1) const
+{
+  return (fabs(locpo[locX]) < m_boundValues[RectangleBounds::bv_halfX] + tol1);
+}
+
+inline bool
+RectangleBounds::insideLoc2(const Amg::Vector2D& locpo, double tol2) const
+{
+  return (fabs(locpo[locY]) < m_boundValues[RectangleBounds::bv_halfY] + tol2);
+}
+
+inline double
+RectangleBounds::halflengthPhi() const
+{
+  return this->halflengthX();
+}
+
+inline double
+RectangleBounds::halflengthEta() const
+{
+  return this->halflengthY();
+}
+
+inline double
+RectangleBounds::halflengthX() const
+{
+  return m_boundValues[RectangleBounds::bv_halfX];
+}
+
+inline double
+RectangleBounds::halflengthY() const
+{
+  return m_boundValues[RectangleBounds::bv_halfY];
+}
+
+inline double
+RectangleBounds::r() const
+{
+  return sqrt(m_boundValues[RectangleBounds::bv_halfX] *
+                m_boundValues[RectangleBounds::bv_halfX] +
+              m_boundValues[RectangleBounds::bv_halfY] *
+                m_boundValues[RectangleBounds::bv_halfY]);
+}
+
+} // end of namespace
+
diff --git a/Tracking/TrkDetDescr/TrkSurfaces/TrkSurfaces/RotatedDiamondBounds.h b/Tracking/TrkDetDescr/TrkSurfaces/TrkSurfaces/RotatedDiamondBounds.h
index 9db6926a50647f0d3dc886e8f1a5c58da771b492..6956160a5a520f8ce15be3744792b81760e3a879 100644
--- a/Tracking/TrkDetDescr/TrkSurfaces/TrkSurfaces/RotatedDiamondBounds.h
+++ b/Tracking/TrkDetDescr/TrkSurfaces/TrkSurfaces/RotatedDiamondBounds.h
@@ -38,7 +38,7 @@ in the constructor.)
 Michiel.Jan.Veen@cern.ch
 */
 
-class RotatedDiamondBounds : public SurfaceBounds
+class RotatedDiamondBounds final: public SurfaceBounds
 {
 
 public:
@@ -141,142 +141,7 @@ private:
   TDD_real_t m_alpha2;
 };
 
-inline RotatedDiamondBounds*
-RotatedDiamondBounds::clone() const
-{
-  return new RotatedDiamondBounds(*this);
-}
-
-inline double
-RotatedDiamondBounds::minHalflengthX() const
-{
-  return m_boundValues[RotatedDiamondBounds::bv_minHalfX];
-}
-
-inline double
-RotatedDiamondBounds::medHalflengthX() const
-{
-  return m_boundValues[RotatedDiamondBounds::bv_medHalfX];
-}
-
-inline double
-RotatedDiamondBounds::maxHalflengthX() const
-{
-  return m_boundValues[RotatedDiamondBounds::bv_maxHalfX];
-}
-
-inline double
-RotatedDiamondBounds::halflengthY1() const
-{
-  return m_boundValues[RotatedDiamondBounds::bv_halfY1];
-}
-
-inline double
-RotatedDiamondBounds::halflengthY2() const
-{
-  return m_boundValues[RotatedDiamondBounds::bv_halfY2];
-}
-
-inline double
-RotatedDiamondBounds::r() const
-{
-  return sqrt(m_boundValues[RotatedDiamondBounds::bv_medHalfX] * m_boundValues[RotatedDiamondBounds::bv_medHalfX] +
-              m_boundValues[RotatedDiamondBounds::bv_halfY1] * m_boundValues[RotatedDiamondBounds::bv_halfY1]);
-}
-
-inline bool
-RotatedDiamondBounds::inside(const Amg::Vector2D& locpo, const BoundaryCheck& bchk) const
-{
-  // locX and locY are interchanged wrt DiamondBounds
-  if (bchk.bcType == 0)
-    return RotatedDiamondBounds::inside(locpo, bchk.toleranceLoc1, bchk.toleranceLoc2);
-  // a fast FALSE
-  double max_ell = bchk.lCovariance(0, 0) > bchk.lCovariance(1, 1) ? bchk.lCovariance(0, 0) : bchk.lCovariance(1, 1);
-  double limit = bchk.nSigmas * sqrt(max_ell);
-  if (locpo[Trk::locX] < -2 * m_boundValues[RotatedDiamondBounds::bv_halfY1] - limit)
-    return false;
-  if (locpo[Trk::locX] > 2 * m_boundValues[RotatedDiamondBounds::bv_halfY2] + limit)
-    return false;
-  // a fast FALSE
-  double fabsX = fabs(locpo[Trk::locY]);
-  if (fabsX > (m_boundValues[RotatedDiamondBounds::bv_medHalfX] + limit))
-    return false;
-  // a fast TRUE
-  double min_ell = bchk.lCovariance(0, 0) < bchk.lCovariance(1, 1) ? bchk.lCovariance(0, 0) : bchk.lCovariance(1, 1);
-  limit = bchk.nSigmas * sqrt(min_ell);
-  if (fabsX <
-      (fmin(m_boundValues[RotatedDiamondBounds::bv_minHalfX], m_boundValues[RotatedDiamondBounds::bv_maxHalfX]) -
-       limit))
-    return true;
-  // a fast TRUE
-  if (fabs(locpo[Trk::locX]) <
-      (fmin(m_boundValues[RotatedDiamondBounds::bv_halfY1], m_boundValues[RotatedDiamondBounds::bv_halfY2]) - limit))
-    return true;
-
-  // compute KDOP and axes for surface polygon
-  std::vector<KDOP> elementKDOP(5);
-  std::vector<Amg::Vector2D> elementP(6);
-  float theta = (bchk.lCovariance(1, 0) != 0 && (bchk.lCovariance(1, 1) - bchk.lCovariance(0, 0)) != 0)
-                  ? .5 * bchk.FastArcTan(2 * bchk.lCovariance(1, 0) / (bchk.lCovariance(1, 1) - bchk.lCovariance(0, 0)))
-                  : 0.;
-  sincosCache scResult = bchk.FastSinCos(theta);
-  AmgMatrix(2, 2) rotMatrix;
-  rotMatrix << scResult.cosC, scResult.sinC, -scResult.sinC, scResult.cosC;
-  AmgMatrix(2, 2) normal;
-  normal << 0, -1, 1, 0;
-  // ellipse is always at (0,0), surface is moved to ellipse position and then rotated
-  // exchange locX and locY
-  Amg::Vector2D locpoF;
-  locpoF[0] = locpo[Trk::locY];
-  locpoF[1] = locpo[Trk::locX];
-  Amg::Vector2D p;
-  p << -m_boundValues[RotatedDiamondBounds::bv_minHalfX], -2. * m_boundValues[RotatedDiamondBounds::bv_halfY1];
-  elementP[0] = (rotMatrix * (p - locpoF));
-  p << -m_boundValues[RotatedDiamondBounds::bv_medHalfX], 0.;
-  elementP[1] = (rotMatrix * (p - locpoF));
-  p << -m_boundValues[RotatedDiamondBounds::bv_maxHalfX], 2. * m_boundValues[RotatedDiamondBounds::bv_halfY2];
-  elementP[2] = (rotMatrix * (p - locpoF));
-  p << m_boundValues[RotatedDiamondBounds::bv_maxHalfX], 2. * m_boundValues[RotatedDiamondBounds::bv_halfY2];
-  elementP[3] = (rotMatrix * (p - locpoF));
-  p << m_boundValues[RotatedDiamondBounds::bv_medHalfX], 0.;
-  elementP[4] = (rotMatrix * (p - locpoF));
-  p << m_boundValues[RotatedDiamondBounds::bv_minHalfX], -2. * m_boundValues[RotatedDiamondBounds::bv_halfY1];
-  elementP[5] = (rotMatrix * (p - locpoF));
-  std::vector<Amg::Vector2D> axis = { normal * (elementP[1] - elementP[0]),
-                                      normal * (elementP[2] - elementP[1]),
-                                      normal * (elementP[3] - elementP[2]),
-                                      normal * (elementP[4] - elementP[3]),
-                                      normal * (elementP[5] - elementP[4]) };
-  bchk.ComputeKDOP(elementP, axis, elementKDOP);
-  // compute KDOP for error ellipse
-  std::vector<KDOP> errelipseKDOP(5);
-  bchk.ComputeKDOP(bchk.EllipseToPoly(3), axis, errelipseKDOP);
-  // check if KDOPs overlap and return result
-  return bchk.TestKDOPKDOP(elementKDOP, errelipseKDOP);
-}
-
-inline bool
-RotatedDiamondBounds::insideLoc1(const Amg::Vector2D& locpo, double tol1) const
-{
-  return (fabs(locpo[locY]) < m_boundValues[RotatedDiamondBounds::bv_medHalfX] + tol1);
-}
-
-inline bool
-RotatedDiamondBounds::insideLoc2(const Amg::Vector2D& locpo, double tol2) const
-{
-  return ((locpo[locX] > -2. * m_boundValues[RotatedDiamondBounds::bv_halfY1] - tol2) &&
-          (locpo[locY] < 2. * m_boundValues[RotatedDiamondBounds::bv_halfY2] + tol2));
-}
-
-inline void
-RotatedDiamondBounds::initCache()
-{
-  m_alpha1 = atan2(m_boundValues[RotatedDiamondBounds::bv_medHalfX] - m_boundValues[RotatedDiamondBounds::bv_minHalfX],
-                   2. * m_boundValues[RotatedDiamondBounds::bv_halfY1]);
-  m_alpha2 = atan2(m_boundValues[RotatedDiamondBounds::bv_medHalfX] - m_boundValues[RotatedDiamondBounds::bv_maxHalfX],
-                   2. * m_boundValues[RotatedDiamondBounds::bv_halfY2]);
-}
-
 } // end of namespace
 
+#include "TrkSurfaces/RotatedDiamondBounds.icc"
 #endif // TRKSURFACES_DIAMONDBOUNDS_H
diff --git a/Tracking/TrkDetDescr/TrkSurfaces/TrkSurfaces/RotatedDiamondBounds.icc b/Tracking/TrkDetDescr/TrkSurfaces/TrkSurfaces/RotatedDiamondBounds.icc
new file mode 100644
index 0000000000000000000000000000000000000000..4352db7577d9066a4e504adb08d710fb13a3679d
--- /dev/null
+++ b/Tracking/TrkDetDescr/TrkSurfaces/TrkSurfaces/RotatedDiamondBounds.icc
@@ -0,0 +1,163 @@
+namespace Trk {
+
+inline RotatedDiamondBounds*
+RotatedDiamondBounds::clone() const
+{
+  return new RotatedDiamondBounds(*this);
+}
+
+inline double
+RotatedDiamondBounds::minHalflengthX() const
+{
+  return m_boundValues[RotatedDiamondBounds::bv_minHalfX];
+}
+
+inline double
+RotatedDiamondBounds::medHalflengthX() const
+{
+  return m_boundValues[RotatedDiamondBounds::bv_medHalfX];
+}
+
+inline double
+RotatedDiamondBounds::maxHalflengthX() const
+{
+  return m_boundValues[RotatedDiamondBounds::bv_maxHalfX];
+}
+
+inline double
+RotatedDiamondBounds::halflengthY1() const
+{
+  return m_boundValues[RotatedDiamondBounds::bv_halfY1];
+}
+
+inline double
+RotatedDiamondBounds::halflengthY2() const
+{
+  return m_boundValues[RotatedDiamondBounds::bv_halfY2];
+}
+
+inline double
+RotatedDiamondBounds::r() const
+{
+  return sqrt(m_boundValues[RotatedDiamondBounds::bv_medHalfX] *
+                m_boundValues[RotatedDiamondBounds::bv_medHalfX] +
+              m_boundValues[RotatedDiamondBounds::bv_halfY1] *
+                m_boundValues[RotatedDiamondBounds::bv_halfY1]);
+}
+
+inline bool
+RotatedDiamondBounds::inside(const Amg::Vector2D& locpo,
+                             const BoundaryCheck& bchk) const
+{
+  // locX and locY are interchanged wrt DiamondBounds
+  if (bchk.bcType == 0)
+    return RotatedDiamondBounds::inside(
+      locpo, bchk.toleranceLoc1, bchk.toleranceLoc2);
+  // a fast FALSE
+  double max_ell = bchk.lCovariance(0, 0) > bchk.lCovariance(1, 1)
+                     ? bchk.lCovariance(0, 0)
+                     : bchk.lCovariance(1, 1);
+  double limit = bchk.nSigmas * sqrt(max_ell);
+  if (locpo[Trk::locX] <
+      -2 * m_boundValues[RotatedDiamondBounds::bv_halfY1] - limit)
+    return false;
+  if (locpo[Trk::locX] >
+      2 * m_boundValues[RotatedDiamondBounds::bv_halfY2] + limit)
+    return false;
+  // a fast FALSE
+  double fabsX = fabs(locpo[Trk::locY]);
+  if (fabsX > (m_boundValues[RotatedDiamondBounds::bv_medHalfX] + limit))
+    return false;
+  // a fast TRUE
+  double min_ell = bchk.lCovariance(0, 0) < bchk.lCovariance(1, 1)
+                     ? bchk.lCovariance(0, 0)
+                     : bchk.lCovariance(1, 1);
+  limit = bchk.nSigmas * sqrt(min_ell);
+  if (fabsX < (fmin(m_boundValues[RotatedDiamondBounds::bv_minHalfX],
+                    m_boundValues[RotatedDiamondBounds::bv_maxHalfX]) -
+               limit))
+    return true;
+  // a fast TRUE
+  if (fabs(locpo[Trk::locX]) <
+      (fmin(m_boundValues[RotatedDiamondBounds::bv_halfY1],
+            m_boundValues[RotatedDiamondBounds::bv_halfY2]) -
+       limit))
+    return true;
+
+  // compute KDOP and axes for surface polygon
+  std::vector<KDOP> elementKDOP(5);
+  std::vector<Amg::Vector2D> elementP(6);
+  float theta =
+    (bchk.lCovariance(1, 0) != 0 &&
+     (bchk.lCovariance(1, 1) - bchk.lCovariance(0, 0)) != 0)
+      ? .5 * bchk.FastArcTan(2 * bchk.lCovariance(1, 0) /
+                             (bchk.lCovariance(1, 1) - bchk.lCovariance(0, 0)))
+      : 0.;
+  sincosCache scResult = bchk.FastSinCos(theta);
+  AmgMatrix(2, 2) rotMatrix;
+  rotMatrix << scResult.cosC, scResult.sinC, -scResult.sinC, scResult.cosC;
+  AmgMatrix(2, 2) normal;
+  normal << 0, -1, 1, 0;
+  // ellipse is always at (0,0), surface is moved to ellipse position and then
+  // rotated exchange locX and locY
+  Amg::Vector2D locpoF;
+  locpoF[0] = locpo[Trk::locY];
+  locpoF[1] = locpo[Trk::locX];
+  Amg::Vector2D p;
+  p << -m_boundValues[RotatedDiamondBounds::bv_minHalfX],
+    -2. * m_boundValues[RotatedDiamondBounds::bv_halfY1];
+  elementP[0] = (rotMatrix * (p - locpoF));
+  p << -m_boundValues[RotatedDiamondBounds::bv_medHalfX], 0.;
+  elementP[1] = (rotMatrix * (p - locpoF));
+  p << -m_boundValues[RotatedDiamondBounds::bv_maxHalfX],
+    2. * m_boundValues[RotatedDiamondBounds::bv_halfY2];
+  elementP[2] = (rotMatrix * (p - locpoF));
+  p << m_boundValues[RotatedDiamondBounds::bv_maxHalfX],
+    2. * m_boundValues[RotatedDiamondBounds::bv_halfY2];
+  elementP[3] = (rotMatrix * (p - locpoF));
+  p << m_boundValues[RotatedDiamondBounds::bv_medHalfX], 0.;
+  elementP[4] = (rotMatrix * (p - locpoF));
+  p << m_boundValues[RotatedDiamondBounds::bv_minHalfX],
+    -2. * m_boundValues[RotatedDiamondBounds::bv_halfY1];
+  elementP[5] = (rotMatrix * (p - locpoF));
+  std::vector<Amg::Vector2D> axis = { normal * (elementP[1] - elementP[0]),
+                                      normal * (elementP[2] - elementP[1]),
+                                      normal * (elementP[3] - elementP[2]),
+                                      normal * (elementP[4] - elementP[3]),
+                                      normal * (elementP[5] - elementP[4]) };
+  bchk.ComputeKDOP(elementP, axis, elementKDOP);
+  // compute KDOP for error ellipse
+  std::vector<KDOP> errelipseKDOP(5);
+  bchk.ComputeKDOP(bchk.EllipseToPoly(3), axis, errelipseKDOP);
+  // check if KDOPs overlap and return result
+  return bchk.TestKDOPKDOP(elementKDOP, errelipseKDOP);
+}
+
+inline bool
+RotatedDiamondBounds::insideLoc1(const Amg::Vector2D& locpo, double tol1) const
+{
+  return (fabs(locpo[locY]) <
+          m_boundValues[RotatedDiamondBounds::bv_medHalfX] + tol1);
+}
+
+inline bool
+RotatedDiamondBounds::insideLoc2(const Amg::Vector2D& locpo, double tol2) const
+{
+  return (
+    (locpo[locX] >
+     -2. * m_boundValues[RotatedDiamondBounds::bv_halfY1] - tol2) &&
+    (locpo[locY] < 2. * m_boundValues[RotatedDiamondBounds::bv_halfY2] + tol2));
+}
+
+inline void
+RotatedDiamondBounds::initCache()
+{
+  m_alpha1 = atan2(m_boundValues[RotatedDiamondBounds::bv_medHalfX] -
+                     m_boundValues[RotatedDiamondBounds::bv_minHalfX],
+                   2. * m_boundValues[RotatedDiamondBounds::bv_halfY1]);
+  m_alpha2 = atan2(m_boundValues[RotatedDiamondBounds::bv_medHalfX] -
+                     m_boundValues[RotatedDiamondBounds::bv_maxHalfX],
+                   2. * m_boundValues[RotatedDiamondBounds::bv_halfY2]);
+}
+
+} // end of namespace
diff --git a/Tracking/TrkDetDescr/TrkSurfaces/TrkSurfaces/RotatedTrapezoidBounds.h b/Tracking/TrkDetDescr/TrkSurfaces/TrkSurfaces/RotatedTrapezoidBounds.h
index 76bba568e297e6825ac7b8cc882d9a665f62348e..0b91054b395a5f0daaac3ce96ff4c72cd9cf66f9 100644
--- a/Tracking/TrkDetDescr/TrkSurfaces/TrkSurfaces/RotatedTrapezoidBounds.h
+++ b/Tracking/TrkDetDescr/TrkSurfaces/TrkSurfaces/RotatedTrapezoidBounds.h
@@ -41,7 +41,7 @@ namespace Trk {
  @author Andreas.Salzburger@cern.ch
  */
 
-class RotatedTrapezoidBounds : public SurfaceBounds
+class RotatedTrapezoidBounds final: public SurfaceBounds
 {
 
 public:
@@ -153,114 +153,7 @@ private:
   TDD_real_t m_kappa;
   TDD_real_t m_delta;
 };
-
-inline RotatedTrapezoidBounds*
-RotatedTrapezoidBounds::clone() const
-{
-  return new RotatedTrapezoidBounds(*this);
-}
-
-inline void
-RotatedTrapezoidBounds::initCache()
-{
-  m_kappa = 0.5 *
-            (m_boundValues[RotatedTrapezoidBounds::bv_maxHalfY] - m_boundValues[RotatedTrapezoidBounds::bv_minHalfY]) /
-            m_boundValues[RotatedTrapezoidBounds::bv_halfX];
-  m_delta =
-    0.5 * (m_boundValues[RotatedTrapezoidBounds::bv_minHalfY] + m_boundValues[RotatedTrapezoidBounds::bv_maxHalfY]);
-}
-
-inline double
-RotatedTrapezoidBounds::halflengthX() const
-{
-  return m_boundValues[RotatedTrapezoidBounds::bv_halfX];
-}
-
-inline double
-RotatedTrapezoidBounds::minHalflengthY() const
-{
-  return m_boundValues[RotatedTrapezoidBounds::bv_minHalfY];
-}
-
-inline double
-RotatedTrapezoidBounds::maxHalflengthY() const
-{
-  return m_boundValues[RotatedTrapezoidBounds::bv_maxHalfY];
-}
-
-inline double
-RotatedTrapezoidBounds::r() const
-{
-  return sqrt(m_boundValues[RotatedTrapezoidBounds::bv_halfX] * m_boundValues[RotatedTrapezoidBounds::bv_halfX] +
-              m_boundValues[RotatedTrapezoidBounds::bv_maxHalfY] * m_boundValues[RotatedTrapezoidBounds::bv_maxHalfY]);
-}
-
-inline bool
-RotatedTrapezoidBounds::inside(const Amg::Vector2D& locpo, const BoundaryCheck& bchk) const
-{
-  if (bchk.bcType == 0)
-    return RotatedTrapezoidBounds::inside(locpo, bchk.toleranceLoc1, bchk.toleranceLoc2);
-
-  // a fast FALSE
-  double fabsX = fabs(locpo[Trk::locX]);
-  double max_ell = bchk.lCovariance(0, 0) > bchk.lCovariance(1, 1) ? bchk.lCovariance(0, 0) : bchk.lCovariance(1, 1);
-  double limit = bchk.nSigmas * sqrt(max_ell);
-  if (fabsX > (m_boundValues[RotatedTrapezoidBounds::bv_halfX] + limit))
-    return false;
-  // a fast FALSE
-  double fabsY = fabs(locpo[Trk::locY]);
-  if (fabsY > (m_boundValues[RotatedTrapezoidBounds::bv_maxHalfY] + limit))
-    return false;
-  // a fast TRUE
-  double min_ell = bchk.lCovariance(0, 0) < bchk.lCovariance(1, 1) ? bchk.lCovariance(0, 0) : bchk.lCovariance(1, 1);
-  limit = bchk.nSigmas * sqrt(min_ell);
-  if (fabsY < (m_boundValues[RotatedTrapezoidBounds::bv_minHalfY] + limit) &&
-      fabsX < (m_boundValues[RotatedTrapezoidBounds::bv_halfX] + limit))
-    return true;
-
-  // compute KDOP and axes for surface polygon
-  std::vector<KDOP> elementKDOP(3);
-  std::vector<Amg::Vector2D> elementP(4);
-  float theta = (bchk.lCovariance(1, 0) != 0 && (bchk.lCovariance(1, 1) - bchk.lCovariance(0, 0)) != 0)
-                  ? .5 * bchk.FastArcTan(2 * bchk.lCovariance(1, 0) / (bchk.lCovariance(1, 1) - bchk.lCovariance(0, 0)))
-                  : 0.;
-  sincosCache scResult = bchk.FastSinCos(theta);
-  AmgMatrix(2, 2) rotMatrix;
-  rotMatrix << scResult.cosC, scResult.sinC, -scResult.sinC, scResult.cosC;
-  AmgMatrix(2, 2) normal;
-  normal << 0, -1, 1, 0;
-  // ellipse is always at (0,0), surface is moved to ellipse position and then rotated
-  Amg::Vector2D p;
-  p << -m_boundValues[RotatedTrapezoidBounds::bv_halfX], m_boundValues[RotatedTrapezoidBounds::bv_minHalfY];
-  elementP[0] = (rotMatrix * (p - locpo));
-  p << -m_boundValues[RotatedTrapezoidBounds::bv_halfX], -m_boundValues[RotatedTrapezoidBounds::bv_minHalfY];
-  elementP[1] = (rotMatrix * (p - locpo));
-  p << m_boundValues[RotatedTrapezoidBounds::bv_halfX], m_boundValues[RotatedTrapezoidBounds::bv_maxHalfY];
-  elementP[2] = (rotMatrix * (p - locpo));
-  p << m_boundValues[RotatedTrapezoidBounds::bv_halfX], -m_boundValues[RotatedTrapezoidBounds::bv_maxHalfY];
-  elementP[3] = (rotMatrix * (p - locpo));
-  std::vector<Amg::Vector2D> axis = { normal * (elementP[1] - elementP[0]),
-                                      normal * (elementP[3] - elementP[1]),
-                                      normal * (elementP[2] - elementP[0]) };
-  bchk.ComputeKDOP(elementP, axis, elementKDOP);
-  // compute KDOP for error ellipse
-  std::vector<KDOP> errelipseKDOP(3);
-  bchk.ComputeKDOP(bchk.EllipseToPoly(3), axis, errelipseKDOP);
-  // check if KDOPs overlap and return result
-  return bchk.TestKDOPKDOP(elementKDOP, errelipseKDOP);
-}
-
-inline bool
-RotatedTrapezoidBounds::insideLoc1(const Amg::Vector2D& locpo, double tol1) const
-{
-  return (fabs(locpo[locX]) < m_boundValues[RotatedTrapezoidBounds::bv_halfX] + tol1);
-}
-
-inline bool
-RotatedTrapezoidBounds::insideLoc2(const Amg::Vector2D& locpo, double tol2) const
-{
-  return (fabs(locpo[locY]) < m_boundValues[RotatedTrapezoidBounds::bv_maxHalfY] + tol2);
-}
 } // end of namespace
 
+#include "TrkSurfaces/RotatedTrapezoidBounds.icc"
 #endif // TRKSURFACES_TRAPEZOIDBOUNDS_H
diff --git a/Tracking/TrkDetDescr/TrkSurfaces/TrkSurfaces/RotatedTrapezoidBounds.icc b/Tracking/TrkDetDescr/TrkSurfaces/TrkSurfaces/RotatedTrapezoidBounds.icc
new file mode 100644
index 0000000000000000000000000000000000000000..5ad28fa833b8ff3a8d925565526c2d388b44e58b
--- /dev/null
+++ b/Tracking/TrkDetDescr/TrkSurfaces/TrkSurfaces/RotatedTrapezoidBounds.icc
@@ -0,0 +1,135 @@
+/*
+  Copyright (C) 2002-2020 CERN for the benefit of the ATLAS collaboration
+*/
+
+namespace Trk {
+inline RotatedTrapezoidBounds*
+RotatedTrapezoidBounds::clone() const
+{
+  return new RotatedTrapezoidBounds(*this);
+}
+
+inline void
+RotatedTrapezoidBounds::initCache()
+{
+  m_kappa = 0.5 *
+            (m_boundValues[RotatedTrapezoidBounds::bv_maxHalfY] -
+             m_boundValues[RotatedTrapezoidBounds::bv_minHalfY]) /
+            m_boundValues[RotatedTrapezoidBounds::bv_halfX];
+  m_delta = 0.5 * (m_boundValues[RotatedTrapezoidBounds::bv_minHalfY] +
+                   m_boundValues[RotatedTrapezoidBounds::bv_maxHalfY]);
+}
+
+inline double
+RotatedTrapezoidBounds::halflengthX() const
+{
+  return m_boundValues[RotatedTrapezoidBounds::bv_halfX];
+}
+
+inline double
+RotatedTrapezoidBounds::minHalflengthY() const
+{
+  return m_boundValues[RotatedTrapezoidBounds::bv_minHalfY];
+}
+
+inline double
+RotatedTrapezoidBounds::maxHalflengthY() const
+{
+  return m_boundValues[RotatedTrapezoidBounds::bv_maxHalfY];
+}
+
+inline double
+RotatedTrapezoidBounds::r() const
+{
+  return sqrt(m_boundValues[RotatedTrapezoidBounds::bv_halfX] *
+                m_boundValues[RotatedTrapezoidBounds::bv_halfX] +
+              m_boundValues[RotatedTrapezoidBounds::bv_maxHalfY] *
+                m_boundValues[RotatedTrapezoidBounds::bv_maxHalfY]);
+}
+
+inline bool
+RotatedTrapezoidBounds::inside(const Amg::Vector2D& locpo,
+                               const BoundaryCheck& bchk) const
+{
+  if (bchk.bcType == 0)
+    return RotatedTrapezoidBounds::inside(
+      locpo, bchk.toleranceLoc1, bchk.toleranceLoc2);
+
+  // a fast FALSE
+  double fabsX = fabs(locpo[Trk::locX]);
+  double max_ell = bchk.lCovariance(0, 0) > bchk.lCovariance(1, 1)
+                     ? bchk.lCovariance(0, 0)
+                     : bchk.lCovariance(1, 1);
+  double limit = bchk.nSigmas * sqrt(max_ell);
+  if (fabsX > (m_boundValues[RotatedTrapezoidBounds::bv_halfX] + limit))
+    return false;
+  // a fast FALSE
+  double fabsY = fabs(locpo[Trk::locY]);
+  if (fabsY > (m_boundValues[RotatedTrapezoidBounds::bv_maxHalfY] + limit))
+    return false;
+  // a fast TRUE
+  double min_ell = bchk.lCovariance(0, 0) < bchk.lCovariance(1, 1)
+                     ? bchk.lCovariance(0, 0)
+                     : bchk.lCovariance(1, 1);
+  limit = bchk.nSigmas * sqrt(min_ell);
+  if (fabsY < (m_boundValues[RotatedTrapezoidBounds::bv_minHalfY] + limit) &&
+      fabsX < (m_boundValues[RotatedTrapezoidBounds::bv_halfX] + limit))
+    return true;
+
+  // compute KDOP and axes for surface polygon
+  std::vector<KDOP> elementKDOP(3);
+  std::vector<Amg::Vector2D> elementP(4);
+  float theta =
+    (bchk.lCovariance(1, 0) != 0 &&
+     (bchk.lCovariance(1, 1) - bchk.lCovariance(0, 0)) != 0)
+      ? .5 * bchk.FastArcTan(2 * bchk.lCovariance(1, 0) /
+                             (bchk.lCovariance(1, 1) - bchk.lCovariance(0, 0)))
+      : 0.;
+  sincosCache scResult = bchk.FastSinCos(theta);
+  AmgMatrix(2, 2) rotMatrix;
+  rotMatrix << scResult.cosC, scResult.sinC, -scResult.sinC, scResult.cosC;
+  AmgMatrix(2, 2) normal;
+  normal << 0, -1, 1, 0;
+  // ellipse is always at (0,0), surface is moved to ellipse position and then
+  // rotated
+  Amg::Vector2D p;
+  p << -m_boundValues[RotatedTrapezoidBounds::bv_halfX],
+    m_boundValues[RotatedTrapezoidBounds::bv_minHalfY];
+  elementP[0] = (rotMatrix * (p - locpo));
+  p << -m_boundValues[RotatedTrapezoidBounds::bv_halfX],
+    -m_boundValues[RotatedTrapezoidBounds::bv_minHalfY];
+  elementP[1] = (rotMatrix * (p - locpo));
+  p << m_boundValues[RotatedTrapezoidBounds::bv_halfX],
+    m_boundValues[RotatedTrapezoidBounds::bv_maxHalfY];
+  elementP[2] = (rotMatrix * (p - locpo));
+  p << m_boundValues[RotatedTrapezoidBounds::bv_halfX],
+    -m_boundValues[RotatedTrapezoidBounds::bv_maxHalfY];
+  elementP[3] = (rotMatrix * (p - locpo));
+  std::vector<Amg::Vector2D> axis = { normal * (elementP[1] - elementP[0]),
+                                      normal * (elementP[3] - elementP[1]),
+                                      normal * (elementP[2] - elementP[0]) };
+  bchk.ComputeKDOP(elementP, axis, elementKDOP);
+  // compute KDOP for error ellipse
+  std::vector<KDOP> errelipseKDOP(3);
+  bchk.ComputeKDOP(bchk.EllipseToPoly(3), axis, errelipseKDOP);
+  // check if KDOPs overlap and return result
+  return bchk.TestKDOPKDOP(elementKDOP, errelipseKDOP);
+}
+
+inline bool
+RotatedTrapezoidBounds::insideLoc1(const Amg::Vector2D& locpo,
+                                   double tol1) const
+{
+  return (fabs(locpo[locX]) <
+          m_boundValues[RotatedTrapezoidBounds::bv_halfX] + tol1);
+}
+
+inline bool
+RotatedTrapezoidBounds::insideLoc2(const Amg::Vector2D& locpo,
+                                   double tol2) const
+{
+  return (fabs(locpo[locY]) <
+          m_boundValues[RotatedTrapezoidBounds::bv_maxHalfY] + tol2);
+}
+} // end of namespace
+