diff --git a/src/copick_utils/cli/seg2fil.py b/src/copick_utils/cli/seg2fil.py index d7d5c95..c6ca9f9 100644 --- a/src/copick_utils/cli/seg2fil.py +++ b/src/copick_utils/cli/seg2fil.py @@ -59,6 +59,14 @@ def volume_in_angstrom3(value, unit, input_uri): description="Also write the instance segmentation of the traced filaments (each voxel holds the ID of its " "filament, the same IDs as the filaments). Smart defaults as for -o; ?instance=true is implied.", ) +@optgroup.option( + "--curve", + type=click.Choice(["catmull-rom", "bspline"]), + default="catmull-rom", + help="Curve stored for each filament: 'catmull-rom' control points through the fitted spline (within half a " + "voxel of it), which ChimeraX-copick and napari-copick edit freely; 'bspline' the exact fit, whose control points " + "can be moved but not added or removed.", +) @add_debug_option def seg2fil( config, @@ -80,15 +88,16 @@ def seg2fil( workers, output_uri, instances_uri, + curve, debug, ): """ Trace filaments in segmentations. Traces the centreline of every filament (e.g. microtubules or actin) in a segmentation and - stores it as a copick Filaments entry: the filament's fitted B-spline (an exact `bspline` - curve that editors can reopen) and the centreline points copick regenerates from it. - Filaments are numbered 1, 2, ... by length, longest first. + stores it as a copick Filaments entry: an editable Catmull-Rom curve through the filament's + fitted spline (`--curve bspline` stores the fit itself) and the centreline points copick + regenerates from it. Filaments are numbered 1, 2, ... by length, longest first. Each connected component is skeletonized and split into branches between ends and junctions. Holes up to `--fill-lumen` are filled first, so a tube labelled by its wall alone traces as one @@ -198,6 +207,7 @@ def seg2fil( smoothing=smoothing, extend_ends=extend_ends, label=label, + curve_kind=curve, **instance_params, ) diff --git a/src/copick_utils/converters/filaments_from_segmentation.py b/src/copick_utils/converters/filaments_from_segmentation.py index a46dda4..240cb45 100644 --- a/src/copick_utils/converters/filaments_from_segmentation.py +++ b/src/copick_utils/converters/filaments_from_segmentation.py @@ -67,13 +67,15 @@ def filaments_from_segmentation( instances_object_name: Optional[str] = None, instances_user_id: Optional[str] = None, instances_session_id: Optional[str] = None, + curve_kind: str = "catmull-rom", **kwargs, ) -> Optional[Tuple["CopickFilaments", Dict[str, float]]]: """ Trace the filaments of a segmentation and store them as copick Filaments. - Each filament is stored with its fitted spline as a ``bspline`` curve (copick regenerates the filament's points - from it) and the fit's settings in ``metadata["fit"]``. IDs are 1..K by length (longest first) for a binary, + Each filament is stored with an editable ``catmull-rom`` curve through its fitted spline, within half a voxel of + it (``curve_kind="bspline"`` stores the fit itself, whose control points editors can move but not add or remove); + copick regenerates the filament's points from the curve, and the fit's settings go in ``metadata["fit"]``. IDs are 1..K by length (longest first) for a binary, multilabel or panoptic-region input; an instance segmentation (or a panoptic segmentation's instances) keeps its IDs. With ``instances_*`` set, the instance segmentation of the traced filaments is stored as well: each label voxel holds the ID of the nearest filament of its connected component, with the same IDs as the @@ -101,12 +103,14 @@ def filaments_from_segmentation( instances_object_name: Object name of the instance segmentation to store, or None for none. instances_user_id: User ID of the instance segmentation. instances_session_id: Session ID of the instance segmentation ({input_session_id} is replaced). + curve_kind: Curve stored per filament: ``catmull-rom`` (editable, derived from the fit) or ``bspline`` (the + exact fit). **kwargs: Additional keyword arguments from the lazy converter. Returns: Tuple of (CopickFilaments, stats dict) or None if the operation failed. """ - from copick.models import CopickFilamentCurve + from copick.models import CopickFilament, CopickFilamentCurve import copick_utils from copick_utils.util.segmentations import new_segmentation @@ -139,16 +143,29 @@ def filaments_from_segmentation( source = f"{segmentation.name}:{segmentation.user_id}/{segmentation.session_id}@{voxel_size}" curves, metadata = [], [] for centreline in centrelines: - curves.append( - CopickFilamentCurve.from_tck( - centreline.tck, - step=voxel_size, - smoothing=centreline.smoothing, - scale=1.0, - ), + curve = CopickFilamentCurve.from_tck( + centreline.tck, + step=voxel_size, + smoothing=centreline.smoothing, + scale=1.0, ) + if curve_kind == "catmull-rom": + # The editors' own conversion: points a voxel apart, within half a voxel of the fit. + curve = CopickFilament.from_curve(centreline.instance_id, curve).editable_curve( + tolerance=voxel_size / 2, + step=voxel_size, + kinds=("catmull-rom",), + ) + curves.append(curve) fit = dict(centreline.metadata) - fit.update({"tool": "copick-utils seg2fil", "version": copick_utils.__version__, "source": source}) + fit.update( + { + "tool": "copick-utils seg2fil", + "version": copick_utils.__version__, + "source": source, + "stored_curve": curve_kind, + }, + ) metadata.append({"fit": fit}) filaments = run.new_filaments(object_name, session_id, user_id, exist_ok=True) diff --git a/tests/test_filaments.py b/tests/test_filaments.py index 896bc74..75de5ac 100644 --- a/tests/test_filaments.py +++ b/tests/test_filaments.py @@ -167,7 +167,7 @@ def test_bspline_tangent_matches_the_curve_direction(): # --- the commands ------------------------------------------------------------------------------------------------------ -def _traced(run, config_path, volume): +def _traced(run, config_path, volume, extra=()): _write_seg(run, volume) _invoke( seg2fil, @@ -182,6 +182,7 @@ def _traced(run, config_path, volume): "microtubule:trace/s@10.0?instance=true", "--min-length", "150", + *extra, ], ) run.refresh() @@ -190,15 +191,13 @@ def _traced(run, config_path, volume): return filaments, instances -def test_ids_are_shared_by_instances_filaments_and_picks(run, config_path): - volume = tube(tube(np.zeros(SHAPE, np.uint8), (20, 40, 5), (20, 40, 75), 3), (20, 5, 40), (20, 75, 40), 3) - filaments, instances = _traced(run, config_path, volume) - ids = sorted(f.instance_id for f in filaments.filaments) - assert ids == [1, 2] - assert sorted(np.unique(instances.numpy()[volume > 0]).tolist()) == ids - - # Each stored curve is the exact fit, and copick's points are its evaluation - for f in filaments.filaments: +def test_seg2fil_stores_the_exact_fit_on_request(run, config_path): + volume = tube(np.zeros(SHAPE, np.uint8), (20, 40, 5), (20, 40, 75), 3) + editable, _ = _traced(run, config_path, volume) + editable_points = {f.instance_id: np.asarray(f.points) for f in editable.filaments} + exact, _ = _traced(run, config_path, volume, extra=("--curve", "bspline")) + for f in exact.filaments: + # the stored curve is the fit, and copick's points are its evaluation assert f.curve.kind == "bspline" and f.curve_is_current() regenerated = evaluate_curve( f.curve.control_points, @@ -209,7 +208,24 @@ def test_ids_are_shared_by_instances_filaments_and_picks(run, config_path): f.curve.knots, ) np.testing.assert_allclose(np.asarray(f.points), regenerated, atol=1e-6) + assert f.metadata["fit"]["stored_curve"] == "bspline" + # the default Catmull-Rom curve stays within half a voxel (5 Å) of the fit + distance = np.linalg.norm(editable_points[f.instance_id][:, None] - np.asarray(f.points)[None], axis=2) + assert distance.min(axis=1).max() <= 5.0 + 1e-6 + + +def test_ids_are_shared_by_instances_filaments_and_picks(run, config_path): + volume = tube(tube(np.zeros(SHAPE, np.uint8), (20, 40, 5), (20, 40, 75), 3), (20, 5, 40), (20, 75, 40), 3) + filaments, instances = _traced(run, config_path, volume) + ids = sorted(f.instance_id for f in filaments.filaments) + assert ids == [1, 2] + assert sorted(np.unique(instances.numpy()[volume > 0]).tolist()) == ids + + # Each stored curve is an editable Catmull-Rom curve, and copick's points are its evaluation + for f in filaments.filaments: + assert f.curve.kind == "catmull-rom" and f.curve_is_current() assert f.metadata["fit"]["method"] == "splprep" and f.metadata["fit"]["tool"] == "copick-utils seg2fil" + assert f.metadata["fit"]["stored_curve"] == "catmull-rom" _invoke(fil2picks, ["-c", config_path, "-i", "microtubule:trace/s", "-o", "microtubule:trace/s", "--spacing", "82"]) run.refresh()