Merge branch 'develop' of https://github.com/boostorg/geometry into feature/envelope_and_expand_for_spherical_cs

doc/release_notes.qbk

@@ -27,6 +27,8 @@
 * Added rtree const_iterator, begin(), end() and the support for Boost.Range.
 * The support for C++11 `std::initializer_list` in geometries models.
 * Disjoint and intersects support the following geometry combinations: multipoint/linestring, multipoint/multilinestring
+* Intersection has been implemented for combinations of pointlike and linear geometries
+* Added implementation for difference(pointlike, linear)
 
 [*Improvements]
 

extensions/test/gis/projections/projections_static.cpp

@@ -225,7 +225,6 @@ void test_all()
     test_forward<bg::projections::eck5_spheroid>(amsterdam, utrecht, "+ellps=WGS84 +units=m", 2);
 
     test_forward<bg::projections::eck6_spheroid>(amsterdam, utrecht,  "+ellps=WGS84 +units=m");
-    test_forward<bg::projections::eck6_ellipsoid>(amsterdam, utrecht,  "+ellps=WGS84 +units=m", 18);
 
     test_forward<bg::projections::eqc_spheroid>(amsterdam, utrecht,  "+ellps=WGS84 +units=m", 5);
     test_forward<bg::projections::eqdc_ellipsoid>(amsterdam, utrecht,  "+ellps=WGS84 +units=m +lat_1=60 +lat_2=0");
@@ -241,7 +240,6 @@ void test_all()
     test_forward<bg::projections::gins8_spheroid>(amsterdam, utrecht, "+ellps=WGS84 +units=m", 7);
 
     test_forward<bg::projections::gn_sinu_spheroid>(amsterdam, utrecht, "+ellps=WGS84 +units=m +m=0.5 +n=1.785");
-    test_forward<bg::projections::gn_sinu_ellipsoid>(amsterdam, utrecht, "+ellps=WGS84 +units=m +m=0.5 +n=1.785", 18);
 
     test_forward<bg::projections::gnom_spheroid>(amsterdam, utrecht, "+ellps=WGS84 +units=m", 50);
     test_forward<bg::projections::goode_spheroid>(amsterdam, utrecht, "+ellps=WGS84 +units=m");
@@ -267,7 +265,6 @@ void test_all()
     test_forward<bg::projections::mbtfpp_spheroid>(amsterdam, utrecht, "+ellps=WGS84 +units=m");
     test_forward<bg::projections::mbtfpq_spheroid>(amsterdam, utrecht, "+ellps=WGS84 +units=m", 2);
 
-    test_forward<bg::projections::mbtfps_ellipsoid>(amsterdam, utrecht, "+ellps=WGS84 +units=m", 18);
     test_forward<bg::projections::mbtfps_spheroid>(amsterdam, utrecht, "+ellps=WGS84 +units=m");
 
     test_forward<bg::projections::merc_ellipsoid>(amsterdam, utrecht, "+ellps=WGS84 +units=m", 22);

include/boost/geometry/algorithms/detail/buffer/buffer_inserter.hpp

@@ -383,6 +383,7 @@ inline void buffer_point(Point const& point, Collection& collection,
     std::vector<OutputPointType> range_out;
     point_strategy.apply(point, distance_strategy, range_out);
     collection.add_piece(strategy::buffer::buffered_point, range_out, false);
+    collection.set_piece_center(point);
     collection.finish_ring();
 }
 

include/boost/geometry/algorithms/detail/buffer/buffer_policies.hpp

@@ -9,6 +9,7 @@
 #ifndef BOOST_GEOMETRY_ALGORITHMS_DETAIL_BUFFER_BUFFER_POLICIES_HPP
 #define BOOST_GEOMETRY_ALGORITHMS_DETAIL_BUFFER_BUFFER_POLICIES_HPP
 
+#define BOOST_GEOMETRY_BUFFER_USE_SIDE_OF_INTERSECTION
 
 #include <cstddef>
 
@@ -122,6 +123,10 @@ struct buffer_turn_info
 
     intersection_location_type location;
 
+#if defined(BOOST_GEOMETRY_BUFFER_USE_SIDE_OF_INTERSECTION)
+    robust_point_type rob_pi, rob_pj, rob_qi, rob_qj;
+#endif
+
     int count_within;
 
     bool within_original;
@@ -130,7 +135,9 @@ struct buffer_turn_info
 
     int count_on_offsetted;
     int count_on_helper;
+#if ! defined(BOOST_GEOMETRY_BUFFER_USE_SIDE_OF_INTERSECTION)
     int count_within_near_offsetted;
+#endif
 
     bool remove_on_multi;
 
@@ -147,7 +154,9 @@ struct buffer_turn_info
         , count_in_original(0)
         , count_on_offsetted(0)
         , count_on_helper(0)
+#if ! defined(BOOST_GEOMETRY_BUFFER_USE_SIDE_OF_INTERSECTION)
         , count_within_near_offsetted(0)
+#endif
         , remove_on_multi(false)
         , count_on_occupied(0)
         , count_on_multi(0)

include/boost/geometry/algorithms/detail/buffer/buffered_piece_collection.hpp

@@ -19,8 +19,10 @@
 #include <boost/geometry/core/coordinate_type.hpp>
 #include <boost/geometry/core/point_type.hpp>
 
+#include <boost/geometry/algorithms/comparable_distance.hpp>
 #include <boost/geometry/algorithms/covered_by.hpp>
 #include <boost/geometry/algorithms/envelope.hpp>
+#include <boost/geometry/algorithms/is_convex.hpp>
 
 #include <boost/geometry/strategies/buffer.hpp>
 
@@ -126,6 +128,11 @@ struct buffered_piece_collection
     typedef geometry::model::ring<robust_point_type> robust_ring_type;
     typedef geometry::model::box<robust_point_type> robust_box_type;
 
+    typedef typename default_comparable_distance_result
+        <
+            robust_point_type
+        >::type robust_comparable_radius_type;
+
     typedef typename strategy::side::services::default_strategy
         <
             typename cs_tag<point_type>::type
@@ -191,13 +198,13 @@ struct buffered_piece_collection
         std::vector<point_type> helper_points; // 4 points for side, 3 points for join - 0 points for flat-end
 #endif
 
+        bool is_convex;
         bool is_monotonic_increasing[2]; // 0=x, 1=y
         bool is_monotonic_decreasing[2]; // 0=x, 1=y
 
         // Monotonic sections of pieces around points
         std::vector<section_type> sections;
 
-
         // Robust representations
         // 3: complete ring
         robust_ring_type robust_ring;
@@ -206,6 +213,27 @@ struct buffered_piece_collection
         robust_box_type robust_offsetted_envelope;
 
         std::vector<robust_turn> robust_turns; // Used only in insert_rescaled_piece_turns - we might use a map instead
+
+        robust_point_type robust_center;
+        robust_comparable_radius_type robust_min_comparable_radius;
+        robust_comparable_radius_type robust_max_comparable_radius;
+
+        piece()
+            : type(strategy::buffer::piece_type_unknown)
+            , index(-1)
+            , left_index(-1)
+            , right_index(-1)
+            , last_segment_index(-1)
+            , offsetted_count(-1)
+            , is_convex(false)
+            , robust_min_comparable_radius(0)
+            , robust_max_comparable_radius(0)
+        {
+            is_monotonic_increasing[0] = false;
+            is_monotonic_increasing[1] = false;
+            is_monotonic_decreasing[0] = false;
+            is_monotonic_decreasing[1] = false;
+        }
     };
 
     struct robust_original
@@ -444,6 +472,7 @@ struct buffered_piece_collection
             {
                 it->location = inside_buffer;
             }
+#if ! defined(BOOST_GEOMETRY_BUFFER_USE_SIDE_OF_INTERSECTION)
             if (it->count_within_near_offsetted > 0)
             {
                 // Within can have in rare cases a rounding issue. We don't discard this
@@ -451,6 +480,7 @@ struct buffered_piece_collection
                 // will never start a new ring from this type of points.
                 it->selectable_start = false;
             }
+#endif
         }
     }
 
@@ -544,6 +574,7 @@ struct buffered_piece_collection
             }
         }
 
+#if ! defined(BOOST_GEOMETRY_BUFFER_USE_SIDE_OF_INTERSECTION)
         // Insert all rescaled turn-points into these rings, to form a
         // reliable integer-based ring. All turns can be compared (inside) to this
         // rings to see if they are inside.
@@ -582,6 +613,7 @@ struct buffered_piece_collection
         }
 
         BOOST_ASSERT(assert_indices_in_robust_rings());
+#endif
     }
 
     template <std::size_t Dimension>
@@ -606,6 +638,8 @@ struct buffered_piece_collection
         pc.is_monotonic_decreasing[0] = true;
         pc.is_monotonic_decreasing[1] = true;
 
+        pc.is_convex = geometry::is_convex(pc.robust_ring);
+
         if (pc.offsetted_count < 2)
         {
             return;
@@ -621,7 +655,6 @@ struct buffered_piece_collection
             current = next;
             ++next;
         }
-
     }
 
     void determine_properties()
@@ -658,7 +691,31 @@ struct buffered_piece_collection
         geometry::sectionalize<false, dimensions>(pc.robust_ring,
                 detail::no_rescale_policy(), pc.sections);
 
-        // TODO (next phase) determine min/max radius
+        // Determine min/max radius
+        typedef geometry::model::referring_segment<robust_point_type const>
+            robust_segment_type;
+
+        typename robust_ring_type::const_iterator current = pc.robust_ring.begin();
+        typename robust_ring_type::const_iterator next = current + 1;
+
+        for (int i = 1; i < pc.offsetted_count; i++)
+        {
+            robust_segment_type s(*current, *next);
+            robust_comparable_radius_type const d
+                = geometry::comparable_distance(pc.robust_center, s);
+
+            if (i == 1 || d < pc.robust_min_comparable_radius)
+            {
+                pc.robust_min_comparable_radius = d;
+            }
+            if (i == 1 || d > pc.robust_max_comparable_radius)
+            {
+                pc.robust_max_comparable_radius = d;
+            }
+
+            current = next;
+            ++next;
+        }
     }
 
     inline void prepare_buffered_point_pieces()
@@ -772,6 +829,13 @@ struct buffered_piece_collection
         }
     }
 
+    inline void set_piece_center(point_type const& center)
+    {
+        BOOST_ASSERT(! m_pieces.empty());
+        geometry::recalculate(m_pieces.back().robust_center, center,
+                m_robust_policy);
+    }
+
     inline void finish_ring(bool is_interior = false, bool has_interiors = false)
     {
         if (m_first_piece_index == -1)

include/boost/geometry/algorithms/detail/buffer/get_piece_turns.hpp

@@ -17,6 +17,7 @@
 #include <boost/geometry/algorithms/detail/overlay/segment_identifier.hpp>
 #include <boost/geometry/algorithms/detail/overlay/get_turn_info.hpp>
 #include <boost/geometry/algorithms/detail/sections/section_functions.hpp>
+#include <boost/geometry/algorithms/detail/buffer/buffer_policies.hpp>
 
 
 namespace boost { namespace geometry
@@ -28,6 +29,31 @@ namespace detail { namespace buffer
 {
 
 
+#if defined(BOOST_GEOMETRY_BUFFER_USE_SIDE_OF_INTERSECTION)
+struct buffer_assign_turn
+{
+    static bool const include_no_turn = false;
+    static bool const include_degenerate = false;
+    static bool const include_opposite = false;
+
+    template
+    <
+        typename Info,
+        typename Point1,
+        typename Point2,
+        typename IntersectionInfo
+    >
+    static inline void apply(Info& info, Point1 const& p1, Point2 const& p2, IntersectionInfo const& iinfo)
+    {
+        info.rob_pi = iinfo.rpi();
+        info.rob_pj = iinfo.rpj();
+        info.rob_qi = iinfo.rqi();
+        info.rob_qj = iinfo.rqj();
+    }
+
+};
+#endif
+
 template
 <
     typename Pieces,
@@ -204,7 +230,11 @@ class piece_turn_visitor
                 // and iterating in sync with them...
                 typedef detail::overlay::get_turn_info
                     <
+#if defined(BOOST_GEOMETRY_BUFFER_USE_SIDE_OF_INTERSECTION)
+                        buffer_assign_turn
+#else
                         detail::overlay::assign_null_policy
+#endif
                     > turn_policy;
 
                 turn_policy::apply(*prev1, *it1, *next1,

include/boost/geometry/algorithms/detail/buffer/turn_in_piece_visitor.hpp

@@ -25,6 +25,11 @@
 #include <boost/geometry/algorithms/detail/overlay/get_turn_info.hpp>
 #include <boost/geometry/policies/compare.hpp>
 #include <boost/geometry/strategies/buffer.hpp>
+#include <boost/geometry/algorithms/detail/buffer/buffer_policies.hpp>
+
+#if defined(BOOST_GEOMETRY_BUFFER_USE_SIDE_OF_INTERSECTION)
+#include <boost/geometry/strategies/cartesian/side_of_intersection.hpp>
+#endif
 
 
 namespace boost { namespace geometry
@@ -89,8 +94,10 @@ enum analyse_result
     analyse_disjoint,
     analyse_within,
     analyse_on_original_boundary,
-    analyse_on_offsetted,
-    analyse_near_offsetted
+    analyse_on_offsetted
+#if ! defined(BOOST_GEOMETRY_BUFFER_USE_SIDE_OF_INTERSECTION)
+    , analyse_near_offsetted
+#endif
 };
 
 template <typename Point>
@@ -112,6 +119,31 @@ inline analyse_result check_segment(Point const& previous,
         Point const& current, Turn const& turn,
         bool from_monotonic)
 {
+
+#if defined(BOOST_GEOMETRY_BUFFER_USE_SIDE_OF_INTERSECTION)
+    typedef geometry::model::referring_segment<Point const> segment_type;
+    segment_type const p(turn.rob_pi, turn.rob_pj);
+    segment_type const q(turn.rob_qi, turn.rob_qj);
+    segment_type const r(previous, current);
+    int const side = strategy::side::side_of_intersection::apply(p, q, r,
+                turn.robust_point);
+
+    if (side == 0)
+    {
+        return analyse_on_offsetted;
+    }
+    if (side == -1 && from_monotonic)
+    {
+        return analyse_within;
+    }
+    if (side == 1 && from_monotonic)
+    {
+        return analyse_disjoint;
+    }
+    return analyse_continue;
+
+#else
+
     typedef typename strategy::side::services::default_strategy
         <
             typename cs_tag<Point>::type
@@ -156,6 +188,7 @@ inline analyse_result check_segment(Point const& previous,
 
     // Not monotonic, on left or right side: continue analysing
     return analyse_continue;
+#endif
 }
 
 
@@ -169,14 +202,22 @@ public :
         typedef typename Turn::robust_point_type point_type;
         typedef typename geometry::coordinate_type<point_type>::type coordinate_type;
 
-        coordinate_type const point_y = geometry::get<1>(turn.robust_point);
-
+#if defined(BOOST_GEOMETRY_BUFFER_USE_SIDE_OF_INTERSECTION)
+        typedef geometry::model::referring_segment<point_type const> segment_type;
+        segment_type const p(turn.rob_pi, turn.rob_pj);
+        segment_type const q(turn.rob_qi, turn.rob_qj);
+#else
         typedef strategy::within::winding<point_type> strategy_type;
 
         typename strategy_type::state_type state;
         strategy_type strategy;
         boost::ignore_unused(strategy);
-        
+#endif
+
+        BOOST_ASSERT(! piece.sections.empty());
+
+        coordinate_type const point_y = geometry::get<1>(turn.robust_point);
+
         for (std::size_t s = 0; s < piece.sections.size(); s++)
         {
             section_type const& section = piece.sections[s];
@@ -191,6 +232,38 @@ public :
                     point_type const& previous = piece.robust_ring[i - 1];
                     point_type const& current = piece.robust_ring[i];
 
+#if defined(BOOST_GEOMETRY_BUFFER_USE_SIDE_OF_INTERSECTION)
+
+                    // First check if it is in range - if it is not, the
+                    // expensive side_of_intersection does not need to be
+                    // applied
+                    coordinate_type y1 = geometry::get<1>(previous);
+                    coordinate_type y2 = geometry::get<1>(current);
+
+                    if (y1 > y2)
+                    {
+                        std::swap(y1, y2);
+                    }
+
+                    if (point_y >= y1 - 1 && point_y <= y2 + 1)
+                    {
+                        segment_type const r(previous, current);
+                        int const side = strategy::side::side_of_intersection::apply(p, q, r,
+                                    turn.robust_point);
+
+                        // Sections are monotonic in y-dimension
+                        if (side == 1)
+                        {
+                            // Left on segment
+                            return analyse_disjoint;
+                        }
+                        else if (side == 0)
+                        {
+                            // Collinear - TODO: check if really on segment
+                            return analyse_on_offsetted;
+                        }
+                    }
+#else
                     analyse_result code = check_segment(previous, current, turn, false);
                     if (code != analyse_continue)
                     {
@@ -200,10 +273,15 @@ public :
                     // Get the state (to determine it is within), we don't have
                     // to cover the on-segment case (covered above)
                     strategy.apply(turn.robust_point, previous, current, state);
+#endif
                 }
             }
         }
 
+#if defined(BOOST_GEOMETRY_BUFFER_USE_SIDE_OF_INTERSECTION)
+        // It is nowhere outside, and not on segment, so it is within
+        return analyse_within;
+#else
         int const code = strategy.result(state);
         if (code == 1)
         {
@@ -216,6 +294,7 @@ public :
 
         // Should normally not occur - on-segment is covered
         return analyse_unknown;
+#endif
     }
 
 };
@@ -228,6 +307,48 @@ class analyse_turn_wrt_piece
                 bool is_original,
                 Point const& offsetted)
     {
+#if defined(BOOST_GEOMETRY_BUFFER_USE_SIDE_OF_INTERSECTION)
+        typedef geometry::model::referring_segment<Point const> segment_type;
+        segment_type const p(turn.rob_pi, turn.rob_pj);
+        segment_type const q(turn.rob_qi, turn.rob_qj);
+        segment_type const r(s1, s2);
+        int const side = strategy::side::side_of_intersection::apply(p, q, r,
+                    turn.robust_point);
+
+        if (side == 1)
+        {
+            // left of segment
+            return analyse_disjoint;
+        }
+        else if (side == 0)
+        {
+            // If is collinear, either on segment or before/after
+            typedef geometry::model::box<Point> box_type;
+
+            box_type box;
+            geometry::assign_inverse(box);
+            geometry::expand(box, s1);
+            geometry::expand(box, s2);
+
+            if (geometry::covered_by(turn.robust_point, box))
+            {
+                // Points on helper-segments are considered as within
+                // Points on original boundary are processed differently
+                return is_original
+                    ? analyse_on_original_boundary
+                    : analyse_within;
+            }
+
+            // It is collinear but not on the segment. Because these
+            // segments are convex, it is outside
+            // Unless the offsetted ring is collinear or concave w.r.t.
+            // helper-segment but that scenario is not yet supported
+            return analyse_disjoint;
+        }
+
+        // right of segment
+        return analyse_continue;
+#else
         typedef typename strategy::side::services::default_strategy
             <
                 typename cs_tag<Point>::type
@@ -276,6 +397,7 @@ class analyse_turn_wrt_piece
 
         // right of segment
         return analyse_continue;
+#endif
     }
 
     template <typename Turn, typename Piece>
@@ -492,6 +614,60 @@ class turn_in_piece_visitor
         return false;
     }
 
+#if defined(BOOST_GEOMETRY_BUFFER_USE_SIDE_OF_INTERSECTION)
+    template <typename Turn, typename Piece>
+    static inline int turn_in_convex_piece(Turn const& turn,
+            Piece const& piece)
+    {
+        typedef typename Turn::robust_point_type point_type;
+        typedef typename Piece::piece_robust_ring_type ring_type;
+        typedef geometry::model::referring_segment<point_type const> segment;
+
+        segment const p(turn.rob_pi, turn.rob_pj);
+        segment const q(turn.rob_qi, turn.rob_qj);
+
+        typedef typename boost::range_iterator<ring_type const>::type iterator_type;
+        iterator_type it = boost::begin(piece.robust_ring);
+        iterator_type end = boost::end(piece.robust_ring);
+
+        // A robust ring is always closed, and always clockwise
+        for (iterator_type previous = it++; it != end; ++previous, ++it)
+        {
+            geometry::equal_to<point_type> comparator;
+            if (comparator(*previous, *it))
+            {
+                // Points are the same
+                continue;
+            }
+
+            segment r(*previous, *it);
+
+            int const side = strategy::side::side_of_intersection::apply(p, q, r,
+                        turn.robust_point);
+
+            if (side == 1)
+            {
+                // IP is left of segment, so it is outside
+                return -1; // outside
+            }
+            else if (side == 0)
+            {
+                // IP is collinear with segment. TODO: we should analyze this further
+                // For now we use the fallback point
+                if (in_box(*previous, *it, turn.robust_point))
+                {
+                    return 0;
+                }
+                else
+                {
+                    return -1; // outside
+                }
+            }
+        }
+        return 1; // inside
+    }
+#endif
+
 
 public:
 
@@ -530,12 +706,37 @@ public:
         }
 
         // TODO: mutable_piece to make some on-demand preparations in analyse
+        Turn& mutable_turn = m_turns[turn.turn_index];
+
+        if (piece.type == geometry::strategy::buffer::buffered_point)
+        {
+            // Optimization for buffer around points: if distance from center
+            // is not between min/max radius, the result is clear
+            typedef typename default_comparable_distance_result
+                <
+                    typename Turn::robust_point_type
+                >::type distance_type;
+
+            distance_type const cd
+                = geometry::comparable_distance(piece.robust_center,
+                        turn.robust_point);
+
+            if (cd < piece.robust_min_comparable_radius)
+            {
+                mutable_turn.count_within++;
+                return;
+            }
+            if (cd > piece.robust_max_comparable_radius)
+            {
+                return;
+            }
+        }
+
         analyse_result analyse_code =
             piece.type == geometry::strategy::buffer::buffered_point
                 ? analyse_turn_wrt_point_piece::apply(turn, piece)
                 : analyse_turn_wrt_piece::apply(turn, piece);
 
-        Turn& mutable_turn = m_turns[turn.turn_index];
         switch(analyse_code)
         {
             case analyse_disjoint :
@@ -549,16 +750,32 @@ public:
             case analyse_within :
                 mutable_turn.count_within++;
                 return;
+#if ! defined(BOOST_GEOMETRY_BUFFER_USE_SIDE_OF_INTERSECTION)
             case analyse_near_offsetted :
                 mutable_turn.count_within_near_offsetted++;
                 return;
+#endif
             default :
                 break;
         }
 
-        // TODO: this point_in_geometry is a performance-bottleneck here and
-        // will be replaced completely by extending analyse_piece functionality
+#if defined(BOOST_GEOMETRY_BUFFER_USE_SIDE_OF_INTERSECTION)
+        // We don't know (yet)
+        int geometry_code = 0;
+        if (piece.is_convex)
+        {
+            geometry_code = turn_in_convex_piece(turn, piece);
+        }
+        else
+        {
+
+            // TODO: this point_in_geometry is a performance-bottleneck here and
+            // will be replaced completely by extending analyse_piece functionality
+            geometry_code = detail::within::point_in_geometry(turn.robust_point, piece.robust_ring);
+        }
+#else
         int geometry_code = detail::within::point_in_geometry(turn.robust_point, piece.robust_ring);
+#endif
 
         if (geometry_code == 1)
         {

include/boost/geometry/algorithms/detail/disjoint/point_point.hpp

@@ -1,12 +1,12 @@
 // Boost.Geometry (aka GGL, Generic Geometry Library)
 
-// Copyright (c) 2007-2014 Barend Gehrels, Amsterdam, the Netherlands.
-// Copyright (c) 2008-2014 Bruno Lalande, Paris, France.
-// Copyright (c) 2009-2014 Mateusz Loskot, London, UK.
-// Copyright (c) 2013-2014 Adam Wulkiewicz, Lodz, Poland
+// Copyright (c) 2007-2015 Barend Gehrels, Amsterdam, the Netherlands.
+// Copyright (c) 2008-2015 Bruno Lalande, Paris, France.
+// Copyright (c) 2009-2015 Mateusz Loskot, London, UK.
+// Copyright (c) 2013-2015 Adam Wulkiewicz, Lodz, Poland
 
-// This file was modified by Oracle on 2013-2014.
-// Modifications copyright (c) 2013-2014, Oracle and/or its affiliates.
+// This file was modified by Oracle on 2013, 2014, 2015.
+// Modifications copyright (c) 2013-2015, Oracle and/or its affiliates.
 
 // Contributed and/or modified by Adam Wulkiewicz, on behalf of Oracle
 // Contributed and/or modified by Menelaos Karavelas, on behalf of Oracle
@@ -23,11 +23,19 @@
 
 #include <cstddef>
 
+#include <boost/type_traits/is_same.hpp>
+
 #include <boost/geometry/core/access.hpp>
-#include <boost/geometry/core/coordinate_dimension.hpp>
+#include <boost/geometry/core/radian_access.hpp>
+#include <boost/geometry/core/coordinate_system.hpp>
+#include <boost/geometry/core/coordinate_type.hpp>
+#include <boost/geometry/core/cs.hpp>
 #include <boost/geometry/core/tags.hpp>
 
+#include <boost/geometry/util/condition.hpp>
 #include <boost/geometry/util/math.hpp>
+#include <boost/geometry/util/normalize_spheroidal_coordinates.hpp>
+#include <boost/geometry/util/select_most_precise.hpp>
 
 #include <boost/geometry/algorithms/dispatch/disjoint.hpp>
 
@@ -40,12 +48,103 @@ namespace boost { namespace geometry
 namespace detail { namespace disjoint
 {
 
+
+class point_point_on_spheroid
+{
+private:
+    template
+    <
+        typename Units,
+        typename CoordinateType1,
+        typename CoordinateType2
+    >
+    static inline bool apply_same_units(CoordinateType1& lon1,
+                                        CoordinateType1& lat1,
+                                        CoordinateType2& lon2,
+                                        CoordinateType2& lat2)
+    {
+        math::normalize_spheroidal_coordinates<Units>(lon1, lat1);
+        math::normalize_spheroidal_coordinates<Units>(lon2, lat2);
+
+        return ! math::equals(lat1, lat2) || ! math::equals(lon1, lon2);
+    }
+
+public:
+    template <typename Point1, typename Point2>
+    static inline bool apply(Point1 const& p1, Point2 const& p2)
+    {
+        typedef typename coordinate_type<Point1>::type coordinate_type1;
+        typedef typename coordinate_type<Point2>::type coordinate_type2;
+
+        typedef typename coordinate_system<Point1>::type::units units1;
+        typedef typename coordinate_system<Point2>::type::units units2;
+
+        if (BOOST_GEOMETRY_CONDITION((boost::is_same<units1, units2>::value)))
+        {
+            coordinate_type1 lon1 = geometry::get<0>(p1);
+            coordinate_type1 lat1 = geometry::get<1>(p1);
+            coordinate_type2 lon2 = geometry::get<0>(p2);
+            coordinate_type2 lat2 = geometry::get<1>(p2);
+
+            return apply_same_units<units1>(lon1, lat1, lon2, lat2);
+        }
+
+        typedef typename geometry::select_most_precise
+            <
+                typename fp_coordinate_type<Point1>::type,
+                typename fp_coordinate_type<Point2>::type
+            >::type calculation_type;
+
+        calculation_type lon1 = get_as_radian<0>(p1);
+        calculation_type lat1 = get_as_radian<1>(p1);
+
+        calculation_type lon2 = get_as_radian<0>(p2);
+        calculation_type lat2 = get_as_radian<1>(p2);
+
+        return apply_same_units<radian>(lon1, lat1, lon2, lat2);
+    }
+};
+
+
+template
+<
+    typename Point1, typename Point2,
+    std::size_t Dimension, std::size_t DimensionCount,
+    typename CSTag1 = typename cs_tag<Point1>::type,
+    typename CSTag2 = CSTag1
+>
+struct point_point
+    : point_point<Point1, Point2, Dimension, DimensionCount, cartesian_tag>
+{};
+
 template
 <
     typename Point1, typename Point2,
     std::size_t Dimension, std::size_t DimensionCount
 >
 struct point_point
+    <
+        Point1, Point2, Dimension, DimensionCount, spherical_equatorial_tag
+    > : point_point_on_spheroid
+{};
+
+template
+<
+    typename Point1, typename Point2,
+    std::size_t Dimension, std::size_t DimensionCount
+>
+struct point_point
+    <
+        Point1, Point2, Dimension, DimensionCount, geographic_tag
+    > : point_point_on_spheroid
+{};
+
+template
+<
+    typename Point1, typename Point2,
+    std::size_t Dimension, std::size_t DimensionCount
+>
+struct point_point<Point1, Point2, Dimension, DimensionCount, cartesian_tag>
 {
     static inline bool apply(Point1 const& p1, Point2 const& p2)
     {
@@ -61,9 +160,11 @@ struct point_point
     }
 };
 
-
 template <typename Point1, typename Point2, std::size_t DimensionCount>
-struct point_point<Point1, Point2, DimensionCount, DimensionCount>
+struct point_point
+    <
+        Point1, Point2, DimensionCount, DimensionCount, cartesian_tag
+    >
 {
     static inline bool apply(Point1 const& , Point2 const& )
     {

include/boost/geometry/algorithms/detail/overlay/clip_linestring.hpp

@@ -1,6 +1,11 @@
 // Boost.Geometry (aka GGL, Generic Geometry Library)
 
-// Copyright (c) 2007-2012 Barend Gehrels, Amsterdam, the Netherlands.
+// Copyright (c) 2007-2015 Barend Gehrels, Amsterdam, the Netherlands.
+
+// This file was modified by Oracle on 2015.
+// Modifications copyright (c) 2015 Oracle and/or its affiliates.
+
+// Contributed and/or modified by Menelaos Karavelas, on behalf of Oracle
 
 // Use, modification and distribution is subject to the Boost Software License,
 // Version 1.0. (See accompanying file LICENSE_1_0.txt or copy at
@@ -160,10 +165,12 @@ template
     typename OutputLinestring,
     typename OutputIterator,
     typename Range,
+    typename RobustPolicy,
     typename Box,
     typename Strategy
 >
 OutputIterator clip_range_with_box(Box const& b, Range const& range,
+            RobustPolicy const&,
             OutputIterator out, Strategy const& strategy)
 {
     if (boost::begin(range) == boost::end(range))

include/boost/geometry/algorithms/detail/overlay/intersection_insert.hpp

@@ -410,13 +410,13 @@ struct intersection_insert
     template <typename RobustPolicy, typename OutputIterator, typename Strategy>
     static inline OutputIterator apply(Linestring const& linestring,
             Box const& box,
-            RobustPolicy const& ,
+            RobustPolicy const& robust_policy,
             OutputIterator out, Strategy const& )
     {
         typedef typename point_type<GeometryOut>::type point_type;
         strategy::intersection::liang_barsky<Box, point_type> lb_strategy;
         return detail::intersection::clip_range_with_box
-            <GeometryOut>(box, linestring, out, lb_strategy);
+            <GeometryOut>(box, linestring, robust_policy, out, lb_strategy);
     }
 };
 
@@ -488,7 +488,7 @@ struct intersection_insert
     template <typename RobustPolicy, typename OutputIterator, typename Strategy>
     static inline OutputIterator apply(Segment const& segment,
             Box const& box,
-            RobustPolicy const& ,// TODO: propagate to clip_range_with_box
+            RobustPolicy const& robust_policy,
             OutputIterator out, Strategy const& )
     {
         geometry::segment_view<Segment> range(segment);
@@ -496,7 +496,7 @@ struct intersection_insert
         typedef typename point_type<GeometryOut>::type point_type;
         strategy::intersection::liang_barsky<Box, point_type> lb_strategy;
         return detail::intersection::clip_range_with_box
-            <GeometryOut>(box, range, out, lb_strategy);
+            <GeometryOut>(box, range, robust_policy, out, lb_strategy);
     }
 };
 

include/boost/geometry/algorithms/detail/relate/less.hpp

@@ -14,6 +14,10 @@
 #ifndef BOOST_GEOMETRY_ALGORITHMS_DETAIL_RELATE_LESS_HPP
 #define BOOST_GEOMETRY_ALGORITHMS_DETAIL_RELATE_LESS_HPP
 
+#include <boost/geometry/core/coordinate_dimension.hpp>
+#include <boost/geometry/core/coordinate_type.hpp>
+#include <boost/geometry/util/math.hpp>
+
 namespace boost { namespace geometry
 {
 

include/boost/geometry/algorithms/is_convex.hpp

@@ -0,0 +1,160 @@
+// Boost.Geometry (aka GGL, Generic Geometry Library)
+
+// Copyright (c) 2015 Barend Gehrels, Amsterdam, the Netherlands.
+
+// Use, modification and distribution is subject to the Boost Software License,
+// Version 1.0. (See accompanying file LICENSE_1_0.txt or copy at
+// http://www.boost.org/LICENSE_1_0.txt)
+
+#ifndef BOOST_GEOMETRY_ALGORITHMS_IS_CONVEX_HPP
+#define BOOST_GEOMETRY_ALGORITHMS_IS_CONVEX_HPP
+
+
+#include <boost/geometry/core/access.hpp>
+#include <boost/geometry/core/closure.hpp>
+#include <boost/geometry/core/cs.hpp>
+#include <boost/geometry/core/coordinate_dimension.hpp>
+#include <boost/geometry/core/point_type.hpp>
+#include <boost/geometry/algorithms/detail/equals/point_point.hpp>
+#include <boost/geometry/iterators/ever_circling_iterator.hpp>
+#include <boost/geometry/strategies/side.hpp>
+#include <boost/geometry/strategies/cartesian/side_by_triangle.hpp>
+#include <boost/geometry/views/detail/normalized_view.hpp>
+
+namespace boost { namespace geometry
+{
+
+
+#ifndef DOXYGEN_NO_DETAIL
+namespace detail { namespace is_convex
+{
+
+struct ring_is_convex
+{
+    template <typename Ring>
+    static inline bool apply(Ring const& ring)
+    {
+        typedef typename geometry::point_type<Ring>::type point_type;
+        typedef typename strategy::side::services::default_strategy
+        <
+            typename cs_tag<point_type>::type
+        >::type side_strategy_type;
+
+        std::size_t n = boost::size(ring);
+        if (boost::size(ring) < core_detail::closure::minimum_ring_size
+                                    <
+                                        geometry::closure<Ring>::value
+                                    >::value)
+        {
+            // (Too) small rings are considered as non-concave, is convex
+            return true;
+        }
+
+        // Walk in clockwise direction, consider ring as closed
+        // (though closure is not important in this algorithm - any dupped
+        //  point is skipped)
+        typedef detail::normalized_view<Ring const> view_type;
+        view_type view(ring);
+
+        typedef geometry::ever_circling_range_iterator<view_type const> it_type;
+        it_type previous(view);
+        it_type current(view);
+        current++;
+
+        std::size_t index = 1;
+        while (equals::equals_point_point(*current, *previous) && index < n)
+        {
+            current++;
+            index++;
+        }
+
+        if (index == n)
+        {
+            // All points are apparently equal
+            return true;
+        }
+
+        it_type next = current;
+        next++;
+        while (equals::equals_point_point(*current, *next))
+        {
+            next++;
+        }
+
+        // We have now three different points on the ring
+        // Walk through all points, use a counter because of the ever-circling
+        // iterator
+        for (std::size_t i = 0; i < n; i++)
+        {
+            int const side = side_strategy_type::apply(*previous, *current, *next);
+            if (side == 1)
+            {
+                // Next is on the left side of clockwise ring:
+                // the piece is not convex
+                return false;
+            }
+
+            previous = current;
+            current = next;
+
+            // Advance next to next different point
+            // (because there are non-equal points, this loop is not infinite)
+            next++;
+            while (equals::equals_point_point(*current, *next))
+            {
+                next++;
+            }
+        }
+        return true;
+    }
+};
+
+
+}} // namespace detail::is_convex
+#endif // DOXYGEN_NO_DETAIL
+
+
+#ifndef DOXYGEN_NO_DISPATCH
+namespace dispatch
+{
+
+template
+<
+    typename Geometry,
+    typename Tag = typename tag<Geometry>::type
+>
+struct is_convex : not_implemented<Tag>
+{};
+
+template <typename Box>
+struct is_convex<Box, box_tag>
+{
+    static inline bool apply(Box const& )
+    {
+        // Any box is convex (TODO: consider spherical boxes)
+        return true;
+    }
+};
+
+template <typename Box>
+struct is_convex<Box, ring_tag> : detail::is_convex::ring_is_convex
+{};
+
+
+} // namespace dispatch
+#endif // DOXYGEN_NO_DISPATCH
+
+// TODO: variants
+
+// TODO: documentation / qbk
+template<typename Geometry>
+inline bool is_convex(Geometry const& geometry)
+{
+    return dispatch::is_convex<Geometry>::apply(geometry);
+}
+
+
+}} // namespace boost::geometry
+
+
+#endif // BOOST_GEOMETRY_ALGORITHMS_IS_CONVEX_HPP

include/boost/geometry/core/radian_access.hpp

@@ -1,8 +1,13 @@
 // Boost.Geometry (aka GGL, Generic Geometry Library)
 
-// Copyright (c) 2007-2012 Barend Gehrels, Amsterdam, the Netherlands.
-// Copyright (c) 2008-2012 Bruno Lalande, Paris, France.
-// Copyright (c) 2009-2012 Mateusz Loskot, London, UK.
+// Copyright (c) 2007-2015 Barend Gehrels, Amsterdam, the Netherlands.
+// Copyright (c) 2008-2015 Bruno Lalande, Paris, France.
+// Copyright (c) 2009-2015 Mateusz Loskot, London, UK.
+
+// This file was modified by Oracle on 2015.
+// Modifications copyright (c) 2015, Oracle and/or its affiliates.
+
+// Contributed and/or modified by Menelaos Karavelas, on behalf of Oracle
 
 // Parts of Boost.Geometry are redesigned from Geodan's Geographic Library
 // (geolib/GGL), copyright (c) 1995-2010 Geodan, Amsterdam, the Netherlands.
@@ -47,7 +52,8 @@ struct degree_radian_converter
         return boost::numeric_cast
             <
                 coordinate_type
-            >(geometry::get<Dimension>(geometry) * geometry::math::d2r);
+            >(geometry::get<Dimension>(geometry)
+              * math::d2r<coordinate_type>());
     }
 
     static inline void set(Geometry& geometry, coordinate_type const& radians)
@@ -55,7 +61,7 @@ struct degree_radian_converter
         geometry::set<Dimension>(geometry, boost::numeric_cast
             <
                 coordinate_type
-            >(radians * geometry::math::r2d));
+            >(radians * math::r2d<coordinate_type>()));
     }
 
 };

include/boost/geometry/extensions/gis/projections/proj/aeqd.hpp

@@ -419,12 +419,14 @@ namespace boost { namespace geometry { namespace projections
             public :
                 virtual projection<Geographic, Cartesian>* create_new(const Parameters& par) const
                 {
-                    if (! par.es)
-                        return new base_v_fi<aeqd_spheroid<Geographic, Cartesian, Parameters>, Geographic, Cartesian, Parameters>(par);
-                    else if (pj_param(par.params, "bguam").i)
+                    bool const guam = pj_param(par.params, "bguam").i;
+
+                    if (par.es && ! guam)
+                        return new base_v_fi<aeqd_ellipsoid<Geographic, Cartesian, Parameters>, Geographic, Cartesian, Parameters>(par);
+                    else if (par.es && guam)
                         return new base_v_fi<aeqd_guam<Geographic, Cartesian, Parameters>, Geographic, Cartesian, Parameters>(par);
                     else
-                        return new base_v_fi<aeqd_ellipsoid<Geographic, Cartesian, Parameters>, Geographic, Cartesian, Parameters>(par);
+                        return new base_v_fi<aeqd_spheroid<Geographic, Cartesian, Parameters>, Geographic, Cartesian, Parameters>(par);
                 }
         };
 

include/boost/geometry/extensions/gis/projections/proj/gn_sinu.hpp

@@ -213,9 +213,9 @@ namespace boost { namespace geometry { namespace projections
         \image html ex_gn_sinu.gif
     */
     template <typename Geographic, typename Cartesian, typename Parameters = parameters>
-    struct gn_sinu_ellipsoid : public detail::gn_sinu::base_gn_sinu_ellipsoid<Geographic, Cartesian, Parameters>
+    struct gn_sinu_spheroid : public detail::gn_sinu::base_gn_sinu_spheroid<Geographic, Cartesian, Parameters>
     {
-        inline gn_sinu_ellipsoid(const Parameters& par) : detail::gn_sinu::base_gn_sinu_ellipsoid<Geographic, Cartesian, Parameters>(par)
+        inline gn_sinu_spheroid(const Parameters& par) : detail::gn_sinu::base_gn_sinu_spheroid<Geographic, Cartesian, Parameters>(par)
         {
             detail::gn_sinu::setup_gn_sinu(this->m_par, this->m_proj_parm);
         }
@@ -243,70 +243,6 @@ namespace boost { namespace geometry { namespace projections
         }
     };
 
-    /*!
-        \brief Eckert VI projection
-        \ingroup projections
-        \tparam Geographic latlong point type
-        \tparam Cartesian xy point type
-        \tparam Parameters parameter type
-        \par Projection characteristics
-         - Pseudocylindrical
-         - Spheroid
-        \par Example
-        \image html ex_eck6.gif
-    */
-    template <typename Geographic, typename Cartesian, typename Parameters = parameters>
-    struct eck6_ellipsoid : public detail::gn_sinu::base_gn_sinu_ellipsoid<Geographic, Cartesian, Parameters>
-    {
-        inline eck6_ellipsoid(const Parameters& par) : detail::gn_sinu::base_gn_sinu_ellipsoid<Geographic, Cartesian, Parameters>(par)
-        {
-            detail::gn_sinu::setup_eck6(this->m_par, this->m_proj_parm);
-        }
-    };
-
-    /*!
-        \brief McBryde-Thomas Flat-Polar Sinusoidal projection
-        \ingroup projections
-        \tparam Geographic latlong point type
-        \tparam Cartesian xy point type
-        \tparam Parameters parameter type
-        \par Projection characteristics
-         - Pseudocylindrical
-         - Spheroid
-        \par Example
-        \image html ex_mbtfps.gif
-    */
-    template <typename Geographic, typename Cartesian, typename Parameters = parameters>
-    struct mbtfps_ellipsoid : public detail::gn_sinu::base_gn_sinu_ellipsoid<Geographic, Cartesian, Parameters>
-    {
-        inline mbtfps_ellipsoid(const Parameters& par) : detail::gn_sinu::base_gn_sinu_ellipsoid<Geographic, Cartesian, Parameters>(par)
-        {
-            detail::gn_sinu::setup_mbtfps(this->m_par, this->m_proj_parm);
-        }
-    };
-
-    /*!
-        \brief General Sinusoidal Series projection
-        \ingroup projections
-        \tparam Geographic latlong point type
-        \tparam Cartesian xy point type
-        \tparam Parameters parameter type
-        \par Projection characteristics
-         - Pseudocylindrical
-         - Spheroid
-         - m= n=
-        \par Example
-        \image html ex_gn_sinu.gif
-    */
-    template <typename Geographic, typename Cartesian, typename Parameters = parameters>
-    struct gn_sinu_spheroid : public detail::gn_sinu::base_gn_sinu_spheroid<Geographic, Cartesian, Parameters>
-    {
-        inline gn_sinu_spheroid(const Parameters& par) : detail::gn_sinu::base_gn_sinu_spheroid<Geographic, Cartesian, Parameters>(par)
-        {
-            detail::gn_sinu::setup_gn_sinu(this->m_par, this->m_proj_parm);
-        }
-    };
-
     /*!
         \brief Sinusoidal (Sanson-Flamsteed) projection
         \ingroup projections

include/boost/geometry/extensions/gis/projections/proj/isea.hpp

@@ -56,11 +56,11 @@ namespace boost { namespace geometry { namespace projections
             static const double DEG144 = 2.51327412287183459075;
             static const double DEG36 = 0.62831853071795864768;
             static const double DEG108 = 1.88495559215387594306;
-            static const double DEG180 = M_PI;
+            static const double DEG180 = boost::math::constants::pi<double>();
             static const double ISEA_SCALE = 0.8301572857837594396028083;
             static const double V_LAT = 0.46364760899944494524;
-            static const double RAD2DEG = (180.0/M_PI);
-            static const double DEG2RAD = (M_PI/180.0);
+            static const double RAD2DEG = (180.0/boost::math::constants::pi<double>());
+            static const double DEG2RAD = (boost::math::constants::pi<double>()/180.0);
             static const double E_RAD = 0.91843818702186776133;
             static const double F_RAD = 0.18871053072122403508;
             static const double TABLE_G = 0.6615845383;
@@ -116,7 +116,7 @@ namespace boost { namespace geometry { namespace projections
                 int ix, iy, iz, s;
                 struct hex h;
 
-                x = x / cos(30 * M_PI / 180.0); /* rotated X coord */
+                x = x / cos(30 * boost::math::constants::pi<double>() / 180.0); /* rotated X coord */
                 y = y - x / 2.0; /* adjustment for rotated X */
 
                 /* adjust for actual hexwidth */
@@ -215,7 +215,7 @@ namespace boost { namespace geometry { namespace projections
             };
 
 
-            /* sqrt(5)/M_PI */
+            /* sqrt(5)/boost::math::constants::pi<double>() */
 
             /* 26.565051177 degrees */
 
@@ -409,7 +409,7 @@ namespace boost { namespace geometry { namespace projections
                     /* TODO I don't know why we do this.  It's not in snyder */
                     /* maybe because we should have picked a better vertex */
                     if (Az < 0.0) {
-                        Az += 2.0 * M_PI;
+                        Az += 2.0 * boost::math::constants::pi<double>();
                     }
 
                     Az_adjust_multiples = 0;
@@ -447,7 +447,7 @@ namespace boost { namespace geometry { namespace projections
                     H = acos(sin(Az) * sin(G) * cos(g) - cos(Az) * cos(G));
 
                     /* eq 7 */
-                    /* Ag = (Az + G + H - DEG180) * M_PI * R * R / DEG180; */
+                    /* Ag = (Az + G + H - DEG180) * boost::math::constants::pi<double>() * R * R / DEG180; */
                     Ag = Az + G + H - DEG180;
 
                     /* eq 8 */
@@ -523,11 +523,11 @@ namespace boost { namespace geometry { namespace projections
 
                 /* normalize longitude */
                 /* TODO can we just do a modulus ? */
-                lambdap = fmod(lambdap, 2 * M_PI);
-                while (lambdap > M_PI)
-                    lambdap -= 2 * M_PI;
-                while (lambdap < -M_PI)
-                    lambdap += 2 * M_PI;
+                lambdap = fmod(lambdap, 2 * boost::math::constants::pi<double>());
+                while (lambdap > boost::math::constants::pi<double>())
+                    lambdap -= 2 * boost::math::constants::pi<double>();
+                while (lambdap < -boost::math::constants::pi<double>())
+                    lambdap += 2 * boost::math::constants::pi<double>();
 
                 phip = asin(sin_phip);
 
@@ -543,19 +543,19 @@ namespace boost { namespace geometry { namespace projections
             {
                 struct isea_geo npt;
 
-                np->lon += M_PI;
+                np->lon += boost::math::constants::pi<double>();
                 npt = snyder_ctran(np, pt);
-                np->lon -= M_PI;
+                np->lon -= boost::math::constants::pi<double>();
 
-                npt.lon -= (M_PI - lon0 + np->lon);
+                npt.lon -= (boost::math::constants::pi<double>() - lon0 + np->lon);
 
-                npt.lon += M_PI;
+                npt.lon += boost::math::constants::pi<double>();
                 /* normalize longitude */
-                npt.lon = fmod(npt.lon, 2 * M_PI);
-                while (npt.lon > M_PI)
-                    npt.lon -= 2 * M_PI;
-                while (npt.lon < -M_PI)
-                    npt.lon += 2 * M_PI;
+                npt.lon = fmod(npt.lon, 2 * boost::math::constants::pi<double>());
+                while (npt.lon > boost::math::constants::pi<double>())
+                    npt.lon -= 2 * boost::math::constants::pi<double>();
+                while (npt.lon < -boost::math::constants::pi<double>())
+                    npt.lon += 2 * boost::math::constants::pi<double>();
 
                 return npt;
             }
@@ -601,7 +601,7 @@ namespace boost { namespace geometry { namespace projections
             {
                 if (!g)
                     return 0;
-                g->o_lat = M_PI / 2.0;
+                g->o_lat = boost::math::constants::pi<double>() / 2.0;
                 g->o_lon = 0.0;
                 g->o_az = 0;
                 return 1;
@@ -637,9 +637,9 @@ namespace boost { namespace geometry { namespace projections
 
                 double          x, y;
 
-                rad = -degrees * M_PI / 180.0;
-                while (rad >= 2.0 * M_PI) rad -= 2.0 * M_PI;
-                while (rad <= -2.0 * M_PI) rad += 2.0 * M_PI;
+                rad = -degrees * boost::math::constants::pi<double>() / 180.0;
+                while (rad >= 2.0 * boost::math::constants::pi<double>()) rad -= 2.0 * boost::math::constants::pi<double>();
+                while (rad <= -2.0 * boost::math::constants::pi<double>()) rad += 2.0 * boost::math::constants::pi<double>();
 
                 x = pt->x * cos(rad) + pt->y * sin(rad);
                 y = -pt->x * sin(rad) + pt->y * cos(rad);
@@ -676,7 +676,7 @@ namespace boost { namespace geometry { namespace projections
                 isea_rotate(pt, downtri ? 240.0 : 60.0);
                 if (downtri) {
                     pt->x += 0.5;
-                    /* pt->y += cos(30.0 * M_PI / 180.0); */
+                    /* pt->y += cos(30.0 * boost::math::constants::pi<double>() / 180.0); */
                     pt->y += .86602540378443864672;
                 }
                 return quad;
@@ -697,7 +697,7 @@ namespace boost { namespace geometry { namespace projections
                 sidelength = (pow(2.0, g->resolution) + 1.0) / 2.0;
 
                 /* apex to base is cos(30deg) */
-                hexwidth = cos(M_PI / 6.0) / sidelength;
+                hexwidth = cos(boost::math::constants::pi<double>() / 6.0) / sidelength;
 
                 /* TODO I think sidelength is always x.5, so
                  * (int)sidelength * 2 + 1 might be just as good

include/boost/geometry/extensions/gis/projections/proj/ob_tran.hpp

@@ -292,6 +292,7 @@ namespace boost { namespace geometry { namespace projections
                     detail::ob_tran::par_ob_tran<Geographic, Cartesian> proj_parm;
                     Parameters p = par;
                     double phip = setup_ob_tran(p, proj_parm, false);
+
                     if (fabs(phip) > detail::ob_tran::TOL)
                         return new base_v_fi<ob_tran_oblique<Geographic, Cartesian, Parameters>, Geographic, Cartesian, Parameters>(par);
                     else

include/boost/geometry/extensions/gis/projections/proj/stere.hpp

@@ -360,29 +360,6 @@ namespace boost { namespace geometry { namespace projections
         }
     };
 
-    /*!
-        \brief Universal Polar Stereographic projection
-        \ingroup projections
-        \tparam Geographic latlong point type
-        \tparam Cartesian xy point type
-        \tparam Parameters parameter type
-        \par Projection characteristics
-         - Azimuthal
-         - Spheroid
-         - Ellipsoid
-         - south
-        \par Example
-        \image html ex_ups.gif
-    */
-    template <typename Geographic, typename Cartesian, typename Parameters = parameters>
-    struct ups_ellipsoid : public detail::stere::base_stere_ellipsoid<Geographic, Cartesian, Parameters>
-    {
-        inline ups_ellipsoid(const Parameters& par) : detail::stere::base_stere_ellipsoid<Geographic, Cartesian, Parameters>(par)
-        {
-            detail::stere::setup_ups(this->m_par, this->m_proj_parm);
-        }
-    };
-
     /*!
         \brief Stereographic projection
         \ingroup projections
@@ -406,6 +383,29 @@ namespace boost { namespace geometry { namespace projections
         }
     };
 
+    /*!
+        \brief Universal Polar Stereographic projection
+        \ingroup projections
+        \tparam Geographic latlong point type
+        \tparam Cartesian xy point type
+        \tparam Parameters parameter type
+        \par Projection characteristics
+         - Azimuthal
+         - Spheroid
+         - Ellipsoid
+         - south
+        \par Example
+        \image html ex_ups.gif
+    */
+    template <typename Geographic, typename Cartesian, typename Parameters = parameters>
+    struct ups_ellipsoid : public detail::stere::base_stere_ellipsoid<Geographic, Cartesian, Parameters>
+    {
+        inline ups_ellipsoid(const Parameters& par) : detail::stere::base_stere_ellipsoid<Geographic, Cartesian, Parameters>(par)
+        {
+            detail::stere::setup_ups(this->m_par, this->m_proj_parm);
+        }
+    };
+
     /*!
         \brief Universal Polar Stereographic projection
         \ingroup projections
@@ -453,7 +453,10 @@ namespace boost { namespace geometry { namespace projections
             public :
                 virtual projection<Geographic, Cartesian>* create_new(const Parameters& par) const
                 {
-                    return new base_v_fi<ups_ellipsoid<Geographic, Cartesian, Parameters>, Geographic, Cartesian, Parameters>(par);
+                    if (par.es)
+                        return new base_v_fi<ups_ellipsoid<Geographic, Cartesian, Parameters>, Geographic, Cartesian, Parameters>(par);
+                    else
+                        return new base_v_fi<ups_spheroid<Geographic, Cartesian, Parameters>, Geographic, Cartesian, Parameters>(par);
                 }
         };
 

include/boost/geometry/extensions/gis/projections/proj/tmerc.hpp

@@ -284,28 +284,6 @@ namespace boost { namespace geometry { namespace projections
         }
     };
 
-    /*!
-        \brief Universal Transverse Mercator (UTM) projection
-        \ingroup projections
-        \tparam Geographic latlong point type
-        \tparam Cartesian xy point type
-        \tparam Parameters parameter type
-        \par Projection characteristics
-         - Cylindrical
-         - Spheroid
-         - zone= south
-        \par Example
-        \image html ex_utm.gif
-    */
-    template <typename Geographic, typename Cartesian, typename Parameters = parameters>
-    struct utm_ellipsoid : public detail::tmerc::base_tmerc_ellipsoid<Geographic, Cartesian, Parameters>
-    {
-        inline utm_ellipsoid(const Parameters& par) : detail::tmerc::base_tmerc_ellipsoid<Geographic, Cartesian, Parameters>(par)
-        {
-            detail::tmerc::setup_utm(this->m_par, this->m_proj_parm);
-        }
-    };
-
     /*!
         \brief Transverse Mercator projection
         \ingroup projections
@@ -328,6 +306,28 @@ namespace boost { namespace geometry { namespace projections
         }
     };
 
+    /*!
+        \brief Universal Transverse Mercator (UTM) projection
+        \ingroup projections
+        \tparam Geographic latlong point type
+        \tparam Cartesian xy point type
+        \tparam Parameters parameter type
+        \par Projection characteristics
+         - Cylindrical
+         - Spheroid
+         - zone= south
+        \par Example
+        \image html ex_utm.gif
+    */
+    template <typename Geographic, typename Cartesian, typename Parameters = parameters>
+    struct utm_ellipsoid : public detail::tmerc::base_tmerc_ellipsoid<Geographic, Cartesian, Parameters>
+    {
+        inline utm_ellipsoid(const Parameters& par) : detail::tmerc::base_tmerc_ellipsoid<Geographic, Cartesian, Parameters>(par)
+        {
+            detail::tmerc::setup_utm(this->m_par, this->m_proj_parm);
+        }
+    };
+
     /*!
         \brief Universal Transverse Mercator (UTM) projection
         \ingroup projections
@@ -374,7 +374,10 @@ namespace boost { namespace geometry { namespace projections
             public :
                 virtual projection<Geographic, Cartesian>* create_new(const Parameters& par) const
                 {
-                    return new base_v_fi<utm_ellipsoid<Geographic, Cartesian, Parameters>, Geographic, Cartesian, Parameters>(par);
+                    if (par.es)
+                        return new base_v_fi<utm_ellipsoid<Geographic, Cartesian, Parameters>, Geographic, Cartesian, Parameters>(par);
+                    else
+                        return new base_v_fi<utm_spheroid<Geographic, Cartesian, Parameters>, Geographic, Cartesian, Parameters>(par);
                 }
         };
 

include/boost/geometry/extensions/gis/projections/project_inverse_transformer.hpp

@@ -60,8 +60,8 @@ struct project_inverse_transformer
     {
         // Latlong (LL -> XY) will be projected, rest will be copied.
         // So first copy third or higher dimensions
-        geometry::detail::convert::point_to_point<Cartesian, LatLong, 2,
-                geometry::dimension<Cartesian>::value> ::copy(p1, p2);
+        geometry::detail::conversion::point_to_point<Cartesian, LatLong, 2,
+                geometry::dimension<Cartesian>::value> ::apply(p1, p2);
         return m_prj->inverse(p1, p2);
     }
 

include/boost/geometry/strategies/buffer.hpp

@@ -66,7 +66,8 @@ enum piece_type
     buffered_round_end,
     buffered_flat_end,
     buffered_point,
-    buffered_concave   // always on the inside
+    buffered_concave,   // always on the inside
+    piece_type_unknown
 };
 
 

include/boost/geometry/strategies/cartesian/side_of_intersection.hpp

@@ -13,8 +13,13 @@
 #include <boost/geometry/arithmetic/determinant.hpp>
 #include <boost/geometry/core/access.hpp>
 #include <boost/geometry/core/coordinate_type.hpp>
+#include <boost/geometry/algorithms/detail/assign_indexed_point.hpp>
+#include <boost/geometry/strategies/cartesian/side_by_triangle.hpp>
 #include <boost/geometry/util/math.hpp>
 
+#include <boost/math/common_factor_ct.hpp>
+#include <boost/math/common_factor_rt.hpp>
+#include <boost/multiprecision/cpp_int.hpp>
 
 namespace boost { namespace geometry
 {
@@ -29,15 +34,74 @@ namespace strategy { namespace side
 // It can be used for either integer (rescaled) points, and also for FP
 class side_of_intersection
 {
+private :
+    template <typename T, typename U>
+    static inline
+    int sign_of_product(T const& a, U const& b)
+    {
+        return a == 0 || b == 0 ? 0
+            : a > 0 && b > 0 ? 1
+            : a < 0 && b < 0 ? 1
+            : -1;
+    }
+
+    template <typename T>
+    static inline
+    int sign_of_compare(T const& a, T const& b, T const& c, T const& d)
+    {
+        // Both a*b and c*d are positive
+        // We have to judge if a*b > c*d
+
+        using namespace boost::multiprecision;
+        cpp_int const lab = cpp_int(a) * cpp_int(b);
+        cpp_int const lcd = cpp_int(c) * cpp_int(d);
+
+        return lab > lcd ? 1
+            :  lab < lcd ? -1
+            :  0
+            ;
+    }
+
+    template <typename T>
+    static inline
+    int sign_of_addition_of_two_products(T const& a, T const& b, T const& c, T const& d)
+    {
+        // sign of a*b+c*d, 1 if positive, -1 if negative, else 0
+        int const ab = sign_of_product(a, b);
+        int const cd = sign_of_product(c, d);
+        if (ab == 0)
+        {
+            return cd;
+        }
+        if (cd == 0)
+        {
+            return ab;
+        }
+
+        if (ab == cd)
+        {
+            // Both positive or both negative
+            return ab;
+        }
+
+        // One is positive, one is negative, both are non zero
+        // If ab is positive, we have to judge if a*b > -c*d (then 1 because sum is positive)
+        // If ab is negative, we have to judge if c*d > -a*b (idem)
+        return ab == 1
+            ? sign_of_compare(a, b, -c, d)
+            : sign_of_compare(c, d, -a, b);
+    }
+
+
 public :
 
     // Calculates the side of the intersection-point (if any) of
     // of segment a//b w.r.t. segment c
     // This is calculated without (re)calculating the IP itself again and fully
     // based on integer mathematics
-    template <typename T, typename Segment>
+    template <typename T, typename Segment, typename Point>
     static inline T side_value(Segment const& a, Segment const& b,
-                Segment const& c)
+                Segment const& c, Point const& fallback_point)
     {
         // The first point of the three segments is reused several times
         T const ax = get<0, 0>(a);
@@ -67,9 +131,15 @@ public :
         if (d == zero)
         {
             // There is no IP of a//b, they are collinear or parallel
-            // We don't have to divide but we can already conclude the side-value
-            // is meaningless and the resulting determinant will be 0
-            return zero;
+            // Assuming they intersect (this method should be called for
+            // segments known to intersect), they are collinear and overlap.
+            // They have one or two intersection points - we don't know and
+            // have to rely on the fallback intersection point
+
+            Point c1, c2;
+            geometry::detail::assign_point_from_index<0>(c, c1);
+            geometry::detail::assign_point_from_index<1>(c, c2);
+            return side_by_triangle<>::apply(c1, c2, fallback_point);
         }
 
         // Cramer's rule: da (see cart_intersect.hpp)
@@ -82,7 +152,9 @@ public :
         // IP is at (ax + (da/d) * dx_a, ay + (da/d) * dy_a)
         // Side of IP is w.r.t. c is: determinant(dx_c, dy_c, ipx-cx, ipy-cy)
         // We replace ipx by expression above and multiply each term by d
-        T const result = geometry::detail::determinant<T>
+
+#ifdef BOOST_GEOMETRY_SIDE_OF_INTERSECTION_DEBUG
+        T const result1 = geometry::detail::determinant<T>
                     (
                         dx_c * d,                   dy_c * d,
                         d * (ax - cx) + dx_a * da,  d * (ay - cy) + dy_a * da
@@ -93,15 +165,51 @@ public :
         // Therefore, the sign is always the same as that result, and the
         // resulting side (left,right,collinear) is the same
 
+        // The first row we divide again by d because of determinant multiply rule
+        T const result2 = d * geometry::detail::determinant<T>
+                    (
+                        dx_c,                   dy_c,
+                        d * (ax - cx) + dx_a * da,  d * (ay - cy) + dy_a * da
+                    );
+        // Write out:
+        T const result3 = d * (dx_c * (d * (ay - cy) + dy_a * da)
+                             - dy_c * (d * (ax - cx) + dx_a * da));
+        // Write out in braces:
+        T const result4 = d * (dx_c * d * (ay - cy) + dx_c * dy_a * da
+                             - dy_c * d * (ax - cx) - dy_c * dx_a * da);
+        // Write in terms of d * XX + da * YY
+        T const result5 = d * (d * (dx_c * (ay - cy) - dy_c * (ax - cx))
+                             + da * (dx_c * dy_a - dy_c * dx_a));
+
+        //return result;
+#endif
+
+        // We consider the results separately
+        // (in the end we only have to return the side-value 1,0 or -1)
+
+        // To avoid multiplications we judge the product (easy, avoids *d)
+        // and the sign of p*q+r*s (more elaborate)
+        T const result = sign_of_product
+            (
+                d,
+                sign_of_addition_of_two_products
+                    (
+                        d, dx_c * (ay - cy) - dy_c * (ax - cx),
+                        da, dx_c * dy_a - dy_c * dx_a
+                    )
+                );
         return result;
 
+
     }
 
-    template <typename Segment>
-    static inline int apply(Segment const& a, Segment const& b, Segment const& c)
+    template <typename Segment, typename Point>
+    static inline int apply(Segment const& a, Segment const& b,
+            Segment const& c,
+            Point const& fallback_point)
     {
         typedef typename geometry::coordinate_type<Segment>::type coordinate_type;
-        coordinate_type const s = side_value<coordinate_type>(a, b, c);
+        coordinate_type const s = side_value<coordinate_type>(a, b, c, fallback_point);
         coordinate_type const zero = coordinate_type();
         return math::equals(s, zero) ? 0
             : s > zero ? 1

include/boost/geometry/strategies/spherical/distance_cross_track.hpp

@@ -372,11 +372,16 @@ public :
         typedef typename return_type<Point, PointOfSegment>::type return_type;
 
 #ifdef BOOST_GEOMETRY_DEBUG_CROSS_TRACK
-        std::cout << "Course " << dsv(sp1) << " to " << dsv(p) << " " << crs_AD * geometry::math::r2d << std::endl;
-        std::cout << "Course " << dsv(sp1) << " to " << dsv(sp2) << " " << crs_AB * geometry::math::r2d << std::endl;
-        std::cout << "Course " << dsv(sp2) << " to " << dsv(p) << " " << crs_BD * geometry::math::r2d << std::endl;
-        std::cout << "Projection AD-AB " << projection1 << " : " << d_crs1 * geometry::math::r2d << std::endl;
-        std::cout << "Projection BD-BA " << projection2 << " : " << d_crs2 * geometry::math::r2d << std::endl;
+        std::cout << "Course " << dsv(sp1) << " to " << dsv(p) << " "
+                  << crs_AD * geometry::math::r2d<return_type>() << std::endl;
+        std::cout << "Course " << dsv(sp1) << " to " << dsv(sp2) << " "
+                  << crs_AB * geometry::math::r2d<return_type>() << std::endl;
+        std::cout << "Course " << dsv(sp2) << " to " << dsv(p) << " "
+                  << crs_BD * geometry::math::r2d << std::endl;
+        std::cout << "Projection AD-AB " << projection1 << " : "
+                  << d_crs1 * geometry::math::r2d<return_type>() << std::endl;
+        std::cout << "Projection BD-BA " << projection2 << " : "
+                  << d_crs2 * geometry::math::r2d<return_type>() << std::endl;
 #endif
 
         // http://williams.best.vwh.net/avform.htm#XTE

include/boost/geometry/strategies/strategy_transform.hpp

@@ -1,8 +1,13 @@
 // Boost.Geometry (aka GGL, Generic Geometry Library)
 
-// Copyright (c) 2007-2012 Barend Gehrels, Amsterdam, the Netherlands.
-// Copyright (c) 2008-2012 Bruno Lalande, Paris, France.
-// Copyright (c) 2009-2012 Mateusz Loskot, London, UK.
+// Copyright (c) 2007-2015 Barend Gehrels, Amsterdam, the Netherlands.
+// Copyright (c) 2008-2015 Bruno Lalande, Paris, France.
+// Copyright (c) 2009-2015 Mateusz Loskot, London, UK.
+
+// This file was modified by Oracle on 2015.
+// Modifications copyright (c) 2015 Oracle and/or its affiliates.
+
+// Contributed and/or modified by Menelaos Karavelas, on behalf of Oracle
 
 // Parts of Boost.Geometry are redesigned from Geodan's Geographic Library
 // (geolib/GGL), copyright (c) 1995-2010 Geodan, Amsterdam, the Netherlands.
@@ -28,6 +33,7 @@
 #include <boost/geometry/strategies/transform.hpp>
 
 #include <boost/geometry/util/math.hpp>
+#include <boost/geometry/util/promote_floating_point.hpp>
 #include <boost/geometry/util/select_coordinate_type.hpp>
 
 namespace boost { namespace geometry
@@ -130,7 +136,15 @@ struct degree_radian_vv
         assert_dimension<P1, 2>();
         assert_dimension<P2, 2>();
 
-        detail::transform_coordinates<P1, P2, 0, 2, F>::transform(p1, p2, math::d2r);
+        typedef typename promote_floating_point
+            <
+                typename select_coordinate_type<P1, P2>::type
+            >::type calculation_type;
+
+        detail::transform_coordinates
+            <
+                P1, P2, 0, 2, F
+            >::transform(p1, p2, math::d2r<calculation_type>());
         return true;
     }
 };
@@ -143,7 +157,16 @@ struct degree_radian_vv_3
         assert_dimension<P1, 3>();
         assert_dimension<P2, 3>();
 
-        detail::transform_coordinates<P1, P2, 0, 2, F>::transform(p1, p2, math::d2r);
+        typedef typename promote_floating_point
+            <
+                typename select_coordinate_type<P1, P2>::type
+            >::type calculation_type;
+
+        detail::transform_coordinates
+            <
+                P1, P2, 0, 2, F
+            >::transform(p1, p2, math::d2r<calculation_type>());
+
         // Copy height or other third dimension
         set<2>(p2, get<2>(p1));
         return true;

include/boost/geometry/strategies/transform/matrix_transformers.hpp

@@ -1,8 +1,13 @@
 // Boost.Geometry (aka GGL, Generic Geometry Library)
 
-// Copyright (c) 2007-2012 Barend Gehrels, Amsterdam, the Netherlands.
-// Copyright (c) 2008-2012 Bruno Lalande, Paris, France.
-// Copyright (c) 2009-2012 Mateusz Loskot, London, UK.
+// Copyright (c) 2007-2015 Barend Gehrels, Amsterdam, the Netherlands.
+// Copyright (c) 2008-2015 Bruno Lalande, Paris, France.
+// Copyright (c) 2009-2015 Mateusz Loskot, London, UK.
+
+// This file was modified by Oracle on 2015.
+// Modifications copyright (c) 2015 Oracle and/or its affiliates.
+
+// Contributed and/or modified by Menelaos Karavelas, on behalf of Oracle
 
 // Parts of Boost.Geometry are redesigned from Geodan's Geographic Library
 // (geolib/GGL), copyright (c) 1995-2010 Geodan, Amsterdam, the Netherlands.
@@ -40,6 +45,7 @@
 #include <boost/geometry/core/coordinate_dimension.hpp>
 #include <boost/geometry/core/cs.hpp>
 #include <boost/geometry/util/math.hpp>
+#include <boost/geometry/util/promote_floating_point.hpp>
 #include <boost/geometry/util/select_coordinate_type.hpp>
 #include <boost/geometry/util/select_most_precise.hpp>
 
@@ -340,7 +346,8 @@ struct as_radian<degree>
     template <typename T>
     static inline T get(T const& value)
     {
-        return value * math::d2r;
+        typedef typename promote_floating_point<T>::type promoted_type;
+        return value * math::d2r<promoted_type>();
     }
 
 };

include/boost/geometry/util/math.hpp

@@ -493,9 +493,20 @@ inline bool larger(T1 const& a, T2 const& b)
 }
 
 
+template <typename T>
+inline T d2r()
+{
+    static T const conversion_coefficient = geometry::math::pi<T>() / T(180.0);
+    return conversion_coefficient;
+}
+
+template <typename T>
+inline T r2d()
+{
+    static T const conversion_coefficient = T(180.0) / geometry::math::pi<T>();
+    return conversion_coefficient;
+}
 
-double const d2r = geometry::math::pi<double>() / 180.0;
-double const r2d = 1.0 / d2r;
 
 /*!
     \brief Calculates the haversine of an angle

include/boost/geometry/util/normalize_spheroidal_coordinates.hpp

@@ -0,0 +1,254 @@
+// Boost.Geometry (aka GGL, Generic Geometry Library)
+
+// Copyright (c) 2015, 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_UTIL_NORMALIZE_SPHEROIDAL_COORDINATES_HPP
+#define BOOST_GEOMETRY_UTIL_NORMALIZE_SPHEROIDAL_COORDINATES_HPP
+
+#include <boost/assert.hpp>
+
+#include <boost/geometry/core/cs.hpp>
+#include <boost/geometry/util/math.hpp>
+
+
+namespace boost { namespace geometry
+{
+
+namespace math 
+{
+
+#ifndef DOXYGEN_NO_DETAIL
+namespace detail
+{
+
+
+template <typename CoordinateType, typename Units>
+struct constants_on_spheroid
+{
+    static inline CoordinateType period()
+    {
+        return math::two_pi<CoordinateType>();
+    }
+
+    static inline CoordinateType half_period()
+    {
+        return math::pi<CoordinateType>();
+    }
+
+    static inline CoordinateType min_longitude()
+    {
+        static CoordinateType const minus_pi = -math::pi<CoordinateType>();
+        return minus_pi;
+    }
+
+    static inline CoordinateType max_longitude()
+    {
+        return math::pi<CoordinateType>();
+    }
+
+    static inline CoordinateType min_latitude()
+    {
+        static CoordinateType const minus_half_pi
+            = -math::half_pi<CoordinateType>();
+        return minus_half_pi;
+    }
+
+    static inline CoordinateType max_latitude()
+    {
+        return math::half_pi<CoordinateType>();
+    }
+};
+
+template <typename CoordinateType>
+struct constants_on_spheroid<CoordinateType, degree>
+{
+    static inline CoordinateType period()
+    {
+        return CoordinateType(360.0);
+    }
+
+    static inline CoordinateType half_period()
+    {
+        return CoordinateType(180.0);
+    }
+
+    static inline CoordinateType min_longitude()
+    {
+        return CoordinateType(-180.0);
+    }
+
+    static inline CoordinateType max_longitude()
+    {
+        return CoordinateType(180.0);
+    }
+
+    static inline CoordinateType min_latitude()
+    {
+        return CoordinateType(-90.0);
+    }
+
+    static inline CoordinateType max_latitude()
+    {
+        return CoordinateType(90.0);
+    }
+};
+
+
+template <typename Units, typename CoordinateType>
+class normalize_spheroidal_coordinates
+{
+    typedef constants_on_spheroid<CoordinateType, Units> constants;
+
+protected:
+    static inline CoordinateType normalize_up(CoordinateType const& value)
+    {
+        return
+            math::mod(value + constants::half_period(), constants::period())
+            - constants::half_period();            
+    }
+
+    static inline CoordinateType normalize_down(CoordinateType const& value)
+    {
+        return
+            math::mod(value - constants::half_period(), constants::period())
+            + constants::half_period();            
+    }
+
+public:
+    static inline void apply(CoordinateType& longitude,
+                             CoordinateType& latitude,
+                             bool normalize_poles = true)
+    {
+
+        // first normalize latitude
+        if (math::larger(latitude, constants::half_period()))
+        {
+            latitude = normalize_up(latitude);
+        }
+        else if (math::smaller(latitude, -constants::half_period()))
+        {
+            latitude = normalize_down(latitude);
+        }
+
+        // fix latitude range
+        if (latitude < constants::min_latitude())
+        {
+            latitude = -constants::half_period() - latitude;
+            longitude -= constants::half_period();
+        }
+        else if (latitude > constants::max_latitude())
+        {
+            latitude = constants::half_period() - latitude;
+            longitude -= constants::half_period();
+        }
+
+        // now normalize longitude
+        if (math::equals(longitude, -constants::half_period())
+            || math::equals(longitude, constants::half_period()))
+        {
+            longitude = constants::half_period();
+        }
+        else if (longitude > constants::half_period())
+        {
+            longitude = normalize_up(longitude);
+            if (math::equals(longitude, -constants::half_period()))
+            {
+                longitude = constants::half_period();
+            }
+        }
+        else if (longitude < -constants::half_period())
+        {
+            longitude = normalize_down(longitude);
+        }
+
+        // finally normalize poles
+        if (normalize_poles)
+        {
+            if (math::equals(math::abs(latitude), constants::max_latitude()))
+            {
+                // for the north and south pole we set the longitude to 0
+                // (works for both radians and degrees)
+                longitude = CoordinateType(0);
+            }
+        }
+
+        BOOST_ASSERT(! math::larger(constants::min_latitude(), latitude));
+        BOOST_ASSERT(! math::larger(latitude, constants::max_latitude()));
+
+        BOOST_ASSERT(math::smaller(constants::min_longitude(), longitude));
+        BOOST_ASSERT(! math::larger(longitude, constants::max_longitude()));
+    }
+};
+
+
+} // namespace detail
+#endif // DOXYGEN_NO_DETAIL
+
+
+template <typename Units, typename CoordinateType>
+inline CoordinateType min_longitude()
+{
+    return detail::constants_on_spheroid
+        <
+            CoordinateType, Units
+        >::min_longitude();
+}
+
+template <typename Units, typename CoordinateType>
+inline CoordinateType max_longitude()
+{
+    return detail::constants_on_spheroid
+        <
+            CoordinateType, Units
+        >::max_longitude();
+}
+
+template <typename Units, typename CoordinateType>
+inline CoordinateType min_latitude()
+{
+    return detail::constants_on_spheroid
+        <
+            CoordinateType, Units
+        >::min_latitude();
+}
+
+template <typename Units, typename CoordinateType>
+inline CoordinateType max_latitude()
+{
+    return detail::constants_on_spheroid
+        <
+            CoordinateType, Units
+        >::max_latitude();
+}
+
+
+/*!
+\brief Short utility to normalize the coordinates on a spheroid
+\tparam Units The units of the coordindate system in the spheroid
+\tparam CoordinateType The type of the coordinates
+\param longitude Longitude
+\param latitude Latitude
+\ingroup utility
+*/
+template <typename Units, typename CoordinateType>
+inline void normalize_spheroidal_coordinates(CoordinateType& longitude,
+                                             CoordinateType& latitude)
+{
+    detail::normalize_spheroidal_coordinates
+        <
+            Units, CoordinateType
+        >::apply(longitude, latitude);
+}
+
+
+} // namespace math
+
+
+}} // namespace boost::geometry
+
+#endif // BOOST_GEOMETRY_UTIL_NORMALIZE_SPHEROIDAL_COORDINATES_HPP

test/algorithms/Jamfile.v2

@@ -30,6 +30,7 @@ test-suite boost-geometry-algorithms
     [ run expand.cpp ]
     [ run expand_on_spheroid.cpp ]
     [ run for_each.cpp ]
+    [ run is_convex.cpp ]
     [ run is_simple.cpp ]
     [ run is_valid.cpp ]
     [ run is_valid_failure.cpp ]

test/algorithms/buffer/linestring_buffer.cpp

@@ -203,9 +203,10 @@ void test_all()
         double const d15 = 1.5;
         test_one<linestring, polygon>("mysql_report_2015_03_02c_asym1", mysql_report_2015_03_02c, join_round(7), end_round(7), 39.714, d10, d15);
         test_one<linestring, polygon>("mysql_report_2015_03_02c_asym2", mysql_report_2015_03_02c, join_round(7), end_round(7), 46.116, d15, d10);
-#if defined(BOOST_GEOMETRY_BUFFER_INCLUDE_FAILING_TESTS)
+
+#if defined(BOOST_GEOMETRY_BUFFER_USE_SIDE_OF_INTERSECTION)
         double const d100 = 10;
-        test_one<linestring, polygon>("mysql_report_2015_04_01", mysql_report_2015_04_01, join_round(32), end_round(32), 1.0/*NON ZERO*/, d100);
+        test_one<linestring, polygon>("mysql_report_2015_04_01", mysql_report_2015_04_01, join_round(32), end_round(32), 632.234, d100);
 #endif
     }
 

test/algorithms/buffer/multi_linestring_buffer.cpp

@@ -100,7 +100,9 @@ void test_all()
     test_one<multi_linestring_type, polygon>("mikado4_small", mikado4, join_round32, end_flat, 1930.785, 10.0);
 
 #ifdef BOOST_GEOMETRY_BUFFER_INCLUDE_FAILING_TESTS
-    test_one<multi_linestring_type, polygon>("mysql_15_04_10", mysql_15_04_10, join_round32, end_round32, 29151950703.779/*something big*/, 0x98);
+    // Cause of failure: coordinates vary so much that
+    // length = geometry::math::sqrt(dx * dx + dy * dy); returns a value of inf for length
+    test_one<multi_linestring_type, polygon>("mysql_15_04_10", mysql_15_04_10, join_round32, end_round32, 99999.999, 0.98);
 #endif
 }
 

test/algorithms/buffer/multi_point_buffer.cpp

@@ -30,7 +30,7 @@ void test_all()
     typedef bg::model::polygon<P, Clockwise> polygon;
     typedef bg::model::multi_point<P> multi_point_type;
 
-    bg::strategy::buffer::join_miter join_miter;
+    bg::strategy::buffer::join_round join;
     bg::strategy::buffer::end_flat end_flat;
     typedef bg::strategy::buffer::distance_symmetric
     <
@@ -40,41 +40,39 @@ void test_all()
 
     double const pi = boost::geometry::math::pi<double>();
 
-    test_one<multi_point_type, polygon>("simplex1", simplex, join_miter, end_flat, 2.0 * pi, 1.0, 1.0);
-    test_one<multi_point_type, polygon>("simplex2", simplex, join_miter, end_flat, 22.8372, 2.0, 2.0);
-    test_one<multi_point_type, polygon>("simplex3", simplex, join_miter, end_flat, 44.5692, 3.0, 3.0);
+    test_one<multi_point_type, polygon>("simplex1", simplex, join, end_flat, 2.0 * pi, 1.0, 1.0);
+    test_one<multi_point_type, polygon>("simplex2", simplex, join, end_flat, 22.8372, 2.0, 2.0);
+    test_one<multi_point_type, polygon>("simplex3", simplex, join, end_flat, 44.5692, 3.0, 3.0);
 
-    test_one<multi_point_type, polygon>("three1", three, join_miter, end_flat, 3.0 * pi, 1.0, 1.0);
-    test_one<multi_point_type, polygon>("three2", three, join_miter, end_flat, 36.7592, 2.0, 2.0);
-    test_one<multi_point_type, polygon>("three19", three, join_miter, end_flat, 33.6914, 1.9, 1.9);
-    test_one<multi_point_type, polygon>("three21", three, join_miter, end_flat, 39.6394, 2.1, 2.1);
-    test_one<multi_point_type, polygon>("three3", three, join_miter, end_flat, 65.533, 3.0, 3.0);
-
-    test_one<multi_point_type, polygon>("multipoint_a", multipoint_a, join_miter, end_flat, 2049.98, 14.0, 14.0);
-    test_one<multi_point_type, polygon>("multipoint_b", multipoint_b, join_miter, end_flat, 7109.88, 15.0, 15.0);
-    test_one<multi_point_type, polygon>("multipoint_b1", multipoint_b, join_miter, end_flat, 6911.89, 14.7, 14.7);
-    test_one<multi_point_type, polygon>("multipoint_b2", multipoint_b, join_miter, end_flat, 7174.79, 15.1, 15.1);
+    test_one<multi_point_type, polygon>("three1", three, join, end_flat, 3.0 * pi, 1.0, 1.0);
+    test_one<multi_point_type, polygon>("three2", three, join, end_flat, 36.7592, 2.0, 2.0);
+    test_one<multi_point_type, polygon>("three19", three, join, end_flat, 33.6914, 1.9, 1.9);
+    test_one<multi_point_type, polygon>("three21", three, join, end_flat, 39.6394, 2.1, 2.1);
+    test_one<multi_point_type, polygon>("three3", three, join, end_flat, 65.533, 3.0, 3.0);
 
+    test_one<multi_point_type, polygon>("multipoint_a", multipoint_a, join, end_flat, 2049.98, 14.0, 14.0);
+    test_one<multi_point_type, polygon>("multipoint_b", multipoint_b, join, end_flat, 7109.88, 15.0, 15.0);
+    test_one<multi_point_type, polygon>("multipoint_b1", multipoint_b, join, end_flat, 6911.89, 14.7, 14.7);
+    test_one<multi_point_type, polygon>("multipoint_b2", multipoint_b, join, end_flat, 7174.79, 15.1, 15.1);
 
     // Grid tests
     {
         bg::strategy::buffer::point_square point_strategy;
 
-
         test_with_custom_strategies<multi_point_type, polygon>("grid_a50",
-                grid_a, join_miter, end_flat,
+                grid_a, join, end_flat,
                 distance_strategy(0.5), side_strategy, point_strategy, 7.0);
 
-#if defined(BOOST_GEOMETRY_BUFFER_INCLUDE_FAILING_TESTS)
+#if defined(BOOST_GEOMETRY_BUFFER_USE_SIDE_OF_INTERSECTION)
         test_with_custom_strategies<multi_point_type, polygon>("grid_a54",
-                grid_a, join_miter, end_flat,
-                distance_strategy(0.54), side_strategy, point_strategy, 99);
+                grid_a, join, end_flat,
+                distance_strategy(0.54), side_strategy, point_strategy, 7.819);
 #endif
 
     }
 
     test_with_custom_strategies<multi_point_type, polygon>("mysql_report_2015_02_25_1_800",
-            mysql_report_2015_02_25_1, join_miter, end_flat,
+            mysql_report_2015_02_25_1, join, end_flat,
             distance_strategy(6051788), side_strategy,
             bg::strategy::buffer::point_circle(800), 115057490003226.125, 1.0);
 }
@@ -92,7 +90,7 @@ void test_many_points_per_circle()
     typedef bg::model::polygon<P, false> polygon;
     typedef bg::model::multi_point<P> multi_point_type;
 
-    bg::strategy::buffer::join_miter join_miter;
+    bg::strategy::buffer::join_round join;
     bg::strategy::buffer::end_flat end_flat;
     typedef bg::strategy::buffer::distance_symmetric
     <
@@ -104,41 +102,62 @@ void test_many_points_per_circle()
 
     double const tolerance = 1.0;
 
+    // Area should be somewhat larger (~>) than pi*distance^2
+    // 6051788: area ~> 115058122875258
+
     // Strategies with many points, which are (very) slow in debug mode
     test_with_custom_strategies<multi_point_type, polygon>(
             "mysql_report_2015_02_25_1_8000",
-            mysql_report_2015_02_25_1, join_miter, end_flat,
+            mysql_report_2015_02_25_1, join, end_flat,
             distance_strategy(6051788), side_strategy, point_circle(8000),
             115058661065242.812, tolerance);
 
     test_with_custom_strategies<multi_point_type, polygon>(
             "mysql_report_2015_02_25_1",
-            mysql_report_2015_02_25_1, join_miter, end_flat,
+            mysql_report_2015_02_25_1, join, end_flat,
             distance_strategy(6051788), side_strategy, point_circle(83585),
             115058672785611.219, tolerance);
 
-    // Takes about 20 seconds in release mode
+    // Takes about 7 seconds in release mode
     test_with_custom_strategies<multi_point_type, polygon>(
             "mysql_report_2015_02_25_1_250k",
-            mysql_report_2015_02_25_1, join_miter, end_flat,
+            mysql_report_2015_02_25_1, join, end_flat,
             distance_strategy(6051788), side_strategy, point_circle(250000),
             115058672880671.531, tolerance);
 
-#if defined(BOOST_GEOMETRY_BUFFER_INCLUDE_FAILING_TESTS)
-    // Takes too long, TODO improve turn_in_piece_visitor
+#if defined(BOOST_GEOMETRY_BUFFER_INCLUDE_SLOW_TESTS)
+    // Takes about 110 seconds in release mode
     test_with_custom_strategies<multi_point_type, polygon>(
-            "mysql_report_2015_02_25_1",
-            mysql_report_2015_02_25_1, join_miter, end_flat,
+            "mysql_report_2015_02_25_1_800k",
+            mysql_report_2015_02_25_1, join, end_flat,
             distance_strategy(6051788), side_strategy, point_circle(800000),
-            115058672799999.999, tolerance); // area to be determined
+            115058672871849.219, tolerance);
 #endif
 
+    // 5666962: area ~> 100890546298964
     test_with_custom_strategies<multi_point_type, polygon>(
             "mysql_report_2015_02_25_2",
-            mysql_report_2015_02_25_2, join_miter, end_flat,
+            mysql_report_2015_02_25_2, join, end_flat,
             distance_strategy(5666962), side_strategy, point_circle(46641),
             100891031341757.344, tolerance);
 
+    // Multipoint b with large distances/many points
+    // Area ~> pi * 10x
+    test_with_custom_strategies<multi_point_type, polygon>(
+            "multipoint_b_50k",
+            multipoint_b, join, end_flat,
+            distance_strategy(1000000), side_strategy, point_circle(50000),
+            3141871558215.26465, tolerance);
+
+#if defined(BOOST_GEOMETRY_BUFFER_INCLUDE_SLOW_TESTS)
+    // Tests optimization min/max radius
+    // Takes about 55 seconds in release mode
+    test_with_custom_strategies<multi_point_type, polygon>(
+            "multipoint_b_500k",
+            multipoint_b, join, end_flat,
+            distance_strategy(10000000), side_strategy, point_circle(500000),
+            314162054419515.562, tolerance);
+#endif
 }
 
 int test_main(int, char* [])

test/algorithms/buffer/multi_polygon_buffer.cpp

@@ -413,7 +413,11 @@ void test_all()
     test_one<multi_polygon_type, polygon_type>("rt_p22", rt_p22, join_miter, end_flat, 26.5711, 1.0);
 
     test_one<multi_polygon_type, polygon_type>("rt_q1", rt_q1, join_miter, end_flat, 27, 1.0);
+#if ! defined(BOOST_GEOMETRY_BUFFER_USE_SIDE_OF_INTERSECTION)
+    // Fails because all IP's on external borders are now selectable, two are on top of each other
+    // and this should be fixed differently
     test_one<multi_polygon_type, polygon_type>("rt_q2", rt_q2, join_miter, end_flat, 26.4853, 1.0);
+#endif
     test_one<multi_polygon_type, polygon_type>("rt_q2", rt_q2, join_miter, end_flat, 0.9697, -0.25);
 
     test_one<multi_polygon_type, polygon_type>("rt_r", rt_r, join_miter, end_flat, 21.0761, 1.0);
@@ -449,8 +453,8 @@ void test_all()
     test_one<multi_polygon_type, polygon_type>("rt_u11_50", rt_u11, join_miter, end_flat, 0.04289, -0.50);
     test_one<multi_polygon_type, polygon_type>("rt_u11_25", rt_u11, join_miter, end_flat, 10.1449, -0.25);
 
-#if defined(BOOST_GEOMETRY_BUFFER_INCLUDE_FAILING_TESTS)
-    test_one<multi_polygon_type, polygon_type>("rt_u12", rt_u12, join_miter, end_flat, 999, 1.0);
+#if defined(BOOST_GEOMETRY_BUFFER_USE_SIDE_OF_INTERSECTION)
+    test_one<multi_polygon_type, polygon_type>("rt_u12", rt_u12, join_miter, end_flat, 142.1348, 1.0);
     test_one<multi_polygon_type, polygon_type>("rt_u13", rt_u13, join_miter, end_flat, 115.4853, 1.0);
 #endif
 

test/algorithms/buffer/test_buffer.hpp

@@ -1,7 +1,7 @@
 // Boost.Geometry (aka GGL, Generic Geometry Library)
-// Unit Test
+// Unit Test Helper
 
-// Copyright (c) 2010-2014 Barend Gehrels, Amsterdam, the Netherlands.
+// Copyright (c) 2010-2015 Barend Gehrels, Amsterdam, the Netherlands.
 // Use, modification and distribution is subject to the Boost Software License,
 // Version 1.0. (See accompanying file LICENSE_1_0.txt or copy at
 // http://www.boost.org/LICENSE_1_0.txt)
@@ -10,29 +10,18 @@
 #ifndef BOOST_GEOMETRY_TEST_BUFFER_HPP
 #define BOOST_GEOMETRY_TEST_BUFFER_HPP
 
+#if defined(TEST_WITH_SVG)
+    // Define before including any buffer headerfile
+    #define BOOST_GEOMETRY_BUFFER_USE_HELPER_POINTS
+#endif
+
 #include <iostream>
 #include <fstream>
 #include <iomanip>
 
-#if defined(TEST_WITH_SVG)
-#define BOOST_GEOMETRY_BUFFER_USE_HELPER_POINTS
-
-// Uncomment next lines if you want to have a zoomed view
-// #define BOOST_GEOMETRY_BUFFER_TEST_SVG_USE_ALTERNATE_BOX
-
-// If possible define box before including this unit with the right view
-#ifdef BOOST_GEOMETRY_BUFFER_TEST_SVG_USE_ALTERNATE_BOX
-#  ifndef BOOST_GEOMETRY_BUFFER_TEST_SVG_ALTERNATE_BOX
-#    define BOOST_GEOMETRY_BUFFER_TEST_SVG_ALTERNATE_BOX "BOX(0 0,100 100)"
-#  endif
-#endif
-
-#endif // TEST_WITH_SVG
-
 #include <boost/foreach.hpp>
 #include <geometry_test_common.hpp>
 
-
 #include <boost/geometry/algorithms/envelope.hpp>
 #include <boost/geometry/algorithms/area.hpp>
 #include <boost/geometry/algorithms/buffer.hpp>
@@ -57,296 +46,8 @@
 
 
 #if defined(TEST_WITH_SVG)
-
-#include <boost/geometry/io/svg/svg_mapper.hpp>
-
-
-inline char piece_type_char(bg::strategy::buffer::piece_type const& type)
-{
-    using namespace bg::strategy::buffer;
-    switch(type)
-    {
-        case buffered_segment : return 's';
-        case buffered_join : return 'j';
-        case buffered_round_end : return 'r';
-        case buffered_flat_end : return 'f';
-        case buffered_point : return 'p';
-        case buffered_concave : return 'c';
-        default : return '?';
-    }
-}
-
-template <typename Geometry, typename Mapper, typename RescalePolicy>
-void post_map(Geometry const& geometry, Mapper& mapper, RescalePolicy const& rescale_policy)
-{
-    typedef typename bg::point_type<Geometry>::type point_type;
-    typedef bg::detail::overlay::turn_info
-    <
-        point_type,
-        typename bg::segment_ratio_type<point_type, RescalePolicy>::type
-    > turn_info;
-
-    std::vector<turn_info> turns;
-
-    bg::detail::self_get_turn_points::no_interrupt_policy policy;
-    bg::self_turns
-        <
-            bg::detail::overlay::assign_null_policy
-        >(geometry, rescale_policy, turns, policy);
-
-    BOOST_FOREACH(turn_info const& turn, turns)
-    {
-        mapper.map(turn.point, "fill:rgb(255,128,0);stroke:rgb(0,0,100);stroke-width:1", 3);
-    }
-}
-
-template <typename SvgMapper, typename Box, typename Tag>
-struct svg_visitor
-{
-#ifdef BOOST_GEOMETRY_BUFFER_TEST_SVG_USE_ALTERNATE_BOX
-    Box m_alternate_box;
-#endif
-
-    class si
-    {
-    private :
-        bg::segment_identifier m_id;
-
-    public :
-        inline si(bg::segment_identifier const& id)
-            : m_id(id)
-        {}
-
-        template <typename Char, typename Traits>
-        inline friend std::basic_ostream<Char, Traits>& operator<<(
-                std::basic_ostream<Char, Traits>& os,
-                si const& s)
-        {
-            os << s.m_id.multi_index << "." << s.m_id.segment_index;
-            return os;
-        }
-    };
-
-
-    SvgMapper& m_mapper;
-
-    svg_visitor(SvgMapper& mapper)
-        : m_mapper(mapper)
-    {}
-
-    template <typename Turns>
-    inline void map_turns(Turns const& turns, bool label_good_turns, bool label_wrong_turns)
-    {
-        namespace bgdb = boost::geometry::detail::buffer;
-        typedef typename boost::range_value<Turns const>::type turn_type;
-        typedef typename turn_type::point_type point_type;
-        typedef typename turn_type::robust_point_type robust_point_type;
-
-        std::map<robust_point_type, int, bg::less<robust_point_type> > offsets;
-
-        for (typename boost::range_iterator<Turns const>::type it =
-            boost::begin(turns); it != boost::end(turns); ++it)
-        {
-#ifdef BOOST_GEOMETRY_BUFFER_TEST_SVG_USE_ALTERNATE_BOX
-            if (bg::disjoint(it->point, m_alternate_box))
-            {
-                continue;
-            }
-#endif
-
-            bool is_good = true;
-            char color = 'g';
-            std::string fill = "fill:rgb(0,255,0);";
-            switch(it->location)
-            {
-                case bgdb::inside_buffer :
-                    fill = "fill:rgb(255,0,0);";
-                    color = 'r';
-                    is_good = false;
-                    break;
-                case bgdb::location_discard :
-                    fill = "fill:rgb(0,0,255);";
-                    color = 'b';
-                    is_good = false;
-                    break;
-            }
-            if (!it->selectable_start)
-            {
-                fill = "fill:rgb(255,192,0);";
-                color = 'o'; // orange
-            }
-            if (it->blocked())
-            {
-                fill = "fill:rgb(128,128,128);";
-                color = '-';
-                is_good = false;
-            }
-
-            fill += "fill-opacity:0.7;";
-
-            m_mapper.map(it->point, fill, 4);
-
-            if ((label_good_turns && is_good) || (label_wrong_turns && ! is_good))
-            {
-                std::ostringstream out;
-                out << it->turn_index
-                    << " " << it->operations[0].piece_index << "/" << it->operations[1].piece_index
-                    << " " << si(it->operations[0].seg_id) << "/" << si(it->operations[1].seg_id)
-
-    //              If you want to see travel information
-                    << std::endl
-                    << " nxt " << it->operations[0].enriched.travels_to_ip_index
-                    << "/" << it->operations[1].enriched.travels_to_ip_index
-                    << " or " << it->operations[0].enriched.next_ip_index
-                    << "/" << it->operations[1].enriched.next_ip_index
-                    //<< " frac " << it->operations[0].fraction
-
-    //                If you want to see robust-point coordinates (e.g. to find duplicates)
-    //                << std::endl
-    //                << " " << bg::get<0>(it->robust_point) << " , " << bg::get<1>(it->robust_point)
-
-                    << std::endl;
-                out << " " << bg::method_char(it->method)
-                    << ":" << bg::operation_char(it->operations[0].operation)
-                    << "/" << bg::operation_char(it->operations[1].operation);
-                out << " "
-                    << (it->count_on_offsetted > 0 ? "b" : "") // b: offsetted border
-                    << (it->count_within_near_offsetted > 0 ? "n" : "")
-                    << (it->count_within > 0 ? "w" : "")
-                    << (it->count_on_helper > 0 ? "h" : "")
-                    << (it->count_on_multi > 0 ? "m" : "")
-                    ;
-
-                offsets[it->get_robust_point()] += 10;
-                int offset = offsets[it->get_robust_point()];
-
-                m_mapper.text(it->point, out.str(), "fill:rgb(0,0,0);font-family='Arial';font-size:9px;", 5, offset);
-
-                offsets[it->get_robust_point()] += 25;
-            }
-        }
-    }
-
-    template <typename Pieces, typename OffsettedRings>
-    inline void map_pieces(Pieces const& pieces,
-                OffsettedRings const& offsetted_rings,
-                bool do_pieces, bool do_indices)
-    {
-        typedef typename boost::range_value<Pieces const>::type piece_type;
-        typedef typename boost::range_value<OffsettedRings const>::type ring_type;
-
-        for(typename boost::range_iterator<Pieces const>::type it = boost::begin(pieces);
-            it != boost::end(pieces);
-            ++it)
-        {
-            const piece_type& piece = *it;
-            bg::segment_identifier seg_id = piece.first_seg_id;
-            if (seg_id.segment_index < 0)
-            {
-                continue;
-            }
-
-            ring_type corner;
-
-
-            ring_type const& ring = offsetted_rings[seg_id.multi_index];
-
-            std::copy(boost::begin(ring) + seg_id.segment_index,
-                    boost::begin(ring) + piece.last_segment_index,
-                    std::back_inserter(corner));
-            std::copy(boost::begin(piece.helper_points),
-                    boost::end(piece.helper_points),
-                    std::back_inserter(corner));
-
-            if (corner.empty())
-            {
-                continue;
-            }
-#ifdef BOOST_GEOMETRY_BUFFER_TEST_SVG_USE_ALTERNATE_BOX
-            if (bg::disjoint(corner, m_alternate_box))
-            {
-                continue;
-            }
-#endif
-
-            if (do_pieces)
-            {
-                std::string style = "opacity:0.3;stroke:rgb(0,0,0);stroke-width:1;";
-                m_mapper.map(corner,
-                    piece.type == bg::strategy::buffer::buffered_segment
-                    ? style + "fill:rgb(255,128,0);"
-                    : style + "fill:rgb(255,0,0);");
-            }
-
-            if (do_indices)
-            {
-                // Label starting piece_index / segment_index
-                typedef typename bg::point_type<ring_type>::type point_type;
-
-                std::ostringstream out;
-                out << piece.index << " (" << piece_type_char(piece.type) << ") " << piece.first_seg_id.segment_index << ".." << piece.last_segment_index - 1;
-                point_type label_point = bg::return_centroid<point_type>(corner);
-
-                if ((piece.type == bg::strategy::buffer::buffered_concave
-                     || piece.type == bg::strategy::buffer::buffered_flat_end)
-                    && corner.size() >= 2u)
-                {
-                    bg::set<0>(label_point, (bg::get<0>(corner[0]) + bg::get<0>(corner[1])) / 2.0);
-                    bg::set<1>(label_point, (bg::get<1>(corner[0]) + bg::get<1>(corner[1])) / 2.0);
-                }
-                m_mapper.text(label_point, out.str(), "fill:rgb(255,0,0);font-family='Arial';font-size:10px;", 5, 5);
-            }
-        }
-    }
-
-    template <typename TraversedRings>
-    inline void map_traversed_rings(TraversedRings const& traversed_rings)
-    {
-        for(typename boost::range_iterator<TraversedRings const>::type it
-                = boost::begin(traversed_rings); it != boost::end(traversed_rings); ++it)
-        {
-            m_mapper.map(*it, "opacity:0.4;fill:none;stroke:rgb(0,255,0);stroke-width:2");
-        }
-    }
-
-    template <typename OffsettedRings>
-    inline void map_offsetted_rings(OffsettedRings const& offsetted_rings)
-    {
-        for(typename boost::range_iterator<OffsettedRings const>::type it
-                = boost::begin(offsetted_rings); it != boost::end(offsetted_rings); ++it)
-        {
-            if (it->discarded())
-            {
-                m_mapper.map(*it, "opacity:0.4;fill:none;stroke:rgb(255,0,0);stroke-width:2");
-            }
-            else
-            {
-                m_mapper.map(*it, "opacity:0.4;fill:none;stroke:rgb(0,0,255);stroke-width:2");
-            }
-        }
-    }
-
-    template <typename PieceCollection>
-    inline void apply(PieceCollection const& collection, int phase)
-    {
-        // Comment next return if you want to see pieces, turns, etc.
-        return;
-
-        if(phase == 0)
-        {
-            map_pieces(collection.m_pieces, collection.offsetted_rings, true, true);
-            map_turns(collection.m_turns, true, false);
-        }
-#if !defined(BOOST_GEOMETRY_BUFFER_TEST_SVG_USE_ALTERNATE_BOX)
-        if (phase == 1)
-        {
-//        map_traversed_rings(collection.traversed_rings);
-//        map_offsetted_rings(collection.offsetted_rings);
-        }
-#endif
-    }
-};
-
+#  include <test_buffer_svg.hpp>
+#  include <test_buffer_svg_per_turn.hpp>
 #endif
 
 //-----------------------------------------------------------------------------
@@ -472,40 +173,25 @@ void test_buffer(std::string const& caseid, Geometry const& geometry,
 
     //std::cout << complete.str() << std::endl;
 
+#if defined(TEST_WITH_SVG_PER_TURN)
+    save_turns_visitor<point_type> visitor;
+#elif defined(TEST_WITH_SVG)
+
+    buffer_svg_mapper<point_type> buffer_mapper(complete.str());
+
     std::ostringstream filename;
     filename << "buffer_" << complete.str() << ".svg";
-
-#if defined(TEST_WITH_SVG)
     std::ofstream svg(filename.str().c_str());
     typedef bg::svg_mapper<point_type> mapper_type;
-    mapper_type mapper(svg, 1000, 1000);
+    mapper_type mapper(svg, 1000, 800);
 
-    svg_visitor<mapper_type, bg::model::box<point_type>, tag> visitor(mapper);
-
-#ifdef BOOST_GEOMETRY_BUFFER_TEST_SVG_USE_ALTERNATE_BOX
-    // Create a zoomed-in view
-    bg::model::box<point_type> alternate_box;
-    bg::read_wkt(BOOST_GEOMETRY_BUFFER_TEST_SVG_ALTERNATE_BOX, alternate_box);
-    mapper.add(alternate_box);
-
-    // Take care non-visible elements are skipped
-    visitor.m_alternate_box = alternate_box;
-#else
-
-    {
-        bg::model::box<point_type> box = envelope;
-        if (distance_strategy.negative())
-        {
-            bg::buffer(box, box, 1.0);
-        }
-        else
-        {
-            bg::buffer(box, box, 1.1 * distance_strategy.max_distance(join_strategy, end_strategy));
-        }
-        mapper.add(box);
-    }
-#endif
+    svg_visitor<mapper_type, bg::model::box<point_type> > visitor(mapper);
 
+    buffer_mapper.prepare(mapper, visitor, envelope,
+            distance_strategy.negative()
+            ? 1.0
+            : 1.1 * distance_strategy.max_distance(join_strategy, end_strategy)
+        );
 #else
     bg::detail::buffer::visit_pieces_default_policy visitor;
 #endif
@@ -557,14 +243,22 @@ void test_buffer(std::string const& caseid, Geometry const& geometry,
 
         if (check_self_intersections)
         {
-            // Be sure resulting polygon does not contain
-            // self-intersections
-            BOOST_CHECK_MESSAGE
-                (
-                    ! bg::detail::overlay::has_self_intersections(buffered,
-                            rescale_policy, false),
-                    complete.str() << " output is self-intersecting. "
-                );
+            try
+            {
+                bool has_self_ips = bg::detail::overlay::has_self_intersections(buffered,
+                rescale_policy, false);
+                // Be sure resulting polygon does not contain
+                // self-intersections
+                BOOST_CHECK_MESSAGE
+                    (
+                        ! has_self_ips,
+                        complete.str() << " output is self-intersecting. "
+                    );
+            }
+            catch(...)
+            {
+                BOOST_MESSAGE("Exception in checking self-intersections");
+            }
         }
     }
 
@@ -579,48 +273,28 @@ void test_buffer(std::string const& caseid, Geometry const& geometry,
 //    BOOST_CHECK_MESSAGE(bg::intersects(buffered) == false, complete.str() <<  " intersects");
 #endif
 
-#if defined(TEST_WITH_SVG)
-    bool const areal = boost::is_same
-        <
-            typename bg::tag_cast<tag, bg::areal_tag>::type, bg::areal_tag
-        >::type::value;
-
-    // Map input geometry in green
-    if (areal)
+#if defined(TEST_WITH_SVG_PER_TURN)
     {
-#ifdef BOOST_GEOMETRY_BUFFER_TEST_SVG_USE_ALTERNATE_BOX_FOR_INPUT
-        // Assuming input is areal
-        bg::model::multi_polygon<GeometryOut> clipped;
-        bg::intersection(geometry, alternate_box, clipped);
-        mapper.map(clipped, "opacity:0.5;fill:rgb(0,128,0);stroke:rgb(0,64,0);stroke-width:2");
-#else
-        mapper.map(geometry, "opacity:0.5;fill:rgb(0,128,0);stroke:rgb(0,64,0);stroke-width:2");
+        // Create a per turn visitor to map per turn, and buffer again with it
+        per_turn_visitor<point_type> ptv(complete.str(), visitor.get_points());
+        bg::detail::buffer::buffer_inserter<GeometryOut>(geometry,
+                            std::back_inserter(buffered),
+                            distance_strategy,
+                            side_strategy,
+                            join_strategy,
+                            end_strategy,
+                            point_strategy,
+                            rescale_policy,
+                            ptv);
+        ptv.map_input_output(geometry, buffered, distance_strategy.negative());
+        // self_ips NYI here
+    }
+#elif defined(TEST_WITH_SVG)
+    buffer_mapper.map_input_output(mapper, geometry, buffered, distance_strategy.negative());
+    buffer_mapper.map_self_ips(mapper, buffered, rescale_policy);
 #endif
-    }
-    else
-    {
-        // TODO: clip input points/linestring
-        mapper.map(geometry, "opacity:0.5;stroke:rgb(0,128,0);stroke-width:10");
-    }
-
-    {
-        // Map buffer in yellow (inflate) and with orange-dots (deflate)
-        std::string style = distance_strategy.negative()
-            ? "opacity:0.4;fill:rgb(255,255,192);stroke:rgb(255,128,0);stroke-width:3"
-            : "opacity:0.4;fill:rgb(255,255,128);stroke:rgb(0,0,0);stroke-width:3";
-
-#ifdef BOOST_GEOMETRY_BUFFER_TEST_SVG_USE_ALTERNATE_BOX
-        // Clip output multi-polygon with box
-        bg::model::multi_polygon<GeometryOut> clipped;
-        bg::intersection(buffered, alternate_box, clipped);
-        mapper.map(clipped, style);
-#else
-        mapper.map(buffered, style);
-#endif
-    }
-    post_map(buffered, mapper, rescale_policy);
-#endif // TEST_WITH_SVG
 
+    // Check for self-intersections
     if (self_ip_count != NULL)
     {
         std::size_t count = 0;

test/algorithms/buffer/test_buffer_svg.hpp

@@ -0,0 +1,468 @@
+// Boost.Geometry (aka GGL, Generic Geometry Library)
+// Unit Test Helper
+
+// Copyright (c) 2010-2015 Barend Gehrels, Amsterdam, the Netherlands.
+// Use, modification and distribution is subject to the Boost Software License,
+// Version 1.0. (See accompanying file LICENSE_1_0.txt or copy at
+// http://www.boost.org/LICENSE_1_0.txt)
+
+
+#ifndef BOOST_GEOMETRY_TEST_BUFFER_SVG_HPP
+#define BOOST_GEOMETRY_TEST_BUFFER_SVG_HPP
+
+#include <fstream>
+#include <sstream>
+
+// Uncomment next lines if you want to have a zoomed view
+//#define BOOST_GEOMETRY_BUFFER_TEST_SVG_USE_ALTERNATE_BOX
+
+// If possible define box before including this unit with the right view
+#ifdef BOOST_GEOMETRY_BUFFER_TEST_SVG_USE_ALTERNATE_BOX
+#  ifndef BOOST_GEOMETRY_BUFFER_TEST_SVG_ALTERNATE_BOX
+#    define BOOST_GEOMETRY_BUFFER_TEST_SVG_ALTERNATE_BOX "BOX(0 0,100 100)"
+#  endif
+#endif
+
+#include <boost/foreach.hpp>
+#include <boost/geometry/io/svg/svg_mapper.hpp>
+
+
+inline char piece_type_char(bg::strategy::buffer::piece_type const& type)
+{
+    using namespace bg::strategy::buffer;
+    switch(type)
+    {
+        case buffered_segment : return 's';
+        case buffered_join : return 'j';
+        case buffered_round_end : return 'r';
+        case buffered_flat_end : return 'f';
+        case buffered_point : return 'p';
+        case buffered_concave : return 'c';
+        default : return '?';
+    }
+}
+
+template <typename SvgMapper, typename Box>
+class svg_visitor
+{
+public :
+    svg_visitor(SvgMapper& mapper)
+        : m_mapper(mapper)
+        , m_zoom(false)
+    {}
+
+    void set_alternate_box(Box const& box)
+    {
+        m_alternate_box = box;
+        m_zoom = true;
+    }
+
+    template <typename PieceCollection>
+    inline void apply(PieceCollection const& collection, int phase)
+    {
+        // Comment next return if you want to see pieces, turns, etc.
+        return;
+
+        if(phase == 0)
+        {
+            map_pieces(collection.m_pieces, collection.offsetted_rings, true, true);
+            map_turns(collection.m_turns, true, false);
+        }
+        if (phase == 1 && ! m_zoom)
+        {
+//        map_traversed_rings(collection.traversed_rings);
+//        map_offsetted_rings(collection.offsetted_rings);
+        }
+    }
+
+private :
+    class si
+    {
+    private :
+        bg::segment_identifier m_id;
+
+    public :
+        inline si(bg::segment_identifier const& id)
+            : m_id(id)
+        {}
+
+        template <typename Char, typename Traits>
+        inline friend std::basic_ostream<Char, Traits>& operator<<(
+                std::basic_ostream<Char, Traits>& os,
+                si const& s)
+        {
+            os << s.m_id.multi_index << "." << s.m_id.segment_index;
+            return os;
+        }
+    };
+
+    template <typename Turns>
+    inline void map_turns(Turns const& turns, bool label_good_turns, bool label_wrong_turns)
+    {
+        namespace bgdb = boost::geometry::detail::buffer;
+        typedef typename boost::range_value<Turns const>::type turn_type;
+        typedef typename turn_type::point_type point_type;
+        typedef typename turn_type::robust_point_type robust_point_type;
+
+        std::map<robust_point_type, int, bg::less<robust_point_type> > offsets;
+
+        for (typename boost::range_iterator<Turns const>::type it =
+            boost::begin(turns); it != boost::end(turns); ++it)
+        {
+            if (m_zoom && bg::disjoint(it->point, m_alternate_box))
+            {
+                continue;
+            }
+
+            bool is_good = true;
+            char color = 'g';
+            std::string fill = "fill:rgb(0,255,0);";
+            switch(it->location)
+            {
+                case bgdb::inside_buffer :
+                    fill = "fill:rgb(255,0,0);";
+                    color = 'r';
+                    is_good = false;
+                    break;
+                case bgdb::location_discard :
+                    fill = "fill:rgb(0,0,255);";
+                    color = 'b';
+                    is_good = false;
+                    break;
+            }
+            if (!it->selectable_start)
+            {
+                fill = "fill:rgb(255,192,0);";
+                color = 'o'; // orange
+            }
+            if (it->blocked())
+            {
+                fill = "fill:rgb(128,128,128);";
+                color = '-';
+                is_good = false;
+            }
+
+            fill += "fill-opacity:0.7;";
+
+            m_mapper.map(it->point, fill, 4);
+
+            if ((label_good_turns && is_good) || (label_wrong_turns && ! is_good))
+            {
+                std::ostringstream out;
+                out << it->turn_index
+                    << " " << it->operations[0].piece_index << "/" << it->operations[1].piece_index
+                    << " " << si(it->operations[0].seg_id) << "/" << si(it->operations[1].seg_id)
+
+    //              If you want to see travel information
+                    << std::endl
+                    << " nxt " << it->operations[0].enriched.travels_to_ip_index
+                    << "/" << it->operations[1].enriched.travels_to_ip_index
+                    << " or " << it->operations[0].enriched.next_ip_index
+                    << "/" << it->operations[1].enriched.next_ip_index
+                    //<< " frac " << it->operations[0].fraction
+
+    //                If you want to see robust-point coordinates (e.g. to find duplicates)
+    //                << std::endl
+    //                << " " << bg::get<0>(it->robust_point) << " , " << bg::get<1>(it->robust_point)
+
+                    << std::endl;
+                out << " " << bg::method_char(it->method)
+                    << ":" << bg::operation_char(it->operations[0].operation)
+                    << "/" << bg::operation_char(it->operations[1].operation);
+                out << " "
+                    << (it->count_on_offsetted > 0 ? "b" : "") // b: offsetted border
+#if ! defined(BOOST_GEOMETRY_BUFFER_USE_SIDE_OF_INTERSECTION)
+                    << (it->count_within_near_offsetted > 0 ? "n" : "")
+#endif
+                    << (it->count_within > 0 ? "w" : "")
+                    << (it->count_on_helper > 0 ? "h" : "")
+                    << (it->count_on_multi > 0 ? "m" : "")
+                    ;
+
+                offsets[it->get_robust_point()] += 10;
+                int offset = offsets[it->get_robust_point()];
+
+                m_mapper.text(it->point, out.str(), "fill:rgb(0,0,0);font-family='Arial';font-size:9px;", 5, offset);
+
+                offsets[it->get_robust_point()] += 25;
+            }
+        }
+    }
+
+    template <typename Pieces, typename OffsettedRings>
+    inline void map_pieces(Pieces const& pieces,
+                OffsettedRings const& offsetted_rings,
+                bool do_pieces, bool do_indices)
+    {
+        typedef typename boost::range_value<Pieces const>::type piece_type;
+        typedef typename boost::range_value<OffsettedRings const>::type ring_type;
+
+        for(typename boost::range_iterator<Pieces const>::type it = boost::begin(pieces);
+            it != boost::end(pieces);
+            ++it)
+        {
+            const piece_type& piece = *it;
+            bg::segment_identifier seg_id = piece.first_seg_id;
+            if (seg_id.segment_index < 0)
+            {
+                continue;
+            }
+
+            ring_type corner;
+
+
+            ring_type const& ring = offsetted_rings[seg_id.multi_index];
+
+            std::copy(boost::begin(ring) + seg_id.segment_index,
+                    boost::begin(ring) + piece.last_segment_index,
+                    std::back_inserter(corner));
+            std::copy(boost::begin(piece.helper_points),
+                    boost::end(piece.helper_points),
+                    std::back_inserter(corner));
+
+            if (corner.empty())
+            {
+                continue;
+            }
+            if (m_zoom && bg::disjoint(corner, m_alternate_box))
+            {
+                continue;
+            }
+
+            if (m_zoom && do_pieces)
+            {
+                try
+                {
+                    std::string style = "opacity:0.3;stroke:rgb(0,0,0);stroke-width:1;";
+                    typedef typename bg::point_type<Box>::type point_type;
+                    bg::model::multi_polygon<bg::model::polygon<point_type> > clipped;
+                    bg::intersection(ring, m_alternate_box, clipped);
+                    m_mapper.map(clipped,
+                        piece.type == bg::strategy::buffer::buffered_segment
+                        ? style + "fill:rgb(255,128,0);"
+                        : style + "fill:rgb(255,0,0);");
+                }
+                catch (...)
+                {
+                    std::cout << "Error for piece " << piece.index << std::endl;
+                }
+            }
+            else if (do_pieces)
+            {
+                std::string style = "opacity:0.3;stroke:rgb(0,0,0);stroke-width:1;";
+                m_mapper.map(corner,
+                    piece.type == bg::strategy::buffer::buffered_segment
+                    ? style + "fill:rgb(255,128,0);"
+                    : style + "fill:rgb(255,0,0);");
+            }
+
+            if (do_indices)
+            {
+                // Label starting piece_index / segment_index
+                typedef typename bg::point_type<ring_type>::type point_type;
+
+                std::ostringstream out;
+                out << piece.index << " (" << piece_type_char(piece.type) << ") " << piece.first_seg_id.segment_index << ".." << piece.last_segment_index - 1;
+                point_type label_point = bg::return_centroid<point_type>(corner);
+
+                if ((piece.type == bg::strategy::buffer::buffered_concave
+                     || piece.type == bg::strategy::buffer::buffered_flat_end)
+                    && corner.size() >= 2u)
+                {
+                    bg::set<0>(label_point, (bg::get<0>(corner[0]) + bg::get<0>(corner[1])) / 2.0);
+                    bg::set<1>(label_point, (bg::get<1>(corner[0]) + bg::get<1>(corner[1])) / 2.0);
+                }
+                m_mapper.text(label_point, out.str(), "fill:rgb(255,0,0);font-family='Arial';font-size:10px;", 5, 5);
+            }
+        }
+    }
+
+    template <typename TraversedRings>
+    inline void map_traversed_rings(TraversedRings const& traversed_rings)
+    {
+        for(typename boost::range_iterator<TraversedRings const>::type it
+                = boost::begin(traversed_rings); it != boost::end(traversed_rings); ++it)
+        {
+            m_mapper.map(*it, "opacity:0.4;fill:none;stroke:rgb(0,255,0);stroke-width:2");
+        }
+    }
+
+    template <typename OffsettedRings>
+    inline void map_offsetted_rings(OffsettedRings const& offsetted_rings)
+    {
+        for(typename boost::range_iterator<OffsettedRings const>::type it
+                = boost::begin(offsetted_rings); it != boost::end(offsetted_rings); ++it)
+        {
+            if (it->discarded())
+            {
+                m_mapper.map(*it, "opacity:0.4;fill:none;stroke:rgb(255,0,0);stroke-width:2");
+            }
+            else
+            {
+                m_mapper.map(*it, "opacity:0.4;fill:none;stroke:rgb(0,0,255);stroke-width:2");
+            }
+        }
+    }
+
+
+    SvgMapper& m_mapper;
+    Box m_alternate_box;
+    bool m_zoom;
+
+};
+
+template <typename Point>
+class buffer_svg_mapper
+{
+public :
+
+    buffer_svg_mapper(std::string const& casename)
+        : m_casename(casename)
+        , m_zoom(false)
+    {}
+
+    template <typename Mapper, typename Visitor, typename Envelope>
+    void prepare(Mapper& mapper, Visitor& visitor, Envelope const& envelope, double box_buffer_distance)
+    {
+#ifdef BOOST_GEOMETRY_BUFFER_TEST_SVG_USE_ALTERNATE_BOX
+        // Create a zoomed-in view
+        bg::read_wkt(BOOST_GEOMETRY_BUFFER_TEST_SVG_ALTERNATE_BOX, alternate_box);
+        mapper.add(alternate_box);
+
+        // Take care non-visible elements are skipped
+        visitor.m_alternate_box = alternate_box;
+        m_zoom = true;
+#else
+        bg::model::box<Point> box = envelope;
+        bg::buffer(box, box, box_buffer_distance);
+        mapper.add(box);
+#endif
+
+        boost::ignore_unused(visitor);
+    }
+
+    void set_alternate_box(bg::model::box<Point> const& box)
+    {
+        m_alternate_box = box;
+        m_zoom = true;
+    }
+
+    template <typename Mapper, typename Geometry, typename GeometryBuffer>
+    void map_input_output(Mapper& mapper, Geometry const& geometry,
+            GeometryBuffer const& buffered, bool negative)
+    {
+        bool const areal = boost::is_same
+            <
+                typename bg::tag_cast
+                    <
+                        typename bg::tag<Geometry>::type,
+                        bg::areal_tag
+                    >::type, bg::areal_tag
+            >::type::value;
+
+        if (m_zoom)
+        {
+            map_io_zoomed(mapper, geometry, buffered, negative, areal);
+        }
+        else
+        {
+            map_io(mapper, geometry, buffered, negative, areal);
+        }
+    }
+
+    template <typename Mapper, typename Geometry, typename RescalePolicy>
+    void map_self_ips(Mapper& mapper, Geometry const& geometry, RescalePolicy const& rescale_policy)
+    {
+        typedef bg::detail::overlay::turn_info
+        <
+            Point,
+            typename bg::segment_ratio_type<Point, RescalePolicy>::type
+        > turn_info;
+
+        std::vector<turn_info> turns;
+
+        bg::detail::self_get_turn_points::no_interrupt_policy policy;
+        bg::self_turns
+            <
+                bg::detail::overlay::assign_null_policy
+            >(geometry, rescale_policy, turns, policy);
+
+        BOOST_FOREACH(turn_info const& turn, turns)
+        {
+            mapper.map(turn.point, "fill:rgb(255,128,0);stroke:rgb(0,0,100);stroke-width:1", 3);
+        }
+    }
+
+private :
+
+    template <typename Mapper, typename Geometry, typename GeometryBuffer>
+    void map_io(Mapper& mapper, Geometry const& geometry,
+            GeometryBuffer const& buffered, bool negative, bool areal)
+    {
+        // Map input geometry in green
+        if (areal)
+        {
+            mapper.map(geometry, "opacity:0.5;fill:rgb(0,128,0);stroke:rgb(0,64,0);stroke-width:2");
+        }
+        else
+        {
+            // TODO: clip input points/linestring
+            mapper.map(geometry, "opacity:0.5;stroke:rgb(0,128,0);stroke-width:10");
+        }
+
+        {
+            // Map buffer in yellow (inflate) and with orange-dots (deflate)
+            std::string style = negative
+                ? "opacity:0.4;fill:rgb(255,255,192);stroke:rgb(255,128,0);stroke-width:3"
+                : "opacity:0.4;fill:rgb(255,255,128);stroke:rgb(0,0,0);stroke-width:3";
+
+            mapper.map(buffered, style);
+        }
+    }
+
+    template <typename Mapper, typename Geometry, typename GeometryBuffer>
+    void map_io_zoomed(Mapper& mapper, Geometry const& geometry,
+            GeometryBuffer const& buffered, bool negative, bool areal)
+    {
+        // Map input geometry in green
+        if (areal)
+        {
+            // Assuming input is areal
+            GeometryBuffer clipped;
+// TODO: the next line does NOT compile for multi-point, TODO: implement that line
+//            bg::intersection(geometry, m_alternate_box, clipped);
+            mapper.map(clipped, "opacity:0.5;fill:rgb(0,128,0);stroke:rgb(0,64,0);stroke-width:2");
+        }
+        else
+        {
+            // TODO: clip input (multi)point/linestring
+            mapper.map(geometry, "opacity:0.5;stroke:rgb(0,128,0);stroke-width:10");
+        }
+
+        {
+            // Map buffer in yellow (inflate) and with orange-dots (deflate)
+            std::string style = negative
+                ? "opacity:0.4;fill:rgb(255,255,192);stroke:rgb(255,128,0);stroke-width:3"
+                : "opacity:0.4;fill:rgb(255,255,128);stroke:rgb(0,0,0);stroke-width:3";
+
+            try
+            {
+                // Clip output multi-polygon with box
+                GeometryBuffer clipped;
+                bg::intersection(buffered, m_alternate_box, clipped);
+                mapper.map(clipped, style);
+            }
+            catch (...)
+            {
+                std::cout << "Error for buffered output " << m_casename << std::endl;
+            }
+        }
+    }
+
+    bg::model::box<Point> m_alternate_box;
+    bool m_zoom;
+    std::string m_casename;
+};
+
+
+#endif

test/algorithms/buffer/test_buffer_svg_per_turn.hpp

@@ -0,0 +1,164 @@
+// Boost.Geometry (aka GGL, Generic Geometry Library)
+// Unit Test Helper
+
+// Copyright (c) 2010-2015 Barend Gehrels, Amsterdam, the Netherlands.
+// Use, modification and distribution is subject to the Boost Software License,
+// Version 1.0. (See accompanying file LICENSE_1_0.txt or copy at
+// http://www.boost.org/LICENSE_1_0.txt)
+
+
+#ifndef BOOST_GEOMETRY_TEST_BUFFER_SVG_PER_TURN_HPP
+#define BOOST_GEOMETRY_TEST_BUFFER_SVG_PER_TURN_HPP
+
+#include <fstream>
+#include <vector>
+
+#include <test_buffer_svg.hpp>
+#include <boost/ptr_container/ptr_vector.hpp>
+
+template <typename Point>
+class save_turns_visitor
+{
+public :
+    typedef std::vector<std::pair<Point, int> > vector_type;
+
+    template <typename Turns>
+    inline void get_turns(Turns const& turns)
+    {
+        for (typename boost::range_iterator<Turns const>::type it =
+            boost::begin(turns); it != boost::end(turns); ++it)
+        {
+            m_points.push_back(std::make_pair(it->point, it->turn_index));
+        }
+    }
+
+    template <typename PieceCollection>
+    inline void apply(PieceCollection const& collection, int phase)
+    {
+        if (phase == 0)
+        {
+            get_turns(collection.m_turns);
+        }
+    }
+
+    vector_type const& get_points() { return m_points; }
+
+private :
+    vector_type m_points;
+};
+
+
+
+template <typename Point>
+class mapper_visitor
+{
+public :
+    mapper_visitor(std::string const& complete, int index, Point const& point)
+        : m_filename(get_filename(complete, index))
+        , m_svg(m_filename.c_str())
+        , m_mapper(m_svg, 1000, 800)
+        , m_visitor(m_mapper)
+        , m_buffer_mapper(m_filename)
+    {
+        box_type box;
+        double half_size = 75.0; // specific for multi_point buffer
+        bg::set<bg::min_corner, 0>(box, bg::get<0>(point) - half_size);
+        bg::set<bg::min_corner, 1>(box, bg::get<1>(point) - half_size);
+        bg::set<bg::max_corner, 0>(box, bg::get<0>(point) + half_size);
+        bg::set<bg::max_corner, 1>(box, bg::get<1>(point) + half_size);
+
+        m_mapper.add(box);
+        m_visitor.set_alternate_box(box);
+        m_buffer_mapper.set_alternate_box(box);
+    }
+
+    // It is used in a ptr vector
+    virtual ~mapper_visitor()
+    {
+    }
+
+    template <typename PieceCollection>
+    inline void apply(PieceCollection const& collection, int phase)
+    {
+        m_visitor.apply(collection, phase);
+    }
+
+    template <typename Geometry, typename GeometryBuffer>
+    void map_input_output(Geometry const& geometry,
+            GeometryBuffer const& buffered, bool negative)
+    {
+        m_buffer_mapper.map_input_output(m_mapper, geometry, buffered, negative);
+    }
+
+private :
+    std::string get_filename(std::string const& complete, int index)
+    {
+        std::ostringstream filename;
+        filename << "buffer_per_turn_" << complete << "_" << index << ".svg";
+        return filename.str();
+    }
+
+    typedef bg::svg_mapper<Point> mapper_type;
+    typedef bg::model::box<Point> box_type;
+
+    std::string m_filename;
+    std::ofstream m_svg;
+    mapper_type m_mapper;
+    svg_visitor<mapper_type, box_type> m_visitor;
+    buffer_svg_mapper<Point> m_buffer_mapper;
+};
+
+template <typename Point>
+class per_turn_visitor
+{
+    // Both fstreams and svg mappers are non-copyable,
+    // therefore we need to use dynamic memory
+    typedef boost::ptr_vector<mapper_visitor<Point> > container_type;
+    container_type mappers;
+
+public :
+
+    typedef std::pair<Point, int> pair_type;
+    typedef std::vector<pair_type> vector_type;
+
+    per_turn_visitor(std::string const& complete_caseid,
+            vector_type const& points)
+    {
+        namespace bg = boost::geometry;
+
+        if (points.size() > 100u)
+        {
+            // Defensive check. Too much intersections. Don't create anything
+            return;
+        }
+
+        BOOST_FOREACH(pair_type const& p, points)
+        {
+            mappers.push_back(new mapper_visitor<Point>(complete_caseid, p.second, p.first));
+        }
+    }
+
+    template <typename PieceCollection>
+    inline void apply(PieceCollection const& collection, int phase)
+    {
+        for(typename container_type::iterator it = mappers.begin();
+            it != mappers.end(); ++it)
+        {
+            it->apply(collection, phase);
+        }
+    }
+
+    template <typename Geometry, typename GeometryBuffer>
+    void map_input_output(Geometry const& geometry,
+            GeometryBuffer const& buffered, bool negative)
+    {
+        for(typename container_type::iterator it = mappers.begin();
+            it != mappers.end(); ++it)
+        {
+           it->map_input_output(geometry, buffered, negative);
+        }
+    }
+};
+
+
+#endif // BOOST_GEOMETRY_TEST_BUFFER_SVG_PER_TURN_HPP

test/algorithms/distance/distance_se_pl_l.cpp

@@ -122,19 +122,19 @@ void test_distance_point_segment(Strategy const& strategy)
 #endif
     typedef test_distance_of_geometries<point_type, segment_type> tester;
 
+    double const d2r = bg::math::d2r<double>();
+
     tester::apply("p-s-01",
                   "POINT(0 0)",
                   "SEGMENT(2 0,3 0)",
-                  2.0 * bg::math::d2r * strategy.radius(),
-                  to_comparable(strategy,
-                                2.0 * bg::math::d2r * strategy.radius()),
+                  2.0 * d2r * strategy.radius(),
+                  to_comparable(strategy, 2.0 * d2r * strategy.radius()),
                   strategy);
     tester::apply("p-s-02",
                   "POINT(2.5 3)",
                   "SEGMENT(2 0,3 0)",
-                  3.0 * bg::math::d2r * strategy.radius(),
-                  to_comparable(strategy,
-                                3.0 * bg::math::d2r * strategy.radius()),
+                  3.0 * d2r * strategy.radius(),
+                  to_comparable(strategy, 3.0 * d2r * strategy.radius()),
                   strategy);
     tester::apply("p-s-03",
                   "POINT(2 0)",
@@ -163,44 +163,38 @@ void test_distance_point_segment(Strategy const& strategy)
     tester::apply("p-s-07",
                   "POINT(90 1e-3)",
                   "SEGMENT(0.5 0,175.5 0)",
-                  1e-3 * bg::math::d2r * strategy.radius(),
-                  to_comparable(strategy,
-                                1e-3 * bg::math::d2r * strategy.radius()),
+                  1e-3 * d2r * strategy.radius(),
+                  to_comparable(strategy, 1e-3 * d2r * strategy.radius()),
                   strategy);
     tester::apply("p-s-08",
                   "POINT(90 1e-4)",
                   "SEGMENT(0.5 0,175.5 0)",
-                  1e-4 * bg::math::d2r * strategy.radius(),
-                  to_comparable(strategy,
-                                1e-4 * bg::math::d2r * strategy.radius()),
+                  1e-4 * d2r * strategy.radius(),
+                  to_comparable(strategy, 1e-4 * d2r * strategy.radius()),
                   strategy);
     tester::apply("p-s-09",
                   "POINT(90 1e-5)",
                   "SEGMENT(0.5 0,175.5 0)",
-                  1e-5 * bg::math::d2r * strategy.radius(),
-                  to_comparable(strategy,
-                                1e-5 * bg::math::d2r * strategy.radius()),
+                  1e-5 * d2r * strategy.radius(),
+                  to_comparable(strategy, 1e-5 * d2r * strategy.radius()),
                   strategy);
     tester::apply("p-s-10",
                   "POINT(90 1e-6)",
                   "SEGMENT(0.5 0,175.5 0)",
-                  1e-6 * bg::math::d2r * strategy.radius(),
-                  to_comparable(strategy,
-                                1e-6 * bg::math::d2r * strategy.radius()),
+                  1e-6 * d2r * strategy.radius(),
+                  to_comparable(strategy, 1e-6 * d2r * strategy.radius()),
                   strategy);
     tester::apply("p-s-11",
                   "POINT(90 1e-7)",
                   "SEGMENT(0.5 0,175.5 0)",
-                  1e-7 * bg::math::d2r * strategy.radius(),
-                  to_comparable(strategy,
-                                1e-7 * bg::math::d2r * strategy.radius()),
+                  1e-7 * d2r * strategy.radius(),
+                  to_comparable(strategy, 1e-7 * d2r * strategy.radius()),
                   strategy);
     tester::apply("p-s-12",
                   "POINT(90 1e-8)",
                   "SEGMENT(0.5 0,175.5 0)",
-                  1e-8 * bg::math::d2r * strategy.radius(),
-                  to_comparable(strategy,
-                                1e-8 * bg::math::d2r * strategy.radius()),
+                  1e-8 * d2r * strategy.radius(),
+                  to_comparable(strategy, 1e-8 * d2r * strategy.radius()),
                   strategy);
 }
 
@@ -215,26 +209,25 @@ void test_distance_point_linestring(Strategy const& strategy)
 #endif
     typedef test_distance_of_geometries<point_type, linestring_type> tester;
 
+    double const d2r = bg::math::d2r<double>();
+
     tester::apply("p-l-01",
                   "POINT(0 0)",
                   "LINESTRING(2 0,2 0)",
-                  2.0 * bg::math::d2r * strategy.radius(),
-                  to_comparable(strategy,
-                                2.0 * bg::math::d2r * strategy.radius()),
+                  2.0 * d2r * strategy.radius(),
+                  to_comparable(strategy, 2.0 * d2r * strategy.radius()),
                   strategy);
     tester::apply("p-l-02",
                   "POINT(0 0)",
                   "LINESTRING(2 0,3 0)",
-                  2.0 * bg::math::d2r * strategy.radius(),
-                  to_comparable(strategy,
-                                2.0 * bg::math::d2r * strategy.radius()),
+                  2.0 * d2r * strategy.radius(),
+                  to_comparable(strategy, 2.0 * d2r * strategy.radius()),
                   strategy);
     tester::apply("p-l-03",
                   "POINT(2.5 3)",
                   "LINESTRING(2 0,3 0)",
-                  3.0 * bg::math::d2r * strategy.radius(),
-                  to_comparable(strategy,
-                                3.0 * bg::math::d2r * strategy.radius()),
+                  3.0 * d2r * strategy.radius(),
+                  to_comparable(strategy, 3.0 * d2r * strategy.radius()),
                   strategy);
     tester::apply("p-l-04",
                   "POINT(2 0)",
@@ -254,9 +247,8 @@ void test_distance_point_linestring(Strategy const& strategy)
     tester::apply("p-l-07",
                   "POINT(7.5 10)",
                   "LINESTRING(1 0,2 0,3 0,4 0,5 0,6 0,7 0,8 0,9 0)",
-                  10.0 * bg::math::d2r * strategy.radius(),
-                  to_comparable(strategy,
-                                10.0 * bg::math::d2r * strategy.radius()),
+                  10.0 * d2r * strategy.radius(),
+                  to_comparable(strategy, 10.0 * d2r * strategy.radius()),
                   strategy);
     tester::apply("p-l-08",
                   "POINT(7.5 10)",
@@ -282,19 +274,19 @@ void test_distance_point_multilinestring(Strategy const& strategy)
             point_type, multi_linestring_type
         > tester;
 
+    double const d2r = bg::math::d2r<double>();
+
     tester::apply("p-ml-01",
                   "POINT(0 0)",
                   "MULTILINESTRING((-5 0,-3 0),(2 0,3 0))",
-                  2.0 * bg::math::d2r * strategy.radius(),
-                  to_comparable(strategy,
-                                2.0 * bg::math::d2r * strategy.radius()),
+                  2.0 * d2r * strategy.radius(),
+                  to_comparable(strategy, 2.0 * d2r * strategy.radius()),
                   strategy);
     tester::apply("p-ml-02",
                   "POINT(2.5 3)",
                   "MULTILINESTRING((-5 0,-3 0),(2 0,3 0))",
-                  3.0 * bg::math::d2r * strategy.radius(),
-                  to_comparable(strategy,
-                                3.0 * bg::math::d2r * strategy.radius()),
+                  3.0 * d2r * strategy.radius(),
+                  to_comparable(strategy, 3.0 * d2r * strategy.radius()),
                   strategy);
     tester::apply("p-ml-03",
                   "POINT(2 0)",
@@ -431,6 +423,8 @@ void test_distance_multipoint_segment(Strategy const& strategy)
 #endif
     typedef test_distance_of_geometries<multi_point_type, segment_type> tester;
 
+    double d2r = bg::math::d2r<double>();
+
     tester::apply("mp-s-01",
                   "MULTIPOINT(0 0,1 0,0 1,1 1)",
                   "SEGMENT(2 0,0 2)",
@@ -442,9 +436,8 @@ void test_distance_multipoint_segment(Strategy const& strategy)
     tester::apply("mp-s-02",
                   "MULTIPOINT(0 0,1 0,0 1,1 1)",
                   "SEGMENT(0 -3,1 -10)",
-                  3.0 * bg::math::d2r * strategy.radius(),
-                  to_comparable(strategy,
-                                3.0 * bg::math::d2r * strategy.radius()),
+                  3.0 * d2r * strategy.radius(),
+                  to_comparable(strategy, 3.0 * d2r * strategy.radius()),
                   strategy);
     tester::apply("mp-s-03",
                   "MULTIPOINT(0 0,1 0,0 1,1 1)",

test/algorithms/distance/distance_se_pl_pl.cpp

@@ -89,7 +89,7 @@ void test_distance_point_point(Strategy const& strategy,
                   "POINT(180 -10)",
                   (is_comparable_strategy
                    ? 1.0
-                   : (180.0 * bg::math::d2r * strategy.radius())),
+                   : (bg::math::pi<double>() * strategy.radius())),
                   1.0,
                   strategy);
     tester::apply("p-p-04",
@@ -97,7 +97,7 @@ void test_distance_point_point(Strategy const& strategy,
                   "POINT(180 0)",
                   (is_comparable_strategy
                    ? 1.0
-                   : (180.0 * bg::math::d2r * strategy.radius())),
+                   : (bg::math::pi<double>() * strategy.radius())),
                   1.0,
                   strategy);
 }

test/algorithms/is_convex.cpp

@@ -0,0 +1,61 @@
+// Boost.Geometry (aka GGL, Generic Geometry Library)
+// Unit Test
+
+// Copyright (c) 2015 Barend Gehrels, Amsterdam, the Netherlands.
+
+// Use, modification and distribution is subject to the Boost Software License,
+// Version 1.0. (See accompanying file LICENSE_1_0.txt or copy at
+// http://www.boost.org/LICENSE_1_0.txt)
+
+#include <cstddef>
+#include <string>
+
+#include <boost/geometry/algorithms/is_convex.hpp>
+
+#include <geometry_test_common.hpp>
+
+#include <boost/geometry/algorithms/correct.hpp>
+#include <boost/geometry/io/wkt/wkt.hpp>
+#include <boost/geometry/strategies/strategies.hpp>
+#include <boost/geometry/geometries/point_xy.hpp>
+#include <boost/geometry/geometries/ring.hpp>
+
+
+template <typename Geometry>
+void test_one(std::string const& case_id, std::string const& wkt, bool expected)
+{
+    Geometry geometry;
+    bg::read_wkt(wkt, geometry);
+    bg::correct(geometry);
+
+    bool detected = bg::is_convex(geometry);
+    BOOST_CHECK_MESSAGE(detected == expected,
+        "Not as expected, case: " << case_id
+            << " / expected: " << expected
+            << " / detected: " << detected);
+}
+
+
+template <typename P>
+void test_all()
+{
+    // rectangular, with concavity
+    std::string const concave1 = "polygon((1 1, 1 4, 3 4, 3 3, 4 3, 4 4, 5 4, 5 1, 1 1))";
+    std::string const triangle = "polygon((1 1, 1 4, 5 1, 1 1))";
+
+    test_one<bg::model::ring<P> >("triangle", triangle, true);
+    test_one<bg::model::ring<P> >("concave1", concave1, false);
+    test_one<bg::model::ring<P, false, false> >("triangle", triangle, true);
+    test_one<bg::model::ring<P, false, false> >("concave1", concave1, false);
+
+    test_one<bg::model::box<P> >("box", "box(0 0,2 2)", true);
+}
+
+
+int test_main(int, char* [])
+{
+    test_all<bg::model::d2::point_xy<int> >();
+    test_all<bg::model::d2::point_xy<double> >();
+
+    return 0;
+}

test/algorithms/relational_operations/Jamfile.v2

@@ -1,11 +1,11 @@
 # Boost.Geometry (aka GGL, Generic Geometry Library)
 #
-# Copyright (c) 2007-2014 Barend Gehrels, Amsterdam, the Netherlands.
-# Copyright (c) 2008-2014 Bruno Lalande, Paris, France.
-# Copyright (c) 2009-2014 Mateusz Loskot, London, UK.
+# Copyright (c) 2007-2015 Barend Gehrels, Amsterdam, the Netherlands.
+# Copyright (c) 2008-2015 Bruno Lalande, Paris, France.
+# Copyright (c) 2009-2015 Mateusz Loskot, London, UK.
 #
-# This file was modified by Oracle on 2014.
-# Modifications copyright (c) 2014, Oracle and/or its affiliates.
+# This file was modified by Oracle on 2014, 2015.
+# Modifications copyright (c) 2014-2015, Oracle and/or its affiliates.
 #
 # Contributed and/or modified by Menelaos Karavelas, on behalf of Oracle
 # Contributed and/or modified by Adam Wulkiewicz, on behalf of Oracle
@@ -19,6 +19,7 @@ test-suite boost-geometry-algorithms-relational
     [ run covered_by.cpp ]
     [ run crosses.cpp ]
     [ run equals.cpp ]
+    [ run equals_on_spheroid.cpp ]
     [ run intersects.cpp ]
     [ run multi_covered_by.cpp ]
     [ run multi_equals.cpp ]

test/algorithms/relational_operations/equals_on_spheroid.cpp

@@ -0,0 +1,183 @@
+// Boost.Geometry (aka GGL, Generic Geometry Library)
+// Unit test
+
+// Copyright (c) 2015, Oracle and/or its affiliates.
+
+// Licensed under the Boost Software License version 1.0.
+// http://www.boost.org/users/license.html
+
+// Contributed and/or modified by Menelaos Karavelas, on behalf of Oracle
+
+#ifndef BOOST_TEST_MODULE
+#define BOOST_TEST_MODULE test_equals_on_spheroid
+#endif
+
+#include <iostream>
+
+#include <boost/test/included/unit_test.hpp>
+
+#include "test_equals.hpp"
+
+#include <boost/geometry/geometries/geometries.hpp>
+
+#include <boost/geometry/core/cs.hpp>
+
+namespace bgm = bg::model;
+
+template <typename P1, typename P2 = P1>
+struct test_point_point
+{
+    static inline void apply(std::string const& header)
+    {
+        std::string const str = header + "-";
+
+        test_geometry<P1, P2>(str + "pp_01", "POINT(0 0)", "POINT(0 0)", true);
+        test_geometry<P1, P2>(str + "pp_02", "POINT(0 0)", "POINT(10 0)", false);
+
+        // points whose longitudes differ by 360 degrees
+        test_geometry<P1, P2>(str + "pp_03", "POINT(0 0)", "POINT(360 0)", true);
+        test_geometry<P1, P2>(str + "pp_04", "POINT(10 0)", "POINT(370 0)", true);
+        test_geometry<P1, P2>(str + "pp_05", "POINT(10 0)", "POINT(-350 0)", true);
+        test_geometry<P1, P2>(str + "pp_06", "POINT(180 10)", "POINT(-180 10)", true);
+        test_geometry<P1, P2>(str + "pp_06a", "POINT(540 10)", "POINT(-540 10)", true);
+        test_geometry<P1, P2>(str + "pp_06b", "POINT(540 370)", "POINT(-540 -350)", true);
+        test_geometry<P1, P2>(str + "pp_06c", "POINT(1260 370)", "POINT(-1260 -350)", true);
+        test_geometry<P1, P2>(str + "pp_06d", "POINT(2340 370)", "POINT(-2340 -350)", true);
+        test_geometry<P1, P2>(str + "pp_06e", "POINT(-180 10)", "POINT(-540 10)", true);
+        test_geometry<P1, P2>(str + "pp_06f", "POINT(180 10)", "POINT(-540 10)", true);
+
+        // north & south pole
+        test_geometry<P1, P2>(str + "pp_07", "POINT(0 90)", "POINT(0 90)", true);
+        test_geometry<P1, P2>(str + "pp_07a", "POINT(0 450)", "POINT(10 -270)", true);
+        test_geometry<P1, P2>(str + "pp_07b", "POINT(0 270)", "POINT(10 90)", false);
+        test_geometry<P1, P2>(str + "pp_07c", "POINT(0 -450)", "POINT(10 90)", false);
+        test_geometry<P1, P2>(str + "pp_08", "POINT(0 90)", "POINT(10 90)", true);
+        test_geometry<P1, P2>(str + "pp_09", "POINT(0 90)", "POINT(0 -90)", false);
+        test_geometry<P1, P2>(str + "pp_10", "POINT(0 -90)", "POINT(0 -90)", true);
+        test_geometry<P1, P2>(str + "pp_11", "POINT(0 -90)", "POINT(10 -90)", true);
+        test_geometry<P1, P2>(str + "pp_11a", "POINT(0 -90)", "POINT(10 90)", false);
+        test_geometry<P1, P2>(str + "pp_12", "POINT(0 -90)", "POINT(0 -85)", false);
+        test_geometry<P1, P2>(str + "pp_13", "POINT(0 90)", "POINT(0 85)", false);
+        test_geometry<P1, P2>(str + "pp_14", "POINT(0 90)", "POINT(10 85)", false);
+
+        // symmetric wrt prime meridian
+        test_geometry<P1, P2>(str + "pp_15", "POINT(-10 45)", "POINT(10 45)", false);
+        test_geometry<P1, P2>(str + "pp_16", "POINT(-170 45)", "POINT(170 45)", false);
+
+        // other points
+        test_geometry<P1, P2>(str + "pp_17", "POINT(-10 45)", "POINT(10 -45)", false);
+        test_geometry<P1, P2>(str + "pp_18", "POINT(-10 -45)", "POINT(10 45)", false);
+        test_geometry<P1, P2>(str + "pp_19", "POINT(10 -135)", "POINT(10 45)", false);
+        test_geometry<P1, P2>(str + "pp_20", "POINT(190 135)", "POINT(10 45)", true);
+        test_geometry<P1, P2>(str + "pp_21", "POINT(190 150)", "POINT(10 30)", true);
+        test_geometry<P1, P2>(str + "pp_21a", "POINT(-170 150)", "POINT(10 30)", true);
+        test_geometry<P1, P2>(str + "pp_22", "POINT(190 -135)", "POINT(10 -45)", true);
+        test_geometry<P1, P2>(str + "pp_23", "POINT(190 -150)", "POINT(10 -30)", true);
+        test_geometry<P1, P2>(str + "pp_23a", "POINT(-170 -150)", "POINT(10 -30)", true);
+    }
+};
+
+
+template <typename P>
+void test_segment_segment(std::string const& header)
+{
+    typedef bgm::segment<P> seg;
+
+    std::string const str = header + "-";
+
+    test_geometry<seg, seg>(str + "ss_01",
+                            "SEGMENT(10 0,180 0)",
+                            "SEGMENT(10 0,-180 0)",
+                            true);
+    test_geometry<seg, seg>(str + "ss_02",
+                            "SEGMENT(0 90,180 0)",
+                            "SEGMENT(10 90,-180 0)",
+                            true);
+    test_geometry<seg, seg>(str + "ss_03",
+                            "SEGMENT(0 90,0 -90)",
+                            "SEGMENT(10 90,20 -90)",
+                            true);
+    test_geometry<seg, seg>(str + "ss_04",
+                            "SEGMENT(10 80,10 -80)",
+                            "SEGMENT(10 80,20 -80)",
+                            false);
+    test_geometry<seg, seg>(str + "ss_05",
+                            "SEGMENT(170 10,-170 10)",
+                            "SEGMENT(170 10,350 10)",
+                            false);
+}
+
+
+BOOST_AUTO_TEST_CASE( equals_point_point_se )
+{
+    typedef bg::cs::spherical_equatorial<bg::degree> cs_type;
+
+    test_point_point<bgm::point<int, 2, cs_type> >::apply("se");
+    test_point_point<bgm::point<double, 2, cs_type> >::apply("se");
+    test_point_point<bgm::point<long double, 2, cs_type> >::apply("se");
+
+    // mixed point types
+    test_point_point
+        <
+            bgm::point<double, 2, cs_type>, bgm::point<int, 2, cs_type>
+        >::apply("se");
+
+    test_point_point
+        <
+            bgm::point<double, 2, cs_type>, bgm::point<long double, 2, cs_type>
+        >::apply("se");
+
+#if defined(HAVE_TTMATH)
+    test_point_point<bgm::point<ttmath_big, 2, cs_type> >::apply("se");
+#endif
+}
+
+BOOST_AUTO_TEST_CASE( equals_point_point_geo )
+{
+    typedef bg::cs::geographic<bg::degree> cs_type;
+
+    test_point_point<bgm::point<int, 2, cs_type> >::apply("geo");
+    test_point_point<bgm::point<double, 2, cs_type> >::apply("geo");
+    test_point_point<bgm::point<long double, 2, cs_type> >::apply("geo");
+
+    // mixed point types
+    test_point_point
+        <
+            bgm::point<double, 2, cs_type>, bgm::point<int, 2, cs_type>
+        >::apply("se");
+
+    test_point_point
+        <
+            bgm::point<double, 2, cs_type>, bgm::point<long double, 2, cs_type>
+        >::apply("se");
+
+#if defined(HAVE_TTMATH)
+    test_point_point<bgm::point<ttmath_big, 2, cs_type> >::apply("geo");
+#endif
+}
+
+BOOST_AUTO_TEST_CASE( equals_segment_segment_se )
+{
+    typedef bg::cs::spherical_equatorial<bg::degree> cs_type;
+
+    test_segment_segment<bgm::point<int, 2, cs_type> >("se");
+    test_segment_segment<bgm::point<double, 2, cs_type> >("se");
+    test_segment_segment<bgm::point<long double, 2, cs_type> >("se");
+
+#if defined(HAVE_TTMATH)
+    test_segment_segment<bgm::point<ttmath_big, 2, cs_type> >("se");
+#endif
+}
+
+BOOST_AUTO_TEST_CASE( equals_segment_segment_geo )
+{
+    typedef bg::cs::geographic<bg::degree> cs_type;
+
+    test_segment_segment<bgm::point<int, 2, cs_type> >("geo");
+    test_segment_segment<bgm::point<double, 2, cs_type> >("geo");
+    test_segment_segment<bgm::point<long double, 2, cs_type> >("geo");
+
+#if defined(HAVE_TTMATH)
+    test_segment_segment<bgm::point<ttmath_big, 2, cs_type> >("geo");
+#endif
+}

test/algorithms/set_operations/intersection/intersection.cpp

@@ -434,6 +434,38 @@ void test_areal_linear()
 
 }
 
+
+template <typename Linestring, typename Box>
+void test_linear_box()
+{
+    typedef bg::model::multi_linestring<Linestring> multi_linestring_type;
+
+    test_one_lp<Linestring, Box, Linestring>
+        ("case-l-b-01",
+         "BOX(-10 -10,10 10)",
+         "LINESTRING(-20 -20, 0 0,20 20)",
+         1, 3, 20 * sqrt(2.0));
+
+    test_one_lp<Linestring, Box, Linestring>
+        ("case-l-b-02",
+         "BOX(-10 -10,10 10)",
+         "LINESTRING(-20 -20, 20 20)",
+         1, 2, 20.0 * sqrt(2.0));
+
+    test_one_lp<Linestring, Box, Linestring>
+        ("case-l-b-02",
+         "BOX(-10 -10,10 10)",
+         "LINESTRING(-20 -20, 20 20,15 0,0 -15)",
+         2, 4, 25.0 * sqrt(2.0));
+
+    test_one_lp<Linestring, Box, multi_linestring_type>
+        ("case-ml-b-01",
+         "BOX(-10 -10,10 10)",
+         "MULTILINESTRING((-20 -20, 20 20),(0 -15,15 0))",
+         2, 4, 25.0 * sqrt(2.0));
+}
+
+
 template <typename P>
 void test_all()
 {
@@ -456,6 +488,8 @@ void test_all()
     test_areal_linear<polygon_ccw_open, linestring>();
 #endif
 
+    test_linear_box<linestring, box>();
+
     // Test polygons clockwise and counter clockwise
     test_areal<polygon>();
 

test/algorithms/transform.cpp

@@ -130,24 +130,34 @@ int test_main(int, char* [])
     test_all<bg::model::d2::point_xy<int>, bg::model::d2::point_xy<float> >(1.0);
 
     test_all<bg::model::point<double, 2, bg::cs::spherical<bg::degree> >,
-        bg::model::point<double, 2, bg::cs::spherical<bg::radian> > >(bg::math::d2r);
+        bg::model::point<double, 2, bg::cs::spherical<bg::radian> > >(bg::math::d2r<double>());
     test_all<bg::model::point<double, 2, bg::cs::spherical<bg::radian> >,
-        bg::model::point<double, 2, bg::cs::spherical<bg::degree> > >(bg::math::r2d);
+        bg::model::point<double, 2, bg::cs::spherical<bg::degree> > >(bg::math::r2d<double>());
 
     test_all<bg::model::point<int, 2, bg::cs::spherical<bg::degree> >,
-        bg::model::point<float, 2, bg::cs::spherical<bg::radian> > >(bg::math::d2r);
+        bg::model::point<float, 2, bg::cs::spherical<bg::radian> > >(bg::math::d2r<float>());
 
     test_transformations<float, bg::degree>(4, 52, 1);
     test_transformations<double, bg::degree>(4, 52, 1);
 
-    test_transformations<float, bg::radian>(3 * bg::math::d2r, 51 * bg::math::d2r, 1);
-    test_transformations<double, bg::radian>(3 * bg::math::d2r, 51 * bg::math::d2r, 1);
+    test_transformations
+        <
+            float, bg::radian
+        >(3 * bg::math::d2r<float>(), 51 * bg::math::d2r<float>(), 1);
+
+    test_transformations
+        <
+            double, bg::radian
+        >(3 * bg::math::d2r<double>(), 51 * bg::math::d2r<double>(), 1);
 
 #if defined(HAVE_TTMATH)
     typedef bg::model::d2::point_xy<ttmath_big > PT;
     test_all<PT, PT>();
     test_transformations<ttmath_big, bg::degree>(4, 52, 1);
-    test_transformations<ttmath_big, bg::radian>(3 * bg::math::d2r, 51 * bg::math::d2r, 1);
+    test_transformations
+        <
+            ttmath_big, bg::radian
+        >(3 * bg::math::d2r<ttmath_big>(), 51 * bg::math::d2r<ttmath_big>(), 1);
 #endif
 
     return 0;

test/strategies/side_of_intersection.cpp

@@ -19,9 +19,13 @@ namespace bg = boost::geometry;
 
 int test_main(int, char* [])
 {
-    typedef bg::model::d2::point_xy<int> point;
+    typedef bg::model::d2::point_xy<boost::long_long_type> point;
     typedef bg::model::segment<point> segment;
 
+    typedef bg::strategy::side::side_of_intersection side;
+
+    point no_fb(-99, -99);
+
     segment a(point(20, 10), point(10, 20));
 
     segment b1(point(11, 16), point(20, 14));  // IP with a: (14.857, 15.143)
@@ -31,27 +35,70 @@ int test_main(int, char* [])
     segment c2(point(15, 16), point(14, 8));
     segment c3(point(15, 16), point(15, 8));
 
-    typedef bg::strategy::side::side_of_intersection side;
 
-    BOOST_CHECK_EQUAL( 1, side::apply(a, b1, c1));
-    BOOST_CHECK_EQUAL(-1, side::apply(a, b1, c2));
-    BOOST_CHECK_EQUAL(-1, side::apply(a, b1, c3));
+    BOOST_CHECK_EQUAL( 1, side::apply(a, b1, c1, no_fb));
+    BOOST_CHECK_EQUAL(-1, side::apply(a, b1, c2, no_fb));
+    BOOST_CHECK_EQUAL(-1, side::apply(a, b1, c3, no_fb));
 
-    BOOST_CHECK_EQUAL( 1, side::apply(a, b2, c1));
-    BOOST_CHECK_EQUAL( 1, side::apply(a, b2, c2));
-    BOOST_CHECK_EQUAL( 0, side::apply(a, b2, c3));
+    BOOST_CHECK_EQUAL( 1, side::apply(a, b2, c1, no_fb));
+    BOOST_CHECK_EQUAL( 1, side::apply(a, b2, c2, no_fb));
+    BOOST_CHECK_EQUAL( 0, side::apply(a, b2, c3, no_fb));
 
-    // Check internal calculation-method:
-    BOOST_CHECK_EQUAL(-1400, side::side_value<int>(a, b1, c2));
-    BOOST_CHECK_EQUAL( 2800, side::side_value<int>(a, b1, c1));
+    // Collinear cases
+    // 1: segments intersecting are collinear (with a fallback point):
+    {
+        point fb(5, 5);
+        segment col1(point(0, 5), point(5, 5));
+        segment col2(point(5, 5), point(10, 5)); // One IP with col1 at (5,5)
+        segment col3(point(5, 0), point(5, 5));
+        BOOST_CHECK_EQUAL( 0, side::apply(col1, col2, col3, fb));
+    }
+    // 2: segment of side calculation collinear with one of the segments
+    {
+        point fb(5, 5);
+        segment col1(point(0, 5), point(10, 5));
+        segment col2(point(5, 5), point(5, 12)); // IP with col1 at (5,5)
+        segment col3(point(5, 1), point(5, 5)); // collinear with col2
+        BOOST_CHECK_EQUAL( 0, side::apply(col1, col2, col3, fb));
+    }
+    {
+        point fb(5, 5);
+        segment col1(point(10, 5), point(0, 5));
+        segment col2(point(5, 5), point(5, 12)); // IP with col1 at (5,5)
+        segment col3(point(5, 1), point(5, 5)); // collinear with col2
+        BOOST_CHECK_EQUAL( 0, side::apply(col1, col2, col3, fb));
+    }
+    {
+        point fb(5, 5);
+        segment col1(point(0, 5), point(10, 5));
+        segment col2(point(5, 12), point(5, 5)); // IP with col1 at (5,5)
+        segment col3(point(5, 1), point(5, 5)); // collinear with col2
+        BOOST_CHECK_EQUAL( 0, side::apply(col1, col2, col3, fb));
+    }
 
-    BOOST_CHECK_EQUAL (2800, side::side_value<int>(a, b1, c1));
-    BOOST_CHECK_EQUAL(-1400, side::side_value<int>(a, b1, c2));
-    BOOST_CHECK_EQUAL(-5600, side::side_value<int>(a, b1, c3));
+    {
+        point fb(517248, -517236);
+        segment col1(point(-172408, -517236), point(862076, -517236));
+        segment col2(point(517248, -862064), point(517248, -172408));
+        segment col3(point(517248, -172408), point(517248, -517236));
+        BOOST_CHECK_EQUAL( 0, side::apply(col1, col2, col3, fb));
+    }
+    {
+        point fb(-221647, -830336);
+        segment col1(point(-153817, -837972), point(-222457, -830244));
+        segment col2(point(-221139, -833615), point(-290654, -384388));
+        segment col3(point(-255266, -814663), point(-264389, -811197));
+        BOOST_CHECK_EQUAL(1, side::apply(col1, col2, col3, fb)); // Left of segment...
+    }
 
-    BOOST_CHECK_EQUAL(12800, side::side_value<int>(a, b2, c1));
-    BOOST_CHECK_EQUAL( 6400, side::side_value<int>(a, b2, c2));
-    BOOST_CHECK_EQUAL(    0, side::side_value<int>(a, b2, c3));
+
+    {
+        point fb(27671131, 30809240);
+        segment col1(point(27671116, 30809247), point(27675474, 30807351));
+        segment col2(point(27666779, 30811130), point(27671139, 30809237));
+        segment col3(point(27671122, 30809244), point(27675480, 30807348));
+        BOOST_CHECK_EQUAL(1, side::apply(col1, col2, col3, fb)); // Left of segment...
+    }
 
     // TODO: we might add a check calculating the IP, determining the side
     // with the normal side strategy, and verify the results are equal

test/strategies/vincenty.cpp

@@ -121,24 +121,27 @@ void test_vincenty(double lon1, double lat1, double lon2, double lat2,
 
     // formula
     {
-        bg::detail::vincenty_inverse<calc_t> vi(lon1 * bg::math::d2r,
-                                                lat1 * bg::math::d2r,
-                                                lon2 * bg::math::d2r,
-                                                lat2 * bg::math::d2r,
+        double const d2r = bg::math::d2r<double>();
+        double const r2d = bg::math::r2d<double>();
+
+        bg::detail::vincenty_inverse<calc_t> vi(lon1 * d2r,
+                                                lat1 * d2r,
+                                                lon2 * d2r,
+                                                lat2 * d2r,
                                                 spheroid);
         calc_t dist = vi.distance();
         calc_t az12 = vi.azimuth12();
         calc_t az21 = vi.azimuth21();
 
-        calc_t az12_deg = az12 * bg::math::r2d;
-        calc_t az21_deg = az21 * bg::math::r2d;
+        calc_t az12_deg = az12 * r2d;
+        calc_t az21_deg = az21 * r2d;
         
         BOOST_CHECK_CLOSE(dist, calc_t(expected_distance), tolerance);
         check_deg("az12_deg", az12_deg, calc_t(expected_azimuth_12), tolerance, error);
         check_deg("az21_deg", az21_deg, calc_t(expected_azimuth_21), tolerance, error);
 
-        bg::detail::vincenty_direct<calc_t> vd(lon1 * bg::math::d2r,
-                                               lat1 * bg::math::d2r,
+        bg::detail::vincenty_direct<calc_t> vd(lon1 * d2r,
+                                               lat1 * d2r,
                                                dist,
                                                az12,
                                                spheroid);
@@ -146,9 +149,9 @@ void test_vincenty(double lon1, double lat1, double lon2, double lat2,
         calc_t direct_lat2 = vd.lat2();
         calc_t direct_az21 = vd.azimuth21();
 
-        calc_t direct_lon2_deg = direct_lon2 * bg::math::r2d;
-        calc_t direct_lat2_deg = direct_lat2 * bg::math::r2d;
-        calc_t direct_az21_deg = direct_az21 * bg::math::r2d;
+        calc_t direct_lon2_deg = direct_lon2 * r2d;
+        calc_t direct_lat2_deg = direct_lat2 * r2d;
+        calc_t direct_az21_deg = direct_az21 * r2d;
         
         check_deg("direct_lon2_deg", direct_lon2_deg, calc_t(lon2), tolerance, error);
         check_deg("direct_lat2_deg", direct_lat2_deg, calc_t(lat2), tolerance, error);