From 77838a899547b9f3714ea296106ef4ca83bd6be3 Mon Sep 17 00:00:00 2001 From: Barend Gehrels Date: Wed, 2 Dec 2020 10:35:39 +0100 Subject: [PATCH 1/5] [copy_segment_point] change implementation to allow also negative offsets, including unit test --- .../detail/overlay/copy_segment_point.hpp | 26 +++--- .../detail/overlay/segment_identifier.hpp | 2 +- test/algorithms/overlay/Jamfile | 1 + .../algorithms/overlay/copy_segment_point.cpp | 88 +++++++++++++++++++ 4 files changed, 102 insertions(+), 15 deletions(-) create mode 100644 test/algorithms/overlay/copy_segment_point.cpp diff --git a/include/boost/geometry/algorithms/detail/overlay/copy_segment_point.hpp b/include/boost/geometry/algorithms/detail/overlay/copy_segment_point.hpp index 0e36cf004..49901a8b8 100644 --- a/include/boost/geometry/algorithms/detail/overlay/copy_segment_point.hpp +++ b/include/boost/geometry/algorithms/detail/overlay/copy_segment_point.hpp @@ -29,7 +29,6 @@ #include #include #include -#include #include #include #include @@ -48,7 +47,7 @@ template ::type iterator; - geometry::ever_circling_iterator it(boost::begin(view), boost::end(view), - boost::begin(view) + seg_id.segment_index, true); + std::size_t const segment_count = boost::size(view) - 1; - for (signed_size_type i = 0; i < offset; ++i, ++it) - { - } + while (offset < 0) { offset += segment_count; } + + signed_size_type const target = (seg_id.segment_index + offset) % segment_count; + + geometry::convert(*(boost::begin(view) + target), point); - geometry::convert(*it, point); return true; } }; @@ -84,7 +82,7 @@ template inline bool copy_segment_point(Geometry const& geometry, - SegmentIdentifier const& seg_id, int offset, + SegmentIdentifier const& seg_id, signed_size_type offset, PointOut& point_out) { concepts::check(); @@ -316,7 +314,7 @@ template typename PointOut > inline bool copy_segment_point(Geometry1 const& geometry1, Geometry2 const& geometry2, - SegmentIdentifier const& seg_id, int offset, + SegmentIdentifier const& seg_id, signed_size_type offset, PointOut& point_out) { concepts::check(); diff --git a/include/boost/geometry/algorithms/detail/overlay/segment_identifier.hpp b/include/boost/geometry/algorithms/detail/overlay/segment_identifier.hpp index b2d6d7f34..4d6c9fc02 100644 --- a/include/boost/geometry/algorithms/detail/overlay/segment_identifier.hpp +++ b/include/boost/geometry/algorithms/detail/overlay/segment_identifier.hpp @@ -27,11 +27,11 @@ namespace boost { namespace geometry { - // Internal struct to uniquely identify a segment // on a linestring,ring // or polygon (needs ring_index) // or multi-geometry (needs multi_index) +// It is always used for clockwise indication (even if the original is anticlockwise) struct segment_identifier { inline segment_identifier() diff --git a/test/algorithms/overlay/Jamfile b/test/algorithms/overlay/Jamfile index be413b335..b9575595c 100644 --- a/test/algorithms/overlay/Jamfile +++ b/test/algorithms/overlay/Jamfile @@ -16,6 +16,7 @@ test-suite boost-geometry-algorithms-overlay : [ run assemble.cpp : : : : algorithms_assemble ] + [ run copy_segment_point.cpp : : : : algorithms_copy_segment_point ] [ run get_turn_info.cpp : : : : algorithms_get_turn_info ] [ run get_turns.cpp : : : : algorithms_get_turns ] [ run get_turns_areal_areal.cpp : : : : algorithms_get_turns_areal_areal ] diff --git a/test/algorithms/overlay/copy_segment_point.cpp b/test/algorithms/overlay/copy_segment_point.cpp new file mode 100644 index 000000000..101908906 --- /dev/null +++ b/test/algorithms/overlay/copy_segment_point.cpp @@ -0,0 +1,88 @@ +// Boost.Geometry +// Unit Test + +// Copyright (c) 2020 Barend Gehrels, Amsterdam, the Netherlands. + +// Use, modification and distribution is subject to the Boost Software License, +// Version 1.0. (See accompanying file LICENSE_1_0.txt or copy at +// http://www.boost.org/LICENSE_1_0.txt) + +#include + +#include + +#include +#include +#include +#include + +template +void test_basic(std::string const& case_id, int line, + std::string const& wkt, bg::segment_identifier const& id, + int offset, std::size_t expected_index) +{ + using point_type = bg::model::point; + using polygon_type = bg::model::polygon; + + polygon_type polygon; + bg::read_wkt(wkt, polygon); + + // Check the result + auto const& ring = bg::exterior_ring(polygon); + + point_type point; + bg::copy_segment_point(polygon, id, offset, point); + + // Sanity check + bool const expectation_in_range = expected_index < ring.size(); + BOOST_CHECK_MESSAGE(expectation_in_range, "Expectation out of range " << expected_index); + if (! expectation_in_range) + { + return; + } + + point_type const expected_point = ring[expected_index]; + + bool const equal = ! bg::disjoint(point, expected_point); + + BOOST_CHECK_MESSAGE(equal, "copy_segment_point: " << case_id + << " line " << line << " at index " << expected_index + << " not as expected: " + << bg::wkt(point) << " vs " << bg::wkt(expected_point)); +} + +template +void test_all(std::string const& case_id, std::string const& wkt) +{ + // Check zero offset, all segment ids + test_basic(case_id, __LINE__, wkt, {0, 0, -1, 0}, 0, 0); + test_basic(case_id, __LINE__, wkt, {0, 0, -1, 1}, 0, 1); + test_basic(case_id, __LINE__, wkt, {0, 0, -1, 2}, 0, 2); + test_basic(case_id, __LINE__, wkt, {0, 0, -1, 3}, 0, 3); + + // Check positive offsets, it should endlessly loop around, regardless of direction or closure + bg::segment_identifier const start{0, 0, -1, 0}; + test_basic(case_id, __LINE__, wkt, start, 1, 1); + test_basic(case_id, __LINE__, wkt, start, 2, 2); + test_basic(case_id, __LINE__, wkt, start, 3, 3); + test_basic(case_id, __LINE__, wkt, start, 4, 0); + test_basic(case_id, __LINE__, wkt, start, 5, 1); + test_basic(case_id, __LINE__, wkt, start, 6, 2); + test_basic(case_id, __LINE__, wkt, start, 7, 3); + + // Check negative offsets + test_basic(case_id, __LINE__, wkt, start, -1, 3); + test_basic(case_id, __LINE__, wkt, start, -2, 2); + test_basic(case_id, __LINE__, wkt, start, -3, 1); + test_basic(case_id, __LINE__, wkt, start, -4, 0); + test_basic(case_id, __LINE__, wkt, start, -5, 3); + test_basic(case_id, __LINE__, wkt, start, -6, 2); + test_basic(case_id, __LINE__, wkt, start, -7, 1); +} + +int test_main(int, char* []) +{ + test_all("closed", "POLYGON((0 2,1 2,1 1,0 1,0 2))"); + test_all("open", "POLYGON((0 2,1 2,1 1,0 1))"); + return 0; +} From abaa211d3ad890ab245cbd8d571ed3a777512944 Mon Sep 17 00:00:00 2001 From: Barend Gehrels Date: Wed, 2 Dec 2020 14:19:50 +0100 Subject: [PATCH 2/5] [sort_by_side] fix cases where the cluster point is approached by segments, but the last point before is colocated with the turn itself This fixes 50% of the errors currently found by recursive_polygons_buffer (when rescaling is turned off) --- .../detail/overlay/copy_segment_point.hpp | 1 - .../algorithms/detail/overlay/get_ring.hpp | 5 +- .../detail/overlay/handle_colocations.hpp | 2 +- .../detail/overlay/sort_by_side.hpp | 98 ++++++++++++++----- .../algorithms/detail/overlay/traversal.hpp | 6 +- .../buffer/buffer_multi_polygon.cpp | 12 +-- test/algorithms/buffer/test_buffer_svg.hpp | 2 +- test/algorithms/overlay/sort_by_side.cpp | 2 +- .../algorithms/overlay/sort_by_side_basic.cpp | 2 +- .../set_operations/difference/difference.cpp | 2 +- .../difference/difference_multi.cpp | 2 +- .../intersection/intersection_multi.cpp | 2 +- .../set_operations/union/union_multi.cpp | 2 +- 13 files changed, 90 insertions(+), 48 deletions(-) diff --git a/include/boost/geometry/algorithms/detail/overlay/copy_segment_point.hpp b/include/boost/geometry/algorithms/detail/overlay/copy_segment_point.hpp index 49901a8b8..10d12085b 100644 --- a/include/boost/geometry/algorithms/detail/overlay/copy_segment_point.hpp +++ b/include/boost/geometry/algorithms/detail/overlay/copy_segment_point.hpp @@ -397,7 +397,6 @@ inline bool copy_segment_points(Geometry1 const& geometry1, Geometry2 const& geo } - }} // namespace boost::geometry #endif // BOOST_GEOMETRY_ALGORITHMS_DETAIL_OVERLAY_COPY_SEGMENT_POINT_HPP diff --git a/include/boost/geometry/algorithms/detail/overlay/get_ring.hpp b/include/boost/geometry/algorithms/detail/overlay/get_ring.hpp index 1cf5b1e9d..4931633b1 100644 --- a/include/boost/geometry/algorithms/detail/overlay/get_ring.hpp +++ b/include/boost/geometry/algorithms/detail/overlay/get_ring.hpp @@ -119,13 +119,12 @@ struct get_ring template -inline std::size_t segment_count_on_ring(Geometry const& geometry, - segment_identifier const& seg_id) +inline std::size_t segment_count_on_ring(Geometry const& geometry, segment_identifier const& seg_id) { typedef typename geometry::tag::type tag; ring_identifier const rid(0, seg_id.multi_index, seg_id.ring_index); // A closed polygon, a triangle of 4 points, including starting point, - // contains 3 segments. So handle as if closed and subtract one. + // contains 3 segments. So handle as if it is closed, and subtract one. return geometry::num_points(detail::overlay::get_ring::apply(rid, geometry), true) - 1; } diff --git a/include/boost/geometry/algorithms/detail/overlay/handle_colocations.hpp b/include/boost/geometry/algorithms/detail/overlay/handle_colocations.hpp index 960e37032..8d6887d95 100644 --- a/include/boost/geometry/algorithms/detail/overlay/handle_colocations.hpp +++ b/include/boost/geometry/algorithms/detail/overlay/handle_colocations.hpp @@ -831,7 +831,7 @@ inline bool fill_sbs(Sbs& sbs, Point& turn_point, } for (int i = 0; i < 2; i++) { - sbs.add(turn.operations[i], turn_index, i, geometry1, geometry2, first); + sbs.add(turn, turn.operations[i], turn_index, i, geometry1, geometry2, first); first = false; } } diff --git a/include/boost/geometry/algorithms/detail/overlay/sort_by_side.hpp b/include/boost/geometry/algorithms/detail/overlay/sort_by_side.hpp index 484a439e0..9b0726deb 100644 --- a/include/boost/geometry/algorithms/detail/overlay/sort_by_side.hpp +++ b/include/boost/geometry/algorithms/detail/overlay/sort_by_side.hpp @@ -66,6 +66,8 @@ struct ranked_point , seg_id(si) {} + using point_type = Point; + Point point; rank_type rank; signed_size_type zone; // index of closed zone, in uu turn there would be 2 zones @@ -252,12 +254,14 @@ public : , m_strategy(strategy) {} + template void add_segment_from(signed_size_type turn_index, int op_index, Point const& point_from, - operation_type op, segment_identifier const& si, + Operation const& op, bool is_origin) { - m_ranked_points.push_back(rp(point_from, turn_index, op_index, dir_from, op, si)); + m_ranked_points.push_back(rp(point_from, turn_index, op_index, + dir_from, op.operation, op.seg_id)); if (is_origin) { m_origin = point_from; @@ -265,24 +269,46 @@ public : } } + template void add_segment_to(signed_size_type turn_index, int op_index, Point const& point_to, - operation_type op, segment_identifier const& si) + Operation const& op) { - m_ranked_points.push_back(rp(point_to, turn_index, op_index, dir_to, op, si)); + m_ranked_points.push_back(rp(point_to, turn_index, op_index, + dir_to, op.operation, op.seg_id)); } + template void add_segment(signed_size_type turn_index, int op_index, Point const& point_from, Point const& point_to, - operation_type op, segment_identifier const& si, - bool is_origin) + Operation const& op, bool is_origin) { - add_segment_from(turn_index, op_index, point_from, op, si, is_origin); - add_segment_to(turn_index, op_index, point_to, op, si); + add_segment_from(turn_index, op_index, point_from, op, is_origin); + add_segment_to(turn_index, op_index, point_to, op); + } + + template + static inline bool approximately_equals(Point1 const& a, Point2 const& b, + T const& limit) + { + // Including distance would introduce cyclic dependencies. + // This simple code works and is efficient for all coordinate systems. + return std::abs(geometry::get<0>(a) - geometry::get<0>(b)) <= limit + && std::abs(geometry::get<1>(a) - geometry::get<1>(b)) <= limit; } template - Point add(Operation const& op, signed_size_type turn_index, int op_index, + static Point walk_back(Operation const& op, int offset, + Geometry1 const& geometry1, + Geometry2 const& geometry2) + { + Point point; + geometry::copy_segment_point(geometry1, geometry2, op.seg_id, offset, point); + return point; + } + + template + Point add(Turn const& turn, Operation const& op, signed_size_type turn_index, int op_index, Geometry1 const& geometry1, Geometry2 const& geometry2, bool is_origin) @@ -291,20 +317,34 @@ public : geometry::copy_segment_points(geometry1, geometry2, op.seg_id, point1, point2, point3); Point const& point_to = op.fraction.is_one() ? point3 : point2; - add_segment(turn_index, op_index, point1, point_to, op.operation, op.seg_id, is_origin); + + int offset = 0; + + // If the point is in the neighbourhood (the limit itself is not important), + // then take a point (or more) further back. + // The limit of offset avoids theoretical infinite loops. In practice it currently + // walks max 1 point back in all cases. + while (approximately_equals(point1, turn.point, 1.0e-6) && offset > -10) + { + point1 = walk_back(op, --offset, geometry1, geometry2); + } + + add_segment(turn_index, op_index, point1, point_to, op, is_origin); + return point1; } - template - void add(Operation const& op, signed_size_type turn_index, int op_index, + template + void add(Turn const& turn, + Operation const& op, signed_size_type turn_index, int op_index, segment_identifier const& departure_seg_id, Geometry1 const& geometry1, Geometry2 const& geometry2, - bool check_origin) + bool is_departure) { - Point const point1 = add(op, turn_index, op_index, geometry1, geometry2, false); + Point potential_origin = add(turn, op, turn_index, op_index, geometry1, geometry2, false); - if (check_origin) + if (is_departure) { bool const is_origin = op.seg_id.source_index == departure_seg_id.source_index @@ -317,7 +357,7 @@ public : if (m_origin_count == 0 || segment_distance < m_origin_segment_distance) { - m_origin = point1; + m_origin = potential_origin; m_origin_segment_distance = segment_distance; } m_origin_count++; @@ -325,25 +365,33 @@ public : } } + template + static signed_size_type segment_count_on_ring(Operation const& op, + Geometry1 const& geometry1, + Geometry2 const& geometry2) + { + // Take wrap into account + // Suppose point_count=10 (10 points, 9 segments), dep.seg_id=7, op.seg_id=2, + // then distance=9-7+2=4, being segments 7,8,0,1 + return op.seg_id.source_index == 0 + ? detail::overlay::segment_count_on_ring(geometry1, op.seg_id) + : detail::overlay::segment_count_on_ring(geometry2, op.seg_id); + } + template static signed_size_type calculate_segment_distance(Operation const& op, segment_identifier const& departure_seg_id, Geometry1 const& geometry1, Geometry2 const& geometry2) { + BOOST_ASSERT(op.seg_id.source_index == departure_seg_id.source_index); + signed_size_type result = op.seg_id.segment_index - departure_seg_id.segment_index; if (op.seg_id.segment_index >= departure_seg_id.segment_index) { // dep.seg_id=5, op.seg_id=7, distance=2, being segments 5,6 - return op.seg_id.segment_index - departure_seg_id.segment_index; + return result; } - // Take wrap into account - // Suppose point_count=10 (10 points, 9 segments), dep.seg_id=7, op.seg_id=2, - // then distance=9-7+2=4, being segments 7,8,0,1 - std::size_t const segment_count - = op.seg_id.source_index == 0 - ? segment_count_on_ring(geometry1, op.seg_id) - : segment_count_on_ring(geometry2, op.seg_id); - return segment_count - departure_seg_id.segment_index + op.seg_id.segment_index; + return segment_count_on_ring(op, geometry1, geometry2) + result; } void apply(Point const& turn_point) diff --git a/include/boost/geometry/algorithms/detail/overlay/traversal.hpp b/include/boost/geometry/algorithms/detail/overlay/traversal.hpp index a6f8b2a9f..f436639cb 100644 --- a/include/boost/geometry/algorithms/detail/overlay/traversal.hpp +++ b/include/boost/geometry/algorithms/detail/overlay/traversal.hpp @@ -748,7 +748,8 @@ public : for (int i = 0; i < 2; i++) { - sbs.add(cluster_turn.operations[i], + sbs.add(cluster_turn, + cluster_turn.operations[i], cluster_turn_index, i, previous_seg_id, m_geometry1, m_geometry2, departure_turn); @@ -823,7 +824,8 @@ public : // Add this turn to the sort-by-side sorter for (int i = 0; i < 2; i++) { - sbs.add(current_turn.operations[i], + sbs.add(current_turn, + current_turn.operations[i], turn_index, i, previous_seg_id, m_geometry1, m_geometry2, true); diff --git a/test/algorithms/buffer/buffer_multi_polygon.cpp b/test/algorithms/buffer/buffer_multi_polygon.cpp index 032b992f2..f8dbac476 100644 --- a/test/algorithms/buffer/buffer_multi_polygon.cpp +++ b/test/algorithms/buffer/buffer_multi_polygon.cpp @@ -558,23 +558,17 @@ void test_all() test_one("nores_mt_6", nores_mt_6, join_round32, end_flat, 16.9018, 1.0); test_one("nores_et_1", nores_et_1, join_round32, end_flat, 18.9866, 1.0); - -#if defined(BOOST_GEOMETRY_USE_RESCALING) || defined(BOOST_GEOMETRY_TEST_FAILURES) test_one("nores_et_2", nores_et_2, join_round32, end_flat, 23.8389, 1.0); test_one("nores_et_3", nores_et_3, join_round32, end_flat, 26.9030, 1.0); -#endif - test_one("nores_et_4", nores_et_4, join_round32, end_flat, 19.9626, 1.0); test_one("nores_et_5", nores_et_5, join_round32, end_flat, 19.9615, 1.0); + test_one("nores_et_6", nores_et_6, join_round32, end_flat, 96.1795, 1.0); + test_one("nores_et_7", nores_et_7, join_round32, end_flat, 105.9577, 1.0); test_one("nores_wn_1", nores_wn_1, join_round32, end_flat, 23.7659, 1.0); test_one("nores_wn_2", nores_wn_2, join_round32, end_flat, 18.2669, 1.0); -#if defined(BOOST_GEOMETRY_USE_RESCALING) || defined(BOOST_GEOMETRY_TEST_FAILURES) - test_one("nores_et_6", nores_et_6, join_round32, end_flat, 96.1795, 1.0); - test_one("nores_et_7", nores_et_7, join_round32, end_flat, 105.9577, 1.0); test_one("nores_wt_1", nores_wt_1, join_round32, end_flat, 80.1609, 1.0); -#endif test_one("neighbouring_small", neighbouring, @@ -619,7 +613,7 @@ int test_main(int, char* []) #endif #if defined(BOOST_GEOMETRY_TEST_FAILURES) - BoostGeometryWriteExpectedFailures(1, 8, 2, 7); + BoostGeometryWriteExpectedFailures(1, 1, 2, 3); #endif return 0; diff --git a/test/algorithms/buffer/test_buffer_svg.hpp b/test/algorithms/buffer/test_buffer_svg.hpp index 2c7db78ce..6f59fe647 100644 --- a/test/algorithms/buffer/test_buffer_svg.hpp +++ b/test/algorithms/buffer/test_buffer_svg.hpp @@ -207,7 +207,7 @@ private : } #endif - bg::model::ring const& corner = piece.m_piece_border.get_full_ring(); + auto const& corner = piece.m_piece_border.get_full_ring(); if (m_zoom && do_pieces) { diff --git a/test/algorithms/overlay/sort_by_side.cpp b/test/algorithms/overlay/sort_by_side.cpp index 1e0f42258..bcd6b3bae 100644 --- a/test/algorithms/overlay/sort_by_side.cpp +++ b/test/algorithms/overlay/sort_by_side.cpp @@ -101,7 +101,7 @@ std::vector gather_cluster_properties( for (int i = 0; i < 2; i++) { turn_operation_type const& op = turn.operations[i]; - sbs.add(op, turn_index, i, geometry1, geometry2, first); + sbs.add(turn, op, turn_index, i, geometry1, geometry2, first); first = false; } } diff --git a/test/algorithms/overlay/sort_by_side_basic.cpp b/test/algorithms/overlay/sort_by_side_basic.cpp index f37bd38c6..e014e7dff 100644 --- a/test/algorithms/overlay/sort_by_side_basic.cpp +++ b/test/algorithms/overlay/sort_by_side_basic.cpp @@ -120,7 +120,7 @@ std::vector apply_get_turns(std::string const& case_id, has_origin = true; } - sbs.add(turn.operations[i], turn_index, i, + sbs.add(turn, turn.operations[i], turn_index, i, geometry1, geometry2, is_origin); } diff --git a/test/algorithms/set_operations/difference/difference.cpp b/test/algorithms/set_operations/difference/difference.cpp index 6d3d49999..ef6b2e140 100644 --- a/test/algorithms/set_operations/difference/difference.cpp +++ b/test/algorithms/set_operations/difference/difference.cpp @@ -624,7 +624,7 @@ int test_main(int, char* []) #if defined(BOOST_GEOMETRY_TEST_FAILURES) // Not yet fully tested for float and long double. // The difference algorithm can generate (additional) slivers - BoostGeometryWriteExpectedFailures(10, 11, 24, 15); + BoostGeometryWriteExpectedFailures(10, 11, 24, 16); #endif return 0; diff --git a/test/algorithms/set_operations/difference/difference_multi.cpp b/test/algorithms/set_operations/difference/difference_multi.cpp index 2203f071f..91e143dff 100644 --- a/test/algorithms/set_operations/difference/difference_multi.cpp +++ b/test/algorithms/set_operations/difference/difference_multi.cpp @@ -524,7 +524,7 @@ int test_main(int, char* []) #if defined(BOOST_GEOMETRY_TEST_FAILURES) // Not yet fully tested for float. // The difference algorithm can generate (additional) slivers - BoostGeometryWriteExpectedFailures(22, 13, 19, 7); + BoostGeometryWriteExpectedFailures(22, 9, 20, 7); #endif return 0; diff --git a/test/algorithms/set_operations/intersection/intersection_multi.cpp b/test/algorithms/set_operations/intersection/intersection_multi.cpp index 585232273..043512267 100644 --- a/test/algorithms/set_operations/intersection/intersection_multi.cpp +++ b/test/algorithms/set_operations/intersection/intersection_multi.cpp @@ -499,7 +499,7 @@ int test_main(int, char* []) #endif #if defined(BOOST_GEOMETRY_TEST_FAILURES) - BoostGeometryWriteExpectedFailures(9, 3, 2, 1); + BoostGeometryWriteExpectedFailures(9, 1, 2, 1); #endif return 0; diff --git a/test/algorithms/set_operations/union/union_multi.cpp b/test/algorithms/set_operations/union/union_multi.cpp index 53c70844b..6d7d7fbc6 100644 --- a/test/algorithms/set_operations/union/union_multi.cpp +++ b/test/algorithms/set_operations/union/union_multi.cpp @@ -489,7 +489,7 @@ int test_main(int, char* []) #endif #if defined(BOOST_GEOMETRY_TEST_FAILURES) - BoostGeometryWriteExpectedFailures(9, 2, 1, 0); + BoostGeometryWriteExpectedFailures(9, 0, 1, 0); #endif return 0; From 324249bb2dd2f348f2dab95359b25db81da7cf4f Mon Sep 17 00:00:00 2001 From: Barend Gehrels Date: Wed, 9 Dec 2020 12:03:17 +0100 Subject: [PATCH 3/5] [copy_segment_point] change offset with modulo, add to box, update unit test --- .../detail/overlay/copy_segment_point.hpp | 36 ++++--- .../algorithms/overlay/copy_segment_point.cpp | 100 +++++++++++++----- test/count_set.hpp | 4 - 3 files changed, 92 insertions(+), 48 deletions(-) diff --git a/include/boost/geometry/algorithms/detail/overlay/copy_segment_point.hpp b/include/boost/geometry/algorithms/detail/overlay/copy_segment_point.hpp index 10d12085b..e9b9a2f04 100644 --- a/include/boost/geometry/algorithms/detail/overlay/copy_segment_point.hpp +++ b/include/boost/geometry/algorithms/detail/overlay/copy_segment_point.hpp @@ -42,6 +42,16 @@ namespace boost { namespace geometry namespace detail { namespace copy_segments { +inline signed_size_type circular_offset(signed_size_type segment_count, signed_size_type index, + signed_size_type offset) +{ + signed_size_type result = (index + offset) % segment_count; + if (result < 0) + { + result += segment_count; + } + return result; +} template struct copy_segment_point_range @@ -66,12 +76,11 @@ struct copy_segment_point_range rview_type view(cview); std::size_t const segment_count = boost::size(view) - 1; + signed_size_type const target = circular_offset(segment_count, seg_id.segment_index, offset); - while (offset < 0) { offset += segment_count; } - - signed_size_type const target = (seg_id.segment_index + offset) % segment_count; - - geometry::convert(*(boost::begin(view) + target), point); + BOOST_GEOMETRY_ASSERT(target >= 0); + BOOST_GEOMETRY_ASSERT(target < boost::size(view)); + geometry::convert(range::at(view, target), point); return true; } @@ -111,15 +120,14 @@ struct copy_segment_point_box SegmentIdentifier const& seg_id, signed_size_type offset, PointOut& point) { - signed_size_type index = seg_id.segment_index; - for (int i = 0; i < offset; i++) - { - index++; - } - boost::array::type, 4> bp; assign_box_corners_oriented(box, bp); - point = bp[index % 4]; + + signed_size_type const target = circular_offset(4, seg_id.segment_index, offset); + BOOST_GEOMETRY_ASSERT(target >= 0); + BOOST_GEOMETRY_ASSERT(target < bp.size()); + + point = bp[target]; return true; } }; @@ -138,7 +146,6 @@ struct copy_segment_point_multi SegmentIdentifier const& seg_id, signed_size_type offset, PointOut& point) { - BOOST_GEOMETRY_ASSERT ( seg_id.multi_index >= 0 @@ -276,9 +283,6 @@ struct copy_segment_point #endif // DOXYGEN_NO_DISPATCH - - - /*! \brief Helper function, copies a point from a segment \ingroup overlay diff --git a/test/algorithms/overlay/copy_segment_point.cpp b/test/algorithms/overlay/copy_segment_point.cpp index 101908906..33c3f9536 100644 --- a/test/algorithms/overlay/copy_segment_point.cpp +++ b/test/algorithms/overlay/copy_segment_point.cpp @@ -16,22 +16,21 @@ #include #include -template -void test_basic(std::string const& case_id, int line, +template +void test_basic(GetRing get_ring, std::string const& case_id, int line, std::string const& wkt, bg::segment_identifier const& id, int offset, std::size_t expected_index) { - using point_type = bg::model::point; - using polygon_type = bg::model::polygon; + using point_type = typename bg::point_type::type; - polygon_type polygon; - bg::read_wkt(wkt, polygon); + Geometry geometry; + bg::read_wkt(wkt, geometry); // Check the result - auto const& ring = bg::exterior_ring(polygon); + auto ring = get_ring(geometry); point_type point; - bg::copy_segment_point(polygon, id, offset, point); + bg::copy_segment_point(geometry, id, offset, point); // Sanity check bool const expectation_in_range = expected_index < ring.size(); @@ -51,38 +50,83 @@ void test_basic(std::string const& case_id, int line, << bg::wkt(point) << " vs " << bg::wkt(expected_point)); } -template -void test_all(std::string const& case_id, std::string const& wkt) +template +void test_geometry(std::string const& case_id, std::string const& wkt, GetRing get_ring) { // Check zero offset, all segment ids - test_basic(case_id, __LINE__, wkt, {0, 0, -1, 0}, 0, 0); - test_basic(case_id, __LINE__, wkt, {0, 0, -1, 1}, 0, 1); - test_basic(case_id, __LINE__, wkt, {0, 0, -1, 2}, 0, 2); - test_basic(case_id, __LINE__, wkt, {0, 0, -1, 3}, 0, 3); + test_basic(get_ring, case_id, __LINE__, wkt, {0, 0, -1, 0}, 0, 0); + test_basic(get_ring, case_id, __LINE__, wkt, {0, 0, -1, 1}, 0, 1); + test_basic(get_ring, case_id, __LINE__, wkt, {0, 0, -1, 2}, 0, 2); + test_basic(get_ring, case_id, __LINE__, wkt, {0, 0, -1, 3}, 0, 3); // Check positive offsets, it should endlessly loop around, regardless of direction or closure bg::segment_identifier const start{0, 0, -1, 0}; - test_basic(case_id, __LINE__, wkt, start, 1, 1); - test_basic(case_id, __LINE__, wkt, start, 2, 2); - test_basic(case_id, __LINE__, wkt, start, 3, 3); - test_basic(case_id, __LINE__, wkt, start, 4, 0); - test_basic(case_id, __LINE__, wkt, start, 5, 1); - test_basic(case_id, __LINE__, wkt, start, 6, 2); - test_basic(case_id, __LINE__, wkt, start, 7, 3); + test_basic(get_ring, case_id, __LINE__, wkt, start, 1, 1); + test_basic(get_ring, case_id, __LINE__, wkt, start, 2, 2); + test_basic(get_ring, case_id, __LINE__, wkt, start, 3, 3); + test_basic(get_ring, case_id, __LINE__, wkt, start, 4, 0); + test_basic(get_ring, case_id, __LINE__, wkt, start, 5, 1); + test_basic(get_ring, case_id, __LINE__, wkt, start, 6, 2); + test_basic(get_ring, case_id, __LINE__, wkt, start, 7, 3); // Check negative offsets - test_basic(case_id, __LINE__, wkt, start, -1, 3); - test_basic(case_id, __LINE__, wkt, start, -2, 2); - test_basic(case_id, __LINE__, wkt, start, -3, 1); - test_basic(case_id, __LINE__, wkt, start, -4, 0); - test_basic(case_id, __LINE__, wkt, start, -5, 3); - test_basic(case_id, __LINE__, wkt, start, -6, 2); - test_basic(case_id, __LINE__, wkt, start, -7, 1); + test_basic(get_ring, case_id, __LINE__, wkt, start, -1, 3); + test_basic(get_ring, case_id, __LINE__, wkt, start, -2, 2); + test_basic(get_ring, case_id, __LINE__, wkt, start, -3, 1); + test_basic(get_ring, case_id, __LINE__, wkt, start, -4, 0); + test_basic(get_ring, case_id, __LINE__, wkt, start, -5, 3); + test_basic(get_ring, case_id, __LINE__, wkt, start, -6, 2); + test_basic(get_ring, case_id, __LINE__, wkt, start, -7, 1); +} + +template +void test_all(std::string const& case_id, std::string const& wkt) +{ + using point_type = bg::model::point; + using polygon_type = bg::model::polygon; + + test_geometry(case_id, wkt, [](polygon_type const& polygon) + { + return bg::exterior_ring(polygon); + }); +} + +template +void test_box(std::string const& case_id, std::string const& wkt) +{ + using point_type = bg::model::point; + using box_type = bg::model::box; + + test_geometry(case_id, wkt, [](box_type const& box) + { + boost::array ring; + bg::detail::assign_box_corners_oriented(box, ring); + return ring; + }); + +} + +void test_circular_offset() +{ + BOOST_CHECK_EQUAL(3, bg::detail::copy_segments::circular_offset(4, 0, -1)); + BOOST_CHECK_EQUAL(2, bg::detail::copy_segments::circular_offset(4, 0, -2)); + BOOST_CHECK_EQUAL(1, bg::detail::copy_segments::circular_offset(4, 0, -3)); + + BOOST_CHECK_EQUAL(6, bg::detail::copy_segments::circular_offset(10, 5, 1)); + BOOST_CHECK_EQUAL(6, bg::detail::copy_segments::circular_offset(10, 5, 11)); + BOOST_CHECK_EQUAL(6, bg::detail::copy_segments::circular_offset(10, 5, 21)); + + BOOST_CHECK_EQUAL(4, bg::detail::copy_segments::circular_offset(10, 5, -1)); + BOOST_CHECK_EQUAL(4, bg::detail::copy_segments::circular_offset(10, 5, -11)); + BOOST_CHECK_EQUAL(4, bg::detail::copy_segments::circular_offset(10, 5, -21)); } int test_main(int, char* []) { + test_circular_offset(); + test_all("closed", "POLYGON((0 2,1 2,1 1,0 1,0 2))"); test_all("open", "POLYGON((0 2,1 2,1 1,0 1))"); + test_box("box", "BOX(0 0,5 5)"); return 0; } diff --git a/test/count_set.hpp b/test/count_set.hpp index 751649171..3d69b5ca4 100644 --- a/test/count_set.hpp +++ b/test/count_set.hpp @@ -29,10 +29,6 @@ struct count_set { m_values.insert(static_cast(value)); } - else - { - std::cout << "EMPTY" << std::endl; - } } count_set(std::size_t value1, std::size_t value2) From 4e8ff8113114b037d5e382af2df4958cac0e918a Mon Sep 17 00:00:00 2001 From: Barend Gehrels Date: Wed, 9 Dec 2020 13:52:52 +0100 Subject: [PATCH 4/5] [sort_by_side] add epsilon to approximately_equals --- .../detail/overlay/sort_by_side.hpp | 28 +++++++++++++++---- .../set_operations/union/union_multi.cpp | 2 +- 2 files changed, 24 insertions(+), 6 deletions(-) diff --git a/include/boost/geometry/algorithms/detail/overlay/sort_by_side.hpp b/include/boost/geometry/algorithms/detail/overlay/sort_by_side.hpp index 9b0726deb..948b03db5 100644 --- a/include/boost/geometry/algorithms/detail/overlay/sort_by_side.hpp +++ b/include/boost/geometry/algorithms/detail/overlay/sort_by_side.hpp @@ -25,6 +25,9 @@ #include #include +#include +#include +#include namespace boost { namespace geometry { @@ -287,14 +290,29 @@ public : add_segment_to(turn_index, op_index, point_to, op); } + // Returns true if two points are approximately equal, tuned by a giga-epsilon constant + // (if constant is 1.0, for type double, the boundary is about 1.0e-7) template static inline bool approximately_equals(Point1 const& a, Point2 const& b, - T const& limit) + T const& limit_giga_epsilon) { // Including distance would introduce cyclic dependencies. - // This simple code works and is efficient for all coordinate systems. - return std::abs(geometry::get<0>(a) - geometry::get<0>(b)) <= limit - && std::abs(geometry::get<1>(a) - geometry::get<1>(b)) <= limit; + using coor_t = typename select_coordinate_type::type; + using calc_t = typename geometry::select_most_precise ::type; + constexpr calc_t machine_giga_epsilon = 1.0e9 * std::numeric_limits::epsilon(); + + calc_t const& a0 = geometry::get<0>(a); + calc_t const& b0 = geometry::get<0>(b); + calc_t const& a1 = geometry::get<1>(a); + calc_t const& b1 = geometry::get<1>(b); + calc_t const one = 1.0; + calc_t const c = math::detail::greatest(a0, b0, a1, b1, one); + + // The maximum limit is avoid, for floating point, large limits like 400 + // (which are be calculated using eps) + constexpr calc_t maxlimit = 1.0e-3; + auto const limit = (std::min)(maxlimit, limit_giga_epsilon * machine_giga_epsilon * c); + return std::abs(a0 - b0) <= limit && std::abs(a1 - b1) <= limit; } template @@ -324,7 +342,7 @@ public : // then take a point (or more) further back. // The limit of offset avoids theoretical infinite loops. In practice it currently // walks max 1 point back in all cases. - while (approximately_equals(point1, turn.point, 1.0e-6) && offset > -10) + while (approximately_equals(point1, turn.point, 1.0) && offset > -10) { point1 = walk_back(op, --offset, geometry1, geometry2); } diff --git a/test/algorithms/set_operations/union/union_multi.cpp b/test/algorithms/set_operations/union/union_multi.cpp index 6d7d7fbc6..bcc0185a5 100644 --- a/test/algorithms/set_operations/union/union_multi.cpp +++ b/test/algorithms/set_operations/union/union_multi.cpp @@ -489,7 +489,7 @@ int test_main(int, char* []) #endif #if defined(BOOST_GEOMETRY_TEST_FAILURES) - BoostGeometryWriteExpectedFailures(9, 0, 1, 0); + BoostGeometryWriteExpectedFailures(9, 0, 3, 0); #endif return 0; From 59e0840d750cbcc7fd1f563ab6fa5fd7a3914ebc Mon Sep 17 00:00:00 2001 From: Barend Gehrels Date: Wed, 16 Dec 2020 14:29:15 +0100 Subject: [PATCH 5/5] [sort_by_side] walk forward for point_to (similarly to walking backwards for point_from) --- .../detail/overlay/sort_by_side.hpp | 25 ++++++++++++------- .../buffer/buffer_multi_polygon.cpp | 4 ++- 2 files changed, 19 insertions(+), 10 deletions(-) diff --git a/include/boost/geometry/algorithms/detail/overlay/sort_by_side.hpp b/include/boost/geometry/algorithms/detail/overlay/sort_by_side.hpp index 948b03db5..316f53a02 100644 --- a/include/boost/geometry/algorithms/detail/overlay/sort_by_side.hpp +++ b/include/boost/geometry/algorithms/detail/overlay/sort_by_side.hpp @@ -316,7 +316,7 @@ public : } template - static Point walk_back(Operation const& op, int offset, + static Point walk_over_ring(Operation const& op, int offset, Geometry1 const& geometry1, Geometry2 const& geometry2) { @@ -331,25 +331,32 @@ public : Geometry2 const& geometry2, bool is_origin) { - Point point1, point2, point3; + Point point_from, point2, point3; geometry::copy_segment_points(geometry1, geometry2, - op.seg_id, point1, point2, point3); - Point const& point_to = op.fraction.is_one() ? point3 : point2; + op.seg_id, point_from, point2, point3); + Point point_to = op.fraction.is_one() ? point3 : point2; - int offset = 0; // If the point is in the neighbourhood (the limit itself is not important), // then take a point (or more) further back. // The limit of offset avoids theoretical infinite loops. In practice it currently // walks max 1 point back in all cases. - while (approximately_equals(point1, turn.point, 1.0) && offset > -10) + int offset = 0; + while (approximately_equals(point_from, turn.point, 1.0) && offset > -10) { - point1 = walk_back(op, --offset, geometry1, geometry2); + point_from = walk_over_ring(op, --offset, geometry1, geometry2); } - add_segment(turn_index, op_index, point1, point_to, op, is_origin); + // Similarly for the point to, walk forward + offset = 0; + while (approximately_equals(point_to, turn.point, 1.0) && offset < 10) + { + point_to = walk_over_ring(op, ++offset, geometry1, geometry2); + } - return point1; + add_segment(turn_index, op_index, point_from, point_to, op, is_origin); + + return point_from; } template diff --git a/test/algorithms/buffer/buffer_multi_polygon.cpp b/test/algorithms/buffer/buffer_multi_polygon.cpp index f8dbac476..f5421ec44 100644 --- a/test/algorithms/buffer/buffer_multi_polygon.cpp +++ b/test/algorithms/buffer/buffer_multi_polygon.cpp @@ -346,7 +346,8 @@ static std::string const nores_wn_2 // Other cases with wrong turn information static std::string const nores_wt_1 = "MULTIPOLYGON(((0 4,0 5,1 4,0 4)),((9 3,9 4,10 4,9 3)),((9 7,10 8,10 7,9 7)),((6 7,7 8,7 7,6 7)),((0 7,0 8,1 8,0 7)),((3 6,4 6,4 5,3 4,3 6)),((3 7,2 6,2 7,3 7)),((3 7,3 8,4 8,4 7,3 7)),((3 3,4 4,4 3,3 3)),((3 3,3 2,2 2,2 3,3 3)),((2 6,2 5,1 5,1 6,2 6)),((6 4,6 3,5 3,5 4,6 4)),((6 4,7 5,7 4,6 4)),((5 1,4 0,4 1,5 1)),((5 1,5 2,6 2,6 1,5 1)))"; - +static std::string const nores_wt_2 + = "MULTIPOLYGON(((1 1,2 2,2 1,1 1)),((0 2,0 3,1 3,0 2)),((4 3,5 4,5 3,4 3)),((4 3,3 3,4 4,4 3)))"; static std::string const neighbouring = "MULTIPOLYGON(((0 0,0 10,10 10,10 0,0 0)),((10 10,10 20,20 20,20 10,10 10)))"; @@ -569,6 +570,7 @@ void test_all() test_one("nores_wn_2", nores_wn_2, join_round32, end_flat, 18.2669, 1.0); test_one("nores_wt_1", nores_wt_1, join_round32, end_flat, 80.1609, 1.0); + test_one("nores_wt_2", nores_wt_2, join_round32, end_flat, 22.1102, 1.0); test_one("neighbouring_small", neighbouring,