[Feature] New Matcher Class: Voronoi

.github/workflows/tests.yml

@@ -10,7 +10,7 @@ on:
     branches: ["main"]
 
 jobs:
-  build:
+  test:
     runs-on: ${{ matrix.os }}
     strategy:
       fail-fast: false
@@ -32,7 +32,7 @@ jobs:
           python -m pip install poetry
       - name: Install dependencies
         run: |
-          python -m poetry install
+          python -m poetry install --extras test
       - name: Test with pytest and create coverage report
         run: |
           python -m poetry run coverage run --source=panoptica -m pytest
@@ -43,3 +43,35 @@ jobs:
         uses: codecov/codecov-action@v4
         with:
           token: ${{ secrets.CODECOV_TOKEN }}
+
+  test-cuda:
+    runs-on: ubuntu-latest
+    strategy:
+      fail-fast: false
+      matrix:
+        python-version: ["3.10", "3.11", "3.12"]
+
+    steps:
+      - uses: actions/checkout@v4
+      - name: Set up Python ${{ matrix.python-version }}
+        uses: actions/setup-python@v5
+        with:
+          python-version: ${{ matrix.python-version }}
+          cache: "pip"
+
+      - name: Configure poetry
+        run: |
+          python -m pip install --upgrade pip
+          python -m pip install poetry
+      - name: Install dependencies with GPU extras
+        run: |
+          python -m poetry install --extras "gpu test" || python -m poetry install --extras test
+      - name: Test CUDA functionality (CPU fallback)
+        run: |
+          python -m poetry run pytest unit_tests/test_cupy_connected_components.py -v
+      - name: Upload coverage results to Codecov (Only on merge to main)
+        # Only upload to Codecov after a merge to the main branch
+        if: github.ref == 'refs/heads/main' && github.event_name == 'push'
+        uses: codecov/codecov-action@v4
+        with:
+          token: ${{ secrets.CODECOV_TOKEN }}

benchmark/benchmark.py

@@ -7,6 +7,16 @@
 # scipy needs to be installed to run this benchmark, we use cc3d as it is quicker for 3D data
 from scipy import ndimage
 
+# Try to import cupy for GPU acceleration
+try:
+    import cupy as cp
+    from cupyx.scipy import ndimage as cp_ndimage
+
+    CUPY_AVAILABLE = True
+except ImportError:
+    CUPY_AVAILABLE = False
+    print("CuPy not available. GPU benchmarks will be skipped.")
+
 
 def generate_random_binary_mask(size: Tuple[int, int, Union[int, None]]) -> np.ndarray:
     """
@@ -64,6 +74,64 @@ def label_cc3d():
     return cc3d_time
 
 
+def benchmark_cupy(mask: np.ndarray):
+    """
+    Benchmark the performance of cupy.ndimage.label for connected component labeling on GPU.
+
+    Args:
+        mask (np.ndarray): Binary mask to label.
+
+    Returns:
+        float: Time taken to label the mask in seconds, or None if CuPy is not available.
+    """
+    if not CUPY_AVAILABLE:
+        return None
+
+    # Transfer data to GPU
+    mask_gpu = cp.asarray(mask)
+
+    # Warmup phase
+    for _ in range(3):
+        cp_ndimage.label(mask_gpu)
+        cp.cuda.Stream.null.synchronize()
+
+    def label_cupy():
+        cp_ndimage.label(mask_gpu)
+        cp.cuda.Stream.null.synchronize()  # Ensure GPU computation is complete
+
+    cupy_time = timeit.timeit(label_cupy, number=10)
+
+    # Clean up GPU memory
+    del mask_gpu
+    cp.get_default_memory_pool().free_all_blocks()
+
+    return cupy_time
+
+
+def benchmark_panoptica_cupy(mask: np.ndarray):
+    """
+    Benchmark the performance of panoptica's CuPy backend for connected component labeling.
+
+    Args:
+        mask (np.ndarray): Binary mask to label.
+
+    Returns:
+        float: Time taken to label the mask in seconds, or None if CuPy is not available.
+    """
+    if not CUPY_AVAILABLE:
+        return None
+
+    from panoptica._functionals import _connected_components
+    from panoptica.utils.constants import CCABackend
+
+    def label_panoptica_cupy():
+        _connected_components(mask, CCABackend.cupy)
+
+    panoptica_cupy_time = timeit.timeit(label_panoptica_cupy, number=10)
+
+    return panoptica_cupy_time
+
+
 def run_benchmarks(volume_sizes: Tuple[Tuple[int, int, Union[int, None]]]) -> None:
     """
     Run benchmark tests for connected component labeling with different volume sizes.
@@ -80,10 +148,15 @@ def run_benchmarks(volume_sizes: Tuple[Tuple[int, int, Union[int, None]]]) -> No
 
         scipy_time = benchmark_scipy(mask)
         cc3d_time = benchmark_cc3d(mask)
+        cupy_time = benchmark_cupy(mask)
 
         print(f"Volume Size: {size}")
         print(f"Scipy Time: {scipy_time:.4f} seconds")
         print(f"CC3D Time: {cc3d_time:.4f} seconds")
+        if cupy_time is not None:
+            print(f"CuPy Time: {cupy_time:.4f} seconds")
+        else:
+            print("CuPy Time: Not available")
         print()
 
 

panoptica/__init__.py

@@ -2,7 +2,11 @@
     ConnectedComponentsInstanceApproximator,
     CCABackend,
 )
-from panoptica.instance_matcher import NaiveThresholdMatching, MaxBipartiteMatching
+from panoptica.instance_matcher import (
+    NaiveThresholdMatching,
+    MaxBipartiteMatching,
+    RegionBasedMatching,
+)
 from panoptica.panoptica_statistics import Panoptica_Statistic, ValueSummary
 from panoptica.panoptica_aggregator import Panoptica_Aggregator
 from panoptica.panoptica_evaluator import Panoptica_Evaluator

panoptica/_functionals.py

@@ -63,6 +63,19 @@ def _connected_components(
         from scipy.ndimage import label
 
         cc_arr, n_instances = label(array)
+    elif cca_backend == CCABackend.cupy:
+        try:
+            import cupy as cp
+            from cupyx.scipy.ndimage import label as cp_label
+
+            array_gpu = cp.asarray(array)
+            cc_arr, n_instances = cp_label(array_gpu)
+            cc_arr = cp.asnumpy(cc_arr)
+        except ImportError:
+            raise ImportError(
+                "CuPy is not installed. Please install CuPy to use the GPU backend. "
+                "You can install it using: pip install cupy-cuda11x or cupy-cuda12x depending on your CUDA version."
+            )
     else:
         raise NotImplementedError(cca_backend)
 

panoptica/instance_matcher.py

@@ -3,11 +3,13 @@
 from typing import Optional, Tuple, List
 
 import numpy as np
+from scipy.ndimage import distance_transform_edt
 
 from panoptica._functionals import (
     _calc_matching_metric_of_overlapping_labels,
     _calc_matching_metric_of_overlapping_partlabels,
     _map_labels,
+    _connected_components,
 )
 from panoptica.metrics import Metric
 from panoptica.utils.processing_pair import (
@@ -17,6 +19,7 @@
 from panoptica.utils.instancelabelmap import InstanceLabelMap
 from panoptica.utils.config import SupportsConfig
 from panoptica.utils.label_group import LabelGroup, LabelPartGroup
+from panoptica.utils.constants import CCABackend
 
 
 @dataclass
@@ -493,3 +496,118 @@ def _yaml_repr(cls, node) -> dict:
             "matching_metric": node._matching_metric,
             "matching_threshold": node._matching_threshold,
         }
+
+
+class RegionBasedMatching(InstanceMatchingAlgorithm):
+    """
+    Instance matching algorithm that performs region-based matching using spatial distance.
+
+    This method assigns prediction instances to ground truth regions based on spatial proximity
+    rather than traditional overlap-based metrics. It uses connected components and distance
+    transforms to create region assignments.
+
+    Note: This matching method does not produce traditional count metrics (TP/FP/FN) as it
+    assigns all predictions to regions. Count metrics will be set to NaN.
+
+    Attributes:
+        cca_backend (CCABackend): Backend for connected component analysis.
+    """
+
+    def __init__(
+        self,
+        cca_backend: CCABackend = CCABackend.scipy,
+    ) -> None:
+        """
+        Initialize the RegionBasedMatching instance.
+
+        Args:
+            cca_backend (CCABackend): Backend for connected component analysis.
+        """
+        self._cca_backend = cca_backend
+
+    def _get_gt_regions(self, gt: np.ndarray) -> Tuple[np.ndarray, int]:
+        """
+        Get ground truth regions using connected components and distance transforms.
+
+        Args:
+            gt: Ground truth array
+
+        Returns:
+            Tuple of (region_map, num_features) where region_map assigns each pixel
+            to the closest ground truth region.
+        """
+        # Step 1: Connected Components
+        labeled_array, num_features = _connected_components(gt, self._cca_backend)
+
+        # Step 2: Compute distance transform for each region
+        distance_map = np.full(gt.shape, np.inf, dtype=np.float32)
+        region_map = np.zeros(gt.shape, dtype=np.int32)
+
+        for region_label in range(1, num_features + 1):
+            # Create region mask
+            region_mask = labeled_array == region_label
+
+            # Compute distance transform
+            distance = distance_transform_edt(~region_mask)
+
+            # Update pixels where this region is closer
+            update_mask = distance < distance_map
+            distance_map[update_mask] = distance[update_mask]
+            region_map[update_mask] = region_label
+
+        return region_map, num_features
+
+    def _match_instances(
+        self,
+        unmatched_instance_pair: UnmatchedInstancePair,
+        context: Optional[MatchingContext] = None,
+        **kwargs,
+    ) -> InstanceLabelMap:
+        """
+        Perform region-based instance matching.
+
+        Args:
+            unmatched_instance_pair (UnmatchedInstancePair): The unmatched instance pair to be matched.
+            context (Optional[MatchingContext]): The matching context.
+            **kwargs: Additional keyword arguments.
+
+        Returns:
+            InstanceLabelMap: The result of the region-based matching.
+        """
+        pred_arr = unmatched_instance_pair.prediction_arr
+        ref_arr = unmatched_instance_pair.reference_arr
+        pred_labels = unmatched_instance_pair.pred_labels
+
+        labelmap = InstanceLabelMap()
+
+        if len(pred_labels) == 0:
+            return labelmap
+
+        # Get ground truth regions
+        region_map, num_features = self._get_gt_regions(ref_arr)
+
+        # For each prediction instance, find which ground truth region it belongs to
+        for pred_label in pred_labels:
+            pred_mask = pred_arr == pred_label
+
+            # Find the most common region assignment for this prediction instance
+            pred_regions = region_map[pred_mask]
+
+            # Remove background (region 0)
+            pred_regions = pred_regions[pred_regions > 0]
+
+            if len(pred_regions) > 0:
+                # Assign to the most common region
+                unique_regions, counts = np.unique(pred_regions, return_counts=True)
+                most_common_region = unique_regions[np.argmax(counts)]
+
+                # Add to labelmap
+                labelmap.add_labelmap_entry(int(pred_label), int(most_common_region))
+
+        return labelmap
+
+    @classmethod
+    def _yaml_repr(cls, node) -> dict:
+        return {
+            "cca_backend": node._cca_backend,
+        }

panoptica/panoptica_evaluator.py

@@ -442,6 +442,11 @@ def panoptic_evaluate(
             else instance_metadata["original_num_refs"]
         )
 
+        # For region-based matching, set TP to NaN since it doesn't use traditional counting
+        tp_value = processing_pair.tp
+        if instance_matcher.__class__.__name__ == "RegionBasedMatching":
+            tp_value = np.nan
+
         processing_pair = PanopticaResult(
             reference_arr=processing_pair.reference_arr,
             prediction_arr=processing_pair.prediction_arr,
@@ -450,7 +455,7 @@ def panoptic_evaluate(
             num_ref_instances=final_num_ref_instances,
             num_ref_labels=instance_metadata["num_ref_labels"],
             label_group=label_group,
-            tp=processing_pair.tp,
+            tp=tp_value,
             list_metrics=processing_pair.list_metrics,
             global_metrics=global_metrics,
             edge_case_handler=edge_case_handler,