diff --git a/CHANGELOG.md b/CHANGELOG.md index 14cd5cc..6b692da 100644 --- a/CHANGELOG.md +++ b/CHANGELOG.md @@ -21,6 +21,19 @@ removed no sooner than the next major (see `docs/API_STABILITY.md`). ### Added +- **AI-ready image and vision foundation (Epics 75–79)**: mutable ROI views, + planar/interleaved layouts and color metadata in `spatialrust-image`; new + feature-gated `spatialrust-vision` preprocessing, warp, detection, mask/RLE, + and dense spatial-map APIs; camera/point-cloud/pipeline bridges; Python/NumPy + bindings and stubs; property tests, OpenCV comparison, Criterion benchmark, + and an end-to-end vision-to-point-cloud demo. + +- **OpenCV-oriented image/camera foundation**: new `spatialrust-image` typed + packed buffers and strided zero-copy views; new `spatialrust-camera` pinhole + projection/unprojection, Brown–Conrady distortion, aligned depth/RGB-D to + XYZ/XYZRGB conversion, Python binding, MVP integration test, Criterion bench, + synthetic demo, and OpenCV `rgbd.depthTo3d` comparison harness. + - **GPU Euclidean clustering** (`segment-euclidean-gpu`): WGSL uniform-grid label propagation, `GpuEuclideanClusterExtractor`, `EuclideanClusterExtractor::extract_with_policy`, MVP `cluster_policy`, and CLI diff --git a/Cargo.toml b/Cargo.toml index 64482a0..f58f0bd 100644 --- a/Cargo.toml +++ b/Cargo.toml @@ -14,6 +14,9 @@ members = [ "crates/spatialrust-metrics", "crates/spatialrust-transform", "crates/spatialrust-voxelize", + "crates/spatialrust-image", + "crates/spatialrust-camera", + "crates/spatialrust-vision", ] # spatialrust-py is a PyO3 cdylib built with maturin; keep it out of the Rust # workspace so `cargo test --workspace` does not require a Python toolchain. @@ -43,6 +46,9 @@ spatialrust-pipeline = { path = "crates/spatialrust-pipeline", version = "1.0.0" spatialrust-metrics = { path = "crates/spatialrust-metrics", version = "1.0.0", default-features = false } spatialrust-transform = { path = "crates/spatialrust-transform", version = "1.0.0", default-features = false } spatialrust-voxelize = { path = "crates/spatialrust-voxelize", version = "1.0.0", default-features = false } +spatialrust-image = { path = "crates/spatialrust-image", version = "1.0.0" } +spatialrust-camera = { path = "crates/spatialrust-camera", version = "1.0.0" } +spatialrust-vision = { path = "crates/spatialrust-vision", version = "1.0.0", default-features = false } bytemuck = { version = "1", features = ["derive"] } las = "0.9" @@ -56,6 +62,7 @@ serde = { version = "1", features = ["derive"] } thiserror = "2" wgpu = "24" criterion = { version = "0.5", features = ["html_reports"] } +proptest = "1" [profile.release] lto = "thin" diff --git a/README.md b/README.md index f0eeac2..f7e9a02 100644 --- a/README.md +++ b/README.md @@ -145,7 +145,7 @@ MVP pipeline is implemented end-to-end: PCD/PLY/LAS/COPC IO, voxel downsampling ## Workspace crates -One dataflow, eleven crates — each pipeline stage maps to the crate that implements it, all sitting on a small math/core/search foundation: +One dataflow, focused crates — each pipeline stage maps to the crate that implements it, all sitting on a small math/core/search foundation:

SpatialRust architecture: Load → Voxel → Normals → Plane → Cluster → Register → Save dataflow with implementing crates, wgpu voxel acceleration, and the core/math/search foundation @@ -156,6 +156,9 @@ One dataflow, eleven crates — each pipeline stage maps to the crate that imple | `spatialrust` | Meta crate / stable re-exports | | `spatialrust-core` | Point schema, metadata, execution traits | | `spatialrust-math` | Vec/Mat/Pose math primitives | +| `spatialrust-image` | Typed image buffers and zero-copy strided views | +| `spatialrust-camera` | Pinhole/Brown–Conrady camera models and RGB-D conversion | +| `spatialrust-vision` | Resize/preprocess, warps, detection postprocess, masks, and dense spatial maps | | `spatialrust-io` | Point cloud readers/writers (PCD, PLY, LAS, COPC) | | `spatialrust-search` | KD-tree search, k-NN / radius graphs | | `spatialrust-filtering` | Voxel / FPS downsample, outlier removal, crop, MLS | @@ -184,6 +187,39 @@ labels = result.labels() # (N,) int32 cluster ids sr.write("labeled.las", result.output) # LAS/PCD/PLY/COPC by extension ``` +Aligned RGB-D images feed the same point-cloud pipeline without an OpenCV +runtime dependency: + +```python +depth = np.ones((480, 640), dtype=np.float32) +rgb = np.zeros((480, 640, 3), dtype=np.uint8) +cloud = sr.rgbd_to_point_cloud( + depth, rgb, fx=525.0, fy=525.0, cx=319.5, cy=239.5 +) +result = sr.run_pipeline(cloud, leaf_size=0.03) +``` + +Rust users enable `camera-rgbd`; projection/unprojection supports optional +Brown–Conrady radial and tangential distortion. The reproducible numerical and +timing comparison against OpenCV is under `bench/opencv_rgbd_comparison/`. + +The `vision-full` feature adds an AI-ready CPU image path with explicit data +ownership: nearest/bilinear/bicubic/area resize, letterbox and CHW normalization, +color conversion, remap/warps, IoU/NMS/Soft-NMS, connected components, contours, +RLE masks, and depth/confidence/flow/point maps. Dense maps bridge explicitly to +calibrated cameras and point clouds; no API performs a hidden device transfer. + +```python +model_image, transform = sr.letterbox_image(rgb, 640, 640) +chw = sr.normalize_image_chw(model_image) # float32 (3,H,W) +keep = sr.nms(boxes_xyxy, scores, iou_threshold=0.5) +cloud = sr.point_map_to_point_cloud(points, confidence, 0.5) +``` + +The reproducible algorithm comparison is in +`bench/opencv_vision_comparison/`; the complete synthetic demo is +`crates/spatialrust-py/examples/vision_ai_pipeline.py`. +

Top-down view of clusters segmented from the public PCL table_scene_lms400 point cloud via a single Python run_pipeline() call

diff --git a/bench/opencv_rgbd_comparison/README.md b/bench/opencv_rgbd_comparison/README.md new file mode 100644 index 0000000..e62ba38 --- /dev/null +++ b/bench/opencv_rgbd_comparison/README.md @@ -0,0 +1,17 @@ +# OpenCV RGB-D comparison + +This harness compares SpatialRust `rgbd_to_point_cloud` against OpenCV's +`cv.rgbd.depthTo3d` on the same deterministic depth image and intrinsics. + +The OpenCV wheel must include the contrib `rgbd` module: + +```bash +pip install maturin numpy opencv-contrib-python +cd crates/spatialrust-py +maturin develop --release +cd ../.. +python bench/opencv_rgbd_comparison/run.py +``` + +The command exits non-zero when valid-point masks differ or maximum XYZ error +exceeds `1e-5` meters. It also reports median runtime for both implementations. diff --git a/bench/opencv_rgbd_comparison/run.py b/bench/opencv_rgbd_comparison/run.py new file mode 100644 index 0000000..c7ca5e7 --- /dev/null +++ b/bench/opencv_rgbd_comparison/run.py @@ -0,0 +1,58 @@ +"""Numerical and timing comparison with OpenCV rgbd.depthTo3d.""" + +from __future__ import annotations + +import statistics +import time + +import cv2 +import numpy as np +import spatialrust as sr + + +def timed(call, repeats: int = 20): + values = [] + result = None + for _ in range(repeats): + start = time.perf_counter() + result = call() + values.append(time.perf_counter() - start) + return result, statistics.median(values) + + +def main() -> None: + if not hasattr(cv2, "rgbd"): + raise RuntimeError("OpenCV rgbd module missing; install opencv-contrib-python") + + height, width = 240, 320 + yy, xx = np.mgrid[:height, :width] + depth = (1.0 + xx * 0.001 + yy * 0.0005).astype(np.float32) + depth[::31, ::29] = np.nan + color = np.empty((height, width, 3), dtype=np.uint8) + color[..., 0] = xx % 256 + color[..., 1] = yy % 256 + color[..., 2] = 127 + fx, fy, cx, cy = 280.0, 282.0, 159.5, 119.5 + intrinsics = np.array([[fx, 0.0, cx], [0.0, fy, cy], [0.0, 0.0, 1.0]]) + + cv_points, cv_seconds = timed(lambda: cv2.rgbd.depthTo3d(depth, intrinsics)) + sr_cloud, sr_seconds = timed( + lambda: sr.rgbd_to_point_cloud(depth, color, fx, fy, cx, cy) + ) + mask = np.isfinite(depth) & (depth > 0) + expected = cv_points[mask].astype(np.float32) + actual = sr_cloud.xyz() + if expected.shape != actual.shape: + raise AssertionError(f"shape mismatch: OpenCV={expected.shape}, SpatialRust={actual.shape}") + max_error = float(np.max(np.abs(expected - actual))) + if max_error > 1e-5: + raise AssertionError(f"maximum XYZ error {max_error:.3e} exceeds 1e-5 m") + + print(f"points: {len(actual)}") + print(f"max XYZ error: {max_error:.3e} m") + print(f"SpatialRust median: {sr_seconds * 1e3:.3f} ms") + print(f"OpenCV median: {cv_seconds * 1e3:.3f} ms") + + +if __name__ == "__main__": + main() diff --git a/bench/opencv_vision_comparison/README.md b/bench/opencv_vision_comparison/README.md new file mode 100644 index 0000000..8682aa4 --- /dev/null +++ b/bench/opencv_vision_comparison/README.md @@ -0,0 +1,16 @@ +# OpenCV vision comparison + +This deterministic harness compares SpatialRust's Python-visible CPU vision +primitives with OpenCV: four resize filters, RGB-to-gray/HSV conversion, +bilinear remap, NMS, and connected-component areas. + +From the repository root, after installing the editable Python extension: + +```powershell +$env:PYTHONPATH=(Resolve-Path '.\.venv\Lib\site-packages').Path +python bench\opencv_vision_comparison\run.py +``` + +The command prints a JSON report and exits non-zero when a documented numerical +tolerance is exceeded. OpenCV is comparison/test tooling only; it is not a Rust +runtime dependency. diff --git a/bench/opencv_vision_comparison/run.py b/bench/opencv_vision_comparison/run.py new file mode 100644 index 0000000..9c588ae --- /dev/null +++ b/bench/opencv_vision_comparison/run.py @@ -0,0 +1,109 @@ +"""Numerically compare SpatialRust vision primitives with OpenCV. + +Run after `maturin develop` so the native `spatialrust` module is importable. +The script exits non-zero when a compatibility tolerance is exceeded. +""" + +from __future__ import annotations + +import json + +import cv2 +import numpy as np + +import spatialrust as sr + + +def max_abs(a: np.ndarray, b: np.ndarray) -> int: + return int(np.max(np.abs(a.astype(np.int16) - b.astype(np.int16)))) + + +def main() -> None: + rng = np.random.default_rng(75) + image = rng.integers(0, 256, size=(73, 97, 3), dtype=np.uint8) + size = (61, 43) + results: dict[str, object] = {"opencv_version": cv2.__version__} + + resize_cases = { + "nearest": getattr(cv2, "INTER_NEAREST_EXACT", cv2.INTER_NEAREST), + "bilinear": cv2.INTER_LINEAR, + "bicubic": cv2.INTER_CUBIC, + "area": cv2.INTER_AREA, + } + limits = {"nearest": 0, "bilinear": 1, "bicubic": 1, "area": 1} + for name, flag in resize_cases.items(): + actual = sr.resize_image(image, size[0], size[1], interpolation=name) + expected = cv2.resize(image, size, interpolation=flag) + error = max_abs(actual, expected) + results[f"resize_{name}_max_u8_error"] = error + if error > limits[name]: + raise AssertionError(f"{name} resize error {error} > {limits[name]}") + + gray = sr.rgb_to_gray_image(image) + gray_cv = cv2.cvtColor(image, cv2.COLOR_RGB2GRAY) + gray_error = max_abs(gray, gray_cv) + results["rgb_to_gray_max_u8_error"] = gray_error + if gray_error > 1: + raise AssertionError(f"gray error {gray_error} > 1") + + hsv = sr.rgb_to_hsv_image(image) + hsv_cv = cv2.cvtColor(image, cv2.COLOR_RGB2HSV) + hue_delta = np.abs(hsv[..., 0].astype(np.int16) - hsv_cv[..., 0].astype(np.int16)) + hue_error = int(np.max(np.minimum(hue_delta, 180 - hue_delta))) + sv_error = max_abs(hsv[..., 1:], hsv_cv[..., 1:]) + results["rgb_to_hsv_max_hue_error"] = hue_error + results["rgb_to_hsv_max_sv_error"] = sv_error + if hue_error > 1 or sv_error > 1: + raise AssertionError(f"HSV error exceeds tolerance: H={hue_error}, SV={sv_error}") + + grid_x, grid_y = np.meshgrid( + np.arange(image.shape[1], dtype=np.float32), + np.arange(image.shape[0], dtype=np.float32), + ) + map_x = grid_x + np.float32(0.3125) + map_y = grid_y - np.float32(0.1875) + remapped = sr.remap_image(image, map_x, map_y, interpolation="bilinear") + remapped_cv = cv2.remap( + image, + map_x, + map_y, + cv2.INTER_LINEAR, + borderMode=cv2.BORDER_CONSTANT, + borderValue=(0, 0, 0), + ) + remap_error = max_abs(remapped, remapped_cv) + results["remap_bilinear_max_u8_error"] = remap_error + if remap_error > 1: + raise AssertionError(f"remap error {remap_error} > 1") + + boxes_xyxy = np.array( + [[0, 0, 20, 20], [2, 2, 18, 18], [30, 30, 45, 45], [32, 32, 46, 46]], + dtype=np.float32, + ) + scores = np.array([0.95, 0.8, 0.9, 0.7], dtype=np.float32) + actual_indices = sr.nms(boxes_xyxy, scores, 0.1, 0.5).tolist() + boxes_xywh = boxes_xyxy.copy() + boxes_xywh[:, 2:] -= boxes_xywh[:, :2] + cv_indices = cv2.dnn.NMSBoxes(boxes_xywh.tolist(), scores.tolist(), 0.1, 0.5) + expected_indices = np.asarray(cv_indices).reshape(-1).astype(int).tolist() + results["nms_indices"] = actual_indices + if actual_indices != expected_indices: + raise AssertionError(f"NMS mismatch: {actual_indices} != {expected_indices}") + + mask = np.zeros((32, 40), dtype=np.uint8) + mask[2:10, 3:12] = 1 + mask[16:30, 20:37] = 1 + _, stats = sr.connected_components_image(mask, connectivity=8) + count_cv, _, stats_cv, _ = cv2.connectedComponentsWithStats(mask, connectivity=8) + areas = sorted(stat[1] for stat in stats) + areas_cv = sorted(int(value) for value in stats_cv[1:count_cv, cv2.CC_STAT_AREA]) + results["connected_component_areas"] = areas + if areas != areas_cv: + raise AssertionError(f"component areas mismatch: {areas} != {areas_cv}") + + results["status"] = "pass" + print(json.dumps(results, indent=2, sort_keys=True)) + + +if __name__ == "__main__": + main() diff --git a/crates/spatialrust-camera/Cargo.toml b/crates/spatialrust-camera/Cargo.toml new file mode 100644 index 0000000..f57735a --- /dev/null +++ b/crates/spatialrust-camera/Cargo.toml @@ -0,0 +1,25 @@ +[package] +name = "spatialrust-camera" +version.workspace = true +edition.workspace = true +license.workspace = true +authors.workspace = true +repository.workspace = true +rust-version.workspace = true +description = "Camera models and RGB-D point cloud conversion for SpatialRust" + +[features] +default = [] + +[dependencies] +spatialrust-core.workspace = true +spatialrust-image.workspace = true +spatialrust-math.workspace = true +thiserror.workspace = true + +[dev-dependencies] +criterion.workspace = true + +[[bench]] +name = "rgbd" +harness = false diff --git a/crates/spatialrust-camera/benches/rgbd.rs b/crates/spatialrust-camera/benches/rgbd.rs new file mode 100644 index 0000000..ac007fc --- /dev/null +++ b/crates/spatialrust-camera/benches/rgbd.rs @@ -0,0 +1,28 @@ +use criterion::{black_box, criterion_group, criterion_main, Criterion}; +use spatialrust_camera::{rgbd_to_point_cloud, CameraIntrinsics, PinholeCamera}; +use spatialrust_image::Image; + +fn benchmark_rgbd(c: &mut Criterion) { + let width = 640; + let height = 480; + let depth = Image::::try_new(width, height, vec![2.0; width * height]).unwrap(); + let color = Image::::try_new(width, height, vec![127; width * height * 3]).unwrap(); + let camera = PinholeCamera::new( + CameraIntrinsics::try_new(525.0, 525.0, 319.5, 239.5, width, height).unwrap(), + ); + + c.bench_function("rgbd_to_point_cloud_640x480", |b| { + b.iter(|| { + rgbd_to_point_cloud( + black_box(depth.view()), + black_box(color.view()), + black_box(&camera), + Default::default(), + ) + .unwrap() + }); + }); +} + +criterion_group!(benches, benchmark_rgbd); +criterion_main!(benches); diff --git a/crates/spatialrust-camera/src/distortion.rs b/crates/spatialrust-camera/src/distortion.rs new file mode 100644 index 0000000..7166b8f --- /dev/null +++ b/crates/spatialrust-camera/src/distortion.rs @@ -0,0 +1,88 @@ +use spatialrust_math::Vec2; + +/// Brown–Conrady radial and tangential lens-distortion coefficients. +#[derive(Clone, Copy, Debug, Default, PartialEq)] +pub struct BrownConrady { + /// First radial coefficient. + pub k1: f64, + /// Second radial coefficient. + pub k2: f64, + /// First tangential coefficient. + pub p1: f64, + /// Second tangential coefficient. + pub p2: f64, + /// Third radial coefficient. + pub k3: f64, +} + +impl BrownConrady { + /// Returns whether all coefficients are zero. + #[must_use] + pub fn is_identity(self) -> bool { + self == Self::default() + } + + /// Distorts normalized pinhole coordinates. + #[must_use] + pub fn distort(self, point: Vec2) -> Vec2 { + let x = point.x; + let y = point.y; + let r2 = x.mul_add(x, y * y); + let radial = 1.0 + r2 * (self.k1 + r2 * (self.k2 + r2 * self.k3)); + Vec2 { + x: x * radial + 2.0 * self.p1 * x * y + self.p2 * (r2 + 2.0 * x * x), + y: y * radial + self.p1 * (r2 + 2.0 * y * y) + 2.0 * self.p2 * x * y, + } + } + + /// Iteratively removes distortion from normalized coordinates. + /// + /// Newton iterations use a numerical 2x2 Jacobian and stop once normalized + /// reprojection error reaches machine precision. + #[must_use] + pub fn undistort(self, distorted: Vec2) -> Vec2 { + if self.is_identity() { + return distorted; + } + let mut estimate = distorted; + const STEP: f64 = 1e-7; + for _ in 0..12 { + let observed = self.distort(estimate); + let error = Vec2 { x: observed.x - distorted.x, y: observed.y - distorted.y }; + if error.x.abs().max(error.y.abs()) < 1e-14 { + break; + } + let dx = self.distort(Vec2 { x: estimate.x + STEP, y: estimate.y }); + let dy = self.distort(Vec2 { x: estimate.x, y: estimate.y + STEP }); + let j00 = (dx.x - observed.x) / STEP; + let j10 = (dx.y - observed.y) / STEP; + let j01 = (dy.x - observed.x) / STEP; + let j11 = (dy.y - observed.y) / STEP; + let determinant = j00.mul_add(j11, -(j01 * j10)); + if determinant.abs() < f64::EPSILON { + break; + } + let delta_x = (j11 * error.x - j01 * error.y) / determinant; + let delta_y = (-j10 * error.x + j00 * error.y) / determinant; + estimate.x -= delta_x; + estimate.y -= delta_y; + } + estimate + } +} + +#[cfg(test)] +mod tests { + use super::BrownConrady; + use spatialrust_math::Vec2; + + #[test] + fn distortion_roundtrip() { + let model = BrownConrady { k1: -0.2, k2: 0.03, p1: 0.001, p2: -0.002, k3: 0.0 }; + for point in [Vec2 { x: -0.4, y: 0.3 }, Vec2 { x: 0.2, y: -0.1 }] { + let recovered = model.undistort(model.distort(point)); + assert!((recovered.x - point.x).abs() < 1e-9); + assert!((recovered.y - point.y).abs() < 1e-9); + } + } +} diff --git a/crates/spatialrust-camera/src/lib.rs b/crates/spatialrust-camera/src/lib.rs new file mode 100644 index 0000000..871e3e5 --- /dev/null +++ b/crates/spatialrust-camera/src/lib.rs @@ -0,0 +1,12 @@ +//! Camera models, projection geometry, and RGB-D conversion. + +#![deny(unsafe_code)] +#![warn(missing_docs)] + +mod distortion; +mod model; +mod rgbd; + +pub use distortion::BrownConrady; +pub use model::{CameraError, CameraIntrinsics, PinholeCamera}; +pub use rgbd::{depth_to_point_cloud, rgbd_to_point_cloud, DepthConversionOptions, RgbdError}; diff --git a/crates/spatialrust-camera/src/model.rs b/crates/spatialrust-camera/src/model.rs new file mode 100644 index 0000000..72a4464 --- /dev/null +++ b/crates/spatialrust-camera/src/model.rs @@ -0,0 +1,137 @@ +use crate::BrownConrady; +use spatialrust_math::{Vec2, Vec3}; + +/// Camera model errors. +#[derive(Clone, Debug, PartialEq, thiserror::Error)] +pub enum CameraError { + /// A focal length was zero, negative, or non-finite. + #[error("camera focal lengths must be finite and positive")] + InvalidFocalLength, + /// The principal point was non-finite. + #[error("camera principal point must be finite")] + InvalidPrincipalPoint, + /// Projection was requested for a point outside the positive camera half-space. + #[error("point depth must be finite and positive, found {0}")] + InvalidDepth(f64), + /// Pixel coordinates were non-finite. + #[error("pixel coordinates must be finite")] + InvalidPixel, +} + +/// Pinhole camera intrinsic parameters and image dimensions. +#[derive(Clone, Copy, Debug, PartialEq)] +pub struct CameraIntrinsics { + /// Horizontal focal length in pixels. + pub fx: f64, + /// Vertical focal length in pixels. + pub fy: f64, + /// Principal point x coordinate in pixels. + pub cx: f64, + /// Principal point y coordinate in pixels. + pub cy: f64, + /// Calibrated image width. + pub width: usize, + /// Calibrated image height. + pub height: usize, +} + +impl CameraIntrinsics { + /// Creates and validates camera intrinsics. + pub fn try_new( + fx: f64, + fy: f64, + cx: f64, + cy: f64, + width: usize, + height: usize, + ) -> Result { + if !fx.is_finite() || !fy.is_finite() || fx <= 0.0 || fy <= 0.0 { + return Err(CameraError::InvalidFocalLength); + } + if !cx.is_finite() || !cy.is_finite() { + return Err(CameraError::InvalidPrincipalPoint); + } + Ok(Self { fx, fy, cx, cy, width, height }) + } +} + +/// A pinhole camera with optional Brown–Conrady lens distortion. +#[derive(Clone, Copy, Debug, PartialEq)] +pub struct PinholeCamera { + /// Intrinsic calibration. + pub intrinsics: CameraIntrinsics, + /// Lens distortion coefficients. + pub distortion: BrownConrady, +} + +impl PinholeCamera { + /// Creates a camera with no lens distortion. + #[must_use] + pub const fn new(intrinsics: CameraIntrinsics) -> Self { + Self { + intrinsics, + distortion: BrownConrady { k1: 0.0, k2: 0.0, p1: 0.0, p2: 0.0, k3: 0.0 }, + } + } + + /// Attaches a Brown–Conrady distortion model. + #[must_use] + pub const fn with_distortion(mut self, distortion: BrownConrady) -> Self { + self.distortion = distortion; + self + } + + /// Projects a camera-space point into distorted pixel coordinates. + pub fn project(&self, point: Vec3) -> Result, CameraError> { + if !point.z.is_finite() || point.z <= 0.0 { + return Err(CameraError::InvalidDepth(point.z)); + } + let normalized = Vec2 { x: point.x / point.z, y: point.y / point.z }; + let distorted = self.distortion.distort(normalized); + Ok(Vec2 { + x: self.intrinsics.fx.mul_add(distorted.x, self.intrinsics.cx), + y: self.intrinsics.fy.mul_add(distorted.y, self.intrinsics.cy), + }) + } + + /// Unprojects a distorted pixel and metric depth into camera space. + pub fn unproject(&self, pixel: Vec2, depth: f64) -> Result, CameraError> { + if !depth.is_finite() || depth <= 0.0 { + return Err(CameraError::InvalidDepth(depth)); + } + if !pixel.x.is_finite() || !pixel.y.is_finite() { + return Err(CameraError::InvalidPixel); + } + let distorted = Vec2 { + x: (pixel.x - self.intrinsics.cx) / self.intrinsics.fx, + y: (pixel.y - self.intrinsics.cy) / self.intrinsics.fy, + }; + let normalized = self.distortion.undistort(distorted); + Ok(Vec3::new(normalized.x * depth, normalized.y * depth, depth)) + } +} + +#[cfg(test)] +mod tests { + use super::{CameraIntrinsics, PinholeCamera}; + use crate::BrownConrady; + use spatialrust_math::Vec3; + + #[test] + fn project_unproject_roundtrip_with_distortion() { + let intrinsics = CameraIntrinsics::try_new(525.0, 520.0, 319.5, 239.5, 640, 480).unwrap(); + let camera = PinholeCamera::new(intrinsics).with_distortion(BrownConrady { + k1: -0.15, + k2: 0.02, + p1: 0.001, + p2: -0.001, + k3: 0.0, + }); + let point = Vec3::new(0.4, -0.2, 2.5); + let pixel = camera.project(point).unwrap(); + let recovered = camera.unproject(pixel, point.z).unwrap(); + assert!((recovered.x - point.x).abs() < 1e-8); + assert!((recovered.y - point.y).abs() < 1e-8); + assert_eq!(recovered.z, point.z); + } +} diff --git a/crates/spatialrust-camera/src/rgbd.rs b/crates/spatialrust-camera/src/rgbd.rs new file mode 100644 index 0000000..bd042b3 --- /dev/null +++ b/crates/spatialrust-camera/src/rgbd.rs @@ -0,0 +1,226 @@ +use crate::{CameraError, PinholeCamera}; +use spatialrust_core::{ + PointBuffer, PointBufferSet, PointCloud, SpatialError, SpatialMetadata, StandardSchemas, +}; +use spatialrust_image::ImageView; +use spatialrust_math::Vec2; + +/// Errors from RGB-D conversion. +#[derive(Debug, thiserror::Error)] +pub enum RgbdError { + /// Camera dimensions and depth image dimensions differ. + #[error("depth dimensions {image_width}x{image_height} do not match camera dimensions {camera_width}x{camera_height}")] + DepthDimensionMismatch { + /// Depth image width. + image_width: usize, + /// Depth image height. + image_height: usize, + /// Calibrated camera width. + camera_width: usize, + /// Calibrated camera height. + camera_height: usize, + }, + /// RGB and depth image dimensions differ. + #[error("color dimensions {color_width}x{color_height} do not match depth dimensions {depth_width}x{depth_height}")] + ColorDimensionMismatch { + /// Color image width. + color_width: usize, + /// Color image height. + color_height: usize, + /// Depth image width. + depth_width: usize, + /// Depth image height. + depth_height: usize, + }, + /// Conversion settings were invalid. + #[error("{0}")] + InvalidOptions(String), + /// Camera geometry failed. + #[error(transparent)] + Camera(#[from] CameraError), + /// Point cloud construction failed. + #[error(transparent)] + Spatial(#[from] SpatialError), +} + +/// Controls conversion from stored depth values to metric camera coordinates. +#[derive(Clone, Copy, Debug, PartialEq)] +pub struct DepthConversionOptions { + /// Multiplier converting each stored depth value to meters. + pub depth_scale: f32, + /// Inclusive minimum accepted metric depth. + pub min_depth: f32, + /// Inclusive maximum accepted metric depth. + pub max_depth: f32, +} + +impl Default for DepthConversionOptions { + fn default() -> Self { + Self { depth_scale: 1.0, min_depth: f32::EPSILON, max_depth: f32::INFINITY } + } +} + +impl DepthConversionOptions { + fn validate(self) -> Result<(), RgbdError> { + if !self.depth_scale.is_finite() || self.depth_scale <= 0.0 { + return Err(RgbdError::InvalidOptions( + "depth_scale must be finite and positive".to_owned(), + )); + } + if !self.min_depth.is_finite() + || self.min_depth < 0.0 + || self.max_depth.is_nan() + || self.max_depth < self.min_depth + { + return Err(RgbdError::InvalidOptions( + "depth range must be ordered, non-negative, and not NaN".to_owned(), + )); + } + Ok(()) + } +} + +fn validate_depth(depth: ImageView<'_, f32, 1>, camera: &PinholeCamera) -> Result<(), RgbdError> { + let intrinsics = camera.intrinsics; + if depth.width() != intrinsics.width || depth.height() != intrinsics.height { + return Err(RgbdError::DepthDimensionMismatch { + image_width: depth.width(), + image_height: depth.height(), + camera_width: intrinsics.width, + camera_height: intrinsics.height, + }); + } + Ok(()) +} + +/// Converts an aligned depth image into an XYZ point cloud. +/// +/// Zero, non-finite, and out-of-range depths are omitted. +pub fn depth_to_point_cloud( + depth: ImageView<'_, f32, 1>, + camera: &PinholeCamera, + options: DepthConversionOptions, +) -> Result { + validate_depth(depth, camera)?; + options.validate()?; + let capacity = depth.width().saturating_mul(depth.height()); + let mut xs = Vec::with_capacity(capacity); + let mut ys = Vec::with_capacity(capacity); + let mut zs = Vec::with_capacity(capacity); + for y in 0..depth.height() { + for x in 0..depth.width() { + let meters = + depth.get(x, y).expect("validated image coordinates")[0] * options.depth_scale; + if !meters.is_finite() || meters < options.min_depth || meters > options.max_depth { + continue; + } + let point = camera.unproject(Vec2 { x: x as f64, y: y as f64 }, meters as f64)?; + xs.push(point.x as f32); + ys.push(point.y as f32); + zs.push(point.z as f32); + } + } + let mut buffers = PointBufferSet::new(); + buffers.insert("x", PointBuffer::from_f32(xs)); + buffers.insert("y", PointBuffer::from_f32(ys)); + buffers.insert("z", PointBuffer::from_f32(zs)); + Ok(PointCloud::try_from_parts( + StandardSchemas::point_xyz(), + buffers, + SpatialMetadata::default(), + )?) +} + +/// Converts aligned RGB and depth images into an XYZRGB point cloud. +/// +/// RGB channel order is preserved as `r`, `g`, and `b` `u8` fields. +pub fn rgbd_to_point_cloud( + depth: ImageView<'_, f32, 1>, + color: ImageView<'_, u8, 3>, + camera: &PinholeCamera, + options: DepthConversionOptions, +) -> Result { + validate_depth(depth, camera)?; + options.validate()?; + if color.width() != depth.width() || color.height() != depth.height() { + return Err(RgbdError::ColorDimensionMismatch { + color_width: color.width(), + color_height: color.height(), + depth_width: depth.width(), + depth_height: depth.height(), + }); + } + + let capacity = depth.width().saturating_mul(depth.height()); + let mut xs = Vec::with_capacity(capacity); + let mut ys = Vec::with_capacity(capacity); + let mut zs = Vec::with_capacity(capacity); + let mut rs = Vec::with_capacity(capacity); + let mut gs = Vec::with_capacity(capacity); + let mut bs = Vec::with_capacity(capacity); + for y in 0..depth.height() { + for x in 0..depth.width() { + let meters = + depth.get(x, y).expect("validated image coordinates")[0] * options.depth_scale; + if !meters.is_finite() || meters < options.min_depth || meters > options.max_depth { + continue; + } + let point = camera.unproject(Vec2 { x: x as f64, y: y as f64 }, meters as f64)?; + let rgb = color.get(x, y).expect("validated image coordinates"); + xs.push(point.x as f32); + ys.push(point.y as f32); + zs.push(point.z as f32); + rs.push(rgb[0]); + gs.push(rgb[1]); + bs.push(rgb[2]); + } + } + let mut buffers = PointBufferSet::new(); + buffers.insert("x", PointBuffer::from_f32(xs)); + buffers.insert("y", PointBuffer::from_f32(ys)); + buffers.insert("z", PointBuffer::from_f32(zs)); + buffers.insert("r", PointBuffer::U8(rs)); + buffers.insert("g", PointBuffer::U8(gs)); + buffers.insert("b", PointBuffer::U8(bs)); + Ok(PointCloud::try_from_parts( + StandardSchemas::point_xyzrgb(), + buffers, + SpatialMetadata::default(), + )?) +} + +#[cfg(test)] +mod tests { + use super::{depth_to_point_cloud, rgbd_to_point_cloud, DepthConversionOptions}; + use crate::{CameraIntrinsics, PinholeCamera}; + use spatialrust_core::PointBuffer; + use spatialrust_image::Image; + + fn camera() -> PinholeCamera { + PinholeCamera::new(CameraIntrinsics::try_new(2.0, 2.0, 0.0, 0.0, 2, 2).unwrap()) + } + + #[test] + fn skips_invalid_depth_and_unprojects() { + let depth = Image::::try_new(2, 2, vec![1.0, 0.0, f32::NAN, 2.0]).unwrap(); + let cloud = depth_to_point_cloud(depth.view(), &camera(), Default::default()).unwrap(); + assert_eq!(cloud.len(), 2); + assert_eq!(cloud.field("x").unwrap().as_f32().unwrap(), &[0.0, 1.0]); + assert_eq!(cloud.field("y").unwrap().as_f32().unwrap(), &[0.0, 1.0]); + assert_eq!(cloud.field("z").unwrap().as_f32().unwrap(), &[1.0, 2.0]); + } + + #[test] + fn rgb_fields_follow_valid_depths() { + let depth = Image::::try_new(2, 2, vec![1.0, 0.0, 3.0, 2.0]).unwrap(); + let color = + Image::::try_new(2, 2, vec![10, 11, 12, 20, 21, 22, 30, 31, 32, 40, 41, 42]) + .unwrap(); + let options = DepthConversionOptions { max_depth: 2.0, ..Default::default() }; + let cloud = rgbd_to_point_cloud(depth.view(), color.view(), &camera(), options).unwrap(); + assert_eq!(cloud.len(), 2); + assert_eq!(cloud.field("r").unwrap(), &PointBuffer::U8(vec![10, 40])); + assert_eq!(cloud.field("g").unwrap(), &PointBuffer::U8(vec![11, 41])); + assert_eq!(cloud.field("b").unwrap(), &PointBuffer::U8(vec![12, 42])); + } +} diff --git a/crates/spatialrust-image/Cargo.toml b/crates/spatialrust-image/Cargo.toml new file mode 100644 index 0000000..3851d30 --- /dev/null +++ b/crates/spatialrust-image/Cargo.toml @@ -0,0 +1,15 @@ +[package] +name = "spatialrust-image" +version.workspace = true +edition.workspace = true +license.workspace = true +authors.workspace = true +repository.workspace = true +rust-version.workspace = true +description = "Typed image buffers and zero-copy image views for SpatialRust" + +[features] +default = [] + +[dependencies] +thiserror.workspace = true diff --git a/crates/spatialrust-image/src/lib.rs b/crates/spatialrust-image/src/lib.rs new file mode 100644 index 0000000..c39d07a --- /dev/null +++ b/crates/spatialrust-image/src/lib.rs @@ -0,0 +1,1031 @@ +//! Typed, CPU-resident image buffers and zero-copy strided views. +//! +//! Channel count is part of the type. Packed interleaved ownership is the +//! default; planar ownership and explicitly-strided views are also available. +//! Device-backed images belong in dedicated GPU crates so transfers remain +//! explicit. + +#![deny(unsafe_code)] +#![warn(missing_docs)] + +use std::ops::{Index, IndexMut}; + +/// Errors raised while constructing or indexing images. +#[derive(Clone, Debug, PartialEq, Eq, thiserror::Error)] +pub enum ImageError { + /// A zero channel image type was requested. + #[error("image channel count must be greater than zero")] + ZeroChannels, + /// Image dimensions overflowed `usize` arithmetic. + #[error("image dimensions overflow addressable memory")] + DimensionOverflow, + /// The provided storage does not match the image layout. + #[error("image storage is too short: need at least {required} elements, found {found}")] + StorageTooShort { + /// Minimum required element count. + required: usize, + /// Provided element count. + found: usize, + }, + /// A row stride cannot hold one row of pixels. + #[error("row stride {stride} is smaller than packed row width {minimum}")] + InvalidStride { + /// Provided stride in scalar elements. + stride: usize, + /// Packed row size in scalar elements. + minimum: usize, + }, + /// A planar channel stride overlaps the preceding channel plane. + #[error("plane stride {stride} is smaller than one plane span {minimum}")] + InvalidPlaneStride { + /// Provided plane stride in scalar elements. + stride: usize, + /// Minimum non-overlapping plane span. + minimum: usize, + }, + /// Color metadata is incompatible with the compile-time channel count. + #[error("color space {color_space:?} requires {expected} channels, image has {found}")] + MetadataChannelMismatch { + /// Declared color space. + color_space: ColorSpace, + /// Required channel count. + expected: usize, + /// Image channel count. + found: usize, + }, + /// A requested image region lies outside its parent image. + #[error("region ({x}, {y}, {region_width}, {region_height}) exceeds image bounds {image_width}x{image_height}")] + InvalidRegion { + /// Region x origin. + x: usize, + /// Region y origin. + y: usize, + /// Region width. + region_width: usize, + /// Region height. + region_height: usize, + /// Parent image width. + image_width: usize, + /// Parent image height. + image_height: usize, + }, + /// Pixel coordinates were outside the image. + #[error("pixel ({x}, {y}) is outside image bounds {width}x{height}")] + OutOfBounds { + /// Pixel x coordinate. + x: usize, + /// Pixel y coordinate. + y: usize, + /// Image width. + width: usize, + /// Image height. + height: usize, + }, +} + +/// Physical channel arrangement in CPU storage. +#[derive(Clone, Copy, Debug, Default, PartialEq, Eq, Hash)] +pub enum ImageLayout { + /// Pixel channels are adjacent (`RGBRGB...`). + #[default] + Interleaved, + /// Each channel occupies a separate plane (`RR...GG...BB...`). + Planar, +} + +/// Semantic interpretation of image channels. +#[derive(Clone, Copy, Debug, Default, PartialEq, Eq, Hash)] +pub enum ColorSpace { + /// Channel semantics are not specified. + #[default] + Unknown, + /// One-channel luminance. + Gray, + /// Nonlinear red, green, blue. + Rgb, + /// Nonlinear blue, green, red. + Bgr, + /// Nonlinear red, green, blue, alpha. + Rgba, + /// Nonlinear blue, green, red, alpha. + Bgra, + /// Linear-light red, green, blue. + LinearRgb, + /// Hue, saturation, value. + Hsv, + /// Metric or sensor depth values. + Depth, + /// Integer semantic or instance labels. + Label, +} + +impl ColorSpace { + /// Returns the required channel count when the color space fixes one. + #[must_use] + pub const fn required_channels(self) -> Option { + match self { + Self::Unknown => None, + Self::Gray | Self::Depth | Self::Label => Some(1), + Self::Rgb | Self::Bgr | Self::LinearRgb | Self::Hsv => Some(3), + Self::Rgba | Self::Bgra => Some(4), + } + } +} + +/// Numeric range convention associated with image channels. +#[derive(Clone, Copy, Debug, Default, PartialEq, Eq, Hash)] +pub enum ColorRange { + /// Range is not specified or is naturally unbounded (for example depth). + #[default] + Unspecified, + /// Full dtype or normalized range. + Full, + /// Video-range encoding such as limited-range YUV. + Limited, +} + +/// Alpha-channel interpretation. +#[derive(Clone, Copy, Debug, Default, PartialEq, Eq, Hash)] +pub enum AlphaMode { + /// No alpha channel is declared. + #[default] + None, + /// Straight (unassociated) alpha. + Straight, + /// RGB values are premultiplied by alpha. + Premultiplied, +} + +/// Lightweight semantic metadata carried by images and borrowed views. +#[derive(Clone, Copy, Debug, Default, PartialEq, Eq, Hash)] +pub struct ImageMetadata { + /// Channel color/depth/label interpretation. + pub color_space: ColorSpace, + /// Numeric range convention. + pub color_range: ColorRange, + /// Alpha convention. + pub alpha_mode: AlphaMode, +} + +impl ImageMetadata { + /// Validates metadata against a compile-time channel count. + pub fn validate(self) -> Result<(), ImageError> { + if CHANNELS == 0 { + return Err(ImageError::ZeroChannels); + } + if let Some(expected) = self.color_space.required_channels() { + if expected != CHANNELS { + return Err(ImageError::MetadataChannelMismatch { + color_space: self.color_space, + expected, + found: CHANNELS, + }); + } + } + Ok(()) + } +} + +/// A checked rectangular image region. +#[derive(Clone, Copy, Debug, Default, PartialEq, Eq, Hash)] +pub struct ImageRegion { + /// Horizontal origin. + pub x: usize, + /// Vertical origin. + pub y: usize, + /// Region width. + pub width: usize, + /// Region height. + pub height: usize, +} + +impl ImageRegion { + /// Creates a rectangular region. + #[must_use] + pub const fn new(x: usize, y: usize, width: usize, height: usize) -> Self { + Self { x, y, width, height } + } + + fn validate(self, image_width: usize, image_height: usize) -> Result<(), ImageError> { + let end_x = self.x.checked_add(self.width).ok_or(ImageError::DimensionOverflow)?; + let end_y = self.y.checked_add(self.height).ok_or(ImageError::DimensionOverflow)?; + if end_x > image_width || end_y > image_height { + return Err(ImageError::InvalidRegion { + x: self.x, + y: self.y, + region_width: self.width, + region_height: self.height, + image_width, + image_height, + }); + } + Ok(()) + } +} + +fn packed_len(width: usize, height: usize) -> Result { + if CHANNELS == 0 { + return Err(ImageError::ZeroChannels); + } + width + .checked_mul(height) + .and_then(|value| value.checked_mul(CHANNELS)) + .ok_or(ImageError::DimensionOverflow) +} + +fn strided_span(height: usize, row_stride: usize, packed_row: usize) -> Result { + if height == 0 { + return Ok(0); + } + row_stride + .checked_mul(height - 1) + .and_then(|offset| offset.checked_add(packed_row)) + .ok_or(ImageError::DimensionOverflow) +} + +/// An owning, densely packed, interleaved image. +#[derive(Clone, Debug, PartialEq, Eq)] +pub struct Image { + width: usize, + height: usize, + data: Vec, + metadata: ImageMetadata, +} + +impl Image { + /// Creates an image from densely packed, interleaved scalar elements. + pub fn try_new(width: usize, height: usize, data: Vec) -> Result { + let required = packed_len::(width, height)?; + if data.len() != required { + return Err(ImageError::StorageTooShort { required, found: data.len() }); + } + Ok(Self { width, height, data, metadata: ImageMetadata::default() }) + } + + /// Creates a packed image and validates its semantic metadata. + pub fn try_new_with_metadata( + width: usize, + height: usize, + data: Vec, + metadata: ImageMetadata, + ) -> Result { + metadata.validate::()?; + let mut image = Self::try_new(width, height, data)?; + image.metadata = metadata; + Ok(image) + } + + /// Returns image width in pixels. + #[must_use] + pub const fn width(&self) -> usize { + self.width + } + + /// Returns image height in pixels. + #[must_use] + pub const fn height(&self) -> usize { + self.height + } + + /// Returns the packed row stride in scalar elements. + #[must_use] + pub const fn row_stride(&self) -> usize { + self.width * CHANNELS + } + + /// Returns physical channel layout. + #[must_use] + pub const fn layout(&self) -> ImageLayout { + ImageLayout::Interleaved + } + + /// Returns semantic image metadata. + #[must_use] + pub const fn metadata(&self) -> ImageMetadata { + self.metadata + } + + /// Replaces semantic metadata after validating the channel count. + pub fn set_metadata(&mut self, metadata: ImageMetadata) -> Result<(), ImageError> { + metadata.validate::()?; + self.metadata = metadata; + Ok(()) + } + + /// Returns the packed scalar storage. + #[must_use] + pub fn as_slice(&self) -> &[T] { + &self.data + } + + /// Returns mutable packed scalar storage. + #[must_use] + pub fn as_mut_slice(&mut self) -> &mut [T] { + &mut self.data + } + + /// Borrows this image as a zero-copy view. + #[must_use] + pub fn view(&self) -> ImageView<'_, T, CHANNELS> { + ImageView { + width: self.width, + height: self.height, + row_stride: self.row_stride(), + data: &self.data, + metadata: self.metadata, + } + } + + /// Borrows this image as a mutable zero-copy view. + #[must_use] + pub fn view_mut(&mut self) -> ImageViewMut<'_, T, CHANNELS> { + ImageViewMut { + width: self.width, + height: self.height, + row_stride: self.width * CHANNELS, + data: &mut self.data, + metadata: self.metadata, + } + } + + /// Returns one pixel, or `None` outside image bounds. + #[must_use] + pub fn get(&self, x: usize, y: usize) -> Option<&[T; CHANNELS]> { + self.view().get(x, y) + } + + /// Returns one mutable pixel, or `None` outside image bounds. + #[must_use] + pub fn get_mut(&mut self, x: usize, y: usize) -> Option<&mut [T; CHANNELS]> { + if x >= self.width || y >= self.height { + return None; + } + let offset = (y * self.width + x) * CHANNELS; + self.data.get_mut(offset..offset + CHANNELS)?.try_into().ok() + } + + /// Consumes the image and returns its packed scalar storage. + #[must_use] + pub fn into_vec(self) -> Vec { + self.data + } +} + +impl Image { + /// Creates a densely packed image filled with one pixel value. + pub fn from_pixel( + width: usize, + height: usize, + pixel: [T; CHANNELS], + ) -> Result { + let required = packed_len::(width, height)?; + let mut data = Vec::with_capacity(required); + for _ in 0..width.saturating_mul(height) { + data.extend_from_slice(&pixel); + } + Self::try_new(width, height, data) + } +} + +impl Index<(usize, usize)> for Image { + type Output = [T; CHANNELS]; + + fn index(&self, (x, y): (usize, usize)) -> &Self::Output { + self.get(x, y).expect("image index out of bounds") + } +} + +impl IndexMut<(usize, usize)> for Image { + fn index_mut(&mut self, (x, y): (usize, usize)) -> &mut Self::Output { + self.get_mut(x, y).expect("image index out of bounds") + } +} + +/// A read-only, zero-copy image view with an explicit row stride. +#[derive(Clone, Copy, Debug)] +pub struct ImageView<'a, T, const CHANNELS: usize> { + width: usize, + height: usize, + row_stride: usize, + data: &'a [T], + metadata: ImageMetadata, +} + +impl<'a, T, const CHANNELS: usize> ImageView<'a, T, CHANNELS> { + /// Creates a view over interleaved storage. + /// + /// `row_stride` is measured in scalar elements, not bytes. + pub fn new( + width: usize, + height: usize, + row_stride: usize, + data: &'a [T], + ) -> Result { + Self::new_with_metadata(width, height, row_stride, data, ImageMetadata::default()) + } + + /// Creates a view with validated semantic metadata. + pub fn new_with_metadata( + width: usize, + height: usize, + row_stride: usize, + data: &'a [T], + metadata: ImageMetadata, + ) -> Result { + metadata.validate::()?; + let packed_row = width.checked_mul(CHANNELS).ok_or(ImageError::DimensionOverflow)?; + if row_stride < packed_row { + return Err(ImageError::InvalidStride { stride: row_stride, minimum: packed_row }); + } + let required = if height == 0 { + 0 + } else { + row_stride + .checked_mul(height - 1) + .and_then(|offset| offset.checked_add(packed_row)) + .ok_or(ImageError::DimensionOverflow)? + }; + if data.len() < required { + return Err(ImageError::StorageTooShort { required, found: data.len() }); + } + Ok(Self { width, height, row_stride, data, metadata }) + } + + /// Returns image width in pixels. + #[must_use] + pub const fn width(self) -> usize { + self.width + } + + /// Returns image height in pixels. + #[must_use] + pub const fn height(self) -> usize { + self.height + } + + /// Returns row stride in scalar elements. + #[must_use] + pub const fn row_stride(self) -> usize { + self.row_stride + } + + /// Returns physical channel layout. + #[must_use] + pub const fn layout(self) -> ImageLayout { + ImageLayout::Interleaved + } + + /// Returns semantic image metadata. + #[must_use] + pub const fn metadata(self) -> ImageMetadata { + self.metadata + } + + /// Returns one pixel, or `None` outside image bounds. + #[must_use] + pub fn get(self, x: usize, y: usize) -> Option<&'a [T; CHANNELS]> { + if x >= self.width || y >= self.height { + return None; + } + let offset = y * self.row_stride + x * CHANNELS; + self.data.get(offset..offset + CHANNELS)?.try_into().ok() + } + + /// Returns a packed row without its trailing padding. + #[must_use] + pub fn row(self, y: usize) -> Option<&'a [T]> { + if y >= self.height { + return None; + } + let start = y * self.row_stride; + self.data.get(start..start + self.width * CHANNELS) + } + + /// Creates a checked zero-copy subview. + pub fn subview(self, region: ImageRegion) -> Result { + region.validate(self.width, self.height)?; + if region.width == 0 || region.height == 0 { + return Ok(Self { + width: region.width, + height: region.height, + row_stride: self.row_stride, + data: &self.data[..0], + metadata: self.metadata, + }); + } + let start = region.y * self.row_stride + region.x * CHANNELS; + let span = (region.height - 1) * self.row_stride + region.width * CHANNELS; + Ok(Self { + width: region.width, + height: region.height, + row_stride: self.row_stride, + data: &self.data[start..start + span], + metadata: self.metadata, + }) + } +} + +/// A mutable, zero-copy interleaved image view with an explicit row stride. +#[derive(Debug)] +pub struct ImageViewMut<'a, T, const CHANNELS: usize> { + width: usize, + height: usize, + row_stride: usize, + data: &'a mut [T], + metadata: ImageMetadata, +} + +impl<'a, T, const CHANNELS: usize> ImageViewMut<'a, T, CHANNELS> { + /// Creates a mutable view over interleaved storage. + pub fn new( + width: usize, + height: usize, + row_stride: usize, + data: &'a mut [T], + ) -> Result { + Self::new_with_metadata(width, height, row_stride, data, ImageMetadata::default()) + } + + /// Creates a mutable view with validated semantic metadata. + pub fn new_with_metadata( + width: usize, + height: usize, + row_stride: usize, + data: &'a mut [T], + metadata: ImageMetadata, + ) -> Result { + metadata.validate::()?; + let packed_row = width.checked_mul(CHANNELS).ok_or(ImageError::DimensionOverflow)?; + if row_stride < packed_row { + return Err(ImageError::InvalidStride { stride: row_stride, minimum: packed_row }); + } + let required = strided_span(height, row_stride, packed_row)?; + if data.len() < required { + return Err(ImageError::StorageTooShort { required, found: data.len() }); + } + Ok(Self { width, height, row_stride, data, metadata }) + } + + /// Returns image width in pixels. + #[must_use] + pub const fn width(&self) -> usize { + self.width + } + + /// Returns image height in pixels. + #[must_use] + pub const fn height(&self) -> usize { + self.height + } + + /// Returns row stride in scalar elements. + #[must_use] + pub const fn row_stride(&self) -> usize { + self.row_stride + } + + /// Returns semantic image metadata. + #[must_use] + pub const fn metadata(&self) -> ImageMetadata { + self.metadata + } + + /// Reborrows this mutable view as read-only. + #[must_use] + pub fn as_view(&self) -> ImageView<'_, T, CHANNELS> { + ImageView { + width: self.width, + height: self.height, + row_stride: self.row_stride, + data: self.data, + metadata: self.metadata, + } + } + + /// Returns one read-only pixel, or `None` outside image bounds. + #[must_use] + pub fn get(&self, x: usize, y: usize) -> Option<&[T; CHANNELS]> { + self.as_view().get(x, y) + } + + /// Returns one mutable pixel, or `None` outside image bounds. + #[must_use] + pub fn get_mut(&mut self, x: usize, y: usize) -> Option<&mut [T; CHANNELS]> { + if x >= self.width || y >= self.height { + return None; + } + let offset = y * self.row_stride + x * CHANNELS; + self.data.get_mut(offset..offset + CHANNELS)?.try_into().ok() + } + + /// Returns a mutable packed row without trailing padding. + #[must_use] + pub fn row_mut(&mut self, y: usize) -> Option<&mut [T]> { + if y >= self.height { + return None; + } + let start = y * self.row_stride; + self.data.get_mut(start..start + self.width * CHANNELS) + } + + /// Creates a checked mutable zero-copy subview. + pub fn subview( + &mut self, + region: ImageRegion, + ) -> Result, ImageError> { + region.validate(self.width, self.height)?; + if region.width == 0 || region.height == 0 { + return Ok(ImageViewMut { + width: region.width, + height: region.height, + row_stride: self.row_stride, + data: &mut self.data[..0], + metadata: self.metadata, + }); + } + let start = region.y * self.row_stride + region.x * CHANNELS; + let span = (region.height - 1) * self.row_stride + region.width * CHANNELS; + Ok(ImageViewMut { + width: region.width, + height: region.height, + row_stride: self.row_stride, + data: &mut self.data[start..start + span], + metadata: self.metadata, + }) + } +} + +/// An owning, densely packed planar image. +/// +/// Storage order is all values for channel 0, followed by channel 1, and so on. +#[derive(Clone, Debug, PartialEq, Eq)] +pub struct PlanarImage { + width: usize, + height: usize, + data: Vec, + metadata: ImageMetadata, +} + +impl PlanarImage { + /// Creates an image from densely packed planar scalar elements. + pub fn try_new(width: usize, height: usize, data: Vec) -> Result { + let required = packed_len::(width, height)?; + if data.len() != required { + return Err(ImageError::StorageTooShort { required, found: data.len() }); + } + Ok(Self { width, height, data, metadata: ImageMetadata::default() }) + } + + /// Creates a planar image with validated semantic metadata. + pub fn try_new_with_metadata( + width: usize, + height: usize, + data: Vec, + metadata: ImageMetadata, + ) -> Result { + metadata.validate::()?; + let mut image = Self::try_new(width, height, data)?; + image.metadata = metadata; + Ok(image) + } + + /// Returns image width in pixels. + #[must_use] + pub const fn width(&self) -> usize { + self.width + } + + /// Returns image height in pixels. + #[must_use] + pub const fn height(&self) -> usize { + self.height + } + + /// Returns physical channel layout. + #[must_use] + pub const fn layout(&self) -> ImageLayout { + ImageLayout::Planar + } + + /// Returns packed row stride within each plane. + #[must_use] + pub const fn row_stride(&self) -> usize { + self.width + } + + /// Returns the scalar distance between channel-plane origins. + #[must_use] + pub const fn plane_stride(&self) -> usize { + self.width * self.height + } + + /// Returns semantic image metadata. + #[must_use] + pub const fn metadata(&self) -> ImageMetadata { + self.metadata + } + + /// Replaces semantic metadata after validating the channel count. + pub fn set_metadata(&mut self, metadata: ImageMetadata) -> Result<(), ImageError> { + metadata.validate::()?; + self.metadata = metadata; + Ok(()) + } + + /// Returns planar scalar storage. + #[must_use] + pub fn as_slice(&self) -> &[T] { + &self.data + } + + /// Returns mutable planar scalar storage. + #[must_use] + pub fn as_mut_slice(&mut self) -> &mut [T] { + &mut self.data + } + + /// Borrows this image as a planar view. + #[must_use] + pub fn view(&self) -> PlanarImageView<'_, T, CHANNELS> { + PlanarImageView { + width: self.width, + height: self.height, + row_stride: self.width, + plane_stride: self.width * self.height, + data: &self.data, + metadata: self.metadata, + } + } + + /// Returns one channel value, or `None` outside image bounds. + #[must_use] + pub fn get(&self, channel: usize, x: usize, y: usize) -> Option<&T> { + self.view().get(channel, x, y) + } + + /// Returns one mutable channel value, or `None` outside image bounds. + #[must_use] + pub fn get_mut(&mut self, channel: usize, x: usize, y: usize) -> Option<&mut T> { + if channel >= CHANNELS || x >= self.width || y >= self.height { + return None; + } + let offset = channel * self.width * self.height + y * self.width + x; + self.data.get_mut(offset) + } + + /// Consumes the image and returns planar scalar storage. + #[must_use] + pub fn into_vec(self) -> Vec { + self.data + } +} + +/// A read-only, zero-copy planar image view with explicit strides. +#[derive(Clone, Copy, Debug)] +pub struct PlanarImageView<'a, T, const CHANNELS: usize> { + width: usize, + height: usize, + row_stride: usize, + plane_stride: usize, + data: &'a [T], + metadata: ImageMetadata, +} + +impl<'a, T, const CHANNELS: usize> PlanarImageView<'a, T, CHANNELS> { + /// Creates a planar view. Strides are measured in scalar elements. + pub fn new( + width: usize, + height: usize, + row_stride: usize, + plane_stride: usize, + data: &'a [T], + ) -> Result { + Self::new_with_metadata( + width, + height, + row_stride, + plane_stride, + data, + ImageMetadata::default(), + ) + } + + /// Creates a planar view with validated semantic metadata. + pub fn new_with_metadata( + width: usize, + height: usize, + row_stride: usize, + plane_stride: usize, + data: &'a [T], + metadata: ImageMetadata, + ) -> Result { + metadata.validate::()?; + if row_stride < width { + return Err(ImageError::InvalidStride { stride: row_stride, minimum: width }); + } + let plane_span = strided_span(height, row_stride, width)?; + if plane_stride < plane_span { + return Err(ImageError::InvalidPlaneStride { + stride: plane_stride, + minimum: plane_span, + }); + } + let required = if width == 0 || height == 0 { + 0 + } else { + plane_stride + .checked_mul(CHANNELS - 1) + .and_then(|offset| offset.checked_add(plane_span)) + .ok_or(ImageError::DimensionOverflow)? + }; + if data.len() < required { + return Err(ImageError::StorageTooShort { required, found: data.len() }); + } + Ok(Self { width, height, row_stride, plane_stride, data, metadata }) + } + + /// Returns image width in pixels. + #[must_use] + pub const fn width(self) -> usize { + self.width + } + + /// Returns image height in pixels. + #[must_use] + pub const fn height(self) -> usize { + self.height + } + + /// Returns row stride within a plane. + #[must_use] + pub const fn row_stride(self) -> usize { + self.row_stride + } + + /// Returns scalar distance between channel-plane origins. + #[must_use] + pub const fn plane_stride(self) -> usize { + self.plane_stride + } + + /// Returns physical channel layout. + #[must_use] + pub const fn layout(self) -> ImageLayout { + ImageLayout::Planar + } + + /// Returns semantic image metadata. + #[must_use] + pub const fn metadata(self) -> ImageMetadata { + self.metadata + } + + /// Returns one channel value, or `None` outside image bounds. + #[must_use] + pub fn get(self, channel: usize, x: usize, y: usize) -> Option<&'a T> { + if channel >= CHANNELS || x >= self.width || y >= self.height { + return None; + } + self.data.get(channel * self.plane_stride + y * self.row_stride + x) + } + + /// Copies one pixel from its channel planes. + #[must_use] + pub fn pixel(self, x: usize, y: usize) -> Option<[T; CHANNELS]> + where + T: Copy, + { + if x >= self.width || y >= self.height { + return None; + } + Some(std::array::from_fn(|channel| { + *self.get(channel, x, y).expect("validated planar coordinate") + })) + } + + /// Creates a checked zero-copy planar subview. + pub fn subview(self, region: ImageRegion) -> Result { + region.validate(self.width, self.height)?; + if region.width == 0 || region.height == 0 { + return Ok(Self { + width: region.width, + height: region.height, + row_stride: self.row_stride, + plane_stride: self.plane_stride, + data: &self.data[..0], + metadata: self.metadata, + }); + } + let start = region.y * self.row_stride + region.x; + let span = (CHANNELS - 1) * self.plane_stride + + (region.height - 1) * self.row_stride + + region.width; + Ok(Self { + width: region.width, + height: region.height, + row_stride: self.row_stride, + plane_stride: self.plane_stride, + data: &self.data[start..start + span], + metadata: self.metadata, + }) + } +} + +/// A one-channel image. +pub type GrayImage = Image; +/// A three-channel RGB image. +pub type RgbImage = Image; + +#[cfg(test)] +mod tests { + use super::{ + AlphaMode, ColorRange, ColorSpace, Image, ImageError, ImageMetadata, ImageRegion, + ImageView, ImageViewMut, PlanarImage, PlanarImageView, + }; + + #[test] + fn packed_image_indexes_pixels() { + let mut image = Image::::try_new(2, 1, vec![1, 2, 3, 4, 5, 6]).unwrap(); + assert_eq!(image[(1, 0)], [4, 5, 6]); + image[(0, 0)] = [7, 8, 9]; + assert_eq!(image.as_slice(), &[7, 8, 9, 4, 5, 6]); + } + + #[test] + fn strided_view_skips_padding() { + let data = [1_u16, 2, 99, 3, 4]; + let view = ImageView::::new(2, 2, 3, &data).unwrap(); + assert_eq!(view.get(0, 1), Some(&[3])); + assert_eq!(view.row(0), Some(&[1, 2][..])); + } + + #[test] + fn rejects_short_storage() { + assert_eq!( + Image::::try_new(2, 2, vec![0; 3]).unwrap_err(), + ImageError::StorageTooShort { required: 4, found: 3 } + ); + } + + #[test] + fn mutable_roi_updates_only_selected_pixels() { + let mut data = [0_u8, 1, 2, 99, 3, 4, 5]; + let mut view = ImageViewMut::::new(3, 2, 4, &mut data).unwrap(); + { + let mut roi = view.subview(ImageRegion::new(1, 0, 2, 2)).unwrap(); + *roi.get_mut(0, 0).unwrap() = [10]; + *roi.get_mut(1, 1).unwrap() = [20]; + } + assert_eq!(data, [0, 10, 2, 99, 3, 4, 20]); + } + + #[test] + fn immutable_roi_preserves_parent_stride() { + let data = [0_u8, 1, 2, 99, 3, 4, 5]; + let view = ImageView::::new(3, 2, 4, &data).unwrap(); + let roi = view.subview(ImageRegion::new(1, 0, 2, 2)).unwrap(); + assert_eq!(roi.row_stride(), 4); + assert_eq!(roi.get(0, 0), Some(&[1])); + assert_eq!(roi.get(1, 1), Some(&[5])); + assert!(view.subview(ImageRegion::new(2, 1, 2, 1)).is_err()); + } + + #[test] + fn planar_image_reads_channels_and_subviews() { + // R plane, G plane, B plane. + let mut image = PlanarImage::::try_new( + 2, + 2, + vec![1, 2, 3, 4, 10, 20, 30, 40, 100, 110, 120, 130], + ) + .unwrap(); + assert_eq!(image.view().pixel(1, 1), Some([4, 40, 130])); + *image.get_mut(1, 0, 1).unwrap() = 31; + let roi = image.view().subview(ImageRegion::new(0, 1, 2, 1)).unwrap(); + assert_eq!(roi.pixel(0, 0), Some([3, 31, 120])); + assert_eq!(roi.pixel(1, 0), Some([4, 40, 130])); + } + + #[test] + fn planar_view_honors_row_and_plane_padding() { + let data = [1_u8, 2, 99, 3, 4, 88, 77, 10, 20, 99, 30, 40]; + let view = PlanarImageView::::new(2, 2, 3, 7, &data).unwrap(); + assert_eq!(view.pixel(0, 1), Some([3, 30])); + assert_eq!(view.pixel(1, 1), Some([4, 40])); + } + + #[test] + fn validates_color_metadata_channel_count() { + let metadata = ImageMetadata { + color_space: ColorSpace::Rgb, + color_range: ColorRange::Full, + alpha_mode: AlphaMode::None, + }; + let image = Image::::try_new_with_metadata(1, 1, vec![1, 2, 3], metadata).unwrap(); + assert_eq!(image.metadata(), metadata); + assert!(matches!( + Image::::try_new_with_metadata(1, 1, vec![1], metadata), + Err(ImageError::MetadataChannelMismatch { .. }) + )); + } +} diff --git a/crates/spatialrust-py/Cargo.toml b/crates/spatialrust-py/Cargo.toml index 341957e..e0c1838 100644 --- a/crates/spatialrust-py/Cargo.toml +++ b/crates/spatialrust-py/Cargo.toml @@ -40,6 +40,8 @@ spatialrust = { path = "../spatialrust", features = [ "register-gicp", "register-ndt", "register-fpfh", + "camera-rgbd", + "vision-full", ] } # Keep this crate out of the main Rust workspace so `cargo test --workspace` diff --git a/crates/spatialrust-py/README.md b/crates/spatialrust-py/README.md index c5f053f..0752a40 100644 --- a/crates/spatialrust-py/README.md +++ b/crates/spatialrust-py/README.md @@ -84,6 +84,13 @@ reloaded = sr.read("labeled.las") | `farthest_point_sampling(cloud, sample_size, seed_index=0)` | Even FPS downsampling to a target count | | `voxelize(cloud, voxel_size=0.1, mode="occupancy")` | Dense 3D occupancy/count grid `(nz, ny, nx)` for ML | | `range_image(cloud, width=1024, height=64, fov_up_deg=3.0, fov_down_deg=-25.0)` | Spherical LiDAR range image `(height, width)` | +| `rgbd_to_point_cloud(depth, color, fx, fy, cx, cy, ...)` | Aligned `(H,W)` depth + `(H,W,3)` RGB to an XYZRGB cloud | +| `resize_image` / `letterbox_image` / `normalize_image_chw` | Model-ready RGB resize, padding, and float32 CHW packing | +| `rgb_to_gray_image` / `rgb_to_hsv_image` / `remap_image` | CPU color conversion and coordinate-map resampling | +| `nms` / `soft_nms` | Detection post-processing for `(N,4)` xyxy boxes | +| `connected_components_image` / `find_mask_contours` | Binary-mask labeling and contour extraction | +| `encode_mask_rle` / `decode_mask_rle` | Row-major or COCO column-major binary-mask RLE | +| `point_map_to_point_cloud` | Filter a dense `(H,W,3)` point map into a native point cloud | | `knn_graph(cloud, k)` / `radius_graph(cloud, radius)` | PyG-style `(2, E)` `edge_index` for GNNs | | `statistical_outlier_removal(cloud, k_neighbors=16, std_mul=1.0)` | Drop points far from their k-NN (SOR) | | `radius_outlier_removal(cloud, radius=0.5, min_neighbors=4)` | Drop points with too few neighbors in radius (ROR) | @@ -122,6 +129,8 @@ python examples/end_to_end.py --png demo.png # full clean->cluster->re python examples/make_gifs.py # rotating cluster + voxel GIFs python examples/ml_preprocess.py --png ml.png # point cloud -> ML tensors python examples/pyg_pointnet_demo.py --input scan.pcd # SpatialRust -> PyG Data -> tiny GCN +python examples/rgbd_pipeline.py # RGB-D -> colored cloud -> MVP +python examples/vision_ai_pipeline.py # image AI -> masks/points -> MVP ``` `segment_room.py` loads a real scan, runs the pipeline, and writes a labeled diff --git a/crates/spatialrust-py/examples/rgbd_pipeline.py b/crates/spatialrust-py/examples/rgbd_pipeline.py new file mode 100644 index 0000000..747c9e1 --- /dev/null +++ b/crates/spatialrust-py/examples/rgbd_pipeline.py @@ -0,0 +1,36 @@ +"""Synthetic RGB-D -> colored cloud -> SpatialRust MVP pipeline demo.""" + +import numpy as np +import spatialrust as sr + + +def main() -> None: + height, width = 48, 64 + depth = np.ones((height, width), dtype=np.float32) + depth[18:30, 26:38] = 0.7 + color = np.zeros((height, width, 3), dtype=np.uint8) + color[..., 1] = 160 + color[18:30, 26:38] = (230, 80, 40) + + cloud = sr.rgbd_to_point_cloud( + depth, + color, + fx=60.0, + fy=60.0, + cx=(width - 1) / 2, + cy=(height - 1) / 2, + ) + result = sr.run_pipeline( + cloud, + leaf_size=0.025, + plane_distance=0.02, + cluster_tolerance=0.08, + min_cluster_size=4, + ) + print(f"RGB-D points: {len(cloud)}") + print(f"plane inliers: {result.plane_inliers}") + print(f"clusters: {result.cluster_count}") + + +if __name__ == "__main__": + main() diff --git a/crates/spatialrust-py/examples/vision_ai_pipeline.py b/crates/spatialrust-py/examples/vision_ai_pipeline.py new file mode 100644 index 0000000..96a068d --- /dev/null +++ b/crates/spatialrust-py/examples/vision_ai_pipeline.py @@ -0,0 +1,70 @@ +"""End-to-end image preprocessing, post-processing, and spatial AI demo.""" + +from __future__ import annotations + +import numpy as np + +import spatialrust as sr + + +def main() -> None: + height, width = 48, 64 + yy, xx = np.mgrid[:height, :width] + image = np.stack( + [ + (xx * 4).clip(0, 255), + (yy * 5).clip(0, 255), + np.full_like(xx, 96), + ], + axis=-1, + ).astype(np.uint8) + + model_image, transform = sr.letterbox_image(image, 64, 64) + model_tensor = sr.normalize_image_chw( + model_image, + mean=(0.485, 0.456, 0.406), + std=(0.229, 0.224, 0.225), + ) + + boxes = np.array([[8, 8, 30, 28], [10, 9, 29, 27], [38, 30, 58, 52]], np.float32) + scores = np.array([0.94, 0.83, 0.88], np.float32) + kept = sr.nms(boxes, scores, score_threshold=0.25, iou_threshold=0.5) + + mask = np.zeros((height, width), dtype=np.uint8) + mask[8:24, 7:26] = 1 + mask[29:43, 39:58] = 1 + labels, components = sr.connected_components_image(mask) + runs = sr.encode_mask_rle(mask) + assert np.array_equal(sr.decode_mask_rle(width, height, runs), mask) + + z = np.ones((height, width), dtype=np.float32) + z[29:43, 39:58] = 0.75 + points = np.stack( + [ + (xx.astype(np.float32) - width / 2) * z / 60.0, + (yy.astype(np.float32) - height / 2) * z / 60.0, + z, + ], + axis=-1, + ) + confidence = np.where(mask != 0, 0.95, 0.8).astype(np.float32) + cloud = sr.point_map_to_point_cloud(points, confidence, min_confidence=0.5) + pipeline = sr.run_pipeline( + cloud, + leaf_size=0.025, + cluster_tolerance=0.1, + min_cluster_size=1, + plane_distance=0.02, + ) + + print(f"model tensor: {model_tensor.shape}, letterbox={transform}") + print(f"detections kept: {kept.tolist()}") + print(f"mask components: {len(components)}, labels={labels.max()}, RLE runs={len(runs)}") + print( + f"point cloud: {len(cloud)} points, plane inliers={pipeline.plane_inliers}, " + f"output={len(pipeline.output)}" + ) + + +if __name__ == "__main__": + main() diff --git a/crates/spatialrust-py/spatialrust.pyi b/crates/spatialrust-py/spatialrust.pyi index b91d152..68f2e24 100644 --- a/crates/spatialrust-py/spatialrust.pyi +++ b/crates/spatialrust-py/spatialrust.pyi @@ -15,7 +15,9 @@ __version__: str # Convenient aliases for the array shapes the bindings exchange. _F32Array = NDArray[np.float32] # positions, grids, range images, transforms _I32Array = NDArray[np.int32] # labels, edge_index +_U32Array = NDArray[np.uint32] _Vec3 = tuple[float, float, float] +_U8Array = NDArray[np.uint8] @final class PointCloud: @@ -41,6 +43,83 @@ class PointCloud: def __len__(self) -> int: ... def __repr__(self) -> str: ... +def rgbd_to_point_cloud( + depth: _F32Array, + color: _U8Array, + fx: float, + fy: float, + cx: float, + cy: float, + depth_scale: float = ..., + min_depth: float = ..., + max_depth: float = ..., + distortion: Optional[tuple[float, float, float, float, float]] = ..., +) -> PointCloud: + """Convert aligned depth and RGB images to an XYZRGB point cloud.""" + ... + +# --------------------------------------------------------------------------- # +# Image preprocessing, detection, masks, and dense spatial data +# --------------------------------------------------------------------------- # +def resize_image( + image: _U8Array, + width: int, + height: int, + interpolation: str = ..., +) -> _U8Array: ... +def letterbox_image( + image: _U8Array, + width: int, + height: int, + interpolation: str = ..., + fill: Optional[tuple[int, int, int]] = ..., +) -> tuple[_U8Array, tuple[float, int, int, int, int]]: ... +def normalize_image_chw( + image: _U8Array, + scale: float = ..., + mean: Optional[tuple[float, float, float]] = ..., + std: Optional[tuple[float, float, float]] = ..., +) -> _F32Array: ... +def rgb_to_gray_image(image: _U8Array) -> _U8Array: ... +def rgb_to_hsv_image(image: _U8Array) -> _U8Array: ... +def remap_image( + image: _U8Array, + map_x: _F32Array, + map_y: _F32Array, + interpolation: str = ..., + fill: Optional[tuple[int, int, int]] = ..., +) -> _U8Array: ... +def nms( + boxes: _F32Array, + scores: _F32Array, + score_threshold: float = ..., + iou_threshold: float = ..., +) -> NDArray[np.int64]: ... +def soft_nms( + boxes: _F32Array, + scores: _F32Array, + score_threshold: float = ..., + iou_threshold: float = ..., + method: str = ..., + sigma: float = ..., +) -> tuple[list[int], list[float]]: ... +def connected_components_image( + mask: _U8Array, + connectivity: int = ..., +) -> tuple[_U32Array, list[tuple[int, int, tuple[float, float, float, float]]]]: ... +def find_mask_contours( + mask: _U8Array, epsilon: float = ... +) -> list[list[tuple[int, int]]]: ... +def encode_mask_rle(mask: _U8Array, coco: bool = ...) -> list[int]: ... +def decode_mask_rle( + width: int, height: int, counts: Sequence[int], coco: bool = ... +) -> _U8Array: ... +def point_map_to_point_cloud( + points: _F32Array, + confidence: Optional[_F32Array] = ..., + min_confidence: float = ..., +) -> PointCloud: ... + @final class PipelineResult: """Result of running the MVP pipeline.""" diff --git a/crates/spatialrust-py/src/lib.rs b/crates/spatialrust-py/src/lib.rs index b44c328..8fd4d69 100644 --- a/crates/spatialrust-py/src/lib.rs +++ b/crates/spatialrust-py/src/lib.rs @@ -8,7 +8,9 @@ #![allow(clippy::useless_conversion)] use numpy::ndarray::{Array2, Array3}; -use numpy::{IntoPyArray, PyArray1, PyArray2, PyArray3, PyReadonlyArray2}; +use numpy::{ + IntoPyArray, PyArray1, PyArray2, PyArray3, PyReadonlyArray1, PyReadonlyArray2, PyReadonlyArray3, +}; use pyo3::exceptions::PyValueError; use pyo3::prelude::*; @@ -42,6 +44,15 @@ use spatialrust::transform::{ normalize_unit_sphere as normalize_unit, oriented_bounding_box as obb, recenter as recenter_op, scale_cloud, }; +use spatialrust::vision::{ + approximate_polygon as approximate_contour, connected_components as label_components, + decode_rle as decode_mask_runs, encode_rle as encode_mask_runs, + find_contours as trace_contours, letterbox as letterbox_op, nms as nms_op, + pack_chw as pack_chw_op, point_map_to_point_cloud as point_map_to_cloud, remap as remap_op, + resize as resize_op, rgb_to_gray as rgb_to_gray_op, rgb_to_hsv as rgb_to_hsv_op, + soft_nms as soft_nms_op, BinaryMask, BorderMode, BoundingBox2, ConfidenceMap, Connectivity, + Interpolation, MaskRle, PointMap, RleOrder, SoftNmsMethod, +}; use spatialrust::voxelize::{ range_image as range_image_proj, voxelize as voxelize_grid, RangeImageConfig, VoxelFill, VoxelGridConfig, @@ -53,6 +64,14 @@ use spatialrust::{ read_point_cloud_file, write_point_cloud_file, ExecutionPolicy, HasPositions3, PointCloud, StandardSchemas, }; +use spatialrust::{ + rgbd_to_point_cloud as rgbd_to_cloud, BrownConrady, CameraIntrinsics, DepthConversionOptions, + Image, PinholeCamera, +}; + +type Vec3Tuple = (f32, f32, f32); +type OrientedBoundingBoxTuple = (Vec3Tuple, Vec3Tuple, Vec); +type ComponentStats = Vec<(u32, usize, (f32, f32, f32, f32))>; fn to_py_err(err: E) -> PyErr { PyValueError::new_err(err.to_string()) @@ -69,6 +88,33 @@ fn parse_policy(policy: &str) -> PyResult { } } +fn parse_interpolation(interpolation: &str) -> PyResult { + match interpolation.to_lowercase().as_str() { + "nearest" => Ok(Interpolation::Nearest), + "bilinear" | "linear" => Ok(Interpolation::Bilinear), + "bicubic" | "cubic" => Ok(Interpolation::Bicubic), + "area" => Ok(Interpolation::Area), + other => Err(PyValueError::new_err(format!( + "unknown interpolation `{other}` (expected: nearest, bilinear, bicubic, area)" + ))), + } +} + +fn rgb_image_from_numpy(array: PyReadonlyArray3<'_, u8>) -> PyResult> { + let view = array.as_array(); + let shape = view.shape(); + if shape.len() != 3 || shape[2] != 3 { + return Err(PyValueError::new_err("expected an (H, W, 3) uint8 RGB array")); + } + Image::try_new(shape[1], shape[0], view.iter().copied().collect()).map_err(to_py_err) +} + +fn gray_u8_image_from_numpy(array: PyReadonlyArray2<'_, u8>) -> PyResult> { + let view = array.as_array(); + let shape = view.shape(); + Image::try_new(shape[1], shape[0], view.iter().copied().collect()).map_err(to_py_err) +} + fn cloud_from_xyz(arr: PyReadonlyArray2<'_, f32>) -> PyResult { let view = arr.as_array(); let shape = view.shape(); @@ -583,7 +629,7 @@ fn farthest_point_sampling( /// then flags points whose tangent-plane neighbors leave a large angular gap. /// Returns a sparse sub-cloud of the boundary points. #[pyfunction] -#[pyo3(signature = (cloud, search_radius=0.1, angle_threshold=1.5708, min_neighbors=5, k_neighbors=20))] +#[pyo3(signature = (cloud, search_radius=0.1, angle_threshold=std::f32::consts::FRAC_PI_2, min_neighbors=5, k_neighbors=20))] fn detect_boundary( cloud: &PyPointCloud, search_radius: f32, @@ -996,9 +1042,7 @@ fn bounding_box(cloud: &PyPointCloud) -> PyResult<((f32, f32, f32), (f32, f32, f /// Oriented (PCA) bounding box as `(center, half_extents, axes_3x3)`. The axes /// are returned principal-first; column `k` of `axes_3x3` is the k-th box axis. #[pyfunction] -fn oriented_bounding_box( - cloud: &PyPointCloud, -) -> PyResult<((f32, f32, f32), (f32, f32, f32), Vec<(f32, f32, f32)>)> { +fn oriented_bounding_box(cloud: &PyPointCloud) -> PyResult { let o = obb(&cloud.inner).map_err(to_py_err)?; let axis = |k: usize| (o.axes.m[0][k], o.axes.m[1][k], o.axes.m[2][k]); Ok(( @@ -1102,6 +1146,388 @@ fn range_image<'py>( Ok(arr.into_pyarray_bound(py)) } +/// Converts aligned `(H, W)` float32 depth and `(H, W, 3)` uint8 RGB images +/// into a colored point cloud. `depth_scale` converts stored values to meters. +#[pyfunction] +#[pyo3(signature = ( + depth, + color, + fx, + fy, + cx, + cy, + depth_scale=1.0, + min_depth=f32::EPSILON, + max_depth=f32::INFINITY, + distortion=None +))] +#[allow(clippy::too_many_arguments)] +fn rgbd_to_point_cloud( + depth: PyReadonlyArray2<'_, f32>, + color: PyReadonlyArray3<'_, u8>, + fx: f64, + fy: f64, + cx: f64, + cy: f64, + depth_scale: f32, + min_depth: f32, + max_depth: f32, + distortion: Option<(f64, f64, f64, f64, f64)>, +) -> PyResult { + let depth_view = depth.as_array(); + let color_view = color.as_array(); + let depth_shape = depth_view.shape(); + let color_shape = color_view.shape(); + let height = depth_shape[0]; + let width = depth_shape[1]; + if color_shape != [height, width, 3] { + return Err(PyValueError::new_err(format!( + "expected color shape ({height}, {width}, 3), found {:?}", + color_shape + ))); + } + + // Iteration follows logical ndarray order, so non-contiguous NumPy views + // are packed explicitly at the Python/native boundary. + let depth_image = Image::::try_new(width, height, depth_view.iter().copied().collect()) + .map_err(to_py_err)?; + let color_image = Image::::try_new(width, height, color_view.iter().copied().collect()) + .map_err(to_py_err)?; + let intrinsics = CameraIntrinsics::try_new(fx, fy, cx, cy, width, height).map_err(to_py_err)?; + let mut camera = PinholeCamera::new(intrinsics); + if let Some((k1, k2, p1, p2, k3)) = distortion { + camera = camera.with_distortion(BrownConrady { k1, k2, p1, p2, k3 }); + } + let options = DepthConversionOptions { depth_scale, min_depth, max_depth }; + let inner = rgbd_to_cloud(depth_image.view(), color_image.view(), &camera, options) + .map_err(to_py_err)?; + Ok(PyPointCloud { inner }) +} + +/// Resizes an `(H, W, 3)` uint8 RGB image. +#[pyfunction] +#[pyo3(signature = (image, width, height, interpolation="bilinear"))] +fn resize_image<'py>( + py: Python<'py>, + image: PyReadonlyArray3<'_, u8>, + width: usize, + height: usize, + interpolation: &str, +) -> PyResult>> { + let image = rgb_image_from_numpy(image)?; + let output = resize_op(image.view(), width, height, parse_interpolation(interpolation)?) + .map_err(to_py_err)?; + let array = Array3::from_shape_vec((height, width, 3), output.into_vec()).map_err(to_py_err)?; + Ok(array.into_pyarray_bound(py)) +} + +/// Letterboxes an RGB image and returns `(image, transform)`, where transform +/// is `(scale, pad_left, pad_top, content_width, content_height)`. +#[pyfunction] +#[pyo3(signature = (image, width, height, interpolation="bilinear", fill=None))] +fn letterbox_image<'py>( + py: Python<'py>, + image: PyReadonlyArray3<'_, u8>, + width: usize, + height: usize, + interpolation: &str, + fill: Option<(u8, u8, u8)>, +) -> PyResult<(Bound<'py, PyArray3>, (f64, usize, usize, usize, usize))> { + let image = rgb_image_from_numpy(image)?; + let (output, transform) = letterbox_op( + image.view(), + width, + height, + parse_interpolation(interpolation)?, + fill.map_or([114; 3], |(r, g, b)| [r, g, b]), + ) + .map_err(to_py_err)?; + let array = Array3::from_shape_vec((height, width, 3), output.into_vec()).map_err(to_py_err)?; + Ok(( + array.into_pyarray_bound(py), + ( + transform.scale, + transform.pad_left, + transform.pad_top, + transform.content_width, + transform.content_height, + ), + )) +} + +/// Normalizes RGB and packs it into a float32 `(3, H, W)` CHW tensor. +#[pyfunction] +#[pyo3(signature = (image, scale=1.0/255.0, mean=None, std=None))] +fn normalize_image_chw<'py>( + py: Python<'py>, + image: PyReadonlyArray3<'_, u8>, + scale: f32, + mean: Option<(f32, f32, f32)>, + std: Option<(f32, f32, f32)>, +) -> PyResult>> { + let image = rgb_image_from_numpy(image)?; + let mean = mean.map_or([0.0; 3], |(r, g, b)| [r, g, b]); + let std = std.map_or([1.0; 3], |(r, g, b)| [r, g, b]); + let output = pack_chw_op(image.view(), scale, mean, std).map_err(to_py_err)?; + let array = Array3::from_shape_vec((3, image.height(), image.width()), output.into_vec()) + .map_err(to_py_err)?; + Ok(array.into_pyarray_bound(py)) +} + +/// Converts an RGB image to an `(H, W)` grayscale image. +#[pyfunction] +fn rgb_to_gray_image<'py>( + py: Python<'py>, + image: PyReadonlyArray3<'_, u8>, +) -> PyResult>> { + let image = rgb_image_from_numpy(image)?; + let output = rgb_to_gray_op(image.view()).map_err(to_py_err)?; + let array = Array2::from_shape_vec((image.height(), image.width()), output.into_vec()) + .map_err(to_py_err)?; + Ok(array.into_pyarray_bound(py)) +} + +/// Converts RGB to OpenCV-style uint8 HSV. +#[pyfunction] +fn rgb_to_hsv_image<'py>( + py: Python<'py>, + image: PyReadonlyArray3<'_, u8>, +) -> PyResult>> { + let image = rgb_image_from_numpy(image)?; + let output = rgb_to_hsv_op(image.view()).map_err(to_py_err)?; + let array = Array3::from_shape_vec((image.height(), image.width(), 3), output.into_vec()) + .map_err(to_py_err)?; + Ok(array.into_pyarray_bound(py)) +} + +/// Remaps an RGB image with absolute float32 source-coordinate maps. +#[pyfunction] +#[pyo3(signature = (image, map_x, map_y, interpolation="bilinear", fill=None))] +fn remap_image<'py>( + py: Python<'py>, + image: PyReadonlyArray3<'_, u8>, + map_x: PyReadonlyArray2<'_, f32>, + map_y: PyReadonlyArray2<'_, f32>, + interpolation: &str, + fill: Option<(u8, u8, u8)>, +) -> PyResult>> { + let image = rgb_image_from_numpy(image)?; + let mx = map_x.as_array(); + let my = map_y.as_array(); + if mx.shape() != my.shape() { + return Err(PyValueError::new_err("map_x and map_y shapes must match")); + } + let height = mx.shape()[0]; + let width = mx.shape()[1]; + let map_x = + Image::::try_new(width, height, mx.iter().copied().collect()).map_err(to_py_err)?; + let map_y = + Image::::try_new(width, height, my.iter().copied().collect()).map_err(to_py_err)?; + let output = remap_op( + image.view(), + map_x.view(), + map_y.view(), + parse_interpolation(interpolation)?, + BorderMode::Constant(fill.map_or([0; 3], |(r, g, b)| [r, g, b])), + ) + .map_err(to_py_err)?; + let array = Array3::from_shape_vec((height, width, 3), output.into_vec()).map_err(to_py_err)?; + Ok(array.into_pyarray_bound(py)) +} + +/// Greedy non-maximum suppression over `(N, 4)` xyxy boxes. +#[pyfunction] +#[pyo3(signature = (boxes, scores, score_threshold=0.0, iou_threshold=0.5))] +fn nms<'py>( + py: Python<'py>, + boxes: PyReadonlyArray2<'_, f32>, + scores: PyReadonlyArray1<'_, f32>, + score_threshold: f32, + iou_threshold: f32, +) -> PyResult>> { + let boxes_view = boxes.as_array(); + if boxes_view.shape().len() != 2 || boxes_view.shape()[1] != 4 { + return Err(PyValueError::new_err("expected boxes with shape (N, 4)")); + } + let mut native_boxes = Vec::with_capacity(boxes_view.shape()[0]); + for row in boxes_view.rows() { + native_boxes + .push(BoundingBox2::try_new(row[0], row[1], row[2], row[3]).map_err(to_py_err)?); + } + let scores: Vec = scores.as_array().iter().copied().collect(); + let indices = nms_op(&native_boxes, &scores, score_threshold, iou_threshold) + .map_err(to_py_err)? + .into_iter() + .map(|index| index as i64) + .collect::>(); + Ok(indices.into_pyarray_bound(py)) +} + +/// Soft-NMS returning `(indices, updated_scores)`. +#[pyfunction] +#[pyo3(signature = (boxes, scores, score_threshold=0.001, iou_threshold=0.5, method="linear", sigma=0.5))] +fn soft_nms( + boxes: PyReadonlyArray2<'_, f32>, + scores: PyReadonlyArray1<'_, f32>, + score_threshold: f32, + iou_threshold: f32, + method: &str, + sigma: f32, +) -> PyResult<(Vec, Vec)> { + let boxes_view = boxes.as_array(); + if boxes_view.shape().len() != 2 || boxes_view.shape()[1] != 4 { + return Err(PyValueError::new_err("expected boxes with shape (N, 4)")); + } + let mut native_boxes = Vec::with_capacity(boxes_view.shape()[0]); + for row in boxes_view.rows() { + native_boxes + .push(BoundingBox2::try_new(row[0], row[1], row[2], row[3]).map_err(to_py_err)?); + } + let method = match method.to_lowercase().as_str() { + "hard" => SoftNmsMethod::Hard, + "linear" => SoftNmsMethod::Linear, + "gaussian" => SoftNmsMethod::Gaussian { sigma }, + other => return Err(PyValueError::new_err(format!("unknown Soft-NMS method `{other}`"))), + }; + let scores: Vec = scores.as_array().iter().copied().collect(); + let result = soft_nms_op(&native_boxes, &scores, score_threshold, iou_threshold, method) + .map_err(to_py_err)?; + Ok(( + result.iter().map(|value| value.index).collect(), + result.iter().map(|value| value.score).collect(), + )) +} + +/// Labels connected foreground regions in a uint8 binary mask. +#[pyfunction] +#[pyo3(signature = (mask, connectivity=8))] +fn connected_components_image<'py>( + py: Python<'py>, + mask: PyReadonlyArray2<'_, u8>, + connectivity: u8, +) -> PyResult<(Bound<'py, PyArray2>, ComponentStats)> { + let image = gray_u8_image_from_numpy(mask)?; + let mask = + BinaryMask::try_new(image.width(), image.height(), image.into_vec()).map_err(to_py_err)?; + let connectivity = match connectivity { + 4 => Connectivity::Four, + 8 => Connectivity::Eight, + _ => return Err(PyValueError::new_err("connectivity must be 4 or 8")), + }; + let result = label_components(&mask, connectivity).map_err(to_py_err)?; + let stats = result + .components + .iter() + .map(|component| { + ( + component.label, + component.area, + ( + component.bbox.x_min, + component.bbox.y_min, + component.bbox.x_max, + component.bbox.y_max, + ), + ) + }) + .collect(); + let labels = Array2::from_shape_vec( + (result.labels.height(), result.labels.width()), + result.labels.as_slice().to_vec(), + ) + .map_err(to_py_err)?; + Ok((labels.into_pyarray_bound(py), stats)) +} + +/// Extracts and optionally simplifies mask contours. +#[pyfunction] +#[pyo3(signature = (mask, epsilon=0.0))] +fn find_mask_contours( + mask: PyReadonlyArray2<'_, u8>, + epsilon: f64, +) -> PyResult>> { + let image = gray_u8_image_from_numpy(mask)?; + let mask = + BinaryMask::try_new(image.width(), image.height(), image.into_vec()).map_err(to_py_err)?; + trace_contours(&mask) + .into_iter() + .map(|contour| { + let contour = if epsilon > 0.0 { + approximate_contour(&contour, epsilon).map_err(to_py_err)? + } else { + contour + }; + Ok(contour.points.into_iter().map(|[x, y]| (x, y)).collect()) + }) + .collect() +} + +/// Encodes a binary mask into alternating run lengths. +#[pyfunction] +#[pyo3(signature = (mask, coco=true))] +fn encode_mask_rle(mask: PyReadonlyArray2<'_, u8>, coco: bool) -> PyResult> { + let image = gray_u8_image_from_numpy(mask)?; + let mask = + BinaryMask::try_new(image.width(), image.height(), image.into_vec()).map_err(to_py_err)?; + Ok(encode_mask_runs(&mask, if coco { RleOrder::CocoColumnMajor } else { RleOrder::RowMajor }) + .counts) +} + +/// Decodes alternating mask run lengths. +#[pyfunction] +#[pyo3(signature = (width, height, counts, coco=true))] +fn decode_mask_rle<'py>( + py: Python<'py>, + width: usize, + height: usize, + counts: Vec, + coco: bool, +) -> PyResult>> { + let rle = MaskRle { + width, + height, + order: if coco { RleOrder::CocoColumnMajor } else { RleOrder::RowMajor }, + counts, + }; + let mask = decode_mask_runs(&rle).map_err(to_py_err)?; + let array = + Array2::from_shape_vec((height, width), mask.into_image().into_vec()).map_err(to_py_err)?; + Ok(array.into_pyarray_bound(py)) +} + +/// Converts an `(H, W, 3)` float32 point map into a cloud. +#[pyfunction] +#[pyo3(signature = (points, confidence=None, min_confidence=0.0))] +fn point_map_to_point_cloud( + points: PyReadonlyArray3<'_, f32>, + confidence: Option>, + min_confidence: f32, +) -> PyResult { + let view = points.as_array(); + let shape = view.shape(); + if shape.len() != 3 || shape[2] != 3 { + return Err(PyValueError::new_err("expected point map shape (H, W, 3)")); + } + let point_map = + PointMap::try_new(shape[1], shape[0], view.iter().copied().collect()).map_err(to_py_err)?; + let confidence_map = if let Some(confidence) = confidence { + let confidence = confidence.as_array(); + Some( + ConfidenceMap::try_new( + confidence.shape()[1], + confidence.shape()[0], + confidence.iter().copied().collect(), + ) + .map_err(to_py_err)?, + ) + } else { + None + }; + let inner = point_map_to_cloud(&point_map, confidence_map.as_ref(), min_confidence) + .map_err(to_py_err)?; + Ok(PyPointCloud { inner }) +} + /// SpatialRust Python bindings. #[pymodule] #[pyo3(name = "spatialrust")] @@ -1147,6 +1573,20 @@ fn spatialrust_module(m: &Bound<'_, PyModule>) -> PyResult<()> { m.add_function(wrap_pyfunction!(oriented_bounding_box, m)?)?; m.add_function(wrap_pyfunction!(voxelize, m)?)?; m.add_function(wrap_pyfunction!(range_image, m)?)?; + m.add_function(wrap_pyfunction!(rgbd_to_point_cloud, m)?)?; + m.add_function(wrap_pyfunction!(resize_image, m)?)?; + m.add_function(wrap_pyfunction!(letterbox_image, m)?)?; + m.add_function(wrap_pyfunction!(normalize_image_chw, m)?)?; + m.add_function(wrap_pyfunction!(rgb_to_gray_image, m)?)?; + m.add_function(wrap_pyfunction!(rgb_to_hsv_image, m)?)?; + m.add_function(wrap_pyfunction!(remap_image, m)?)?; + 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!(find_mask_contours, m)?)?; + m.add_function(wrap_pyfunction!(encode_mask_rle, m)?)?; + m.add_function(wrap_pyfunction!(decode_mask_rle, m)?)?; + m.add_function(wrap_pyfunction!(point_map_to_point_cloud, m)?)?; m.add_function(wrap_pyfunction!(knn_graph, m)?)?; m.add_function(wrap_pyfunction!(radius_graph, m)?)?; m.add_function(wrap_pyfunction!(register_icp, m)?)?; diff --git a/crates/spatialrust-py/tests/test_bindings.py b/crates/spatialrust-py/tests/test_bindings.py index 42b15c0..540a66f 100644 --- a/crates/spatialrust-py/tests/test_bindings.py +++ b/crates/spatialrust-py/tests/test_bindings.py @@ -63,6 +63,12 @@ def test_exports_present(): for name in ( "PointCloud", "voxel_downsample", "dbscan", "register_icp", "voxelize", "knn_graph", "chamfer_distance", "oriented_bounding_box", + "rgbd_to_point_cloud", + "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", "encode_mask_rle", "decode_mask_rle", + "point_map_to_point_cloud", ): assert hasattr(sr, name), f"missing export: {name}" @@ -91,6 +97,95 @@ def test_unlabeled_cloud_has_no_labels(plane): assert plane.labels() is None +def test_rgbd_to_point_cloud(): + depth = np.array([[1.0, 0.0], [np.nan, 2.0]], dtype=np.float32) + color = np.array( + [[[10, 11, 12], [20, 21, 22]], [[30, 31, 32], [40, 41, 42]]], + dtype=np.uint8, + ) + cloud = sr.rgbd_to_point_cloud(depth, color, 2.0, 2.0, 0.0, 0.0) + assert len(cloud) == 2 + assert {"x", "y", "z", "r", "g", "b"}.issubset(set(cloud.field_names())) + np.testing.assert_allclose( + cloud.xyz(), np.array([[0.0, 0.0, 1.0], [1.0, 1.0, 2.0]], dtype=np.float32) + ) + + +def test_image_resize_letterbox_and_normalize(): + image = np.array( + [[[255, 0, 0], [0, 255, 0]], [[0, 0, 255], [255, 255, 255]]], + dtype=np.uint8, + ) + resized = sr.resize_image(image, 4, 4, interpolation="nearest") + assert resized.shape == (4, 4, 3) + np.testing.assert_array_equal(resized[0, 0], image[0, 0]) + + letterboxed, transform = sr.letterbox_image(image, 4, 6, fill=(7, 8, 9)) + assert letterboxed.shape == (6, 4, 3) + assert transform == (2.0, 0, 1, 4, 4) + np.testing.assert_array_equal(letterboxed[0, 0], [7, 8, 9]) + + chw = sr.normalize_image_chw(image) + assert chw.shape == (3, 2, 2) + assert chw.dtype == np.float32 + np.testing.assert_allclose(chw[:, 0, 0], [1.0, 0.0, 0.0], atol=1e-6) + + +def test_image_color_and_remap(): + image = np.array([[[255, 0, 0], [0, 255, 0]]], dtype=np.uint8) + gray = sr.rgb_to_gray_image(image) + assert gray.shape == (1, 2) + np.testing.assert_allclose(gray, [[76, 150]], atol=1) + hsv = sr.rgb_to_hsv_image(image) + np.testing.assert_array_equal(hsv[0, 0], [0, 255, 255]) + np.testing.assert_array_equal(hsv[0, 1], [60, 255, 255]) + + map_x = np.array([[0.0, 1.0]], dtype=np.float32) + map_y = np.zeros((1, 2), dtype=np.float32) + remapped = sr.remap_image(image, map_x, map_y, interpolation="nearest") + np.testing.assert_array_equal(remapped, image) + + +def test_detection_nms_and_soft_nms(): + boxes = np.array( + [[0, 0, 10, 10], [1, 1, 9, 9], [20, 20, 30, 30]], dtype=np.float32 + ) + scores = np.array([0.9, 0.8, 0.7], dtype=np.float32) + np.testing.assert_array_equal(sr.nms(boxes, scores), [0, 2]) + indices, updated = sr.soft_nms(boxes, scores, method="linear") + assert indices[0] == 0 + assert len(indices) == len(updated) == 3 + assert updated[-1] < 0.8 + + +def test_mask_components_contours_and_rle(): + mask = np.zeros((5, 7), dtype=np.uint8) + mask[1:3, 1:3] = 1 + mask[2:4, 5:7] = 1 + labels, stats = sr.connected_components_image(mask, connectivity=4) + assert labels.shape == mask.shape + assert labels.dtype == np.uint32 + assert sorted(stat[1] for stat in stats) == [4, 4] + contours = sr.find_mask_contours(mask) + assert len(contours) == 2 + + for coco in (False, True): + counts = sr.encode_mask_rle(mask, coco=coco) + decoded = sr.decode_mask_rle(7, 5, counts, coco=coco) + np.testing.assert_array_equal(decoded, mask) + + +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 + ) + confidence = np.array([[0.9, 0.2], [1.0, 0.8]], dtype=np.float32) + cloud = sr.point_map_to_point_cloud(points, confidence, min_confidence=0.5) + np.testing.assert_allclose( + cloud.xyz(), np.array([[0, 0, 1], [1, 1, 2]], dtype=np.float32) + ) + + # --------------------------------------------------------------------------- # # Filters # --------------------------------------------------------------------------- # diff --git a/crates/spatialrust-vision/Cargo.toml b/crates/spatialrust-vision/Cargo.toml new file mode 100644 index 0000000..0769aef --- /dev/null +++ b/crates/spatialrust-vision/Cargo.toml @@ -0,0 +1,35 @@ +[package] +name = "spatialrust-vision" +version.workspace = true +edition.workspace = true +license.workspace = true +authors.workspace = true +repository.workspace = true +rust-version.workspace = true +description = "AI-ready image processing and vision algorithms for SpatialRust" + +[features] +default = [] +resize = [] +preprocess = ["resize"] +warp = ["resize"] +detection = [] +dense = ["detection"] +spatial = ["dense", "dep:spatialrust-core", "dep:spatialrust-camera"] +full = ["preprocess", "warp", "detection", "dense", "spatial"] + +[dependencies] +spatialrust-image.workspace = true +spatialrust-math.workspace = true +spatialrust-core = { workspace = true, optional = true } +spatialrust-camera = { workspace = true, optional = true } +thiserror.workspace = true + +[dev-dependencies] +criterion.workspace = true +proptest.workspace = true + +[[bench]] +name = "preprocess" +harness = false +required-features = ["preprocess"] diff --git a/crates/spatialrust-vision/benches/preprocess.rs b/crates/spatialrust-vision/benches/preprocess.rs new file mode 100644 index 0000000..0a0258d --- /dev/null +++ b/crates/spatialrust-vision/benches/preprocess.rs @@ -0,0 +1,18 @@ +use criterion::{black_box, criterion_group, criterion_main, Criterion}; +use spatialrust_image::Image; +use spatialrust_vision::{letterbox, pack_chw, Interpolation}; + +fn benchmark_preprocess(c: &mut Criterion) { + let input = Image::::try_new(1280, 720, vec![127; 1280 * 720 * 3]).unwrap(); + c.bench_function("letterbox_normalize_chw_1280x720_to_640", |b| { + b.iter(|| { + let (resized, _) = + letterbox(black_box(input.view()), 640, 640, Interpolation::Bilinear, [114; 3]) + .unwrap(); + pack_chw(resized.view(), 1.0 / 255.0, [0.0; 3], [1.0; 3]).unwrap() + }); + }); +} + +criterion_group!(benches, benchmark_preprocess); +criterion_main!(benches); diff --git a/crates/spatialrust-vision/src/dense.rs b/crates/spatialrust-vision/src/dense.rs new file mode 100644 index 0000000..02fde7f --- /dev/null +++ b/crates/spatialrust-vision/src/dense.rs @@ -0,0 +1,685 @@ +//! Dense mask, depth, flow, confidence, and point-map primitives. + +use std::collections::{BTreeMap, VecDeque}; + +use spatialrust_image::{ColorSpace, Image, ImageMetadata, ImageView}; + +use crate::{BoundingBox2, VisionError, VisionResult}; + +/// Pixel connectivity used by binary-mask algorithms. +#[derive(Clone, Copy, Debug, Default, PartialEq, Eq, Hash)] +pub enum Connectivity { + /// Horizontal and vertical neighbors. + Four, + /// Horizontal, vertical, and diagonal neighbors. + #[default] + Eight, +} + +/// Validated binary mask (`0` background, `1` foreground). +#[derive(Clone, Debug, PartialEq, Eq)] +pub struct BinaryMask { + image: Image, +} + +impl BinaryMask { + /// Creates a mask and rejects values other than zero and one. + pub fn try_new(width: usize, height: usize, data: Vec) -> VisionResult { + if data.iter().any(|&value| value > 1) { + return Err(VisionError::InvalidParameter( + "binary mask values must be zero or one".to_owned(), + )); + } + let metadata = ImageMetadata { color_space: ColorSpace::Label, ..Default::default() }; + Ok(Self { image: Image::try_new_with_metadata(width, height, data, metadata)? }) + } + + /// Thresholds a float score image into a binary mask. + pub fn from_threshold(input: ImageView<'_, f32, 1>, threshold: f32) -> VisionResult { + if !threshold.is_finite() { + return Err(VisionError::InvalidParameter("mask threshold must be finite".to_owned())); + } + let mut data = Vec::with_capacity(input.width() * input.height()); + for y in 0..input.height() { + for x in 0..input.width() { + let value = input.get(x, y).expect("coordinate in bounds")[0]; + data.push(u8::from(value.is_finite() && value >= threshold)); + } + } + Self::try_new(input.width(), input.height(), data) + } + + /// Returns mask width. + #[must_use] + pub const fn width(&self) -> usize { + self.image.width() + } + + /// Returns mask height. + #[must_use] + pub const fn height(&self) -> usize { + self.image.height() + } + + /// Borrows the underlying image. + #[must_use] + pub fn image(&self) -> &Image { + &self.image + } + + /// Borrows a mask view. + #[must_use] + pub fn view(&self) -> ImageView<'_, u8, 1> { + self.image.view() + } + + /// Returns whether one pixel is foreground. + #[must_use] + pub fn contains(&self, x: usize, y: usize) -> bool { + self.image.get(x, y).is_some_and(|pixel| pixel[0] != 0) + } + + /// Counts foreground pixels. + #[must_use] + pub fn area(&self) -> usize { + self.image.as_slice().iter().filter(|&&value| value != 0).count() + } + + /// Consumes the wrapper and returns its image. + #[must_use] + pub fn into_image(self) -> Image { + self.image + } +} + +/// Connected-component label image (`0` is background). +#[derive(Clone, Debug, PartialEq, Eq)] +pub struct LabelImage { + image: Image, +} + +impl LabelImage { + /// Returns label image width. + #[must_use] + pub const fn width(&self) -> usize { + self.image.width() + } + + /// Returns label image height. + #[must_use] + pub const fn height(&self) -> usize { + self.image.height() + } + + /// Returns packed labels. + #[must_use] + pub fn as_slice(&self) -> &[u32] { + self.image.as_slice() + } + + /// Returns one label. + #[must_use] + pub fn get(&self, x: usize, y: usize) -> Option { + self.image.get(x, y).map(|pixel| pixel[0]) + } + + /// Borrows the underlying image. + #[must_use] + pub fn image(&self) -> &Image { + &self.image + } +} + +/// Per-component statistics. +#[derive(Clone, Copy, Debug, PartialEq)] +pub struct ComponentStats { + /// Positive component label. + pub label: u32, + /// Foreground pixel count. + pub area: usize, + /// Half-open pixel bounding box. + pub bbox: BoundingBox2, + /// Mean pixel-center coordinate. + pub centroid: [f64; 2], +} + +/// Connected-component labeling output. +#[derive(Clone, Debug, PartialEq)] +pub struct ConnectedComponents { + /// Per-pixel labels. + pub labels: LabelImage, + /// Statistics ordered by positive label. + pub components: Vec, +} + +/// Labels foreground components and computes bounding boxes and centroids. +pub fn connected_components( + mask: &BinaryMask, + connectivity: Connectivity, +) -> VisionResult { + let width = mask.width(); + let height = mask.height(); + let mut labels = vec![0_u32; width * height]; + let mut components = Vec::new(); + let mut queue = VecDeque::new(); + let mut next_label = 1_u32; + for seed_y in 0..height { + for seed_x in 0..width { + let seed = seed_y * width + seed_x; + if !mask.contains(seed_x, seed_y) || labels[seed] != 0 { + continue; + } + labels[seed] = next_label; + queue.push_back((seed_x, seed_y)); + let mut area = 0_usize; + let mut min_x = seed_x; + let mut min_y = seed_y; + let mut max_x = seed_x; + let mut max_y = seed_y; + let mut sum_x = 0.0_f64; + let mut sum_y = 0.0_f64; + while let Some((x, y)) = queue.pop_front() { + area += 1; + min_x = min_x.min(x); + min_y = min_y.min(y); + max_x = max_x.max(x); + max_y = max_y.max(y); + sum_x += x as f64 + 0.5; + sum_y += y as f64 + 0.5; + for (nx, ny) in neighbors(x, y, width, height, connectivity) { + let index = ny * width + nx; + if mask.contains(nx, ny) && labels[index] == 0 { + labels[index] = next_label; + queue.push_back((nx, ny)); + } + } + } + components.push(ComponentStats { + label: next_label, + area, + bbox: BoundingBox2 { + x_min: min_x as f32, + y_min: min_y as f32, + x_max: (max_x + 1) as f32, + y_max: (max_y + 1) as f32, + }, + centroid: [sum_x / area as f64, sum_y / area as f64], + }); + next_label = next_label + .checked_add(1) + .ok_or_else(|| VisionError::InvalidDimensions("too many components".to_owned()))?; + } + } + let metadata = ImageMetadata { color_space: ColorSpace::Label, ..Default::default() }; + let image = Image::try_new_with_metadata(width, height, labels, metadata)?; + Ok(ConnectedComponents { labels: LabelImage { image }, components }) +} + +fn neighbors( + x: usize, + y: usize, + width: usize, + height: usize, + connectivity: Connectivity, +) -> impl Iterator { + let offsets: &[(isize, isize)] = match connectivity { + Connectivity::Four => &[(0, -1), (-1, 0), (1, 0), (0, 1)], + Connectivity::Eight => { + &[(-1, -1), (0, -1), (1, -1), (-1, 0), (1, 0), (-1, 1), (0, 1), (1, 1)] + } + }; + offsets.iter().filter_map(move |&(dx, dy)| { + let nx = x.checked_add_signed(dx)?; + let ny = y.checked_add_signed(dy)?; + (nx < width && ny < height).then_some((nx, ny)) + }) +} + +/// One closed polygonal contour on pixel-grid corner coordinates. +#[derive(Clone, Debug, PartialEq, Eq)] +pub struct Contour { + /// Ordered vertices. The first vertex is not repeated at the end. + pub points: Vec<[i32; 2]>, +} + +/// Extracts oriented boundary loops, including hole contours. +pub fn find_contours(mask: &BinaryMask) -> Vec { + type Point = (i32, i32); + let mut edges: BTreeMap> = BTreeMap::new(); + let foreground = |x: isize, y: isize| { + x >= 0 + && y >= 0 + && (x as usize) < mask.width() + && (y as usize) < mask.height() + && mask.contains(x as usize, y as usize) + }; + for y in 0..mask.height() as isize { + for x in 0..mask.width() as isize { + if !foreground(x, y) { + continue; + } + let x = x as i32; + let y = y as i32; + if !foreground(x as isize, y as isize - 1) { + add_edge(&mut edges, (x, y), (x + 1, y)); + } + if !foreground(x as isize + 1, y as isize) { + add_edge(&mut edges, (x + 1, y), (x + 1, y + 1)); + } + if !foreground(x as isize, y as isize + 1) { + add_edge(&mut edges, (x + 1, y + 1), (x, y + 1)); + } + if !foreground(x as isize - 1, y as isize) { + add_edge(&mut edges, (x, y + 1), (x, y)); + } + } + } + + let mut contours = Vec::new(); + while let Some((&start, _)) = edges.iter().next() { + let mut current = start; + let mut points = Vec::new(); + loop { + points.push([current.0, current.1]); + let Some(next) = take_edge(&mut edges, current) else { + break; + }; + current = next; + if current == start { + break; + } + } + if points.len() >= 3 { + contours.push(Contour { points }); + } + } + contours +} + +fn add_edge(edges: &mut BTreeMap<(i32, i32), Vec<(i32, i32)>>, start: (i32, i32), end: (i32, i32)) { + edges.entry(start).or_default().push(end); +} + +fn take_edge( + edges: &mut BTreeMap<(i32, i32), Vec<(i32, i32)>>, + start: (i32, i32), +) -> Option<(i32, i32)> { + let values = edges.get_mut(&start)?; + values.sort_unstable(); + let end = values.remove(0); + if values.is_empty() { + edges.remove(&start); + } + Some(end) +} + +/// Simplifies a closed contour with Ramer–Douglas–Peucker approximation. +pub fn approximate_polygon(contour: &Contour, epsilon: f64) -> VisionResult { + if !epsilon.is_finite() || epsilon < 0.0 { + return Err(VisionError::InvalidParameter( + "polygon epsilon must be finite and non-negative".to_owned(), + )); + } + if contour.points.len() <= 3 || epsilon == 0.0 { + return Ok(contour.clone()); + } + // Rotate at the point farthest from vertex 0, then simplify the two open + // arcs independently so a closed loop has stable endpoints. + let first = contour.points[0]; + let split = contour + .points + .iter() + .enumerate() + .max_by_key(|(_, point)| squared_distance(**point, first)) + .map_or(0, |(index, _)| index); + let mut left = rdp(&contour.points[..=split], epsilon); + let mut wrapped = contour.points[split..].to_vec(); + wrapped.push(first); + let right = rdp(&wrapped, epsilon); + left.extend(right.into_iter().skip(1).take_while(|point| *point != first)); + Ok(Contour { points: left }) +} + +fn rdp(points: &[[i32; 2]], epsilon: f64) -> Vec<[i32; 2]> { + if points.len() <= 2 { + return points.to_vec(); + } + let first = points[0]; + let last = points[points.len() - 1]; + let mut max_distance = 0.0; + let mut split = 0; + for (index, &point) in points.iter().enumerate().take(points.len() - 1).skip(1) { + let distance = point_segment_distance(point, first, last); + if distance > max_distance { + max_distance = distance; + split = index; + } + } + if max_distance <= epsilon { + return vec![first, last]; + } + let mut left = rdp(&points[..=split], epsilon); + let right = rdp(&points[split..], epsilon); + left.extend(right.into_iter().skip(1)); + left +} + +fn squared_distance(a: [i32; 2], b: [i32; 2]) -> i64 { + let dx = i64::from(a[0]) - i64::from(b[0]); + let dy = i64::from(a[1]) - i64::from(b[1]); + dx * dx + dy * dy +} + +fn point_segment_distance(point: [i32; 2], start: [i32; 2], end: [i32; 2]) -> f64 { + let px = f64::from(point[0]); + let py = f64::from(point[1]); + let sx = f64::from(start[0]); + let sy = f64::from(start[1]); + let dx = f64::from(end[0] - start[0]); + let dy = f64::from(end[1] - start[1]); + let length_squared = dx.mul_add(dx, dy * dy); + if length_squared == 0.0 { + return (px - sx).hypot(py - sy); + } + let t = ((px - sx).mul_add(dx, (py - sy) * dy) / length_squared).clamp(0.0, 1.0); + (px - sx - t * dx).hypot(py - sy - t * dy) +} + +/// Linearization order used by run-length encoded masks. +#[derive(Clone, Copy, Debug, Default, PartialEq, Eq, Hash)] +pub enum RleOrder { + /// Conventional row-major traversal. + #[default] + RowMajor, + /// COCO-compatible column-major traversal. + CocoColumnMajor, +} + +/// Alternating zero/one run lengths, always beginning with a zero run. +#[derive(Clone, Debug, PartialEq, Eq)] +pub struct MaskRle { + /// Mask width. + pub width: usize, + /// Mask height. + pub height: usize, + /// Linearization order. + pub order: RleOrder, + /// Alternating zero and one run lengths. + pub counts: Vec, +} + +/// Encodes a binary mask into alternating runs. +#[must_use] +pub fn encode_rle(mask: &BinaryMask, order: RleOrder) -> MaskRle { + let mut counts = vec![0_usize]; + let mut current = 0_u8; + for index in 0..mask.width() * mask.height() { + let (x, y) = linear_coordinate(index, mask.width(), mask.height(), order); + let value = u8::from(mask.contains(x, y)); + if value == current { + *counts.last_mut().expect("initial run exists") += 1; + } else { + counts.push(1); + current = value; + } + } + MaskRle { width: mask.width(), height: mask.height(), order, counts } +} + +/// Decodes a run-length mask and validates its total length. +pub fn decode_rle(rle: &MaskRle) -> VisionResult { + let total = rle + .width + .checked_mul(rle.height) + .ok_or_else(|| VisionError::InvalidDimensions("RLE dimensions overflow".to_owned()))?; + let encoded = rle.counts.iter().try_fold(0_usize, |sum, &count| sum.checked_add(count)); + if encoded != Some(total) { + return Err(VisionError::InvalidParameter( + "RLE run lengths do not match mask dimensions".to_owned(), + )); + } + let mut data = vec![0_u8; total]; + let mut linear = 0; + let mut value = 0_u8; + for &count in &rle.counts { + for _ in 0..count { + let (x, y) = linear_coordinate(linear, rle.width, rle.height, rle.order); + data[y * rle.width + x] = value; + linear += 1; + } + value ^= 1; + } + BinaryMask::try_new(rle.width, rle.height, data) +} + +fn linear_coordinate(index: usize, width: usize, height: usize, order: RleOrder) -> (usize, usize) { + match order { + RleOrder::RowMajor => (index % width, index / width), + RleOrder::CocoColumnMajor => (index / height, index % height), + } +} + +/// Metric or relative dense depth map. Non-positive and non-finite values are invalid. +#[derive(Clone, Debug, PartialEq)] +pub struct DepthMap { + image: Image, +} + +impl DepthMap { + /// Wraps a depth image and marks its semantic metadata. + pub fn try_new(width: usize, height: usize, data: Vec) -> VisionResult { + if data.iter().any(|value| value.is_infinite()) { + return Err(VisionError::InvalidParameter( + "depth values may be finite or NaN, but not infinite".to_owned(), + )); + } + let metadata = ImageMetadata { color_space: ColorSpace::Depth, ..Default::default() }; + Ok(Self { image: Image::try_new_with_metadata(width, height, data, metadata)? }) + } + + /// Borrows the depth image. + #[must_use] + pub fn image(&self) -> &Image { + &self.image + } + + /// Builds a validity mask for an inclusive depth range. + pub fn valid_mask(&self, min_depth: f32, max_depth: f32) -> VisionResult { + if !min_depth.is_finite() || min_depth < 0.0 || max_depth.is_nan() || max_depth < min_depth + { + return Err(VisionError::InvalidParameter("invalid depth range".to_owned())); + } + BinaryMask::try_new( + self.image.width(), + self.image.height(), + self.image + .as_slice() + .iter() + .map(|&value| { + u8::from(value.is_finite() && value >= min_depth && value <= max_depth) + }) + .collect(), + ) + } +} + +/// Dense confidence values constrained to `[0, 1]`. +#[derive(Clone, Debug, PartialEq)] +pub struct ConfidenceMap { + image: Image, +} + +impl ConfidenceMap { + /// Creates a validated confidence map. + pub fn try_new(width: usize, height: usize, data: Vec) -> VisionResult { + if data.iter().any(|value| !value.is_finite() || !(0.0..=1.0).contains(value)) { + return Err(VisionError::InvalidParameter( + "confidence values must be finite and in [0, 1]".to_owned(), + )); + } + Ok(Self { image: Image::try_new(width, height, data)? }) + } + + /// Borrows the confidence image. + #[must_use] + pub fn image(&self) -> &Image { + &self.image + } + + /// Returns confidence map width. + #[must_use] + pub const fn width(&self) -> usize { + self.image.width() + } + + /// Returns confidence map height. + #[must_use] + pub const fn height(&self) -> usize { + self.image.height() + } +} + +/// Dense `(dx, dy)` optical-flow field in pixel units. +#[derive(Clone, Debug, PartialEq)] +pub struct FlowField { + image: Image, +} + +impl FlowField { + /// Creates a flow field. NaN marks invalid flow; infinity is rejected. + pub fn try_new(width: usize, height: usize, data: Vec) -> VisionResult { + if data.iter().any(|value| value.is_infinite()) { + return Err(VisionError::InvalidParameter( + "flow values may be finite or NaN, but not infinite".to_owned(), + )); + } + Ok(Self { image: Image::try_new(width, height, data)? }) + } + + /// Borrows the flow image. + #[must_use] + pub fn image(&self) -> &Image { + &self.image + } + + /// Converts displacement vectors to absolute remap coordinates. + pub fn to_remap(&self) -> VisionResult<(Image, Image)> { + let mut map_x = Vec::with_capacity(self.image.width() * self.image.height()); + let mut map_y = Vec::with_capacity(self.image.width() * self.image.height()); + for y in 0..self.image.height() { + for x in 0..self.image.width() { + let flow = self.image.get(x, y).expect("coordinate in bounds"); + map_x.push(x as f32 + flow[0]); + map_y.push(y as f32 + flow[1]); + } + } + Ok(( + Image::try_new(self.image.width(), self.image.height(), map_x)?, + Image::try_new(self.image.width(), self.image.height(), map_y)?, + )) + } +} + +/// Dense per-pixel XYZ point map. Non-finite triples represent invalid points. +#[derive(Clone, Debug, PartialEq)] +pub struct PointMap { + image: Image, +} + +impl PointMap { + /// Creates a point map and rejects infinite components. + pub fn try_new(width: usize, height: usize, data: Vec) -> VisionResult { + if data.iter().any(|value| value.is_infinite()) { + return Err(VisionError::InvalidParameter( + "point-map values may be finite or NaN, but not infinite".to_owned(), + )); + } + Ok(Self { image: Image::try_new(width, height, data)? }) + } + + /// Borrows the XYZ image. + #[must_use] + pub fn image(&self) -> &Image { + &self.image + } + + /// Returns point-map width. + #[must_use] + pub const fn width(&self) -> usize { + self.image.width() + } + + /// Returns point-map height. + #[must_use] + pub const fn height(&self) -> usize { + self.image.height() + } + + /// Returns a mask where all XYZ components are finite. + pub fn valid_mask(&self) -> VisionResult { + let mut data = Vec::with_capacity(self.image.width() * self.image.height()); + for y in 0..self.image.height() { + for x in 0..self.image.width() { + let point = self.image.get(x, y).expect("coordinate in bounds"); + data.push(u8::from(point.iter().all(|value| value.is_finite()))); + } + } + BinaryMask::try_new(self.image.width(), self.image.height(), data) + } +} + +#[cfg(test)] +mod tests { + use super::{ + approximate_polygon, connected_components, decode_rle, encode_rle, find_contours, + BinaryMask, Connectivity, DepthMap, FlowField, PointMap, RleOrder, + }; + use spatialrust_image::Image; + + #[test] + fn threshold_and_components_find_two_regions() { + let scores = Image::::try_new( + 4, + 3, + vec![1.0, 1.0, 0.0, 0.0, 1.0, 1.0, 0.0, 1.0, 0.0, 0.0, 0.0, 1.0], + ) + .unwrap(); + let mask = BinaryMask::from_threshold(scores.view(), 0.5).unwrap(); + let result = connected_components(&mask, Connectivity::Four).unwrap(); + assert_eq!(result.components.len(), 2); + assert_eq!(result.components[0].area, 4); + assert_eq!(result.components[1].area, 2); + assert_eq!(result.components[0].centroid, [1.0, 1.0]); + } + + #[test] + fn contours_trace_outer_and_hole_loops() { + let mask = BinaryMask::try_new(3, 3, vec![1, 1, 1, 1, 0, 1, 1, 1, 1]).unwrap(); + let contours = find_contours(&mask); + assert_eq!(contours.len(), 2); + assert!(contours.iter().all(|contour| contour.points.len() >= 4)); + let simplified = approximate_polygon(&contours[0], 0.1).unwrap(); + assert!(simplified.points.len() <= contours[0].points.len()); + } + + #[test] + fn rle_roundtrips_both_orders() { + let mask = BinaryMask::try_new(3, 2, vec![0, 1, 1, 1, 0, 1]).unwrap(); + for order in [RleOrder::RowMajor, RleOrder::CocoColumnMajor] { + let encoded = encode_rle(&mask, order); + assert_eq!(decode_rle(&encoded).unwrap(), mask); + } + } + + #[test] + fn depth_flow_and_point_maps_expose_validity() { + let depth = DepthMap::try_new(2, 1, vec![1.0, f32::NAN]).unwrap(); + assert_eq!(depth.valid_mask(0.1, 10.0).unwrap().image().as_slice(), &[1, 0]); + let flow = FlowField::try_new(2, 1, vec![1.0, 2.0, -1.0, 0.0]).unwrap(); + let (mx, my) = flow.to_remap().unwrap(); + assert_eq!(mx.as_slice(), &[1.0, 0.0]); + assert_eq!(my.as_slice(), &[2.0, 0.0]); + let points = PointMap::try_new(2, 1, vec![0.0, 0.0, 1.0, f32::NAN, 0.0, 1.0]).unwrap(); + assert_eq!(points.valid_mask().unwrap().image().as_slice(), &[1, 0]); + } +} diff --git a/crates/spatialrust-vision/src/detection.rs b/crates/spatialrust-vision/src/detection.rs new file mode 100644 index 0000000..a9860b0 --- /dev/null +++ b/crates/spatialrust-vision/src/detection.rs @@ -0,0 +1,364 @@ +//! Detection post-processing primitives. + +use std::cmp::Ordering; + +use crate::{VisionError, VisionResult}; + +/// Axis-aligned bounding box using half-open `(min, max)` coordinates. +#[derive(Clone, Copy, Debug, Default, PartialEq)] +pub struct BoundingBox2 { + /// Minimum x coordinate. + pub x_min: f32, + /// Minimum y coordinate. + pub y_min: f32, + /// Maximum x coordinate. + pub x_max: f32, + /// Maximum y coordinate. + pub y_max: f32, +} + +impl BoundingBox2 { + /// Creates a validated axis-aligned box. + pub fn try_new(x_min: f32, y_min: f32, x_max: f32, y_max: f32) -> VisionResult { + let bbox = Self { x_min, y_min, x_max, y_max }; + if !bbox.is_valid() { + return Err(VisionError::InvalidParameter( + "box coordinates must be finite and max >= min".to_owned(), + )); + } + Ok(bbox) + } + + /// Returns whether coordinates are finite and ordered. + #[must_use] + pub fn is_valid(self) -> bool { + self.x_min.is_finite() + && self.y_min.is_finite() + && self.x_max.is_finite() + && self.y_max.is_finite() + && self.x_max >= self.x_min + && self.y_max >= self.y_min + } + + /// Box width. + #[must_use] + pub fn width(self) -> f32 { + (self.x_max - self.x_min).max(0.0) + } + + /// Box height. + #[must_use] + pub fn height(self) -> f32 { + (self.y_max - self.y_min).max(0.0) + } + + /// Box area. + #[must_use] + pub fn area(self) -> f32 { + self.width() * self.height() + } + + /// Intersection box, including zero-area edge contact. + #[must_use] + pub fn intersection(self, other: Self) -> Option { + let result = Self { + x_min: self.x_min.max(other.x_min), + y_min: self.y_min.max(other.y_min), + x_max: self.x_max.min(other.x_max), + y_max: self.y_max.min(other.y_max), + }; + result.is_valid().then_some(result) + } + + /// Intersection over union, returning zero for two zero-area boxes. + #[must_use] + pub fn iou(self, other: Self) -> f32 { + let intersection = self.intersection(other).map_or(0.0, Self::area); + let union = self.area() + other.area() - intersection; + if union > 0.0 { + intersection / union + } else { + 0.0 + } + } + + /// Generalized intersection over union (GIoU). + #[must_use] + pub fn generalized_iou(self, other: Self) -> f32 { + let iou = self.iou(other); + let enclosing = Self { + x_min: self.x_min.min(other.x_min), + y_min: self.y_min.min(other.y_min), + x_max: self.x_max.max(other.x_max), + y_max: self.y_max.max(other.y_max), + }; + let enclosing_area = enclosing.area(); + if enclosing_area == 0.0 { + return iou; + } + let intersection = self.intersection(other).map_or(0.0, Self::area); + let union = self.area() + other.area() - intersection; + iou - (enclosing_area - union) / enclosing_area + } + + /// Clips coordinates into an image rectangle. + #[must_use] + pub fn clip(self, width: f32, height: f32) -> Self { + Self { + x_min: self.x_min.clamp(0.0, width), + y_min: self.y_min.clamp(0.0, height), + x_max: self.x_max.clamp(0.0, width), + y_max: self.y_max.clamp(0.0, height), + } + } + + /// Applies uniform scaling and translation, useful for letterbox mapping. + #[must_use] + pub fn scale_translate(self, scale: f32, tx: f32, ty: f32) -> Self { + Self { + x_min: self.x_min.mul_add(scale, tx), + y_min: self.y_min.mul_add(scale, ty), + x_max: self.x_max.mul_add(scale, tx), + y_max: self.y_max.mul_add(scale, ty), + } + } +} + +/// One scored detection used by class-aware post-processing. +#[derive(Clone, Copy, Debug, PartialEq)] +pub struct Detection { + /// Detected box. + pub bbox: BoundingBox2, + /// Confidence score. + pub score: f32, + /// Model-defined class identifier. + pub class_id: i64, +} + +/// Soft-NMS score-decay strategy. +#[derive(Clone, Copy, Debug, PartialEq)] +pub enum SoftNmsMethod { + /// Hard suppression above the IoU threshold. + Hard, + /// Linear score decay above the IoU threshold. + Linear, + /// Gaussian score decay at every overlap. + Gaussian { + /// Positive Gaussian variance parameter. + sigma: f32, + }, +} + +/// Index and possibly updated score returned by Soft-NMS. +#[derive(Clone, Copy, Debug, PartialEq)] +pub struct ScoredIndex { + /// Index into the original input arrays. + pub index: usize, + /// Score after overlap decay. + pub score: f32, +} + +/// Greedy non-maximum suppression. Returned indices are score-descending. +pub fn nms( + boxes: &[BoundingBox2], + scores: &[f32], + score_threshold: f32, + iou_threshold: f32, +) -> VisionResult> { + validate_nms_inputs(boxes, scores, score_threshold, iou_threshold)?; + let mut order: Vec = scores + .iter() + .enumerate() + .filter_map(|(index, &score)| (score >= score_threshold).then_some(index)) + .collect(); + sort_indices_by_score(&mut order, scores); + let mut keep: Vec = Vec::with_capacity(order.len()); + 'candidate: for index in order { + for &selected in &keep { + if boxes[index].iou(boxes[selected]) > iou_threshold { + continue 'candidate; + } + } + keep.push(index); + } + Ok(keep) +} + +/// Class-aware greedy NMS over detection records. +pub fn batched_nms( + detections: &[Detection], + score_threshold: f32, + iou_threshold: f32, +) -> VisionResult> { + if !score_threshold.is_finite() + || !iou_threshold.is_finite() + || !(0.0..=1.0).contains(&iou_threshold) + { + return Err(VisionError::InvalidParameter( + "NMS thresholds must be finite and IoU in [0, 1]".to_owned(), + )); + } + if detections.iter().any(|detection| !detection.bbox.is_valid() || !detection.score.is_finite()) + { + return Err(VisionError::InvalidParameter( + "detections must contain valid boxes and finite scores".to_owned(), + )); + } + let scores: Vec = detections.iter().map(|detection| detection.score).collect(); + let mut order: Vec = scores + .iter() + .enumerate() + .filter_map(|(index, &score)| (score >= score_threshold).then_some(index)) + .collect(); + sort_indices_by_score(&mut order, &scores); + let mut keep: Vec = Vec::with_capacity(order.len()); + 'candidate: for index in order { + for &selected in &keep { + if detections[index].class_id == detections[selected].class_id + && detections[index].bbox.iou(detections[selected].bbox) > iou_threshold + { + continue 'candidate; + } + } + keep.push(index); + } + Ok(keep) +} + +/// Soft non-maximum suppression with deterministic score ordering. +pub fn soft_nms( + boxes: &[BoundingBox2], + scores: &[f32], + score_threshold: f32, + iou_threshold: f32, + method: SoftNmsMethod, +) -> VisionResult> { + validate_nms_inputs(boxes, scores, score_threshold, iou_threshold)?; + if matches!(method, SoftNmsMethod::Gaussian { sigma } if !sigma.is_finite() || sigma <= 0.0) { + return Err(VisionError::InvalidParameter( + "Soft-NMS Gaussian sigma must be finite and positive".to_owned(), + )); + } + let mut candidates: Vec = scores + .iter() + .copied() + .enumerate() + .map(|(index, score)| ScoredIndex { index, score }) + .collect(); + let mut output = Vec::new(); + while !candidates.is_empty() { + candidates.sort_by(score_order); + let selected = candidates.remove(0); + if selected.score < score_threshold { + break; + } + output.push(selected); + for candidate in &mut candidates { + let overlap = boxes[selected.index].iou(boxes[candidate.index]); + let weight = match method { + SoftNmsMethod::Hard => { + if overlap <= iou_threshold { + 1.0 + } else { + 0.0 + } + } + SoftNmsMethod::Linear => { + if overlap > iou_threshold { + 1.0 - overlap + } else { + 1.0 + } + } + SoftNmsMethod::Gaussian { sigma } => (-(overlap * overlap) / sigma).exp(), + }; + candidate.score *= weight; + } + candidates.retain(|candidate| candidate.score >= score_threshold); + } + Ok(output) +} + +fn validate_nms_inputs( + boxes: &[BoundingBox2], + scores: &[f32], + score_threshold: f32, + iou_threshold: f32, +) -> VisionResult<()> { + if boxes.len() != scores.len() { + return Err(VisionError::ShapeMismatch( + "boxes and scores must have equal lengths".to_owned(), + )); + } + if !score_threshold.is_finite() + || !iou_threshold.is_finite() + || !(0.0..=1.0).contains(&iou_threshold) + { + return Err(VisionError::InvalidParameter( + "NMS thresholds must be finite and IoU in [0, 1]".to_owned(), + )); + } + if boxes.iter().any(|bbox| !bbox.is_valid()) || scores.iter().any(|score| !score.is_finite()) { + return Err(VisionError::InvalidParameter( + "boxes must be valid and scores finite".to_owned(), + )); + } + Ok(()) +} + +fn sort_indices_by_score(indices: &mut [usize], scores: &[f32]) { + indices.sort_by(|&left, &right| { + scores[right] + .partial_cmp(&scores[left]) + .unwrap_or(Ordering::Equal) + .then_with(|| left.cmp(&right)) + }); +} + +fn score_order(left: &ScoredIndex, right: &ScoredIndex) -> Ordering { + right + .score + .partial_cmp(&left.score) + .unwrap_or(Ordering::Equal) + .then_with(|| left.index.cmp(&right.index)) +} + +#[cfg(test)] +mod tests { + use super::{batched_nms, nms, soft_nms, BoundingBox2, Detection, SoftNmsMethod}; + + fn bbox(x0: f32, y0: f32, x1: f32, y1: f32) -> BoundingBox2 { + BoundingBox2::try_new(x0, y0, x1, y1).unwrap() + } + + #[test] + fn iou_matches_known_overlap() { + let a = bbox(0.0, 0.0, 2.0, 2.0); + let b = bbox(1.0, 1.0, 3.0, 3.0); + assert!((a.iou(b) - 1.0 / 7.0).abs() < 1e-6); + assert!(a.generalized_iou(b) < a.iou(b)); + } + + #[test] + fn nms_suppresses_lower_scored_overlap() { + let boxes = [bbox(0.0, 0.0, 2.0, 2.0), bbox(0.1, 0.1, 2.1, 2.1), bbox(5.0, 5.0, 6.0, 6.0)]; + assert_eq!(nms(&boxes, &[0.9, 0.8, 0.7], 0.0, 0.5).unwrap(), vec![0, 2]); + } + + #[test] + fn batched_nms_keeps_overlapping_different_classes() { + let detections = [ + Detection { bbox: bbox(0.0, 0.0, 2.0, 2.0), score: 0.9, class_id: 1 }, + Detection { bbox: bbox(0.0, 0.0, 2.0, 2.0), score: 0.8, class_id: 2 }, + ]; + assert_eq!(batched_nms(&detections, 0.0, 0.5).unwrap(), vec![0, 1]); + } + + #[test] + fn soft_nms_decays_overlapping_score() { + let boxes = [bbox(0.0, 0.0, 2.0, 2.0), bbox(0.1, 0.1, 2.1, 2.1)]; + let result = soft_nms(&boxes, &[0.9, 0.8], 0.01, 0.5, SoftNmsMethod::Linear).unwrap(); + assert_eq!(result[0].index, 0); + assert!(result[1].score < 0.8); + } +} diff --git a/crates/spatialrust-vision/src/error.rs b/crates/spatialrust-vision/src/error.rs new file mode 100644 index 0000000..38407a8 --- /dev/null +++ b/crates/spatialrust-vision/src/error.rs @@ -0,0 +1,24 @@ +use spatialrust_image::ImageError; + +/// Result type for image and vision operations. +pub type VisionResult = Result; + +/// Errors raised by image and vision algorithms. +#[derive(Clone, Debug, PartialEq, thiserror::Error)] +pub enum VisionError { + /// Image construction or layout validation failed. + #[error(transparent)] + Image(#[from] ImageError), + /// An operation received unusable image dimensions. + #[error("invalid image dimensions: {0}")] + InvalidDimensions(String), + /// A numeric parameter was invalid. + #[error("invalid parameter: {0}")] + InvalidParameter(String), + /// Input collections or maps had incompatible shapes. + #[error("shape mismatch: {0}")] + ShapeMismatch(String), + /// A geometric transform could not be inverted. + #[error("transform is singular")] + SingularTransform, +} diff --git a/crates/spatialrust-vision/src/lib.rs b/crates/spatialrust-vision/src/lib.rs new file mode 100644 index 0000000..0669424 --- /dev/null +++ b/crates/spatialrust-vision/src/lib.rs @@ -0,0 +1,39 @@ +//! AI-ready CPU image processing and vision algorithms. +//! +//! Algorithms are feature-gated by area. GPU implementations belong in an +//! explicit backend and must not introduce hidden device transfers. + +#![deny(unsafe_code)] +#![warn(missing_docs)] + +mod error; +mod pixel; + +#[cfg(feature = "dense")] +mod dense; +#[cfg(feature = "detection")] +mod detection; +#[cfg(feature = "preprocess")] +mod preprocess; +#[cfg(feature = "resize")] +mod resize; +#[cfg(feature = "spatial")] +mod spatial; +#[cfg(feature = "warp")] +mod warp; + +pub use error::{VisionError, VisionResult}; +pub use pixel::PixelComponent; + +#[cfg(feature = "dense")] +pub use dense::*; +#[cfg(feature = "detection")] +pub use detection::*; +#[cfg(feature = "preprocess")] +pub use preprocess::*; +#[cfg(feature = "resize")] +pub use resize::*; +#[cfg(feature = "spatial")] +pub use spatial::*; +#[cfg(feature = "warp")] +pub use warp::*; diff --git a/crates/spatialrust-vision/src/pixel.rs b/crates/spatialrust-vision/src/pixel.rs new file mode 100644 index 0000000..10d27d9 --- /dev/null +++ b/crates/spatialrust-vision/src/pixel.rs @@ -0,0 +1,47 @@ +/// Scalar component supported by generic CPU image kernels. +pub trait PixelComponent: Copy + Send + Sync + 'static { + /// Converts a scalar into the kernel accumulator representation. + fn to_f64(self) -> f64; + /// Converts a finite accumulator value back into the scalar dtype. + fn from_f64(value: f64) -> Self; +} + +impl PixelComponent for u8 { + fn to_f64(self) -> f64 { + f64::from(self) + } + + fn from_f64(value: f64) -> Self { + value.round().clamp(0.0, 255.0) as Self + } +} + +impl PixelComponent for u16 { + fn to_f64(self) -> f64 { + f64::from(self) + } + + fn from_f64(value: f64) -> Self { + value.round().clamp(0.0, f64::from(u16::MAX)) as Self + } +} + +impl PixelComponent for f32 { + fn to_f64(self) -> f64 { + f64::from(self) + } + + fn from_f64(value: f64) -> Self { + value as Self + } +} + +impl PixelComponent for f64 { + fn to_f64(self) -> f64 { + self + } + + fn from_f64(value: f64) -> Self { + value + } +} diff --git a/crates/spatialrust-vision/src/preprocess.rs b/crates/spatialrust-vision/src/preprocess.rs new file mode 100644 index 0000000..0c8530a --- /dev/null +++ b/crates/spatialrust-vision/src/preprocess.rs @@ -0,0 +1,332 @@ +use spatialrust_image::{ColorSpace, Image, ImageMetadata, ImageRegion, ImageView, PlanarImage}; + +use crate::{resize, Interpolation, PixelComponent, VisionError, VisionResult}; + +/// Padding applied around an image. +#[derive(Clone, Copy, Debug, Default, PartialEq, Eq, Hash)] +pub struct Padding { + /// Columns before the source image. + pub left: usize, + /// Columns after the source image. + pub right: usize, + /// Rows before the source image. + pub top: usize, + /// Rows after the source image. + pub bottom: usize, +} + +/// Mapping between a source image and its letterboxed output. +#[derive(Clone, Copy, Debug, PartialEq)] +pub struct LetterboxTransform { + /// Uniform source-to-output scale. + pub scale: f64, + /// Left padding in output pixels. + pub pad_left: usize, + /// Top padding in output pixels. + pub pad_top: usize, + /// Resized content width. + pub content_width: usize, + /// Resized content height. + pub content_height: usize, + /// Final output width. + pub output_width: usize, + /// Final output height. + pub output_height: usize, +} + +impl LetterboxTransform { + /// Maps a source pixel coordinate into letterboxed output coordinates. + #[must_use] + pub fn map_point(self, x: f64, y: f64) -> (f64, f64) { + (x.mul_add(self.scale, self.pad_left as f64), y.mul_add(self.scale, self.pad_top as f64)) + } + + /// Maps an output coordinate back into the source image. + #[must_use] + pub fn unmap_point(self, x: f64, y: f64) -> (f64, f64) { + ((x - self.pad_left as f64) / self.scale, (y - self.pad_top as f64) / self.scale) + } +} + +/// Copies a checked rectangular region into a packed image. +pub fn crop( + input: ImageView<'_, T, CHANNELS>, + region: ImageRegion, +) -> VisionResult> { + let view = input.subview(region)?; + let mut data = Vec::with_capacity(region.width * region.height * CHANNELS); + for y in 0..region.height { + data.extend_from_slice(view.row(y).expect("validated crop row")); + } + Ok(Image::try_new_with_metadata(region.width, region.height, data, input.metadata())?) +} + +/// Pads an image with a constant pixel value. +pub fn pad( + input: ImageView<'_, T, CHANNELS>, + padding: Padding, + value: [T; CHANNELS], +) -> VisionResult> { + let width = input + .width() + .checked_add(padding.left) + .and_then(|value| value.checked_add(padding.right)) + .ok_or_else(|| VisionError::InvalidDimensions("padding width overflow".to_owned()))?; + let height = input + .height() + .checked_add(padding.top) + .and_then(|value| value.checked_add(padding.bottom)) + .ok_or_else(|| VisionError::InvalidDimensions("padding height overflow".to_owned()))?; + let mut output = Vec::with_capacity(width * height * CHANNELS); + for y in 0..height { + for x in 0..width { + if x >= padding.left + && x < padding.left + input.width() + && y >= padding.top + && y < padding.top + input.height() + { + output.extend_from_slice( + input.get(x - padding.left, y - padding.top).expect("validated pad coordinate"), + ); + } else { + output.extend_from_slice(&value); + } + } + } + Ok(Image::try_new_with_metadata(width, height, output, input.metadata())?) +} + +/// Aspect-preserving resize followed by centered constant padding. +pub fn letterbox( + input: ImageView<'_, T, CHANNELS>, + output_width: usize, + output_height: usize, + interpolation: Interpolation, + value: [T; CHANNELS], +) -> VisionResult<(Image, LetterboxTransform)> { + if input.width() == 0 || input.height() == 0 || output_width == 0 || output_height == 0 { + return Err(VisionError::InvalidDimensions( + "letterbox requires non-empty input and output".to_owned(), + )); + } + let scale = (output_width as f64 / input.width() as f64) + .min(output_height as f64 / input.height() as f64); + let content_width = ((input.width() as f64 * scale).round() as usize).clamp(1, output_width); + let content_height = ((input.height() as f64 * scale).round() as usize).clamp(1, output_height); + let resized = resize(input, content_width, content_height, interpolation)?; + let remaining_x = output_width - content_width; + let remaining_y = output_height - content_height; + let padding = Padding { + left: remaining_x / 2, + right: remaining_x - remaining_x / 2, + top: remaining_y / 2, + bottom: remaining_y - remaining_y / 2, + }; + let transform = LetterboxTransform { + scale, + pad_left: padding.left, + pad_top: padding.top, + content_width, + content_height, + output_width, + output_height, + }; + Ok((pad(resized.view(), padding, value)?, transform)) +} + +/// Converts packed scalar values to `f32`, then applies `value * scale`, mean +/// subtraction, and per-channel standard-deviation normalization. +pub fn normalize( + input: ImageView<'_, T, CHANNELS>, + scale: f32, + mean: [f32; CHANNELS], + std: [f32; CHANNELS], +) -> VisionResult> { + if !scale.is_finite() || std.iter().any(|value| !value.is_finite() || *value == 0.0) { + return Err(VisionError::InvalidParameter( + "normalization scale/std must be finite and std non-zero".to_owned(), + )); + } + let mut output = Vec::with_capacity(input.width() * input.height() * CHANNELS); + for y in 0..input.height() { + for x in 0..input.width() { + let pixel = input.get(x, y).expect("image coordinate in bounds"); + for channel in 0..CHANNELS { + output + .push((pixel[channel].to_f64() as f32 * scale - mean[channel]) / std[channel]); + } + } + } + Ok(Image::try_new_with_metadata(input.width(), input.height(), output, input.metadata())?) +} + +/// Packs and normalizes an interleaved image into planar CHW storage. +pub fn pack_chw( + input: ImageView<'_, T, CHANNELS>, + scale: f32, + mean: [f32; CHANNELS], + std: [f32; CHANNELS], +) -> VisionResult> { + if !scale.is_finite() || std.iter().any(|value| !value.is_finite() || *value == 0.0) { + return Err(VisionError::InvalidParameter( + "normalization scale/std must be finite and std non-zero".to_owned(), + )); + } + let plane_len = input.width() * input.height(); + let mut output = vec![0.0_f32; plane_len * CHANNELS]; + for y in 0..input.height() { + for x in 0..input.width() { + let pixel = input.get(x, y).expect("image coordinate in bounds"); + let index = y * input.width() + x; + for channel in 0..CHANNELS { + output[channel * plane_len + index] = + (pixel[channel].to_f64() as f32 * scale - mean[channel]) / std[channel]; + } + } + } + Ok(PlanarImage::try_new_with_metadata(input.width(), input.height(), output, input.metadata())?) +} + +/// Swaps the red and blue channels of a three-channel image. +pub fn swap_red_blue(input: ImageView<'_, T, 3>) -> VisionResult> { + let mut data = Vec::with_capacity(input.width() * input.height() * 3); + for y in 0..input.height() { + for x in 0..input.width() { + let pixel = input.get(x, y).expect("image coordinate in bounds"); + data.extend_from_slice(&[pixel[2], pixel[1], pixel[0]]); + } + } + let mut metadata = input.metadata(); + metadata.color_space = match metadata.color_space { + ColorSpace::Rgb => ColorSpace::Bgr, + ColorSpace::Bgr => ColorSpace::Rgb, + other => other, + }; + Ok(Image::try_new_with_metadata(input.width(), input.height(), data, metadata)?) +} + +/// Converts RGB `u8` pixels to BT.601 luma using fixed-point coefficients. +pub fn rgb_to_gray(input: ImageView<'_, u8, 3>) -> VisionResult> { + let mut data = Vec::with_capacity(input.width() * input.height()); + for y in 0..input.height() { + for x in 0..input.width() { + let pixel = input.get(x, y).expect("image coordinate in bounds"); + let value = (77_u32 * u32::from(pixel[0]) + + 150_u32 * u32::from(pixel[1]) + + 29_u32 * u32::from(pixel[2]) + + 128) + >> 8; + data.push(value as u8); + } + } + let metadata = ImageMetadata { color_space: ColorSpace::Gray, ..input.metadata() }; + Ok(Image::try_new_with_metadata(input.width(), input.height(), data, metadata)?) +} + +/// Replicates a gray channel into RGB. +pub fn gray_to_rgb(input: ImageView<'_, T, 1>) -> VisionResult> { + let mut data = Vec::with_capacity(input.width() * input.height() * 3); + for y in 0..input.height() { + for x in 0..input.width() { + let value = input.get(x, y).expect("image coordinate in bounds")[0]; + data.extend_from_slice(&[value, value, value]); + } + } + let metadata = ImageMetadata { color_space: ColorSpace::Rgb, ..input.metadata() }; + Ok(Image::try_new_with_metadata(input.width(), input.height(), data, metadata)?) +} + +/// Converts RGB `u8` to OpenCV-style HSV (`H` in `0..=179`, `S/V` in `0..=255`). +pub fn rgb_to_hsv(input: ImageView<'_, u8, 3>) -> VisionResult> { + let mut data = Vec::with_capacity(input.width() * input.height() * 3); + for y in 0..input.height() { + for x in 0..input.width() { + let [r, g, b] = *input.get(x, y).expect("image coordinate in bounds"); + let rf = f64::from(r) / 255.0; + let gf = f64::from(g) / 255.0; + let bf = f64::from(b) / 255.0; + let max = rf.max(gf).max(bf); + let min = rf.min(gf).min(bf); + let delta = max - min; + let mut hue = if delta == 0.0 { + 0.0 + } else if max == rf { + 60.0 * ((gf - bf) / delta).rem_euclid(6.0) + } else if max == gf { + 60.0 * ((bf - rf) / delta + 2.0) + } else { + 60.0 * ((rf - gf) / delta + 4.0) + }; + if hue >= 360.0 { + hue = 0.0; + } + let saturation = if max == 0.0 { 0.0 } else { delta / max }; + data.extend_from_slice(&[ + (hue / 2.0).round().clamp(0.0, 179.0) as u8, + (saturation * 255.0).round() as u8, + (max * 255.0).round() as u8, + ]); + } + } + let metadata = ImageMetadata { color_space: ColorSpace::Hsv, ..input.metadata() }; + Ok(Image::try_new_with_metadata(input.width(), input.height(), data, metadata)?) +} + +#[cfg(test)] +mod tests { + use super::{ + crop, gray_to_rgb, letterbox, normalize, pack_chw, pad, rgb_to_gray, rgb_to_hsv, Padding, + }; + use crate::Interpolation; + use spatialrust_image::{ColorSpace, Image, ImageMetadata, ImageRegion}; + + #[test] + fn crop_and_pad_roundtrip_center() { + let image = Image::::try_new(3, 2, vec![1, 2, 3, 4, 5, 6]).unwrap(); + let cropped = crop(image.view(), ImageRegion::new(1, 0, 2, 2)).unwrap(); + assert_eq!(cropped.as_slice(), &[2, 3, 5, 6]); + let padded = + pad(cropped.view(), Padding { left: 1, right: 0, top: 1, bottom: 0 }, [0]).unwrap(); + assert_eq!(padded.as_slice(), &[0, 0, 0, 0, 2, 3, 0, 5, 6]); + } + + #[test] + fn letterbox_preserves_aspect_and_maps_points() { + let image = Image::::try_new(4, 2, vec![1; 8]).unwrap(); + let (output, transform) = + letterbox(image.view(), 8, 8, Interpolation::Nearest, [0]).unwrap(); + assert_eq!((output.width(), output.height()), (8, 8)); + assert_eq!((transform.content_width, transform.content_height), (8, 4)); + assert_eq!(transform.pad_top, 2); + let mapped = transform.map_point(1.0, 1.0); + assert_eq!(transform.unmap_point(mapped.0, mapped.1), (1.0, 1.0)); + } + + #[test] + fn normalize_and_chw_have_expected_layout() { + let image = Image::::try_new(2, 1, vec![0, 10, 20, 30, 40, 50]).unwrap(); + let normalized = normalize(image.view(), 0.1, [0.0; 3], [1.0; 3]).unwrap(); + assert_eq!(normalized.as_slice(), &[0.0, 1.0, 2.0, 3.0, 4.0, 5.0]); + let chw = pack_chw(image.view(), 0.1, [0.0; 3], [1.0; 3]).unwrap(); + assert_eq!(chw.as_slice(), &[0.0, 3.0, 1.0, 4.0, 2.0, 5.0]); + } + + #[test] + fn color_conversions_match_known_primaries() { + let metadata = ImageMetadata { color_space: ColorSpace::Rgb, ..Default::default() }; + let image = Image::::try_new_with_metadata( + 3, + 1, + vec![255, 0, 0, 0, 255, 0, 0, 0, 255], + metadata, + ) + .unwrap(); + assert_eq!(rgb_to_gray(image.view()).unwrap().as_slice(), &[77, 149, 29]); + assert_eq!( + rgb_to_hsv(image.view()).unwrap().as_slice(), + &[0, 255, 255, 60, 255, 255, 120, 255, 255] + ); + let gray = Image::::try_new(1, 1, vec![12]).unwrap(); + assert_eq!(gray_to_rgb(gray.view()).unwrap().as_slice(), &[12, 12, 12]); + } +} diff --git a/crates/spatialrust-vision/src/resize.rs b/crates/spatialrust-vision/src/resize.rs new file mode 100644 index 0000000..0a5fa27 --- /dev/null +++ b/crates/spatialrust-vision/src/resize.rs @@ -0,0 +1,212 @@ +use spatialrust_image::{Image, ImageView}; + +use crate::{PixelComponent, VisionError, VisionResult}; + +/// Sampling filter used by image resampling and geometric warps. +#[derive(Clone, Copy, Debug, Default, PartialEq, Eq, Hash)] +pub enum Interpolation { + /// Closest source pixel using half-pixel coordinate mapping. + Nearest, + /// Bilinear interpolation using half-pixel coordinate mapping. + #[default] + Bilinear, + /// Bicubic interpolation with OpenCV-compatible cubic coefficient `-0.75`. + Bicubic, + /// Pixel-area integration when shrinking; bilinear when enlarging. + Area, +} + +/// Resizes an interleaved image while preserving semantic metadata. +pub fn resize( + input: ImageView<'_, T, CHANNELS>, + output_width: usize, + output_height: usize, + interpolation: Interpolation, +) -> VisionResult> { + if output_width == 0 || output_height == 0 { + return Image::try_new_with_metadata( + output_width, + output_height, + Vec::new(), + input.metadata(), + ) + .map_err(Into::into); + } + if input.width() == 0 || input.height() == 0 { + return Err(VisionError::InvalidDimensions( + "cannot resize an empty input to a non-empty output".to_owned(), + )); + } + + let mut output = Vec::with_capacity(output_width * output_height * CHANNELS); + let area_downsample = interpolation == Interpolation::Area + && (output_width < input.width() || output_height < input.height()); + for y in 0..output_height { + for x in 0..output_width { + let pixel = if area_downsample { + sample_area(input, x, y, output_width, output_height) + } else { + let sx = half_pixel_coordinate(x, input.width(), output_width); + let sy = half_pixel_coordinate(y, input.height(), output_height); + match interpolation { + Interpolation::Nearest => sample_nearest(input, sx, sy), + Interpolation::Bilinear | Interpolation::Area => sample_bilinear(input, sx, sy), + Interpolation::Bicubic => sample_bicubic(input, sx, sy), + } + }; + output.extend_from_slice(&pixel); + } + } + Ok(Image::try_new_with_metadata(output_width, output_height, output, input.metadata())?) +} + +fn half_pixel_coordinate(output: usize, input_len: usize, output_len: usize) -> f64 { + (output as f64 + 0.5) * input_len as f64 / output_len as f64 - 0.5 +} + +fn clamped_index(value: isize, length: usize) -> usize { + value.clamp(0, length.saturating_sub(1) as isize) as usize +} + +pub(crate) fn sample_nearest( + input: ImageView<'_, T, CHANNELS>, + x: f64, + y: f64, +) -> [T; CHANNELS] { + let ix = clamped_index(x.round() as isize, input.width()); + let iy = clamped_index(y.round() as isize, input.height()); + *input.get(ix, iy).expect("clamped resize coordinate") +} + +pub(crate) fn sample_bilinear( + input: ImageView<'_, T, CHANNELS>, + x: f64, + y: f64, +) -> [T; CHANNELS] { + let x0_raw = x.floor() as isize; + let y0_raw = y.floor() as isize; + let x1_raw = x0_raw + 1; + let y1_raw = y0_raw + 1; + let wx = x - x.floor(); + let wy = y - y.floor(); + let x0 = clamped_index(x0_raw, input.width()); + let x1 = clamped_index(x1_raw, input.width()); + let y0 = clamped_index(y0_raw, input.height()); + let y1 = clamped_index(y1_raw, input.height()); + let p00 = input.get(x0, y0).expect("clamped resize coordinate"); + let p10 = input.get(x1, y0).expect("clamped resize coordinate"); + let p01 = input.get(x0, y1).expect("clamped resize coordinate"); + let p11 = input.get(x1, y1).expect("clamped resize coordinate"); + std::array::from_fn(|channel| { + let top = p00[channel].to_f64().mul_add(1.0 - wx, p10[channel].to_f64() * wx); + let bottom = p01[channel].to_f64().mul_add(1.0 - wx, p11[channel].to_f64() * wx); + T::from_f64(top.mul_add(1.0 - wy, bottom * wy)) + }) +} + +fn cubic_weight(distance: f64) -> f64 { + let x = distance.abs(); + const A: f64 = -0.75; + if x <= 1.0 { + (A + 2.0) * x * x * x - (A + 3.0) * x * x + 1.0 + } else if x < 2.0 { + A * x * x * x - 5.0 * A * x * x + 8.0 * A * x - 4.0 * A + } else { + 0.0 + } +} + +fn sample_bicubic( + input: ImageView<'_, T, CHANNELS>, + x: f64, + y: f64, +) -> [T; CHANNELS] { + let base_x = x.floor() as isize; + let base_y = y.floor() as isize; + let mut sums = [0.0_f64; CHANNELS]; + let mut total_weight = 0.0; + for dy in -1..=2 { + let wy = cubic_weight(y - (base_y + dy) as f64); + let iy = clamped_index(base_y + dy, input.height()); + for dx in -1..=2 { + let weight = wy * cubic_weight(x - (base_x + dx) as f64); + let ix = clamped_index(base_x + dx, input.width()); + let pixel = input.get(ix, iy).expect("clamped resize coordinate"); + for channel in 0..CHANNELS { + sums[channel] += pixel[channel].to_f64() * weight; + } + total_weight += weight; + } + } + std::array::from_fn(|channel| T::from_f64(sums[channel] / total_weight)) +} + +fn sample_area( + input: ImageView<'_, T, CHANNELS>, + output_x: usize, + output_y: usize, + output_width: usize, + output_height: usize, +) -> [T; CHANNELS] { + let scale_x = input.width() as f64 / output_width as f64; + let scale_y = input.height() as f64 / output_height as f64; + let start_x = output_x as f64 * scale_x; + let end_x = (output_x + 1) as f64 * scale_x; + let start_y = output_y as f64 * scale_y; + let end_y = (output_y + 1) as f64 * scale_y; + let mut sums = [0.0_f64; CHANNELS]; + let mut total_weight = 0.0; + for iy in start_y.floor() as usize..end_y.ceil().min(input.height() as f64) as usize { + let overlap_y = (end_y.min(iy as f64 + 1.0) - start_y.max(iy as f64)).max(0.0); + for ix in start_x.floor() as usize..end_x.ceil().min(input.width() as f64) as usize { + let overlap_x = (end_x.min(ix as f64 + 1.0) - start_x.max(ix as f64)).max(0.0); + let weight = overlap_x * overlap_y; + let pixel = input.get(ix, iy).expect("area sample within image"); + for channel in 0..CHANNELS { + sums[channel] += pixel[channel].to_f64() * weight; + } + total_weight += weight; + } + } + std::array::from_fn(|channel| T::from_f64(sums[channel] / total_weight)) +} + +#[cfg(test)] +mod tests { + use super::{resize, Interpolation}; + use spatialrust_image::Image; + + #[test] + fn nearest_repeats_pixels() { + let input = Image::::try_new(2, 1, vec![10, 20]).unwrap(); + let output = resize(input.view(), 4, 1, Interpolation::Nearest).unwrap(); + assert_eq!(output.as_slice(), &[10, 10, 20, 20]); + } + + #[test] + fn bilinear_interpolates_center() { + let input = Image::::try_new(2, 2, vec![0.0, 10.0, 20.0, 30.0]).unwrap(); + let output = resize(input.view(), 3, 3, Interpolation::Bilinear).unwrap(); + assert!((output[(1, 1)][0] - 15.0).abs() < 1e-6); + } + + #[test] + fn area_computes_block_average() { + let input = Image::::try_new(2, 2, vec![0, 10, 20, 30]).unwrap(); + let output = resize(input.view(), 1, 1, Interpolation::Area).unwrap(); + assert_eq!(output.as_slice(), &[15]); + } + + #[test] + fn identity_resize_is_exact_for_all_filters() { + let input = Image::::try_new(2, 2, (0..12).collect()).unwrap(); + for filter in [ + Interpolation::Nearest, + Interpolation::Bilinear, + Interpolation::Bicubic, + Interpolation::Area, + ] { + assert_eq!(resize(input.view(), 2, 2, filter).unwrap(), input); + } + } +} diff --git a/crates/spatialrust-vision/src/spatial.rs b/crates/spatialrust-vision/src/spatial.rs new file mode 100644 index 0000000..6d9090c --- /dev/null +++ b/crates/spatialrust-vision/src/spatial.rs @@ -0,0 +1,94 @@ +//! Bridges dense vision outputs into SpatialRust camera and point-cloud APIs. + +use spatialrust_camera::{depth_to_point_cloud, DepthConversionOptions, PinholeCamera}; +use spatialrust_core::{PointBuffer, PointBufferSet, PointCloud, SpatialMetadata, StandardSchemas}; + +use crate::{ConfidenceMap, DepthMap, PointMap, VisionError, VisionResult}; + +/// Unprojects a depth map with a calibrated camera into an XYZ point cloud. +pub fn depth_map_to_point_cloud( + depth: &DepthMap, + camera: &PinholeCamera, + options: DepthConversionOptions, +) -> VisionResult { + depth_to_point_cloud(depth.image().view(), camera, options) + .map_err(|error| VisionError::InvalidParameter(error.to_string())) +} + +/// Flattens valid point-map pixels into an XYZ cloud, optionally filtering by confidence. +pub fn point_map_to_point_cloud( + point_map: &PointMap, + confidence: Option<&ConfidenceMap>, + min_confidence: f32, +) -> VisionResult { + if !min_confidence.is_finite() || !(0.0..=1.0).contains(&min_confidence) { + return Err(VisionError::InvalidParameter( + "minimum confidence must be finite and in [0, 1]".to_owned(), + )); + } + if let Some(confidence) = confidence { + if confidence.width() != point_map.width() || confidence.height() != point_map.height() { + return Err(VisionError::ShapeMismatch( + "point map and confidence map dimensions must match".to_owned(), + )); + } + } + + let capacity = point_map.width().saturating_mul(point_map.height()); + let mut xs = Vec::with_capacity(capacity); + let mut ys = Vec::with_capacity(capacity); + let mut zs = Vec::with_capacity(capacity); + for y in 0..point_map.height() { + for x in 0..point_map.width() { + let point = point_map.image().get(x, y).expect("point-map coordinate in bounds"); + if !point.iter().all(|value| value.is_finite()) { + continue; + } + if let Some(confidence) = confidence { + let score = + confidence.image().get(x, y).expect("confidence coordinate in bounds")[0]; + if score < min_confidence { + continue; + } + } + xs.push(point[0]); + ys.push(point[1]); + zs.push(point[2]); + } + } + + let mut buffers = PointBufferSet::new(); + buffers.insert("x", PointBuffer::from_f32(xs)); + buffers.insert("y", PointBuffer::from_f32(ys)); + buffers.insert("z", PointBuffer::from_f32(zs)); + PointCloud::try_from_parts(StandardSchemas::point_xyz(), buffers, SpatialMetadata::default()) + .map_err(|error| VisionError::InvalidParameter(error.to_string())) +} + +#[cfg(test)] +mod tests { + use super::{depth_map_to_point_cloud, point_map_to_point_cloud}; + use crate::{ConfidenceMap, DepthMap, PointMap}; + use spatialrust_camera::{CameraIntrinsics, PinholeCamera}; + + #[test] + fn depth_map_uses_camera_unprojection() { + let depth = DepthMap::try_new(2, 1, vec![1.0, 2.0]).unwrap(); + let camera = + PinholeCamera::new(CameraIntrinsics::try_new(2.0, 2.0, 0.0, 0.0, 2, 1).unwrap()); + let cloud = depth_map_to_point_cloud(&depth, &camera, Default::default()).unwrap(); + assert_eq!(cloud.field("x").unwrap().as_f32().unwrap(), &[0.0, 1.0]); + assert_eq!(cloud.field("z").unwrap().as_f32().unwrap(), &[1.0, 2.0]); + } + + #[test] + fn point_map_filters_invalid_and_low_confidence_points() { + let points = + PointMap::try_new(3, 1, vec![0.0, 0.0, 1.0, 1.0, 0.0, 1.0, f32::NAN, 0.0, 1.0]) + .unwrap(); + let confidence = ConfidenceMap::try_new(3, 1, vec![0.9, 0.1, 1.0]).unwrap(); + let cloud = point_map_to_point_cloud(&points, Some(&confidence), 0.5).unwrap(); + assert_eq!(cloud.len(), 1); + assert_eq!(cloud.field("x").unwrap().as_f32().unwrap(), &[0.0]); + } +} diff --git a/crates/spatialrust-vision/src/warp.rs b/crates/spatialrust-vision/src/warp.rs new file mode 100644 index 0000000..d77cec6 --- /dev/null +++ b/crates/spatialrust-vision/src/warp.rs @@ -0,0 +1,420 @@ +//! Image remapping and geometric warp primitives. + +use spatialrust_image::{Image, ImageView}; + +use crate::{Interpolation, PixelComponent, VisionError, VisionResult}; + +/// Out-of-bounds sampling behavior for geometric image operations. +#[derive(Clone, Copy, Debug, PartialEq)] +pub enum BorderMode { + /// Returns a fixed pixel outside the source image. + Constant([T; CHANNELS]), + /// Repeats the closest edge pixel. + Replicate, + /// Reflects including the edge pixel (`fedcba|abcdefgh|hgfedc`). + Reflect, + /// Reflects without repeating the edge (`gfedcb|abcdefgh|gfedcb`). + Reflect101, + /// Periodically wraps source coordinates. + Wrap, +} + +/// A source-to-destination 2D affine transform. +#[derive(Clone, Copy, Debug, PartialEq)] +pub struct AffineTransform { + /// First two rows of a homogeneous 3x3 transform. + pub matrix: [[f64; 3]; 2], +} + +impl AffineTransform { + /// Identity affine transform. + #[must_use] + pub const fn identity() -> Self { + Self { matrix: [[1.0, 0.0, 0.0], [0.0, 1.0, 0.0]] } + } + + /// Applies the source-to-destination transform. + #[must_use] + pub fn map_point(self, x: f64, y: f64) -> (f64, f64) { + ( + self.matrix[0][0].mul_add(x, self.matrix[0][1].mul_add(y, self.matrix[0][2])), + self.matrix[1][0].mul_add(x, self.matrix[1][1].mul_add(y, self.matrix[1][2])), + ) + } + + /// Computes the inverse affine transform. + pub fn inverse(self) -> VisionResult { + let a = self.matrix[0][0]; + let b = self.matrix[0][1]; + let c = self.matrix[0][2]; + let d = self.matrix[1][0]; + let e = self.matrix[1][1]; + let f = self.matrix[1][2]; + let determinant = a.mul_add(e, -(b * d)); + if !determinant.is_finite() || determinant.abs() <= f64::EPSILON { + return Err(VisionError::SingularTransform); + } + let inv = 1.0 / determinant; + Ok(Self { + matrix: [ + [e * inv, -b * inv, (b * f - e * c) * inv], + [-d * inv, a * inv, (d * c - a * f) * inv], + ], + }) + } +} + +/// A source-to-destination projective 2D transform. +#[derive(Clone, Copy, Debug, PartialEq)] +pub struct PerspectiveTransform { + /// Homogeneous 3x3 transform matrix. + pub matrix: [[f64; 3]; 3], +} + +impl PerspectiveTransform { + /// Identity projective transform. + #[must_use] + pub const fn identity() -> Self { + Self { matrix: [[1.0, 0.0, 0.0], [0.0, 1.0, 0.0], [0.0, 0.0, 1.0]] } + } + + /// Applies the source-to-destination homogeneous transform. + #[must_use] + pub fn map_point(self, x: f64, y: f64) -> Option<(f64, f64)> { + let denominator = + self.matrix[2][0].mul_add(x, self.matrix[2][1].mul_add(y, self.matrix[2][2])); + if !denominator.is_finite() || denominator.abs() <= f64::EPSILON { + return None; + } + Some(( + self.matrix[0][0].mul_add(x, self.matrix[0][1].mul_add(y, self.matrix[0][2])) + / denominator, + self.matrix[1][0].mul_add(x, self.matrix[1][1].mul_add(y, self.matrix[1][2])) + / denominator, + )) + } + + /// Computes the inverse projective transform. + pub fn inverse(self) -> VisionResult { + let m = self.matrix; + let c00 = m[1][1].mul_add(m[2][2], -(m[1][2] * m[2][1])); + let c01 = -(m[1][0].mul_add(m[2][2], -(m[1][2] * m[2][0]))); + let c02 = m[1][0].mul_add(m[2][1], -(m[1][1] * m[2][0])); + let c10 = -(m[0][1].mul_add(m[2][2], -(m[0][2] * m[2][1]))); + let c11 = m[0][0].mul_add(m[2][2], -(m[0][2] * m[2][0])); + let c12 = -(m[0][0].mul_add(m[2][1], -(m[0][1] * m[2][0]))); + let c20 = m[0][1].mul_add(m[1][2], -(m[0][2] * m[1][1])); + let c21 = -(m[0][0].mul_add(m[1][2], -(m[0][2] * m[1][0]))); + let c22 = m[0][0].mul_add(m[1][1], -(m[0][1] * m[1][0])); + let determinant = m[0][0].mul_add(c00, m[0][1].mul_add(c01, m[0][2] * c02)); + if !determinant.is_finite() || determinant.abs() <= f64::EPSILON { + return Err(VisionError::SingularTransform); + } + let inv = 1.0 / determinant; + Ok(Self { + matrix: [ + [c00 * inv, c10 * inv, c20 * inv], + [c01 * inv, c11 * inv, c21 * inv], + [c02 * inv, c12 * inv, c22 * inv], + ], + }) + } +} + +/// Samples an image at absolute coordinates supplied by two single-channel maps. +pub fn remap( + input: ImageView<'_, T, CHANNELS>, + map_x: ImageView<'_, f32, 1>, + map_y: ImageView<'_, f32, 1>, + interpolation: Interpolation, + border: BorderMode, +) -> VisionResult> { + if map_x.width() != map_y.width() || map_x.height() != map_y.height() { + return Err(VisionError::ShapeMismatch("map_x and map_y dimensions must match".to_owned())); + } + if interpolation == Interpolation::Area { + return Err(VisionError::InvalidParameter( + "area interpolation is not defined for arbitrary remap".to_owned(), + )); + } + let mut output = Vec::with_capacity(map_x.width() * map_x.height() * CHANNELS); + for y in 0..map_x.height() { + for x in 0..map_x.width() { + let sx = f64::from(map_x.get(x, y).expect("map coordinate in bounds")[0]); + let sy = f64::from(map_y.get(x, y).expect("map coordinate in bounds")[0]); + let pixel = sample(input, sx, sy, interpolation, border); + output.extend_from_slice(&pixel); + } + } + Ok(Image::try_new_with_metadata(map_x.width(), map_x.height(), output, input.metadata())?) +} + +/// Warps an image with a source-to-destination affine transform. +pub fn warp_affine( + input: ImageView<'_, T, CHANNELS>, + transform: AffineTransform, + output_width: usize, + output_height: usize, + interpolation: Interpolation, + border: BorderMode, +) -> VisionResult> { + let inverse = transform.inverse()?; + warp_with_mapping(input, output_width, output_height, interpolation, border, |x, y| { + Some(inverse.map_point(x, y)) + }) +} + +/// Warps an image with a source-to-destination projective transform. +pub fn warp_perspective( + input: ImageView<'_, T, CHANNELS>, + transform: PerspectiveTransform, + output_width: usize, + output_height: usize, + interpolation: Interpolation, + border: BorderMode, +) -> VisionResult> { + let inverse = transform.inverse()?; + warp_with_mapping(input, output_width, output_height, interpolation, border, |x, y| { + inverse.map_point(x, y) + }) +} + +fn warp_with_mapping( + input: ImageView<'_, T, CHANNELS>, + output_width: usize, + output_height: usize, + interpolation: Interpolation, + border: BorderMode, + mut mapping: impl FnMut(f64, f64) -> Option<(f64, f64)>, +) -> VisionResult> { + if interpolation == Interpolation::Area { + return Err(VisionError::InvalidParameter( + "area interpolation is not defined for geometric warp".to_owned(), + )); + } + let mut output = Vec::with_capacity(output_width * output_height * CHANNELS); + for y in 0..output_height { + for x in 0..output_width { + let pixel = mapping(x as f64, y as f64).map_or_else( + || constant_pixel(border), + |(sx, sy)| sample(input, sx, sy, interpolation, border), + ); + output.extend_from_slice(&pixel); + } + } + Ok(Image::try_new_with_metadata(output_width, output_height, output, input.metadata())?) +} + +fn constant_pixel( + border: BorderMode, +) -> [T; CHANNELS] { + match border { + BorderMode::Constant(pixel) => pixel, + BorderMode::Replicate | BorderMode::Reflect | BorderMode::Reflect101 | BorderMode::Wrap => { + std::array::from_fn(|_| T::from_f64(0.0)) + } + } +} + +fn sample( + input: ImageView<'_, T, CHANNELS>, + x: f64, + y: f64, + interpolation: Interpolation, + border: BorderMode, +) -> [T; CHANNELS] { + if !x.is_finite() || !y.is_finite() || input.width() == 0 || input.height() == 0 { + return constant_pixel(border); + } + match interpolation { + Interpolation::Nearest => fetch(input, x.round() as isize, y.round() as isize, border), + Interpolation::Bilinear => sample_bilinear(input, x, y, border), + Interpolation::Bicubic => sample_bicubic(input, x, y, border), + Interpolation::Area => unreachable!("area rejected before sampling"), + } +} + +fn sample_bilinear( + input: ImageView<'_, T, CHANNELS>, + x: f64, + y: f64, + border: BorderMode, +) -> [T; CHANNELS] { + let x0 = x.floor() as isize; + let y0 = y.floor() as isize; + let wx = x - x.floor(); + let wy = y - y.floor(); + let p00 = fetch(input, x0, y0, border); + let p10 = fetch(input, x0 + 1, y0, border); + let p01 = fetch(input, x0, y0 + 1, border); + let p11 = fetch(input, x0 + 1, y0 + 1, border); + std::array::from_fn(|channel| { + let top = p00[channel].to_f64().mul_add(1.0 - wx, p10[channel].to_f64() * wx); + let bottom = p01[channel].to_f64().mul_add(1.0 - wx, p11[channel].to_f64() * wx); + T::from_f64(top.mul_add(1.0 - wy, bottom * wy)) + }) +} + +fn cubic_weight(distance: f64) -> f64 { + let x = distance.abs(); + const A: f64 = -0.75; + if x <= 1.0 { + (A + 2.0) * x * x * x - (A + 3.0) * x * x + 1.0 + } else if x < 2.0 { + A * x * x * x - 5.0 * A * x * x + 8.0 * A * x - 4.0 * A + } else { + 0.0 + } +} + +fn sample_bicubic( + input: ImageView<'_, T, CHANNELS>, + x: f64, + y: f64, + border: BorderMode, +) -> [T; CHANNELS] { + let base_x = x.floor() as isize; + let base_y = y.floor() as isize; + let mut sums = [0.0_f64; CHANNELS]; + let mut total_weight = 0.0; + for dy in -1..=2 { + let wy = cubic_weight(y - (base_y + dy) as f64); + for dx in -1..=2 { + let weight = wy * cubic_weight(x - (base_x + dx) as f64); + let pixel = fetch(input, base_x + dx, base_y + dy, border); + for channel in 0..CHANNELS { + sums[channel] += pixel[channel].to_f64() * weight; + } + total_weight += weight; + } + } + std::array::from_fn(|channel| T::from_f64(sums[channel] / total_weight)) +} + +fn fetch( + input: ImageView<'_, T, CHANNELS>, + x: isize, + y: isize, + border: BorderMode, +) -> [T; CHANNELS] { + if x >= 0 && y >= 0 && x < input.width() as isize && y < input.height() as isize { + return *input.get(x as usize, y as usize).expect("coordinate checked"); + } + match border { + BorderMode::Constant(pixel) => pixel, + BorderMode::Replicate => { + let ix = x.clamp(0, input.width().saturating_sub(1) as isize) as usize; + let iy = y.clamp(0, input.height().saturating_sub(1) as isize) as usize; + *input.get(ix, iy).expect("replicated coordinate") + } + BorderMode::Reflect => { + let ix = border_index(x, input.width(), false); + let iy = border_index(y, input.height(), false); + *input.get(ix, iy).expect("reflected coordinate") + } + BorderMode::Reflect101 => { + let ix = border_index(x, input.width(), true); + let iy = border_index(y, input.height(), true); + *input.get(ix, iy).expect("reflected coordinate") + } + BorderMode::Wrap => { + let ix = x.rem_euclid(input.width() as isize) as usize; + let iy = y.rem_euclid(input.height() as isize) as usize; + *input.get(ix, iy).expect("wrapped coordinate") + } + } +} + +fn border_index(mut index: isize, length: usize, reflect101: bool) -> usize { + if length <= 1 { + return 0; + } + let length = length as isize; + while index < 0 || index >= length { + index = if index < 0 { + if reflect101 { + -index + } else { + -index - 1 + } + } else if reflect101 { + 2 * length - index - 2 + } else { + 2 * length - index - 1 + }; + } + index as usize +} + +#[cfg(test)] +mod tests { + use super::{ + remap, warp_affine, warp_perspective, AffineTransform, BorderMode, PerspectiveTransform, + }; + use crate::Interpolation; + use spatialrust_image::Image; + + #[test] + fn identity_remap_is_exact() { + let input = Image::::try_new(2, 2, vec![1, 2, 3, 4]).unwrap(); + let mx = Image::::try_new(2, 2, vec![0.0, 1.0, 0.0, 1.0]).unwrap(); + let my = Image::::try_new(2, 2, vec![0.0, 0.0, 1.0, 1.0]).unwrap(); + let output = remap( + input.view(), + mx.view(), + my.view(), + Interpolation::Bilinear, + BorderMode::Constant([0]), + ) + .unwrap(); + assert_eq!(output, input); + } + + #[test] + fn affine_translation_uses_constant_border() { + let input = Image::::try_new(3, 1, vec![1, 2, 3]).unwrap(); + let transform = AffineTransform { matrix: [[1.0, 0.0, 1.0], [0.0, 1.0, 0.0]] }; + let output = warp_affine( + input.view(), + transform, + 3, + 1, + Interpolation::Nearest, + BorderMode::Constant([9]), + ) + .unwrap(); + assert_eq!(output.as_slice(), &[9, 1, 2]); + } + + #[test] + fn perspective_identity_is_exact() { + let input = Image::::try_new(2, 2, vec![1, 2, 3, 4]).unwrap(); + let output = warp_perspective( + input.view(), + PerspectiveTransform::identity(), + 2, + 2, + Interpolation::Nearest, + BorderMode::Replicate, + ) + .unwrap(); + assert_eq!(output, input); + } + + #[test] + fn border_modes_are_distinct() { + let input = Image::::try_new(3, 1, vec![10, 20, 30]).unwrap(); + let mx = Image::::try_new(1, 1, vec![-1.0]).unwrap(); + let my = Image::::try_new(1, 1, vec![0.0]).unwrap(); + let run = |border| { + remap(input.view(), mx.view(), my.view(), Interpolation::Nearest, border) + .unwrap() + .as_slice()[0] + }; + assert_eq!(run(BorderMode::Constant([5])), 5); + assert_eq!(run(BorderMode::Replicate), 10); + assert_eq!(run(BorderMode::Reflect), 10); + assert_eq!(run(BorderMode::Reflect101), 20); + assert_eq!(run(BorderMode::Wrap), 30); + } +} diff --git a/crates/spatialrust-vision/tests/properties.rs b/crates/spatialrust-vision/tests/properties.rs new file mode 100644 index 0000000..8994658 --- /dev/null +++ b/crates/spatialrust-vision/tests/properties.rs @@ -0,0 +1,74 @@ +//! Cross-module properties over generated image and geometry inputs. + +#![cfg(feature = "full")] + +use proptest::prelude::*; +use spatialrust_image::Image; +use spatialrust_vision::{ + decode_rle, encode_rle, resize, BinaryMask, BoundingBox2, Interpolation, RleOrder, +}; + +proptest! { + #[test] + fn identity_resize_preserves_u8_rgb( + width in 1usize..24, + height in 1usize..24, + seed in any::(), + ) { + let len = width * height * 3; + let data = (0..len) + .map(|index| seed.wrapping_add((index as u8).wrapping_mul(37))) + .collect::>(); + let image = Image::::try_new(width, height, data.clone()).unwrap(); + for interpolation in [ + Interpolation::Nearest, + Interpolation::Bilinear, + Interpolation::Bicubic, + Interpolation::Area, + ] { + let output = resize(image.view(), width, height, interpolation).unwrap(); + prop_assert_eq!(output.as_slice(), data.as_slice()); + } + } + + #[test] + fn mask_rle_round_trips_both_orders( + width in 1usize..32, + height in 1usize..32, + seed in any::(), + ) { + let mut state = seed; + let data = (0..width * height) + .map(|_| { + state ^= state << 13; + state ^= state >> 7; + state ^= state << 17; + (state & 1) as u8 + }) + .collect::>(); + let mask = BinaryMask::try_new(width, height, data.clone()).unwrap(); + for order in [RleOrder::RowMajor, RleOrder::CocoColumnMajor] { + let decoded = decode_rle(&encode_rle(&mask, order)).unwrap(); + prop_assert_eq!(decoded.image().as_slice(), data.as_slice()); + } + } + + #[test] + fn iou_is_symmetric_and_bounded( + ax in -100.0f32..100.0, + ay in -100.0f32..100.0, + aw in 0.0f32..100.0, + ah in 0.0f32..100.0, + bx in -100.0f32..100.0, + by in -100.0f32..100.0, + bw in 0.0f32..100.0, + bh in 0.0f32..100.0, + ) { + let a = BoundingBox2::try_new(ax, ay, ax + aw, ay + ah).unwrap(); + let b = BoundingBox2::try_new(bx, by, bx + bw, by + bh).unwrap(); + let ab = a.iou(b); + let ba = b.iou(a); + prop_assert!((ab - ba).abs() <= f32::EPSILON); + prop_assert!((0.0..=1.0).contains(&ab)); + } +} diff --git a/crates/spatialrust/Cargo.toml b/crates/spatialrust/Cargo.toml index dd2697a..563eb97 100644 --- a/crates/spatialrust/Cargo.toml +++ b/crates/spatialrust/Cargo.toml @@ -126,6 +126,22 @@ metrics-distance = ["spatialrust-metrics/metrics-distance"] transform-ops = ["spatialrust-transform/transform-ops"] voxelize-occupancy = ["spatialrust-voxelize/voxelize-occupancy"] voxelize-range-image = ["spatialrust-voxelize/voxelize-range-image"] +image = ["dep:spatialrust-image"] +camera-rgbd = ["image", "dep:spatialrust-camera"] +vision = ["image", "dep:spatialrust-vision"] +vision-resize = ["vision", "spatialrust-vision/resize"] +vision-preprocess = ["vision-resize", "spatialrust-vision/preprocess"] +vision-warp = ["vision-resize", "spatialrust-vision/warp"] +vision-detection = ["vision", "spatialrust-vision/detection"] +vision-dense = ["vision-detection", "spatialrust-vision/dense"] +vision-spatial = ["vision-dense", "camera-rgbd", "spatialrust-vision/spatial"] +vision-full = [ + "vision-preprocess", + "vision-warp", + "vision-detection", + "vision-dense", + "vision-spatial", +] [dependencies] spatialrust-core = { workspace = true } @@ -141,6 +157,9 @@ spatialrust-pipeline = { workspace = true } spatialrust-metrics = { workspace = true } spatialrust-transform = { workspace = true } spatialrust-voxelize = { workspace = true } +spatialrust-image = { workspace = true, optional = true } +spatialrust-camera = { workspace = true, optional = true } +spatialrust-vision = { workspace = true, optional = true } [dev-dependencies] diff --git a/crates/spatialrust/src/lib.rs b/crates/spatialrust/src/lib.rs index 21c0bae..c78e306 100644 --- a/crates/spatialrust/src/lib.rs +++ b/crates/spatialrust/src/lib.rs @@ -20,6 +20,13 @@ pub use spatialrust_segmentation as segmentation; pub use spatialrust_transform as transform; pub use spatialrust_voxelize as voxelize; +#[cfg(feature = "camera-rgbd")] +pub use spatialrust_camera as camera; +#[cfg(feature = "image")] +pub use spatialrust_image as image; +#[cfg(feature = "vision")] +pub use spatialrust_vision as vision; + pub use spatialrust_core::{ CpuDevice, DType, Device, DeviceKind, ExecutionPolicy, FieldSemantic, FrameId, HasIntensity, HasNormals3, HasPositions3, PointBuffer, PointCloud, PointCloudBuilder, PointField, @@ -197,6 +204,21 @@ pub use spatialrust_voxelize::{voxelize, OccupancyGrid, VoxelFill, VoxelGridConf #[cfg(feature = "voxelize-range-image")] pub use spatialrust_voxelize::{range_image, RangeImage, RangeImageConfig}; +#[cfg(feature = "image")] +pub use spatialrust_image::{ + AlphaMode, ColorRange, ColorSpace, GrayImage, Image, ImageError, ImageLayout, ImageMetadata, + ImageRegion, ImageView, ImageViewMut, PlanarImage, PlanarImageView, RgbImage, +}; + +#[cfg(feature = "camera-rgbd")] +pub use spatialrust_camera::{ + depth_to_point_cloud, rgbd_to_point_cloud, BrownConrady, CameraError, CameraIntrinsics, + DepthConversionOptions, PinholeCamera, RgbdError, +}; + +#[cfg(feature = "vision")] +pub use spatialrust_vision::*; + #[cfg(feature = "pipeline-mvp")] pub use spatialrust_pipeline::{ MvpIcpConfig, MvpPipeline, MvpPipelineConfig, MvpPipelineResult, MvpRegistrationMethod, diff --git a/crates/spatialrust/tests/rgbd_pipeline.rs b/crates/spatialrust/tests/rgbd_pipeline.rs new file mode 100644 index 0000000..e6c1248 --- /dev/null +++ b/crates/spatialrust/tests/rgbd_pipeline.rs @@ -0,0 +1,56 @@ +#![cfg(all(feature = "camera-rgbd", feature = "mvp"))] + +use spatialrust::{ + rgbd_to_point_cloud, CameraIntrinsics, DepthConversionOptions, EuclideanClusterConfig, Image, + MvpPipeline, MvpPipelineConfig, NormalEstimationConfig, PinholeCamera, RansacPlaneConfig, Vec3, + VoxelGridDownsampleConfig, +}; + +#[test] +fn rgbd_cloud_runs_through_mvp_pipeline() { + let width = 16; + let height = 16; + let mut depths = vec![1.0_f32; width * height]; + for y in 6..10 { + for x in 6..10 { + depths[y * width + x] = 0.7; + } + } + let depth = Image::::try_new(width, height, depths).unwrap(); + let color = Image::::try_new(width, height, vec![128; width * height * 3]).unwrap(); + let camera = + PinholeCamera::new(CameraIntrinsics::try_new(20.0, 20.0, 7.5, 7.5, width, height).unwrap()); + let cloud = + rgbd_to_point_cloud(depth.view(), color.view(), &camera, DepthConversionOptions::default()) + .unwrap(); + + let pipeline = MvpPipeline::new(MvpPipelineConfig { + voxel: VoxelGridDownsampleConfig::centroid(0.025), + normals: NormalEstimationConfig { + k_neighbors: 8, + min_neighbors: 3, + viewpoint: Some(Vec3::new(0.0, 0.0, 0.0)), + ..Default::default() + }, + plane: RansacPlaneConfig { + distance_threshold: 0.02, + max_iterations: 200, + min_inliers: 50, + seed: 42, + ..Default::default() + }, + cluster: EuclideanClusterConfig { + cluster_tolerance: 0.1, + min_cluster_size: 1, + max_cluster_size: usize::MAX, + ..Default::default() + }, + icp: None, + ..Default::default() + }); + let result = pipeline.run(&cloud).unwrap(); + assert_eq!(cloud.len(), width * height); + assert!(result.plane.inlier_count >= 50); + assert!(!result.output.is_empty()); + assert!(result.output.schema().find_semantic(spatialrust::FieldSemantic::ColorR).is_some()); +} diff --git a/crates/spatialrust/tests/vision_ai_pipeline.rs b/crates/spatialrust/tests/vision_ai_pipeline.rs new file mode 100644 index 0000000..46be684 --- /dev/null +++ b/crates/spatialrust/tests/vision_ai_pipeline.rs @@ -0,0 +1,58 @@ +#![cfg(all(feature = "vision-full", feature = "mvp"))] + +use spatialrust::{ + point_map_to_point_cloud, ConfidenceMap, EuclideanClusterConfig, MvpPipeline, + MvpPipelineConfig, NormalEstimationConfig, PointMap, RansacPlaneConfig, Vec3, + VoxelGridDownsampleConfig, +}; + +#[test] +fn dense_ai_point_map_runs_through_spatial_pipeline() { + let width = 16; + let height = 16; + let mut points = Vec::with_capacity(width * height * 3); + let mut confidence = vec![1.0_f32; width * height]; + for y in 0..height { + for x in 0..width { + let z = if (6..10).contains(&x) && (6..10).contains(&y) { 0.7 } else { 1.0 }; + points.extend_from_slice(&[ + (x as f32 - 7.5) * z / 20.0, + (y as f32 - 7.5) * z / 20.0, + z, + ]); + } + } + confidence[0] = 0.1; + let point_map = PointMap::try_new(width, height, points).unwrap(); + let confidence = ConfidenceMap::try_new(width, height, confidence).unwrap(); + let cloud = point_map_to_point_cloud(&point_map, Some(&confidence), 0.5).unwrap(); + + let pipeline = MvpPipeline::new(MvpPipelineConfig { + voxel: VoxelGridDownsampleConfig::centroid(0.025), + normals: NormalEstimationConfig { + k_neighbors: 8, + min_neighbors: 3, + viewpoint: Some(Vec3::new(0.0, 0.0, 0.0)), + ..Default::default() + }, + plane: RansacPlaneConfig { + distance_threshold: 0.02, + max_iterations: 200, + min_inliers: 50, + seed: 75, + ..Default::default() + }, + cluster: EuclideanClusterConfig { + cluster_tolerance: 0.1, + min_cluster_size: 1, + max_cluster_size: usize::MAX, + ..Default::default() + }, + icp: None, + ..Default::default() + }); + let result = pipeline.run(&cloud).unwrap(); + assert_eq!(cloud.len(), width * height - 1); + assert!(result.plane.inlier_count >= 50); + assert!(!result.output.is_empty()); +} diff --git a/docs/API_STABILITY.md b/docs/API_STABILITY.md index 343fdf0..ad5666e 100644 --- a/docs/API_STABILITY.md +++ b/docs/API_STABILITY.md @@ -59,6 +59,8 @@ until their individual 1.0 milestones. | --- | --- | | MVP CLI flags | `--bounds`, `--resolution`, `--repeat` may gain aliases | | HTTP COPC (`mvp-http`) | URL IO is stable; timeout/retry policy may change | +| Image/camera (`image`, `camera-rgbd`) | Typed image, calibration, distortion, and RGB-D APIs are provisional | +| Vision (`vision-*`) | CPU preprocessing, warp, detection, masks, and dense spatial bridges are provisional | ## Algorithm crates @@ -71,13 +73,16 @@ spatialrust- / feature- | Crate | 1.0 status | Notes | | --- | --- | --- | | `spatialrust-math` | Stable primitives | `Vec3`, `Mat4`, `Isometry3` | +| `spatialrust-image` | Provisional | Packed ownership and strided CPU views; no hidden device transfers | +| `spatialrust-camera` | Provisional | Pinhole/Brown–Conrady and RGB-D conversion | +| `spatialrust-vision` | Provisional | Feature-gated CPU image algorithms and explicit point-cloud bridges | | `spatialrust-search` | Stable with features | KD-tree behind `search-kdtree`; **chunked query traits** and **`search-parallel`** provisional | | `spatialrust-filtering` | Provisional | GPU thresholds may move | | `spatialrust-features` | Provisional | Normal GPU path still tuning | | `spatialrust-segmentation` | Provisional | RANSAC configs may extend; **GPU plane scoring** behind `segment-ransac-plane-gpu` | | `spatialrust-registration` | Provisional | New backends (TEASER++, etc.) expected | | `spatialrust-gpu` | Provisional | `WgpuRuntime`, `GpuBufferPool` upload/recycle API stable; kernel APIs still tuning | -| `spatialrust-py` | Stable user surface | Stubs enforced by `mypy.stubtest` in CI | +| `spatialrust-py` | Stable user surface | Stubs enforced by `mypy.stubtest`; new vision functions remain provisional with the Rust APIs | ## Explicitly out of 1.0 scope diff --git a/docs/ARCHITECTURE.md b/docs/ARCHITECTURE.md index fc2a1f6..73359d9 100644 --- a/docs/ARCHITECTURE.md +++ b/docs/ARCHITECTURE.md @@ -22,6 +22,9 @@ North star: **Rust-native spatial computing** - `spatialrust-math` — Vec/Mat/Pose math - `spatialrust-io` — readers/writers (feature-gated formats) - `spatialrust-gpu` — device buffers and GPU runtime +- `spatialrust-image` — typed CPU image buffers and strided zero-copy views +- `spatialrust-camera` — camera models, distortion, and RGB-D/point-cloud bridge +- `spatialrust-vision` — feature-gated CPU preprocessing, warps, detection, masks, and dense maps ## MVP scope @@ -47,6 +50,20 @@ math -> core -> search/geometry/io/gpu -> algorithms -> integration Forbidden: `core -> io`, `core -> gpu impl`, `core -> ros2`, `core -> ai`. +Image and camera dependency direction: + +``` +math -> image +math -> image -> vision +math + image + core -> camera -> vision::spatial/rgbd/odometry +``` + +`spatialrust-image` remains independent of `spatialrust-core`. GPU image storage +must use a dedicated backend and explicit upload/readback APIs. +`spatialrust-vision` keeps preprocessing, warp, detection, dense-map, and spatial +bridges in separate additive features. CPU APIs never perform implicit device +copies; future GPU/CUDA implementations belong behind explicit backend features. + ## Roadmap epics | Year | Focus | diff --git a/notes/2026-07-14_ai_vision_foundation.md b/notes/2026-07-14_ai_vision_foundation.md new file mode 100644 index 0000000..f69301e --- /dev/null +++ b/notes/2026-07-14_ai_vision_foundation.md @@ -0,0 +1,34 @@ +# AI-ready image and vision foundation (Epics 75–79) + +Date: 2026-07-14 + +## Delivered + +- `C:\Users\rsasa\Workspace\SpatialRust\crates\spatialrust-image` now provides + checked mutable strided views and ROI/subviews, planar and interleaved layouts, + and explicit color/range/alpha metadata. +- `C:\Users\rsasa\Workspace\SpatialRust\crates\spatialrust-vision` contains + independently gated resize, preprocess, warp, detection, dense-map, and + spatial bridge modules. CPU ownership and camera/point-cloud conversion are + explicit; there are no hidden device copies. +- `C:\Users\rsasa\Workspace\SpatialRust\crates\spatialrust-py\src\lib.rs` and + `C:\Users\rsasa\Workspace\SpatialRust\crates\spatialrust-py\spatialrust.pyi` + expose the AI preprocessing/post-processing and PointMap bridge to NumPy. +- `C:\Users\rsasa\Workspace\SpatialRust\crates\spatialrust-py\examples\vision_ai_pipeline.py` + runs letterbox/CHW, NMS, mask components/RLE, PointMap conversion, and the MVP + point-cloud pipeline in one executable example. + +## Verification + +- Rust vision unit tests: 22 passed. +- Generated property tests: 3 passed (identity resize, both RLE orders, IoU laws). +- Python binding tests: 48 passed; `mypy.stubtest` passed. +- All seven `spatialrust-vision` features and all eight meta-crate `vision-*` + features built separately with default features disabled. +- `C:\Users\rsasa\Workspace\SpatialRust\bench\opencv_vision_comparison\run.py` + passed against OpenCV 4.13.0: resize max uint8 error 0/1/1/0 for + nearest/bilinear/bicubic/area, gray/HSV max error 1, remap error 0, and exact + NMS indices/component areas. +- `C:\Users\rsasa\Workspace\SpatialRust\crates\spatialrust-vision\benches\preprocess.rs` + measured 1280x720 RGB letterbox plus normalization/CHW packing to 640x640 at + 12.075–13.621 ms (Criterion interval on this development machine). diff --git a/notes/2026-07-14_opencv_rgbd_foundation.md b/notes/2026-07-14_opencv_rgbd_foundation.md new file mode 100644 index 0000000..ecea426 --- /dev/null +++ b/notes/2026-07-14_opencv_rgbd_foundation.md @@ -0,0 +1,39 @@ +# OpenCV-oriented RGB-D foundation + +Date: 2026-07-14 + +## Delivered + +- `spatialrust-image`: typed packed images and validated strided views +- `spatialrust-camera`: pinhole projection/unprojection, Brown–Conrady + distortion, depth to XYZ, and aligned RGB-D to XYZRGB +- `camera-rgbd` meta-crate feature and full MVP integration test +- NumPy/Python `rgbd_to_point_cloud` API, type stub, contract test, and demo +- OpenCV `cv2.rgbd.depthTo3d` numerical/timing comparison harness +- Criterion CPU benchmark for a 640x480 RGB-D frame + +## Local validation + +Windows, release build, Python 3.12.10: + +| Check | Result | +| --- | --- | +| OpenCV comparison points | 76,704 | +| Maximum XYZ difference | `5.960e-08 m` | +| SpatialRust Python median | `2.065 ms` | +| OpenCV Python median | `0.264 ms` | +| Native SpatialRust 640x480 RGB-D | `6.47–6.86 ms` | + +The native benchmark processes 307,200 colored points. The OpenCV comparison +uses a 320x240 depth frame and includes Python API boundary costs for both +libraries. Results are machine-specific and should be rerun before publishing. + +## Reproduction + +```powershell +cargo test -p spatialrust-camera +cargo test -p spatialrust --features mvp,camera-rgbd --test rgbd_pipeline +cargo bench -p spatialrust-camera --bench rgbd +python bench/opencv_rgbd_comparison/run.py +python crates/spatialrust-py/examples/rgbd_pipeline.py +```