From 076d1077c5ecb0fd3fc482fee4599fb9e2db42d0 Mon Sep 17 00:00:00 2001 From: Adam Wulkiewicz Date: Fri, 13 Feb 2015 21:25:30 +0100 Subject: [PATCH] [strategies][policies] Increase robustness of collinear segments intersection. Do not use the ratios when checking the relation of endpoints and the other segment. The ratios depend on segment lengths and if one of the segments was a lot longer than the other one the direction and intersection results were inconsistent. E.g. the endpoints of one segment was detected inside/outside the longer segment and in the same time (using different check) both endpoints was detected equal to one of the endpoints of the longer segment. Then depending on the order of the segments 2 intersection points were generated or 1 which could cause an assertion failure in turn handler. --- .../geometry/policies/relate/direction.hpp | 56 ++++++++++++++--- .../policies/relate/intersection_points.hpp | 9 +-- .../boost/geometry/policies/relate/tupled.hpp | 14 +++-- .../strategies/cartesian/cart_intersect.hpp | 61 ++++++++++++++----- 4 files changed, 107 insertions(+), 33 deletions(-) diff --git a/include/boost/geometry/policies/relate/direction.hpp b/include/boost/geometry/policies/relate/direction.hpp index 02fed94b1..a0eb55533 100644 --- a/include/boost/geometry/policies/relate/direction.hpp +++ b/include/boost/geometry/policies/relate/direction.hpp @@ -206,22 +206,60 @@ struct segments_direction } } + static inline int arrival_from_position_value(int /*v_from*/, int v_to) + { + return v_to == 2 ? 1 + : v_to == 1 || v_to == 3 ? 0 + //: v_from >= 1 && v_from <= 3 ? -1 + : -1; + + // NOTE: this should be an equivalent of the above for the other order + /* (v_from < 3 && v_to > 3) || (v_from > 3 && v_to < 3) ? 1 + : v_from == 3 || v_to == 3 ? 0 + : -1;*/ + } + + static inline void analyse_position_value(int pos_val, + int & in_segment_count, + int & on_end_count, + int & outside_segment_count) + { + if ( pos_val == 1 || pos_val == 3 ) + { + on_end_count++; + } + else if ( pos_val == 2 ) + { + in_segment_count++; + } + else + { + outside_segment_count++; + } + } + template static inline return_type segments_collinear( - Segment1 const& , Segment2 const&, + Segment1 const& , Segment2 const& , + int a1_wrt_b, int a2_wrt_b, int b1_wrt_a, int b2_wrt_a, Ratio const& ra_from_wrt_b, Ratio const& ra_to_wrt_b, - Ratio const& rb_from_wrt_a, Ratio const& rb_to_wrt_a) + Ratio const& /*rb_from_wrt_a*/, Ratio const& /*rb_to_wrt_a*/) { // If segments are opposite, the ratio of the FROM w.r.t. the other // is larger than the ratio of the TO w.r.t. the other bool const opposite = ra_to_wrt_b < ra_from_wrt_b; - + // NOTE: an alternative version could compare ratios only + // if both ax_wrt_b positions were equal, otherwise use int positions + /*bool const opposite = a1_wrt_b == a2_wrt_b + ? ra_to_wrt_b < ra_from_wrt_b + : a2_wrt_b < a1_wrt_b;*/ + return_type r('c', opposite); // IMPORTANT: the order of conditions is different as in intersection_points.hpp // We assign A in 0 and B in 1 - r.arrival[0] = arrival_value(ra_from_wrt_b, ra_to_wrt_b); - r.arrival[1] = arrival_value(rb_from_wrt_a, rb_to_wrt_a); + r.arrival[0] = arrival_from_position_value(a1_wrt_b, a2_wrt_b); + r.arrival[1] = arrival_from_position_value(b1_wrt_a, b2_wrt_a); // Analyse them int a_in_segment_count = 0; @@ -230,13 +268,13 @@ struct segments_direction int b_in_segment_count = 0; int b_on_end_count = 0; int b_outside_segment_count = 0; - analyze(ra_from_wrt_b, + analyse_position_value(a1_wrt_b, a_in_segment_count, a_on_end_count, a_outside_segment_count); - analyze(ra_to_wrt_b, + analyse_position_value(a2_wrt_b, a_in_segment_count, a_on_end_count, a_outside_segment_count); - analyze(rb_from_wrt_a, + analyse_position_value(b1_wrt_a, b_in_segment_count, b_on_end_count, b_outside_segment_count); - analyze(rb_to_wrt_a, + analyse_position_value(b2_wrt_a, b_in_segment_count, b_on_end_count, b_outside_segment_count); if (a_on_end_count == 1 diff --git a/include/boost/geometry/policies/relate/intersection_points.hpp b/include/boost/geometry/policies/relate/intersection_points.hpp index 14248e771..9881d9882 100644 --- a/include/boost/geometry/policies/relate/intersection_points.hpp +++ b/include/boost/geometry/policies/relate/intersection_points.hpp @@ -101,6 +101,7 @@ struct segments_intersection_points template static inline return_type segments_collinear( Segment1 const& a, Segment2 const& b, + int a1_wrt_b, int a2_wrt_b, int b1_wrt_a, int b2_wrt_a, Ratio const& ra_from_wrt_b, Ratio const& ra_to_wrt_b, Ratio const& rb_from_wrt_a, Ratio const& rb_to_wrt_a) { @@ -112,7 +113,7 @@ struct segments_intersection_points // if index would be 2 this indicate an (currently uncatched) error // IMPORTANT: the order of conditions is different as in direction.hpp - if (ra_from_wrt_b.on_segment() + if (a1_wrt_b >= 1 && a1_wrt_b <= 3 // ra_from_wrt_b.on_segment() && index < 2) { // a1--------->a2 @@ -126,7 +127,7 @@ struct segments_intersection_points index++; count_a++; } - if (rb_from_wrt_a.in_segment() + if (b1_wrt_a == 2 //rb_from_wrt_a.in_segment() && index < 2) { // We take the first intersection point of B @@ -143,7 +144,7 @@ struct segments_intersection_points count_b++; } - if (ra_to_wrt_b.on_segment() + if (a2_wrt_b >= 1 && a2_wrt_b <= 3 //ra_to_wrt_b.on_segment() && index < 2) { // Similarly, second IP (here a2) @@ -155,7 +156,7 @@ struct segments_intersection_points index++; count_a++; } - if (rb_to_wrt_a.in_segment() + if (b2_wrt_a == 2 // rb_to_wrt_a.in_segment() && index < 2) { detail::assign_point_from_index<1>(b, result.intersections[index]); diff --git a/include/boost/geometry/policies/relate/tupled.hpp b/include/boost/geometry/policies/relate/tupled.hpp index 6da9095c4..eadf5aa6c 100644 --- a/include/boost/geometry/policies/relate/tupled.hpp +++ b/include/boost/geometry/policies/relate/tupled.hpp @@ -47,13 +47,19 @@ struct segments_tupled template static inline return_type segments_collinear( - Segment1 const& segment1, Segment2 const& segment2, - Ratio const& ra1, Ratio const& ra2, Ratio const& rb1, Ratio const& rb2) + Segment1 const& segment1, Segment2 const& segment2, + int pa1, int pa2, int pb1, int pb2, + Ratio const& ra1, Ratio const& ra2, + Ratio const& rb1, Ratio const& rb2) { return boost::make_tuple ( - Policy1::segments_collinear(segment1, segment2, ra1, ra2, rb1, rb2), - Policy2::segments_collinear(segment1, segment2, ra1, ra2, rb1, rb2) + Policy1::segments_collinear(segment1, segment2, + pa1, pa2, pb1, pb2, + ra1, ra2, rb1, rb2), + Policy2::segments_collinear(segment1, segment2, + pa1, pa2, pb1, pb2, + ra1, ra2, rb1, rb2) ); } diff --git a/include/boost/geometry/strategies/cartesian/cart_intersect.hpp b/include/boost/geometry/strategies/cartesian/cart_intersect.hpp index 4c2c075bc..4e0d1f1a8 100644 --- a/include/boost/geometry/strategies/cartesian/cart_intersect.hpp +++ b/include/boost/geometry/strategies/cartesian/cart_intersect.hpp @@ -371,39 +371,50 @@ private: RatioType rb_to(ob_2 - oa_1, length_a); // use absolute measure to detect endpoints intersection - // NOTE: some of those conditions shouldn't be met in the same time, - // e.g. a1 == b1 && a2 == b1 but they're checked anyway - // CONSIDER: below the ratio is modified if it should indicate that - // the endpoints intersect (but maybe it doesn't) - // should also the opposite case be handled explicitly? That is, - // if ratio indicates that the endpoints intersect but they aren't? - if (math::equals(oa_1, ob_1)) + // NOTE: it'd be possible to calculate bx_wrt_a using ax_wrt_b values + int const a1_wrt_b = position_value(oa_1, ob_1, ob_2); + int const a2_wrt_b = position_value(oa_2, ob_1, ob_2); + int const b1_wrt_a = position_value(ob_1, oa_1, oa_2); + int const b2_wrt_a = position_value(ob_2, oa_1, oa_2); + + // fix the ratios if necessary + // CONSIDER: fixing ratios also in other cases, if they're inconsistent + // e.g. if ratio == 1 or 0 (so IP at the endpoint) + // but position value indicates that the IP is in the middle of the segment + // because one of the segments is very long + // In such case the ratios could be moved into the middle direction + // by some small value (e.g. EPS+1ULP) + if (a1_wrt_b == 1) { ra_from.assign(0, 1); rb_from.assign(0, 1); } - if (math::equals(oa_2, ob_1)) + else if (a1_wrt_b == 3) + { + ra_from.assign(1, 1); + rb_to.assign(0, 1); + } + + if (a2_wrt_b == 1) { ra_to.assign(0, 1); rb_from.assign(1, 1); } - if (math::equals(oa_1, ob_2)) - { - ra_from.assign(1, 1); - rb_to.assign(0, 1); - } - if (math::equals(oa_2, ob_2)) + else if (a2_wrt_b == 3) { ra_to.assign(1, 1); rb_to.assign(1, 1); } - if ((ra_from.left() && ra_to.left()) || (ra_from.right() && ra_to.right())) + if ((a1_wrt_b < 1 && a2_wrt_b < 1) || (a1_wrt_b > 3 && a2_wrt_b > 3)) + //if ((ra_from.left() && ra_to.left()) || (ra_from.right() && ra_to.right())) { return Policy::disjoint(); } - return Policy::segments_collinear(a, b, ra_from, ra_to, rb_from, rb_to); + return Policy::segments_collinear(a, b, + a1_wrt_b, a2_wrt_b, b1_wrt_a, b2_wrt_a, + ra_from, ra_to, rb_from, rb_to); } /// Relate segments where one is degenerate @@ -433,6 +444,24 @@ private: return Policy::one_degenerate(degenerate_segment, ratio, a_degenerate); } + + template + static inline int position_value(ProjCoord1 const& ca1, + ProjCoord2 const& cb1, + ProjCoord2 const& cb2) + { + // S1x 0 1 2 3 4 + // S2 |----------> + return math::equals(ca1, cb1) ? 1 + : math::equals(ca1, cb2) ? 3 + : cb1 < cb2 ? + ( ca1 < cb1 ? 0 + : ca1 > cb2 ? 4 + : 2 ) + : ( ca1 > cb1 ? 0 + : ca1 < cb2 ? 4 + : 2 ); + } };