From 58bc22dfb2d9f909ea4f4c4147a370bb4d5d90d7 Mon Sep 17 00:00:00 2001 From: Menelaos Karavelas Date: Thu, 27 Aug 2015 16:32:48 +0300 Subject: [PATCH 01/22] [algorithms][overlay][union] The rings produced during the traversal phase for the union operation for areal geometries may exhibit self intersections (t:u/u and m:u/u turns); this commit adds a post-processing step to traversal, during which step the rings produced during the traversal are split into simple rings; this is done by inserting m:u/u turns into the rings produced by the traversal, and then splitting the rings into smaller ones using their duplicate points (which correspond to t:u/u turns); --- .../geometry/algorithms/detail/overlay/overlay.hpp | 4 + .../algorithms/detail/overlay/split_rings.hpp | 417 +++++++++++++++++++++ 2 files changed, 421 insertions(+) create mode 100644 include/boost/geometry/algorithms/detail/overlay/split_rings.hpp diff --git a/include/boost/geometry/algorithms/detail/overlay/overlay.hpp b/include/boost/geometry/algorithms/detail/overlay/overlay.hpp index baf9d4777d..57eb25146f 100644 --- a/include/boost/geometry/algorithms/detail/overlay/overlay.hpp +++ b/include/boost/geometry/algorithms/detail/overlay/overlay.hpp @@ -40,6 +40,7 @@ #include #include #include +#include #include #include @@ -249,6 +250,9 @@ std::cout << "traverse" << std::endl; select_rings(geometry1, geometry2, turn_info_per_ring, selected_ring_properties); + // split the rings into simple rings + split_rings::apply(rings, robust_policy); + // Add rings created during traversal { ring_identifier id(2, 0, -1); diff --git a/include/boost/geometry/algorithms/detail/overlay/split_rings.hpp b/include/boost/geometry/algorithms/detail/overlay/split_rings.hpp new file mode 100644 index 0000000000..32520b37b2 --- /dev/null +++ b/include/boost/geometry/algorithms/detail/overlay/split_rings.hpp @@ -0,0 +1,417 @@ +// Boost.Geometry (aka GGL, Generic Geometry Library) + +// 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_GEOMETRY_ALGORITHMS_DETAIL_OVERLAY_SPLIT_RINGS_HPP +#define BOOST_GEOMETRY_ALGORITHMS_DETAIL_OVERLAY_SPLIT_RINGS_HPP + +#include +#include +#include +#include + +#include + +#include +#include +#include + +#include +#include + +#include + +#include + +#include +#include +#include +#include + + +namespace boost { namespace geometry +{ + +namespace detail { namespace overlay +{ + + +// the ring_as_dcl class implements a ring using a doubly-connected list; +// it is a model of a bidirectional-traversal Boost.Range and +// supports: +// (1) constant-time push_back(), size() and empty() operations; +// (2) splitting at two vertices (in order to create two rings out of +// the initial ring) also in constant time. +// +// The second operation cannot be done in constant time using typical +// random access ranges, which is precisely the reason why this class +// has been introduced and used. +template +< + typename Point, + closure_selector Closure, + typename List = std::list +> +class ring_as_dcl +{ +public: + typedef Point point_type; + typedef List list_type; + typedef typename List::size_type size_type; + typedef typename List::iterator iterator; + typedef typename List::const_iterator const_iterator; + + ring_as_dcl() + : m_list() + {} + + inline void push_back(Point const& point) + { + m_list.push_back(point); + } + + inline void split_at(iterator pos1, iterator pos2, ring_as_dcl& other) + { + // for this method to work, I think that the iterator pos1 + // must preceed iterator pos2 in the original list node sequence + other.m_list.splice(other.m_list.end(), m_list, pos1, pos2); + + if (BOOST_GEOMETRY_CONDITION(Closure == closed)) + { + other.m_list.push_back(*other.m_list.begin()); + } + } + + iterator begin() { return m_list.begin(); } + const_iterator begin() const { return m_list.begin(); } + + iterator end() { return m_list.end(); } + const_iterator end() const { return m_list.end(); } + +private: + List m_list; +}; + + +template +class find_duplicate_points +{ + static bool const is_closed = Closure == closed; + + struct point_iterator_less + { + template + inline bool operator()(PointIterator it1, PointIterator it2) const + { + return geometry::less + < + typename std::iterator_traits::value_type + >()(*it1, *it2); + } + }; + +public: + template + static inline bool apply(Ring& ring, Iterator& pos1, Iterator& pos2) + { + typedef typename Ring::iterator iterator_type; + typedef std::set point_set_type; + point_set_type point_set; + + iterator_type last = is_closed ? --boost::end(ring) : boost::end(ring); + for (iterator_type it = boost::begin(ring); it != last; ++it) + { + typename point_set_type::const_iterator pit = point_set.find(it); + if (pit != point_set.end()) + { + pos1 = *pit; + pos2 = it; + return true; + } + else + { + point_set.insert(it); + } + } + + // initialize pos1 and pos2 with something + pos1 = boost::begin(ring); + pos2 = pos1; + return false; + } +}; + + +template +struct split_ring +{ + template + static inline void apply(Ring const&, RingCollection&, RobustPolicy const&) + { + } +}; + +// specialization for union +// TODO: add another specialization for intersection once implemented +template +class split_ring +{ + typedef turn_info + < + typename point_type::type, + typename geometry::segment_ratio_type + < + typename point_type::type, RobustPolicy + >::type + > turn_type; + + typedef std::deque turns_container_type; + + struct no_interrupt_policy + { + static bool const enabled = false; + static bool const has_intersections = false; + + template + static inline bool apply(Range const&) + { + return false; + } + }; + + template + static inline + typename Turn::turn_operation_type get_correct_op(Turn const& t) + { + return + (t.operations[0].fraction.is_zero() + || t.operations[0].fraction.is_one()) + ? + t.operations[1] + : + t.operations[0] + ; + } + + template + struct muu_turn_less + { + bool operator()(MUU_Turn const& t1, MUU_Turn const& t2) const + { + BOOST_GEOMETRY_ASSERT(t1.method == method_touch_interior); + BOOST_GEOMETRY_ASSERT(t1.operations[0].operation + == operation_union); + BOOST_GEOMETRY_ASSERT(t1.operations[1].operation + == operation_union); + + BOOST_GEOMETRY_ASSERT(t2.method == method_touch_interior); + BOOST_GEOMETRY_ASSERT(t2.operations[0].operation + == operation_union); + BOOST_GEOMETRY_ASSERT(t2.operations[1].operation + == operation_union); + + typename MUU_Turn::turn_operation_type op1 = get_correct_op(t1); + typename MUU_Turn::turn_operation_type op2 = get_correct_op(t2); + + BOOST_GEOMETRY_ASSERT(! op1.fraction.is_zero() + && ! op1.fraction.is_one()); + BOOST_GEOMETRY_ASSERT(! op2.fraction.is_zero() + && ! op2.fraction.is_one()); + + + if (op1.seg_id.segment_index != op2.seg_id.segment_index) + { + return op1.seg_id.segment_index < op2.seg_id.segment_index; + } + return op1.fraction < op2.fraction; + } + }; + + template + static inline void get_self_turns(Ring const& ring, + turns_container_type& turns, + RobustPolicy const& robust_policy, + InterruptPolicy const& policy) + { + geometry::self_turns + < + get_turn_info + >(ring, robust_policy, turns, policy); + } + + static inline void get_self_turns(Ring const& ring, + turns_container_type& turns, + RobustPolicy const& robust_policy) + { + get_self_turns(ring, turns, robust_policy, no_interrupt_policy()); + } + + template + static inline void insert_muu_turns(Ring const& ring, + MUU_Turns const& muu_turns, + RingOut& ring_out) + { + typedef typename boost::range_iterator + < + MUU_Turns const + >::type iterator_type; + + typename boost::range_size::type point_index = 0; + for (iterator_type it = muu_turns.begin(); it != muu_turns.end(); ++it) + { + signed_size_type segment_index + = get_correct_op(*it).seg_id.segment_index; + + while (point_index <= segment_index) + { + ring_out.push_back(ring[point_index]); + ++point_index; + } + + ring_out.push_back(it->point); + } + + while (point_index < ring.size()) + { + ring_out.push_back(ring[point_index]); + ++point_index; + } + } + + template + static inline void copy_to_collection(RingIn const& ring, + Collection& collection) + { + typedef typename boost::range_value::type ring_out_type; + + ring_out_type tmp; + for (typename boost::range_iterator::type + it = ring.begin(); + it != ring.end(); + ++it) + { + geometry::traits::push_back::apply(tmp, *it); + } + collection.push_back(tmp); + } + + template + static inline void split_one_ring(RingType& ring, Collection& collection) + { + // we implement split_one_ring recursively (instead of using a + // stack or a queue) so that rings do not get copied + + typename boost::range_iterator::type pos1, pos2; + bool duplicate_found = find_duplicate_points + < + Closure + >::apply(ring, pos1, pos2); + + if (duplicate_found) + { + RingType other_ring; + ring.split_at(pos1, pos2, other_ring); + split_one_ring(ring, collection); + split_one_ring(other_ring, collection); + } + else + { + copy_to_collection(ring, collection); + } + } + +public: + template + static inline void apply(Ring const& ring, + RingCollection& collection, + RobustPolicy const& robust_policy) + { + typedef std::set > muu_turn_set; + typedef ring_as_dcl + < + typename point_type::type, closure::value + > ring_dcl_type; + + + // compute the ring's self turns + turns_container_type turns; + get_self_turns(ring, turns, robust_policy); + + // collect the ring's m:u/u turns; + // notice the use of std::set; we want to record coinciding + // m:u/u turns only once + muu_turn_set muu_turns; + for (typename turns_container_type::const_iterator it = turns.begin(); + it != turns.end(); + ++it) + { + if (it->method == method_touch_interior) + { + BOOST_GEOMETRY_ASSERT(it->operations[0].operation + == operation_union); + BOOST_GEOMETRY_ASSERT(it->operations[1].operation + == operation_union); + muu_turns.insert(*it); + } + } + + // insert the m:u/u turns as points in the original ring + ring_dcl_type output; + insert_muu_turns(ring, muu_turns, output); + + // split the ring into simple rings + split_one_ring::value>(output, collection); + } +}; + + +template +struct split_rings +{ + template + static inline void apply(RingCollection&, RobustPolicy const&) + { + } +}; + +// specialization for union +// TODO: add another specialization for intersection once implemented +template <> +struct split_rings +{ + template + static inline void apply(RingCollection& collection, + RobustPolicy const& robust_policy) + { + typedef typename boost::range_iterator + < + RingCollection + >::type ring_iterator_type; + + RingCollection new_collection; + for (ring_iterator_type rit = boost::begin(collection); + rit != boost::end(collection); + ++rit) + { + split_ring + < + overlay_union, + typename boost::range_value::type, + RobustPolicy + >::apply(*rit, new_collection, robust_policy); + } + collection.swap(new_collection); + } +}; + + +}} // namespace detail::overlay + +}} // namespace boost::geometry + +#endif // BOOST_GEOMETRY_ALGORITHMS_DETAIL_OVERLAY_SPLIT_RINGS_HPP From ac30dd488c5beb9e88da45eff4cce2a187922639 Mon Sep 17 00:00:00 2001 From: Menelaos Karavelas Date: Thu, 27 Aug 2015 17:22:13 +0300 Subject: [PATCH 02/22] [test][algorithms][union] add more test cases; update results for a few of the old ones; --- test/algorithms/overlay/multi_overlay_cases.hpp | 11 +++++ test/algorithms/overlay/overlay_cases.hpp | 51 ++++++++++++++++++++++ test/algorithms/set_operations/union/union.cpp | 16 ++++++- .../set_operations/union/union_multi.cpp | 10 ++++- 4 files changed, 85 insertions(+), 3 deletions(-) diff --git a/test/algorithms/overlay/multi_overlay_cases.hpp b/test/algorithms/overlay/multi_overlay_cases.hpp index c44b17adb5..5130904355 100644 --- a/test/algorithms/overlay/multi_overlay_cases.hpp +++ b/test/algorithms/overlay/multi_overlay_cases.hpp @@ -345,6 +345,17 @@ static std::string case_107_multi[2] = "MULTIPOLYGON(((5 7,6 8,6 10,7 9,8 10,8 8,7 8,6 7,6 6,5 7)))" }; +static std::string case_108_multi[2] = + { + "MULTIPOLYGON(((0 0,0 40,40 40,40 0,0 0),(10 10,30 10,30 30,10 30,10 10)))", + "MULTIPOLYGON(((10 10,10 20,20 10,10 10)),((20 10,30 20,30 10,20 10)),((10 20,10 30,20 20,10 20)),((20 20,30 30,30 20,20 20)))" + }; + +static std::string case_109_multi[2] = + { + "MULTIPOLYGON(((0 0,0 40,40 40,40 0,0 0),(10 10,30 10,30 30,10 30,10 10)))", + "MULTIPOLYGON(((15 10,10 15,10 17,15 10)),((15 10,10 20,10 22,15 10)),((15 10,10 25,10 27,15 10)),((25 10,30 17,30 15,25 10)),((25 10,30 22,30 20,25 10)),((25 10,30 27,30 25,25 10)),((18 10,20 30,19 10,18 10)),((21 10,20 30,22 10,21 10)))" + }; static std::string case_recursive_boxes_1[2] = { diff --git a/test/algorithms/overlay/overlay_cases.hpp b/test/algorithms/overlay/overlay_cases.hpp index 1a28e0c067..097d6a529d 100644 --- a/test/algorithms/overlay/overlay_cases.hpp +++ b/test/algorithms/overlay/overlay_cases.hpp @@ -180,6 +180,7 @@ static std::string case_35[2] = { "POLYGON((2 2,4 2,4 1,2 2))" }; static std::string case_36[2] = { + // union as one polygon with a hole touching the exterior ring "POLYGON((1 0,0 3,4 2,1 0))", "POLYGON((1 5,5 5,4 2,3 3,2 1,1 2,1 5))" }; @@ -353,6 +354,56 @@ static std::string case_79[2] = { "POLYGON((0 0,0 5,5 5,5 0,2 0,2 2,1 2,1 0,0 0))" }; +static std::string case_80[2] = + { + // union has one polygon with two holes; one of them is + // touching the exterior ring + // reported by MySQL QA on Aug 19, 2015 + "POLYGON((0 6,-11 -6,6 0,0 6),(3 1,5 0,-2 0,3 1))", + "POLYGON((5 4,6 0,9 12,-7 -12,5 -19,5 4))" + }; + +static std::string case_81[2] = + { + // union has a polygon with one hole touching the exterior ring + "POLYGON((0 0,10 10,20 0,0 0))", + "POLYGON((10 5,30 10,20 0,20 5,10 5))" + }; + +static std::string case_82[2] = + { + "POLYGON((0 0,10 10,20 0,0 0))", + "POLYGON((10 10,30 10,20 0,20 5,10 10))" + }; + +static std::string case_83[2] = + { + // union as a single polygon and two holes both touching the + // exterior ring at vertices + "POLYGON((0 0,10 10,20 0,0 0))", + "POLYGON((10 5,20 7,10 10,30 10,20 0,20 5,10 5))" + }; + +static std::string case_84[2] = + { + "POLYGON((0 0,10 10,20 0,0 0))", + "POLYGON((15 5,20 7,10 10,30 10,20 0,20 5,15 5))" + }; + +static std::string case_85[2] = + { + // union has a single polygon and two holes that touch each + // other at a vertex + "POLYGON((0 0,0 40,40 40,40 0,0 0),(10 10,30 10,30 30,10 30,10 10))", + "POLYGON((5 15,5 30,30 15,5 15))" + }; + +static std::string case_86[2] = + { + "POLYGON((0 0,0 40,40 40,40 0,20 0,0 0),(10 10,20 0,30 10,30 30,10 30,10 10))", + "POLYGON((10 10,10 30,30 30,30 10,10 10))" + }; + static std::string case_many_situations[2] = { "POLYGON((2 6,2 14,10 18,18 14,18 6,16 5,14 4,12 3,10 2,8 3,6 4,4 5,2 6))", "POLYGON((2 6,2 7,2 8,2 9,2 10,2 11,2 12,1 14" diff --git a/test/algorithms/set_operations/union/union.cpp b/test/algorithms/set_operations/union/union.cpp index 09c993db8a..b2b3622d21 100644 --- a/test/algorithms/set_operations/union/union.cpp +++ b/test/algorithms/set_operations/union/union.cpp @@ -169,7 +169,7 @@ void test_areal() test_one("33", case_33[0], case_33[1], 2, 0, 8, 4.5); test_one("36", - case_36[0], case_36[1], 1, 0, 10, 14.375); + case_36[0], case_36[1], 1, 1, 10, 14.375); test_one("40", case_40[0], case_40[1], 2, 0, 18, 11); @@ -181,6 +181,20 @@ void test_areal() test_one("59_iet", case_59[0], case_59[2], 1, 1, 14, 17.20833); + test_one("80", + case_80[0], case_80[1], 1, 2, 18, 221.369); + test_one("81", + case_81[0], case_81[1], 1, 1, 10, 147.5); + test_one("82", + case_82[0], case_82[1], 2, 0, 9, 175); + test_one("83", + case_83[0], case_83[1], 1, 2, 13, 172.917); + test_one("84", + case_84[0], case_84[1], 2, 0, 11, 170); + test_one("85", + case_85[0], case_85[1], 1, 2, 15, 1320); + test_one("86", + case_86[0], case_86[1], 1, 1, 10, 1500); /* test_one(102, diff --git a/test/algorithms/set_operations/union/union_multi.cpp b/test/algorithms/set_operations/union/union_multi.cpp index a80a249cae..9641f8b169 100644 --- a/test/algorithms/set_operations/union/union_multi.cpp +++ b/test/algorithms/set_operations/union/union_multi.cpp @@ -92,7 +92,7 @@ void test_areal() 1, 0, 13, 6); test_one("case_101_multi", case_101_multi[0], case_101_multi[1], - 1, 0, 32, 22.25); + 1, 3, 35, 22.25); test_one("case_103_multi", case_103_multi[0], case_103_multi[1], 1, 0, 7, 25); @@ -102,6 +102,12 @@ void test_areal() test_one("case_105_multi", case_105_multi[0], case_105_multi[1], 1, 0, 5, 25); + test_one("case_108_multi", + case_108_multi[0], case_108_multi[1], + 1, 2, 14, 1400); + test_one("case_109_multi", + case_109_multi[0], case_109_multi[1], + 1, 9, 45, 1250); test_one("case_recursive_boxes_1", case_recursive_boxes_1[0], case_recursive_boxes_1[1], @@ -111,7 +117,7 @@ void test_areal() 1, 0, 14, 100.0); // Area from SQL Server test_one("case_recursive_boxes_3", case_recursive_boxes_3[0], case_recursive_boxes_3[1], - 17, 0, 159, 56.5); // Area from SQL Server + 17, 6, 165, 56.5); // Area from SQL Server test_one("ggl_list_20120915_h2_a", ggl_list_20120915_h2[0], ggl_list_20120915_h2[1], From 19d2684efee4c9de58b6630b5bf18725e7dedd59 Mon Sep 17 00:00:00 2001 From: Menelaos Karavelas Date: Fri, 28 Aug 2015 13:05:52 +0300 Subject: [PATCH 03/22] [algorithm][overlay][union] in split_rings<>: implement the split_one_ring method in a non-recursive way using a stack; --- .../algorithms/detail/overlay/split_rings.hpp | 58 ++++++++++++++-------- 1 file changed, 37 insertions(+), 21 deletions(-) diff --git a/include/boost/geometry/algorithms/detail/overlay/split_rings.hpp b/include/boost/geometry/algorithms/detail/overlay/split_rings.hpp index 32520b37b2..988ae45341 100644 --- a/include/boost/geometry/algorithms/detail/overlay/split_rings.hpp +++ b/include/boost/geometry/algorithms/detail/overlay/split_rings.hpp @@ -14,6 +14,7 @@ #include #include #include +#include #include @@ -87,12 +88,16 @@ class ring_as_dcl } } - iterator begin() { return m_list.begin(); } - const_iterator begin() const { return m_list.begin(); } + inline iterator begin() { return m_list.begin(); } + inline const_iterator begin() const { return m_list.begin(); } - iterator end() { return m_list.end(); } - const_iterator end() const { return m_list.end(); } + inline iterator end() { return m_list.end(); } + inline const_iterator end() const { return m_list.end(); } + inline void swap(ring_as_dcl& other) + { + m_list.swap(other.m_list); + } private: List m_list; }; @@ -303,25 +308,36 @@ class split_ring template static inline void split_one_ring(RingType& ring, Collection& collection) { - // we implement split_one_ring recursively (instead of using a - // stack or a queue) so that rings do not get copied + // create and initialize a stack with the input ring + std::stack stack; + stack.push(RingType()); + stack.top().swap(ring); + + // while the stack is not empty: + // look for duplicates, split and push to stack; + // otherwise, copy to output collection + while (! stack.empty()) + { + RingType& top_ring = stack.top(); - typename boost::range_iterator::type pos1, pos2; - bool duplicate_found = find_duplicate_points - < - Closure - >::apply(ring, pos1, pos2); + typename boost::range_iterator::type pos1, pos2; + bool duplicate_found = find_duplicate_points + < + Closure + >::apply(top_ring, pos1, pos2); - if (duplicate_found) - { - RingType other_ring; - ring.split_at(pos1, pos2, other_ring); - split_one_ring(ring, collection); - split_one_ring(other_ring, collection); - } - else - { - copy_to_collection(ring, collection); + if (duplicate_found) + { + RingType other_ring; + stack.push(other_ring); + top_ring.split_at(pos1, pos2, other_ring); + stack.top().swap(other_ring); + } + else + { + copy_to_collection(top_ring, collection); + stack.pop(); + } } } From 9b50df119f52f552ffbaeef061f0da1dc7944bbb Mon Sep 17 00:00:00 2001 From: Menelaos Karavelas Date: Mon, 31 Aug 2015 00:40:38 +0300 Subject: [PATCH 04/22] [algorithms][overlay][union][split_rings] optimize the search for duplicate points by calling std::set<>::insert() instead of first std::set<>::find() and then possibly std::set<>::insert(); encapsulate moving an element to the top of a stack inside the move_to_top() method; --- .../algorithms/detail/overlay/split_rings.hpp | 27 +++++++++++++--------- 1 file changed, 16 insertions(+), 11 deletions(-) diff --git a/include/boost/geometry/algorithms/detail/overlay/split_rings.hpp b/include/boost/geometry/algorithms/detail/overlay/split_rings.hpp index 988ae45341..ce9bde4bfd 100644 --- a/include/boost/geometry/algorithms/detail/overlay/split_rings.hpp +++ b/include/boost/geometry/algorithms/detail/overlay/split_rings.hpp @@ -131,17 +131,15 @@ class find_duplicate_points iterator_type last = is_closed ? --boost::end(ring) : boost::end(ring); for (iterator_type it = boost::begin(ring); it != last; ++it) { - typename point_set_type::const_iterator pit = point_set.find(it); - if (pit != point_set.end()) + std::pair res + = point_set.insert(it); + + if (! res.second) { - pos1 = *pit; + pos1 = *res.first; pos2 = it; return true; } - else - { - point_set.insert(it); - } } // initialize pos1 and pos2 with something @@ -305,13 +303,21 @@ class split_ring collection.push_back(tmp); } + template + static inline void move_to_top(Stack& stack, + typename Stack::value_type& value) + { + typedef typename Stack::value_type value_type; + stack.push(value_type()); + stack.top().swap(value); + } + template static inline void split_one_ring(RingType& ring, Collection& collection) { // create and initialize a stack with the input ring std::stack stack; - stack.push(RingType()); - stack.top().swap(ring); + move_to_top(stack, ring); // while the stack is not empty: // look for duplicates, split and push to stack; @@ -329,9 +335,8 @@ class split_ring if (duplicate_found) { RingType other_ring; - stack.push(other_ring); top_ring.split_at(pos1, pos2, other_ring); - stack.top().swap(other_ring); + move_to_top(stack, other_ring); } else { From e542d0d9ce709fd3922d7e6d7a6e92638c1e5930 Mon Sep 17 00:00:00 2001 From: Menelaos Karavelas Date: Mon, 31 Aug 2015 01:21:50 +0300 Subject: [PATCH 05/22] [algorithms][overlay][union][split_rings] in find_duplicate_points, make the ring a Ring const&; --- .../boost/geometry/algorithms/detail/overlay/split_rings.hpp | 10 ++++++---- 1 file changed, 6 insertions(+), 4 deletions(-) diff --git a/include/boost/geometry/algorithms/detail/overlay/split_rings.hpp b/include/boost/geometry/algorithms/detail/overlay/split_rings.hpp index ce9bde4bfd..158600a143 100644 --- a/include/boost/geometry/algorithms/detail/overlay/split_rings.hpp +++ b/include/boost/geometry/algorithms/detail/overlay/split_rings.hpp @@ -122,14 +122,16 @@ class find_duplicate_points public: template - static inline bool apply(Ring& ring, Iterator& pos1, Iterator& pos2) + static inline bool apply(Ring const& ring, Iterator& pos1, Iterator& pos2) { typedef typename Ring::iterator iterator_type; typedef std::set point_set_type; point_set_type point_set; - iterator_type last = is_closed ? --boost::end(ring) : boost::end(ring); - for (iterator_type it = boost::begin(ring); it != last; ++it) + Ring& nc_ring = const_cast(ring); + iterator_type last + = is_closed ? --boost::end(nc_ring) : boost::end(nc_ring); + for (iterator_type it = boost::begin(nc_ring); it != last; ++it) { std::pair res = point_set.insert(it); @@ -143,7 +145,7 @@ class find_duplicate_points } // initialize pos1 and pos2 with something - pos1 = boost::begin(ring); + pos1 = boost::begin(nc_ring); pos2 = pos1; return false; } From ae36c7afaab05db627ea38447257e14ea36e68c4 Mon Sep 17 00:00:00 2001 From: Menelaos Karavelas Date: Mon, 31 Aug 2015 13:09:37 +0300 Subject: [PATCH 06/22] [algorithms][overlay][union][split_rings] avoid signed/unsigned integer comparison warning --- include/boost/geometry/algorithms/detail/overlay/split_rings.hpp | 5 +++-- 1 file changed, 3 insertions(+), 2 deletions(-) diff --git a/include/boost/geometry/algorithms/detail/overlay/split_rings.hpp b/include/boost/geometry/algorithms/detail/overlay/split_rings.hpp index 158600a143..5483611367 100644 --- a/include/boost/geometry/algorithms/detail/overlay/split_rings.hpp +++ b/include/boost/geometry/algorithms/detail/overlay/split_rings.hpp @@ -261,18 +261,19 @@ class split_ring MUU_Turns const& muu_turns, RingOut& ring_out) { + typedef typename boost::range_size::type size_type; typedef typename boost::range_iterator < MUU_Turns const >::type iterator_type; - typename boost::range_size::type point_index = 0; + size_type point_index = 0; for (iterator_type it = muu_turns.begin(); it != muu_turns.end(); ++it) { signed_size_type segment_index = get_correct_op(*it).seg_id.segment_index; - while (point_index <= segment_index) + while (point_index <= static_cast(segment_index)) { ring_out.push_back(ring[point_index]); ++point_index; From f8a0d9302f5bdf0c3fddf200cb0ba4e1c6e93505 Mon Sep 17 00:00:00 2001 From: Menelaos Karavelas Date: Sat, 5 Sep 2015 01:49:27 +0300 Subject: [PATCH 07/22] [algorithms][overlay][intersection] split rings in very much the same way as done for bg::union_ --- .../algorithms/detail/overlay/split_rings.hpp | 29 ++++++++++++++++++++-- 1 file changed, 27 insertions(+), 2 deletions(-) diff --git a/include/boost/geometry/algorithms/detail/overlay/split_rings.hpp b/include/boost/geometry/algorithms/detail/overlay/split_rings.hpp index 5483611367..20e7241374 100644 --- a/include/boost/geometry/algorithms/detail/overlay/split_rings.hpp +++ b/include/boost/geometry/algorithms/detail/overlay/split_rings.hpp @@ -156,8 +156,14 @@ template struct split_ring { template - static inline void apply(Ring const&, RingCollection&, RobustPolicy const&) + static inline void apply(Ring const& ring, + RingCollection& collection, + RobustPolicy const& robust_policy) { + split_ring + < + overlay_union, Ring, RobustPolicy + >::apply(ring, collection, robust_policy); } }; @@ -398,8 +404,27 @@ template struct split_rings { template - static inline void apply(RingCollection&, RobustPolicy const&) + static inline void apply(RingCollection& collection, + RobustPolicy const& robust_policy) { + typedef typename boost::range_iterator + < + RingCollection + >::type ring_iterator_type; + + RingCollection new_collection; + for (ring_iterator_type rit = boost::begin(collection); + rit != boost::end(collection); + ++rit) + { + split_ring + < + OverlayType, + typename boost::range_value::type, + RobustPolicy + >::apply(*rit, new_collection, robust_policy); + } + collection.swap(new_collection); } }; From 1f2bbedb28f4feb23ecb70254ee5e64987b5a72f Mon Sep 17 00:00:00 2001 From: Menelaos Karavelas Date: Sat, 5 Sep 2015 01:50:54 +0300 Subject: [PATCH 08/22] [test][algorithms][intersection][overlay] add more test cases for overlays; enhance testing of intersection(A,A) by allowing to check on the number of holes expected in the output; add more test cases for intersection(A,A); --- test/algorithms/overlay/overlay_cases.hpp | 18 +++++ .../set_operations/intersection/intersection.cpp | 14 ++++ .../intersection/test_intersection.hpp | 81 +++++++++++++++++++--- 3 files changed, 102 insertions(+), 11 deletions(-) diff --git a/test/algorithms/overlay/overlay_cases.hpp b/test/algorithms/overlay/overlay_cases.hpp index 097d6a529d..5d5af388fd 100644 --- a/test/algorithms/overlay/overlay_cases.hpp +++ b/test/algorithms/overlay/overlay_cases.hpp @@ -404,6 +404,24 @@ static std::string case_86[2] = "POLYGON((10 10,10 30,30 30,30 10,10 10))" }; +static std::string case_87[2] = + { + "POLYGON((0 5,-6 -17,12 17,0 5),(4 6,5 5,0 1,4 6))", + "POLYGON((3 9,-15 -5,13 -11,3 9))" + }; + +static std::string case_88[2] = + { + "POLYGON((5 6,-15 -13,1 -8,5 6))", + "POLYGON((0 8,-19 6,18 -17,20 8,11 17,0 8),(3 2,3 -1,1 0,3 2),(1 3,4 4,0 -1,1 3))" + }; + +static std::string case_89[2] = + { + "POLYGON((0 0,0 40,40 40,40 0,0 0),(10 10,20 19,20 20,10 10),(20 20,30 30,20 21,20 20))", + "POLYGON((10 10,10 30,30 30,30 10,10 10))" + }; + static std::string case_many_situations[2] = { "POLYGON((2 6,2 14,10 18,18 14,18 6,16 5,14 4,12 3,10 2,8 3,6 4,4 5,2 6))", "POLYGON((2 6,2 7,2 8,2 9,2 10,2 11,2 12,1 14" diff --git a/test/algorithms/set_operations/intersection/intersection.cpp b/test/algorithms/set_operations/intersection/intersection.cpp index 4fd54c9b0f..ddf171f459 100644 --- a/test/algorithms/set_operations/intersection/intersection.cpp +++ b/test/algorithms/set_operations/intersection/intersection.cpp @@ -324,6 +324,20 @@ void test_areal() test_one("buffer_mp2", buffer_mp2[0], buffer_mp2[1], 1, 29, 0.457126); + test_one_with_holes("case87", + case_87[0], case_87[1], + 1, 1, 11, 54.70134); + + test_one_with_holes("case88", + case_88[0], case_88[1], + 1, 1, 13, 35.933385); + +#ifdef BOOST_GEOMETRY_INCLUDE_FAILING_TESTS + test_one_with_holes("case89", + case_89[0], case_89[1], + 2, 0, 13, 390); +#endif // BOOST_GEOMETRY_INCLUDE_FAILING_TESTS + return; diff --git a/test/algorithms/set_operations/intersection/test_intersection.hpp b/test/algorithms/set_operations/intersection/test_intersection.hpp index 1c45032a98..2dd8fb644e 100644 --- a/test/algorithms/set_operations/intersection/test_intersection.hpp +++ b/test/algorithms/set_operations/intersection/test_intersection.hpp @@ -39,7 +39,9 @@ typename bg::default_area_result::type check_result( std::vector const& intersection_output, std::string const& caseid, - std::size_t expected_count = 0, int expected_point_count = 0, + std::size_t expected_count = 0, + int expected_hole_count = -1, + int expected_point_count = 0, double expected_length_or_area = 0, double percentage = 0.0001, bool debug = false) @@ -48,6 +50,7 @@ check_result( typename bg::default_area_result::type length_or_area = 0; int n = 0; + int holes = 0; for (typename std::vector::const_iterator it = intersection_output.begin(); it != intersection_output.end(); ++it) @@ -56,6 +59,7 @@ check_result( { // here n should rather be of type std::size_t, but expected_point_count // is set to -1 in some test cases so type int was left for now + holes += static_cast(bg::num_interior_rings(*it)); n += static_cast(bg::num_points(*it, true)); } @@ -72,6 +76,14 @@ check_result( #if ! defined(BOOST_GEOMETRY_NO_BOOST_TEST) #if ! defined(BOOST_GEOMETRY_NO_ROBUSTNESS) + + BOOST_CHECK_MESSAGE(expected_hole_count < 0 || holes == expected_hole_count, + "intersection: " << caseid + << " #holes expected: " << expected_hole_count + << " detected: " << holes + << " type: " << (type_for_assert_message()) + ); + if (expected_point_count > 0) { BOOST_CHECK_MESSAGE(bg::math::abs(n - expected_point_count) < 3, @@ -111,9 +123,11 @@ check_result( template -typename bg::default_area_result::type test_intersection(std::string const& caseid, +typename bg::default_area_result::type test_intersection_with_holes(std::string const& caseid, G1 const& g1, G2 const& g2, - std::size_t expected_count = 0, int expected_point_count = 0, + std::size_t expected_count = 0, + int expected_hole_count = -1, + int expected_point_count = 0, double expected_length_or_area = 0, double percentage = 0.0001, bool debug = false) @@ -147,26 +161,30 @@ typename bg::default_area_result::type test_intersection(std::string const& std::vector intersection_output; bg::intersection(g1, g2, intersection_output); - check_result(intersection_output, caseid, expected_count, expected_point_count, + check_result(intersection_output, caseid, expected_count, + expected_hole_count, expected_point_count, expected_length_or_area, percentage, debug); // Check variant behaviour intersection_output.clear(); bg::intersection(boost::variant(g1), g2, intersection_output); - check_result(intersection_output, caseid, expected_count, expected_point_count, + check_result(intersection_output, caseid, expected_count, + expected_hole_count, expected_point_count, expected_length_or_area, percentage, debug); intersection_output.clear(); bg::intersection(g1, boost::variant(g2), intersection_output); - check_result(intersection_output, caseid, expected_count, expected_point_count, + check_result(intersection_output, caseid, expected_count, + expected_hole_count, expected_point_count, expected_length_or_area, percentage, debug); intersection_output.clear(); bg::intersection(boost::variant(g1), boost::variant(g2), intersection_output); - check_result(intersection_output, caseid, expected_count, expected_point_count, + check_result(intersection_output, caseid, expected_count, + expected_hole_count, expected_point_count, expected_length_or_area, percentage, debug); #if defined(TEST_WITH_SVG) @@ -225,10 +243,29 @@ typename bg::default_area_result::type test_intersection(std::string const& return length_or_area; } +template +typename bg::default_area_result::type test_intersection(std::string const& caseid, + G1 const& g1, G2 const& g2, + std::size_t expected_count = 0, + int expected_point_count = 0, + double expected_length_or_area = 0, + double percentage = 0.0001, + bool debug = false) +{ + return test_intersection_with_holes + < + OutputType, CalculationType, G1, G2 + >(caseid, g1, g2, expected_count, -1, expected_point_count, + expected_length_or_area, percentage, debug); +} + + template -typename bg::default_area_result::type test_one(std::string const& caseid, +typename bg::default_area_result::type test_one_with_holes(std::string const& caseid, std::string const& wkt1, std::string const& wkt2, - std::size_t expected_count = 0, int expected_point_count = 0, + std::size_t expected_count = 0, + int expected_hole_count = -1, + int expected_point_count = 0, double expected_length_or_area = 0, double percentage = 0.0001, bool debug = false) @@ -243,12 +280,28 @@ typename bg::default_area_result::type test_one(std::string const& caseid, bg::correct(g1); bg::correct(g2); - return test_intersection(caseid, g1, g2, - expected_count, expected_point_count, + return test_intersection_with_holes(caseid, g1, g2, + expected_count, expected_hole_count, expected_point_count, expected_length_or_area, percentage, debug); } +template +typename bg::default_area_result::type test_one(std::string const& caseid, + std::string const& wkt1, std::string const& wkt2, + std::size_t expected_count = 0, + int expected_point_count = 0, + double expected_length_or_area = 0, + double percentage = 0.0001, + bool debug = false) +{ + return test_one_with_holes + < + OutputType, G1, G2 + >(caseid, wkt1, wkt2, expected_count, -1, expected_point_count, + expected_length_or_area, percentage, debug); +} + template void test_one_lp(std::string const& caseid, std::string const& wkt_areal, std::string const& wkt_linear, @@ -257,6 +310,9 @@ void test_one_lp(std::string const& caseid, double percentage = 0.0001, bool debug1 = false, bool debug2 = false) { +#ifdef BOOST_GEOMETRY_TEST_DEBUG + std::cout << caseid << " -- start" << std::endl; +#endif Areal areal; bg::read_wkt(wkt_areal, areal); bg::correct(areal); @@ -274,6 +330,9 @@ void test_one_lp(std::string const& caseid, test_intersection(caseid + "_rev", areal, linear, expected_count, expected_point_count, expected_length, percentage, debug2); +#ifdef BOOST_GEOMETRY_TEST_DEBUG + std::cout << caseid << " -- end" << std::endl; +#endif } template From 4ae6e4935081fe1954ba568121c7e8072f315f89 Mon Sep 17 00:00:00 2001 From: Menelaos Karavelas Date: Sun, 6 Sep 2015 02:05:44 +0300 Subject: [PATCH 09/22] [test][algorithms][difference] update expected results (due to changes related to the splitting of rings produced during traversal) --- test/algorithms/set_operations/difference/difference_multi.cpp | 2 +- .../set_operations/difference/difference_multi_spike.cpp | 8 ++++---- 2 files changed, 5 insertions(+), 5 deletions(-) diff --git a/test/algorithms/set_operations/difference/difference_multi.cpp b/test/algorithms/set_operations/difference/difference_multi.cpp index 09ba7f74e3..335a086708 100644 --- a/test/algorithms/set_operations/difference/difference_multi.cpp +++ b/test/algorithms/set_operations/difference/difference_multi.cpp @@ -136,7 +136,7 @@ void test_areal() test_one("ggl_list_20120221_volker", ggl_list_20120221_volker[0], ggl_list_20120221_volker[1], - 2, 12, 7962.66, 1, 18, 2775258.93, + 2, 12, 7962.66, 2, 18, 2775258.93, 0.001); #if ! defined(BOOST_GEOMETRY_NO_ROBUSTNESS) diff --git a/test/algorithms/set_operations/difference/difference_multi_spike.cpp b/test/algorithms/set_operations/difference/difference_multi_spike.cpp index 7ec3d2ce0a..0fbf3dcf60 100644 --- a/test/algorithms/set_operations/difference/difference_multi_spike.cpp +++ b/test/algorithms/set_operations/difference/difference_multi_spike.cpp @@ -38,26 +38,26 @@ void test_spikes_in_ticket_8364() test_one("ticket_8364_step3", "MULTIPOLYGON(((3232 2532,2136 2790,1032 1764,1032 1458,1032 1212,2136 2328,3232 2220,3232 1056,1031 1056,1031 2856,3232 2856,3232 2532)))", "MULTIPOLYGON(((1032 2130,2052 2712,1032 1764,1032 2130)),((3234 2580,3234 2532,2558 2690,3234 2580)),((2558 2690,2136 2760,2052 2712,2136 2790,2558 2690)))", - 2, + if_typed(2, 3), if_typed(19, 22), if_typed(2775595.5, 2775256.487954), // SQL Server: 2775256.47588724 3, -1, // don't check point-count if_typed(7893.0, 7810.487954), // SQL Server: 7810.48711165739 - if_typed(1, 5), + if_typed(1, 6), -1, if_typed(2783349.5, 2775256.487954 + 7810.487954)); test_one("ticket_8364_step4", "MULTIPOLYGON(((2567 2688,2136 2790,2052 2712,1032 2130,1032 1764,1032 1458,1032 1212,2136 2328,3232 2220,3232 1056,1031 1056,1031 2856,3232 2856,3232 2580,2567 2688)))", "MULTIPOLYGON(((1032 2556,1778 2556,1032 2130,1032 2556)),((3234 2580,3234 2556,1778 2556,2136 2760,3234 2580)))", - 1, + if_typed(1, 2), if_typed(17, 20), if_typed(2615783.5, 2616029.559567), // SQL Server: 2616029.55616044 1, if_typed(9, 11), if_typed(161133.5, 161054.559567), // SQL Server: 161054.560110092 - if_typed(1, 2), + if_typed(1, 3), if_typed(25, 31), if_typed(2776875.5, 2616029.559567 + 161054.559567)); } From 21779bccc4eafb1d7a449856d0d70089b923e009 Mon Sep 17 00:00:00 2001 From: Menelaos Karavelas Date: Sun, 6 Sep 2015 02:06:47 +0300 Subject: [PATCH 10/22] [test][algorithms][intersection] update expected results (due to changes related to the splitting of rings produced during traversal) --- .../set_operations/intersection/intersection_multi.cpp | 11 +++++++---- 1 file changed, 7 insertions(+), 4 deletions(-) diff --git a/test/algorithms/set_operations/intersection/intersection_multi.cpp b/test/algorithms/set_operations/intersection/intersection_multi.cpp index 2399935cb3..f20f46ca77 100644 --- a/test/algorithms/set_operations/intersection/intersection_multi.cpp +++ b/test/algorithms/set_operations/intersection/intersection_multi.cpp @@ -93,17 +93,20 @@ void test_areal() test_one("case_recursive_boxes_1", case_recursive_boxes_1[0], case_recursive_boxes_1[1], 10, 97, 47.0); - test_one("case_recursive_boxes_2", + + test_one_with_holes( + "case_recursive_boxes_2", case_recursive_boxes_2[0], case_recursive_boxes_2[1], - 1, 47, 90.0); // Area from SQL Server + 1, 6, 52, 90.0); // Area from SQL Server test_one("case_recursive_boxes_3", case_recursive_boxes_3[0], case_recursive_boxes_3[1], 19, 87, 12.5); // Area from SQL Server - test_one("case_recursive_boxes_4", + test_one_with_holes( + "case_recursive_boxes_4", case_recursive_boxes_4[0], case_recursive_boxes_4[1], - 13, 157, 67.0); // Area from SQL Server + 13, 7, 162, 67.0); // Area from SQL Server test_one("ggl_list_20120915_h2_a", ggl_list_20120915_h2[0], ggl_list_20120915_h2[1], From 9550d8b0d1e9ee79301faae069e72398c65c3c5a Mon Sep 17 00:00:00 2001 From: Menelaos Karavelas Date: Tue, 22 Sep 2015 09:22:22 +0300 Subject: [PATCH 11/22] [test][algorithms][overlay] add more cases for set operations --- test/algorithms/overlay/overlay_cases.hpp | 43 +++++++++++++++++++++++++++++++ 1 file changed, 43 insertions(+) diff --git a/test/algorithms/overlay/overlay_cases.hpp b/test/algorithms/overlay/overlay_cases.hpp index 5d5af388fd..a255ea381c 100644 --- a/test/algorithms/overlay/overlay_cases.hpp +++ b/test/algorithms/overlay/overlay_cases.hpp @@ -422,6 +422,49 @@ static std::string case_89[2] = "POLYGON((10 10,10 30,30 30,30 10,10 10))" }; + +static std::string case_90[2] = + { + "POLYGON((1 8,-17 -13,19 -9,1 8))", + "POLYGON((8 6,5 7,-1 4,-8 -7,0 -17,8 6),(3 6,5 5,0 -2,3 6))" + }; + +static std::string case_91[2] = + { + "POLYGON((0 0,50 0,50 50,25 50,0 50,0 0),(25 50,40 40,10 40,25 50))", + "POLYGON((20 20,20 60,30 60,30 20,20 20))" + }; + +static std::string case_92[2] = + { + "POLYGON((0 0,50 0,50 50,25 50,0 50,0 0),(25 50,40 40,10 40,25 50))", + "POLYGON((5 20,5 60,30 60,30 20,5 20))" + }; + +static std::string case_93[2] = + { + "POLYGON((0 0,50 0,50 50,24 50,26 50,0 50,0 0),(24 50,24 40,10 40,24 50),(26 50,40 40,26 40,26 50))", + "POLYGON((20 20,20 60,30 60,30 20,20 20))" + }; + +static std::string case_94[2] = + { + "POLYGON((0 0,50 0,50 50,24 50,26 50,0 50,0 0),(24 50,24 40,10 40,24 50),(26 50,40 40,26 40,26 50))", + "POLYGON((5 20,5 60,30 60,30 20,5 20))" + }; + +static std::string case_95[2] = + { + "POLYGON((0 0,50 0,50 50,25 50,0 50,0 0),(25 50,24 40,10 40,25 50),(25 50,40 40,26 40,25 50))", + "POLYGON((20 20,20 60,30 60,30 20,20 20))" + }; + +static std::string case_96[2] = + { + "POLYGON((0 0,50 0,50 50,25 50,0 50,0 0),(24 50,25 40,10 40,24 50),(26 50,40 40,25 40,26 50))", + "POLYGON((20 20,20 60,30 60,30 20,20 20))" + }; + static std::string case_many_situations[2] = { "POLYGON((2 6,2 14,10 18,18 14,18 6,16 5,14 4,12 3,10 2,8 3,6 4,4 5,2 6))", "POLYGON((2 6,2 7,2 8,2 9,2 10,2 11,2 12,1 14" From 7d8641d5aaab3f1f129e48b6d7d324b97982be4e Mon Sep 17 00:00:00 2001 From: Menelaos Karavelas Date: Tue, 22 Sep 2015 14:41:58 +0300 Subject: [PATCH 12/22] [test][algorithms][intersection] add more test cases; add function to test the validity of the intersection outcome; --- .../set_operations/intersection/intersection.cpp | 35 ++++++++++++++++++++++ .../intersection/test_intersection.hpp | 15 ++++++++++ 2 files changed, 50 insertions(+) diff --git a/test/algorithms/set_operations/intersection/intersection.cpp b/test/algorithms/set_operations/intersection/intersection.cpp index ddf171f459..b13e2335bc 100644 --- a/test/algorithms/set_operations/intersection/intersection.cpp +++ b/test/algorithms/set_operations/intersection/intersection.cpp @@ -338,6 +338,41 @@ void test_areal() 2, 0, 13, 390); #endif // BOOST_GEOMETRY_INCLUDE_FAILING_TESTS + test_one_with_holes("case90", + case_90[0], case_90[1], + 2, 0, 13, 143.067882); + test_validity(case_90[0], case_90[1]); + + test_one_with_holes("case91", + case_91[0], case_91[1], + 3, 0, 13, 216.6667); + test_validity(case_91[0], case_91[1]); + + test_one_with_holes("case92", + case_92[0], case_92[1], + 2, 0, 13, 633.3333); + test_validity(case_92[0], case_92[1]); + + test_one_with_holes("case93", + case_93[0], case_93[1], + 3, 0, 17, 231.42857); + test_validity(case_93[0], case_93[1]); + + test_one_with_holes("case94", + case_94[0], case_94[1], + 2, 1, 16, 645.71429); + test_validity(case_94[0], case_94[1]); + + test_one_with_holes("case95", + case_95[0], case_95[1], + 3, 0, 16, 226.66667); + test_validity(case_95[0], case_95[1]); + + test_one_with_holes("case96", + case_96[0], case_96[1], + 4, 0, 19, 221.42857); + test_validity(case_96[0], case_96[1]); + return; diff --git a/test/algorithms/set_operations/intersection/test_intersection.hpp b/test/algorithms/set_operations/intersection/test_intersection.hpp index 2dd8fb644e..c2a2943385 100644 --- a/test/algorithms/set_operations/intersection/test_intersection.hpp +++ b/test/algorithms/set_operations/intersection/test_intersection.hpp @@ -18,6 +18,7 @@ #include #include #include +#include #include #include @@ -352,4 +353,18 @@ void test_point_output(std::string const& wkt1, std::string const& wkt2, unsigne } +template +inline void test_validity(std::string const& wkt1, + std::string const& wkt2) +{ + Areal1 a1; + Areal2 a2; + bg::read_wkt(wkt1, a1); + bg::read_wkt(wkt2, a2); + bg::model::multi_polygon out; + bg::intersection(a1, a2, out); + BOOST_CHECK(bg::is_valid(out)); +} + + #endif From a7e5f8cf795e8941a7882e7241bea934cb4151c5 Mon Sep 17 00:00:00 2001 From: Menelaos Karavelas Date: Tue, 22 Sep 2015 14:42:57 +0300 Subject: [PATCH 13/22] [algorithms][overlay][split_rings] rings produced during the traversal may have not only m:u/u turns but also m:i/i turns (induced by self-intersection points belonging to different polygons in the set of output rings); insert not only m:u/u turns but also m:i/i turns present in the rings to be split; --- .../algorithms/detail/overlay/split_rings.hpp | 57 +++++++++++----------- 1 file changed, 29 insertions(+), 28 deletions(-) diff --git a/include/boost/geometry/algorithms/detail/overlay/split_rings.hpp b/include/boost/geometry/algorithms/detail/overlay/split_rings.hpp index 20e7241374..671cc5d33e 100644 --- a/include/boost/geometry/algorithms/detail/overlay/split_rings.hpp +++ b/include/boost/geometry/algorithms/detail/overlay/split_rings.hpp @@ -209,25 +209,23 @@ class split_ring ; } - template - struct muu_turn_less + template + struct maa_turn_less { - bool operator()(MUU_Turn const& t1, MUU_Turn const& t2) const + bool operator()(MAA_Turn const& t1, MAA_Turn const& t2) const { BOOST_GEOMETRY_ASSERT(t1.method == method_touch_interior); - BOOST_GEOMETRY_ASSERT(t1.operations[0].operation - == operation_union); - BOOST_GEOMETRY_ASSERT(t1.operations[1].operation - == operation_union); + BOOST_GEOMETRY_ASSERT(t1.both(operation_union) + || + t1.both(operation_intersection)); BOOST_GEOMETRY_ASSERT(t2.method == method_touch_interior); - BOOST_GEOMETRY_ASSERT(t2.operations[0].operation - == operation_union); - BOOST_GEOMETRY_ASSERT(t2.operations[1].operation - == operation_union); + BOOST_GEOMETRY_ASSERT(t2.both(operation_union) + || + t2.both(operation_intersection)); - typename MUU_Turn::turn_operation_type op1 = get_correct_op(t1); - typename MUU_Turn::turn_operation_type op2 = get_correct_op(t2); + typename MAA_Turn::turn_operation_type op1 = get_correct_op(t1); + typename MAA_Turn::turn_operation_type op2 = get_correct_op(t2); BOOST_GEOMETRY_ASSERT(! op1.fraction.is_zero() && ! op1.fraction.is_one()); @@ -262,19 +260,19 @@ class split_ring get_self_turns(ring, turns, robust_policy, no_interrupt_policy()); } - template - static inline void insert_muu_turns(Ring const& ring, - MUU_Turns const& muu_turns, + template + static inline void insert_maa_turns(Ring const& ring, + MAA_Turns const& maa_turns, RingOut& ring_out) { typedef typename boost::range_size::type size_type; typedef typename boost::range_iterator < - MUU_Turns const + MAA_Turns const >::type iterator_type; size_type point_index = 0; - for (iterator_type it = muu_turns.begin(); it != muu_turns.end(); ++it) + for (iterator_type it = maa_turns.begin(); it != maa_turns.end(); ++it) { signed_size_type segment_index = get_correct_op(*it).seg_id.segment_index; @@ -361,7 +359,7 @@ class split_ring RingCollection& collection, RobustPolicy const& robust_policy) { - typedef std::set > muu_turn_set; + typedef std::set > maa_turn_set; typedef ring_as_dcl < typename point_type::type, closure::value @@ -372,27 +370,30 @@ class split_ring turns_container_type turns; get_self_turns(ring, turns, robust_policy); - // collect the ring's m:u/u turns; + // collect the ring's m:u/u and m:i/i turns (the latter can + // appear when we perform an intersection and the intersection + // result consists of a multipolygon whose polygons touch each + // other); // notice the use of std::set; we want to record coinciding - // m:u/u turns only once - muu_turn_set muu_turns; + // m:u/u and m:i/i turns only once + maa_turn_set maa_turns; for (typename turns_container_type::const_iterator it = turns.begin(); it != turns.end(); ++it) { if (it->method == method_touch_interior) { - BOOST_GEOMETRY_ASSERT(it->operations[0].operation - == operation_union); - BOOST_GEOMETRY_ASSERT(it->operations[1].operation - == operation_union); - muu_turns.insert(*it); + BOOST_GEOMETRY_ASSERT(it->both(operation_union) + || + it->both(operation_intersection)); + + maa_turns.insert(*it); } } // insert the m:u/u turns as points in the original ring ring_dcl_type output; - insert_muu_turns(ring, muu_turns, output); + insert_maa_turns(ring, maa_turns, output); // split the ring into simple rings split_one_ring::value>(output, collection); From 70103bd6258e007b5f7f660a5cbe096098877dab Mon Sep 17 00:00:00 2001 From: Menelaos Karavelas Date: Mon, 5 Oct 2015 11:42:00 +0300 Subject: [PATCH 14/22] [algorithms][buffer] the rings produced during the union operations performed during the buffer computation may have self-intersection points; these rings need to be split in order to produce valid buffers; fix: apply the split_rings routine to the rings produced as above; --- .../geometry/algorithms/detail/buffer/buffered_piece_collection.hpp | 6 ++++++ 1 file changed, 6 insertions(+) diff --git a/include/boost/geometry/algorithms/detail/buffer/buffered_piece_collection.hpp b/include/boost/geometry/algorithms/detail/buffer/buffered_piece_collection.hpp index 545d89cb9b..7519cb56df 100644 --- a/include/boost/geometry/algorithms/detail/buffer/buffered_piece_collection.hpp +++ b/include/boost/geometry/algorithms/detail/buffer/buffered_piece_collection.hpp @@ -41,6 +41,7 @@ #include #include #include +#include #include #include #include @@ -1331,6 +1332,11 @@ struct buffered_piece_collection traverser::apply(offsetted_rings, offsetted_rings, detail::overlay::operation_union, m_robust_policy, m_turns, traversed_rings); + + overlay::split_rings + < + overlay_union + >::apply(traversed_rings, m_robust_policy); } inline void reverse() From 2c7d4486b8ecc9932dc737b4ac1844f061d194b5 Mon Sep 17 00:00:00 2001 From: Menelaos Karavelas Date: Wed, 7 Oct 2015 10:38:27 +0300 Subject: [PATCH 15/22] [algorithms][overlay] instead of checking an assertion, through a proper exception instead; the previous behavior can be turned on by defining the appropriate macro; --- .../algorithms/detail/overlay/split_rings.hpp | 28 ++++++++++++++++++++-- 1 file changed, 26 insertions(+), 2 deletions(-) diff --git a/include/boost/geometry/algorithms/detail/overlay/split_rings.hpp b/include/boost/geometry/algorithms/detail/overlay/split_rings.hpp index 671cc5d33e..5833665ba7 100644 --- a/include/boost/geometry/algorithms/detail/overlay/split_rings.hpp +++ b/include/boost/geometry/algorithms/detail/overlay/split_rings.hpp @@ -30,6 +30,7 @@ #include #include +#include #include #include #include @@ -214,15 +215,30 @@ class split_ring { bool operator()(MAA_Turn const& t1, MAA_Turn const& t2) const { +#if ! defined(BOOST_GEOMETRY_OVERLAY_NO_THROW) + if (t1.method != method_touch_interior + || + (! t1.both(operation_union) + && ! t1.both(operation_intersection)) + || + t2.method != method_touch_interior + || + (! t2.both(operation_union) + && ! t2.both(operation_intersection)) + ) + { + throw inconsistent_turns_exception(); + } +#else BOOST_GEOMETRY_ASSERT(t1.method == method_touch_interior); BOOST_GEOMETRY_ASSERT(t1.both(operation_union) || t1.both(operation_intersection)); - BOOST_GEOMETRY_ASSERT(t2.method == method_touch_interior); BOOST_GEOMETRY_ASSERT(t2.both(operation_union) || t2.both(operation_intersection)); +#endif typename MAA_Turn::turn_operation_type op1 = get_correct_op(t1); typename MAA_Turn::turn_operation_type op2 = get_correct_op(t2); @@ -383,10 +399,18 @@ class split_ring { if (it->method == method_touch_interior) { +#if ! defined(BOOST_GEOMETRY_OVERLAY_NO_THROW) + if (! it->both(operation_union) + && + ! it->both(operation_intersection)) + { + throw inconsistent_turns_exception(); + } +#else BOOST_GEOMETRY_ASSERT(it->both(operation_union) || it->both(operation_intersection)); - +#endif maa_turns.insert(*it); } } From 477aef136578abbd73bb9b618eea8f00d5839bd0 Mon Sep 17 00:00:00 2001 From: Menelaos Karavelas Date: Wed, 14 Oct 2015 09:23:08 +0300 Subject: [PATCH 16/22] [test][algorithms][union] update test results --- test/algorithms/set_operations/union/union_multi.cpp | 4 ++-- 1 file changed, 2 insertions(+), 2 deletions(-) diff --git a/test/algorithms/set_operations/union/union_multi.cpp b/test/algorithms/set_operations/union/union_multi.cpp index 80064bc925..38ea938d92 100644 --- a/test/algorithms/set_operations/union/union_multi.cpp +++ b/test/algorithms/set_operations/union/union_multi.cpp @@ -121,10 +121,10 @@ void test_areal() test_one("case_recursive_boxes_4", case_recursive_boxes_4[0], case_recursive_boxes_4[1], - 1, 1, 42, 96.75); + 1, 2, 43, 96.75); test_one("case_recursive_boxes_5", case_recursive_boxes_5[0], case_recursive_boxes_5[1], - 3, 2, 110, 70.0); + 3, 10, 118, 70.0); test_one("ggl_list_20120915_h2_a", ggl_list_20120915_h2[0], ggl_list_20120915_h2[1], From d97a38ee578948666166cf1c77347af46f42ce44 Mon Sep 17 00:00:00 2001 From: Menelaos Karavelas Date: Fri, 16 Oct 2015 12:09:34 +0300 Subject: [PATCH 17/22] [test][algorithms][buffer] add test case regarding the buffer of a multilinestring --- test/algorithms/buffer/buffer_multi_linestring.cpp | 2 ++ 1 file changed, 2 insertions(+) diff --git a/test/algorithms/buffer/buffer_multi_linestring.cpp b/test/algorithms/buffer/buffer_multi_linestring.cpp index 6501d131ad..20fcb18d85 100644 --- a/test/algorithms/buffer/buffer_multi_linestring.cpp +++ b/test/algorithms/buffer/buffer_multi_linestring.cpp @@ -139,6 +139,8 @@ void test_all() mysql_2015_09_08b, join_round32, end_round32, 1.32832149026508268e+19, 2061380362.0, same_distance, true, 1.0e12); + + test_one("mysql_report_2015_09_21", "MULTILINESTRING((-5 15,7 15,19 -10,-11 -2),(2 13,2 -9))", join_round32, end_round32, 186.550431076137, 1.0); } From 24b09b110a7288b478a9d3441eb61fce7d059353 Mon Sep 17 00:00:00 2001 From: Menelaos Karavelas Date: Fri, 16 Oct 2015 12:10:33 +0300 Subject: [PATCH 18/22] [test][algorithms][overlay] add some test cases --- test/algorithms/overlay/multi_overlay_cases.hpp | 6 ++++ test/algorithms/overlay/overlay_cases.hpp | 42 +++++++++++++++++++++++++ 2 files changed, 48 insertions(+) diff --git a/test/algorithms/overlay/multi_overlay_cases.hpp b/test/algorithms/overlay/multi_overlay_cases.hpp index 90787a383c..e040338dfa 100644 --- a/test/algorithms/overlay/multi_overlay_cases.hpp +++ b/test/algorithms/overlay/multi_overlay_cases.hpp @@ -357,6 +357,12 @@ static std::string case_109_multi[2] = "MULTIPOLYGON(((15 10,10 15,10 17,15 10)),((15 10,10 20,10 22,15 10)),((15 10,10 25,10 27,15 10)),((25 10,30 17,30 15,25 10)),((25 10,30 22,30 20,25 10)),((25 10,30 27,30 25,25 10)),((18 10,20 30,19 10,18 10)),((21 10,20 30,22 10,21 10)))" }; +static std::string case_110_multi[2] = + { + "MULTIPOLYGON(((4 5,12 11,-12 -3,4 5)))", + "MULTIPOLYGON(((5 4,-14 0,1 0,5 4)),((1 6,13 0,10 12,1 6)))" + }; + static std::string case_recursive_boxes_1[2] = { // == 70 diff --git a/test/algorithms/overlay/overlay_cases.hpp b/test/algorithms/overlay/overlay_cases.hpp index a255ea381c..a82a2165ed 100644 --- a/test/algorithms/overlay/overlay_cases.hpp +++ b/test/algorithms/overlay/overlay_cases.hpp @@ -465,6 +465,48 @@ static std::string case_96[2] = "POLYGON((20 20,20 60,30 60,30 20,20 20))" }; +static std::string case_97[2] = + { + "POLYGON((6 7,18 14,-8 1,0 0,18 -8,6 7),(6 0,-4 3,5 3,6 0))", + "POLYGON((2 3,-3 5,-10 -1,2 3))" + }; + +static std::string case_98[2] = + { + "POLYGON((6 7,18 14,-8 1,0 0,18 -8,6 7),(6 0,-4 3,5 3,6 0))", + "POLYGON((0 7,-5 6,11 -13,0 7))" + }; + +static std::string case_99[2] = + { + "POLYGON((1 1,1 2398046010000,30000 2398046010000,30000 1,1 1))", + "POLYGON((8000 10000,8000 12000,120115188075850000 12000,120115188075850000 10000,8000 10000))" + }; + +static std::string case_100[2] = + { + "POLYGON((-28011 -32449,-28011 4652.2100,2667 4652.2100,2667 -32449,-28011 -32449))", + "POLYGON((-8560 68,-8560 2305843009213693948,-2087 2305843009213693948,-2087 68,-8560 68))" + }; + +static std::string case_101[2] = + { + "POLYGON((2 5514.6191,2 28951,17237 28951,17237 5514.6191,2 5514.6191))", + "POLYGON((8128 -20054,8128 144115188075855874,11192 144115188075855874,11192 -20054,8128 -20054))" + }; + +static std::string case_102[2] = + { + "POLYGON((-15032 -2251799813685244,-15032 8388611,16781 8388611,16781 -2251799813685244,-15032 -2251799813685244))", + "POLYGON((-21087 20243,-21087 31851,18514 31851,18514 20243,-21087 20243))" + }; + +static std::string case_103[2] = + { + "POLYGON((8 8,-14 1,-6 3,8 8))", + "POLYGON((6 6,3 14,-9 7,6 6))" + }; + static std::string case_many_situations[2] = { "POLYGON((2 6,2 14,10 18,18 14,18 6,16 5,14 4,12 3,10 2,8 3,6 4,4 5,2 6))", "POLYGON((2 6,2 7,2 8,2 9,2 10,2 11,2 12,1 14" From 5d54fdccb8231061f7617334060bfe81c7ffd17d Mon Sep 17 00:00:00 2001 From: Menelaos Karavelas Date: Fri, 16 Oct 2015 12:13:05 +0300 Subject: [PATCH 19/22] [test][algorithms][intersection] add testing of additional test cases; update some test results; update code regarding the validity check of the intersection output; --- .../set_operations/intersection/intersection.cpp | 23 +++++++++++++++------- .../intersection/intersection_multi.cpp | 9 ++++++++- .../intersection/test_intersection.hpp | 16 ++++++++++++--- 3 files changed, 37 insertions(+), 11 deletions(-) diff --git a/test/algorithms/set_operations/intersection/intersection.cpp b/test/algorithms/set_operations/intersection/intersection.cpp index e21b774e72..d805384b69 100644 --- a/test/algorithms/set_operations/intersection/intersection.cpp +++ b/test/algorithms/set_operations/intersection/intersection.cpp @@ -341,38 +341,47 @@ void test_areal() test_one_with_holes("case90", case_90[0], case_90[1], 2, 0, 13, 143.067882); - test_validity(case_90[0], case_90[1]); + test_validity("case90", case_90[0], case_90[1]); test_one_with_holes("case91", case_91[0], case_91[1], 3, 0, 13, 216.6667); - test_validity(case_91[0], case_91[1]); + test_validity("case91", case_91[0], case_91[1]); test_one_with_holes("case92", case_92[0], case_92[1], 2, 0, 13, 633.3333); - test_validity(case_92[0], case_92[1]); + test_validity("case92", case_92[0], case_92[1]); test_one_with_holes("case93", case_93[0], case_93[1], 3, 0, 17, 231.42857); - test_validity(case_93[0], case_93[1]); + test_validity("case93", case_93[0], case_93[1]); test_one_with_holes("case94", case_94[0], case_94[1], 2, 1, 16, 645.71429); - test_validity(case_94[0], case_94[1]); + test_validity("case94", case_94[0], case_94[1]); test_one_with_holes("case95", case_95[0], case_95[1], 3, 0, 16, 226.66667); - test_validity(case_95[0], case_95[1]); + test_validity("case95", case_95[0], case_95[1]); test_one_with_holes("case96", case_96[0], case_96[1], 4, 0, 19, 221.42857); - test_validity(case_96[0], case_96[1]); + test_validity("case96", case_96[0], case_96[1]); + test_one_with_holes("case97", + case_97[0], case_97[1], + 2, 0, 10, 11.81244); + test_validity("case97", case_97[0], case_97[1]); + + test_one_with_holes("case98", + case_98[0], case_98[1], + 2, 0, 10, 17.754744); + test_validity("case98", case_98[0], case_98[1]); return; diff --git a/test/algorithms/set_operations/intersection/intersection_multi.cpp b/test/algorithms/set_operations/intersection/intersection_multi.cpp index f20f46ca77..8c9668e0d0 100644 --- a/test/algorithms/set_operations/intersection/intersection_multi.cpp +++ b/test/algorithms/set_operations/intersection/intersection_multi.cpp @@ -90,6 +90,13 @@ void test_areal() test_one("case_107_multi", case_107_multi[0], case_107_multi[1], 2, 10, 1.5); + + test_one("case_110_multi", + case_110_multi[0], case_110_multi[1], + 2, 11, 9.8050578678287668); + test_validity("case_110_multi", + case_110_multi[0], case_110_multi[1]); + test_one("case_recursive_boxes_1", case_recursive_boxes_1[0], case_recursive_boxes_1[1], 10, 97, 47.0); @@ -106,7 +113,7 @@ void test_areal() test_one_with_holes( "case_recursive_boxes_4", case_recursive_boxes_4[0], case_recursive_boxes_4[1], - 13, 7, 162, 67.0); // Area from SQL Server + 13, 7, 169, 67.0); // Area from SQL Server test_one("ggl_list_20120915_h2_a", ggl_list_20120915_h2[0], ggl_list_20120915_h2[1], diff --git a/test/algorithms/set_operations/intersection/test_intersection.hpp b/test/algorithms/set_operations/intersection/test_intersection.hpp index c2a2943385..1184dd416f 100644 --- a/test/algorithms/set_operations/intersection/test_intersection.hpp +++ b/test/algorithms/set_operations/intersection/test_intersection.hpp @@ -354,17 +354,27 @@ void test_point_output(std::string const& wkt1, std::string const& wkt2, unsigne template -inline void test_validity(std::string const& wkt1, +inline void test_validity(std::string const& caseid, + std::string const& wkt1, std::string const& wkt2) { Areal1 a1; Areal2 a2; bg::read_wkt(wkt1, a1); bg::read_wkt(wkt2, a2); + bg::correct(a1); + bg::correct(a2); + bg::model::multi_polygon out; bg::intersection(a1, a2, out); - BOOST_CHECK(bg::is_valid(out)); -} + std::string reason; + bool b = bg::is_valid(out, reason); + BOOST_CHECK_MESSAGE(b, + "caseid: " << caseid << "; g1: " << bg::wkt(a1) + << "; g2: " << bg::wkt(a2) + << "; i: " << bg::wkt(out) + << "; reason: " << reason); +} #endif From 976f5829e17f51c071da3d9573480ff4701eb301 Mon Sep 17 00:00:00 2001 From: Menelaos Karavelas Date: Fri, 16 Oct 2015 12:15:38 +0300 Subject: [PATCH 20/22] [test][algorithms][union] add function to check the validity of the union set operation; add more test cases; --- .../algorithms/set_operations/union/test_union.hpp | 24 +++++++++++++++++++++- test/algorithms/set_operations/union/union.cpp | 12 +++++++++++ .../set_operations/union/union_multi.cpp | 7 +++++++ 3 files changed, 42 insertions(+), 1 deletion(-) diff --git a/test/algorithms/set_operations/union/test_union.hpp b/test/algorithms/set_operations/union/test_union.hpp index 3d582d2ca6..53c779d011 100644 --- a/test/algorithms/set_operations/union/test_union.hpp +++ b/test/algorithms/set_operations/union/test_union.hpp @@ -28,6 +28,7 @@ #include #include #include +#include #include #include @@ -202,6 +203,27 @@ void test_one(std::string const& caseid, std::string const& wkt1, std::string co expected_area, percentage); } - +template +inline void test_validity(std::string const& caseid, + std::string const& wkt1, + std::string const& wkt2) +{ + Areal1 a1; + Areal2 a2; + bg::read_wkt(wkt1, a1); + bg::read_wkt(wkt2, a2); + bg::correct(a1); + bg::correct(a2); + + bg::model::multi_polygon out; + bg::union_(a1, a2, out); + + std::string reason; + bool b = bg::is_valid(out, reason); + BOOST_CHECK_MESSAGE(b, + "caseid: " << caseid + << "; geometry: " << bg::wkt(out) + << "; reason: " << reason); +} #endif diff --git a/test/algorithms/set_operations/union/union.cpp b/test/algorithms/set_operations/union/union.cpp index be7884ac8b..651ce0a79f 100644 --- a/test/algorithms/set_operations/union/union.cpp +++ b/test/algorithms/set_operations/union/union.cpp @@ -196,6 +196,18 @@ void test_areal() test_one("86", case_86[0], case_86[1], 1, 1, 10, 1500); + test_one("97", + case_97[0], case_97[1], 1, 1, 17, 127.6875598); + test_validity("97", case_97[0], case_97[1]); + + test_one("98", + case_98[0], case_98[1], 1, 2, 23, 155.2452558); + test_validity("98", case_98[0], case_98[1]); + + test_one("103", + case_103[0], case_103[1], 1, 0, 10, 64.072499); + test_validity("103", case_103[0], case_103[1]); + /* test_one(102, simplex_normal[0], simplex_reversed[1], diff --git a/test/algorithms/set_operations/union/union_multi.cpp b/test/algorithms/set_operations/union/union_multi.cpp index 38ea938d92..dbc435b1cc 100644 --- a/test/algorithms/set_operations/union/union_multi.cpp +++ b/test/algorithms/set_operations/union/union_multi.cpp @@ -102,6 +102,7 @@ void test_areal() test_one("case_105_multi", case_105_multi[0], case_105_multi[1], 1, 0, 5, 25); + test_one("case_108_multi", case_108_multi[0], case_108_multi[1], 1, 2, 14, 1400); @@ -109,6 +110,12 @@ void test_areal() case_109_multi[0], case_109_multi[1], 1, 9, 45, 1250); + test_one("case_110_multi", + case_110_multi[0], case_110_multi[1], + 1, 1, 19, 99.194942); + test_validity("case_110_multi", + case_110_multi[0], case_110_multi[1]); + test_one("case_recursive_boxes_1", case_recursive_boxes_1[0], case_recursive_boxes_1[1], 1, 1, 36, 97.0); From 6949f5f8d49a3a50614dffd45f7368c841a5a584 Mon Sep 17 00:00:00 2001 From: Menelaos Karavelas Date: Fri, 16 Oct 2015 12:37:36 +0300 Subject: [PATCH 21/22] [test][algorithms][difference][sym_difference] add function to check the validity of the difference and the symmetric difference set operations; add more test cases; --- .../set_operations/difference/difference.cpp | 44 ++++++++++++++++++++ .../set_operations/difference/difference_multi.cpp | 6 +++ .../set_operations/difference/test_difference.hpp | 48 +++++++++++++++++++++- 3 files changed, 97 insertions(+), 1 deletion(-) diff --git a/test/algorithms/set_operations/difference/difference.cpp b/test/algorithms/set_operations/difference/difference.cpp index dd20682310..d6d233baa3 100644 --- a/test/algorithms/set_operations/difference/difference.cpp +++ b/test/algorithms/set_operations/difference/difference.cpp @@ -411,6 +411,50 @@ void test_all() 8, 36, 2.43452380952381, 7, 33, 3.18452380952381); + test_one( + "case97", case_97[0], case_97[1], + 2, 15, 105.6875598, + 3, 14, 10.18755981); + test_validity("case97", case_97[0], case_97[1]); + +#if ! defined(BOOST_GEOMETRY_NO_ROBUSTNESS) + test_one( + "case98", case_98[0], case_98[1], + 3, 20, 99.7452558, + 3, 14, 37.7452558); +#else + test_one( + "case98", case_98[0], case_98[1], + 3, 20, 99.7452558, + 3, 14, 37.7452558, + 5, 33, 99.7452558 + 37.7452558); +#endif + test_validity("case98", case_98[0], case_98[1]); + +#if defined(BOOST_GEOMETRY_NO_ROBUSTNESS) + test_one( + "case99", case_99[0], case_99[1], + 1, 9, 71938982209960000, + 1, 5, 2.40230376e+20); + test_one( + "case100", case_100[0], case_100[1], + 1, 9, 1108087364.695, + 1, 5, 1.49257218e+22, + 1, 10, 1108087364.695); + test_one( + "case101", case_101[0], case_101[1], + 2, 10, 332066650.1372, + 2, 10, 4.41568936e+20, + 2, 15, 410392630.5984); + test_one( + "case102", case_102[0], case_102[1], + 2, 10, 7.16365077392662e19, + 2, 10, 90403104.0, + 2, 18, 7.16365077392662e19); + test_validity( + "case102", case_102[0], case_102[1]); +#endif // BOOST_GEOMETRY_NO_ROBUSTNESS + test_one("winded", winded[0], winded[1], 3, 37, 61, diff --git a/test/algorithms/set_operations/difference/difference_multi.cpp b/test/algorithms/set_operations/difference/difference_multi.cpp index 335a086708..9c95430777 100644 --- a/test/algorithms/set_operations/difference/difference_multi.cpp +++ b/test/algorithms/set_operations/difference/difference_multi.cpp @@ -91,6 +91,12 @@ void test_areal() case_78_multi[0], case_78_multi[1], 1, 5, 1.0, 1, 5, 1.0); + test_one("case_110_multi", + case_110_multi[0], case_110_multi[1], + 3, 13, 6.194942132, 4, 20, 83.194942); + test_validity("case_110_multi", + case_110_multi[0], case_110_multi[1]); + // Ticket on GGL list 2011/10/25 // to mix polygon/multipolygon in call to difference test_one("ggl_list_20111025_vd_pp", diff --git a/test/algorithms/set_operations/difference/test_difference.hpp b/test/algorithms/set_operations/difference/test_difference.hpp index 72eb218a63..29894483cb 100644 --- a/test/algorithms/set_operations/difference/test_difference.hpp +++ b/test/algorithms/set_operations/difference/test_difference.hpp @@ -21,6 +21,7 @@ #include #include +#include #include #include @@ -346,6 +347,51 @@ void test_one_lp(std::string const& caseid, difference_output(lp + caseid, g1, g2, pieces); } - +template +inline void test_validity(std::string const& caseid, + std::string const& wkt1, + std::string const& wkt2, + bool check_sym_difference = true) +{ + Areal1 a1; + Areal2 a2; + bg::read_wkt(wkt1, a1); + bg::read_wkt(wkt2, a2); + bg::correct(a1); + bg::correct(a2); + + bg::model::multi_polygon df12; + bg::difference(a1, a2, df12); + std::string reason12; + bool b12 = bg::is_valid(df12, reason12); + BOOST_CHECK_MESSAGE(b12, + caseid << "_a; g1: " << bg::wkt(a1) + << "; g2: " << bg::wkt(a2) + << "; df12: " << bg::wkt(df12) + << "; reason: " << reason12); + + bg::model::multi_polygon df21; + bg::difference(a1, a2, df21); + std::string reason21; + bool b21 = bg::is_valid(df21, reason21); + BOOST_CHECK_MESSAGE(b21, + caseid << "_b; g2: " << bg::wkt(a2) + << "; g1: " << bg::wkt(a1) + << "; df21: " << bg::wkt(df21) + << "; reason: " << reason21); + + if (check_sym_difference) + { + bg::model::multi_polygon sdf; + bg::sym_difference(a1, a2, sdf); + std::string reason_s; + bool b_s = bg::is_valid(sdf, reason_s); + BOOST_CHECK_MESSAGE(b_s, + caseid << "_s; g1: " << bg::wkt(a1) + << "; g2: " << bg::wkt(a2) + << "; sdf: " << bg::wkt(sdf) + << "; reason: " << reason_s); + } +} #endif From 0770716f391e5ddb05eac33bb1b340f25620847c Mon Sep 17 00:00:00 2001 From: Menelaos Karavelas Date: Fri, 16 Oct 2015 12:39:08 +0300 Subject: [PATCH 22/22] [algorithms][overlay] add preprocessing step in overlay algorithm: before computing the overlay, compute the self turns of the two geometries and insert the points corresponding to any touch interior self turns they have; This preprocessing step is important for robustness reasons: suppose we have a geometry that has a touch interior self turn, the touch point T of which lies in the interior of some segment AB of the geometry. During a set operation the segment AB may be clipped to another segment CD containting T, where at least one of C and D lies in the interior of AB. By construction T must be a touch point that lies in the interior of CD, but numerically (i.e., since C and D may be computed only approximately) T may be considered as not a point on CD, which can then indicate that the computed geometry is invalid. The solution proposed in this commit basically inserts T in the geometry that has T as a self turns, and this enforces T to be present on CD when the set operation is performed. --- .../detail/overlay/insert_touch_interior_turns.hpp | 375 +++++++++++++++++++++ .../geometry/algorithms/detail/overlay/overlay.hpp | 80 ++++- 2 files changed, 437 insertions(+), 18 deletions(-) create mode 100644 include/boost/geometry/algorithms/detail/overlay/insert_touch_interior_turns.hpp diff --git a/include/boost/geometry/algorithms/detail/overlay/insert_touch_interior_turns.hpp b/include/boost/geometry/algorithms/detail/overlay/insert_touch_interior_turns.hpp new file mode 100644 index 0000000000..bf1becf659 --- /dev/null +++ b/include/boost/geometry/algorithms/detail/overlay/insert_touch_interior_turns.hpp @@ -0,0 +1,375 @@ +// Boost.Geometry (aka GGL, Generic Geometry Library) + +// 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_GEOMETRY_ALGORITHMS_DETAIL_OVERLAY_INSERT_TOUCH_INTERIOR_TURNS_HPP +#define BOOST_GEOMETRY_ALGORITHMS_DETAIL_OVERLAY_INSERT_TOUCH_INTERIOR_TURNS_HPP + +#include +#include +#include + +#include + +#include +#include +#include +#include +#include + +#include + +#include +#include +#include + +#include +#include +#include + + +namespace boost { namespace geometry +{ + +namespace detail { namespace overlay +{ + + +struct self_turns_no_interrupt_policy +{ + static bool const enabled = false; + static bool const has_intersections = false; + + template + static inline bool apply(Range const&) + { + return false; + } +}; + + +template +static inline +typename Turn::turn_operation_type get_correct_op(Turn const& t) +{ + return + (t.operations[0].fraction.is_zero() + || t.operations[0].fraction.is_one()) + ? + t.operations[1] + : + t.operations[0] + ; +} + + +template +struct maa_turn_less +{ + bool operator()(MAA_Turn const& t1, MAA_Turn const& t2) const + { + BOOST_GEOMETRY_ASSERT(t1.method == method_touch_interior); + BOOST_GEOMETRY_ASSERT(t2.method == method_touch_interior); + + typename MAA_Turn::turn_operation_type op1 = get_correct_op(t1); + typename MAA_Turn::turn_operation_type op2 = get_correct_op(t2); + + BOOST_GEOMETRY_ASSERT(! op1.fraction.is_zero() + && ! op1.fraction.is_one()); + BOOST_GEOMETRY_ASSERT(! op2.fraction.is_zero() + && ! op2.fraction.is_one()); + + + if (op1.seg_id.multi_index != op2.seg_id.multi_index) + { + return op1.seg_id.multi_index < op2.seg_id.multi_index; + } + if (op1.seg_id.ring_index != op2.seg_id.ring_index) + { + return op1.seg_id.ring_index < op2.seg_id.ring_index; + } + if (op1.seg_id.segment_index != op2.seg_id.segment_index) + { + return op1.seg_id.segment_index < op2.seg_id.segment_index; + } + return op1.fraction < op2.fraction; + } +}; + + +template ::type> +struct insert_maa_turns +{}; + +template +struct insert_maa_turns +{ + template + static inline TurnIterator apply(Box const& box, + TurnIterator first, + TurnIterator, + BoxOut& box_out, + int = -1, + int = -1) + { + geometry::convert(box, box_out); + return first; + } +}; + + +template +struct insert_maa_turns +{ + template + static inline TurnIterator apply(Ring const& ring, + TurnIterator first, + TurnIterator last, + RingOut& ring_out, + int ring_index = -1, + int multi_index = -1) + { + typedef typename boost::range_iterator + < + Ring const + >::type iterator_type; + + typedef typename std::iterator_traits + < + TurnIterator + >::value_type::turn_operation_type operation_type; + + typename boost::range_size::type point_index = 0; + for (iterator_type it = boost::begin(ring); + it != boost::end(ring); + ++it, ++point_index) + { + geometry::append(ring_out, ring[point_index]); + while (first != last) + { + operation_type op = get_correct_op(*first); + if (op.seg_id.multi_index == multi_index + && + op.seg_id.ring_index == ring_index + && + op.seg_id.segment_index + == static_cast(point_index)) + { + geometry::append(ring_out, first->point); + ++first; + } + else + { + break; + } + } + } + return first; + } +}; + + +template +struct insert_maa_turns +{ + template + < + typename Ring, + typename RingIterator, + typename PolygonOut, + typename TurnIterator + > + static inline + TurnIterator do_interior_rings(RingIterator rings_first, + RingIterator rings_last, + TurnIterator first, + TurnIterator last, + PolygonOut& polygon_out, + int multi_index) + { + typename geometry::ring_type::type ring_out; + int ring_index = 0; + for (RingIterator rit = rings_first; + rit != rings_last; + ++rit, ++ring_index) + { + geometry::clear(ring_out); + first = insert_maa_turns::apply(*rit, + first, + last, + ring_out, + ring_index, + multi_index); + geometry::traits::push_back + < + typename boost::remove_reference + < + typename traits::interior_mutable_type + < + PolygonOut + >::type + >::type + >::apply(geometry::interior_rings(polygon_out), ring_out); + } + return first; + } + + template + < + typename Ring, + typename InteriorRings, + typename PolygonOut, + typename TurnIterator + > + static inline + TurnIterator do_interior_rings(InteriorRings const& irings, + TurnIterator first, + TurnIterator last, + PolygonOut& polygon_out, + int multi_index) + { + return do_interior_rings(boost::begin(irings), + boost::end(irings), + first, + last, + polygon_out, + multi_index); + } + + template + static inline TurnIterator apply(Polygon const& polygon, + TurnIterator first, + TurnIterator last, + PolygonOut& polygon_out, + int multi_index = -1) + { + typedef typename geometry::ring_type::type ring_type; + + typename geometry::ring_type::type ring_out; + + // do exterior ring + first = insert_maa_turns + < + ring_type + >::apply(exterior_ring(polygon), + first, + last, + ring_out, + -1, + multi_index); + geometry::append(polygon_out, ring_out, -1); + + return do_interior_rings(interior_rings(polygon), + first, + last, + polygon_out, + multi_index); + } +}; + +template +struct insert_maa_turns +{ + template + static inline TurnIterator apply(MultiPolygon const& multi_polygon, + TurnIterator first, + TurnIterator last, + MultiPolygonOut& multi_polygon_out) + { + typedef typename boost::range_iterator + < + MultiPolygon const + >::type iterator_type; + + typename boost::range_value::type polygon_out; + + int polygon_index = 0; + for (iterator_type it = boost::begin(multi_polygon); + it != boost::end(multi_polygon); + ++it, ++polygon_index) + { + geometry::clear(polygon_out); + + first = insert_maa_turns + < + typename boost::range_value::type + >::apply(*it, first, last, polygon_out, polygon_index); + + geometry::traits::push_back + < + MultiPolygonOut + >::apply(multi_polygon_out, polygon_out); + } + return first; + } +}; + + +// returns true if the input geometry has been modified (in which case +// the modified geometry is stored in geometry_out), false otherwise +template +inline bool insert_touch_interior_turns(GeometryIn const& geometry_in, + GeometryOut& geometry_out, + RobustPolicy const& robust_policy) +{ + typedef turn_info + < + typename point_type::type, + typename geometry::segment_ratio_type + < + typename point_type::type, RobustPolicy + >::type + > turn_type; + + typedef std::deque turns_container_type; + + self_turns_no_interrupt_policy interrupt_policy; + + // compute self turns + turns_container_type turns; + geometry::self_turns + < + get_turn_info + >(geometry_in, robust_policy, turns, interrupt_policy); + + // select touch interior turns + typedef std::set > maa_turn_set; + + maa_turn_set maa_turns; + for (typename turns_container_type::const_iterator it = turns.begin(); + it != turns.end(); + ++it) + { + if (it->method == method_touch_interior) + { + maa_turns.insert(*it); + } + } + + // if not such turns indicate that the input geometry need not change + if (maa_turns.empty()) + { + return false; + } + + // insert the touch interior turns + insert_maa_turns + < + GeometryIn + >::apply(geometry_in, maa_turns.begin(), maa_turns.end(), geometry_out); + + // return that the input geometry needed change + return true; +} + + +}} // detail::overlay + +}} // boost::geometry + +#endif // BOOST_GEOMETRY_ALGORITHMS_DETAIL_OVERLAY_INSERT_TOUCH_INTERIOR_TURNS_HPP diff --git a/include/boost/geometry/algorithms/detail/overlay/overlay.hpp b/include/boost/geometry/algorithms/detail/overlay/overlay.hpp index 57eb25146f..0a923ed380 100644 --- a/include/boost/geometry/algorithms/detail/overlay/overlay.hpp +++ b/include/boost/geometry/algorithms/detail/overlay/overlay.hpp @@ -22,6 +22,7 @@ #include #include +#include #include #include @@ -38,8 +39,10 @@ #include #include +#include #include #include +#include #include #include @@ -158,31 +161,16 @@ template typename GeometryOut, overlay_type Direction > -struct overlay +class overlay { +private: template - static inline OutputIterator apply( + static inline OutputIterator do_overlay( Geometry1 const& geometry1, Geometry2 const& geometry2, RobustPolicy const& robust_policy, OutputIterator out, Strategy const& ) { - bool const is_empty1 = geometry::is_empty(geometry1); - bool const is_empty2 = geometry::is_empty(geometry2); - - if (is_empty1 && is_empty2) - { - return out; - } - - if (is_empty1 || is_empty2) - { - return return_if_one_input_is_empty - < - GeometryOut, Direction, ReverseOut - >(geometry1, geometry2, out); - } - typedef typename geometry::point_type::type point_type; typedef detail::overlay::traversal_turn_info < @@ -271,6 +259,62 @@ std::cout << "traverse" << std::endl; return add_rings(selected_ring_properties, geometry1, geometry2, rings, out); } + +public: + template + static inline OutputIterator apply( + Geometry1 const& geometry1, Geometry2 const& geometry2, + RobustPolicy const& robust_policy, + OutputIterator out, + Strategy const& strategy) + { + bool const is_empty1 = geometry::is_empty(geometry1); + bool const is_empty2 = geometry::is_empty(geometry2); + + if (is_empty1 && is_empty2) + { + return out; + } + + if (is_empty1 || is_empty2) + { + return return_if_one_input_is_empty + < + GeometryOut, Direction, ReverseOut + >(geometry1, geometry2, out); + } + + Geometry1 modified_geometry1; + bool modified1 = insert_touch_interior_turns(geometry1, + modified_geometry1, + robust_policy); + + Geometry2 modified_geometry2; + bool modified2 = insert_touch_interior_turns(geometry2, + modified_geometry2, + robust_policy); + + if (modified1 && modified2) + { + return do_overlay(modified_geometry1, modified_geometry2, + robust_policy, out, strategy); + } + else if (! modified1 && modified2) + { + return do_overlay(geometry1, modified_geometry2, + robust_policy, out, strategy); + } + else if (modified1 && ! modified2) + { + return do_overlay(modified_geometry1, geometry2, + robust_policy, out, strategy); + } + else + { + return do_overlay(geometry1, geometry2, + robust_policy, out, strategy); + } + } };