From 06ebf17a8d0eb2a23291ab584accd7b917474aba Mon Sep 17 00:00:00 2001
From: Menelaos Karavelas <menelaos.karavelas@gmail.com>
Date: Thu, 8 May 2014 22:24:01 +0300
Subject: [PATCH] [disjoint] add alternative (more robust) implementation of
 segment-box disjointness test; the alternative implementation does not use
 any divisions

---
 .../detail/disjoint/segment_box.hpp           | 263 ++++++++++++++++++
 1 file changed, 263 insertions(+)
 create mode 100644 include/boost/geometry/algorithms/detail/disjoint/segment_box.hpp

diff --git a/include/boost/geometry/algorithms/detail/disjoint/segment_box.hpp b/include/boost/geometry/algorithms/detail/disjoint/segment_box.hpp
new file mode 100644
index 000000000..f4b150826
--- /dev/null
+++ b/include/boost/geometry/algorithms/detail/disjoint/segment_box.hpp
@@ -0,0 +1,263 @@
+// Boost.Geometry (aka GGL, Generic Geometry Library)
+
+// Copyright (c) 2014, Oracle and/or its affiliates.
+
+// Contributed and/or modified by Menelaos Karavelas, on behalf of Oracle
+
+// Licensed under the Boost Software License version 1.0.
+// http://www.boost.org/users/license.html
+
+#ifndef BOOST_GEOMETRY_ALGORITHMS_DETAIL_DISJOINT_SEGMENT_BOX_HPP
+#define BOOST_GEOMETRY_ALGORITHMS_DETAIL_DISJOINT_SEGMENT_BOX_HPP
+
+#include <cstddef>
+#include <utility>
+
+#include <boost/geometry/util/calculation_type.hpp>
+
+#include <boost/geometry/core/access.hpp>
+#include <boost/geometry/core/coordinate_dimension.hpp>
+#include <boost/geometry/core/point_type.hpp>
+
+#include <boost/geometry/algorithms/detail/assign_indexed_point.hpp>
+
+
+namespace boost { namespace geometry
+{
+
+
+#ifndef DOXYGEN_NO_DETAIL
+namespace detail { namespace disjoint
+{
+
+
+
+template
+<
+    typename RelativeDistance,
+    typename SegmentPoint,
+    typename Box,
+    std::size_t I,
+    std::size_t Dimension
+>
+struct disjoint_segment_box_impl
+{
+    template <typename RelativeDistancePair>
+    static inline bool apply(SegmentPoint const& p0,
+                             SegmentPoint const& p1,
+                             Box const& box,
+                             RelativeDistancePair& t_min,
+                             RelativeDistancePair& t_max)
+    {
+        typedef typename coordinate_type<Box>::type box_coordinate_type;
+        typedef typename coordinate_type
+            <
+                SegmentPoint
+            >::type point_coordinate_type;
+
+        RelativeDistance c_p0 = boost::numeric_cast
+            <
+                point_coordinate_type
+            >( geometry::get<I>(p0) );
+
+        RelativeDistance c_p1 = boost::numeric_cast
+            <
+                point_coordinate_type
+            >( geometry::get<I>(p1) );
+
+        RelativeDistance c_b_min = boost::numeric_cast
+            <
+                box_coordinate_type
+            >( geometry::get<geometry::min_corner, I>(box) );
+
+        RelativeDistance c_b_max = boost::numeric_cast
+            <
+                box_coordinate_type
+            >( geometry::get<geometry::max_corner, I>(box) );
+
+        RelativeDistance ti_min, ti_max, diff;
+
+        if ( geometry::get<I>(p1) >= geometry::get<I>(p0) )
+        {
+            diff = c_p1 - c_p0;
+            ti_min = c_b_min - c_p0;
+            ti_max = c_b_max - c_p0;
+        }
+        else
+        {
+            diff = c_p0 - c_p1;
+            ti_min = c_p0 - c_b_max;
+            ti_max = c_p0 - c_b_min;
+        }
+
+        RelativeDistance t_min_x_diff = t_min.first * diff;
+        RelativeDistance t_max_x_diff = t_max.first * diff;
+
+        if ( t_min_x_diff > ti_max * t_min.second
+             || t_max_x_diff < ti_min * t_max.second )
+        {
+            return true;
+        }
+
+        if ( ti_min * t_min.second > t_min_x_diff )
+        {
+            t_min.first = ti_min;
+            t_min.second = diff;
+        }
+        if ( ti_max * t_max.second < t_max_x_diff )
+        {
+            t_max.first = ti_max;
+            t_max.second = diff;
+        }
+
+        if ( t_min.first > t_min.second || t_max.first < 0 )
+        {
+            return true;
+        }
+
+        return disjoint_segment_box_impl
+            <
+                RelativeDistance,
+                SegmentPoint,
+                Box, 
+                I + 1,
+                Dimension
+            >::apply(p0, p1, box, t_min, t_max);
+    }
+};
+
+
+template
+<
+    typename RelativeDistance,
+    typename SegmentPoint,
+    typename Box,
+    std::size_t Dimension
+>
+struct disjoint_segment_box_impl
+    <
+        RelativeDistance, SegmentPoint, Box, 0, Dimension
+    >
+{
+    static inline bool apply(SegmentPoint const& p0,
+                             SegmentPoint const& p1,
+                             Box const& box)
+    {
+        typedef typename coordinate_type<Box>::type box_coordinate_type;
+        typedef typename coordinate_type
+            <
+                SegmentPoint
+            >::type point_coordinate_type;
+
+        RelativeDistance c_p0 = boost::numeric_cast
+            <
+                point_coordinate_type
+            >( geometry::get<0>(p0) );
+
+        RelativeDistance c_p1 = boost::numeric_cast
+            <
+                point_coordinate_type
+            >( geometry::get<0>(p1) );
+
+        RelativeDistance c_b_min = boost::numeric_cast
+            <
+                box_coordinate_type
+            >( geometry::get<geometry::min_corner, 0>(box) );
+
+        RelativeDistance c_b_max = boost::numeric_cast
+            <
+                box_coordinate_type
+            >( geometry::get<geometry::max_corner, 0>(box) );
+
+        std::pair<RelativeDistance, RelativeDistance> t_min, t_max;
+
+        RelativeDistance diff;
+
+        if ( geometry::get<0>(p1) >= geometry::get<0>(p0) )
+        {
+            diff = c_p1 - c_p0;
+            t_min.first = c_b_min - c_p0;
+            t_max.first = c_b_max - c_p0;
+        }
+        else
+        {
+            diff = c_p0 - c_p1;
+            t_min.first = c_p0 - c_b_max;
+            t_max.first = c_p0 - c_b_min;
+        }
+
+        if ( t_min.first > diff || t_max.first < 0 )
+        {
+            return true;
+        }
+
+        t_min.second = t_max.second = diff;
+
+        return disjoint_segment_box_impl
+            <
+                RelativeDistance, SegmentPoint, Box, 1, Dimension
+            >::apply(p0, p1, box, t_min, t_max);
+    }
+};
+
+
+template
+<
+    typename RelativeDistance,
+    typename SegmentPoint,
+    typename Box,
+    std::size_t Dimension
+>
+struct disjoint_segment_box_impl
+    <
+        RelativeDistance, SegmentPoint, Box, Dimension, Dimension
+    >
+{
+    template <typename RelativeDistancePair>
+    static inline bool apply(SegmentPoint const&, SegmentPoint const&,
+                             Box const&,
+                             RelativeDistancePair&, RelativeDistancePair&)
+    {
+        return false;
+    }
+};
+
+
+//=========================================================================
+
+
+template <typename Segment, typename Box>
+struct disjoint_segment_box_no_numeric_limits
+{ 
+   static inline bool apply(Segment const& segment, Box const& box)
+    {
+        assert_dimension_equal<Segment, Box>();
+
+        typedef typename util::calculation_type::geometric::binary
+            <
+                Segment, Box, void
+            >::type relative_distance_type;
+
+        typedef typename point_type<Segment>::type segment_point_type;
+        segment_point_type p0, p1;
+        geometry::detail::assign_point_from_index<0>(segment, p0);
+        geometry::detail::assign_point_from_index<1>(segment, p1);
+
+        return disjoint_segment_box_impl
+            <
+                relative_distance_type, segment_point_type, Box,
+                0, dimension<Box>::value
+            >::apply(p0, p1, box);
+    }
+};
+
+
+}} // namespace detail::disjoint
+#endif // DOXYGEN_NO_DETAIL
+
+
+
+}} // namespace boost::geometry
+
+
+#endif // BOOST_GEOMETRY_ALGORITHMS_DETAIL_DISJOINT_SEGMENT_BOX_HPP