From 06be9cf705f05205f6c8fc6e2bf84829ec3d67d9 Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?S=C3=A9bastien=20Loriot?= Date: Mon, 11 May 2015 22:37:29 +0200 Subject: [PATCH] bug-fix: fix the orientation of the convex polyhedron --- .../CGAL/Convex_hull_3/dual/halfspace_intersection_3.h | 6 +++--- .../dual/halfspace_intersection_with_constructions_3.h | 2 +- .../test/Convex_hull_3/test_halfspace_intersections.cpp | 9 +++++++++ 3 files changed, 13 insertions(+), 4 deletions(-) diff --git a/Convex_hull_3/include/CGAL/Convex_hull_3/dual/halfspace_intersection_3.h b/Convex_hull_3/include/CGAL/Convex_hull_3/dual/halfspace_intersection_3.h index 302e6c791c5..1f9061d44c2 100644 --- a/Convex_hull_3/include/CGAL/Convex_hull_3/dual/halfspace_intersection_3.h +++ b/Convex_hull_3/include/CGAL/Convex_hull_3/dual/halfspace_intersection_3.h @@ -142,7 +142,7 @@ namespace CGAL B.begin_facet(); do { B.add_vertex_to_facet(primal_vertices[hf->facet()]); - } while (++hf != h0); + } while (--hf != h0); B.end_facet(); } @@ -269,7 +269,7 @@ namespace CGAL P.delegate(build_primal); // Posterior check if the origin is inside the computed polyhedron - CGAL_assertion_msg(!Convex_hull_3::internal::point_inside_convex_polyhedron(P, p_origin), "halfspace_intersection_3: origin not in the polyhedron"); + CGAL_assertion_msg(Convex_hull_3::internal::point_inside_convex_polyhedron(P, p_origin), "halfspace_intersection_3: origin not in the polyhedron"); } else { // choose exact integral type #ifdef CGAL_USE_GMP @@ -293,7 +293,7 @@ namespace CGAL P.delegate(build_primal); // Posterior check if the origin is inside the computed polyhedron - CGAL_assertion_msg(!Convex_hull_3::internal::point_inside_convex_polyhedron(P, origin), "halfspace_intersection_3: origin not in the polyhedron"); + CGAL_assertion_msg(Convex_hull_3::internal::point_inside_convex_polyhedron(P, origin), "halfspace_intersection_3: origin not in the polyhedron"); } } diff --git a/Convex_hull_3/include/CGAL/Convex_hull_3/dual/halfspace_intersection_with_constructions_3.h b/Convex_hull_3/include/CGAL/Convex_hull_3/dual/halfspace_intersection_with_constructions_3.h index c91a5a8f198..628367860dd 100644 --- a/Convex_hull_3/include/CGAL/Convex_hull_3/dual/halfspace_intersection_with_constructions_3.h +++ b/Convex_hull_3/include/CGAL/Convex_hull_3/dual/halfspace_intersection_with_constructions_3.h @@ -100,7 +100,7 @@ namespace CGAL do { B.add_vertex_to_facet(extreme_points[hf->facet()]); - } while (++hf != h0); + } while (--hf != h0); B.end_facet(); } diff --git a/Convex_hull_3/test/Convex_hull_3/test_halfspace_intersections.cpp b/Convex_hull_3/test/Convex_hull_3/test_halfspace_intersections.cpp index 3f12ee09e40..4a851b6fe29 100644 --- a/Convex_hull_3/test/Convex_hull_3/test_halfspace_intersections.cpp +++ b/Convex_hull_3/test/Convex_hull_3/test_halfspace_intersections.cpp @@ -4,6 +4,8 @@ #include #include +#include + typedef CGAL::Exact_predicates_inexact_constructions_kernel K; typedef K::Plane_3 Plane; @@ -59,6 +61,13 @@ int main (void) { assert(P4.size_of_vertices()==8 && P4.size_of_facets()==6); assert(P5.size_of_vertices()==8 && P5.size_of_facets()==6); + using CGAL::Convex_hull_3::internal::point_inside_convex_polyhedron; + assert(point_inside_convex_polyhedron(P1, Point(0,0,0))); + assert(point_inside_convex_polyhedron(P2, Point(0,0,0))); + assert(point_inside_convex_polyhedron(P3, Point(0,0,0))); + assert(point_inside_convex_polyhedron(P4, Point(0,0,0))); + assert(point_inside_convex_polyhedron(P5, Point(0,0,0))); + return 0; }