diff --git a/CHANGELOG.md b/CHANGELOG.md index a805765..f841c45 100644 --- a/CHANGELOG.md +++ b/CHANGELOG.md @@ -9,6 +9,13 @@ and this project adheres to [Semantic Versioning](https://semver.org/spec/v2.0.0 ### Added +Add missing igl functions to enable compas_slicer to drop the `libigl` dependency entirely: + +- `cotmatrix` module: `trimesh_cotmatrix`, `trimesh_cotmatrix_entries` +- `grad` module: `trimesh_grad` (gradient operator) +- `meshing` module: `trimesh_cut_mesh`, `trimesh_face_components` +- `simplify` module: `ramer_douglas_peucker` (polyline simplification) + ### Changed ### Removed diff --git a/CMakeLists.txt b/CMakeLists.txt index 8910f00..b58dec0 100644 --- a/CMakeLists.txt +++ b/CMakeLists.txt @@ -151,8 +151,10 @@ endfunction() # === Add your modules === add_nanobind_module(_boundaries src/boundaries.cpp) +add_nanobind_module(_cotmatrix src/cotmatrix.cpp) add_nanobind_module(_curvature src/curvature.cpp) add_nanobind_module(_geodistance src/geodistance.cpp) +add_nanobind_module(_grad src/grad.cpp) add_nanobind_module(_intersections src/intersections.cpp) add_nanobind_module(_isolines src/isolines.cpp) add_nanobind_module(_mapping src/mapping.cpp) @@ -160,4 +162,5 @@ add_nanobind_module(_meshing src/meshing.cpp) add_nanobind_module(_massmatrix src/massmatrix.cpp) add_nanobind_module(_parametrisation src/parametrisation.cpp) add_nanobind_module(_planarize src/planarize.cpp) +add_nanobind_module(_simplify src/simplify.cpp) diff --git a/src/compas_libigl/__init__.py b/src/compas_libigl/__init__.py index c434ba4..8e7997a 100644 --- a/src/compas_libigl/__init__.py +++ b/src/compas_libigl/__init__.py @@ -6,11 +6,14 @@ __all_plugins__ = [ + "compas_libigl.cotmatrix", "compas_libigl.geodistance", + "compas_libigl.grad", "compas_libigl.intersections", "compas_libigl.isolines", "compas_libigl.massmatrix", + "compas_libigl.meshing", "compas_libigl.parametrisation", "compas_libigl.planarize", - "compas_libigl.meshing", + "compas_libigl.simplify", ] diff --git a/src/compas_libigl/cotmatrix.py b/src/compas_libigl/cotmatrix.py new file mode 100644 index 0000000..f405621 --- /dev/null +++ b/src/compas_libigl/cotmatrix.py @@ -0,0 +1,47 @@ +import numpy as np +from compas.plugins import plugin + +from compas_libigl import _cotmatrix + + +@plugin(category="trimesh") +def trimesh_cotmatrix(M): + """Compute the cotangent Laplacian matrix of a triangle mesh. + + Parameters + ---------- + M : tuple[list[list[float]], list[list[int]]] + A mesh represented by a tuple of (vertices, faces) + where vertices are 3D points and faces are triangles + + Returns + ------- + scipy.sparse.csc_matrix + The cotangent Laplacian matrix in sparse format. + """ + V, F = M + V = np.asarray(V, dtype=np.float64) + F = np.asarray(F, dtype=np.int32) + return _cotmatrix.trimesh_cotmatrix(V, F) + + +@plugin(category="trimesh") +def trimesh_cotmatrix_entries(M): + """Compute cotangent values for each edge in each triangle. + + Parameters + ---------- + M : tuple[list[list[float]], list[list[int]]] + A mesh represented by a tuple of (vertices, faces) + where vertices are 3D points and faces are triangles + + Returns + ------- + numpy.ndarray + A matrix of shape (F, 3) containing cotangent values per edge. + For each face, contains cotan of angles opposite to edges (1,2), (2,0), (0,1). + """ + V, F = M + V = np.asarray(V, dtype=np.float64) + F = np.asarray(F, dtype=np.int32) + return _cotmatrix.trimesh_cotmatrix_entries(V, F) diff --git a/src/compas_libigl/grad.py b/src/compas_libigl/grad.py new file mode 100644 index 0000000..14d2a6b --- /dev/null +++ b/src/compas_libigl/grad.py @@ -0,0 +1,27 @@ +import numpy as np +from compas.plugins import plugin + +from compas_libigl import _grad + + +@plugin(category="trimesh") +def trimesh_grad(M): + """Compute the gradient operator for a triangle mesh. + + Parameters + ---------- + M : tuple[list[list[float]], list[list[int]]] + A mesh represented by a tuple of (vertices, faces) + where vertices are 3D points and faces are triangles + + Returns + ------- + scipy.sparse.csc_matrix + The gradient operator as a sparse matrix of size (3*F, V). + When multiplied by a scalar field on vertices, produces gradient + vectors per face (stacked as Gx, Gy, Gz). + """ + V, F = M + V = np.asarray(V, dtype=np.float64) + F = np.asarray(F, dtype=np.int32) + return _grad.trimesh_grad(V, F) diff --git a/src/compas_libigl/meshing.py b/src/compas_libigl/meshing.py index f25b7a8..7134ac4 100644 --- a/src/compas_libigl/meshing.py +++ b/src/compas_libigl/meshing.py @@ -65,3 +65,53 @@ def trimesh_remesh_along_isolines(M, scalars, isovalues): ISO = np.asarray(isovalues, dtype=np.float64) V2, F2, S2, G2 = _meshing.trimesh_remesh_along_isolines(V, F, S, ISO) return V2.tolist(), F2.tolist(), S2.tolist(), G2.tolist() + + +@plugin(category="trimesh") +def trimesh_cut_mesh(M, cuts): + """Cut a mesh along specified edges. + + Parameters + ---------- + M : tuple[list[list[float]], list[list[int]]] + A mesh represented by a tuple of (vertices, faces) + where vertices are 3D points and faces are triangles + cuts : list[list[int]] + A matrix of shape (F, 3) with boolean flags indicating + which edges to cut. For each face, the three values + correspond to edges (v0,v1), (v1,v2), (v2,v0). + + Returns + ------- + tuple[list[list[float]], list[list[int]]] + A tuple containing + * the vertices of the cut mesh (with duplicated vertices along cuts), + * the faces of the cut mesh. + """ + V, F = M + V = np.asarray(V, dtype=np.float64) + F = np.asarray(F, dtype=np.int32) + C = np.asarray(cuts, dtype=np.int32) + Vn, Fn = _meshing.trimesh_cut_mesh(V, F, C) + return Vn.tolist(), Fn.tolist() + + +@plugin(category="trimesh") +def trimesh_face_components(M): + """Compute connected components of mesh faces. + + Parameters + ---------- + M : tuple[list[list[float]], list[list[int]]] + A mesh represented by a tuple of (vertices, faces) + where vertices are 3D points and faces are triangles + + Returns + ------- + list[int] + Component ID per face. Faces sharing edges belong to the same component. + """ + V, F = M + F = np.asarray(F, dtype=np.int32) + C = _meshing.trimesh_face_components(F) + return C.tolist() diff --git a/src/compas_libigl/simplify.py b/src/compas_libigl/simplify.py new file mode 100644 index 0000000..62494d6 --- /dev/null +++ b/src/compas_libigl/simplify.py @@ -0,0 +1,29 @@ +import numpy as np +from compas.plugins import plugin + +from compas_libigl import _simplify + + +@plugin(category="polyline") +def ramer_douglas_peucker(points, threshold): + """Simplify a polyline using Ramer-Douglas-Peucker algorithm. + + Parameters + ---------- + points : list[list[float]] + A list of 3D points representing the polyline. + threshold : float + Maximum distance threshold for simplification. + Points that deviate less than this from the simplified line are removed. + + Returns + ------- + tuple[list[list[float]], list[int], list[list[float]]] + A tuple containing + * the simplified polyline points (S), + * indices in original polyline that were kept (J), + * for each original point, the corresponding point on the simplified curve (Q). + """ + P = np.asarray(points, dtype=np.float64) + S, J, Q = _simplify.ramer_douglas_peucker(P, float(threshold)) + return S.tolist(), J.tolist(), Q.tolist() diff --git a/src/cotmatrix.cpp b/src/cotmatrix.cpp new file mode 100644 index 0000000..008709f --- /dev/null +++ b/src/cotmatrix.cpp @@ -0,0 +1,39 @@ +#include "cotmatrix.hpp" + +Eigen::SparseMatrix +trimesh_cotmatrix( + Eigen::Ref V, + Eigen::Ref F +) { + Eigen::SparseMatrix L; + igl::cotmatrix(V, F, L); + return L; +} + +compas::RowMatrixXd +trimesh_cotmatrix_entries( + Eigen::Ref V, + Eigen::Ref F +) { + Eigen::MatrixXd C; + igl::cotmatrix_entries(V, F, C); + compas::RowMatrixXd C_row = C; + return C_row; +} + +NB_MODULE(_cotmatrix, m) { + + m.def( + "trimesh_cotmatrix", + &trimesh_cotmatrix, + "Compute the cotangent Laplacian matrix for a triangle mesh.", + "V"_a, "F"_a + ); + + m.def( + "trimesh_cotmatrix_entries", + &trimesh_cotmatrix_entries, + "Compute cotangent values for each edge in each triangle.", + "V"_a, "F"_a + ); +} diff --git a/src/cotmatrix.hpp b/src/cotmatrix.hpp new file mode 100644 index 0000000..24e6387 --- /dev/null +++ b/src/cotmatrix.hpp @@ -0,0 +1,31 @@ +#pragma once + +#include "compas.hpp" +#include +#include +#include +#include + +/** + * Compute the cotangent Laplacian matrix of a triangle mesh. + * + * @param V The vertex positions of the mesh. + * @param F The face indices of the mesh. + * @return The cotangent Laplacian as a sparse matrix. + */ +Eigen::SparseMatrix trimesh_cotmatrix( + Eigen::Ref V, + Eigen::Ref F +); + +/** + * Compute the cotangent values for each edge in each triangle of a mesh. + * + * @param V The vertex positions of the mesh. + * @param F The face indices of the mesh. + * @return A matrix of size F x 3 containing cotangent values per edge. + */ +compas::RowMatrixXd trimesh_cotmatrix_entries( + Eigen::Ref V, + Eigen::Ref F +); diff --git a/src/grad.cpp b/src/grad.cpp new file mode 100644 index 0000000..571115e --- /dev/null +++ b/src/grad.cpp @@ -0,0 +1,21 @@ +#include "grad.hpp" + +Eigen::SparseMatrix +trimesh_grad( + Eigen::Ref V, + Eigen::Ref F +) { + Eigen::SparseMatrix G; + igl::grad(V, F, G); + return G; +} + +NB_MODULE(_grad, m) { + + m.def( + "trimesh_grad", + &trimesh_grad, + "Compute the gradient operator for a triangle mesh.", + "V"_a, "F"_a + ); +} diff --git a/src/grad.hpp b/src/grad.hpp new file mode 100644 index 0000000..0035d3a --- /dev/null +++ b/src/grad.hpp @@ -0,0 +1,18 @@ +#pragma once + +#include "compas.hpp" +#include +#include +#include + +/** + * Compute the gradient operator for a triangle mesh. + * + * @param V The vertex positions of the mesh. + * @param F The face indices of the mesh. + * @return The gradient operator as a sparse matrix of size 3F x V. + */ +Eigen::SparseMatrix trimesh_grad( + Eigen::Ref V, + Eigen::Ref F +); diff --git a/src/meshing.cpp b/src/meshing.cpp index 9f656d8..f0d4b2e 100644 --- a/src/meshing.cpp +++ b/src/meshing.cpp @@ -89,6 +89,28 @@ trimesh_remesh_along_isolines( return std::make_tuple(V, F, S, face_groups); } +std::tuple +trimesh_cut_mesh( + Eigen::Ref V, + Eigen::Ref F, + Eigen::Ref cuts +) { + compas::RowMatrixXd Vn; + compas::RowMatrixXi Fn; + Eigen::VectorXi I; // birth vertices (original vertex indices) + igl::cut_mesh(V, F, cuts, Vn, Fn, I); + return std::make_tuple(Vn, Fn); +} + +Eigen::VectorXi +trimesh_face_components( + Eigen::Ref F +) { + Eigen::VectorXi C; + igl::facet_components(F, C); + return C; +} + NB_MODULE(_meshing, m) { m.doc() = "Mesh remeshing functions using libigl"; @@ -100,7 +122,7 @@ NB_MODULE(_meshing, m) { "F1"_a, "S1"_a, "s"_a); - + m.def( "trimesh_remesh_along_isolines", &trimesh_remesh_along_isolines, @@ -109,4 +131,18 @@ NB_MODULE(_meshing, m) { "F1"_a, "S1"_a, "values"_a); + + m.def( + "trimesh_cut_mesh", + &trimesh_cut_mesh, + "Cut a mesh along specified edges, duplicating vertices.", + "V"_a, + "F"_a, + "cuts"_a); + + m.def( + "trimesh_face_components", + &trimesh_face_components, + "Compute connected components of mesh faces.", + "F"_a); } diff --git a/src/meshing.hpp b/src/meshing.hpp index 0e34ec5..28e6d45 100644 --- a/src/meshing.hpp +++ b/src/meshing.hpp @@ -3,6 +3,7 @@ #include "compas.hpp" #include #include +#include #include #include @@ -37,3 +38,27 @@ trimesh_remesh_along_isolines( Eigen::Ref F_initial, Eigen::Ref S_initial, Eigen::Ref values); + +/** + * Cut a mesh along specified edges. + * + * @param V #V x 3 matrix of vertex coordinates + * @param F #F x 3 matrix of triangle indices + * @param cuts #F x 3 matrix of boolean cut flags per edge + * @return Tuple of (vertices, faces) with vertices duplicated along cuts + */ +std::tuple +trimesh_cut_mesh( + Eigen::Ref V, + Eigen::Ref F, + Eigen::Ref cuts); + +/** + * Compute connected components of mesh faces. + * + * @param F #F x 3 matrix of triangle indices + * @return Vector of component IDs per face + */ +Eigen::VectorXi +trimesh_face_components( + Eigen::Ref F); diff --git a/src/simplify.cpp b/src/simplify.cpp new file mode 100644 index 0000000..66f1f15 --- /dev/null +++ b/src/simplify.cpp @@ -0,0 +1,28 @@ +#include "simplify.hpp" + +std::tuple +ramer_douglas_peucker( + Eigen::Ref P, + double threshold +) { + Eigen::MatrixXd S; // Simplified polyline + Eigen::VectorXi J; // Indices of kept points in original polyline + Eigen::MatrixXd Q; // For each point in original P, the corresponding point on simplified curve + + igl::ramer_douglas_peucker(P, threshold, S, J, Q); + + compas::RowMatrixXd S_row = S; + compas::RowMatrixXd Q_row = Q; + return std::make_tuple(S_row, J, Q_row); +} + +NB_MODULE(_simplify, m) { + m.doc() = "Polyline simplification functions using libigl"; + + m.def( + "ramer_douglas_peucker", + &ramer_douglas_peucker, + "Simplify a polyline using Ramer-Douglas-Peucker algorithm.", + "P"_a, + "threshold"_a); +} diff --git a/src/simplify.hpp b/src/simplify.hpp new file mode 100644 index 0000000..6655dd2 --- /dev/null +++ b/src/simplify.hpp @@ -0,0 +1,17 @@ +#pragma once + +#include "compas.hpp" +#include +#include + +/** + * Simplify a polyline using Ramer-Douglas-Peucker algorithm. + * + * @param P #P x 3 matrix of polyline points + * @param threshold Maximum distance threshold for simplification + * @return Tuple of (simplified points, indices of kept points, mapping from original to simplified) + */ +std::tuple +ramer_douglas_peucker( + Eigen::Ref P, + double threshold); diff --git a/tests/test_cotmatrix.py b/tests/test_cotmatrix.py new file mode 100644 index 0000000..a21ad85 --- /dev/null +++ b/tests/test_cotmatrix.py @@ -0,0 +1,23 @@ +import compas +from compas_libigl.cotmatrix import trimesh_cotmatrix, trimesh_cotmatrix_entries +from compas.datastructures import Mesh + + +def test_trimesh_cotmatrix(): + mesh = Mesh.from_off(compas.get("tubemesh.off")) + mesh.quads_to_triangles() + M = mesh.to_vertices_and_faces() + L = trimesh_cotmatrix(M) + # Cotmatrix is V x V sparse matrix + assert L.shape[0] == mesh.number_of_vertices() + assert L.shape[1] == mesh.number_of_vertices() + + +def test_trimesh_cotmatrix_entries(): + mesh = Mesh.from_off(compas.get("tubemesh.off")) + mesh.quads_to_triangles() + M = mesh.to_vertices_and_faces() + C = trimesh_cotmatrix_entries(M) + # Cotmatrix entries is F x 3 (one cotan per edge per face) + assert C.shape[0] == mesh.number_of_faces() + assert C.shape[1] == 3 diff --git a/tests/test_grad.py b/tests/test_grad.py new file mode 100644 index 0000000..ba1b234 --- /dev/null +++ b/tests/test_grad.py @@ -0,0 +1,14 @@ +import compas +from compas_libigl.grad import trimesh_grad +from compas.datastructures import Mesh + + +def test_trimesh_grad(): + mesh = Mesh.from_off(compas.get("tubemesh.off")) + mesh.quads_to_triangles() + M = mesh.to_vertices_and_faces() + G = trimesh_grad(M) + # Gradient operator is (3*F) x V sparse matrix + # When multiplied by V-length scalar field, produces 3F gradient components + assert G.shape[0] == 3 * mesh.number_of_faces() + assert G.shape[1] == mesh.number_of_vertices() diff --git a/tests/test_meshing.py b/tests/test_meshing.py index 696f7b2..1e7230b 100644 --- a/tests/test_meshing.py +++ b/tests/test_meshing.py @@ -1,5 +1,10 @@ import compas -from compas_libigl.meshing import trimesh_remesh_along_isoline, trimesh_remesh_along_isolines +from compas_libigl.meshing import ( + trimesh_remesh_along_isoline, + trimesh_remesh_along_isolines, + trimesh_cut_mesh, + trimesh_face_components, +) from compas.datastructures import Mesh @@ -29,3 +34,27 @@ def test_trimesh_remesh_along_isolines(): assert len(F2) > 0 assert len(S2) == len(V2) assert len(G2) == len(F2) + + +def test_trimesh_cut_mesh(): + mesh = Mesh.from_off(compas.get("tubemesh.off")) + mesh.quads_to_triangles() + M = mesh.to_vertices_and_faces() + V, F = M + # Create cut flags: all ones means cut all edges + cuts = [[1, 1, 1] for _ in F] + Vn, Fn = trimesh_cut_mesh(M, cuts) + # After cutting all edges, we should have more vertices + assert len(Vn) >= len(V) + assert len(Fn) == len(F) + + +def test_trimesh_face_components(): + mesh = Mesh.from_off(compas.get("tubemesh.off")) + mesh.quads_to_triangles() + M = mesh.to_vertices_and_faces() + C = trimesh_face_components(M) + # Component IDs, one per face + assert len(C) == mesh.number_of_faces() + # Connected mesh should have single component (component 0) + assert max(C) == 0 diff --git a/tests/test_simplify.py b/tests/test_simplify.py new file mode 100644 index 0000000..81e616b --- /dev/null +++ b/tests/test_simplify.py @@ -0,0 +1,37 @@ +from compas_libigl.simplify import ramer_douglas_peucker + + +def test_ramer_douglas_peucker(): + # Simple polyline: a zigzag that can be simplified + points = [ + [0.0, 0.0, 0.0], + [1.0, 0.1, 0.0], # slightly off the line + [2.0, 0.0, 0.0], + [3.0, 0.05, 0.0], # slightly off the line + [4.0, 0.0, 0.0], + ] + # With high threshold, should simplify to just endpoints + S, J, Q = ramer_douglas_peucker(points, threshold=0.5) + assert len(S) <= len(points) + assert len(S) >= 2 # At least start and end points + # Q maps each original point to simplified curve - same length as input + assert len(Q) == len(points) + assert len(Q[0]) == 3 # Each Q entry is a 3D point + # With zero threshold, should keep all points + S2, J2, Q2 = ramer_douglas_peucker(points, threshold=0.0) + assert len(S2) == len(points) + + +def test_ramer_douglas_peucker_straight_line(): + # Perfectly straight line should simplify to 2 points + points = [ + [0.0, 0.0, 0.0], + [1.0, 0.0, 0.0], + [2.0, 0.0, 0.0], + [3.0, 0.0, 0.0], + [4.0, 0.0, 0.0], + ] + S, J, Q = ramer_douglas_peucker(points, threshold=0.01) + assert len(S) == 2 + assert S[0] == [0.0, 0.0, 0.0] + assert S[1] == [4.0, 0.0, 0.0]