From 16782d1e8b4539365d5dbea87fab3f2d025217d7 Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?S=C3=A9bastien=20Loriot?= Date: Mon, 8 Apr 2024 15:46:48 +0200 Subject: [PATCH 1/4] do not depend on vertex order when opposite angles are identical --- .../triangulate_faces.h | 40 ++++++++++++++----- 1 file changed, 29 insertions(+), 11 deletions(-) diff --git a/Polygon_mesh_processing/include/CGAL/Polygon_mesh_processing/triangulate_faces.h b/Polygon_mesh_processing/include/CGAL/Polygon_mesh_processing/triangulate_faces.h index 028e728e89e4..7d71c2e63e8b 100644 --- a/Polygon_mesh_processing/include/CGAL/Polygon_mesh_processing/triangulate_faces.h +++ b/Polygon_mesh_processing/include/CGAL/Polygon_mesh_processing/triangulate_faces.h @@ -196,15 +196,15 @@ class Triangulate_polygon_mesh_modifier if(original_size == 4) { - halfedge_descriptor v0, v1, v2, v3; - v0 = halfedge(f, pmesh); - Point_ref p0 = get(vpm, target(v0, pmesh)); - v1 = next(v0, pmesh); - Point_ref p1 = get(vpm, target(v1, pmesh)); - v2 = next(v1, pmesh); - Point_ref p2 = get(vpm, target(v2, pmesh)); - v3 = next(v2, pmesh); - Point_ref p3 = get(vpm, target(v3, pmesh)); + std::array verts; + verts[0] = halfedge(f, pmesh); + verts[1] = next(verts[0], pmesh); + verts[2] = next(verts[1], pmesh); + verts[3] = next(verts[2], pmesh); + Point_ref p0 = get(vpm, target(verts[0], pmesh)); + Point_ref p1 = get(vpm, target(verts[1], pmesh)); + Point_ref p2 = get(vpm, target(verts[2], pmesh)); + Point_ref p3 = get(vpm, target(verts[3], pmesh)); /* Chooses the diagonal that will split the quad in two triangles that maximize * the scalar product of the un-normalized normals of the two triangles. @@ -219,8 +219,26 @@ class Triangulate_polygon_mesh_modifier const FT p1p3 = cross_product(p2-p1, p3-p2) * cross_product(p0-p3, p1-p0); const FT p0p2 = cross_product(p1-p0, p1-p2) * cross_product(p3-p2, p3-p0); - halfedge_descriptor res = (p0p2>p1p3) ? CGAL::Euler::split_face(v0, v2, pmesh) - : CGAL::Euler::split_face(v1, v3, pmesh); + + halfedge_descriptor res = boost::graph_traits::null_halfedge(); + + if (p1p3!=p0p2) + res = (p0p2>p1p3) + ? CGAL::Euler::split_face(verts[0], verts[2], pmesh) + : CGAL::Euler::split_face(verts[1], verts[3], pmesh); + else + { + halfedge_descriptor m = + *(std::min_element)(verts.begin(), verts.end(), + [&pmesh,vpm](halfedge_descriptor v0, halfedge_descriptor v1) + { + return lexicographically_xyz_smaller(get(vpm, target(v0, pmesh)), + get(vpm, target(v1, pmesh))); + }); + res = (m==verts[0] || m==verts[2]) + ? CGAL::Euler::split_face(verts[0], verts[2], pmesh) + : CGAL::Euler::split_face(verts[1], verts[3], pmesh); + } visitor.after_subface_created(face(res, pmesh)); visitor.after_subface_created(face(opposite(res, pmesh), pmesh)); From ebe6d6719cb3b9941c6eb800c83f2928427511e5 Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?S=C3=A9bastien=20Loriot?= Date: Thu, 11 Apr 2024 08:55:10 +0200 Subject: [PATCH 2/4] precompute vectors --- .../Polygon_mesh_processing/triangulate_faces.h | 13 +++++++++---- 1 file changed, 9 insertions(+), 4 deletions(-) diff --git a/Polygon_mesh_processing/include/CGAL/Polygon_mesh_processing/triangulate_faces.h b/Polygon_mesh_processing/include/CGAL/Polygon_mesh_processing/triangulate_faces.h index 7d71c2e63e8b..230e2bf57445 100644 --- a/Polygon_mesh_processing/include/CGAL/Polygon_mesh_processing/triangulate_faces.h +++ b/Polygon_mesh_processing/include/CGAL/Polygon_mesh_processing/triangulate_faces.h @@ -217,13 +217,18 @@ class Triangulate_polygon_mesh_modifier */ visitor.before_subface_creations(f); - const FT p1p3 = cross_product(p2-p1, p3-p2) * cross_product(p0-p3, p1-p0); - const FT p0p2 = cross_product(p1-p0, p1-p2) * cross_product(p3-p2, p3-p0); + typename Traits::Vector_3 p0p1=p1-p0; + typename Traits::Vector_3 p1p2=p2-p1; + typename Traits::Vector_3 p2p3=p3-p2; + typename Traits::Vector_3 p0p3=p3-p0; + + const FT delta1 = cross_product(p1p2, p2p3) * cross_product(-p0p3, p0p1); + const FT delta2 = cross_product(p0p1, -p1p2) * cross_product(p2p3, p0p3); halfedge_descriptor res = boost::graph_traits::null_halfedge(); - if (p1p3!=p0p2) - res = (p0p2>p1p3) + if (delta1!=delta2) + res = (delta2>delta1) ? CGAL::Euler::split_face(verts[0], verts[2], pmesh) : CGAL::Euler::split_face(verts[1], verts[3], pmesh); else From 147b313cc860d1e6c1b8cee65f95a0e3bc2f7611 Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?S=C3=A9bastien=20Loriot?= Date: Thu, 11 Apr 2024 09:36:16 +0200 Subject: [PATCH 3/4] also apply deterministic fix to polygon soup function --- .../triangulate_faces.h | 44 ++++++++++++++++--- 1 file changed, 37 insertions(+), 7 deletions(-) diff --git a/Polygon_mesh_processing/include/CGAL/Polygon_mesh_processing/triangulate_faces.h b/Polygon_mesh_processing/include/CGAL/Polygon_mesh_processing/triangulate_faces.h index 230e2bf57445..a2d234a563e7 100644 --- a/Polygon_mesh_processing/include/CGAL/Polygon_mesh_processing/triangulate_faces.h +++ b/Polygon_mesh_processing/include/CGAL/Polygon_mesh_processing/triangulate_faces.h @@ -549,6 +549,7 @@ class Triangulate_polygon_soup_modifier NamedParameters, Def_Kernel>::type; using FT = typename Traits::FT; + using PID = typename std::iterator_traits::value_type; // Visitor using Visitor = typename internal_np::Lookup_named_param_def< @@ -587,17 +588,46 @@ class Triangulate_polygon_soup_modifier */ visitor.before_subface_creations(polygon); - const FT p1p3 = cross_product(p2-p1, p3-p2) * cross_product(p0-p3, p1-p0); - const FT p0p2 = cross_product(p1-p0, p1-p2) * cross_product(p3-p2, p3-p0); - if(p0p2 > p1p3) + + typename Traits::Vector_3 p0p1=p1-p0; + typename Traits::Vector_3 p1p2=p2-p1; + typename Traits::Vector_3 p2p3=p3-p2; + typename Traits::Vector_3 p0p3=p3-p0; + + const FT delta1 = cross_product(p1p2, p2p3) * cross_product(-p0p3, p0p1); + const FT delta2 = cross_product(p0p1, -p1p2) * cross_product(p2p3, p0p3); + if (delta1!=delta2) { - triangulated_polygons.push_back({polygon[0], polygon[1], polygon[2]}); - triangulated_polygons.push_back({polygon[0], polygon[2], polygon[3]}); + if(delta2 > delta1) + { + triangulated_polygons.push_back({polygon[0], polygon[1], polygon[2]}); + triangulated_polygons.push_back({polygon[0], polygon[2], polygon[3]}); + } + else + { + triangulated_polygons.push_back({polygon[0], polygon[1], polygon[3]}); + triangulated_polygons.push_back({polygon[1], polygon[2], polygon[3]}); + } } else { - triangulated_polygons.push_back({polygon[0], polygon[1], polygon[3]}); - triangulated_polygons.push_back({polygon[1], polygon[2], polygon[3]}); + PID mid = + *(std::min_element)(polygon.begin(), polygon.end(), + [&points,pm](PID id1 , PID id2) + { + return lexicographically_xyz_smaller(get(pm, points[id1]), + get(pm, points[id2])); + }); + if (mid==0|| mid==2) + { + triangulated_polygons.push_back({polygon[0], polygon[1], polygon[2]}); + triangulated_polygons.push_back({polygon[0], polygon[2], polygon[3]}); + } + else + { + triangulated_polygons.push_back({polygon[0], polygon[1], polygon[3]}); + triangulated_polygons.push_back({polygon[1], polygon[2], polygon[3]}); + } } visitor.after_subface_created(triangulated_polygons[triangulated_polygons.size()-2]); From cfec3bb0cceebed48231370839af8991ea36e051 Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?S=C3=A9bastien=20Loriot?= Date: Thu, 11 Apr 2024 11:44:35 +0200 Subject: [PATCH 4/4] fix values due to the diagonal change in degenerate triangulation --- ...meshing_with_isolated_constraints_test.cpp | 2 +- .../Polygon_mesh_processing/test_pmp_clip.cpp | 28 +++++++++---------- 2 files changed, 15 insertions(+), 15 deletions(-) diff --git a/Polygon_mesh_processing/test/Polygon_mesh_processing/remeshing_with_isolated_constraints_test.cpp b/Polygon_mesh_processing/test/Polygon_mesh_processing/remeshing_with_isolated_constraints_test.cpp index 4e7ebe960925..c36a1da56969 100644 --- a/Polygon_mesh_processing/test/Polygon_mesh_processing/remeshing_with_isolated_constraints_test.cpp +++ b/Polygon_mesh_processing/test/Polygon_mesh_processing/remeshing_with_isolated_constraints_test.cpp @@ -30,7 +30,7 @@ int main(int, char**) std::set fs; auto selected_faces = make_boolean_property_map(fs); - fs.insert(*(faces(sm).begin())); + fs.insert(Polygon_mesh::Face_index(190)); CGAL::expand_face_selection(fs, sm, 1, selected_faces, CGAL::Emptyset_iterator()); std::cout << fs.size() << " faces in the range" << std::endl; assert(fs.size() == 4); diff --git a/Polygon_mesh_processing/test/Polygon_mesh_processing/test_pmp_clip.cpp b/Polygon_mesh_processing/test/Polygon_mesh_processing/test_pmp_clip.cpp index c564403312a0..c6c1bb1aeea6 100644 --- a/Polygon_mesh_processing/test/Polygon_mesh_processing/test_pmp_clip.cpp +++ b/Polygon_mesh_processing/test/Polygon_mesh_processing/test_pmp_clip.cpp @@ -393,31 +393,31 @@ void test() // -> closed mesh, true/true PMP::clip(tm1, K::Plane_3(-1,0,0,0), params::use_compact_clipper(true).clip_volume(true)); - assert(faces(tm1).size() == 12); + assert(faces(tm1).size() == 14); assert(CGAL::is_closed(tm1)); // -> closed mesh, false/true PMP::clip(tm1, K::Plane_3(-1,0,0,0), params::use_compact_clipper(false).clip_volume(true)); - assert(faces(tm1).size() == 12); + assert(faces(tm1).size() == 14); assert(CGAL::is_closed(tm1)); // -> closed mesh, true/false PMP::clip(tm1, K::Plane_3(-1,0,0,0), params::use_compact_clipper(true).clip_volume(false)); - assert(faces(tm1).size() == 12); + assert(faces(tm1).size() == 14); assert(CGAL::is_closed(tm1)); // -> closed mesh, false/false PMP::clip(tm1, K::Plane_3(1,0,0,-1), params::use_compact_clipper(false).clip_volume(false)); - assert(faces(tm1).size() == 10); + assert(faces(tm1).size() == 12); assert(!CGAL::is_closed(tm1)); // -> open mesh true/true PMP::clip(tm1, K::Plane_3(-1,0,0,0), params::use_compact_clipper(true).clip_volume(true)); - assert(faces(tm1).size() == 10); + assert(faces(tm1).size() == 12); // -> open mesh true/false PMP::clip(tm1, K::Plane_3(-1,0,0,0), params::use_compact_clipper(true).clip_volume(false)); - assert(faces(tm1).size() == 10); + assert(faces(tm1).size() == 12); // -> open mesh false/false PMP::clip(tm1, K::Plane_3(-1,0,0,0), params::use_compact_clipper(false).clip_volume(false)); @@ -539,7 +539,7 @@ void test() TriangleMesh tm1; std::ifstream("data-coref/open_large_cube.off") >> tm1; PMP::clip(tm1, K::Plane_3(0,0,-1,1), CGAL::parameters::use_compact_clipper(true)); - assert(vertices(tm1).size()==176); + assert(vertices(tm1).size()==178); } { @@ -553,7 +553,7 @@ void test() TriangleMesh tm1; std::ifstream("data-coref/open_large_cube.off") >> tm1; PMP::clip(tm1, K::Plane_3(0,0,-1,1), CGAL::parameters::use_compact_clipper(true).allow_self_intersections(true)); - assert(vertices(tm1).size()==176); + assert(vertices(tm1).size()==178); } } @@ -802,11 +802,11 @@ void test_isocuboid() assert(meshes.size() == 10); //if the order is not deterministc, put the num_vertices in a list and check //if the list does contain all those numbers. - assert(num_vertices(meshes[0]) == 2657); + assert(num_vertices(meshes[0]) == 2663); assert(num_vertices(meshes[1]) == 131 ); assert(num_vertices(meshes[2]) == 32 ); - assert(num_vertices(meshes[3]) == 123 ); - assert(num_vertices(meshes[4]) == 220 ); + assert(num_vertices(meshes[3]) == 125 ); + assert(num_vertices(meshes[4]) == 224 ); assert(num_vertices(meshes[5]) == 107 ); assert(num_vertices(meshes[6]) == 121 ); assert(num_vertices(meshes[7]) == 56 ); @@ -836,8 +836,8 @@ void test_isocuboid() for (int i=0; i<4; ++i) sizes.insert(vertices(meshes[i]).size()); - assert(sizes.count(22)==1); - assert(sizes.count(23)==1); + assert(sizes.count(20)==1); + assert(sizes.count(21)==1); assert(sizes.count(7)==1); assert(sizes.count(4)==1); @@ -852,7 +852,7 @@ void test_isocuboid() assert(meshes.size() == 2); //if the order is not deterministc, put the num_vertices in a list and check //if the list does contain all those numbers. - assert(vertices(meshes[0]).size() == 22); + assert(vertices(meshes[0]).size() == 20); assert(vertices(meshes[1]).size() == 4); } int main()