Merge pull request #57 from jf---/jf/polylines-module

add polylines module: simplify + closest point queries

Author: petrasvestartas

CHANGELOG.md

@@ -10,6 +10,10 @@ and this project adheres to [Semantic Versioning](https://semver.org/spec/v2.0.0
 ### Added
 * Added  `compas_cgal.skeletonization.mesh_skeleton_with_mapping`.
 
+- Added `simplify_polyline` and `simplify_polylines` functions for polyline simplification using Douglas-Peucker algorithm
+- Added `closest_points_on_polyline` function for batch closest point queries on polylines
+
+
 ### Changed
 
 ### Removed

CMakeLists.txt

@@ -236,5 +236,5 @@ add_nanobind_module(_slicer src/slicer.cpp)
 add_nanobind_module(_straight_skeleton_2 src/straight_skeleton_2.cpp)
 add_nanobind_module(_triangulation src/triangulation.cpp)
 add_nanobind_module(_subdivision src/subdivision.cpp)
-
+add_nanobind_module(_polylines src/polylines.cpp)
 

docs/_images/example_simplify_polylines.png

@@ -0,0 +1,42 @@
+import numpy as np
+from compas.geometry import Polyline
+from compas_viewer import Viewer
+from compas_viewer.config import Config
+
+from compas_cgal.polylines import simplify_polyline
+
+
+# Create a complex polyline: spiral with noise
+n_points = 200
+t = np.linspace(0, 6 * np.pi, n_points)
+radius = 1 + t / 10
+noise = np.random.uniform(-0.05, 0.05, n_points)
+
+x = radius * np.cos(t) + noise
+y = radius * np.sin(t) + noise
+z = t / 5
+
+original_points = np.column_stack([x, y, z])
+
+# Simplify with different thresholds
+simplified_high = simplify_polyline(original_points, threshold=0.5)
+
+# Create polylines offset in X for visualization
+offset = 6.0
+original_polyline = Polyline(original_points.tolist())
+simplified_high_polyline = Polyline((simplified_high + [offset, 0, 0]).tolist())
+
+# ==============================================================================
+# Visualize
+# ==============================================================================
+
+config = Config()
+config.camera.target = [5, 0, 5]
+config.camera.position = [5, -20, 10]
+
+viewer = Viewer(config=config)
+
+viewer.scene.add(original_polyline, linewidth=2, show_points=False)
+viewer.scene.add(simplified_high_polyline, linewidth=3, show_points=False)
+
+viewer.show()

docs/examples/example_simplify_polylines.py

@@ -0,0 +1,17 @@
+Polyline Simplification
+=======================
+
+This example demonstrates how to simplify polylines using the Douglas-Peucker algorithm in COMPAS CGAL.
+
+Key Features:
+
+* Douglas-Peucker polyline simplification
+* XY-plane simplification with Z-coordinate preservation
+* Side-by-side visualization of original and simplified polylines
+
+.. figure:: /_images/example_simplify_polylines.png
+    :figclass: figure
+    :class: figure-img img-fluid
+
+.. literalinclude:: example_simplify_polylines.py
+   :language: python

docs/examples/example_simplify_polylines.rst

@@ -24,6 +24,7 @@
     "compas_cgal.triangulation",
     "compas_cgal.slicer",
     "compas_cgal.subdivision",
+    "compas_cgal.polylines",
 ]
 
 __all__ = ["HOME", "DATA", "DOCS", "TEMP"]

src/compas_cgal/__init__.py

@@ -0,0 +1,104 @@
+"""Polyline utilities using CGAL."""
+
+from typing import List
+from typing import Union
+
+import numpy as np
+from numpy.typing import NDArray
+
+from compas_cgal import _types_std  # noqa: F401  # Load vector type bindings
+from compas_cgal._polylines import closest_points_on_polyline as _closest_points
+from compas_cgal._polylines import simplify_polylines as _simplify
+
+PointsList = Union[List[List[float]], NDArray]
+
+
+__all__ = ["simplify_polylines", "simplify_polyline", "closest_points_on_polyline"]
+
+
+def simplify_polylines(polylines: List[PointsList], threshold: float) -> List[NDArray]:
+    """Simplify multiple polylines using Douglas-Peucker algorithm.
+
+    Simplification is performed in the XY plane only. For 3D polylines,
+    Z coordinates are preserved but not considered in distance calculations.
+
+    Parameters
+    ----------
+    polylines : list of array-like
+        List of polylines. Each polyline is a sequence of 2D or 3D points.
+    threshold : float
+        Distance threshold for simplification. Higher values remove more points.
+
+    Returns
+    -------
+    list of ndarray
+        Simplified polylines as numpy arrays.
+
+    Examples
+    --------
+    >>> polylines = [[[0, 0], [1, 0.01], [2, 0]], [[0, 0], [0, 1], [1, 1]]]
+    >>> simplified = simplify_polylines(polylines, threshold=0.1)
+    >>> len(simplified[0])  # First polyline simplified
+    2
+    >>> len(simplified[1])  # Second has corner, preserved
+    3
+
+    """
+    if threshold < 0:
+        raise ValueError("threshold must be non-negative")
+    arrays = [np.asarray(p, dtype=np.float64) for p in polylines]
+    return _simplify(arrays, threshold)
+
+
+def simplify_polyline(polyline: PointsList, threshold: float) -> NDArray:
+    """Simplify a single polyline using Douglas-Peucker algorithm.
+
+    Simplification is performed in the XY plane only. For 3D polylines,
+    Z coordinates are preserved but not considered in distance calculations.
+
+    Parameters
+    ----------
+    polyline : array-like
+        Sequence of 2D or 3D points.
+    threshold : float
+        Distance threshold for simplification.
+
+    Returns
+    -------
+    ndarray
+        Simplified polyline as numpy array.
+
+    """
+    result = simplify_polylines([polyline], threshold)
+    return result[0]
+
+
+def closest_points_on_polyline(query_points: PointsList, polyline: PointsList) -> NDArray:
+    """Find closest points on a polyline for a batch of query points.
+
+    Uses CGAL's AABB tree for efficient batch queries.
+
+    Parameters
+    ----------
+    query_points : array-like
+        Query points as Mx2 or Mx3 array.
+    polyline : array-like
+        Polyline as Nx2 or Nx3 array.
+
+    Returns
+    -------
+    ndarray
+        Closest points on the polyline (same shape as query_points).
+
+    Examples
+    --------
+    >>> polyline = [[0, 0], [10, 0]]
+    >>> queries = [[5, 5], [3, -2]]
+    >>> closest = closest_points_on_polyline(queries, polyline)
+    >>> closest[0]  # Closest to (5, 5) on horizontal line
+    array([5., 0.])
+
+    """
+    queries = np.asarray(query_points, dtype=np.float64)
+    poly = np.asarray(polyline, dtype=np.float64)
+    return _closest_points(queries, poly)