Boost logo

Boost-Commit :

Subject: [Boost-commit] svn:boost r77296 - in trunk: boost/geometry/algorithms/detail/overlay boost/geometry/policies/relate boost/geometry/strategies boost/geometry/strategies/cartesian boost/geometry/util libs/geometry/test/algorithms libs/geometry/test/algorithms/overlay
From: barend.gehrels_at_[hidden]
Date: 2012-03-10 14:11:19


Author: barendgehrels
Date: 2012-03-10 14:11:17 EST (Sat, 10 Mar 2012)
New Revision: 77296
URL: http://svn.boost.org/trac/boost/changeset/77296

Log:
[geometry] Fixed several robustness issues: non-valid polygons/rings are not added anymore;
  collinear is now symmetric (it could happen that A was collinear w.r.t. B but not vice versa, that is now resolved);
  vertical/horizontal (nearly collinear) segments are now checked later (this and previous bug were found by buffer-high-volume tests).

Added testcases (buffer_rt_f did cause problems with previous implementation, fixed now). Updated testcases (some cases are quite sensitive to implementation while output is still valid)

Text files modified:
   trunk/boost/geometry/algorithms/detail/overlay/add_rings.hpp | 27 +++++
   trunk/boost/geometry/algorithms/detail/overlay/assign_parents.hpp | 7
   trunk/boost/geometry/policies/relate/direction.hpp | 2
   trunk/boost/geometry/policies/relate/intersection_points.hpp | 41 +------
   trunk/boost/geometry/policies/relate/tupled.hpp | 6
   trunk/boost/geometry/strategies/cartesian/cart_intersect.hpp | 206 +++++++++++++++++++--------------------
   trunk/boost/geometry/strategies/side_info.hpp | 21 ++-
   trunk/boost/geometry/util/math.hpp | 44 ++++++++
   trunk/libs/geometry/test/algorithms/difference.cpp | 31 ++++-
   trunk/libs/geometry/test/algorithms/intersection.cpp | 33 ++++-
   trunk/libs/geometry/test/algorithms/overlay/overlay_cases.hpp | 20 +++
   trunk/libs/geometry/test/algorithms/union.cpp | 49 ++++++---
   12 files changed, 296 insertions(+), 191 deletions(-)

Modified: trunk/boost/geometry/algorithms/detail/overlay/add_rings.hpp
==============================================================================
--- trunk/boost/geometry/algorithms/detail/overlay/add_rings.hpp (original)
+++ trunk/boost/geometry/algorithms/detail/overlay/add_rings.hpp 2012-03-10 14:11:17 EST (Sat, 10 Mar 2012)
@@ -9,6 +9,7 @@
 #ifndef BOOST_GEOMETRY_ALGORITHMS_DETAIL_OVERLAY_ADD_RINGS_HPP
 #define BOOST_GEOMETRY_ALGORITHMS_DETAIL_OVERLAY_ADD_RINGS_HPP
 
+#include <boost/geometry/core/closure.hpp>
 #include <boost/geometry/algorithms/detail/overlay/convert_ring.hpp>
 #include <boost/geometry/algorithms/detail/overlay/get_ring.hpp>
 
@@ -73,6 +74,21 @@
             OutputIterator out)
 {
     typedef typename SelectionMap::const_iterator iterator;
+ typedef typename SelectionMap::mapped_type property_type;
+ typedef typename property_type::area_type area_type;
+
+ area_type const zero = 0;
+ int const min_num_points = core_detail::closure::minimum_ring_size
+ <
+ geometry::closure
+ <
+ typename boost::range_value
+ <
+ RingCollection const
+ >::type
+ >::value
+ >::value;
+
 
     for (iterator it = boost::begin(map);
         it != boost::end(map);
@@ -99,7 +115,16 @@
                             *child_it, mit->second.reversed, true);
                 }
             }
- *out++ = result;
+
+ // Only add rings if they satisfy minimal requirements.
+ // This cannot be done earlier (during traversal), not
+ // everything is figured out yet (sum of positive/negative rings)
+ // TODO: individual rings can still contain less than 3 points.
+ if (geometry::num_points(result) >= min_num_points
+ && math::larger(geometry::area(result), zero))
+ {
+ *out++ = result;
+ }
         }
     }
     return out;

Modified: trunk/boost/geometry/algorithms/detail/overlay/assign_parents.hpp
==============================================================================
--- trunk/boost/geometry/algorithms/detail/overlay/assign_parents.hpp (original)
+++ trunk/boost/geometry/algorithms/detail/overlay/assign_parents.hpp 2012-03-10 14:11:17 EST (Sat, 10 Mar 2012)
@@ -130,7 +130,7 @@
             return;
         }
 
- if (outer.real_area > 0)
+ if (math::larger(outer.real_area, 0))
         {
             if (inner.real_area < 0 || m_check_for_orientation)
             {
@@ -317,13 +317,14 @@
>
 inline void assign_parents(Geometry const& geometry,
             RingCollection const& collection,
- RingMap& ring_map)
+ RingMap& ring_map,
+ bool check_for_orientation)
 {
     // Call it with an empty geometry
     // (ring_map should be empty for source_id==1)
 
     Geometry empty;
- assign_parents(geometry, empty, collection, ring_map, true);
+ assign_parents(geometry, empty, collection, ring_map, check_for_orientation);
 }
 
 

Modified: trunk/boost/geometry/policies/relate/direction.hpp
==============================================================================
--- trunk/boost/geometry/policies/relate/direction.hpp (original)
+++ trunk/boost/geometry/policies/relate/direction.hpp 2012-03-10 14:11:17 EST (Sat, 10 Mar 2012)
@@ -126,7 +126,9 @@
     typedef typename select_most_precise<coordinate_type, double>::type rtype;
 
 
+ template <typename R>
     static inline return_type segments_intersect(side_info const& sides,
+ R const&,
                     coordinate_type const& dx1, coordinate_type const& dy1,
                     coordinate_type const& dx2, coordinate_type const& dy2,
                     S1 const& s1, S2 const& s2)

Modified: trunk/boost/geometry/policies/relate/intersection_points.hpp
==============================================================================
--- trunk/boost/geometry/policies/relate/intersection_points.hpp (original)
+++ trunk/boost/geometry/policies/relate/intersection_points.hpp 2012-03-10 14:11:17 EST (Sat, 10 Mar 2012)
@@ -36,56 +36,33 @@
     typedef ReturnType return_type;
     typedef S1 segment_type1;
     typedef S2 segment_type2;
+
     typedef typename select_calculation_type
         <
             S1, S2, CalculationType
>::type coordinate_type;
 
+ template <typename R>
     static inline return_type segments_intersect(side_info const&,
+ R const& r,
                     coordinate_type const& dx1, coordinate_type const& dy1,
                     coordinate_type const& dx2, coordinate_type const& dy2,
                     S1 const& s1, S2 const& s2)
     {
- return_type result;
         typedef typename geometry::coordinate_type
             <
                 typename return_type::point_type
- >::type coordinate_type;
-
- // Get the same type, but at least a double (also used for divisions)
- typedef typename select_most_precise
- <
- coordinate_type, double
- >::type promoted_type;
-
- promoted_type const s1x = get<0, 0>(s1);
- promoted_type const s1y = get<0, 1>(s1);
+ >::type return_coordinate_type;
 
- // Calculate other determinants - Cramers rule
- promoted_type const wx = get<0, 0>(s1) - get<0, 0>(s2);
- promoted_type const wy = get<0, 1>(s1) - get<0, 1>(s2);
- promoted_type const d = detail::determinant<promoted_type>(dx1, dy1, dx2, dy2);
- promoted_type const da = detail::determinant<promoted_type>(dx2, dy2, wx, wy);
-
- // r: ratio 0-1 where intersection divides A/B
- promoted_type r = da / d;
- promoted_type const zero = promoted_type();
- promoted_type const one = 1;
- // Handle robustness issues
- if (r < zero)
- {
- r = zero;
- }
- else if (r > one)
- {
- r = one;
- }
+ coordinate_type const s1x = get<0, 0>(s1);
+ coordinate_type const s1y = get<0, 1>(s1);
 
+ return_type result;
         result.count = 1;
         set<0>(result.intersections[0],
- boost::numeric_cast<coordinate_type>(s1x + r * promoted_type(dx1)));
+ boost::numeric_cast<return_coordinate_type>(R(s1x) + r * R(dx1)));
         set<1>(result.intersections[0],
- boost::numeric_cast<coordinate_type>(s1y + r * promoted_type(dy1)));
+ boost::numeric_cast<return_coordinate_type>(R(s1y) + r * R(dy1)));
 
         return result;
     }

Modified: trunk/boost/geometry/policies/relate/tupled.hpp
==============================================================================
--- trunk/boost/geometry/policies/relate/tupled.hpp (original)
+++ trunk/boost/geometry/policies/relate/tupled.hpp 2012-03-10 14:11:17 EST (Sat, 10 Mar 2012)
@@ -49,16 +49,18 @@
     // Get the same type, but at least a double
     typedef typename select_most_precise<coordinate_type, double>::type rtype;
 
+ template <typename R>
     static inline return_type segments_intersect(side_info const& sides,
+ R const& r,
                     coordinate_type const& dx1, coordinate_type const& dy1,
                     coordinate_type const& dx2, coordinate_type const& dy2,
                     segment_type1 const& s1, segment_type2 const& s2)
     {
         return boost::make_tuple
             (
- Policy1::segments_intersect(sides,
+ Policy1::segments_intersect(sides, r,
                     dx1, dy1, dx2, dy2, s1, s2),
- Policy2::segments_intersect(sides,
+ Policy2::segments_intersect(sides, r,
                     dx1, dy1, dx2, dy2, s1, s2)
             );
     }

Modified: trunk/boost/geometry/strategies/cartesian/cart_intersect.hpp
==============================================================================
--- trunk/boost/geometry/strategies/cartesian/cart_intersect.hpp (original)
+++ trunk/boost/geometry/strategies/cartesian/cart_intersect.hpp 2012-03-10 14:11:17 EST (Sat, 10 Mar 2012)
@@ -41,21 +41,17 @@
 namespace detail
 {
 
-template <typename Segment, size_t Dimension>
-struct segment_arrange
+template <std::size_t Dimension, typename Segment, typename T>
+static inline void segment_arrange(Segment const& s, T& s_1, T& s_2, bool& swapped)
 {
- template <typename T>
- static inline void apply(Segment const& s, T& s_1, T& s_2, bool& swapped)
+ s_1 = get<0, Dimension>(s);
+ s_2 = get<1, Dimension>(s);
+ if (s_1 > s_2)
     {
- s_1 = get<0, Dimension>(s);
- s_2 = get<1, Dimension>(s);
- if (s_1 > s_2)
- {
- std::swap(s_1, s_2);
- swapped = true;
- }
+ std::swap(s_1, s_2);
+ swapped = true;
     }
-};
+}
 
 template <std::size_t Index, typename Segment>
 inline typename geometry::point_type<Segment>::type get_from_index(
@@ -121,34 +117,18 @@
             coordinate_type const& dx_a, coordinate_type const& dy_a,
             coordinate_type const& dx_b, coordinate_type const& dy_b)
     {
- // 1) Handle "disjoint", common case.
- // per dimension, 2 cases: a_1----------a_2 b_1-------b_2 or B left of A
- coordinate_type ax_1, ax_2, bx_1, bx_2;
- bool ax_swapped = false, bx_swapped = false;
- detail::segment_arrange<segment_type1, 0>::apply(a, ax_1, ax_2, ax_swapped);
- detail::segment_arrange<segment_type2, 0>::apply(b, bx_1, bx_2, bx_swapped);
- if (ax_2 < bx_1 || ax_1 > bx_2)
- {
- return Policy::disjoint();
- }
-
- // 1b) In Y-dimension
- coordinate_type ay_1, ay_2, by_1, by_2;
- bool ay_swapped = false, by_swapped = false;
- detail::segment_arrange<segment_type1, 1>::apply(a, ay_1, ay_2, ay_swapped);
- detail::segment_arrange<segment_type2, 1>::apply(b, by_1, by_2, by_swapped);
- if (ay_2 < by_1 || ay_1 > by_2)
- {
- return Policy::disjoint();
- }
-
         typedef side::side_by_triangle<coordinate_type> side;
         side_info sides;
 
- // 2) Calculate sides
- // Note: Do NOT yet calculate the determinant here, but use the SIDE strategy.
- // Determinant calculation is not robust; side (orient) can be made robust
- // (and is much robuster even without measures)
+ sides.set<0>
+ (
+ side::apply(detail::get_from_index<0>(b)
+ , detail::get_from_index<1>(b)
+ , detail::get_from_index<0>(a)),
+ side::apply(detail::get_from_index<0>(b)
+ , detail::get_from_index<1>(b)
+ , detail::get_from_index<1>(a))
+ );
         sides.set<1>
             (
                 side::apply(detail::get_from_index<0>(a)
@@ -159,25 +139,22 @@
                     , detail::get_from_index<1>(b))
             );
 
- if (sides.same<1>())
- {
- // Both points are at same side of other segment, we can leave
- return Policy::disjoint();
- }
+ bool collinear = sides.collinear();
 
- // 2b) For other segment
- sides.set<0>
- (
- side::apply(detail::get_from_index<0>(b)
- , detail::get_from_index<1>(b)
- , detail::get_from_index<0>(a)),
- side::apply(detail::get_from_index<0>(b)
- , detail::get_from_index<1>(b)
- , detail::get_from_index<1>(a))
- );
+ if ((sides.zero<0>() && ! sides.zero<1>()) || (sides.zero<1>() && ! sides.zero<0>()))
+ {
+ // If one of the segments is collinear, the other must be as well.
+ // So handle it as collinear.
+ // (In float/double epsilon margins it can easily occur that one or two of them are -1/1)
+ // sides.debug();
+ sides.set<0>(0,0);
+ sides.set<1>(0,0);
+ collinear = true;
+ }
 
- if (sides.same<0>())
+ if (sides.same<0>() || sides.same<1>())
         {
+ // Both points are at same side of other segment, we can leave
             return Policy::disjoint();
         }
 
@@ -192,82 +169,95 @@
             return Policy::degenerate(b, false);
         }
 
- bool collinear = sides.collinear();
-
- // Get the same type, but at least a double (also used for divisions)
         typedef typename select_most_precise
             <
                 coordinate_type, double
>::type promoted_type;
 
- // Calculate the determinant/2D cross product
- // (Note, because we only check on zero,
- // the order a/b does not care)
- promoted_type const d = geometry::detail::determinant
- <
- promoted_type
- >(dx_a, dy_a, dx_b, dy_b);
-
- // Determinant d should be nonzero.
- // If it is zero, we have an robustness issue here,
- // (and besides that we cannot divide by it)
- promoted_type const pt_zero = promoted_type();
- if(math::equals(d, pt_zero) && ! collinear)
- {
-#ifdef BOOST_GEOMETRY_DEBUG_INTERSECTION
- std::cout << "Determinant zero? Type : "
- << typeid(coordinate_type).name()
- << std::endl;
-
- std::cout << " dx_a : " << dx_a << std::endl;
- std::cout << " dy_a : " << dy_a << std::endl;
- std::cout << " dx_b : " << dx_b << std::endl;
- std::cout << " dy_b : " << dy_b << std::endl;
-
- std::cout << " side a <-> b.first : " << sides.get<0,0>() << std::endl;
- std::cout << " side a <-> b.second: " << sides.get<0,1>() << std::endl;
- std::cout << " side b <-> a.first : " << sides.get<1,0>() << std::endl;
- std::cout << " side b <-> a.second: " << sides.get<1,1>() << std::endl;
-#endif
-
- if (sides.as_collinear())
- {
- sides.set<0>(0,0);
- sides.set<1>(0,0);
- collinear = true;
- }
- else
- {
- return Policy::error("Determinant zero!");
- }
- }
+ // r: ratio 0-1 where intersection divides A/B
+ // (only calculated for non-collinear segments)
+ promoted_type r;
+ if (! collinear)
+ {
+ // Calculate determinants - Cramers rule
+ coordinate_type const wx = get<0, 0>(a) - get<0, 0>(b);
+ coordinate_type const wy = get<0, 1>(a) - get<0, 1>(b);
+ coordinate_type const d = geometry::detail::determinant<coordinate_type>(dx_a, dy_a, dx_b, dy_b);
+ coordinate_type const da = geometry::detail::determinant<coordinate_type>(dx_b, dy_b, wx, wy);
+
+ coordinate_type const zero = coordinate_type();
+ if (math::equals(d, zero))
+ {
+ // This is still a collinear case (because of FP imprecision this can occur here)
+ // sides.debug();
+ sides.set<0>(0,0);
+ sides.set<1>(0,0);
+ collinear = true;
+ }
+ else
+ {
+ r = da / d;
+
+ // Ratio should lie between 0 and 1
+ // Also these three conditions might be of FP imprecision, the segments were actually (nearly) collinear
+ promoted_type const zero = 0;
+ promoted_type const one = 1;
+ promoted_type const epsilon = std::numeric_limits<double>::epsilon();
+ if (r < zero)
+ {
+ if (r < -epsilon)
+ {
+ return Policy::disjoint();
+ }
+ r = zero;
+ }
+ else if (r > one)
+ {
+ if (r > one + epsilon)
+ {
+ return Policy::disjoint();
+ }
+ r = one;
+ }
+ }
+ }
 
         if(collinear)
         {
- // Segments are collinear. We'll find out how.
- if (math::equals(dx_b, zero))
+ if (math::abs(dx_a) + math::abs(dx_b) < math::abs(dy_a) + math::abs(dy_a))
             {
- // Vertical -> Check y-direction
- return relate_collinear(a, b,
- ay_1, ay_2, by_1, by_2,
- ay_swapped, by_swapped);
+ // Y direction contains larger segments (maybe dx is zero)
+ return relate_collinear<1>(a, b);
             }
             else
             {
- // Check x-direction
- return relate_collinear(a, b,
- ax_1, ax_2, bx_1, bx_2,
- ax_swapped, bx_swapped);
+ return relate_collinear<0>(a, b);
             }
         }
 
- return Policy::segments_intersect(sides,
+ return Policy::segments_intersect(sides, r,
             dx_a, dy_a, dx_b, dy_b,
             a, b);
     }
 
 private :
 
+ template <std::size_t Dimension>
+ static inline return_type relate_collinear(segment_type1 const& a,
+ segment_type2 const& b)
+ {
+ coordinate_type a_1, a_2, b_1, b_2;
+ bool a_swapped = false, b_swapped = false;
+ detail::segment_arrange<Dimension>(a, a_1, a_2, a_swapped);
+ detail::segment_arrange<Dimension>(b, b_1, b_2, b_swapped);
+ if (math::smaller(a_2, b_1) || math::larger(a_1, b_2))
+ //if (a_2 < b_1 || a_1 > b_2)
+ {
+ return Policy::disjoint();
+ }
+ return relate_collinear(a, b, a_1, a_2, b_1, b_2, a_swapped, b_swapped);
+ }
+
     /// Relate segments known collinear
     static inline return_type relate_collinear(segment_type1 const& a
             , segment_type2 const& b

Modified: trunk/boost/geometry/strategies/side_info.hpp
==============================================================================
--- trunk/boost/geometry/strategies/side_info.hpp (original)
+++ trunk/boost/geometry/strategies/side_info.hpp 2012-03-10 14:11:17 EST (Sat, 10 Mar 2012)
@@ -67,21 +67,28 @@
             && sides[1].second == 0;
     }
 
- // If one of the segments is collinear, the other must be as well.
- // So handle it as collinear.
- // (In floating point margins it can occur that one of them is 1!)
- inline bool as_collinear() const
+ template <int Which>
+ inline bool zero() const
     {
- return sides[0].first * sides[0].second == 0
- || sides[1].first * sides[1].second == 0;
+ return sides[Which].first == 0 && sides[Which].second == 0;
     }
 
+ inline void debug() const
+ {
+ std::cout << sides[0].first << " "
+ << sides[0].second << " "
+ << sides[1].first << " "
+ << sides[1].second
+ << std::endl;
+ }
+
+
     inline void reverse()
     {
         std::swap(sides[0], sides[1]);
     }
 
-private :
+//private :
     std::pair<int, int> sides[2];
 
 };

Modified: trunk/boost/geometry/util/math.hpp
==============================================================================
--- trunk/boost/geometry/util/math.hpp (original)
+++ trunk/boost/geometry/util/math.hpp 2012-03-10 14:11:17 EST (Sat, 10 Mar 2012)
@@ -62,6 +62,28 @@
     }
 };
 
+template <typename Type, bool IsFloatingPoint>
+struct smaller
+{
+ static inline bool apply(Type const& a, Type const& b)
+ {
+ return a < b;
+ }
+};
+
+template <typename Type>
+struct smaller<Type, true>
+{
+ static inline bool apply(Type const& a, Type const& b)
+ {
+ if (equals<Type, true>::apply(a, b))
+ {
+ return false;
+ }
+ return a < b;
+ }
+};
+
 
 template <typename Type, bool IsFloatingPoint>
 struct equals_with_epsilon : public equals<Type, IsFloatingPoint> {};
@@ -126,6 +148,28 @@
>::apply(a, b);
 }
 
+template <typename T1, typename T2>
+inline bool smaller(T1 const& a, T2 const& b)
+{
+ typedef typename select_most_precise<T1, T2>::type select_type;
+ return detail::smaller
+ <
+ select_type,
+ boost::is_floating_point<select_type>::type::value
+ >::apply(a, b);
+}
+
+template <typename T1, typename T2>
+inline bool larger(T1 const& a, T2 const& b)
+{
+ typedef typename select_most_precise<T1, T2>::type select_type;
+ return detail::smaller
+ <
+ select_type,
+ boost::is_floating_point<select_type>::type::value
+ >::apply(b, a);
+}
+
 
 
 double const d2r = geometry::math::pi<double>() / 180.0;

Modified: trunk/libs/geometry/test/algorithms/difference.cpp
==============================================================================
--- trunk/libs/geometry/test/algorithms/difference.cpp (original)
+++ trunk/libs/geometry/test/algorithms/difference.cpp 2012-03-10 14:11:17 EST (Sat, 10 Mar 2012)
@@ -7,6 +7,8 @@
 // Version 1.0. (See accompanying file LICENSE_1_0.txt or copy at
 // http://www.boost.org/LICENSE_1_0.txt)
 
+#define TEST_ISOVIST
+
 //#define BOOST_GEOMETRY_CHECK_WITH_POSTGIS
 
 //#define BOOST_GEOMETRY_DEBUG_SEGMENT_IDENTIFIER
@@ -258,12 +260,17 @@
     ***/
 
 #ifdef _MSC_VER
- // Isovist (submitted by Brandon during Formal Review)
+#ifdef TEST_ISOVIST
     test_one<polygon, polygon, polygon>("isovist",
         isovist1[0], isovist1[1],
- 4, 0, 0.279121891701124,
- 4, 0, 224.889211358929,
- 0.01);
+ if_typed_tt<ct>(4, 2), 0, 0.279121891701124,
+ if_typed_tt<ct>(4, 3), 0, if_typed_tt<ct>(224.889211358929, 223.777),
+ if_typed_tt<ct>(0.001, 0.2));
+
+ // SQL Server gives: 0.279121891701124 and 224.889211358929
+ // PostGIS gives: 0.279121991127244 and 224.889205853156
+
+#endif
 
     test_one<polygon, polygon, polygon>("ggl_list_20110306_javier",
         ggl_list_20110306_javier[0], ggl_list_20110306_javier[1],
@@ -277,11 +284,14 @@
         1, 0, 3200.4,
         0.01);
 
- test_one<polygon, polygon, polygon>("ggl_list_20110716_enrico",
- ggl_list_20110716_enrico[0], ggl_list_20110716_enrico[1],
- 3, 0, 35723.8506317139,
- 1, 0, 58456.4964294434
- );
+ if (! boost::is_same<ct, float>::value)
+ {
+ test_one<polygon, polygon, polygon>("ggl_list_20110716_enrico",
+ ggl_list_20110716_enrico[0], ggl_list_20110716_enrico[1],
+ 3, 0, 35723.8506317139,
+ 1, 0, 58456.4964294434
+ );
+ }
 
     test_one<polygon, polygon, polygon>("ggl_list_20110820_christophe",
         ggl_list_20110820_christophe[0], ggl_list_20110820_christophe[1],
@@ -300,7 +310,7 @@
     // Because we cannot predict this, we only test for MSVC
     test_one<polygon, polygon, polygon>("ggl_list_20110627_phillip",
         ggl_list_20110627_phillip[0], ggl_list_20110627_phillip[1],
- if_typed<ct, double>(0, 1), 0,
+ if_typed_tt<ct>(1, 0), 0,
             if_typed_tt<ct>(0.0000000000001105367, 0.0),
         1, 0, 3577.40960816756,
         0.01
@@ -465,6 +475,7 @@
     test_all<bg::model::d2::point_xy<float> >();
 
 #ifdef HAVE_TTMATH
+ std::cout << "Testing TTMATH" << std::endl;
     test_all<bg::model::d2::point_xy<ttmath_big> >();
     //test_difference_parcel_precision<ttmath_big>();
 #endif

Modified: trunk/libs/geometry/test/algorithms/intersection.cpp
==============================================================================
--- trunk/libs/geometry/test/algorithms/intersection.cpp (original)
+++ trunk/libs/geometry/test/algorithms/intersection.cpp 2012-03-10 14:11:17 EST (Sat, 10 Mar 2012)
@@ -15,6 +15,8 @@
 #include <iostream>
 #include <string>
 
+#define TEST_ISOVIST
+
 //#define BOOST_GEOMETRY_DEBUG_SEGMENT_IDENTIFIER
 //#define BOOST_GEOMETRY_DEBUG_INTERSECTION
 //#define BOOST_GEOMETRY_DEBUG_TRAVERSE
@@ -164,14 +166,17 @@
 
     typedef typename bg::coordinate_type<Polygon>::type ct;
 
+#ifdef TEST_ISOVIST
 #ifdef _MSC_VER
- // Isovist (submitted by Brandon during Formal Review)
     test_one<Polygon, Polygon, Polygon>("isovist",
         isovist1[0], isovist1[1],
- 1,
- if_typed<ct, float>(19, if_typed<ct, double>(20, 20)),
- 88.19203,
- if_typed<ct, float>(0.5, if_typed<ct, double>(0.1, 0.01)));
+ 1, 20, 88.19203,
+ if_typed_tt<ct>(0.01, 0.1));
+
+ // SQL Server gives: 88.1920416352664
+ // PostGIS gives: 88.19203677911
+
+#endif
 #endif
 
     //std::cout << typeid(ct).name() << std::endl;
@@ -191,13 +196,20 @@
         1, if_typed_tt<ct>(6, 5), 11151.6618);
 
 #ifdef _MSC_VER // gcc/linux behaves differently
- test_one<Polygon, Polygon, Polygon>("ggl_list_20110716_enrico",
- ggl_list_20110716_enrico[0], ggl_list_20110716_enrico[1],
- 3,
- if_typed<ct, float>(19, if_typed<ct, double>(22, 21)),
- 35723.8506317139);
+ if (! boost::is_same<ct, float>::value)
+ {
+ test_one<Polygon, Polygon, Polygon>("ggl_list_20110716_enrico",
+ ggl_list_20110716_enrico[0], ggl_list_20110716_enrico[1],
+ 3,
+ if_typed<ct, float>(19, if_typed<ct, double>(22, 21)),
+ 35723.8506317139);
+ }
 #endif
 
+ test_one<Polygon, Polygon, Polygon>("buffer_rt_f", buffer_rt_f[0], buffer_rt_f[1],
+ 1, 4, 0.00029437899183903937, 0.01);
+
+
     return;
 
 
@@ -496,6 +508,7 @@
     test_all<bg::model::d2::point_xy<float> >();
 
 #if defined(HAVE_TTMATH)
+ std::cout << "Testing TTMATH" << std::endl;
     test_all<bg::model::d2::point_xy<ttmath_big> >();
 #endif
 

Modified: trunk/libs/geometry/test/algorithms/overlay/overlay_cases.hpp
==============================================================================
--- trunk/libs/geometry/test/algorithms/overlay/overlay_cases.hpp (original)
+++ trunk/libs/geometry/test/algorithms/overlay/overlay_cases.hpp 2012-03-10 14:11:17 EST (Sat, 10 Mar 2012)
@@ -480,6 +480,7 @@
     "POLYGON((184861.118 606901.158,184893.787 606898.483,184925.043 606913.4,184927.174 606951.759,184912.9 606987.146,184877.87 606989.232,184885.103 607023.774,184899.058 607022.743,184906.008 607044.948,184966.465 607025.02,184968.442 606961.3,185024.768 606947.402,185024.544 606941.355,185027.007 606937.323,185030.366 606934.187,185035.516 606933.963,185040.442 606935.531,185042.905 606939.115,185088.364 606931.385,185089.139 607015.509,185095.2 607011.3,185118.827 606995.545,185126.813 606991.995,185177.727 606973.799,185181.482 606966.676,185193.571 606977.795,185193.711 606960.3,185189.352 606779.02,185167.515 606783.844,185086.96 606801.241,185011.707 606817.809,185000 606819.304,184994.034 606819.794,184976.398 606819.572,184956.654 606817.131,184934.913 606813.137,184893.097 606804.927,184884.445 606831.555,184866.919 606883.481,184861.118 606901.158),(184907.556 607013.301,184905.782 607009.972,184906.004 607005.978,184908.444 606998.877,184912.215 606994.218,184919.314 606993.996,184922.42 6069
95.771,184925.747 606998.877,184926.413 607002.872,184925.747 607007.753,184922.42 607012.191,184917.096 607015.298,184911.771 607015.298,184907.556 607013.301))"};
 
 
+// Isovist (submitted by Brandon during Formal Review)
 static std::string isovist[2] =
     {
     "POLYGON((37.29449462890625 1.7902572154998779, 46.296027072709599 -2.4984308554828116, 45.389434814453125 -4.5143837928771973, 47.585065917176543 -6.1314922196594779, 46.523914387974358 -8.5152102535033496, 42.699958801269531 -4.4278755187988281, 42.577877044677734 -4.4900407791137695, 42.577911376953125 -4.4901103973388672, 40.758884429931641 -5.418975830078125, 40.6978759765625 -5.4500408172607422, 41.590042114257813 -7.2021245956420898, 57.297810222148939 -37.546793343968417, 50.974888957147442 -30.277285722290763, 37.140213012695313 1.3446992635726929, 37.000419616699219 1.664225697517395, 37.29449462890625 1.7902572154998779))",
@@ -559,5 +560,24 @@
     "POLYGON((-95269304 222758,-95260668 419862,-95234760 615696,-95192088 808228,-95132906 996442,-95057214 1178814,-94966028 1354074,-94860110 1520444,-94739968 1676908,-94606618 1822450,-94999048 2214880,-95165164 2033778,-95314770 1838706,-95446850 1631442,-95560388 1413510,-95654368 1186434,-95728282 951992,-95781368 711962,-95813626 468376,-95824294 222758,-95269304 222758))"
     };
 
+static std::string ggl_list_20120229_volker[3] =
+ {
+ "POLYGON((1716 1554,2076 2250,2436 2352,2796 1248,3156 2484,3516 2688,3516 2688,3156 2484,2796 1248,2436 2352,2076 2250, 1716 1554))",
+ "POLYGON((2500 1600,2500 2300,3200 2300,3200 1600,2500 1600))",
+ "POLYGON((1716 1554,2076 2250,2436 2352,2796 1248,3156 2484,3516 2688,3156 2483,2796 1247,2436 2351,2076 2249, 1716 1554))",
+ };
+
+static std::string buffer_rt_a[2] =
+ {
+ "POLYGON((1 7,1 8,1.0012 8.04907,1.00482 8.09802,1.01082 8.14673,1.01921 8.19509,1.02997 8.24298,1.04306 8.29028,1.05846 8.33689,1.07612 8.38268,1.09601 8.42756,1.11808 8.4714,1.14227 8.5141,1.16853 8.55557,1.19679 8.5957,1.22699 8.63439,1.25905 8.67156,1.29289 8.70711,1.32844 8.74095,1.36561 8.77301,1.4043 8.80321,1.44443 8.83147,1.4859 8.85773,1.5286 8.88192,1.57244 8.90399,1.61732 8.92388,1.66311 8.94154,1.70972 8.95694,1.75702 8.97003,1.80491 8.98079,1.85327 8.98918,1.90198 8.99518,1.95093 8.9988,2 9,3 9,3.04907 8.9988,3.09802 8.99518,3.14673 8.98918,3.19509 8.98079,3.24298 8.97003,3.29028 8.95694,3.33689 8.94154,3.38268 8.92388,3.42756 8.90399,3.4714 8.88192,3.5141 8.85773,3.55557 8.83147,3.5957 8.80321,3.63439 8.77301,3.67156 8.74095,3.70711 8.70711,3.74095 8.67156,3.77301 8.63439,3.80321 8.5957,3.83147 8.55557,3.85773 8.5141,3.88192 8.4714,3.90399 8.42756,3.92388 8.38268,3.94154 8.33689,3.95694 8.29028,3.97003 8.24298,3.98079 8.19509,3.98918 8.14673,3.99518 8.09802,3.9988 8.04907,4 8,4 7,3.9988 6
.95093,3.99518 6.90198,3.98918 6.85327,3.98079 6.80491,3.97003 6.75702,3.95694 6.70972,3.94154 6.66311,3.92388 6.61732,3.90399 6.57244,3.88192 6.5286,3.85773 6.4859,3.83147 6.44443,3.80321 6.4043,3.77301 6.36561,3.74095 6.32844,3.70711 6.29289,3.67156 6.25905,3.63439 6.22699,3.5957 6.19679,3.55557 6.16853,3.5141 6.14227,3.4714 6.11808,3.42756 6.09601,3.38268 6.07612,3.33689 6.05846,3.29028 6.04306,3.24298 6.02997,3.19509 6.01921,3.14673 6.01082,3.09802 6.00482,3.04907 6.0012,3 6,2 6,1.95093 6.0012,1.90198 6.00482,1.85327 6.01082,1.80491 6.01921,1.75702 6.02997,1.70972 6.04306,1.66311 6.05846,1.61732 6.07612,1.57244 6.09601,1.5286 6.11808,1.4859 6.14227,1.44443 6.16853,1.4043 6.19679,1.36561 6.22699,1.32844 6.25905,1.29289 6.29289,1.25905 6.32844,1.22699 6.36561,1.19679 6.4043,1.16853 6.44443,1.14227 6.4859,1.11808 6.5286,1.09601 6.57244,1.07612 6.61732,1.05846 6.66311,1.04306 6.70972,1.02997 6.75702,1.01921 6.80491,1.01082 6.85327,1.00482 6.90198,1.0012 6.95093,1 7))",
+ "POLYGON((3 6,4 6,4.04907 5.9988,4.09802 5.99518,4.14673 5.98918,4.19509 5.98079,4.24298 5.97003,4.29028 5.95694,4.33689 5.94154,4.38268 5.92388,4.42756 5.90399,4.4714 5.88192,4.5141 5.85773,4.55557 5.83147,4.5957 5.80321,4.63439 5.77301,4.67156 5.74095,4.70711 5.70711,4.74095 5.67156,4.77301 5.63439,4.80321 5.5957,4.83147 5.55557,4.85773 5.5141,4.88192 5.4714,4.90399 5.42756,4.92388 5.38268,4.94154 5.33689,4.95694 5.29028,4.97003 5.24298,4.98079 5.19509,4.98918 5.14673,4.99518 5.09802,4.9988 5.04907,5 5,5 4,4.9988 3.95093,4.99518 3.90198,4.98918 3.85327,4.98079 3.80491,4.97003 3.75702,4.95694 3.70972,4.94154 3.66311,4.92388 3.61732,4.90399 3.57244,4.88192 3.5286,4.85773 3.4859,4.83147 3.44443,4.80321 3.4043,4.77301 3.36561,4.74095 3.32844,4.70711 3.29289,4.67156 3.25905,4.63439 3.22699,4.5957 3.19679,4.55557 3.16853,4.5141 3.14227,4.4714 3.11808,4.42756 3.09601,4.38268 3.07612,4.33689 3.05846,4.29028 3.04306,4.24298 3.02997,4.19509 3.01921,4.14673 3.01082,4.09802 3.00482,4.04907 3.0012,4 3,3 3,3 3,2 3,
1.95093 3.0012,1.90198 3.00482,1.85327 3.01082,1.80491 3.01921,1.75702 3.02997,1.70972 3.04306,1.66311 3.05846,1.61732 3.07612,1.57244 3.09601,1.5286 3.11808,1.4859 3.14227,1.44443 3.16853,1.4043 3.19679,1.36561 3.22699,1.32844 3.25905,1.29289 3.29289,1.25905 3.32844,1.22699 3.36561,1.19679 3.4043,1.16853 3.44443,1.14227 3.4859,1.11808 3.5286,1.09601 3.57244,1.07612 3.61732,1.05846 3.66311,1.04306 3.70972,1.02997 3.75702,1.01921 3.80491,1.01082 3.85327,1.00482 3.90198,1.0012 3.95093,1 4,1 5,1.0012 5.04907,1.00482 5.09802,1.01082 5.14673,1.01921 5.19509,1.02997 5.24298,1.04306 5.29028,1.05846 5.33689,1.07612 5.38268,1.09601 5.42756,1.11808 5.4714,1.14227 5.5141,1.16853 5.55557,1.19679 5.5957,1.22699 5.63439,1.25905 5.67156,1.29289 5.70711,1.32844 5.74095,1.36561 5.77301,1.4043 5.80321,1.44443 5.83147,1.4859 5.85773,1.5286 5.88192,1.57244 5.90399,1.61732 5.92388,1.66311 5.94154,1.70972 5.95694,1.75702 5.97003,1.80491 5.98079,1.85327 5.98918,1.90198 5.99518,1.95093 5.9988,2 6,3 6,3 6,3 6,3 6,3 6,3 6,3 6,3 6,3 6
,3 6,3 6,3 6,3 6,3 6,3 6,3 6,3 6,3 6,3 6,3 6,3 6,3 6,3 6,3 6,3 6,3 6,3 6,3 6,3 6,3 6,3 6,3 6,3 6))"
+ };
+
+static std::string buffer_rt_f[2] =
+ {
+ "POLYGON((-0.29999999999999999 6.0000000000000000,-0.29999999999999999 7.0000000000000000,-0.30000000000000027 7.7242640687119302,0.21213203435596423 7.2121320343559638,1.2121320343559643 6.2121320343559638,1.7242640687119293 5.7000000000000002,1.0000000000000000 5.7000000000000002,0.00000000000000000 5.7000000000000002,-0.30000000000000027 5.7000000000000002,-0.29999999999999999 6.0000000000000000))",
+ "POLYGON((1.3000000000000000 9.0000000000000000,1.3000000000000000 8.0000000000000000,1.3000000000000007 7.7000000000000002,1.0000000000000000 7.7000000000000002,0.00000000000000000 7.7000000000000002,-0.29999999999999982 7.7000000000000002,-0.29999999999999999 8.0000000000000000,-0.29999999999999999 9.0000000000000000,-0.29999999999999982 9.3000000000000007,0.00000000000000000 9.3000000000000007,1.0000000000000000 9.3000000000000007,1.3000000000000007 9.3000000000000007,1.3000000000000000 9.0000000000000000))"
+ };
+
 
 #endif // BOOST_GEOMETRY_TEST_OVERLAY_CASES_HPP

Modified: trunk/libs/geometry/test/algorithms/union.cpp
==============================================================================
--- trunk/libs/geometry/test/algorithms/union.cpp (original)
+++ trunk/libs/geometry/test/algorithms/union.cpp 2012-03-10 14:11:17 EST (Sat, 10 Mar 2012)
@@ -15,6 +15,8 @@
 #include <iostream>
 #include <string>
 
+#define TEST_ISOVIST
+
 //#define BOOST_GEOMETRY_DEBUG_ASSEMBLE
 //#define BOOST_GEOMETRY_DEBUG_IDENTIFIER
 
@@ -237,16 +239,40 @@
         129904.197692871);
 #endif
 
+#ifdef TEST_ISOVIST
 #ifdef _MSC_VER
- // Isovist (submitted by Brandon during Formal Review)
     test_one<Polygon, Polygon, Polygon>("isovist",
         isovist1[0], isovist1[1],
         1,
         0,
- if_typed<ct, float>(71,
- if_typed<ct, double>(70, 73)),
- 313.36036462);
+ if_typed<ct, float>(71, if_typed<ct, double>(70, 73)),
+ 313.36036462, 0.1);
+
+ // SQL Server gives: 313.360374193241
+ // PostGIS gives: 313.360364623393
+
+#endif
 #endif
+
+ // Ticket 5103 https://svn.boost.org/trac/boost/ticket/5103
+ // This ticket was actually reported for Boost.Polygon
+ // We check it for Boost.Geometry as well.
+ // SQL Server gives: 2515271331437.69
+ // PostGIS gives: 2515271327070.52
+ // Boost.Geometry gives: 2515271327070.5237746891 (ttmath)
+ // 2515271327070.5156 (double)
+ // 2515271320603.0000 (int)
+ // Note the int-test was tested outside of this unit test. It is in two points 0.37 off (logical for an int).
+ // Because of the width of the polygon (400000 meter) this causes a substantial difference.
+
+ test_one<Polygon, Polygon, Polygon>("ticket_5103", ticket_5103[0], ticket_5103[1],
+ 1, 0, 25, 2515271327070.5);
+
+ test_one<Polygon, Polygon, Polygon>("buffer_rt_a", buffer_rt_a[0], buffer_rt_a[1],
+ 1, 0, 265, 19.280667);
+
+ test_one<Polygon, Polygon, Polygon>("buffer_rt_f", buffer_rt_f[0], buffer_rt_f[1],
+ 1, 0, if_typed<ct, double>(22, 23), 4.60853);
 }
 
 template <typename P>
@@ -302,20 +328,6 @@
     test_one<polygon, box, polygon>("box_poly8", "box(0 0, 3 3)",
             "POLYGON((2 2, 1 4, 2 4, 3 3, 2 2))",
                 1, 0, 8, 10.25);
-
- // Ticket 5103 https://svn.boost.org/trac/boost/ticket/5103
- // This ticket was actually reported for Boost.Polygon
- // but it is apparently a difficult case so we check it for Boost.Geometry as well.
- // SQL Server gives: 2515271331437.69
- // PostGIS gives: 2515271327070.52
- // Boost.Geometry gives: 2515271327070.5237746891 (ttmath)
- // 2515271327070.5156 (double)
- // 2515271320603.0000 (int)
- // Note the int-test was tested externally - it is in two points 0.37 off (makes sense).
- // Because of the width of the polygon (400000 meter) this might indeed cause a substantial difference.
-
- test_one<polygon, polygon, polygon>("ticket_5103", ticket_5103[0], ticket_5103[1],
- 1, 0, 25, 2515271327070.5);
 }
 
 
@@ -328,6 +340,7 @@
     //test_all<bg::model::d2::point_xy<long double> >();
 
 #if defined(HAVE_TTMATH)
+ std::cout << "Testing TTMATH" << std::endl;
     test_all<bg::model::d2::point_xy<ttmath_big> >();
 #endif
 #endif


Boost-Commit list run by bdawes at acm.org, david.abrahams at rcn.com, gregod at cs.rpi.edu, cpdaniel at pacbell.net, john at johnmaddock.co.uk