diff --git a/include/boost/geometry/algorithms/detail/overlay/discard_duplicate_turns.hpp b/include/boost/geometry/algorithms/detail/overlay/discard_duplicate_turns.hpp index 604dd9fd36..f85483a1d5 100644 --- a/include/boost/geometry/algorithms/detail/overlay/discard_duplicate_turns.hpp +++ b/include/boost/geometry/algorithms/detail/overlay/discard_duplicate_turns.hpp @@ -14,9 +14,11 @@ #define BOOST_GEOMETRY_ALGORITHMS_DETAIL_OVERLAY_DISCARD_DUPLICATE_TURNS_HPP #include +#include #include #include +#include namespace boost { namespace geometry { @@ -88,72 +90,154 @@ inline bool corresponding_turn(Turn const& turn, Turn const& start_turn, return count == 2; } -template -inline void discard_duplicate_start_turns(Turns& turns, - Geometry0 const& geometry0, - Geometry1 const& geometry1) +// Verify turns (other than start, and cross) if they are present in the map, and if so, +// if they have the other turns as corresponding, discard the start turn. +template +void discard_duplicate_start_turns(Turns& turns, + TurnBySegmentMap const& start_turns_by_segment, + Geometry0 const& geometry0, + Geometry1 const& geometry1) { - // Start turns are generated, in case the previous turn is missed. - // But often it is not missed, and then it should be deleted. - // This is how it can be - // (in float, collinear, points far apart due to floating point precision) - // [m, i s:0, v:6 1/1 (1) // u s:1, v:5 pnt (2.54044, 3.12623)] - // [s, i s:0, v:7 0/1 (0) // u s:1, v:5 pnt (2.70711, 3.29289)] - using multi_and_ring_id_type = std::pair; - auto adapt_id = [](segment_identifier const& seg_id) { return multi_and_ring_id_type{seg_id.multi_index, seg_id.ring_index}; }; - // 1 Build map of start turns (multi/ring-id -> turn indices) - std::map> start_turns_per_segment; - std::size_t index = 0; for (auto const& turn : turns) { - if (turn.method == method_start) + // Any turn which "crosses" does not have a corresponding turn. + // Also avoid comparing "start" with itself + if (turn.method == method_crosses || turn.method == method_start) + { + continue; + } + for (auto const& op : turn.operations) { - for (auto const& op : turn.operations) + auto it = start_turns_by_segment.find(adapt_id(op.seg_id)); + if (it != start_turns_by_segment.end()) { - start_turns_per_segment[adapt_id(op.seg_id)].push_back(index); + for (std::size_t const& i : it->second) + { + if (turns[i].cluster_id == turn.cluster_id + && corresponding_turn(turn, turns[i], geometry0, geometry1)) + { + turns[i].discarded = true; + } + } } } - index++; } +} - // 2: Verify all other turns if they are present in the map, and if so, - // if they have the other turns as corresponding - for (auto const& turn : turns) +// Discard turns for the following (rare) case: +// - they are consecutive +// - the first has a touch, the second a touch_interior +// And then one of the segments touches the others next in the middle. +// This is geometrically not possible, and caused by floating point precision. +// Discard the second (touch interior) +template +void discard_touch_touch_interior_turns(Turns& turns, + Geometry0 const& geometry0, + Geometry1 const& geometry1) +{ + for (auto& current_turn : turns) { - // Any turn which "crosses" does not have a corresponding turn. - // Also avoid comparing "start" with itself. - if (turn.method != method_crosses && turn.method != method_start) + if (current_turn.method != method_touch_interior) { - for (auto const& op : turn.operations) + // Because touch_interior is a rarer case, it is more efficient to start with that + continue; + } + for (auto const& previous_turn : turns) + { + if (previous_turn.method != method_touch) { - auto it = start_turns_per_segment.find(adapt_id(op.seg_id)); - if (it != start_turns_per_segment.end()) - { - for (std::size_t const& i : it->second) - { - if (turns[i].cluster_id != turn.cluster_id) - { - // The turns are not part of the same cluster, - // or one is clustered and the other is not. - // This is not corresponding. - continue; - } - if (corresponding_turn(turn, turns[i], - geometry0, geometry1)) - { - turns[i].discarded = true; - } - } - } + continue; + } + + // Compare 0 with 0 and 1 with 1 + // Note that 0 has always source 0 and 1 has always source 1 + // (not in buffer). Therefore this comparison is OK. + // MAYBE we need to check for buffer. + bool const common0 = current_turn.operations[0].seg_id == previous_turn.operations[0].seg_id; + bool const common1 = current_turn.operations[1].seg_id == previous_turn.operations[1].seg_id; + + // If one of the operations is common, and the other is not, then there is one comment segment. + bool const has_common_segment = common0 != common1; + + if (! has_common_segment) + { + continue; } + + // If the second index (1) is common, we need to check consecutivity of the first index (0) + // and vice versa. + bool const consecutive = + common1 ? is_consecutive(previous_turn.operations[0].seg_id, current_turn.operations[0].seg_id, geometry0, geometry1) + : is_consecutive(previous_turn.operations[1].seg_id, current_turn.operations[1].seg_id, geometry0, geometry1); + + if (consecutive) + { + current_turn.discarded = true; + } + } + } +} + +template +void discard_duplicate_turns(Turns& turns, + Geometry0 const& geometry0, + Geometry1 const& geometry1) +{ + // Start turns are generated, in case the previous turn is missed. + // But often it is not missed, and then it should be deleted. + // This is how it can be + // (in float, collinear, points far apart due to floating point precision) + // [m, i s:0, v:6 1/1 (1) // u s:1, v:5 pnt (2.54044, 3.12623)] + // [s, i s:0, v:7 0/1 (0) // u s:1, v:5 pnt (2.70711, 3.29289)] + // + // Also, if two turns are consecutive, and one is touch and the other touch_interior, + // the touch_interior is discarded. + + using multi_and_ring_id_type = std::pair; + + auto add_to_map = [](auto const& turn, auto& map, std::size_t index) + { + auto adapt_id = [](segment_identifier const& seg_id) + { + return multi_and_ring_id_type{seg_id.multi_index, seg_id.ring_index}; + }; + for (auto const& op : turn.operations) + { + map[adapt_id(op.seg_id)].insert(index); + } + }; + + // Build map of start turns (multi/ring-id -> turn indices) + // and count touch and touch_interior turns (to verify if later checks are needed) + std::map> start_turns_by_segment; + std::size_t touch_count = 0; + std::size_t touch_interior_count = 0; + for (auto const& item : util::enumerate(turns)) + { + auto const& turn = item.value; + switch(turn.method) + { + case method_start: add_to_map(turn, start_turns_by_segment, item.index); break; + case method_touch: touch_count++; break; + case method_touch_interior: touch_interior_count++; break; + default: break; } - index++; + } + + if (!start_turns_by_segment.empty()) + { + discard_duplicate_start_turns(turns, start_turns_by_segment, geometry0, geometry1); + } + + if (touch_count > 0 && touch_interior_count > 0) + { + discard_touch_touch_interior_turns(turns, geometry0, geometry1); } } diff --git a/include/boost/geometry/algorithms/detail/overlay/enrich_intersection_points.hpp b/include/boost/geometry/algorithms/detail/overlay/enrich_intersection_points.hpp index 40eb0ed5a5..215a53298f 100644 --- a/include/boost/geometry/algorithms/detail/overlay/enrich_intersection_points.hpp +++ b/include/boost/geometry/algorithms/detail/overlay/enrich_intersection_points.hpp @@ -442,7 +442,7 @@ inline void enrich_intersection_points(Turns& turns, has_colocations = ! clusters.empty(); } - discard_duplicate_start_turns(turns, geometry1, geometry2); + discard_duplicate_turns(turns, geometry1, geometry2); // Discard turns not part of target overlay for (auto& turn : turns) diff --git a/test/algorithms/overlay/overlay_cases.hpp b/test/algorithms/overlay/overlay_cases.hpp index 3ded5aa154..f1b1bb6134 100644 --- a/test/algorithms/overlay/overlay_cases.hpp +++ b/test/algorithms/overlay/overlay_cases.hpp @@ -1176,6 +1176,12 @@ static std::string issue_1186[2] = "POLYGON((-13848.1446527556 6710443.1496919869,-13848.2559722463 6710440.2884572418,-13847.8106942832 6710440.1096301023,-13847.6993747924 6710443.1496919869,-13847.3654163201 6710442.9708647905,-13846.0295824308 6710442.9708647905,-13846.4748603939 6710435.1024718173,-13847.8106942832 6710435.1024718173,-13848.1446527556 6710435.1024718173,-13849.8144451172 6710443.1496919869,-13848.1446527556 6710443.1496919869),(-13847.4767358109 6710440.1096301023,-13847.8106942832 6710440.1096301023,-13847.9220137740 6710439.9308029665,-13847.5880553017 6710439.7519758362,-13847.4767358109 6710440.1096301023))" }; +static std::string issue_1226[2] = +{ + "POLYGON((-0.90478776881879274807 0.51756843862589896332,-0.91 0.48,-1.2 0.4,-1.4 1.9,-0.90478776881879274807 0.51756843862589896332))", + "POLYGON((-0.91943242964602156508 0.55292377741135378955,-0.90478776881879174887 0.51756843862590162786,-0.91 0.48,-0.91943242964602156508 0.55292377741135378955))" +}; + // Triangle, nearly a line static std::string issue_1229[2] = { @@ -1207,6 +1213,12 @@ static std::string issue_1295[2] = "POLYGON((5.7885861030591501 1.2258588412299519e-16, 4.9808473496052885 -0.91020845951432006, 3.0646473238074286e-16 -0.91020843804041685, 0.0 0.0, 5.7885861030591501 1.2258588412299519e-16))" }; +static std::string issue_1326[2] = +{ + "POLYGON((-2.89115098058377 -0.25,-9.34535427093506 -0.25,-9.34535427093506 1.75,-8.741949140272411 1.75,-8.75132236364928 1.78066623424071,-8.72379640157045 1.98517649674012,-8.697157657744199 2.08892604820013,-8.46516432879152 2.10694502899749,-8.450862779791599 2.11134809476174,-8.38919349162909 2.11284570322019,-8.123031147354849 2.13351860011817,-6.54486857587017 2.15952214498774,-6.33512517561985 2.16974517238479,-6.09011821087143 2.17207715057343,-5.97563428310822 2.17637335012843,-5.57133260872244 2.18630252700583,-5.07887848661571 2.18416804606516,-5.0359022373011 2.18529905658276,-1.19310678637991 2.2589664950919,-0.345028190203889 2.01010408632087,-0.0299072068264135 1.77159468023247,-0.02022234569985 1.75,0 1.75,0 1.7049095146299,0.246630961988226 1.15498766340636,0.268832205492555 0.460222803594299,0.180040345346232 0.355084664215225,0.180036363994124 0.284803921497104,0.08623260647096399 0.174496413890264,0.0780023740039386 0.150962829485067,0.0280089059206103 0.0970044986024092,0 0.07806844106142739,0 -0.25,-0.452569566324256 -0.25,-0.959885872828952 -0.5758450741921251,-0.984492757676798 -0.593037244494381,-0.986053937110457 -0.592652605781279,-1.03380615828627 -0.623323463099221,-1.05433104239722 -0.636707963526427,-1.05455277432398 -0.636648843252325,-1.08230325522736 -0.654472747969937,-1.10154815200623 -0.653773778128255,-1.27688251591406 -0.61035431255586,-2.30803967424193 -0.423843667604139,-2.89115098058377 -0.25))", + "POLYGON((-8.78625207542057 1.47256621042848,-8.7713880762347 1.56833897112359,-8.54870763952499 1.98231380430185,-8.38536606998748 2.07941797642276,-8.33146477996655 2.0901167844798,-6.33512517561985 2.16974517238479,-5.57210759837464 2.18803744441122,-1.10154617615509 2.18080908730283,-0.182335452668298 1.77473544450602,-0.101072642349846 1.70875977436368,-0.00289004633636303 1.59446412805576,0.123812278264287 1.4079409294325,0.130584672613603 1.37357814410482,0.17646634130118 0.897752998467649,0.184490648597164 0.600369907975512,0.180036363994123 0.284803921497103,0.07694307247317141 0.151908777171867,0.0269420690240159 0.09795742893367709,-0.984492757676803 -0.593037244494384,-1.07964346095732 -0.633124388260264,-8.500679967807759 1.25877844922534,-8.78625207542057 1.47256621042848))" +}; + static std::string ggl_list_20120229_volker[3] = { "POLYGON((1716 1554,2076 2250,2436 2352,2796 1248,3156 2484,3516 2688,3516 2688,3156 2484,2796 1248,2436 2352,2076 2250, 1716 1554))", diff --git a/test/algorithms/set_operations/difference/difference.cpp b/test/algorithms/set_operations/difference/difference.cpp index 819c8c8523..0ceec3e35c 100644 --- a/test/algorithms/set_operations/difference/difference.cpp +++ b/test/algorithms/set_operations/difference/difference.cpp @@ -613,6 +613,13 @@ void test_all() TEST_DIFFERENCE(issue_1138, 1, 203161.751, 2, 1237551.0171, 1); + { + ut_settings settings; + settings.set_test_validity(false); + settings.validity_of_sym = false; + TEST_DIFFERENCE_WITH(issue_1226, 1, 0.238037722, 0, 0.0, 1, settings); + } + TEST_DIFFERENCE(issue_1231, 2, 36.798659456837477, 3, 195.2986, 5); TEST_DIFFERENCE(issue_1244, 3, 8, 3, 2, 6); @@ -626,6 +633,8 @@ void test_all() TEST_DIFFERENCE(issue_1295, 1, 9.999, 1, 9.999, 1); #endif + TEST_DIFFERENCE(issue_1326, 3, 6.7128537626409130468, 6, 0.00372806966532758478, 9); + TEST_DIFFERENCE(mysql_21977775, 2, 160.856568913, 2, 92.3565689126, 4); TEST_DIFFERENCE(mysql_21965285, 1, 92.0, 1, 14.0, 1); TEST_DIFFERENCE(mysql_23023665_1, 1, 92.0, 1, 142.5, 2); diff --git a/test/algorithms/set_operations/intersection/intersection.cpp b/test/algorithms/set_operations/intersection/intersection.cpp index 29e4da6144..309ea7d3b2 100644 --- a/test/algorithms/set_operations/intersection/intersection.cpp +++ b/test/algorithms/set_operations/intersection/intersection.cpp @@ -301,6 +301,7 @@ void test_areal() TEST_INTERSECTION(issue_861, 1, -1, 1.4715007684573677693e-10); + TEST_INTERSECTION(issue_1226, 1, -1, 0.00036722862); TEST_INTERSECTION(issue_1229, 0, -1, 0); TEST_INTERSECTION(issue_1231, 1, -1, 54.701340543162516); @@ -309,6 +310,7 @@ void test_areal() TEST_INTERSECTION(issue_1293, 1, -1, 1.49123); TEST_INTERSECTION(issue_1295, 1, -1, 4.90121); + TEST_INTERSECTION(issue_1326, 1, -1, 16.4844); test_one("buffer_mp1", buffer_mp1[0], buffer_mp1[1], 1, 31, 2.271707796); diff --git a/test/algorithms/set_operations/union/union.cpp b/test/algorithms/set_operations/union/union.cpp index e4ea757b82..6e7fb1be48 100644 --- a/test/algorithms/set_operations/union/union.cpp +++ b/test/algorithms/set_operations/union/union.cpp @@ -466,6 +466,9 @@ void test_areal() TEST_UNION(issue_1186, 1, 1, -1, 21.6189); TEST_UNION_REV(issue_1186, 1, 1, -1, 21.6189); + TEST_UNION(issue_1226, 1, 0, -1, 0.238405); + TEST_UNION_REV(issue_1226, 1, 0, -1, 0.238405); + TEST_UNION(issue_1229, 1, 0, -1, 384869.166); TEST_UNION_REV(issue_1229, 1, 0, -1, 384869.166); @@ -481,6 +484,9 @@ void test_areal() TEST_UNION(issue_1295, 1, 0, -1, 17.56273); TEST_UNION_REV(issue_1295, 1, 0, -1, 17.56273); + TEST_UNION(issue_1326, 1, 0, -1, 23.201); + TEST_UNION_REV(issue_1326, 1, 0, -1, 23.201); + TEST_UNION(geos_1, 1, 0, -1, expectation_limits(3458.0, 3461.3203125)); TEST_UNION(geos_2, 1, 0, -1, expectation_limits(349.0625, 350.55102539)); TEST_UNION(geos_3, 1, 0, -1, 29391548.4998779);