@@ -1073,6 +1073,7 @@ namespace OrthoTree
10731073
10741074 // Adaptors
10751075
1076+ // Provides basic accessor and mutator methods for generic geometric types (points, boxes, rays, planes).
10761077 template <dim_t DIMENSION_NO, typename TVector, typename TBox, typename TRay, typename TPlane, typename TGeometry = double >
10771078 struct AdaptorGeneralBasics
10781079 {
@@ -1091,7 +1092,7 @@ namespace OrthoTree
10911092 static inline constexpr TGeometry GetPlaneOrigoDistance (TPlane const & plane) noexcept { return plane.OrigoDistance ; }
10921093 };
10931094
1094-
1095+ // Provides general vector/box/ray/plane operations based on a basic adaptor interface. If the geometric types are connected to an BLAS, it is recommended to implement a custom AdaptorGeneral.
10951096 template <dim_t DIMENSION_NO, typename TVector, typename TBox, typename TRay, typename TPlane, typename TGeometry, typename TAdaptorBasics>
10961097 struct AdaptorGeneralBase : TAdaptorBasics
10971098 {
@@ -1237,65 +1238,52 @@ namespace OrthoTree
12371238 }
12381239 }
12391240
1240-
1241- static constexpr std::optional<double > GetRayBoxDistance (TBox const & box, TVector const & rayBasePoint, TVector const & rayHeading, TGeometry tolerance) noexcept
1241+ static constexpr std::optional<double > GetRayBoxDistance (TBox const & box, TVector const & rayOrigin, TVector const & rayDirection, TGeometry tolerance) noexcept
12421242 {
1243- if (DoesBoxContainPoint (box, rayBasePoint, tolerance))
1243+ assert (tolerance >= 0 && " Tolerance cannot be negative!" );
1244+
1245+ if (DoesBoxContainPoint (box, rayOrigin, tolerance))
12441246 return 0.0 ;
12451247
12461248 auto constexpr inf = std::numeric_limits<double >::max ();
12471249
1248- auto minBoxDistances = std::array< double , DIMENSION_NO>{} ;
1249- auto maxBoxDistances = std::array< double , DIMENSION_NO>{} ;
1250+ double minBoxDistance = -inf ;
1251+ double maxBoxDistance = +inf ;
12501252 for (dim_t dimensionID = 0 ; dimensionID < DIMENSION_NO; ++dimensionID)
12511253 {
1252- auto const dirComp = Base::GetPointC (rayHeading, dimensionID);
1253- if (dirComp == 0 )
1254+ auto const origin = Base::GetPointC (rayOrigin, dimensionID);
1255+ auto const direction = Base::GetPointC (rayDirection, dimensionID);
1256+ auto const boxMin = Base::GetBoxMinC (box, dimensionID) - tolerance;
1257+ auto const boxMax = Base::GetBoxMaxC (box, dimensionID) + tolerance;
1258+
1259+ if (direction == 0 )
12541260 {
12551261 if (tolerance != 0.0 )
12561262 {
12571263 // Box should be within tolerance (<, not <=)
1258-
1259- assert (tolerance > 0 );
1260- if (Base::GetBoxMaxC (box, dimensionID) + tolerance <= Base::GetPointC (rayBasePoint, dimensionID))
1261- return std::nullopt ;
1262-
1263- if (Base::GetBoxMinC (box, dimensionID) - tolerance >= Base::GetPointC (rayBasePoint, dimensionID))
1264+ if (origin <= boxMin || boxMax <= origin)
12641265 return std::nullopt ;
12651266 }
12661267 else
12671268 {
1268- if (Base::GetBoxMaxC (box, dimensionID) < Base::GetPointC (rayBasePoint, dimensionID))
1269- return std::nullopt ;
1270-
1271- if (Base::GetBoxMinC (box, dimensionID) > Base::GetPointC (rayBasePoint, dimensionID))
1269+ if (origin < boxMin || boxMax < origin)
12721270 return std::nullopt ;
12731271 }
1274-
1275- minBoxDistances[dimensionID] = -inf;
1276- maxBoxDistances[dimensionID] = +inf;
12771272 }
12781273 else
12791274 {
1280- auto const minBox = Base::GetBoxMinC (box, dimensionID) - tolerance;
1281- auto const maxBox = Base::GetBoxMaxC (box, dimensionID) + tolerance;
1282- auto const pointComp = Base::GetPointC (rayBasePoint, dimensionID);
1283- auto const dirCompRecip = 1.0 / dirComp;
1284- if (dirComp < 0.0 )
1285- {
1286- minBoxDistances[dimensionID] = (maxBox - pointComp) * dirCompRecip;
1287- maxBoxDistances[dimensionID] = (minBox - pointComp) * dirCompRecip;
1288- }
1289- else
1290- {
1291- minBoxDistances[dimensionID] = (minBox - pointComp) * dirCompRecip;
1292- maxBoxDistances[dimensionID] = (maxBox - pointComp) * dirCompRecip;
1293- }
1275+ double const directionReciprocal = 1.0 / direction;
1276+ double t1 = (boxMin - origin) * directionReciprocal;
1277+ double t2 = (boxMax - origin) * directionReciprocal;
1278+ if (t1 > t2)
1279+ std::swap (t1, t2);
1280+
1281+ minBoxDistance = std::max (minBoxDistance, t1);
1282+ maxBoxDistance = std::min (maxBoxDistance, t2);
12941283 }
12951284 }
12961285
1297- auto const minBoxDistance = *std::max_element (minBoxDistances.begin (), minBoxDistances.end ());
1298- auto const maxBoxDistance = *std::min_element (maxBoxDistances.begin (), maxBoxDistances.end ());
1286+ assert (maxBoxDistance != inf && " rayDirection is a zero vector!" );
12991287 if (minBoxDistance > maxBoxDistance || maxBoxDistance < 0.0 )
13001288 return std::nullopt ;
13011289 else
@@ -1751,64 +1739,51 @@ namespace OrthoTree
17511739 }
17521740
17531741 static inline constexpr std::optional<Geometry> GetRayBoxDistanceAD (
1754- Vector const & center, Vector const & halfSizes, TVector const & rayBasePoint , TVector const & rayHeading , TGeometry tolerance) noexcept
1742+ Vector const & center, Vector const & halfSizes, TVector const & rayOrigin , TVector const & rayDirection , TGeometry tolerance) noexcept
17551743 {
1756- if (DoesBoxContainPointAD (center, halfSizes, rayBasePoint, tolerance))
1744+ assert (tolerance >= 0 && " Tolerance cannot be negative!" );
1745+ if (DoesBoxContainPointAD (center, halfSizes, rayOrigin, tolerance))
17571746 return Geometry{};
17581747
17591748 auto constexpr inf = std::numeric_limits<Geometry>::max ();
1749+ auto minBoxDistance = -inf;
1750+ auto maxBoxDistance = +inf;
1751+ auto const tolerance_ = Geometry (tolerance);
17601752
1761- auto minBoxDistances = Vector{};
1762- auto maxBoxDistances = Vector{};
17631753 for (dim_t dimensionID = 0 ; dimensionID < DIMENSION_NO; ++dimensionID)
17641754 {
1765- auto const dirComp = Geometry (AD::GetPointC (rayHeading, dimensionID));
1766- auto const minBox = center[dimensionID] - halfSizes[dimensionID] - Geometry (tolerance);
1767- auto const maxBox = center[dimensionID] + halfSizes[dimensionID] + Geometry (tolerance);
1768- if (dirComp == 0 )
1755+ auto const origin = Geometry (AD::GetPointC (rayOrigin, dimensionID));
1756+ auto const direction = Geometry (AD::GetPointC (rayDirection, dimensionID));
1757+ auto const boxMin = center[dimensionID] - halfSizes[dimensionID] - tolerance_;
1758+ auto const boxMax = center[dimensionID] + halfSizes[dimensionID] + tolerance_;
1759+ if (direction == 0 )
17691760 {
17701761 if (tolerance != 0.0 )
17711762 {
17721763 // Box should be within tolerance (<, not <=)
1773-
1774- assert (tolerance > 0 );
1775- if (maxBox <= AD::GetPointC (rayBasePoint, dimensionID))
1776- return std::nullopt ;
1777-
1778- if (minBox >= AD::GetPointC (rayBasePoint, dimensionID))
1764+ if (origin <= boxMin || boxMax <= origin)
17791765 return std::nullopt ;
17801766 }
17811767 else
17821768 {
1783- if (maxBox < AD::GetPointC (rayBasePoint, dimensionID))
1784- return std::nullopt ;
1785-
1786- if (minBox > AD::GetPointC (rayBasePoint, dimensionID))
1769+ if (origin < boxMin || boxMax < origin)
17871770 return std::nullopt ;
17881771 }
1789-
1790- minBoxDistances[dimensionID] = -inf;
1791- maxBoxDistances[dimensionID] = +inf;
17921772 }
17931773 else
17941774 {
1795- auto const pointComp = Geometry (AD::GetPointC (rayBasePoint, dimensionID));
1796- auto const dirCompRecip = Geometry (1.0 ) / dirComp;
1797- if (dirComp < Geometry{})
1798- {
1799- minBoxDistances[dimensionID] = (maxBox - pointComp) * dirCompRecip;
1800- maxBoxDistances[dimensionID] = (minBox - pointComp) * dirCompRecip;
1801- }
1802- else
1803- {
1804- minBoxDistances[dimensionID] = (minBox - pointComp) * dirCompRecip;
1805- maxBoxDistances[dimensionID] = (maxBox - pointComp) * dirCompRecip;
1806- }
1775+ auto const directionReciprocal = Geometry (1 ) / direction;
1776+ auto t1 = (boxMin - origin) * directionReciprocal;
1777+ auto t2 = (boxMax - origin) * directionReciprocal;
1778+ if (t1 > t2)
1779+ std::swap (t1, t2);
1780+
1781+ minBoxDistance = std::max (minBoxDistance, t1);
1782+ maxBoxDistance = std::min (maxBoxDistance, t2);
18071783 }
18081784 }
18091785
1810- auto const minBoxDistance = *std::max_element (minBoxDistances.begin (), minBoxDistances.end ());
1811- auto const maxBoxDistance = *std::min_element (maxBoxDistances.begin (), maxBoxDistances.end ());
1786+ assert (maxBoxDistance != inf && " rayDirection is a zero vector!" );
18121787 if (minBoxDistance > maxBoxDistance || maxBoxDistance < 0.0 )
18131788 return std::nullopt ;
18141789 else
0 commit comments