Merge pull request #167 from computational-cell-analytics/bioimage-cpp-migration

Migrate to bioimage-cpp

Author: constantinpape

CLAUDE.md

@@ -0,0 +1,93 @@
+# CLAUDE.md
+
+This file provides guidance to coding agents when working in this repository.
+
+## Overview
+
+SynapseNet (`synapse_net`) is a Python package for deep-learning segmentation and analysis of synapses in electron microscopy, mainly electron tomography. It ships U-Net models for segmenting synaptic vesicles, active zones, compartments, mitochondria, cristae, and ribbon-synapse structures, plus analysis tools (distance measurement, morphometry, vesicle pooling) and IMOD import/export. It is exposed three ways: a **napari plugin**, a **CLI** (console scripts), and a **Python library**.
+
+## Environment & Install
+
+The dependencies for this package are available via conda. We are in the process of removing dependency on nifty and vigra in favor of bioimage-cpp.
+This will enable pure pip installation.
+
+Here is how to set-up an environment with a clean installation.
+
+```bash
+conda env create -f environment.yaml   # creates env "synapse-net"
+conda activate synapse-net
+pip install -e .                        # or: pip install --no-deps -e .  (as CI does)
+```
+
+IMOD export commands additionally require the external IMOD suite installed.
+
+## Testing
+
+```bash
+python -m unittest discover -s test -v                       # full suite (what CI runs)
+python -m unittest test.test_inference.TestInference         # one test class
+python -m unittest test.test_inference.TestInference.test_run_segmentation  # one test
+```
+
+Tests are **network-dependent**: they download pretrained models (via `pooch` from GWDG ownCloud) and sample data the first time they run, into the OS cache dir (`pooch.os_cache("synapse-net")`). A first run is slow and requires internet; offline runs will fail on model/sample download, not on a code bug.
+
+## Common commands
+
+```bash
+# Build API docs (pdoc, google docstring format) — output to tmp/
+python build_doc.py -o
+
+# CLI entry points (defined in setup.py console_scripts):
+synapse_net.run_segmentation -i <mrc-or-dir> -o <out-dir> -m vesicles_3d
+synapse_net.export_to_imod_points  -i <mrc-dir> -s <seg-dir> -o <mod-dir>
+synapse_net.export_to_imod_objects -i <mrc-dir> -s <seg-dir> -o <mod-dir>
+synapse_net.run_supervised_training -h
+synapse_net.run_domain_adaptation   -h
+synapse_net.visualize_vesicle_pools -h
+```
+
+## Architecture
+
+The package is organized by capability. The big picture spans several modules:
+
+### Models & the model registry (`synapse_net/inference/inference.py`)
+Models are identified by a **`model_type` string** (`vesicles_3d`, `vesicles_2d`, `vesicles_cryo`, `active_zone`, `compartments`, `mitochondria`/`mitochondria2`, `cristae`/`cristae2`/`cristae3`, `ribbon`, plus CLI-only `vesicles_*` variants). `_get_model_registry()` maps each name to a sha256 + GWDG ownCloud download URL and fetches via `pooch`. `get_model(model_type)` downloads (if needed) and `torch.load`s the checkpoint.
+
+**Voxel-size scaling is central.** Each model was trained at a specific resolution (`get_model_training_resolution`, in nm/axis). `compute_scale_from_voxel_size(voxel_size, model_type)` produces a zyx scale factor so input data is resized to match training resolution before inference and resized back after. Always pair a model with the right `scale` rather than feeding raw-resolution data.
+
+### Inference dispatch
+`run_segmentation(image, model, model_type, tiling, scale, **kwargs)` is the top-level entry. It branches on `model_type` to a per-structure function in its own module: `inference/vesicles.py`, `active_zone.py`, `compartments.py`, `mitochondria.py`, `cristae.py`, `ribbon_synapse.py`, `actin.py`. The `ribbon` path is special — it runs multi-output prediction then optional post-processing that requires a vesicle segmentation passed via `extra_segmentation`.
+
+`inference/util.py` holds the shared machinery:
+- `get_prediction` / `get_prediction_torch_em` — tiled prediction via `torch_em.util.prediction.predict_with_halo`. **Tiling is a dict** `{"tile": {"z","y","x"}, "halo": {"z","y","x"}}`; note `get_prediction` subtracts `2*halo` from each tile dim before calling the torch_em layer.
+- `_Scaler` — applies/reverts the voxel-size scale (order-0 for segmentations).
+- `inference_helper` — drives file-in/file-out batch segmentation over a folder; used by the CLI.
+- A `bioimageio` prediction path exists in `get_prediction` but currently raises `NotImplementedError` (model format migration in progress; see the `bioimage-cpp-migration` branch).
+
+Per-structure **post-processing** lives in `inference/postprocessing/` (vesicles, compartments, membranes, ribbon, presynaptic_density), kept separate from raw prediction.
+
+### Training (`synapse_net/training/`)
+Built on `torch_em` U-Nets (`UNet2d`, `AnisotropicUNet`). Two regimes:
+- `supervised_training.py` — needs data + manual labels. Can init weights from a pretrained `model_type`.
+- `domain_adaptation.py` — unsupervised student-teacher (mean-teacher) adaptation to a new condition without labels; only works if the source model already partially detects the structure.
+- `semisupervised_training.py`, `transform.py` — supporting pieces.
+
+### Tools (`synapse_net/tools/`)
+- `cli.py` — argparse wrappers behind the console scripts.
+- Napari widgets (`segmentation_widget.py`, `distance_measure_widget.py`, `morphology_widget.py`, `vesicle_pool_widget.py`, `postprocessing_widget.py`, `size_filter_widget.py`), all subclassing `base_widget.py` and registered in `synapse_net/napari.yaml`. **When adding/renaming a widget, command, reader, or sample, update `napari.yaml` and the `setup.py` entry points to match.**
+- `volume_reader.py` — napari reader for `.mrc`/`.rec`/`.h5`.
+
+### Analysis & I/O
+- `distance_measurements.py`, `morphology.py` — quantification of segmentation results.
+- `imod/to_imod.py`, `imod/export.py` — convert segmentations to IMOD point models (spheres per vesicle) or closed-contour object models.
+- `file_utils.py` — I/O. `read_mrc` reads `.mrc`/`.rec`; **voxel sizes are converted Angstrom→nm by dividing by 10**, and data axes are flipped to match Python order. Also reads OME-Zarr and streams from the CryoET Data Portal (`cryoet_data_portal`, `s3fs` — both optional imports).
+- `ground_truth/` — annotation/GT-curation utilities (vesicle matching, shape refinement, edge filtering).
+
+### `scripts/`
+Project- and collaboration-specific scripts (subfolders `cooper`, `inner_ear`, `cryo`, `rizzoli`, `wichmann`, `portal`, `baselines`, …) tied to the publication and ongoing collaborations. These are **not** part of the installable package and not covered by tests; many paths are gitignored. Don't treat them as library API.
+
+## Conventions
+
+- Spatial dicts use explicit axis keys (`{"x":..,"y":..,"z":..}`); arrays are zyx. Be careful that tiling/scale lists are in zyx order while voxel-size dicts are keyed by axis name.
+- Public functions use Google-style docstrings (docs are generated with `pdoc -d google`); match that style when adding API.
+- Optional heavy deps (`zarr`, `s3fs`, `cryoet_data_portal`) are imported defensively and guarded — keep new optional integrations the same way.

environment.yaml

@@ -4,16 +4,16 @@ name:
     synapse-net
 dependencies:
     - bioimageio.core
+    - bioimage-cpp
     - kornia
     - magicgui
-    - napari >=0.5.0,<0.6.0
-    - nifty >=1.2.3
+    - napari >=0.5.0,<0.7.0
     - pip
     - pyqt
-    - python-elf
+    - python-elf >=0.9.0
     - pytorch
     - tensorboard
-    - torch_em >=0.8.1
+    - torch_em >=0.9.0
     - torchvision
     - trimesh
     - zarr <3

scripts/inner_ear/processing/process_manual_segmentation.py

@@ -6,7 +6,7 @@
 import imageio.v3 as imageio
 import mrcfile
 import numpy as np
-import nifty.tools as nt
+from bioimage_cpp.utils import take_dict
 from elf.io import open_file
 from skimage.transform import rescale
 from tqdm import tqdm
@@ -60,7 +60,7 @@ def export_vesicles(input_path, data_path, export_path):
     labels, label_names = process_labels(labels, label_names)
     imageio.imwrite(export_path, segmentation, compression="zlib")
 
-    segmentation = nt.takeDict(labels, segmentation)
+    segmentation = take_dict(labels, segmentation)
     export_pool_path = os.path.join(os.path.split(export_path)[0], "manual_pools.tif")
     imageio.imwrite(export_pool_path, segmentation, compression="zlib")
 

scripts/inner_ear/processing/run_analyis.py

@@ -2,11 +2,11 @@
 import warnings
 from pathlib import Path
 
+import bioimage_cpp as bic
 import imageio.v3 as imageio
 import mrcfile
 import numpy as np
 import pandas
-import vigra
 
 from synapse_net.file_utils import get_data_path
 from synapse_net.distance_measurements import (
@@ -302,8 +302,8 @@ def analyze_distances(
 def _relabel_vesicles(path):
     print("Relabel vesicles at", path)
     seg = _load_segmentation(path, None)
-    seg = vigra.analysis.labelVolumeWithBackground(seg.astype("uint32"))
-    seg, _, _ = vigra.analysis.relabelConsecutive(seg, start_label=1, keep_zeros=True)
+    seg = bic.segmentation.label(seg.astype("uint32"), connectivity=1)
+    seg, _, _ = bic.segmentation.relabel_sequential(seg, offset=1)
     imageio.imwrite(path, seg, compression="zlib")
 
 

setup.py

@@ -11,6 +11,32 @@
     author="Constantin Pape; Sarah Muth; Luca Freckmann",
     url="https://github.com/computational-cell-analytics/synapse-net",
     license="MIT",
+    install_requires=[
+        "bioimage-cpp",
+        "python-elf>=0.9.0",
+        "torch_em>=0.9.0",
+        "numpy",
+        "scipy",
+        "scikit-image",
+        "scikit-learn",
+        "h5py",
+        "imageio",
+        "mrcfile",
+        "pandas",
+        "pooch",
+        "tqdm",
+        "trimesh",
+    ],
+    extras_require={
+        # Dependencies for the napari plugin GUI. Install with `pip install synapse_net[napari]`.
+        "napari": [
+            "napari>=0.5.0,<0.7.0",
+            "magicgui",
+            "qtpy",
+            "superqt",
+            "napari-skimage-regionprops",
+        ],
+    },
     entry_points={
         "console_scripts": [
             "synapse_net.run_segmentation = synapse_net.tools.cli:segmentation_cli",

synapse_net/ground_truth/az_evaluation.py

@@ -4,13 +4,13 @@
 import h5py
 import pandas as pd
 import numpy as np
-import vigra
 
 from elf.evaluation.matching import _compute_scores, _compute_tps
 from elf.evaluation import dice_score
 from elf.segmentation.workflows import simple_multicut_workflow
 from scipy.ndimage import binary_dilation, binary_closing, distance_transform_edt, binary_opening
 from skimage.measure import label, regionprops, regionprops_table
+from skimage.morphology import local_maxima
 from skimage.segmentation import relabel_sequential, watershed
 from tqdm import tqdm
 
@@ -143,8 +143,7 @@ def _get_presynaptic_mask(boundary_map, vesicles):
 
     def _compute_mask_2d(z):
         distances = distance_transform_edt(boundary_map[z] < 0.25).astype("float32")
-        seeds = vigra.analysis.localMaxima(distances, marker=np.nan, allowAtBorder=True, allowPlateaus=True)
-        seeds = label(np.isnan(seeds))
+        seeds = label(local_maxima(distances, allow_borders=True))
         overseg = watershed(boundary_map[z], markers=seeds)
         seg = simple_multicut_workflow(
             boundary_map[z], use_2dws=False, watershed=overseg, n_threads=1, beta=0.6

synapse_net/ground_truth/shape_refinement.py

@@ -5,6 +5,7 @@
 
 import numpy as np
 
+import bioimage_cpp as bic
 from scipy.ndimage import binary_erosion, binary_dilation
 from skimage import img_as_ubyte
 from skimage.filters import gaussian, rank, sato, sobel
@@ -13,11 +14,6 @@
 from skimage.segmentation import watershed
 from tqdm import tqdm
 
-try:
-    import vigra
-except ImportError:
-    vigra = None
-
 FILTERS = ("sobel", "laplace", "ggm", "structure-tensor", "sato")
 
 
@@ -48,7 +44,7 @@ def edge_filter(
             - "sobel": Edges are found by smoothing the data and then applying a sobel filter.
             - "laplace": Edges are found with a laplacian of gaussian filter.
             - "ggm": Edges are found with a gaussian gradient magnitude filter.
-            - "structure-tensor": Edges are found based on the 2nd eigenvalue of the structure tensor.
+            - "structure-tensor": Edges are found based on the 1st eigenvalue of the structure tensor.
             - "sato": Edges are found with a sato-filter, followed by smoothing and leveling.
         per_slice: Compute the filter per slice instead of for the whole volume.
         n_threads: Number of threads for parallel computation over the slices.
@@ -58,8 +54,6 @@ def edge_filter(
     """
     if method not in FILTERS:
         raise ValueError(f"Invalid edge filter method: {method}. Expect one of {FILTERS}.")
-    if method in FILTERS[1:] and vigra is None:
-        raise ValueError(f"Filter {method} requires vigra.")
 
     if per_slice and data.ndim == 3:
         n_threads = mp.cpu_count() if n_threads is None else n_threads
@@ -73,14 +67,14 @@ def edge_filter(
         edge_map = gaussian(data, sigma=sigma)
         edge_map = sobel(edge_map)
     elif method == "laplace":
-        edge_map = vigra.filters.laplacianOfGaussian(data.astype("float32"), sigma)
+        edge_map = bic.filters.laplacian_of_gaussian(data.astype("float32"), sigma)
     elif method == "ggm":
-        edge_map = vigra.filters.gaussianGradientMagnitude(data.astype("float32"), sigma)
+        edge_map = bic.filters.gaussian_gradient_magnitude(data.astype("float32"), sigma)
     elif method == "structure-tensor":
         inner_scale, outer_scale = sigma, sigma * 0.5
-        edge_map = vigra.filters.structureTensorEigenvalues(
-            data.astype("float32"), innerScale=inner_scale, outerScale=outer_scale
-        )[..., 1]
+        edge_map = bic.filters.structure_tensor_eigenvalues(
+            data.astype("float32"), inner_scale, outer_scale
+        )[..., 0]
     elif method == "sato":
         edge_map = _sato_filter(data, sigma)
 

synapse_net/ground_truth/vesicles.py

@@ -14,13 +14,13 @@
 
 def _check_volume(raw, vesicles, labels, title=None, **extra_segmentations):
     import napari
-    from nifty.tools import takeDict
+    from bioimage_cpp.utils import take_dict
 
     if labels is None:
         label_vol = None
     else:
         labels[0] = 0
-        label_vol = takeDict(labels, vesicles)
+        label_vol = take_dict(labels, vesicles)
 
     v = napari.Viewer()
     if raw is not None:

synapse_net/inference/compartments.py

@@ -2,16 +2,15 @@
 from typing import Dict, List, Optional, Tuple, Union
 
 import numpy as np
-import vigra
 import torch
 
+import bioimage_cpp as bic
 import elf.segmentation as eseg
-import nifty
 from elf.tracking.tracking_utils import compute_edges_from_overlap
 from scipy.ndimage import distance_transform_edt, binary_closing
 from skimage.measure import label, regionprops
 from skimage.segmentation import watershed
-from skimage.morphology import remove_small_holes
+from skimage.morphology import remove_small_holes, local_maxima
 
 from synapse_net.inference.util import get_prediction, _Scaler, _postprocess_seg_3d
 
@@ -41,8 +40,7 @@ def _segment_compartments_2d(
     seeds[np.isin(seeds, remove_ids)] = 0
 
     # Compute the small seeds = local maxima of the in-plane distance map
-    small_seeds = vigra.analysis.localMaxima(distances_z, marker=np.nan, allowAtBorder=True, allowPlateaus=True)
-    small_seeds = label(np.isnan(small_seeds))
+    small_seeds = label(local_maxima(distances_z, allow_borders=True))
 
     # We only keep small seeds that don't intersect with a large seed.
     props = regionprops(small_seeds, seeds)
@@ -87,16 +85,16 @@ def _merge_segmentation_3d(seg_2d, beta=0.5, min_z_extent=10):
     overlaps = np.array([edge["score"] for edge in edges])
 
     n_nodes = int(seg_2d.max() + 1)
-    graph = nifty.graph.undirectedGraph(n_nodes)
-    graph.insertEdges(uv_ids)
+    graph = bic.graph.UndirectedGraph(n_nodes)
+    graph.insert_edges(uv_ids)
 
     costs = eseg.multicut.compute_edge_costs(1.0 - overlaps)
     # set background weights to be maximally repulsive
     bg_edges = (uv_ids == 0).any(axis=1)
     costs[bg_edges] = -8.0
 
     node_labels = eseg.multicut.multicut_decomposition(graph, costs, beta=beta)
-    segmentation = nifty.tools.take(node_labels, seg_2d)
+    segmentation = node_labels[seg_2d]
 
     if min_z_extent is not None and min_z_extent > 0:
         props = regionprops(segmentation)

synapse_net/tools/size_filter_widget.py

@@ -3,7 +3,7 @@
 import napari.viewer
 
 import numpy as np
-import nifty.tools as nt
+from bioimage_cpp.utils import take_dict
 
 from qtpy.QtWidgets import QHBoxLayout, QVBoxLayout, QPushButton, QSpinBox, QLabel
 from skimage.measure import regionprops, label
@@ -70,7 +70,7 @@ def _filter_segmentation(self, segmentation, size_threshold, apply_label):
         if apply_label:
             mapping = {prop.label: int(prop.max_intensity) for prop in props if prop.label not in filter_ids}
             mapping[0] = 0
-            segmentation = nt.takeDict(mapping, segmentation)
+            segmentation = take_dict(mapping, segmentation)
         return segmentation.astype(dtype)
 
     def on_size_filter(self):