Merge pull request #60 from compas-dev/copilot/add-additional-parameter-inputs

Add surface meshing parameters to Poisson surface reconstruction

Author: tomvanmele

CHANGELOG.md

@@ -19,6 +19,7 @@ and this project adheres to [Semantic Versioning](https://semver.org/spec/v2.0.0
 
 
 ### Changed
+* Added optional surface meshing parameters (`sm_angle`, `sm_radius`, `sm_distance`) to `compas_cgal.reconstruction.poisson_surface_reconstruction` for controlling mesh quality and density.
 
 ### Removed
 

_codeql_detected_source_root

@@ -0,0 +1 @@
+.

docs/examples/example_reconstruction_poisson_surface_reconstruction.py

@@ -23,7 +23,16 @@
 # Reconstruction
 # =============================================================================
 
-V, F = poisson_surface_reconstruction(points, normals)
+# Reconstruct surface with custom parameters
+# Using larger sm_radius and sm_distance values reduces mesh complexity
+# and helps filter out vertices that don't belong to the original point cloud
+V, F = poisson_surface_reconstruction(
+    points,
+    normals,
+    sm_angle=20.0,  # Surface meshing angle bound (degrees)
+    sm_radius=30.0,  # Surface meshing radius bound (factor of avg spacing)
+    sm_distance=0.375,  # Surface meshing distance bound (factor of avg spacing)
+)
 mesh = Mesh.from_vertices_and_faces(V, F)
 
 # ==============================================================================

src/compas_cgal/reconstruction.py

@@ -14,6 +14,9 @@
 def poisson_surface_reconstruction(
     points: Union[list[Point], FloatNx3],
     normals: Union[list[Vector], FloatNx3],
+    sm_angle: float = 20.0,
+    sm_radius: float = 30.0,
+    sm_distance: float = 0.375,
 ) -> Tuple[FloatNx3, IntNx3]:
     """Reconstruct a surface from a point cloud using the Poisson surface reconstruction algorithm.
 
@@ -23,6 +26,20 @@ def poisson_surface_reconstruction(
         The points of the point cloud.
     normals
         The normals of the point cloud.
+    sm_angle : float, optional
+        Surface meshing angle bound in degrees.
+        Controls the minimum angle of triangles in the output mesh.
+        Default is 20.0.
+    sm_radius : float, optional
+        Surface meshing radius bound as a factor of average spacing.
+        Controls the size of triangles relative to the point cloud density.
+        Larger values result in coarser meshes with fewer vertices.
+        Default is 30.0.
+    sm_distance : float, optional
+        Surface meshing approximation error bound as a factor of average spacing.
+        Controls how closely the mesh approximates the implicit surface.
+        Larger values result in coarser meshes with fewer vertices that may deviate more from the original point cloud.
+        Default is 0.375.
 
     Returns
     -------
@@ -46,6 +63,18 @@ def poisson_surface_reconstruction(
     2. Well-oriented normals
     3. Points distributed across a meaningful surface
 
+    The surface meshing parameters (sm_angle, sm_radius, sm_distance) control the quality and
+    density of the output mesh. Increasing sm_radius and sm_distance will typically result in
+    fewer mesh vertices, which can help filter out vertices that don't belong to the original
+    point cloud, but may also reduce detail.
+
+    Examples
+    --------
+    >>> points = [[0, 0, 0], [1, 0, 0], [0, 1, 0], [0, 0, 1]]
+    >>> normals = [[0, 0, 1], [0, 0, 1], [0, 0, 1], [0, 0, 1]]
+    >>> V, F = poisson_surface_reconstruction(points, normals)
+    >>> # Use larger sm_radius and sm_distance to reduce mesh complexity
+    >>> V, F = poisson_surface_reconstruction(points, normals, sm_radius=50.0, sm_distance=0.5)
     """
     # Convert input to numpy arrays with proper type and memory layout
     P = np.asarray(points, dtype=np.float64, order="C")
@@ -67,7 +96,7 @@ def poisson_surface_reconstruction(
         N = N / norms[:, np.newaxis]
 
     try:
-        return _reconstruction.poisson_surface_reconstruction(P, N)
+        return _reconstruction.poisson_surface_reconstruction(P, N, sm_angle, sm_radius, sm_distance)
     except RuntimeError as e:
         raise RuntimeError(f"Poisson surface reconstruction failed: {str(e)}")
 

src/reconstruction.cpp

@@ -7,7 +7,10 @@ typedef CGAL::Parallel_if_available_tag ConcurrencyTag;
 std::tuple<compas::RowMatrixXd, compas::RowMatrixXi>
 poisson_surface_reconstruction(
     Eigen::Ref<const compas::RowMatrixXd> points,
-    Eigen::Ref<const compas::RowMatrixXd> normals)
+    Eigen::Ref<const compas::RowMatrixXd> normals,
+    double sm_angle,
+    double sm_radius,
+    double sm_distance)
 {
     compas::Polyhedron mesh;
     std::vector<PointVectorPair> points_with_normals;
@@ -30,7 +33,10 @@ poisson_surface_reconstruction(
         CGAL::First_of_pair_property_map<PointVectorPair>(),
         CGAL::Second_of_pair_property_map<PointVectorPair>(),
         mesh,
-        average_spacing);
+        average_spacing,
+        sm_angle,
+        sm_radius,
+        sm_distance);
 
     return compas::polyhedron_to_vertices_and_faces(mesh);
 }
@@ -238,7 +244,10 @@ NB_MODULE(_reconstruction, m) {
         &poisson_surface_reconstruction,
         "Perform Poisson surface reconstruction on an oriented pointcloud with normals",
         "points"_a,
-        "normals"_a
+        "normals"_a,
+        "sm_angle"_a = 20.0,
+        "sm_radius"_a = 30.0,
+        "sm_distance"_a = 0.375
     );
 
     m.def(

src/reconstruction.h

@@ -19,14 +19,20 @@
  *
  * @param points Matrix of point positions as Nx3 matrix in row-major order (float64)
  * @param normals Matrix of point normals as Nx3 matrix in row-major order (float64)
+ * @param sm_angle Surface meshing angle bound in degrees (default: 20.0)
+ * @param sm_radius Surface meshing radius bound as a factor of average spacing (default: 30.0)
+ * @param sm_distance Surface meshing approximation error bound as a factor of average spacing (default: 0.375)
  * @return std::tuple<RowMatrixXd, RowMatrixXi> containing:
  *         - vertices as Rx3 matrix (float64)
  *         - faces as Sx3 matrix (int32)
  */
 std::tuple<compas::RowMatrixXd, compas::RowMatrixXi>
 poisson_surface_reconstruction(
     Eigen::Ref<const compas::RowMatrixXd> points,
-    Eigen::Ref<const compas::RowMatrixXd> normals);
+    Eigen::Ref<const compas::RowMatrixXd> normals,
+    double sm_angle = 20.0,
+    double sm_radius = 30.0,
+    double sm_distance = 0.375);
 
 /**
  * @brief Remove outliers from a pointcloud.

tests/test_reconstruction_poisson_surface_reconstruction.py

@@ -7,17 +7,35 @@
 from compas.geometry import Scale
 from compas_cgal.reconstruction import poisson_surface_reconstruction
 
+# Test data file path
+TEST_DATA_FILE = Path(__file__).parent.parent / "data" / "oni.xyz"
 
-def test_reconstruction_poisson_surface_reconstruction():
-    FILE = Path(__file__).parent.parent / "data" / "oni.xyz"
 
+def load_test_data(file_path):
+    """Load point cloud data from file.
+
+    Parameters
+    ----------
+    file_path : Path
+        Path to the xyz file containing points and normals.
+
+    Returns
+    -------
+    tuple
+        (points, normals) where each is a list of 3D coordinates.
+    """
     points = []
     normals = []
-    with open(FILE, "r") as f:
+    with open(file_path, "r") as f:
         for line in f:
             x, y, z, nx, ny, nz = line.strip().split()
             points.append([float(x), float(y), float(z)])
             normals.append([float(nx), float(ny), float(nz)])
+    return points, normals
+
+
+def test_reconstruction_poisson_surface_reconstruction():
+    points, normals = load_test_data(TEST_DATA_FILE)
 
     V, F = poisson_surface_reconstruction(points, normals)
     mesh = Mesh.from_vertices_and_faces(V, F)
@@ -33,3 +51,38 @@ def test_reconstruction_poisson_surface_reconstruction():
     assert mesh.is_manifold()
     assert mesh.number_of_vertices() > 0
     assert mesh.number_of_faces() > 0
+
+
+def test_reconstruction_poisson_surface_reconstruction_with_parameters():
+    """Test Poisson surface reconstruction with custom parameters to reduce mesh complexity."""
+    points, normals = load_test_data(TEST_DATA_FILE)
+
+    # Test with larger sm_radius and sm_distance to reduce mesh complexity
+    V1, F1 = poisson_surface_reconstruction(points, normals, sm_radius=50.0, sm_distance=0.5)
+    mesh1 = Mesh.from_vertices_and_faces(V1, F1)
+
+    # Test with default parameters
+    V2, F2 = poisson_surface_reconstruction(points, normals)
+    mesh2 = Mesh.from_vertices_and_faces(V2, F2)
+
+    # Mesh with larger parameters should have fewer or equal vertices
+    # (though this isn't guaranteed in all cases due to CGAL's internal logic)
+    assert mesh1.is_manifold()
+    assert mesh1.number_of_vertices() > 0
+    assert mesh1.number_of_faces() > 0
+    assert mesh2.is_manifold()
+    assert mesh2.number_of_vertices() > 0
+    assert mesh2.number_of_faces() > 0
+
+
+def test_reconstruction_poisson_surface_reconstruction_angle_parameter():
+    """Test Poisson surface reconstruction with custom angle parameter."""
+    points, normals = load_test_data(TEST_DATA_FILE)
+
+    # Test with custom angle parameter
+    V, F = poisson_surface_reconstruction(points, normals, sm_angle=25.0)
+    mesh = Mesh.from_vertices_and_faces(V, F)
+
+    assert mesh.is_manifold()
+    assert mesh.number_of_vertices() > 0
+    assert mesh.number_of_faces() > 0