diff --git a/BGL/include/CGAL/boost/graph/Euler_operations.h b/BGL/include/CGAL/boost/graph/Euler_operations.h index d575de0f17a..8d5268da6e0 100644 --- a/BGL/include/CGAL/boost/graph/Euler_operations.h +++ b/BGL/include/CGAL/boost/graph/Euler_operations.h @@ -13,6 +13,8 @@ #define CGAL_EULER_OPERATIONS_H #include +#include +#include #include #include @@ -558,6 +560,127 @@ add_edge(typename boost::graph_traits::vertex_descriptor s, return e; } + +/** +* checks whether a new face defined by a range of vertices (identified by their descriptors, +* `boost::graph_traits::%vertex_descriptor`) can be added. +*/ +template +bool can_add_face(const VertexRange& vrange, const PMesh& sm) +{ + typedef typename boost::graph_traits::vertex_descriptor vertex_descriptor; + typedef typename boost::graph_traits::halfedge_descriptor halfedge_descriptor; + + std::vector::vertex_descriptor> face(vrange.begin(), vrange.end()); + + std::size_t N = face.size(); + std::vector f2(face); + std::sort(f2.begin(), f2.end()); + + typename std::vector::iterator it = std::unique(f2.begin(),f2.end()); + + if((N > 0) && (it != f2.end())){ + return false; + } + + if(N < 3){ + return false; + } + + face.push_back(face.front()); + + for(std::size_t i=0; i < N; ++i){ + halfedge_descriptor hd; + bool found; + boost::tie(hd,found) = halfedge(face[i],face[i+1],sm); + if(found && (! is_border(hd,sm))){ + return false; + } + } + + for(std::size_t i=0; i < N; ++i){ + if(halfedge(face[i],sm) == boost::graph_traits::null_halfedge()){ + continue; + } + + if(! is_border(face[i],sm)){ + return false; + } + } + //Test if all halfedges of the new face + //are possibly consecutive border halfedges in the HDS. + //Possibly because it may be not directly encoded in the HDS + //(using next() function ). This situation can occur when one or + //more facets share only a vertex: For example, the new facet we try to add + //would make the vertex indices[i] a manifold but this should be forbidden + //if a facet only incident to that vertex has already been inserted. + //We check this for each vertex of the sequence. + for(std::size_t i = 0; i < N; ++i) { + std::size_t prev_index= (i-1+N)%N; + std::size_t next_index= (i+1)%N; + vertex_descriptor previous_vertex = face[ prev_index ]; + vertex_descriptor next_vertex = face[ next_index ]; + + halfedge_descriptor halfedge_around_vertex = halfedge(face[i],sm); + + if ( halfedge_around_vertex == boost::graph_traits::null_halfedge() || + halfedge(previous_vertex,sm) == boost::graph_traits::null_halfedge()|| + halfedge(next_vertex,sm) == boost::graph_traits::null_halfedge() + ) continue; + + halfedge_descriptor start=halfedge_around_vertex; + //halfedges pointing to/running out from vertex indices[i] + //and that need to be possibly consecutive + halfedge_descriptor prev_hd= boost::graph_traits::null_halfedge(),next_hd= boost::graph_traits::null_halfedge(); + + halfedge_around_vertex = opposite(next(halfedge_around_vertex,sm),sm); + //look for a halfedge incident to vertex indices[i] + //and which opposite is incident to previous_vertex + do{ + if(target(opposite(halfedge_around_vertex,sm),sm)==previous_vertex){ + prev_hd=halfedge_around_vertex; + CGAL_precondition(is_border(prev_hd,sm)); + break; + } + halfedge_around_vertex = opposite(next(halfedge_around_vertex,sm),sm); + } + while (halfedge_around_vertex!=start); + + if (prev_hd != boost::graph_traits::null_halfedge()){ + halfedge_around_vertex = opposite(next(halfedge_around_vertex,sm),sm); + //prev_hd and next are already consecutive in the HDS + if (target(opposite(halfedge_around_vertex,sm),sm)==next_vertex) continue; + + //look for a border halfedge which opposite is + //incident to next_vertex: set next halfedge + do + { + if (target(opposite(halfedge_around_vertex,sm),sm)==next_vertex){ + next_hd = opposite(halfedge_around_vertex,sm); + break; + } + halfedge_around_vertex = opposite(next(halfedge_around_vertex,sm),sm); + } + while(halfedge_around_vertex != prev_hd); + if (next_hd==boost::graph_traits::null_halfedge()) continue; + + //check if no constraint prevents + //prev_hd and next_hd to be adjacent: + do{ + halfedge_around_vertex = opposite(next(halfedge_around_vertex, sm),sm); + if ( is_border(opposite(halfedge_around_vertex,sm),sm) ) break; + } + while (halfedge_around_vertex != prev_hd); + if (halfedge_around_vertex == prev_hd) return false; + start = halfedge_around_vertex; + } + } + + return true; +} + + + /** * adds a new face defined by a range of vertices (identified by their descriptors, * `boost::graph_traits::%vertex_descriptor`). diff --git a/BGL/test/BGL/CMakeLists.txt b/BGL/test/BGL/CMakeLists.txt index 875cf10fd5d..dfe324a32af 100644 --- a/BGL/test/BGL/CMakeLists.txt +++ b/BGL/test/BGL/CMakeLists.txt @@ -95,6 +95,8 @@ create_single_source_cgal_program( "test_Face_filtered_graph.cpp" ) create_single_source_cgal_program( "test_Euler_operations.cpp" ) +create_single_source_cgal_program( "test_test_face.cpp" ) + create_single_source_cgal_program( "test_Collapse_edge.cpp" ) create_single_source_cgal_program( "test_graph_traits.cpp" ) diff --git a/BGL/test/BGL/test_test_face.cpp b/BGL/test/BGL/test_test_face.cpp new file mode 100644 index 00000000000..e69f7e3176e --- /dev/null +++ b/BGL/test/BGL/test_test_face.cpp @@ -0,0 +1,60 @@ +#include +#include +#include + +#include +#include + +typedef CGAL::Simple_cartesian K; +typedef K::Point_3 Point_3; +typedef CGAL::Surface_mesh SM; + +typedef boost::graph_traits::vertex_descriptor vertex_descriptor; +typedef boost::graph_traits::vertex_iterator vertex_iterator; +typedef std::vector V; +int main() +{ + { + SM sm; + + vertex_descriptor vp = CGAL::add_vertex(sm); + vertex_descriptor vq = CGAL::add_vertex(sm); + vertex_descriptor vr = CGAL::add_vertex(sm); + vertex_descriptor vs = CGAL::add_vertex(sm); + std::array face0; + assert( ! CGAL::Euler::can_add_face(face0,sm) ); + std::array face1; + assert( ! CGAL::Euler::can_add_face(face1,sm) ); + std::array face2; + assert( ! CGAL::Euler::can_add_face(face2,sm) ); + + std::array face = { vp, vq, vr }; + CGAL::Euler::add_face(face, sm); + + assert( ! CGAL::Euler::can_add_face(face,sm) ); + std::swap(face[0],face[1]); + assert( CGAL::Euler::can_add_face(face,sm) ); + + face[2] = vs; + assert( CGAL::Euler::can_add_face(face,sm) ); + std::swap(face[0],face[1]); + assert( ! CGAL::Euler::can_add_face(face,sm) ); + } + + { + SM sm; + Point_3 p(0,0,0), q(1,0,0), r(0,1,0), s(0,0,1); + CGAL::make_tetrahedron(p, q, r, s, sm); + std::array face; + vertex_iterator it = vertices(sm).first; + face[0] = *it; + ++it; + face[1] = *it; + face[2] = CGAL::add_vertex(sm); + assert( ! CGAL::Euler::can_add_face(face,sm) ); + std::swap(face[0],face[1]); + assert( ! CGAL::Euler::can_add_face(face,sm) ); + } + + return 0; +} diff --git a/Installation/CHANGES.md b/Installation/CHANGES.md index 80607cf638e..14d0ced92a4 100644 --- a/Installation/CHANGES.md +++ b/Installation/CHANGES.md @@ -9,6 +9,7 @@ Release date: December 2020 - Added the convenience header `CGAL/boost/graph/graph_traits_inheritance_macros.h` that allows to easily make any class inheriting from a model of a face graph concept, a model of the same concept. +- Added the function `can_add_face()`, which tests whether a new face defined by a range of vertices can be added. ### [3D Convex Hulls](https://doc.cgal.org/5.2/Manual/packages.html#PkgConvexHull3) - Added the function `CGAL::halfspace_intersection_interior_point_3()` that can be used to retrieve @@ -27,7 +28,7 @@ Release date: December 2020 to the make-x-monotone functions are fixed to use the new return type. - Changed `decompose()` interface to use `boost::variant` instead of legacy [`CGAL::Object`](https://doc.cgal.org/5.1/STL_Extension/classCGAL_1_1Object.html) - As exaplained above, the code is backward compatible. However, it is recommended + As explained above, the code is backward compatible. However, it is recommended that all calls to `decompose()` are fixed to use the new interface. ### [Polygon Mesh Processing](https://doc.cgal.org/5.2/Manual/packages.html#PkgPolygonMeshProcessing) diff --git a/Surface_mesh/include/CGAL/Surface_mesh/Surface_mesh.h b/Surface_mesh/include/CGAL/Surface_mesh/Surface_mesh.h index 4d056ae11e2..89f092155ac 100644 --- a/Surface_mesh/include/CGAL/Surface_mesh/Surface_mesh.h +++ b/Surface_mesh/include/CGAL/Surface_mesh/Surface_mesh.h @@ -9,11 +9,9 @@ // SPDX-License-Identifier: GPL-3.0-or-later OR LicenseRef-Commercial // - #ifndef CGAL_SURFACE_MESH_H #define CGAL_SURFACE_MESH_H - #include #include