Add exact Euclidean distance transform

CHANGELOG.md

@@ -21,6 +21,11 @@ removed no sooner than the next major (see `docs/API_STABILITY.md`).
 
 ### Added
 
+- **Exact Euclidean distance transform**: `spatialrust-vision` now computes
+  foreground-to-nearest-background L2 distances in linear time, supports
+  anisotropic pixel spacing, exposes a NumPy binding, and includes native
+  Criterion plus OpenCV `DIST_MASK_PRECISE` comparison coverage.
+
 - **Vision 2 performance roadmap and documentation site**: Epics 112–120 now
   define cost attribution, reusable workspaces, safe CPU dispatch, resize/color,
   Gaussian/Sobel, morphology, Canny, explicit GPU-chain, and release-gate work.

README.md

@@ -160,6 +160,7 @@ ratio; these are machine-specific measurements, not universal guarantees.
 | Sobel X 3×3 | OpenCV 36.6× | OpenCV 35.6× | OpenCV 38.6× |
 | Morphology open 5×5 | OpenCV 909.0× | OpenCV 1,027.2× | OpenCV 1,082.8× |
 | Canny | OpenCV 13.0× | OpenCV 14.5× | OpenCV 15.3× |
+| Exact Euclidean distance transform | OpenCV 10.60× | OpenCV 12.63× | OpenCV 12.35× |
 
 The current CPU result is deliberately mixed: SpatialRust's fused typed CHW
 path wins, while OpenCV's tuned general-purpose image kernels lead the present
@@ -181,6 +182,7 @@ The same deterministic RGB inputs passed all VGA, 1080p, and 4K gates:
 | Sobel X 3×3 | Exact values (max error 0) |
 | Morphology open 5×5 | Exact pixels (max error 0) |
 | Canny | Precision, recall, F1, and IoU all 1.0 |
+| Exact Euclidean distance transform | Exact values on canonical profiles; separate irregular-mask max float error `9.54e-7` |
 
 The broader correctness harness also checks filters, analysis, keypoints,
 matching, and geometry with documented tolerances (exact pixels where we claim

bench/opencv_comparison/manifest.json

@@ -32,6 +32,7 @@
     { "id": "sobel_x", "domain": "imgproc", "modes": ["allocate"] },
     { "id": "canny", "domain": "imgproc", "modes": ["allocate"] },
     { "id": "morphology_open", "domain": "imgproc", "modes": ["allocate"] },
+    { "id": "distance_transform_edt", "domain": "imgproc", "modes": ["allocate"] },
     { "id": "orb", "domain": "feature2d", "modes": ["allocate"] },
     { "id": "stereo_bm", "domain": "calib3d", "modes": ["allocate"] },
     { "id": "pinhole_calibration", "domain": "calib3d", "modes": ["allocate"] },

bench/opencv_vision_comparison/performance.py

@@ -314,6 +314,31 @@ def main() -> None:
             min_sample_time_ms=MIN_SAMPLE_TIME_MS,
         )
 
+        distance_mask = np.where(gray_cv > 96, 255, 0).astype(np.uint8)
+        distance_cv = cv2.distanceTransform(
+            distance_mask, cv2.DIST_L2, cv2.DIST_MASK_PRECISE
+        )
+        distance_sr = sr.distance_transform_edt(distance_mask)
+        distance_accuracy = numerical_accuracy(
+            distance_cv, distance_sr, float(np.hypot(width, height))
+        )
+        accuracy[f"{profile}_distance_transform_edt"] = distance_accuracy
+        if distance_accuracy["max_absolute_error"] > 1e-5:
+            raise AssertionError(
+                f"{profile} distance-transform max error "
+                f"{distance_accuracy['max_absolute_error']} > 1e-5"
+            )
+        _, _, cv_distance, sr_distance = timed_pair(
+            lambda: cv2.distanceTransform(
+                distance_mask, cv2.DIST_L2, cv2.DIST_MASK_PRECISE
+            ),
+            lambda: sr.distance_transform_edt(distance_mask),
+            warmup=args.warmup,
+            repeats=repeats,
+            seed=113,
+            min_sample_time_ms=MIN_SAMPLE_TIME_MS,
+        )
+
         rows = (
             ("resize_bilinear", "opencv", "allocate", cv_resize_alloc),
             ("resize_bilinear", "spatialrust", "allocate", sr_resize_alloc),
@@ -334,6 +359,8 @@ def main() -> None:
             ("morphology_open", "spatialrust", "allocate", sr_morphology),
             ("canny", "opencv", "allocate", cv_canny),
             ("canny", "spatialrust", "allocate", sr_canny),
+            ("distance_transform_edt", "opencv", "allocate", cv_distance),
+            ("distance_transform_edt", "spatialrust", "allocate", sr_distance),
         )
         measurements.extend(
             measurement(workload, implementation, mode, width, height, timing)
@@ -352,6 +379,7 @@ def main() -> None:
             "sobel_x": speed_comparison(cv_sobel, sr_sobel),
             "morphology_open": speed_comparison(cv_morphology, sr_morphology),
             "canny": speed_comparison(cv_canny, sr_canny),
+            "distance_transform_edt": speed_comparison(cv_distance, sr_distance),
         }
 
     environment_receipt = environment(

bench/opencv_vision_comparison/run.py

@@ -396,6 +396,18 @@ def main() -> None:
     if areas != areas_cv:
         raise AssertionError(f"component areas mismatch: {areas} != {areas_cv}")
 
+    distance_mask = np.ones((67, 89), dtype=np.uint8)
+    distance_mask[::11, ::13] = 0
+    distance_mask[20:28, 31:40] = 0
+    distance_sr = sr.distance_transform_edt(distance_mask)
+    distance_cv = cv2.distanceTransform(
+        distance_mask, cv2.DIST_L2, cv2.DIST_MASK_PRECISE
+    )
+    distance_error = float(np.max(np.abs(distance_sr - distance_cv)))
+    results["distance_transform_edt_max_f32_error"] = distance_error
+    if distance_error > 1e-5:
+        raise AssertionError(f"exact distance-transform error {distance_error} > 1e-5")
+
     # Geometry: planar homography residual agreement (not scale-normalized identity).
     source = np.array(
         [[10.0, 12.0], [70.0, 14.0], [18.0, 55.0], [66.0, 60.0], [40.0, 34.0], [28.0, 22.0]],

crates/spatialrust-platform/src/stability.rs

@@ -121,6 +121,8 @@ impl StabilityRegistry {
             "spatialrust-vision::video",
             "spatialrust-vision::odometry",
             "spatialrust-vision::photography",
+            "spatialrust-vision::distance_transform_edt",
+            "spatialrust-vision::distance_transform_edt_with_spacing",
             "spatialrust-runtime::execution-graph",
             "spatialrust-vision::ai-adapters",
             "spatialrust-gpu::GpuImage",

crates/spatialrust-py/spatialrust.pyi

@@ -36,7 +36,8 @@ __all__: list[str] = [
     "histogram_image", "equalize_histogram_image", "clahe_image",
     "integral_image_u8", "canny_image", "resize_image", "letterbox_image",
     "normalize_image_chw", "rgb_to_gray_image", "rgb_to_hsv_image", "remap_image",
-    "nms", "soft_nms", "connected_components_image", "find_mask_contours",
+    "nms", "soft_nms", "connected_components_image", "distance_transform_edt",
+    "find_mask_contours",
     "encode_mask_rle", "decode_mask_rle", "point_map_to_point_cloud", "knn_graph",
     "radius_graph", "register_icp", "register_point_to_plane", "register_gicp",
     "register_ndt", "register_fpfh_ransac", "register_fpfh_keypoints",
@@ -284,6 +285,9 @@ def connected_components_image(
     mask: _U8Array,
     connectivity: int = ...,
 ) -> tuple[_U32Array, list[tuple[int, int, tuple[float, float, float, float]]]]: ...
+def distance_transform_edt(
+    mask: _U8Array, spacing: tuple[float, float] = ...
+) -> _F32Array: ...
 def find_mask_contours(
     mask: _U8Array, epsilon: float = ...
 ) -> list[list[tuple[int, int]]]: ...

crates/spatialrust-py/src/lib.rs

@@ -75,6 +75,7 @@ use spatialrust::vision::{
     connected_components as label_components, decode_rle as decode_mask_runs,
     detect_and_describe_orb as detect_and_describe_orb_op, detect_fast as detect_fast_op,
     detect_harris as detect_harris_op, detect_shi_tomasi as detect_shi_tomasi_op,
+    distance_transform_edt_with_spacing as distance_transform_edt_op,
     encode_rle as encode_mask_runs, equalize_histogram as equalize_histogram_op,
     estimate_homography_ransac as estimate_homography_ransac_op,
     estimate_rgbd_odometry as estimate_rgbd_odometry_op, filter2d as filter2d_op,
@@ -3121,6 +3122,25 @@ fn connected_components_image<'py>(
     Ok((labels.into_pyarray_bound(py), stats))
 }
 
+/// Computes the exact Euclidean distance to the nearest zero-valued mask pixel.
+#[pyfunction]
+#[pyo3(signature = (mask, spacing=(1.0, 1.0)))]
+fn distance_transform_edt<'py>(
+    py: Python<'py>,
+    mask: PyReadonlyArray2<'_, u8>,
+    spacing: (f32, f32),
+) -> PyResult<Bound<'py, PyArray2<f32>>> {
+    let image = gray_u8_image_from_numpy(mask)?;
+    let (width, height) = (image.width(), image.height());
+    let binary = image.into_vec().into_iter().map(|value| u8::from(value != 0)).collect();
+    let mask = BinaryMask::try_new(width, height, binary).map_err(to_py_err)?;
+    let distances = distance_transform_edt_op(&mask, spacing.0, spacing.1).map_err(to_py_err)?;
+    let array =
+        Array2::from_shape_vec((distances.height(), distances.width()), distances.into_vec())
+            .map_err(to_py_err)?;
+    Ok(array.into_pyarray_bound(py))
+}
+
 /// Extracts and optionally simplifies mask contours.
 #[pyfunction]
 #[pyo3(signature = (mask, epsilon=0.0))]
@@ -3308,6 +3328,7 @@ fn spatialrust_module(m: &Bound<'_, PyModule>) -> PyResult<()> {
     m.add_function(wrap_pyfunction!(nms, m)?)?;
     m.add_function(wrap_pyfunction!(soft_nms, m)?)?;
     m.add_function(wrap_pyfunction!(connected_components_image, m)?)?;
+    m.add_function(wrap_pyfunction!(distance_transform_edt, m)?)?;
     m.add_function(wrap_pyfunction!(find_mask_contours, m)?)?;
     m.add_function(wrap_pyfunction!(encode_mask_rle, m)?)?;
     m.add_function(wrap_pyfunction!(decode_mask_rle, m)?)?;

crates/spatialrust-py/tests/test_bindings.py

@@ -66,7 +66,7 @@ def test_exports_present():
         "rgbd_to_point_cloud", "depth_to_xyz",
         "resize_image", "letterbox_image", "normalize_image_chw",
         "rgb_to_gray_image", "rgb_to_hsv_image", "remap_image",
-        "nms", "soft_nms", "connected_components_image",
+        "nms", "soft_nms", "connected_components_image", "distance_transform_edt",
         "find_mask_contours", "encode_mask_rle", "decode_mask_rle",
         "point_map_to_point_cloud",
     ):
@@ -199,6 +199,26 @@ def test_mask_components_contours_and_rle():
         np.testing.assert_array_equal(decoded, mask)
 
 
+def test_exact_euclidean_distance_transform():
+    mask = np.full((3, 4), 255, dtype=np.uint8)
+    mask[0, 0] = 0
+    actual = sr.distance_transform_edt(mask)
+    expected = np.array(
+        [
+            [0.0, 1.0, 2.0, 3.0],
+            [1.0, np.sqrt(2.0), np.sqrt(5.0), np.sqrt(10.0)],
+            [2.0, np.sqrt(5.0), np.sqrt(8.0), np.sqrt(13.0)],
+        ],
+        dtype=np.float32,
+    )
+    assert actual.dtype == np.float32
+    np.testing.assert_allclose(actual, expected, atol=1e-6)
+
+    anisotropic = sr.distance_transform_edt(mask, spacing=(2.0, 3.0))
+    assert anisotropic[0, 1] == pytest.approx(2.0)
+    assert anisotropic[1, 0] == pytest.approx(3.0)
+
+
 def test_point_map_to_point_cloud_filters_invalid_and_low_confidence():
     points = np.array(
         [[[0, 0, 1], [1, 0, 1]], [[0, 1, np.nan], [1, 1, 2]]], dtype=np.float32

crates/spatialrust-vision/Cargo.toml

@@ -67,6 +67,11 @@ name = "analysis"
 harness = false
 required-features = ["imgproc-analysis"]
 
+[[bench]]
+name = "dense"
+harness = false
+required-features = ["dense"]
+
 [[bench]]
 name = "canny"
 harness = false

crates/spatialrust-vision/benches/dense.rs

@@ -0,0 +1,21 @@
+use criterion::{black_box, criterion_group, criterion_main, BenchmarkId, Criterion, Throughput};
+use spatialrust_vision::{distance_transform_edt, BinaryMask};
+
+fn benchmark_exact_distance_transform(c: &mut Criterion) {
+    let mut group = c.benchmark_group("distance_transform_edt");
+    group.sample_size(10);
+    for &(name, width, height) in &[("vga", 640, 480), ("1080p", 1920, 1080), ("4k", 3840, 2160)] {
+        let data = (0..width * height)
+            .map(|index| u8::from(index % 97 != 0 && index % (width * 11 + 1) != 0))
+            .collect();
+        let mask = BinaryMask::try_new(width, height, data).unwrap();
+        group.throughput(Throughput::Elements((width * height) as u64));
+        group.bench_function(BenchmarkId::from_parameter(name), |b| {
+            b.iter(|| black_box(distance_transform_edt(black_box(&mask)).unwrap()));
+        });
+    }
+    group.finish();
+}
+
+criterion_group!(benches, benchmark_exact_distance_transform);
+criterion_main!(benches);