Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

Sign up for GitHub

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Make triangulate_faces deterministic #8129

Merged
merged 4 commits into from
Apr 22, 2024
Merged

Make triangulate_faces deterministic #8129

Show file tree
Hide file tree
Changes from all commits
Commits
Show all changes
4 commits
Select commit Hold shift + click to select a range
16782d1
do not depend on vertex order when opposite angles are identical
sloriot Apr 8, 2024
ebe6d67
precompute vectors
sloriot Apr 11, 2024
147b313
also apply deterministic fix to polygon soup function
sloriot Apr 11, 2024
cfec3bb
fix values due to the diagonal change in degenerate triangulation
sloriot Apr 11, 2024
File filter

Filter by extension

Filter by extension
Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
93 changes: 73 additions & 20 deletions Polygon_mesh_processing/include/CGAL/Polygon_mesh_processing/triangulate_faces.h
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -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<halfedge_descriptor,4> 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.
Expand All @@ -217,10 +217,33 @@ 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);
halfedge_descriptor res = (p0p2>p1p3) ? CGAL::Euler::split_face(v0, v2, pmesh)
: CGAL::Euler::split_face(v1, v3, pmesh);
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<PolygonMesh>::null_halfedge();

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
{
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));
Expand Down Expand Up @@ -526,6 +549,7 @@ class Triangulate_polygon_soup_modifier
NamedParameters,
Def_Kernel>::type;
using FT = typename Traits::FT;
using PID = typename std::iterator_traits<typename Polygon::const_iterator>::value_type;

// Visitor
using Visitor = typename internal_np::Lookup_named_param_def<
Expand Down Expand Up @@ -564,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]);
Expand Down
2 changes: 1 addition & 1 deletion ...mesh_processing/test/Polygon_mesh_processing/remeshing_with_isolated_constraints_test.cpp
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -30,7 +30,7 @@ int main(int, char**)

std::set<face_descriptor> 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);
Expand Down
28 changes: 14 additions & 14 deletions Polygon_mesh_processing/test/Polygon_mesh_processing/test_pmp_clip.cpp
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -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));
Expand Down Expand Up @@ -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);
}

{
Expand All @@ -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);
}
}

Expand Down Expand Up @@ -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 );
Expand Down Expand Up @@ -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);

Expand All @@ -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()
Expand Down
Loading