Epic 148A/148B: point-cloud conformance harness (PCL/PDAL)

CHANGELOG.md

@@ -21,6 +21,15 @@ removed no sooner than the next major (see `docs/API_STABILITY.md`).
 
 ### Added
 
+- **Epic 148A/148B point-cloud conformance harness**: a versioned point-cloud
+  benchmark manifest (`spatialrust.pointcloud-benchmark-manifest.v1`) and a
+  stdlib-only report contract (`spatialrust.pointcloud-comparison.v1`) with
+  seven contract tests, plus a PDAL comparison runner
+  (`bench/pdal_comparison/`). `bench_ops` now also times `translate_xyz`
+  behind `transform-ops` to match PDAL's `filters.transformation`. A dated,
+  honest SpatialRust-vs-PCL receipt (`bench/pcl_comparison/receipt-2026-08-07.json`)
+  records voxel 1.70×, normals 4.56×, SOR 4.91×, and radius-outlier 6.62× on
+  the identical public 460,400-point cloud.
 - **Epic 147F Python 3D Tiles bindings**: `export_tiles3d` converts a
   `PointCloud` into a tileset directory and `export_copc_tiles3d` converts a
   `.copc.laz` file into a bounded tileset, both returning a

bench/pcl_comparison/README.md

@@ -37,16 +37,17 @@ that fallback also needs NumPy and the SpatialRust Python extension.
 
 ## Indicative results
 
-Measured on one Windows machine (PCL 1.15.1 via vcpkg, MSYS2 g++ 16.1.0,
-release builds, 460,400-point public PCL `table_scene_lms400.pcd`). Throughput
-varies by CPU; run it yourself for numbers on your hardware.
+Measured on one Linux machine (PCL 1.14.x via libpcl-dev, g++ 16, release
+Rust build, 460,400-point public PCL `table_scene_lms400.pcd`, single-threaded
+where that is the default). Throughput varies by CPU and PCL build; run it
+yourself for numbers on your hardware. Dated receipt: `receipt-2026-08-07.json`.
 
 | Operation | SpatialRust | PCL | Speedup |
-| --- | ---: | ---: | ---: |
-| Radius Outlier Removal | 0.0899 s | 1.8784 s | **20.89× faster** |
-| Statistical Outlier Removal | 0.1664 s | 2.0933 s | **12.58× faster** |
-| Normal estimation | 0.1461 s | 1.9750 s | **13.52× faster** |
-| Voxel downsample | 0.0104 s | 0.0181 s | **1.74× faster** |
+| --- | ---: | ---: | :--- |
+| Voxel downsample | 0.0104 s | 0.0177 s | **1.70× faster** |
+| Normal estimation | 0.2171 s | 0.9893 s | **4.56× faster** |
+| Statistical Outlier Removal | 0.2297 s | 1.1272 s | **4.91× faster** |
+| Radius Outlier Removal | 0.1088 s | 0.7200 s | **6.62× faster** |
 
 SpatialRust is faster on neighborhood-statistics and density operations (radius
 outlier removal uses an early-exit density test; normals and SOR win too; voxel

bench/pcl_comparison/manifest.json

@@ -0,0 +1,37 @@
+{
+  "schema_version": "spatialrust.pointcloud-benchmark-manifest.v1",
+  "profiles": {
+    "small": { "points": 460400, "label": "public PCL table_scene_lms400" },
+    "full": { "points": 2000000, "label": "synthetic room scene" }
+  },
+  "required_statistics": [
+    "mean",
+    "median",
+    "p95",
+    "min",
+    "max",
+    "stdev",
+    "coefficient_of_variation",
+    "median_absolute_deviation",
+    "batch_size",
+    "samples"
+  ],
+  "publication_receipt": [
+    "platform",
+    "machine",
+    "logical_cpu_count",
+    "python_version",
+    "pcl_version",
+    "pdal_version",
+    "open3d_version",
+    "spatialrust_version"
+  ],
+  "workloads": [
+    { "id": "voxel_downsample", "domain": "filtering", "libraries": ["pcl", "pdal", "open3d"] },
+    { "id": "normal_estimation", "domain": "features", "libraries": ["pcl", "open3d"] },
+    { "id": "statistical_outlier_removal", "domain": "filtering", "libraries": ["pcl", "open3d"] },
+    { "id": "radius_outlier_removal", "domain": "filtering", "libraries": ["pcl", "open3d"] },
+    { "id": "translate_xyz", "domain": "io-transform", "libraries": ["pdal"] },
+    { "id": "reprojection", "domain": "io-transform", "libraries": ["pdal"] }
+  ]
+}

bench/pcl_comparison/receipt-2026-08-07.json

@@ -0,0 +1,61 @@
+{
+  "environment": {
+    "logical_cpu_count": 8,
+    "machine": "x86_64",
+    "open3d_version": null,
+    "pcl_version": "1.14.0 (libpcl-dev)",
+    "pdal_version": null,
+    "platform": "Linux-7.0.0-28-generic-x86_64-with-glibc2.39",
+    "processor": "x86_64",
+    "python_implementation": "CPython",
+    "python_version": "3.12.3",
+    "spatialrust_version": "1.2.0"
+  },
+  "generated_at": "2026-08-07T08:11:24.924875+00:00",
+  "kind": "performance",
+  "results": {
+    "input": {
+      "name": "table_scene_lms400.pcd",
+      "points": 460400,
+      "source": "https://github.com/PointCloudLibrary/data/blob/master/tutorials/table_scene_lms400.pcd"
+    },
+    "note": "honest single-run numbers; rerun on target hardware before portability claims",
+    "operations": [
+      {
+        "id": "voxel_downsample",
+        "output_points_pcl": 2215,
+        "output_points_sr": 2194,
+        "pcl_seconds": 0.0177,
+        "spatialrust_seconds": 0.0104,
+        "speedup": 1.7
+      },
+      {
+        "id": "normal_estimation",
+        "output_points_pcl": 460400,
+        "output_points_sr": 460400,
+        "pcl_seconds": 0.9893,
+        "spatialrust_seconds": 0.2171,
+        "speedup": 4.56
+      },
+      {
+        "id": "statistical_outlier_removal",
+        "output_points_pcl": 451983,
+        "output_points_sr": 451983,
+        "pcl_seconds": 1.1272,
+        "spatialrust_seconds": 0.2297,
+        "speedup": 4.91
+      },
+      {
+        "id": "radius_outlier_removal",
+        "output_points_pcl": 460379,
+        "output_points_sr": 460379,
+        "pcl_seconds": 0.72,
+        "spatialrust_seconds": 0.1088,
+        "speedup": 6.62
+      }
+    ]
+  },
+  "schema_version": "spatialrust.pointcloud-comparison.v1",
+  "status": "pass",
+  "suite": "pcl_comparison"
+}

bench/pcl_comparison/report.py

@@ -0,0 +1,209 @@
+"""Machine-readable report contract for point-cloud comparison harnesses.
+
+This module deliberately depends only on the Python standard library so its
+contract tests can run without building SpatialRust or installing PCL / PDAL /
+Open3D. It mirrors the OpenCV comparison contract but names point-cloud
+libraries and workloads, and records output point counts for correctness
+equivalence.
+"""
+
+from __future__ import annotations
+
+import gc
+import json
+import math
+import os
+import platform
+import statistics
+import time
+from datetime import datetime, timezone
+from pathlib import Path
+from typing import Callable, Iterable, TypeVar
+
+
+SCHEMA_VERSION = "spatialrust.pointcloud-comparison.v1"
+T = TypeVar("T")
+
+
+def environment(
+    *,
+    pcl_version: str | None,
+    pdal_version: str | None,
+    open3d_version: str | None,
+    spatialrust_version: str,
+) -> dict[str, object]:
+    """Return the minimum environment receipt required for publication."""
+
+    return {
+        "python_version": platform.python_version(),
+        "python_implementation": platform.python_implementation(),
+        "platform": platform.platform(),
+        "machine": platform.machine(),
+        "processor": platform.processor() or None,
+        "logical_cpu_count": os.cpu_count(),
+        "pcl_version": pcl_version,
+        "pdal_version": pdal_version,
+        "open3d_version": open3d_version,
+        "spatialrust_version": spatialrust_version,
+    }
+
+
+def percentile(values: Iterable[float], quantile: float) -> float:
+    samples = sorted(values)
+    if not samples:
+        raise ValueError("at least one timing sample is required")
+    if not 0.0 <= quantile <= 1.0:
+        raise ValueError("quantile must be in [0, 1]")
+    index = (len(samples) - 1) * quantile
+    lower = int(index)
+    upper = min(lower + 1, len(samples) - 1)
+    fraction = index - lower
+    return samples[lower] * (1.0 - fraction) + samples[upper] * fraction
+
+
+def timed(call: Callable[[], T], *, warmup: int, repeats: int) -> tuple[T, dict[str, object]]:
+    """Benchmark a callable and return its last result plus timing statistics."""
+
+    if warmup < 0 or repeats < 1:
+        raise ValueError("warmup must be non-negative and repeats must be positive")
+    for _ in range(warmup):
+        call()
+    samples_ms: list[float] = []
+    result: T | None = None
+    gc_enabled = gc.isenabled()
+    gc.disable()
+    try:
+        for _ in range(repeats):
+            start = time.perf_counter_ns()
+            result = call()
+            samples_ms.append((time.perf_counter_ns() - start) / 1_000_000.0)
+    finally:
+        if gc_enabled:
+            gc.enable()
+    return result, timing_statistics(samples_ms, warmup=warmup)
+
+
+def timing_statistics(
+    samples_ms: list[float], *, warmup: int, batch_size: int = 1
+) -> dict[str, object]:
+    """Summarize raw millisecond samples with robust dispersion statistics."""
+
+    if not samples_ms:
+        raise ValueError("at least one timing sample is required")
+    median = statistics.median(samples_ms)
+    deviations = [abs(value - median) for value in samples_ms]
+    mean = statistics.fmean(samples_ms)
+    stdev = statistics.stdev(samples_ms) if len(samples_ms) > 1 else 0.0
+    return {
+        "unit": "ms",
+        "warmup": warmup,
+        "repeats": len(samples_ms),
+        "batch_size": batch_size,
+        "mean": mean,
+        "median": median,
+        "p95": percentile(samples_ms, 0.95),
+        "min": min(samples_ms),
+        "max": max(samples_ms),
+        "stdev": stdev,
+        "coefficient_of_variation": stdev / mean if mean > 0.0 else 0.0,
+        "median_absolute_deviation": statistics.median(deviations),
+        "samples": samples_ms,
+    }
+
+
+def make_report(
+    *,
+    suite: str,
+    kind: str,
+    status: str,
+    environment_receipt: dict[str, object],
+    results: object,
+) -> dict[str, object]:
+    if kind not in {"correctness", "performance", "aggregate"}:
+        raise ValueError(f"unsupported report kind: {kind}")
+    if status not in {"pass", "fail"}:
+        raise ValueError(f"unsupported report status: {status}")
+    return {
+        "schema_version": SCHEMA_VERSION,
+        "suite": suite,
+        "kind": kind,
+        "status": status,
+        "generated_at": datetime.now(timezone.utc).isoformat(),
+        "environment": environment_receipt,
+        "results": results,
+    }
+
+
+def emit_report(report: dict[str, object], output: Path | None = None) -> None:
+    errors = validate_report(report)
+    if errors:
+        raise ValueError(f"invalid report: {'; '.join(errors)}")
+    encoded = json.dumps(report, indent=2, sort_keys=True, allow_nan=False)
+    if output is not None:
+        output.parent.mkdir(parents=True, exist_ok=True)
+        output.write_text(encoded + "\n", encoding="utf-8")
+    print(encoded)
+
+
+def validate_report(report: object) -> list[str]:
+    """Validate the stable structural contract without third-party packages."""
+
+    errors: list[str] = []
+    if not isinstance(report, dict):
+        return ["report must be a JSON object"]
+
+    def find_non_finite(value: object, path: str) -> None:
+        if isinstance(value, float) and not math.isfinite(value):
+            errors.append(f"{path} must be finite")
+        elif isinstance(value, dict):
+            for key, child in value.items():
+                find_non_finite(child, f"{path}.{key}")
+        elif isinstance(value, list):
+            for index, child in enumerate(value):
+                find_non_finite(child, f"{path}[{index}]")
+
+    find_non_finite(report, "report")
+    required = {
+        "schema_version",
+        "suite",
+        "kind",
+        "status",
+        "generated_at",
+        "environment",
+        "results",
+    }
+    missing = sorted(required - report.keys())
+    if missing:
+        errors.append(f"missing keys: {', '.join(missing)}")
+    if report.get("schema_version") != SCHEMA_VERSION:
+        errors.append(f"schema_version must be {SCHEMA_VERSION}")
+    if report.get("kind") not in {"correctness", "performance", "aggregate"}:
+        errors.append("kind must be correctness, performance, or aggregate")
+    if report.get("status") not in {"pass", "fail"}:
+        errors.append("status must be pass or fail")
+    environment_value = report.get("environment")
+    if not isinstance(environment_value, dict):
+        errors.append("environment must be an object")
+    else:
+        for key in (
+            "platform",
+            "machine",
+            "logical_cpu_count",
+            "python_version",
+            "pcl_version",
+            "pdal_version",
+            "open3d_version",
+            "spatialrust_version",
+        ):
+            if key not in environment_value:
+                errors.append(f"environment missing {key}")
+    return errors
+
+
+def load_report(path: Path) -> dict[str, object]:
+    with path.open("r", encoding="utf-8") as handle:
+        value = json.load(handle)
+    errors = validate_report(value)
+    if errors:
+        raise ValueError(f"invalid report {path}: {'; '.join(errors)}")
+    return value