feat(preprocessing): add binary closing and hole filling option

vesskel/_napari.py

@@ -48,6 +48,10 @@ def _extraction_params(
             include_vessel_radius: bool = False,
             junction_cleanup: bool = False,
             cleanup_threshold_factor: float = 2.5,
+            fill_holes: bool = False,
+            closing_iterations: int = 0,
+            max_hole_size: int = 0,
+            show_preprocessed: bool = False,
         ) -> None:
             return None
 
@@ -67,6 +71,24 @@ def _extraction_params(
                 "max": 10.0,
                 "step": 0.1,
             },
+            fill_holes={"annotation": bool, "value": False},
+            closing_iterations={
+                "annotation": int,
+                "value": 0,
+                "widget_type": "SpinBox",
+                "min": 0,
+                "max": 10,
+                "step": 1,
+            },
+            max_hole_size={
+                "annotation": int,
+                "value": 0,
+                "widget_type": "SpinBox",
+                "min": 0,
+                "max": 100000,
+                "step": 100,
+            },
+            show_preprocessed={"annotation": bool, "value": False},
         )
 
         # ---------- output parameters (magicgui) ----------
@@ -129,6 +151,23 @@ def _output_params(
         advanced_group.append(self.junction_cleanup_widget)
         advanced_group.append(self.cleanup_threshold_widget)
 
+        self.fill_holes_widget = extraction_gui.fill_holes
+        self.fill_holes_widget.label = "Fill holes in segmentation"
+
+        self.closing_iterations_widget = extraction_gui.closing_iterations
+        self.closing_iterations_widget.label = "Closing iterations"
+
+        self.max_hole_size_widget = extraction_gui.max_hole_size
+        self.max_hole_size_widget.label = "Max hole size (pixels)"
+
+        self.show_preprocessed_widget = extraction_gui.show_preprocessed
+        self.show_preprocessed_widget.label = "Show preprocessed binary layer"
+
+        advanced_group.append(self.fill_holes_widget)
+        advanced_group.append(self.closing_iterations_widget)
+        advanced_group.append(self.max_hole_size_widget)
+        advanced_group.append(self.show_preprocessed_widget)
+
         # ============================================================
         # Output Settings (CLI file export options)
         # ============================================================
@@ -202,6 +241,10 @@ def _get_current_pipeline_config(self) -> PipelineConfig:
                 vessel_radius=self.include_vessel_radius_widget.value,
                 junction_cleanup=junction_cleanup,
                 cleanup_threshold_factor=self.cleanup_threshold_widget.value,
+                fill_holes=self.fill_holes_widget.value,
+                closing_iterations=self.closing_iterations_widget.value,
+                max_hole_size=self.max_hole_size_widget.value,
+                show_preprocessed=self.show_preprocessed_widget.value,
             ),
             output=OutputConfig(
                 write_skeleton_npy=self.write_skeleton_npy_widget.value,
@@ -221,6 +264,10 @@ def _set_pipeline_config(self, config: PipelineConfig) -> None:
         self.include_vessel_radius_widget.value = e.vessel_radius
         self.junction_cleanup_widget.value = e.junction_cleanup
         self.cleanup_threshold_widget.value = e.cleanup_threshold_factor
+        self.fill_holes_widget.value = e.fill_holes
+        self.closing_iterations_widget.value = e.closing_iterations
+        self.max_hole_size_widget.value = e.max_hole_size
+        self.show_preprocessed_widget.value = e.show_preprocessed
 
         o = config.output
         self.write_skeleton_npy_widget.value = o.write_skeleton_npy
@@ -299,6 +346,19 @@ def _on_analyze(self) -> None:
                 )
             )
 
+            # -- optional: show preprocessed binary layer ----------------
+            if (
+                pipeline_config.extraction.show_preprocessed
+                and result.preprocessed_binary is not None
+            ):
+                self.viewer.add_layer(
+                    Layer.create(
+                        result.preprocessed_binary,
+                        {"name": f"{img.name}_preprocessed"},
+                        "labels",
+                    )
+                )
+
             for data, meta, layer_type in result.layers:
                 try:
                     layer = Layer.create(data, meta, layer_type)

vesskel/config.py

@@ -31,6 +31,10 @@ class ExtractionConfig:
     vessel_radius: bool = False
     junction_cleanup: bool = False
     cleanup_threshold_factor: float = 2.5
+    closing_iterations: int = 0
+    fill_holes: bool = False
+    max_hole_size: int = 0
+    show_preprocessed: bool = False
 
     def to_dict(self) -> dict[str, Any]:
         return {
@@ -41,6 +45,10 @@ def to_dict(self) -> dict[str, Any]:
             "vessel_radius": self.vessel_radius,
             "junction_cleanup": self.junction_cleanup,
             "cleanup_threshold_factor": self.cleanup_threshold_factor,
+            "closing_iterations": self.closing_iterations,
+            "fill_holes": self.fill_holes,
+            "max_hole_size": self.max_hole_size,
+            "show_preprocessed": self.show_preprocessed,
         }
 
     @classmethod
@@ -54,6 +62,10 @@ def from_dict(cls, data: dict[str, Any]) -> ExtractionConfig:
             vessel_radius=data.get("vessel_radius", False),
             junction_cleanup=data.get("junction_cleanup", False),
             cleanup_threshold_factor=data.get("cleanup_threshold_factor", 2.5),
+            closing_iterations=data.get("closing_iterations", 0),
+            fill_holes=data.get("fill_holes", False),
+            max_hole_size=data.get("max_hole_size", 0),
+            show_preprocessed=data.get("show_preprocessed", False),
         )
 
 

vesskel/pipeline.py

@@ -6,6 +6,7 @@
 from typing import TYPE_CHECKING
 
 import numpy as np
+from scipy import ndimage as ndi
 from skan import summarize
 
 from vesskel._utils import to_binary
@@ -34,6 +35,7 @@ class AnalysisResult:
     summary_features: dict[str, float]
     branch_records: list[dict[str, object]]
     radius_matrix: np.ndarray | None = None
+    preprocessed_binary: np.ndarray | None = None
 
 
 def analyze_binary_image(
@@ -53,6 +55,33 @@ def analyze_binary_image(
         Full pipeline configuration.
     """
     binary = to_binary(image)
+
+    # -- optional: morphological preprocessing -------------------------
+    preprocessed_binary: np.ndarray | None = None
+    if config.extraction.closing_iterations > 0:
+        structure = ndi.generate_binary_structure(binary.ndim, 1)
+        binary = ndi.binary_closing(
+            binary, structure=structure, iterations=config.extraction.closing_iterations
+        ).astype(binary.dtype)
+
+    if config.extraction.fill_holes:
+        before_fill = binary.copy()
+        binary = ndi.binary_fill_holes(binary).astype(binary.dtype)
+
+        if config.extraction.max_hole_size > 0:
+            diff = binary.astype(np.int8) - before_fill.astype(np.int8)
+            filled = diff > 0
+            if filled.any():
+                labels, n = ndi.label(filled)
+                sizes = np.bincount(labels.ravel())
+                big = sizes > config.extraction.max_hole_size
+                big[0] = False
+                revert = big[labels]
+                binary[revert] = 0
+
+    if config.extraction.closing_iterations > 0 or config.extraction.fill_holes:
+        preprocessed_binary = binary
+
     skeleton = lee94_thin(binary)
 
     if not skeleton.any():
@@ -129,4 +158,5 @@ def analyze_binary_image(
         summary_features=summary_features,
         branch_records=branch_records,
         radius_matrix=radius_matrix,
+        preprocessed_binary=preprocessed_binary,
     )