From 3aced7f7e76bcbee31a776893de8b6d3de6429bd Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?S=C3=A9bastien=20Loriot?= Date: Thu, 19 Nov 2020 19:39:40 +0100 Subject: [PATCH 1/3] add projected segment overlapping case + clean up --- ...Triangulation_2_projection_traits_base_3.h | 80 ++++++++++++++++--- 1 file changed, 69 insertions(+), 11 deletions(-) diff --git a/Triangulation_2/include/CGAL/Triangulation_2/internal/Triangulation_2_projection_traits_base_3.h b/Triangulation_2/include/CGAL/Triangulation_2/internal/Triangulation_2_projection_traits_base_3.h index 984d2381035..337c915289e 100644 --- a/Triangulation_2/include/CGAL/Triangulation_2/internal/Triangulation_2_projection_traits_base_3.h +++ b/Triangulation_2/include/CGAL/Triangulation_2/internal/Triangulation_2_projection_traits_base_3.h @@ -202,36 +202,94 @@ public: Object planes_intersection = intersection(plane_1, plane_2); if(planes_intersection.empty()) { +#ifdef CGAL_T2_PTB_3_DEBUG std::cerr << "planes_intersection is empty\n"; +#endif return planes_intersection; } if(const Line* line = object_cast(&planes_intersection)) { + // check if the intersection line intersects both segments by + // checking if a point on the intersection line is between + // the segments endpoints const Point& pi = line->point(0); - if(cross_product(normal, pi - s1.source()) + if( cross_product(normal, pi - s1.source()) * cross_product(normal, pi - s1.target()) > FT(0) || - cross_product(normal, pi - s2.source()) + cross_product(normal, pi - s2.source()) * cross_product(normal, pi - s2.target()) > FT(0) ) { - // the intersection of the lines is not inside the segments + // the intersection of the supporting lines is not inside both segments +#ifdef CGAL_T2_PTB_3_DEBUG std::cerr << "intersection not inside\n"; return Object(); +#endif } else { - // Let the plane passing through s1.source() and with normal - // the cross product of s1.to_vector() and s2.to_vector(). That - // plane should intersect *l, now. - return intersection(*line, Plane_3(s1.source(), - cross_product(s1.to_vector(), - s2.to_vector()))); + // Let the plane passing through s1.source() and with normal the traits' normal. + // The segment intersection point computed as the intersection + // of the *l and that plane. Note that if segments are not coplanar the positive + // of the intersection point is somehow arbitrary. + return intersection(*line, Plane_3(s1.source(), normal)); } } if(object_cast(&planes_intersection)) { - std::cerr << "coplanar lines\n"; - CGAL_error(); +#ifdef CGAL_T2_PTB_3_DEBUG + std::cerr << "coplanar supporting lines\n"; +#endif + auto is_inside_segment = [](const Segment& s, const Point& q) + { + return Vector_3(q,s.source()) * Vector_3(q,s.target()) <=0; + }; + + bool src2_in_s1 = is_inside_segment(s1, s2.source()); + bool tgt2_in_s1 = is_inside_segment(s1, s2.target()); + bool src1_in_s2 = is_inside_segment(s2, s1.source()); + bool tgt1_in_s2 = is_inside_segment(s2, s1.target()); + + if (src1_in_s2 && tgt1_in_s2) return make_object(s1); + if (src2_in_s1 && tgt2_in_s1) return make_object(s2); + + if (src1_in_s2) + { + if (src2_in_s1) + { + if (cross_product(normal, Vector_3(s1.source(), s2.source())) != NULL_VECTOR) + return make_object(Segment(s1.source(), s2.source())); + else + return make_object(s1.source()); + } + if (tgt2_in_s1) + { + if (cross_product(normal, Vector_3(s1.source(), s2.target())) != NULL_VECTOR) + return make_object(Segment(s1.source(), s2.target())); + else + return make_object(s1.source()); + } + // should never get here with a Kernel with exact constructions + return make_object(s1.source()); + } + if (tgt1_in_s2) + { + if (src2_in_s1) + { + if (cross_product(normal, Vector_3(s1.target(), s2.source())) != NULL_VECTOR) + return make_object(Segment(s1.target(), s2.source())); + else + return make_object(s1.target()); + } + if (tgt2_in_s1) + { + if (cross_product(normal, Vector_3(s1.target(), s2.target())) != NULL_VECTOR) + return make_object(Segment(s1.target(), s2.target())); + else + return make_object(s1.target()); + } + // should never get here with a Kernel with exact constructions + return make_object(s1.target()); + } return Object(); } return Object(); From 3082808434e497a85819f944726d6fe4bbffed4d Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?S=C3=A9bastien=20Loriot?= Date: Mon, 1 Feb 2021 10:12:16 +0100 Subject: [PATCH 2/3] restore original code --- .../Triangulation_2_projection_traits_base_3.h | 11 ++++++----- 1 file changed, 6 insertions(+), 5 deletions(-) diff --git a/Triangulation_2/include/CGAL/Triangulation_2/internal/Triangulation_2_projection_traits_base_3.h b/Triangulation_2/include/CGAL/Triangulation_2/internal/Triangulation_2_projection_traits_base_3.h index 337c915289e..3e586a1eb8a 100644 --- a/Triangulation_2/include/CGAL/Triangulation_2/internal/Triangulation_2_projection_traits_base_3.h +++ b/Triangulation_2/include/CGAL/Triangulation_2/internal/Triangulation_2_projection_traits_base_3.h @@ -227,11 +227,12 @@ public: } else { - // Let the plane passing through s1.source() and with normal the traits' normal. - // The segment intersection point computed as the intersection - // of the *l and that plane. Note that if segments are not coplanar the positive - // of the intersection point is somehow arbitrary. - return intersection(*line, Plane_3(s1.source(), normal)); + // Let the plane passing through s1.source() and with normal + // the cross product of s1.to_vector() and s2.to_vector(). That + // plane should intersect *l, now. + return intersection(*line, Plane_3(s1.source(), + cross_product(s1.to_vector(), + s2.to_vector()))); } } if(object_cast(&planes_intersection)) From 2c1b967db15c3b43b04eb034d14ed938faaa0c6c Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?S=C3=A9bastien=20Loriot?= Date: Mon, 1 Feb 2021 11:19:30 +0100 Subject: [PATCH 3/3] add tests for intersection of segments --- ...t_cdt_2_projection_traits_special_case.cpp | 68 +++++++++++++++++++ 1 file changed, 68 insertions(+) diff --git a/Triangulation_2/test/Triangulation_2/test_cdt_2_projection_traits_special_case.cpp b/Triangulation_2/test/Triangulation_2/test_cdt_2_projection_traits_special_case.cpp index 26693b95a28..bf163771380 100644 --- a/Triangulation_2/test/Triangulation_2/test_cdt_2_projection_traits_special_case.cpp +++ b/Triangulation_2/test/Triangulation_2/test_cdt_2_projection_traits_special_case.cpp @@ -60,6 +60,67 @@ bool test(std::string test_name) return cdt.is_valid(); } + +template +bool test_segment_intersections(std::string test_name) +{ + std::cerr << "Testing " << test_name << std::endl; + + CDT_2_traits traits(typename K::Vector_3(0,0,1)); + typename CDT_2_traits::Intersect_2 intersect = traits.intersect_2_object(); + + typedef typename CDT_2_traits::Segment_2 Segment_2; + typedef typename CDT_2_traits::Point_2 Point_2; + + // test point intersection + Segment_2 s1( Point_2(0,0,0), Point_2(1,1,0) ), + s2( Point_2(0,1,0), Point_2(1,0,0) ); + + CGAL::Object o = intersect(s1,s2); + assert( o.is() ); + +// test segment overlap + Point_2 pts[4] = { Point_2(0,0,0), Point_2(1,0,1), Point_2(2,0,2), Point_2(4,0,3) }; + + o = intersect( Segment_2(pts[0], pts[1]), Segment_2(pts[2], pts[3]) ); + assert( o.empty() ); + + // pure overlap + o = intersect( Segment_2(pts[0], pts[2]), Segment_2(pts[1], pts[3]) ); + assert( o.is() ); + o = intersect( Segment_2(pts[0], pts[2]), Segment_2(pts[3], pts[1]) ); + assert( o.is() ); + o = intersect( Segment_2(pts[2], pts[0]), Segment_2(pts[1], pts[3]) ); + assert( o.is() ); + o = intersect( Segment_2(pts[2], pts[0]), Segment_2(pts[3], pts[1]) ); + assert( o.is() ); + // segment fully included + o = intersect( Segment_2(pts[0], pts[3]), Segment_2(pts[1], pts[2]) ); + assert( o.is() ); + assert( CGAL::object_cast(o) == Segment_2(pts[1], pts[2]) ); + // segment fully included with shared vertex + o = intersect( Segment_2(pts[0], pts[1]), Segment_2(pts[0], pts[2]) ); + assert( o.is() ); + assert( CGAL::object_cast(o) == Segment_2(pts[0], pts[1]) ); + // segment sharing a vertex + o = intersect( Segment_2(pts[0], pts[1]), Segment_2(pts[1], pts[2]) ); + assert( o.is() ); + assert( CGAL::object_cast(o) == pts[1]); + + // degenerate segment + Segment_2 sd(Point_2(1,0,6), Point_2(1,0,7)); + o = intersect( Segment_2(pts[0], pts[2]), sd ); + assert( o.is() ); + o = intersect( Segment_2(pts[0], pts[1]), sd ); + assert( o.is() ); + o = intersect( Segment_2(pts[1], pts[0]), sd ); + assert( o.is() ); + o = intersect( Segment_2(pts[2], pts[3]), sd ); + assert( o.empty() ); + + return true; +} + int main() { std::cerr.precision(17); @@ -75,5 +136,12 @@ int main() test > ("CDT_2 in a 3D plane, with Epeck"); + ok = ok && test_segment_intersections > + ("CDT_2 traits intersection with Epick"); + ok = ok && test_segment_intersections > + ("CDT_2 traits intersection with Epeck"); + return ok ? 0 : 1; }