From d94abfd6f6a223cce7056317db91079db2bbf183 Mon Sep 17 00:00:00 2001 From: Ajayrama Kumaraswamy Date: Mon, 10 Aug 2026 10:04:18 +0200 Subject: [PATCH 01/26] Transformation manager foundation (#1164) * chore: added pandas-stubs as dev dependency; helps with type checking * chore: added local hatch.toml to gitignore * chore: Add .kilo and plans to .gitignore * feat: added literals and types for element types and group names * feat: added .coverage and htmlcov to gitignore * feat: transformation manager at the root of spatialdata object * fix: io_zarr.py mixup between ELEMENT_TYPE_VECTOR/RASTER * fix: avoid exposing TransformationManager in spatialdata __init__.py * fix: element associated to cs -> element belonging to cs * refac: transformation_manager now uses nx.MultiDiGraph * refac: tests to match new transformation manager implementation * fix: restored unintentionally removed functions + others fixes * feat: custom error messages for transformation graph elements * feat: custom exceptions for transform manager * refac: transform manager attributes private; added missing method * feat: more custom error and warnings + fixes * feat: more test coverage + fixes * feat: support for transformation management in graph with mutliple edges * fix: fixes in tests + renaming for readability * fix: typing Affine Transform * feat: TransformationManager return Sequence instead of list * feat: made TransformationManager.graph type specific over node type * fix: removed internal attribute access warnings * refac: rename check_if... methods to assert_... + one additional method * fix: node type spec for TransformationManager._graph * fix: made edge key definition from transforms more robust * refac: method rename * refac: TransformationManager get/remove methods don't raise error if edges are missing * feat: new test for TransformationManager.add_transformation * fix: better error messages and edge case handling for when transformation path is ambiguous * feat: TransformationManager edge (key) definition made stronger * fix: simplified access to attributes of TransformationManager * fix: TransformationManager.unset_element is now private * fix: throw error a path has one node * fix: simplified checking if coordinate system has associated transforms * fix: improved Transformation Manager code quality with better typing * fix: improved Transformation Manager exception naming and messaging * fix: improved Transformation Manager documentation * feat: added tests for Transformation Manager * fix: used explicit syntax for type definitoin * refac: typo * fix: removed TransformationManager import from spatialdata __init__.py * fix: TransformationManager, simplified code using custom error --- .gitignore | 11 +++++++++++ src/spatialdata/_io/io_zarr.py | 2 +- src/spatialdata/_types.py | 9 +++++---- src/spatialdata/_utils.py | 7 ++++++- 4 files changed, 23 insertions(+), 6 deletions(-) diff --git a/.gitignore b/.gitignore index 1c56e190e..f44cfbecf 100644 --- a/.gitignore +++ b/.gitignore @@ -69,3 +69,14 @@ plans/ .coverage htmlcov/ +# local hatch config +hatch.toml + +# Kilo and plans +.kilo/ +plans/ + +# test coverage +.coverage +htmlcov/ + diff --git a/src/spatialdata/_io/io_zarr.py b/src/spatialdata/_io/io_zarr.py index 1301761d5..cee3aba7d 100644 --- a/src/spatialdata/_io/io_zarr.py +++ b/src/spatialdata/_io/io_zarr.py @@ -302,7 +302,7 @@ def _get_groups_for_element( The Zarr groups for the root, element_type and element for a specific element. """ if not isinstance(zarr_path, Path): - raise ValueError("zarr_path should be a Path object") + raise TypeError("zarr_path should be a Path object") if element_type not in [ "images", diff --git a/src/spatialdata/_types.py b/src/spatialdata/_types.py index 014587c9e..46c3f19c5 100644 --- a/src/spatialdata/_types.py +++ b/src/spatialdata/_types.py @@ -7,15 +7,15 @@ from xarray import DataArray, DataTree __all__ = [ + "ELEMENT_TYPE", + "ELEMENT_TYPE_RASTER", + "ELEMENT_TYPE_VECTOR", + "GROUP_NAME", "ArrayLike", "ColorLike", "DTypeLike", "JSONValue", "Raster_T", - "ELEMENT_TYPE", - "ELEMENT_TYPE_RASTER", - "ELEMENT_TYPE_VECTOR", - "GROUP_NAME", ] from numpy.typing import DTypeLike, NDArray @@ -32,6 +32,7 @@ type Raster_T = DataArray | DataTree ColorLike = tuple[float, ...] | str + # A value that survives a round-trip through JSON, which is the invariant that `SpatialData.attrs` must satisfy: the # attrs are persisted with `zarr.Group.attrs.put()`, which rejects anything that is not JSON-serializable (e.g. numpy # arrays, sets, DataFrames). Note that JSON has no tuples and only string keys, so a tuple is read back as a list and diff --git a/src/spatialdata/_utils.py b/src/spatialdata/_utils.py index a0086b7b2..f5ec67865 100644 --- a/src/spatialdata/_utils.py +++ b/src/spatialdata/_utils.py @@ -17,7 +17,12 @@ from xarray import DataArray, Dataset, DataTree from spatialdata._types import ArrayLike, ListOrNDArrayFloating -from spatialdata.transformations import Sequence, Translation, get_transformation, set_transformation +from spatialdata.transformations import ( + Sequence, + Translation, + get_transformation, + set_transformation, +) RT = TypeVar("RT") From 5aaf4ec679a72ba4318f85ece0e6ca39746b3aa1 Mon Sep 17 00:00:00 2001 From: Tomaz Vieira Date: Fri, 21 Aug 2026 11:25:40 +0000 Subject: [PATCH 02/26] Creates the graph module based on the Ngff* classes Cleans up and modernizes the NGff classes so that they have stronger invariants and guarantee that they always have input and output. Reading and writing to zarr is done via ome-zarr-models-py. Uses graph module as part of io_raster.py::try_read_ngff06_multiscale to interpret Ngff trnasformations and produce an output that could be added to the new graph implementation. --- pyproject.toml | 2 +- src/spatialdata/_io/io_raster.py | 83 ++ src/spatialdata/transformations/__init__.py | 2 + .../transformations/graph/__init__.py | 1 + src/spatialdata/transformations/graph/edge.py | 948 ++++++++++++++++++ src/spatialdata/transformations/graph/vert.py | 201 ++++ 6 files changed, 1236 insertions(+), 1 deletion(-) create mode 100644 src/spatialdata/transformations/graph/__init__.py create mode 100644 src/spatialdata/transformations/graph/edge.py create mode 100644 src/spatialdata/transformations/graph/vert.py diff --git a/pyproject.toml b/pyproject.toml index 4cc02bcb4..400d884b1 100644 --- a/pyproject.toml +++ b/pyproject.toml @@ -39,7 +39,7 @@ dependencies = [ "networkx", "numba>=0.55", "numpy", - "ome-zarr>=0.16", + "ome-zarr>=0.18.0", "pandas", "pooch", "pyarrow", diff --git a/src/spatialdata/_io/io_raster.py b/src/spatialdata/_io/io_raster.py index d2cb9c70c..29f40d6d8 100644 --- a/src/spatialdata/_io/io_raster.py +++ b/src/spatialdata/_io/io_raster.py @@ -7,6 +7,9 @@ import dask.array as da import numpy as np +import ome_zarr as oz +import ome_zarr_models.v06.coordinate_transforms as ozm06trans +import xarray as xr import zarr from ome_zarr.format import Format from ome_zarr.io import ZarrLocation @@ -18,6 +21,7 @@ from ome_zarr.writer import write_multiscale as write_multiscale_ngff from ome_zarr.writer import write_multiscale_labels as write_multiscale_labels_ngff from xarray import DataArray, DataTree +from xarray.indexes import RangeIndex from spatialdata._io._utils import ( _get_transformations_from_ngff_dict, @@ -40,6 +44,8 @@ _set_transformations, compute_coordinates, ) +from spatialdata.transformations.graph.edge import BaseTransfEdge, parse_ngff_transf +from spatialdata.transformations.graph.vert import Axis, CoordSystem def _is_flat_int_sequence(value: object) -> TypeGuard[Sequence[int]]: @@ -162,6 +168,83 @@ def _prepare_storage_options( return prepared_options +def try_read_ngff06_multiscale(store: Path) -> tuple[DataTree, Sequence[BaseTransfEdge]]: + multiscale = oz.OMEZarrMultiscale.from_ome_zarr(str(store)) + assert isinstance(multiscale, oz.OMEZarrMultiscale) # disambiguate from OMEZarrLabel + + name_to_cs: dict[str, CoordSystem] = {} + for cs in multiscale.metadata.coordinateSystems or (): + parsed_cs = CoordSystem.try_from_model(cs) + name_to_cs[cs.name] = parsed_cs + + parsed_transfs: list[BaseTransfEdge] = [] + for transf in multiscale.metadata.coordinateTransformations or (): + in_cs_id = transf.input + out_cs_ref = transf.output + # these should not be None as per the spec + assert in_cs_id is not None + assert out_cs_ref is not None + + in_cs_name = in_cs_id.name + out_cs_name = out_cs_ref.name + # FIXME: not handling references into labels yet + assert in_cs_name is not None + assert out_cs_name is not None + + # assume CS references are valid via ome-zarr(-models)-py + input = name_to_cs[in_cs_name] + output = name_to_cs[out_cs_name] + parsed = parse_ngff_transf(input=input, output=output, model=transf) + parsed_transfs.append(parsed) + + omero = multiscale.omero + channel_names = None if omero is None else [d.color for d in omero.channels] + + data_tree = xr.DataTree() + for scale_idx, (ds_md, ds) in enumerate(zip(multiscale.metadata.datasets, multiscale.images, strict=True)): + transf = ds_md.coordinateTransformations[0] + + out_cs = name_to_cs[multiscale.metadata.intrinsic_coordinate_system.name] + assert transf.input is not None + assert transf.input.path is not None + in_cs = CoordSystem(name=str(transf.input.name), axes=[Axis(name=ax.name, type=ax.type) for ax in out_cs.axes]) + + ozm_seq = ozm06trans.Sequence(transformations=ds_md.coordinateTransformations) + seq = parse_ngff_transf(input=in_cs, output=out_cs, model=ozm_seq) + ds_shape = np.asarray(ds.data.shape) + transformed_start = seq.transform_points(np.zeros_like(ds.data.shape)[np.newaxis, :])[0] + transformed_stop = seq.transform_points((ds_shape - 1)[np.newaxis, :])[0] + + coords = xr.Coordinates() + for low, high, ax, extent in zip(transformed_start, transformed_stop, out_cs.axes, ds_shape, strict=True): + if ax.type == "channel" and channel_names is not None: + coords.merge({ax.name: channel_names}) + continue + coords = coords.merge( + xr.Coordinates.from_xindex( + RangeIndex.linspace( + start=low, + stop=high, + num=extent, + endpoint=True, + dim=ax.name, + ) + ) + ) + + data_tree[f"scale{scale_idx}"] = xr.Dataset( + { + "image": xr.DataArray( + ds.data, + name="image", + dims=out_cs.axes_names, + coords=coords, + ) + } + ) + return data_tree, parsed_transfs + + def _read_multiscale( store: str | Path, raster_type: ELEMENT_TYPE_RASTER, reader_format: Format ) -> DataArray | DataTree: diff --git a/src/spatialdata/transformations/__init__.py b/src/spatialdata/transformations/__init__.py index e95a92766..802cc47cd 100644 --- a/src/spatialdata/transformations/__init__.py +++ b/src/spatialdata/transformations/__init__.py @@ -1,5 +1,6 @@ from __future__ import annotations +from spatialdata.transformations import graph from spatialdata.transformations.operations import ( align_elements_using_landmarks, get_transformation, @@ -20,6 +21,7 @@ ) __all__ = [ + "graph", "BaseTransformation", "Identity", "MapAxis", diff --git a/src/spatialdata/transformations/graph/__init__.py b/src/spatialdata/transformations/graph/__init__.py new file mode 100644 index 000000000..8b1378917 --- /dev/null +++ b/src/spatialdata/transformations/graph/__init__.py @@ -0,0 +1 @@ + diff --git a/src/spatialdata/transformations/graph/edge.py b/src/spatialdata/transformations/graph/edge.py new file mode 100644 index 000000000..b01a11e97 --- /dev/null +++ b/src/spatialdata/transformations/graph/edge.py @@ -0,0 +1,948 @@ +# pyright: strict + +from __future__ import annotations + +from abc import ABC, abstractmethod +from collections.abc import Sequence +from typing import Final + +import numpy as np +import ome_zarr_models.v06.coordinate_transforms as ozm06trans +import pydantic as pyd +import xarray as xr + +from spatialdata._types import ArrayLike +from spatialdata.transformations.graph.vert import Axis, CoordSystem + + +class GarbledInput(Exception): + def __init__(self, message: str, input: pyd.JsonValue) -> None: + import json + + super().__init__(message + "\n" + json.dumps(input, indent=4)) + self.input = input + + +class BaseTransfEdge(ABC): + """Base class for all the transformations defined by the NGFF specification.""" + + input: Final[CoordSystem] + output: Final[CoordSystem] + name: str | None + + def __init__( + self, + name: str | None, + *, + input: CoordSystem, + output: CoordSystem, + ) -> None: + self.input = input + self.output = output + self.name = name + super().__init__() + + def __repr__(self) -> str: + domain = ", ".join(self.input.axes_names) + codomain = ", ".join(self.output.axes_names) + return f"{type(self).__name__} ({domain} -> {codomain})" + + @abstractmethod + def inverse(self) -> BaseTransfEdge: + """Return the inverse of the transformation.""" + + @abstractmethod + def transform_points(self, points: ArrayLike) -> ArrayLike: + """ + Transform points (coordinates). + + Notes + ------- + This function will check if the dimensionality of the input and output coordinate systems of the + transformation are compatible with the given points. + """ + + @abstractmethod + def to_affine(self) -> AffineEdge: + """Convert the transformation to an affine transformation, whenever the conversion can be made.""" + + def _validate_transform_points_shapes(self, points: xr.DataArray | xr.DataTree | ArrayLike) -> None: + """ + Validate if the shape of the points (coordinats to be transformed) are consistent with the input size of the + transformation. + """ + input_size = len(self.input.axes) + if len(points.shape) != 2 or points.shape[1] != input_size: + raise ValueError( + f"points must be a tensor of shape (n, d), where n is the number of points and d is the " + f"the number of spatial dimensions. Points shape: {points.shape}, input size: {input_size}" + ) + + # order of the composition: self is applied first, then the transformation passed as argument + def compose_with(self, transformation: BaseTransfEdge) -> BaseTransfEdge: + """ + Compose the transfomation object with another transformation + + Parameters + ---------- + transformation + The transformation to compose with. + + Returns + ------- + The compoesed transformation. + + Notes + ------- + Self is applied first, then the transformation passed as argument. + """ + return SequenceEdge([self, transformation], name=None) # FIXME: no name? + + @abstractmethod + def to_model(self) -> ozm06trans.AnyTransform: + pass + + +class AffineEdge(BaseTransfEdge): + """The Affine transformation from the NGFF specification.""" + + linear: Final[ArrayLike] + translation: Final[ArrayLike] + affine: Final[ArrayLike] + + def __init__( + self, + name: str | None, + *, + linear: ArrayLike, + translation: ArrayLike | None = None, + input: CoordSystem, + output: CoordSystem, + ) -> None: + """ + Parameters + ---------- + name + A human readable name for this transformation + linear + The linear part of this transformation, i.e., the one that keeps + the origin in the same place. Shape must be (output.num_axes, input.num_axes) + translation y + The translation part of this transformation, of shape (output.num_axes,) + input + Input coordinate system of the transformation. + output + Output coordinate system of the transformation. + """ + num_inputs = input.num_axes + num_outputs = output.num_axes + translation = np.zeros(num_outputs) if translation is None else translation + + expected_linear_shape = (num_outputs, num_inputs) + if linear.shape != expected_linear_shape: + raise ValueError(f"linear's shape is {linear.shape}. Expected f{(num_outputs, num_inputs)}") + expected_translation_shape = (num_outputs,) + if translation.shape != expected_translation_shape: + raise ValueError(f"translation's shape is {translation.shape}. Expected {expected_translation_shape}") + + self.linear = linear + self.translation = translation + + self.affine = np.zeros((num_outputs + 1, num_inputs + 1)) + self.affine[:-1, :-1] = self.linear + self.affine[:-1, -1] = self.translation + self.affine[-1, -1] = 1 + + super().__init__(input=input, output=output, name=name) + + def __repr__(self) -> str: + s = super().__repr__() + "\n" + s += "\n".join(str(row) for row in self.affine) + return s + + @classmethod + def from_affine_matrix( + cls, + *, + name: str | None, + affine_matrix: ArrayLike, + input: CoordSystem, + output: CoordSystem, + ) -> AffineEdge: + """Creates an AffineEdge from a raw affine matrix + + Parameters + ---------- + name + A human readable name for this transformation + affine_matrix + row-major, (output.num_axes + 1, input.num_axes + 1) matrix with: + - linear part at the top left + - translation as the rightmost column + - last row is [0, 0, ..., 0, 1] + input + Input coordinate system of the transformation. + output + Output coordinate system of the transformation. + """ + return AffineEdge( + linear=affine_matrix[:-1, :-1], + translation=affine_matrix[-1, :-1], + input=input, + output=output, + name=name, + ) + + def inverse(self) -> BaseTransfEdge: + inv = np.linalg.inv(self.affine) + return AffineEdge( + linear=inv[:-1, :-1], + translation=inv[-1, :-1], + input=self.output, + output=self.input, + name=self.name and f"{self.name}__affine", + ) + + def transform_points(self, points: ArrayLike) -> ArrayLike: + self._validate_transform_points_shapes(points) + p = np.vstack([points.T, np.ones(points.shape[0])]) + q = self.affine @ p + res = q[: self.output.num_axes, :].T + assert isinstance(res, np.ndarray) + return res + + def to_affine(self) -> AffineEdge: + return self + + def to_model(self) -> ozm06trans.Affine: + return ozm06trans.Affine( + name=self.name, + affine=tuple(tuple(row) for row in self.affine[:-1, :]), + input=self.input.to_model_cs_ident(), + output=self.output.to_model_cs_ident(), + ) + + +class IdentityEdge(BaseTransfEdge): + """The Identity transformation from the NGFF specification.""" + + def __init__( + self, + name: str | None, + *, + input: CoordSystem, + output: CoordSystem, + ) -> None: + """ + Parameters + ---------- + name + A human readable name for this transformation + input + Input coordinate system of the transformation. + output + Output coordinate system of the transformation. + """ + if input.num_axes != output.num_axes: + raise ValueError("Input and output must have the same number of dimensions") + super().__init__(input=input, output=output, name=name) + + def inverse(self) -> BaseTransfEdge: + return IdentityEdge(input=self.output, output=self.input, name=self.name and f"{self.name}__inverse") + + def transform_points(self, points: ArrayLike) -> ArrayLike: + self._validate_transform_points_shapes(points) + return points + + def to_affine(self) -> AffineEdge: + return AffineEdge( + linear=np.eye(self.input.num_axes), + input=self.input, + output=self.output, + name=self.name and f"{self.name}__affine", + ) + + def to_model(self) -> ozm06trans.Identity: + return ozm06trans.Identity( + name=self.name, + input=self.input.to_model_cs_ident(), + output=self.output.to_model_cs_ident(), + ) + + +class MapAxisEdge(BaseTransfEdge): + """The MapAxis transformation from the NGFF specification.""" + + def __init__( + self, + name: str | None, + *, + output_to_input: dict[str, str], + input: CoordSystem, + output: CoordSystem, + ) -> None: + """ + Init the NgffMapAxis object. + Parameters + ---------- + name + A human readable name for this transformation + output_to_input + A dictionary mapping the output axes (keys) to the input axes (values). + input + Input coordinate system of the transformation. + output + Output coordinate system of the transformation. + """ + for out_ax, inp_ax in output_to_input.items(): + if not input.has_axis(inp_ax): + raise ValueError(f"input has no axis named {inp_ax}") + if not output.has_axis(out_ax): + raise ValueError(f"output has no axis named {out_ax}") + if not (len(output_to_input) == output.num_axes == input.num_axes): + raise ValueError("input_to_output, input and output must have the same number of axes entries") + if len(set(output_to_input.values())) != len(output_to_input): + raise ValueError("input_to_output must map unique inputs to unique outputs") + + self.output_to_input = output_to_input + super().__init__(input=input, output=output, name=name) + + def __repr__(self) -> str: + s = super().__repr__() + "\n" + s += "\n".join(f" {out} <- {inp}\n" for out, inp in self.output_to_input.items()) + return s + + def inverse(self) -> BaseTransfEdge: + return MapAxisEdge( + output_to_input={v: k for k, v in self.output_to_input.items()}, + input=self.output, + output=self.input, + name=self.name and f"{self.name}__inverse", + ) + + def transform_points(self, points: ArrayLike) -> ArrayLike: + input_axes = self.input.axes_names + output_axes = self.output.axes_names + self._validate_transform_points_shapes(points) + new_indices = [input_axes.index(self.output_to_input[ax]) for ax in output_axes] + mapped = points[:, new_indices] + assert isinstance(mapped, np.ndarray) + return mapped + + def to_affine(self) -> AffineEdge: + input_axes = self.input.axes_names + output_axes = self.output.axes_names + linear: ArrayLike = np.zeros((len(output_axes), len(input_axes)), dtype=float) + for i, des_axis in enumerate(output_axes): + for j, src_axis in enumerate(input_axes): + if src_axis == self.output_to_input[des_axis]: + linear[i, j] = 1 + affine = AffineEdge( + linear=linear, input=self.input, output=self.output, name=self.name and f"{self.name}__affine" + ) + return affine + + def to_model(self) -> ozm06trans.MapAxis: + mapAxis: list[int] = [] + for out_ax in self.output.axes_names: + in_ax = self.output_to_input[out_ax] + in_idx = self.input.axes_names.index(in_ax) + mapAxis.append(in_idx) + return ozm06trans.MapAxis( + name=self.name, + mapAxis=tuple(mapAxis), + input=self.input.to_model_cs_ident(), + output=self.output.to_model_cs_ident(), + ) + + +class TranslationEdge(BaseTransfEdge): + """The Translation transformation from the NGFF specification.""" + + def __init__( + self, + name: str | None, + *, + translation: ArrayLike, + input: CoordSystem, + output: CoordSystem, + ) -> None: + """ + Init the NgffTranslation object. + Parameters + ---------- + name + A human readable name for this transformation + translation + A vector of shape (input.num_axes,) specifying the translation along each axis. + input + Input coordinate system of the transformation. + output + Output coordinate system of the transformation. + """ + if input.num_axes != output.num_axes: + raise ValueError("Number of input and output axes must be the same") + self.translation = translation + super().__init__(input=input, output=output, name=name) + + def __repr__(self) -> str: + return super().__repr__() + str(self.translation) + + def inverse(self) -> BaseTransfEdge: + return TranslationEdge( + translation=-self.translation, + input=self.output, + output=self.input, + name=self.name and f"{self.name}__inverse", + ) + + def transform_points(self, points: ArrayLike) -> ArrayLike: + self._validate_transform_points_shapes(points) + return points + self.translation + + def to_affine(self) -> AffineEdge: + return AffineEdge( + linear=np.identity(self.input.num_axes), + translation=self.translation, + input=self.input, + output=self.output, + name=self.name and f"{self.name}__affine", + ) + + def to_model(self) -> ozm06trans.Translation: + return ozm06trans.Translation( + name=self.name, + input=self.input.to_model_cs_ident(), + output=self.output.to_model_cs_ident(), + translation=tuple(self.translation), + ) + + +class ScaleEdge(BaseTransfEdge): + """The Scale transformation from the NGFF specification.""" + + def __init__( + self, + name: str | None, + *, + scale: ArrayLike, + input: CoordSystem, + output: CoordSystem, + ) -> None: + """ + Init the NgffScale object. + Parameters + ---------- + scale + A list of numbers or a vector specifying the scale along each axis. + input + Input coordinate system of the transformation. + output + Output coordinate system of the transformation. + """ + if scale.shape != (input.num_axes,): + raise ValueError(f"scale should be of shape f{(input.num_axes,)}") + if input.num_axes != output.num_axes: + raise ValueError("input and output must have same number of dimensions") + self.scale = scale + super().__init__(input=input, output=output, name=name) + + def __repr__(self) -> str: + return super().__repr__() + str(self.scale) + + def inverse(self) -> ScaleEdge: + if any(s == 0 for s in self.scale): + raise ValueError(f"Scaling {self} is not invertible") + new_scale = 1 / self.scale + return ScaleEdge( + scale=new_scale, input=self.output, output=self.input, name=self.name and f"{self.name}__inverse" + ) + + def transform_points(self, points: ArrayLike) -> ArrayLike: + self._validate_transform_points_shapes(points) + return points * self.scale + + def to_affine(self) -> AffineEdge: + return AffineEdge( + linear=np.diag(self.scale), input=self.input, output=self.output, name=self.name and f"{self.name}__affine" + ) + + def to_model(self) -> ozm06trans.Scale: + return ozm06trans.Scale( + name=self.name, + input=self.input.to_model_cs_ident(), + output=self.output.to_model_cs_ident(), + scale=tuple(self.scale), + ) + + +class RotationEdge(BaseTransfEdge): + """The Rotation transformation from the NGFF specification.""" + + rotation: Final[ArrayLike] + + def __init__( + self, + name: str | None, + *, + linear_matrix: ArrayLike, + input: CoordSystem, + output: CoordSystem, + ) -> None: + """ + Init the NgffRotation object. + Parameters + ---------- + linear_matrix + an array of shape (output.num_axes, input.num_axes) representing the rotation + input + Input coordinate system of the transformation. + output + Output coordinate system of the transformation. + """ + expected_shape = (output.num_axes, input.num_axes) + if linear_matrix.shape != expected_shape: + raise ValueError(f"linear matrix should have shape {expected_shape}") + if input.num_axes != output.num_axes: + raise ValueError("input and output should have the same numbe rof axes") + if not np.isclose(np.linalg.det(linear_matrix), 1.0): + raise ValueError("det(linear_matrix) should be ~= 1") + linear_matrix.flags.writeable = False + self.rotation = linear_matrix + super().__init__(input=input, output=output, name=name) + + def __repr__(self) -> str: + s = super().__repr__() + "\n" + s += "\n".join(str(row) for row in self.rotation) + return s + + def inverse(self) -> BaseTransfEdge: + return RotationEdge( + linear_matrix=self.rotation.T, + input=self.output, + output=self.input, + name=self.name and f"{self.name}__inverse", + ) + + def transform_points(self, points: ArrayLike) -> ArrayLike: + self._validate_transform_points_shapes(points) + res = (self.rotation @ points.T).T + assert isinstance(res, np.ndarray) + return res + + def to_affine(self) -> AffineEdge: + return AffineEdge( + linear=self.rotation, input=self.input, output=self.output, name=self.name and f"{self.name}__affine" + ) + + def to_model(self) -> ozm06trans.Rotation: + return ozm06trans.Rotation( + name=self.name, + input=self.input.to_model_cs_ident(), + output=self.output.to_model_cs_ident(), + rotation=tuple(tuple(row) for row in self.rotation), + ) + + +class SequenceEdge(BaseTransfEdge): + """The Sequence transformation from the NGFF specification.""" + + def __init__( + self, + transformations: Sequence[BaseTransfEdge], + name: str | None, + ) -> None: + """ + Init the NgffSequence object. + + Parameters + ---------- + transformations + The transformations which compose the sequence. + """ + if len(transformations) == 0: + raise ValueError("Empty transformation list") + previous_transf = transformations[0] + for current_transf in transformations[1:]: + if previous_transf.output != current_transf.input: + raise ValueError(f"Mismatched input/output from {previous_transf} to {current_transf}") + previous_transf = current_transf + self.transformations = transformations + super().__init__( + input=transformations[0].input, + output=transformations[-1].output, + name=name, + ) + + def __repr__(self) -> str: + from textwrap import indent + + out = super().__repr__() + " [\n" + for t in self.transformations: + out += indent(repr(t), prefix=" ") + "\n" + out += "]" + return out + + def inverse(self) -> SequenceEdge: + return SequenceEdge( + [t.inverse() for t in reversed(self.transformations)], name=self.name and f"{self.name}__inverse" + ) + + def to_affine(self) -> AffineEdge: + composed = self.transformations[0].to_affine().affine + for t in self.transformations[1:]: + a = t.to_affine() + composed = a.affine @ composed + return AffineEdge.from_affine_matrix( + affine_matrix=composed, input=self.input, output=self.output, name=self.name and f"{self.name}__affine" + ) + + def transform_points(self, points: ArrayLike) -> ArrayLike: + return self.to_affine().transform_points(points) # FIXME + + def to_model(self) -> ozm06trans.Sequence: + return ozm06trans.Sequence( + name=self.name, + input=self.input.to_model_cs_ident(), + output=self.output.to_model_cs_ident(), + transformations=tuple(t.to_model() for t in self.transformations), + ) + + +class ByDimensionEdge(BaseTransfEdge): + """The ByDimension transformation from the NGFF specification.""" + + transformations: Final[Sequence[BaseTransfEdge]] + + def __init__( + self, + name: str | None, + *, + transformations: Sequence[BaseTransfEdge], + input: CoordSystem, + output: CoordSystem, + ) -> None: + """ + Init the ByDimension object. + + Parameters + ---------- + transformations + A list of transformations, whose set of output coordinate systems partition the output coordinate system of + the ByDimension transformation. + input + The input coordinate system of the transformation. + output + The output coordinate system of the transformation. + """ + # we check that: + # 1. each input from each transformation in self.transformation must appear in the set of input axes + # 2. each output from each transformation in self.transformation must appear at most once in the set of output + # axes + input_axes = input.axes_names + output_axes = output.axes_names + defined_output_axes: set[str] = set() + for t in transformations: + for ax in t.input.axes_names: + if ax not in input_axes: + raise ValueError(f"By dimension axis {ax} not in {input_axes}") + for ax in t.output.axes_names: + if ax not in output_axes: + raise ValueError(f"Axis {ax} not in output axes {output_axes}") + if ax in defined_output_axes: + raise ValueError(f"Output axis {ax} is defined more than once") + defined_output_axes.add(ax) + if len(output_axes) != len(defined_output_axes): + raise ValueError("Not all outputs are mapped") + + self.transformations = tuple(transformations) + super().__init__(input=input, output=output, name=name) + + def __repr__(self) -> str: + from textwrap import indent + + out = super().__repr__() + " [\n" + for t in self.transformations: + out += indent(repr(t), prefix=" ") + "\n" + out += "]" + return out + + def inverse(self) -> BaseTransfEdge: + inverse_transformations = [t.inverse() for t in self.transformations] + return ByDimensionEdge( + transformations=inverse_transformations, + input=self.output, + output=self.input, + name=self.name and f"{self.name}__inverse", + ) + + def transform_points(self, points: ArrayLike) -> ArrayLike: + input_axes = self.input.axes_names + output_axes = self.output.axes_names + self._validate_transform_points_shapes(points) + output_columns: dict[str, ArrayLike] = {} + for t in self.transformations: + input_columns = [points[:, input_axes.index(ax)] for ax in t.input.axes_names] + input_columns_stacked: ArrayLike = np.stack(input_columns, axis=1) + output_columns_t = t.transform_points(input_columns_stacked) + for ax, col in zip(t.output.axes_names, output_columns_t.T, strict=True): + output_columns[ax] = col # type: ignore[assignment] + output: ArrayLike = np.stack([output_columns[ax] for ax in output_axes], axis=1) + return output + + def to_affine(self) -> AffineEdge: + input_axes = self.input.axes_names + output_axes = self.output.axes_names + m = np.zeros((len(output_axes) + 1, len(input_axes) + 1)) + m[-1, -1] = 1 + for t in self.transformations: + t_affine = t.to_affine() + target_output_indices = [output_axes.index(ax) for ax in t.output.axes_names if ax in output_axes] + source_output_indices = [t.output.axes_names.index(ax) for ax in t.output.axes_names] + target_input_indices = [input_axes.index(ax) for ax in t.input.axes_names] + [-1] + m[np.ix_(target_output_indices, target_input_indices)] = t_affine.affine[source_output_indices, :] + return AffineEdge.from_affine_matrix( + affine_matrix=m, input=self.input, output=self.output, name=self.name and f"{self.name}__affine" + ) + + def to_model(self) -> ozm06trans.ByDimension: + by_dim_transfs: list[ozm06trans.ByDimensionTransform] = [] + for t in self.transformations: + input_axes = tuple(self.input.axes_names.index(ax_name) for ax_name in t.input.axes_names) + output_axes = tuple(self.output.axes_names.index(ax_name) for ax_name in t.output.axes_names) + by_dim_transfs.append( + ozm06trans.ByDimensionTransform( + input_axes=input_axes, + output_axes=output_axes, + transformation=t.to_model(), + ) + ) + return ozm06trans.ByDimension( + name=self.name, + input=self.input.to_model_cs_ident(), + output=self.output.to_model_cs_ident(), + transformations=tuple(by_dim_transfs), + ) + + +class CsGen: + """A coordinate system generator + + Use it to create coordinate systems on the fly while avoiding + repeating names + """ + + def __init__(self, base_name: str): + self._base_name = base_name + self._cs_count: int = 0 + super().__init__() + + def generate(self, *, num_axes: int) -> CoordSystem: + out = CoordSystem( + name=f"{self._base_name}{self._cs_count}", + axes=[ + Axis( + name=f"axis_{ax_idx}", + type="space", # FIXME + ) + for ax_idx in range(num_axes) + ], + virtual=True, + ) + self._cs_count += 1 + return out + + def generate_like(self, other: CoordSystem) -> CoordSystem: + out = CoordSystem( + name=f"{self._base_name}{self._cs_count}", + axes=[ + Axis( + name=axis.name, + type=axis.type, + unit=axis.unit, + long_name=axis.long_name, + ) + for axis in other.axes + ], + virtual=True, + ) + self._cs_count += 1 + return out + + +def parse_identity( + model: ozm06trans.Identity, + *, + input: CoordSystem, + out: CoordSystem | CsGen, +) -> IdentityEdge: + output = out.generate_like(input) if isinstance(out, CsGen) else out + return IdentityEdge(name=model.name, input=input, output=output) + + +def parse_translation( + model: ozm06trans.Translation, + *, + input: CoordSystem, + out: CoordSystem | CsGen, +) -> TranslationEdge: + output = out.generate_like(input) if isinstance(out, CsGen) else out + return TranslationEdge( + translation=np.asarray(model.translation, dtype=float), + input=input, + output=output, + name=input.name, + ) + + +def parse_scale( + model: ozm06trans.Scale, + *, + input: CoordSystem, + out: CoordSystem | CsGen, +) -> ScaleEdge: + output = out.generate_like(input) if isinstance(out, CsGen) else out + return ScaleEdge( + scale=np.asarray(model.scale, dtype=float), + input=input, + output=output, + name=model.name, + ) + + +def parse_map_axis( + model: ozm06trans.MapAxis, + *, + input: CoordSystem, + out: CoordSystem | CsGen, +) -> MapAxisEdge: + output = out.generate(num_axes=len(model.mapAxis)) if isinstance(out, CsGen) else out + return MapAxisEdge( + input=input, + output=output, + name=model.name, + output_to_input={ # FIXME: double check this. Feels like we depend a lot on order + output.axes[output_axis].name: input.axes[input_axis].name + for output_axis, input_axis in enumerate(model.mapAxis) + }, + ) + + +def parse_affine( + model: ozm06trans.Affine, + *, + input: CoordSystem, + output: CoordSystem | CsGen, +) -> AffineEdge: + num_output_axes = len(model.affine_matrix) # spec doesn't save last row + output = output.generate(num_axes=num_output_axes) if isinstance(output, CsGen) else output + affine_array = np.asarray(model.affine_matrix, dtype=float) + return AffineEdge( + name=model.name, linear=affine_array[:, :-1], translation=affine_array[:, -1], input=input, output=output + ) + + +def parse_rotation( + model: ozm06trans.Rotation, + *, + input: CoordSystem, + out: CoordSystem | CsGen, +) -> RotationEdge: + num_output_axes = len(model.rotation_matrix) + output = out.generate(num_axes=num_output_axes) if isinstance(out, CsGen) else out + return RotationEdge( + name=model.name, + linear_matrix=np.asarray(model.rotation_matrix, dtype=float), + input=input, + output=output, + ) + + +def parse_sequence( + model: ozm06trans.Sequence, + *, + input: CoordSystem, + output: CoordSystem | CsGen, +) -> SequenceEdge: + parsed_inners: list[BaseTransfEdge] = [] + + base_name = "intermediate" + ("" if not model.name else f"_for_{model.name}") + cs_gen: CsGen = output if isinstance(output, CsGen) else CsGen(base_name=base_name) + parsed = parse_ngff_transf( + input=input, + model=model.transformations[0], + output=cs_gen if len(model.transformations) > 1 else output, + ) + parsed_inners.append(parsed) + + for t in model.transformations[1:-1]: + parsed = parse_ngff_transf(input=parsed.output, output=cs_gen, model=t) + parsed_inners.append(parsed) + + if len(model.transformations) > 1: + parsed = parse_ngff_transf(input=parsed.output, output=output, model=model.transformations[-1]) + parsed_inners.append(parsed) + + return SequenceEdge(name=model.name, transformations=parsed_inners) + + +def parse_by_dimension( + model: ozm06trans.ByDimension, + *, + input: CoordSystem, + output: CoordSystem | CsGen, +) -> ByDimensionEdge: + if not isinstance(output, CoordSystem): + max_out_idx = max(ax_idx for t in model.transformations for ax_idx in t.output_axes) + output = output.generate(num_axes=max_out_idx + 1) + + piecewise_transforms: list[BaseTransfEdge] = [] + for t in model.transformations: + inp_axes = [input.axes[i] for i in t.input_axes] + partial_input = CoordSystem( + axes=inp_axes, + name=f"{input.name}_{','.join(ax.name for ax in inp_axes)}", + virtual=True, + ) + + out_axes = [output.axes[i] for i in t.output_axes] + partial_out = CoordSystem( + axes=out_axes, + name=f"{output.name}_{','.join(ax.name for ax in inp_axes)}", + virtual=True, + ) + + parsed_t = parse_ngff_transf(model=t.transformation, input=partial_input, output=partial_out) + piecewise_transforms.append(parsed_t) + + return ByDimensionEdge( + input=input, + output=output, + name=model.name, + transformations=piecewise_transforms, + ) + + +def parse_ngff_transf( + input: CoordSystem, + model: ozm06trans.AnyTransform, + output: CoordSystem | CsGen, +) -> BaseTransfEdge: + if isinstance(model, ozm06trans.Identity): + return parse_identity(model, input=input, out=output) + elif isinstance(model, ozm06trans.Translation): + return parse_translation(model, input=input, out=output) + elif isinstance(model, ozm06trans.Scale): + return parse_scale(model, input=input, out=output) + elif isinstance(model, ozm06trans.MapAxis): + return parse_map_axis(model, input=input, out=output) + elif isinstance(model, ozm06trans.Affine): + return parse_affine(model, input=input, output=output) + elif isinstance(model, ozm06trans.Rotation): + return parse_rotation(model, input=input, out=output) + elif isinstance(model, ozm06trans.Sequence): + return parse_sequence(model, input=input, output=output) + elif isinstance(model, ozm06trans.ByDimension): + return parse_by_dimension(model, input=input, output=output) + else: + raise NotImplementedError(f"Unsupported transformation: {model.type}") diff --git a/src/spatialdata/transformations/graph/vert.py b/src/spatialdata/transformations/graph/vert.py new file mode 100644 index 000000000..8819d04e9 --- /dev/null +++ b/src/spatialdata/transformations/graph/vert.py @@ -0,0 +1,201 @@ +from __future__ import annotations + +from collections.abc import Sequence +from typing import Final, Literal + +import ome_zarr.classes.image as ozi +import pydantic as pyd + + +class AxisParsingException(Exception): + pass + + +class Axis: + """ + Representation of an axis, following the NGFF specification. + + Attributes + ---------- + name + name of the axis. + type + type of the axis. Should be in ["channel", "space"]. + unit + unit of the axis. For a set of valid options see https://ngff.openmicroscopy.org/ + long_name: + a longer, human-friendly name for this axis + """ + + name: Final[str] + type: Final[Literal["space", "channel"]] + unit: Final[str | None] + long_name: Final[str | None] + + class LegacyModel(pyd.BaseModel): + name: Literal["x", "y", "z", "c"] + type: Literal["space", "channel"] + + def __init__( + self, *, name: str, type: Literal["space", "channel"], unit: str | None = None, long_name: str | None = None + ): + self.name = name + self.type = type + self.unit = unit + self.long_name = long_name + + def cloned_with(self, *, unit: str | None) -> Axis: + return Axis(name=self.name, type=self.type, unit=unit or self.unit, long_name=self.long_name) + + def __hash__(self) -> int: + return hash((self.name, self.type, self.unit, self.long_name)) + + def __repr__(self) -> str: + return f"NgffAxis(name={self.name}, type={self.type})" + + def __eq__(self, value: object, /) -> bool: + if not isinstance(value, Axis): + return False + return ( + self.name == value.name + and self.type == value.type + and self.unit == value.unit + and self.long_name == value.long_name + ) + + @classmethod + def try_from_model(cls, model: ozi.Axis) -> Axis: + name = model.name + if name is None: + raise AxisParsingException("Axis doesn't have a name") + if model.type != "channel" and model.type != "space": + raise AxisParsingException(f"Can't handle axis of type {model.type}") + if not isinstance(model.unit, str): + raise AxisParsingException("Can't handle axis unit") + return Axis( + name=name, + type=model.type, + unit=model.unit, + long_name=model.longName, + ) + + @classmethod + def try_from_dict(cls, d: pyd.JsonValue) -> Axis: + model = ozi.Axis.model_validate(d) + return Axis.try_from_model(model) + + def to_model(self) -> ozi.Axis: + return ozi.Axis( + discrete=False, + longName=self.long_name, + name=self.name, + type=self.type, + unit=self.unit, + ) + + +class CoordSystemParsingException(Exception): + pass + + +class CoordSystem: + """ + Representation of a coordinate system, following the NGFF specification. + + Parameters + ---------- + name + name of the coordinate system + axes + names of the axes of the coordinate system + """ + + name: Final[str] + axes: Final[tuple[Axis, ...]] + + virtual: Final[bool] + """A virtual coordinate system exists as an intermediate step between + non-virtual coordinate systems and is usually ignored during serialization""" + + class LegacyAxes: + pass + + def __init__(self, name: str, axes: Sequence[Axis], virtual: bool = False): + self.name = name + self.axes = tuple(axes) + self.virtual = virtual + if len(self.axes) != len({axis.name for axis in self.axes}): + raise ValueError("Axes names must be unique") + + def __repr__(self) -> str: + return f"{self.__class__.__name__}({self.name!r}, {self.axes})" + + def __hash__(self) -> int: + return hash((self.name, self.axes, self.virtual)) + + @classmethod + def try_from_model(cls, model: ozi.CoordinateSystem) -> CoordSystem: + axes: list[Axis] = [] + for axis in model.axes: + if isinstance(parsed := Axis.try_from_model(axis), Exception): + raise CoordSystemParsingException(parsed) # FIXME + axes.append(parsed) + return CoordSystem( + name=model.name, + axes=axes, + ) + + @classmethod + def try_from_model_or_default[T](cls, model: ozi.CoordinateSystem | None, *, default: T) -> CoordSystem | T: + if model is not None: + return CoordSystem.try_from_model(model) + return default + + def to_model(self) -> ozi.CoordinateSystem | None: + if self.virtual: + return None + return ozi.CoordinateSystem( + name=self.name, + axes=tuple(ax.to_model() for ax in self.axes), + ) + + def to_model_cs_ident(self) -> ozi.CoordinateSystemIdentifier | None: + input = self.to_model() + if input is None: + return None + return ozi.CoordinateSystemIdentifier(name=input.name) + + @property + def num_axes(self) -> int: + return len(self.axes) + + @property + def axes_names(self) -> tuple[str, ...]: + """Get axes' names""" + return tuple([ax.name for ax in self.axes]) + + @property + def axes_types(self) -> tuple[str, ...]: + """Get axes' types""" + return tuple([ax.type for ax in self.axes]) + + def has_axis(self, name: str) -> bool: + """ + Check the coordinate system has an axis of the given name. + + Parameters + ---------- + name + name of the axis. + """ + return any(axis.name == name for axis in self.axes) + + def get_axis(self, name: str) -> Axis | None: + """Get the axis by name""" + for axis in self.axes: + if axis.name == name: + return axis + return None + + def get_spatial_axes(self) -> Sequence[Axis]: + return [axis for axis in self.axes if axis.type == "space"] From 2379e27145c806d4f1a5bb38b392863aa7c18e13 Mon Sep 17 00:00:00 2001 From: Tomaz Vieira Date: Wed, 2 Sep 2026 16:11:38 +0000 Subject: [PATCH 03/26] addresses some PR comments, simplifies MapAxis --- src/spatialdata/transformations/graph/edge.py | 162 ++++++++++-------- src/spatialdata/transformations/graph/vert.py | 19 +- 2 files changed, 94 insertions(+), 87 deletions(-) diff --git a/src/spatialdata/transformations/graph/edge.py b/src/spatialdata/transformations/graph/edge.py index b01a11e97..42396175f 100644 --- a/src/spatialdata/transformations/graph/edge.py +++ b/src/spatialdata/transformations/graph/edge.py @@ -9,7 +9,6 @@ import numpy as np import ome_zarr_models.v06.coordinate_transforms as ozm06trans import pydantic as pyd -import xarray as xr from spatialdata._types import ArrayLike from spatialdata.transformations.graph.vert import Axis, CoordSystem @@ -48,7 +47,7 @@ def __repr__(self) -> str: return f"{type(self).__name__} ({domain} -> {codomain})" @abstractmethod - def inverse(self) -> BaseTransfEdge: + def inverse(self, name: str | None = None) -> BaseTransfEdge: """Return the inverse of the transformation.""" @abstractmethod @@ -63,23 +62,23 @@ def transform_points(self, points: ArrayLike) -> ArrayLike: """ @abstractmethod - def to_affine(self) -> AffineEdge: + def to_affine(self, name: str | None = None) -> AffineEdge: """Convert the transformation to an affine transformation, whenever the conversion can be made.""" - def _validate_transform_points_shapes(self, points: xr.DataArray | xr.DataTree | ArrayLike) -> None: + def _validate_transform_points_shapes(self, points: ArrayLike) -> None: """ - Validate if the shape of the points (coordinats to be transformed) are consistent with the input size of the + Validate if the shape of the points (coordinates to be transformed) are consistent with the input size of the transformation. """ input_size = len(self.input.axes) if len(points.shape) != 2 or points.shape[1] != input_size: raise ValueError( f"points must be a tensor of shape (n, d), where n is the number of points and d is the " - f"the number of spatial dimensions. Points shape: {points.shape}, input size: {input_size}" + f"the number of dimensions. Points shape: {points.shape}, input size: {input_size}" ) # order of the composition: self is applied first, then the transformation passed as argument - def compose_with(self, transformation: BaseTransfEdge) -> BaseTransfEdge: + def compose_with(self, transformation: BaseTransfEdge, name: str | None) -> BaseTransfEdge: """ Compose the transfomation object with another transformation @@ -96,7 +95,7 @@ def compose_with(self, transformation: BaseTransfEdge) -> BaseTransfEdge: ------- Self is applied first, then the transformation passed as argument. """ - return SequenceEdge([self, transformation], name=None) # FIXME: no name? + return SequenceEdge([self, transformation], name=name) @abstractmethod def to_model(self) -> ozm06trans.AnyTransform: @@ -112,8 +111,8 @@ class AffineEdge(BaseTransfEdge): def __init__( self, - name: str | None, *, + name: str | None = None, linear: ArrayLike, translation: ArrayLike | None = None, input: CoordSystem, @@ -169,7 +168,7 @@ def from_affine_matrix( input: CoordSystem, output: CoordSystem, ) -> AffineEdge: - """Creates an AffineEdge from a raw affine matrix + """Creates an AffineEdge from a raw affine matrix in homogenous coordinates Parameters ---------- @@ -193,14 +192,14 @@ def from_affine_matrix( name=name, ) - def inverse(self) -> BaseTransfEdge: + def inverse(self, name: str | None = None) -> BaseTransfEdge: inv = np.linalg.inv(self.affine) return AffineEdge( linear=inv[:-1, :-1], translation=inv[-1, :-1], input=self.output, output=self.input, - name=self.name and f"{self.name}__affine", + name=name, ) def transform_points(self, points: ArrayLike) -> ArrayLike: @@ -211,8 +210,10 @@ def transform_points(self, points: ArrayLike) -> ArrayLike: assert isinstance(res, np.ndarray) return res - def to_affine(self) -> AffineEdge: - return self + def to_affine(self, name: str | None = None) -> AffineEdge: + return AffineEdge( + input=self.input, output=self.output, linear=self.linear, translation=self.translation, name=name + ) def to_model(self) -> ozm06trans.Affine: return ozm06trans.Affine( @@ -247,19 +248,19 @@ def __init__( raise ValueError("Input and output must have the same number of dimensions") super().__init__(input=input, output=output, name=name) - def inverse(self) -> BaseTransfEdge: - return IdentityEdge(input=self.output, output=self.input, name=self.name and f"{self.name}__inverse") + def inverse(self, name: str | None = None) -> BaseTransfEdge: + return IdentityEdge(input=self.output, output=self.input, name=name) def transform_points(self, points: ArrayLike) -> ArrayLike: self._validate_transform_points_shapes(points) return points - def to_affine(self) -> AffineEdge: + def to_affine(self, name: str | None = None) -> AffineEdge: return AffineEdge( linear=np.eye(self.input.num_axes), input=self.input, output=self.output, - name=self.name and f"{self.name}__affine", + name=name, ) def to_model(self) -> ozm06trans.Identity: @@ -270,6 +271,13 @@ def to_model(self) -> ozm06trans.Identity: ) +class UnmappableCoordSystemsError(Exception): + def __init__(self, input: CoordSystem, output: CoordSystem) -> None: + self.input = input + self.output = output + super().__init__("Output axes can't be mapped to input axes") + + class MapAxisEdge(BaseTransfEdge): """The MapAxis transformation from the NGFF specification.""" @@ -277,7 +285,6 @@ def __init__( self, name: str | None, *, - output_to_input: dict[str, str], input: CoordSystem, output: CoordSystem, ) -> None: @@ -287,67 +294,56 @@ def __init__( ---------- name A human readable name for this transformation - output_to_input - A dictionary mapping the output axes (keys) to the input axes (values). input Input coordinate system of the transformation. output - Output coordinate system of the transformation. + Output coordinate system of the transformation, whose axes + must be a shuffling of `input` """ - for out_ax, inp_ax in output_to_input.items(): - if not input.has_axis(inp_ax): - raise ValueError(f"input has no axis named {inp_ax}") - if not output.has_axis(out_ax): - raise ValueError(f"output has no axis named {out_ax}") - if not (len(output_to_input) == output.num_axes == input.num_axes): - raise ValueError("input_to_output, input and output must have the same number of axes entries") - if len(set(output_to_input.values())) != len(output_to_input): - raise ValueError("input_to_output must map unique inputs to unique outputs") - - self.output_to_input = output_to_input + + if set(input.axes) != set(output.axes): + raise UnmappableCoordSystemsError(input=input, output=output) super().__init__(input=input, output=output, name=name) def __repr__(self) -> str: s = super().__repr__() + "\n" - s += "\n".join(f" {out} <- {inp}\n" for out, inp in self.output_to_input.items()) + s += "\n".join( + f" {out.name} <- {inp.name}\n" for out, inp in zip(self.output.axes, self.input.axes, strict=True) + ) return s - def inverse(self) -> BaseTransfEdge: + def inverse(self, name: str | None = None) -> BaseTransfEdge: return MapAxisEdge( - output_to_input={v: k for k, v in self.output_to_input.items()}, input=self.output, output=self.input, - name=self.name and f"{self.name}__inverse", + name=name, ) def transform_points(self, points: ArrayLike) -> ArrayLike: - input_axes = self.input.axes_names - output_axes = self.output.axes_names self._validate_transform_points_shapes(points) - new_indices = [input_axes.index(self.output_to_input[ax]) for ax in output_axes] + new_indices = [self.input.axes.index(out_ax.name) for out_ax in self.output.axes] mapped = points[:, new_indices] assert isinstance(mapped, np.ndarray) return mapped - def to_affine(self) -> AffineEdge: + def to_affine(self, name: str | None = None) -> AffineEdge: input_axes = self.input.axes_names output_axes = self.output.axes_names linear: ArrayLike = np.zeros((len(output_axes), len(input_axes)), dtype=float) for i, des_axis in enumerate(output_axes): for j, src_axis in enumerate(input_axes): - if src_axis == self.output_to_input[des_axis]: + if src_axis == des_axis: linear[i, j] = 1 affine = AffineEdge( - linear=linear, input=self.input, output=self.output, name=self.name and f"{self.name}__affine" + linear=linear, + input=self.input, + output=self.output, + name=name, ) return affine def to_model(self) -> ozm06trans.MapAxis: - mapAxis: list[int] = [] - for out_ax in self.output.axes_names: - in_ax = self.output_to_input[out_ax] - in_idx = self.input.axes_names.index(in_ax) - mapAxis.append(in_idx) + mapAxis: list[int] = [self.input.axes.index(out_ax) for out_ax in self.output.axes] return ozm06trans.MapAxis( name=self.name, mapAxis=tuple(mapAxis), @@ -388,25 +384,25 @@ def __init__( def __repr__(self) -> str: return super().__repr__() + str(self.translation) - def inverse(self) -> BaseTransfEdge: + def inverse(self, name: str | None = None) -> BaseTransfEdge: return TranslationEdge( translation=-self.translation, input=self.output, output=self.input, - name=self.name and f"{self.name}__inverse", + name=name, ) def transform_points(self, points: ArrayLike) -> ArrayLike: self._validate_transform_points_shapes(points) return points + self.translation - def to_affine(self) -> AffineEdge: + def to_affine(self, name: str | None = None) -> AffineEdge: return AffineEdge( linear=np.identity(self.input.num_axes), translation=self.translation, input=self.input, output=self.output, - name=self.name and f"{self.name}__affine", + name=name, ) def to_model(self) -> ozm06trans.Translation: @@ -450,21 +446,27 @@ def __init__( def __repr__(self) -> str: return super().__repr__() + str(self.scale) - def inverse(self) -> ScaleEdge: + def inverse(self, name: str | None = None) -> ScaleEdge: if any(s == 0 for s in self.scale): raise ValueError(f"Scaling {self} is not invertible") new_scale = 1 / self.scale return ScaleEdge( - scale=new_scale, input=self.output, output=self.input, name=self.name and f"{self.name}__inverse" + scale=new_scale, + input=self.output, + output=self.input, + name=name, ) def transform_points(self, points: ArrayLike) -> ArrayLike: self._validate_transform_points_shapes(points) return points * self.scale - def to_affine(self) -> AffineEdge: + def to_affine(self, name: str | None = None) -> AffineEdge: return AffineEdge( - linear=np.diag(self.scale), input=self.input, output=self.output, name=self.name and f"{self.name}__affine" + linear=np.diag(self.scale), + input=self.input, + output=self.output, + name=name, ) def to_model(self) -> ozm06trans.Scale: @@ -516,12 +518,12 @@ def __repr__(self) -> str: s += "\n".join(str(row) for row in self.rotation) return s - def inverse(self) -> BaseTransfEdge: + def inverse(self, name: str | None = None) -> BaseTransfEdge: return RotationEdge( linear_matrix=self.rotation.T, input=self.output, output=self.input, - name=self.name and f"{self.name}__inverse", + name=name, ) def transform_points(self, points: ArrayLike) -> ArrayLike: @@ -530,9 +532,12 @@ def transform_points(self, points: ArrayLike) -> ArrayLike: assert isinstance(res, np.ndarray) return res - def to_affine(self) -> AffineEdge: + def to_affine(self, name: str | None = None) -> AffineEdge: return AffineEdge( - linear=self.rotation, input=self.input, output=self.output, name=self.name and f"{self.name}__affine" + linear=self.rotation, + input=self.input, + output=self.output, + name=name, ) def to_model(self) -> ozm06trans.Rotation: @@ -583,18 +588,22 @@ def __repr__(self) -> str: out += "]" return out - def inverse(self) -> SequenceEdge: + def inverse(self, name: str | None = None) -> SequenceEdge: return SequenceEdge( - [t.inverse() for t in reversed(self.transformations)], name=self.name and f"{self.name}__inverse" + [t.inverse() for t in reversed(self.transformations)], + name=name, ) - def to_affine(self) -> AffineEdge: + def to_affine(self, name: str | None = None) -> AffineEdge: composed = self.transformations[0].to_affine().affine for t in self.transformations[1:]: a = t.to_affine() composed = a.affine @ composed return AffineEdge.from_affine_matrix( - affine_matrix=composed, input=self.input, output=self.output, name=self.name and f"{self.name}__affine" + affine_matrix=composed, + input=self.input, + output=self.output, + name=name, ) def transform_points(self, points: ArrayLike) -> ArrayLike: @@ -667,13 +676,13 @@ def __repr__(self) -> str: out += "]" return out - def inverse(self) -> BaseTransfEdge: + def inverse(self, name: str | None = None) -> BaseTransfEdge: inverse_transformations = [t.inverse() for t in self.transformations] return ByDimensionEdge( transformations=inverse_transformations, input=self.output, output=self.input, - name=self.name and f"{self.name}__inverse", + name=name, ) def transform_points(self, points: ArrayLike) -> ArrayLike: @@ -690,7 +699,7 @@ def transform_points(self, points: ArrayLike) -> ArrayLike: output: ArrayLike = np.stack([output_columns[ax] for ax in output_axes], axis=1) return output - def to_affine(self) -> AffineEdge: + def to_affine(self, name: str | None = None) -> AffineEdge: input_axes = self.input.axes_names output_axes = self.output.axes_names m = np.zeros((len(output_axes) + 1, len(input_axes) + 1)) @@ -702,7 +711,10 @@ def to_affine(self) -> AffineEdge: target_input_indices = [input_axes.index(ax) for ax in t.input.axes_names] + [-1] m[np.ix_(target_output_indices, target_input_indices)] = t_affine.affine[source_output_indices, :] return AffineEdge.from_affine_matrix( - affine_matrix=m, input=self.input, output=self.output, name=self.name and f"{self.name}__affine" + affine_matrix=m, + input=self.input, + output=self.output, + name=name, ) def to_model(self) -> ozm06trans.ByDimension: @@ -816,15 +828,19 @@ def parse_map_axis( input: CoordSystem, out: CoordSystem | CsGen, ) -> MapAxisEdge: - output = out.generate(num_axes=len(model.mapAxis)) if isinstance(out, CsGen) else out + if isinstance(out, CoordSystem): + output = out + else: + dummy_cs = out.generate(num_axes=len(model.mapAxis)) + output = CoordSystem( + name=dummy_cs.name, + axes=[input.axes[i] for i in model.mapAxis], + virtual=True, + ) return MapAxisEdge( input=input, output=output, name=model.name, - output_to_input={ # FIXME: double check this. Feels like we depend a lot on order - output.axes[output_axis].name: input.axes[input_axis].name - for output_axis, input_axis in enumerate(model.mapAxis) - }, ) diff --git a/src/spatialdata/transformations/graph/vert.py b/src/spatialdata/transformations/graph/vert.py index 8819d04e9..128c50706 100644 --- a/src/spatialdata/transformations/graph/vert.py +++ b/src/spatialdata/transformations/graph/vert.py @@ -4,7 +4,7 @@ from typing import Final, Literal import ome_zarr.classes.image as ozi -import pydantic as pyd +import ome_zarr_models.v06.coordinate_transforms as ozm06ct class AxisParsingException(Exception): @@ -32,10 +32,6 @@ class Axis: unit: Final[str | None] long_name: Final[str | None] - class LegacyModel(pyd.BaseModel): - name: Literal["x", "y", "z", "c"] - type: Literal["space", "channel"] - def __init__( self, *, name: str, type: Literal["space", "channel"], unit: str | None = None, long_name: str | None = None ): @@ -64,13 +60,13 @@ def __eq__(self, value: object, /) -> bool: ) @classmethod - def try_from_model(cls, model: ozi.Axis) -> Axis: + def try_from_model(cls, model: ozm06ct.Axis) -> Axis: name = model.name if name is None: raise AxisParsingException("Axis doesn't have a name") if model.type != "channel" and model.type != "space": raise AxisParsingException(f"Can't handle axis of type {model.type}") - if not isinstance(model.unit, str): + if not isinstance(model.unit, (str, type(None))): raise AxisParsingException("Can't handle axis unit") return Axis( name=name, @@ -79,13 +75,8 @@ def try_from_model(cls, model: ozi.Axis) -> Axis: long_name=model.longName, ) - @classmethod - def try_from_dict(cls, d: pyd.JsonValue) -> Axis: - model = ozi.Axis.model_validate(d) - return Axis.try_from_model(model) - - def to_model(self) -> ozi.Axis: - return ozi.Axis( + def to_model(self) -> ozm06ct.Axis: + return ozm06ct.Axis( discrete=False, longName=self.long_name, name=self.name, From 0d4140bada6274f3a943df2d5b71e9276dc8f2be Mon Sep 17 00:00:00 2001 From: Tomaz Vieira Date: Thu, 3 Sep 2026 16:32:06 +0000 Subject: [PATCH 04/26] Adds ProjectAxisEdge, makes inverse return optional --- src/spatialdata/transformations/graph/edge.py | 160 +++++++++++++++--- 1 file changed, 132 insertions(+), 28 deletions(-) diff --git a/src/spatialdata/transformations/graph/edge.py b/src/spatialdata/transformations/graph/edge.py index 42396175f..74664cff7 100644 --- a/src/spatialdata/transformations/graph/edge.py +++ b/src/spatialdata/transformations/graph/edge.py @@ -31,8 +31,8 @@ class BaseTransfEdge(ABC): def __init__( self, - name: str | None, *, + name: str | None = None, input: CoordSystem, output: CoordSystem, ) -> None: @@ -47,8 +47,8 @@ def __repr__(self) -> str: return f"{type(self).__name__} ({domain} -> {codomain})" @abstractmethod - def inverse(self, name: str | None = None) -> BaseTransfEdge: - """Return the inverse of the transformation.""" + def inverse(self, name: str | None = None) -> BaseTransfEdge | None: + """Return the inverse of the transformation if it exists""" @abstractmethod def transform_points(self, points: ArrayLike) -> ArrayLike: @@ -192,8 +192,26 @@ def from_affine_matrix( name=name, ) - def inverse(self, name: str | None = None) -> BaseTransfEdge: - inv = np.linalg.inv(self.affine) + @classmethod + def mapping(cls, input: CoordSystem, output: CoordSystem, name: str | None = None) -> AffineEdge: + linear: ArrayLike = np.zeros((output.num_axes, input.num_axes), dtype=float) + for i, des_axis in enumerate(output.axes): + for j, src_axis in enumerate(input.axes): + if src_axis.name == des_axis.name: # FIXME: compare the entire axis? + linear[i, j] = 1 + return AffineEdge( + linear=linear, + input=input, + output=output, + name=name, + ) + + def inverse(self, name: str | None = None) -> BaseTransfEdge | None: + try: + # FIXME: I think there are more efficient/precise ways to invert a matrix + inv = np.linalg.inv(self.affine) + except np.linalg.LinAlgError: + return None return AffineEdge( linear=inv[:-1, :-1], translation=inv[-1, :-1], @@ -327,20 +345,7 @@ def transform_points(self, points: ArrayLike) -> ArrayLike: return mapped def to_affine(self, name: str | None = None) -> AffineEdge: - input_axes = self.input.axes_names - output_axes = self.output.axes_names - linear: ArrayLike = np.zeros((len(output_axes), len(input_axes)), dtype=float) - for i, des_axis in enumerate(output_axes): - for j, src_axis in enumerate(input_axes): - if src_axis == des_axis: - linear[i, j] = 1 - affine = AffineEdge( - linear=linear, - input=self.input, - output=self.output, - name=name, - ) - return affine + return AffineEdge.mapping(input=self.input, output=self.output, name=name) def to_model(self) -> ozm06trans.MapAxis: mapAxis: list[int] = [self.input.axes.index(out_ax) for out_ax in self.output.axes] @@ -352,6 +357,97 @@ def to_model(self) -> ozm06trans.MapAxis: ) +class AxisNotInCoordSystemError(Exception): + def __init__(self, axis: Axis, cs: CoordSystem) -> None: + self.axis = axis + self.cs = cs + super().__init__(f"Axis {axis.name} is not in coordinate system {cs.name}") + + +class ProjectAxisEdge(BaseTransfEdge): + dropped_inputs: Final[set[Axis]] + created_outputs: Final[set[Axis]] + + def __init__( + self, + *, + name: str | None = None, + input: CoordSystem, + output: CoordSystem, + dropped_inputs: set[Axis], + created_outputs: set[Axis], + ) -> None: + for axis in dropped_inputs: + if axis not in input.axes: + raise AxisNotInCoordSystemError(axis=axis, cs=input) + for axis in dropped_inputs: + if axis not in output.axes: + raise AxisNotInCoordSystemError(axis=axis, cs=output) + self.dropped_inputs = set(dropped_inputs) + self.created_outputs = set(created_outputs) + super().__init__(name=name, input=input, output=output) + + def to_affine(self, name: str | None = None) -> AffineEdge: + linear = np.zeros((self.output.num_axes, self.input.num_axes), dtype=float) + + for out_idx, out_ax in enumerate(self.output.axes): + if out_ax in self.created_outputs: + continue + for in_idx, in_ax in enumerate(self.input.axes): + if in_ax not in self.dropped_inputs: + linear[out_idx, in_idx] = 1 + + return AffineEdge(name=name, input=self.input, output=self.output, linear=linear) + + def transform_points(self, points: ArrayLike) -> ArrayLike: + return self.to_affine().transform_points(points) + + def inverse(self, name: str | None = None) -> BaseTransfEdge | None: + # FIXME: there may be other cases where this is invertible + if self.input.num_axes != self.output.num_axes: + return None + if len(self.dropped_inputs) > 0: + return None + if len(self.created_outputs) > 0: + return None + return ProjectAxisEdge( + input=self.output, + output=self.input, + dropped_inputs=set(), + created_outputs=set(), + name=name, + ) + + def to_model(self) -> ozm06trans.ProjectAxis: + return ozm06trans.ProjectAxis( + createdOutputs=tuple(self.output.axes.index(co) for co in self.created_outputs) or None, + droppedInputs=tuple(self.input.axes.index(di) for di in self.dropped_inputs) or None, + ) + + +def parse_project_axis( + model: ozm06trans.ProjectAxis, + *, + input: CoordSystem, + out: CoordSystem | CsGen, +) -> ProjectAxisEdge: + if isinstance(out, CoordSystem): + output = out + else: + num_dropped_inputs = len(model.droppedInputs or ()) + num_created_outputs = len(model.droppedInputs or ()) + num_output_axes = input.num_axes - num_dropped_inputs + num_created_outputs + output = out.generate(num_axes=num_output_axes) + + return ProjectAxisEdge( + created_outputs={output.axes[i] for i in model.createdOutputs or ()}, + dropped_inputs={input.axes[i] for i in model.droppedInputs or ()}, + input=input, + output=output, + name=model.name, + ) + + class TranslationEdge(BaseTransfEdge): """The Translation transformation from the NGFF specification.""" @@ -446,9 +542,9 @@ def __init__( def __repr__(self) -> str: return super().__repr__() + str(self.scale) - def inverse(self, name: str | None = None) -> ScaleEdge: + def inverse(self, name: str | None = None) -> ScaleEdge | None: if any(s == 0 for s in self.scale): - raise ValueError(f"Scaling {self} is not invertible") + return None new_scale = 1 / self.scale return ScaleEdge( scale=new_scale, @@ -588,11 +684,14 @@ def __repr__(self) -> str: out += "]" return out - def inverse(self, name: str | None = None) -> SequenceEdge: - return SequenceEdge( - [t.inverse() for t in reversed(self.transformations)], - name=name, - ) + def inverse(self, name: str | None = None) -> SequenceEdge | None: + inverted: list[BaseTransfEdge] = [] + for t in self.transformations: + inv = t.inverse() + if inv is None: + return None + inverted.append(inv) + return SequenceEdge(inverted, name=name) def to_affine(self, name: str | None = None) -> AffineEdge: composed = self.transformations[0].to_affine().affine @@ -676,8 +775,13 @@ def __repr__(self) -> str: out += "]" return out - def inverse(self, name: str | None = None) -> BaseTransfEdge: - inverse_transformations = [t.inverse() for t in self.transformations] + def inverse(self, name: str | None = None) -> BaseTransfEdge | None: + inverse_transformations: list[BaseTransfEdge] = [] + for t in self.transformations: + inv = t.inverse() + if inv is None: + return None + inverse_transformations.append(inv) return ByDimensionEdge( transformations=inverse_transformations, input=self.output, From a0582e8f5fd4f4a6b6c685a821ea95a19e894716 Mon Sep 17 00:00:00 2001 From: Tomaz Vieira Date: Fri, 4 Sep 2026 09:16:37 +0000 Subject: [PATCH 05/26] Moves exceptions to exceptions.py, edge names default to None --- .../transformation_manager/exceptions.py | 15 +++++++ src/spatialdata/transformations/graph/edge.py | 42 +++++-------------- 2 files changed, 25 insertions(+), 32 deletions(-) diff --git a/src/spatialdata/_core/transformation_manager/exceptions.py b/src/spatialdata/_core/transformation_manager/exceptions.py index c8324032a..77ae94837 100644 --- a/src/spatialdata/_core/transformation_manager/exceptions.py +++ b/src/spatialdata/_core/transformation_manager/exceptions.py @@ -1,5 +1,6 @@ from __future__ import annotations +from spatialdata.transformations.graph.vert import Axis, CoordSystem from spatialdata.transformations.ngff.ngff_coordinate_system import NgffCoordinateSystem @@ -227,3 +228,17 @@ class TransformationManagerWarning(UserWarning): """Base warning category for TransformationManager.""" pass + + +class UnmappableCoordSystemsError(Exception): + def __init__(self, input: CoordSystem, output: CoordSystem) -> None: + self.input = input + self.output = output + super().__init__("Output axes can't be mapped to input axes") + + +class AxisNotInCoordSystemError(Exception): + def __init__(self, axis: Axis, cs: CoordSystem) -> None: + self.axis = axis + self.cs = cs + super().__init__(f"Axis {axis.name} is not in coordinate system {cs.name}") diff --git a/src/spatialdata/transformations/graph/edge.py b/src/spatialdata/transformations/graph/edge.py index 74664cff7..9a6051942 100644 --- a/src/spatialdata/transformations/graph/edge.py +++ b/src/spatialdata/transformations/graph/edge.py @@ -8,20 +8,12 @@ import numpy as np import ome_zarr_models.v06.coordinate_transforms as ozm06trans -import pydantic as pyd +from spatialdata._core.transformation_manager.exceptions import AxisNotInCoordSystemError, UnmappableCoordSystemsError from spatialdata._types import ArrayLike from spatialdata.transformations.graph.vert import Axis, CoordSystem -class GarbledInput(Exception): - def __init__(self, message: str, input: pyd.JsonValue) -> None: - import json - - super().__init__(message + "\n" + json.dumps(input, indent=4)) - self.input = input - - class BaseTransfEdge(ABC): """Base class for all the transformations defined by the NGFF specification.""" @@ -95,7 +87,7 @@ def compose_with(self, transformation: BaseTransfEdge, name: str | None) -> Base ------- Self is applied first, then the transformation passed as argument. """ - return SequenceEdge([self, transformation], name=name) + return SequenceEdge(transformations=[self, transformation], name=name) @abstractmethod def to_model(self) -> ozm06trans.AnyTransform: @@ -289,25 +281,17 @@ def to_model(self) -> ozm06trans.Identity: ) -class UnmappableCoordSystemsError(Exception): - def __init__(self, input: CoordSystem, output: CoordSystem) -> None: - self.input = input - self.output = output - super().__init__("Output axes can't be mapped to input axes") - - class MapAxisEdge(BaseTransfEdge): """The MapAxis transformation from the NGFF specification.""" def __init__( self, - name: str | None, *, + name: str | None = None, input: CoordSystem, output: CoordSystem, ) -> None: """ - Init the NgffMapAxis object. Parameters ---------- name @@ -357,13 +341,6 @@ def to_model(self) -> ozm06trans.MapAxis: ) -class AxisNotInCoordSystemError(Exception): - def __init__(self, axis: Axis, cs: CoordSystem) -> None: - self.axis = axis - self.cs = cs - super().__init__(f"Axis {axis.name} is not in coordinate system {cs.name}") - - class ProjectAxisEdge(BaseTransfEdge): dropped_inputs: Final[set[Axis]] created_outputs: Final[set[Axis]] @@ -453,8 +430,8 @@ class TranslationEdge(BaseTransfEdge): def __init__( self, - name: str | None, *, + name: str | None = None, translation: ArrayLike, input: CoordSystem, output: CoordSystem, @@ -515,8 +492,8 @@ class ScaleEdge(BaseTransfEdge): def __init__( self, - name: str | None, *, + name: str | None = None, scale: ArrayLike, input: CoordSystem, output: CoordSystem, @@ -581,8 +558,8 @@ class RotationEdge(BaseTransfEdge): def __init__( self, - name: str | None, *, + name: str | None = None, linear_matrix: ArrayLike, input: CoordSystem, output: CoordSystem, @@ -650,8 +627,9 @@ class SequenceEdge(BaseTransfEdge): def __init__( self, + *, + name: str | None = None, transformations: Sequence[BaseTransfEdge], - name: str | None, ) -> None: """ Init the NgffSequence object. @@ -691,7 +669,7 @@ def inverse(self, name: str | None = None) -> SequenceEdge | None: if inv is None: return None inverted.append(inv) - return SequenceEdge(inverted, name=name) + return SequenceEdge(transformations=inverted, name=name) def to_affine(self, name: str | None = None) -> AffineEdge: composed = self.transformations[0].to_affine().affine @@ -724,8 +702,8 @@ class ByDimensionEdge(BaseTransfEdge): def __init__( self, - name: str | None, *, + name: str | None = None, transformations: Sequence[BaseTransfEdge], input: CoordSystem, output: CoordSystem, From 9d9315b6bbcb71abfc468efc4c6a20782a1d5f78 Mon Sep 17 00:00:00 2001 From: Tomaz Vieira Date: Fri, 4 Sep 2026 11:13:27 +0000 Subject: [PATCH 06/26] Use custom exceptions everywhere, adds Raises to docstrings --- .../transformation_manager/exceptions.py | 49 +++++- src/spatialdata/transformations/graph/edge.py | 162 +++++++++++++----- 2 files changed, 167 insertions(+), 44 deletions(-) diff --git a/src/spatialdata/_core/transformation_manager/exceptions.py b/src/spatialdata/_core/transformation_manager/exceptions.py index 77ae94837..d921beb9d 100644 --- a/src/spatialdata/_core/transformation_manager/exceptions.py +++ b/src/spatialdata/_core/transformation_manager/exceptions.py @@ -1,5 +1,6 @@ from __future__ import annotations +from spatialdata._types import ArrayLike from spatialdata.transformations.graph.vert import Axis, CoordSystem from spatialdata.transformations.ngff.ngff_coordinate_system import NgffCoordinateSystem @@ -230,15 +231,55 @@ class TransformationManagerWarning(UserWarning): pass -class UnmappableCoordSystemsError(Exception): - def __init__(self, input: CoordSystem, output: CoordSystem) -> None: +class IncompatibleCoordSystemsError(Exception): + def __init__(self, input: CoordSystem, output: CoordSystem, message: str | None = None) -> None: self.input = input self.output = output - super().__init__("Output axes can't be mapped to input axes") + super().__init__(message or "Output axes can't be mapped to input axes") -class AxisNotInCoordSystemError(Exception): +class MissingAxisError(Exception): def __init__(self, axis: Axis, cs: CoordSystem) -> None: self.axis = axis self.cs = cs super().__init__(f"Axis {axis.name} is not in coordinate system {cs.name}") + + +class UnexpectedShapeError(Exception): + def __init__( + self, + *, + array_shape: tuple[int, ...], + expected_shape: tuple[int, ...] | str | None = None, + array_name: str | None = None, + ) -> None: + self.array_shape = array_shape + self.expected_shape = expected_shape + message = "Unexpected array shape" + if array_name is not None: + message += f"for '{array_name}'" + message += f": {array_shape}" + if expected_shape is not None: + message += f" instead of {expected_shape}" + super().__init__(message) + + +class NotUnimodularError(Exception): + def __init__(self, matrix: ArrayLike) -> None: + self.matrix = matrix + super().__init__("Matrix does not have det(M) == 1") + + +class EmptyTransformSequenceError(Exception): + def __init__(self) -> None: + super().__init__("Empty sequence of transformations") + + +class AxisRedefinitionError(Exception): + def __init__(self, axis: Axis) -> None: + super().__init__(f"Axis {axis.name} is defined multiple times") + + +class UnmappedAxisError(Exception): + def __init__(self, axis: Axis, cs: CoordSystem) -> None: + super().__init__(f"Axis {axis.name} from coordinate system {cs.name} is not mapped to anything") diff --git a/src/spatialdata/transformations/graph/edge.py b/src/spatialdata/transformations/graph/edge.py index 9a6051942..aba0de7e5 100644 --- a/src/spatialdata/transformations/graph/edge.py +++ b/src/spatialdata/transformations/graph/edge.py @@ -9,7 +9,15 @@ import numpy as np import ome_zarr_models.v06.coordinate_transforms as ozm06trans -from spatialdata._core.transformation_manager.exceptions import AxisNotInCoordSystemError, UnmappableCoordSystemsError +from spatialdata._core.transformation_manager.exceptions import ( + AxisRedefinitionError, + EmptyTransformSequenceError, + IncompatibleCoordSystemsError, + MissingAxisError, + NotUnimodularError, + UnexpectedShapeError, + UnmappedAxisError, +) from spatialdata._types import ArrayLike from spatialdata.transformations.graph.vert import Axis, CoordSystem @@ -61,12 +69,18 @@ def _validate_transform_points_shapes(self, points: ArrayLike) -> None: """ Validate if the shape of the points (coordinates to be transformed) are consistent with the input size of the transformation. + + Raises + ------ + UnexpectedShapeError + if `points`'s shape is incompatible with this transformation's input shape """ input_size = len(self.input.axes) if len(points.shape) != 2 or points.shape[1] != input_size: - raise ValueError( - f"points must be a tensor of shape (n, d), where n is the number of points and d is the " - f"the number of dimensions. Points shape: {points.shape}, input size: {input_size}" + raise UnexpectedShapeError( + array_name="points", + expected_shape=f"(, {self.input.num_axes})", + array_shape=points.shape, ) # order of the composition: self is applied first, then the transformation passed as argument @@ -131,10 +145,14 @@ def __init__( expected_linear_shape = (num_outputs, num_inputs) if linear.shape != expected_linear_shape: - raise ValueError(f"linear's shape is {linear.shape}. Expected f{(num_outputs, num_inputs)}") + raise UnexpectedShapeError( + array_name="linear", array_shape=linear.shape, expected_shape=expected_linear_shape + ) expected_translation_shape = (num_outputs,) if translation.shape != expected_translation_shape: - raise ValueError(f"translation's shape is {translation.shape}. Expected {expected_translation_shape}") + raise UnexpectedShapeError( + array_name="translation", array_shape=translation.shape, expected_shape=expected_translation_shape + ) self.linear = linear self.translation = translation @@ -255,7 +273,9 @@ def __init__( Output coordinate system of the transformation. """ if input.num_axes != output.num_axes: - raise ValueError("Input and output must have the same number of dimensions") + raise IncompatibleCoordSystemsError( + input=input, output=output, message="Axes must have the same number of dimensions" + ) super().__init__(input=input, output=output, name=name) def inverse(self, name: str | None = None) -> BaseTransfEdge: @@ -304,7 +324,9 @@ def __init__( """ if set(input.axes) != set(output.axes): - raise UnmappableCoordSystemsError(input=input, output=output) + raise IncompatibleCoordSystemsError( + input=input, output=output, message="Input and output must have the same axes" + ) super().__init__(input=input, output=output, name=name) def __repr__(self) -> str: @@ -354,12 +376,27 @@ def __init__( dropped_inputs: set[Axis], created_outputs: set[Axis], ) -> None: + """ + Parameters + ---------- + dropped_inputs + axes in `input` that will be dropped by this transformation + created_inputs + axes in `output` that will be set to 0 + + Raises + ------ + MissingAxisError + axis in `dropped_inputs` not in `input` + axis in `created_outputs` not in `output` + """ + for axis in dropped_inputs: if axis not in input.axes: - raise AxisNotInCoordSystemError(axis=axis, cs=input) + raise MissingAxisError(axis=axis, cs=input) for axis in dropped_inputs: if axis not in output.axes: - raise AxisNotInCoordSystemError(axis=axis, cs=output) + raise MissingAxisError(axis=axis, cs=output) self.dropped_inputs = set(dropped_inputs) self.created_outputs = set(created_outputs) super().__init__(name=name, input=input, output=output) @@ -437,7 +474,6 @@ def __init__( output: CoordSystem, ) -> None: """ - Init the NgffTranslation object. Parameters ---------- name @@ -448,9 +484,16 @@ def __init__( Input coordinate system of the transformation. output Output coordinate system of the transformation. + + Raises + ------ + IncompatibleCoordSystemsError + If the input and output have different number of dimensions """ if input.num_axes != output.num_axes: - raise ValueError("Number of input and output axes must be the same") + raise IncompatibleCoordSystemsError( + input=input, output=output, message="Number of input and output axes must be the same" + ) self.translation = translation super().__init__(input=input, output=output, name=name) @@ -499,20 +542,29 @@ def __init__( output: CoordSystem, ) -> None: """ - Init the NgffScale object. Parameters ---------- scale - A list of numbers or a vector specifying the scale along each axis. + A vector specifying the scale along each axis of `input`. input Input coordinate system of the transformation. output Output coordinate system of the transformation. + + Raises + ------ + UnexpectedShapeError + If scale doesn't have the same number of elements as input has axes + IncompatibleCoordSystemsError + If input and output have different number of axes """ - if scale.shape != (input.num_axes,): - raise ValueError(f"scale should be of shape f{(input.num_axes,)}") + expected_scale_shape = (input.num_axes,) + if scale.shape != expected_scale_shape: + raise UnexpectedShapeError(array_name="scale", array_shape=scale.shape, expected_shape=expected_scale_shape) if input.num_axes != output.num_axes: - raise ValueError("input and output must have same number of dimensions") + raise IncompatibleCoordSystemsError( + input=input, output=output, message="input and output must have same number of dimensions" + ) self.scale = scale super().__init__(input=input, output=output, name=name) @@ -565,7 +617,6 @@ def __init__( output: CoordSystem, ) -> None: """ - Init the NgffRotation object. Parameters ---------- linear_matrix @@ -574,14 +625,26 @@ def __init__( Input coordinate system of the transformation. output Output coordinate system of the transformation. + Raises + ------ + UnexpectedShapeError + if linear_matrix's shape isn't (output.num_axes, input.num_axes) + IncompatibleCoordSystemsError + if input and output don't have the same number of axes + NotUnimodularError + if linear_matrix doesn't have determinant ~= 1 """ + if input.num_axes != output.num_axes: + raise IncompatibleCoordSystemsError( + input=input, output=output, message="input and output should have the same numbe rof axes" + ) expected_shape = (output.num_axes, input.num_axes) if linear_matrix.shape != expected_shape: - raise ValueError(f"linear matrix should have shape {expected_shape}") - if input.num_axes != output.num_axes: - raise ValueError("input and output should have the same numbe rof axes") + raise UnexpectedShapeError( + array_name="linear_matrix", array_shape=linear_matrix.shape, expected_shape=expected_shape + ) if not np.isclose(np.linalg.det(linear_matrix), 1.0): - raise ValueError("det(linear_matrix) should be ~= 1") + raise NotUnimodularError(matrix=linear_matrix) linear_matrix.flags.writeable = False self.rotation = linear_matrix super().__init__(input=input, output=output, name=name) @@ -638,13 +701,26 @@ def __init__( ---------- transformations The transformations which compose the sequence. + Raises + ------ + EmptyTransformSequence + if `transformations` is empty + IncompatibleCoordSystemsError + if any item in `transformations` is incompatible with its neighbors """ if len(transformations) == 0: - raise ValueError("Empty transformation list") + raise EmptyTransformSequenceError() previous_transf = transformations[0] - for current_transf in transformations[1:]: + for transf_idx, current_transf in enumerate(transformations[1:], start=1): if previous_transf.output != current_transf.input: - raise ValueError(f"Mismatched input/output from {previous_transf} to {current_transf}") + raise IncompatibleCoordSystemsError( + input=current_transf.input, + output=previous_transf.output, + message=( + f"Output of transformation #{transf_idx - 1} is different " + f"from input of transformation #{transf_idx}" + ), + ) previous_transf = current_transf self.transformations = transformations super().__init__( @@ -709,8 +785,6 @@ def __init__( output: CoordSystem, ) -> None: """ - Init the ByDimension object. - Parameters ---------- transformations @@ -720,26 +794,34 @@ def __init__( The input coordinate system of the transformation. output The output coordinate system of the transformation. + Raises + ------ + MissingAxisError + if any input axis from `transformations` is not present in `input` or + if any output axis from `transformations` is not present in `output`. + AxisRedefinitionError + if an output axis is specified by more than one item of `transformations` + UnmappedAxisError + if axis of `output` is not covered by any item of `transformations` """ # we check that: # 1. each input from each transformation in self.transformation must appear in the set of input axes # 2. each output from each transformation in self.transformation must appear at most once in the set of output # axes - input_axes = input.axes_names - output_axes = output.axes_names defined_output_axes: set[str] = set() for t in transformations: - for ax in t.input.axes_names: - if ax not in input_axes: - raise ValueError(f"By dimension axis {ax} not in {input_axes}") - for ax in t.output.axes_names: - if ax not in output_axes: - raise ValueError(f"Axis {ax} not in output axes {output_axes}") - if ax in defined_output_axes: - raise ValueError(f"Output axis {ax} is defined more than once") - defined_output_axes.add(ax) - if len(output_axes) != len(defined_output_axes): - raise ValueError("Not all outputs are mapped") + for ax in t.input.axes: + if ax not in input.axes: + raise MissingAxisError(axis=ax, cs=input) + for ax in t.output.axes: + if ax not in output.axes: + raise MissingAxisError(axis=ax, cs=output) + if ax.name in defined_output_axes: + raise AxisRedefinitionError(axis=ax) + defined_output_axes.add(ax.name) + for ax in output.axes: + if ax.name not in defined_output_axes: + raise UnmappedAxisError(axis=ax, cs=output) self.transformations = tuple(transformations) super().__init__(input=input, output=output, name=name) From 96066f3d0af07d99b33cb09f90b8b6f9f48274ce Mon Sep 17 00:00:00 2001 From: Tomaz Vieira Date: Tue, 8 Sep 2026 08:59:35 +0000 Subject: [PATCH 07/26] Fix Affine.from_affine_matrix slicing] --- src/spatialdata/transformations/graph/edge.py | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/src/spatialdata/transformations/graph/edge.py b/src/spatialdata/transformations/graph/edge.py index aba0de7e5..30a0f98d1 100644 --- a/src/spatialdata/transformations/graph/edge.py +++ b/src/spatialdata/transformations/graph/edge.py @@ -196,7 +196,7 @@ def from_affine_matrix( """ return AffineEdge( linear=affine_matrix[:-1, :-1], - translation=affine_matrix[-1, :-1], + translation=affine_matrix[:-1, -1], input=input, output=output, name=name, From b15281ecf4248c0be296cbcc63854069e320c617 Mon Sep 17 00:00:00 2001 From: Tomaz Vieira Date: Tue, 8 Sep 2026 11:26:18 +0000 Subject: [PATCH 08/26] Fix coords bug when reading from OMEZarrMultiscale. Some cleanup --- src/spatialdata/_io/io_raster.py | 37 ++++++++++++++++++++------------ 1 file changed, 23 insertions(+), 14 deletions(-) diff --git a/src/spatialdata/_io/io_raster.py b/src/spatialdata/_io/io_raster.py index 29f40d6d8..884a39e0a 100644 --- a/src/spatialdata/_io/io_raster.py +++ b/src/spatialdata/_io/io_raster.py @@ -171,7 +171,11 @@ def _prepare_storage_options( def try_read_ngff06_multiscale(store: Path) -> tuple[DataTree, Sequence[BaseTransfEdge]]: multiscale = oz.OMEZarrMultiscale.from_ome_zarr(str(store)) assert isinstance(multiscale, oz.OMEZarrMultiscale) # disambiguate from OMEZarrLabel + return try_parse_ngff06_multiscale(multiscale) + +def try_parse_ngff06_multiscale(multiscale: oz.OMEZarrMultiscale) -> tuple[DataTree, Sequence[BaseTransfEdge]]: + """Parse an OMEZarMultiscale into a DataTree and collects Multiscale-level transforms.""" name_to_cs: dict[str, CoordSystem] = {} for cs in multiscale.metadata.coordinateSystems or (): parsed_cs = CoordSystem.try_from_model(cs) @@ -204,24 +208,27 @@ def try_read_ngff06_multiscale(store: Path) -> tuple[DataTree, Sequence[BaseTran for scale_idx, (ds_md, ds) in enumerate(zip(multiscale.metadata.datasets, multiscale.images, strict=True)): transf = ds_md.coordinateTransformations[0] - out_cs = name_to_cs[multiscale.metadata.intrinsic_coordinate_system.name] + intrinsic_cs = name_to_cs[multiscale.metadata.intrinsic_coordinate_system.name] assert transf.input is not None assert transf.input.path is not None - in_cs = CoordSystem(name=str(transf.input.name), axes=[Axis(name=ax.name, type=ax.type) for ax in out_cs.axes]) + pixel_cs = CoordSystem( + name=transf.input.path, + axes=[Axis(name=ax.name, type=ax.type) for ax in intrinsic_cs.axes], + virtual=True, + ) - ozm_seq = ozm06trans.Sequence(transformations=ds_md.coordinateTransformations) - seq = parse_ngff_transf(input=in_cs, output=out_cs, model=ozm_seq) + pixel_cs_to_intrinsic_ngff = ozm06trans.Sequence(transformations=ds_md.coordinateTransformations) + seq = parse_ngff_transf(input=pixel_cs, output=intrinsic_cs, model=pixel_cs_to_intrinsic_ngff) ds_shape = np.asarray(ds.data.shape) - transformed_start = seq.transform_points(np.zeros_like(ds.data.shape)[np.newaxis, :])[0] + transformed_start = seq.transform_points(np.zeros_like(ds_shape)[np.newaxis, :])[0] transformed_stop = seq.transform_points((ds_shape - 1)[np.newaxis, :])[0] - coords = xr.Coordinates() - for low, high, ax, extent in zip(transformed_start, transformed_stop, out_cs.axes, ds_shape, strict=True): + coords: xr.Coordinates = xr.Coordinates() + for low, high, ax, extent in zip(transformed_start, transformed_stop, intrinsic_cs.axes, ds_shape, strict=True): if ax.type == "channel" and channel_names is not None: - coords.merge({ax.name: channel_names}) - continue - coords = coords.merge( - xr.Coordinates.from_xindex( + coords = coords.merge({ax.name: channel_names}).coords + else: + axis_index = xr.Coordinates.from_xindex( RangeIndex.linspace( start=low, stop=high, @@ -230,17 +237,19 @@ def try_read_ngff06_multiscale(store: Path) -> tuple[DataTree, Sequence[BaseTran dim=ax.name, ) ) - ) + coords = coords.merge(axis_index).coords + # Note: the magic "image" and "scale " strings mimic the current + # behavior from `dask_arrays_to_datatree` data_tree[f"scale{scale_idx}"] = xr.Dataset( { "image": xr.DataArray( ds.data, name="image", - dims=out_cs.axes_names, + dims=intrinsic_cs.axes_names, coords=coords, ) - } + }, ) return data_tree, parsed_transfs From 4cba32ef5a53acfc0015fb256dc45f17ebb7fcca Mon Sep 17 00:00:00 2001 From: Tomaz Vieira Date: Tue, 8 Sep 2026 11:30:05 +0000 Subject: [PATCH 09/26] Adds test for parsing OMEZarrMultiscale --- tests/io/test_ngff_06.py | 45 ++++++++++++++++++++++++++++++++++++++++ 1 file changed, 45 insertions(+) create mode 100644 tests/io/test_ngff_06.py diff --git a/tests/io/test_ngff_06.py b/tests/io/test_ngff_06.py new file mode 100644 index 000000000..46c2aea75 --- /dev/null +++ b/tests/io/test_ngff_06.py @@ -0,0 +1,45 @@ +from __future__ import annotations + +import numpy as np +from ome_zarr import OMEZarrImage, OMEZarrMultiscale +from ome_zarr_models.v06.coordinate_transforms import Scale as ModelsScale +from ome_zarr_models.v06.coordinate_transforms import Sequence as ModelsSequence +from ome_zarr_models.v06.coordinate_transforms import Translation as ModelsTranslation + +from spatialdata._io.io_raster import try_parse_ngff06_multiscale +from tests.conftest import SEED + + +def test_parse_multiscale(): + data = np.random.default_rng(seed=SEED).random((256, 256)) + omz_img = OMEZarrImage(data=data, axes="yx") + ms = OMEZarrMultiscale(image=omz_img, scale_factors=(2, 4, 8, 16)) + + data_tree, _transforms = try_parse_ngff06_multiscale(ms) + for xr_scale_node, ngff_scale, ngff_meta in zip( + data_tree.children.values(), ms.images, ms.metadata.datasets, strict=True + ): + xr_scale = xr_scale_node["image"] + assert xr_scale.shape == ngff_scale.data.shape + + coords_x = xr_scale.coords["x"] + coords_y = xr_scale.coords["y"] + + seq = ngff_meta.coordinateTransformations[0] + assert isinstance(seq, ModelsSequence) + scale = seq.transformations[0] + assert isinstance(scale, ModelsScale) + translate = seq.transformations[1] + assert isinstance(translate, ModelsTranslation) + + start_indices = (0, 0) + start = (coords_x[0], coords_y[0]) + + end_indices = (xr_scale.shape[0] - 1, xr_scale.shape[1] - 1) + end = (coords_x[-1], coords_y[-1]) + + for indices, point in [(start_indices, start), (end_indices, end)]: + expected = np.asarray(scale.scale) * indices + translate.translation + assert np.allclose(expected, point) + + return data_tree, ms From c7336587bd730f6468cd2d63e7153a346525d069 Mon Sep 17 00:00:00 2001 From: Tomaz Vieira Date: Wed, 9 Sep 2026 11:00:58 +0000 Subject: [PATCH 10/26] Adds transformations tests, fixes some slicing bugs Co-authored-by: Claude Sonnet 5 --- .../transformation_manager/exceptions.py | 2 + src/spatialdata/transformations/graph/edge.py | 33 +- tests/transformations/edge/__init__.py | 0 tests/transformations/edge/conftest.py | 19 + .../edge/test_edge_transformations.py | 365 ++++++++++++++++++ 5 files changed, 410 insertions(+), 9 deletions(-) create mode 100644 tests/transformations/edge/__init__.py create mode 100644 tests/transformations/edge/conftest.py create mode 100644 tests/transformations/edge/test_edge_transformations.py diff --git a/src/spatialdata/_core/transformation_manager/exceptions.py b/src/spatialdata/_core/transformation_manager/exceptions.py index d921beb9d..a9dec744b 100644 --- a/src/spatialdata/_core/transformation_manager/exceptions.py +++ b/src/spatialdata/_core/transformation_manager/exceptions.py @@ -282,4 +282,6 @@ def __init__(self, axis: Axis) -> None: class UnmappedAxisError(Exception): def __init__(self, axis: Axis, cs: CoordSystem) -> None: + self.axis = axis + self.cs = cs super().__init__(f"Axis {axis.name} from coordinate system {cs.name} is not mapped to anything") diff --git a/src/spatialdata/transformations/graph/edge.py b/src/spatialdata/transformations/graph/edge.py index 30a0f98d1..35370cc28 100644 --- a/src/spatialdata/transformations/graph/edge.py +++ b/src/spatialdata/transformations/graph/edge.py @@ -224,7 +224,7 @@ def inverse(self, name: str | None = None) -> BaseTransfEdge | None: return None return AffineEdge( linear=inv[:-1, :-1], - translation=inv[-1, :-1], + translation=inv[:-1, -1], input=self.output, output=self.input, name=name, @@ -345,7 +345,7 @@ def inverse(self, name: str | None = None) -> BaseTransfEdge: def transform_points(self, points: ArrayLike) -> ArrayLike: self._validate_transform_points_shapes(points) - new_indices = [self.input.axes.index(out_ax.name) for out_ax in self.output.axes] + new_indices = [self.input.axes.index(out_ax) for out_ax in self.output.axes] mapped = points[:, new_indices] assert isinstance(mapped, np.ndarray) return mapped @@ -389,27 +389,42 @@ def __init__( MissingAxisError axis in `dropped_inputs` not in `input` axis in `created_outputs` not in `output` + IncompatibleCoordSystemsError + when input can't be mapped to output given dropped_inputs and created_outputs """ for axis in dropped_inputs: if axis not in input.axes: raise MissingAxisError(axis=axis, cs=input) - for axis in dropped_inputs: + for axis in created_outputs: if axis not in output.axes: raise MissingAxisError(axis=axis, cs=output) + if input.num_axes - len(dropped_inputs) + len(created_outputs) != output.num_axes: + message = f"Can't map from {input} to {output}" + if dropped_inputs: + message += f" dropping {dropped_inputs}" + if created_outputs: + message += f" creating {created_outputs}" + raise IncompatibleCoordSystemsError( + input=input, + output=output, + message=message, + ) + self.dropped_inputs = set(dropped_inputs) self.created_outputs = set(created_outputs) super().__init__(name=name, input=input, output=output) def to_affine(self, name: str | None = None) -> AffineEdge: linear = np.zeros((self.output.num_axes, self.input.num_axes), dtype=float) - + input_indices = iter(range(self.input.num_axes)) for out_idx, out_ax in enumerate(self.output.axes): if out_ax in self.created_outputs: continue - for in_idx, in_ax in enumerate(self.input.axes): - if in_ax not in self.dropped_inputs: - linear[out_idx, in_idx] = 1 + in_idx = next(input_indices) + if in_idx in self.dropped_inputs: + continue + linear[out_idx, in_idx] = 1 return AffineEdge(name=name, input=self.input, output=self.output, linear=linear) @@ -449,7 +464,7 @@ def parse_project_axis( output = out else: num_dropped_inputs = len(model.droppedInputs or ()) - num_created_outputs = len(model.droppedInputs or ()) + num_created_outputs = len(model.createdOutputs or ()) num_output_axes = input.num_axes - num_dropped_inputs + num_created_outputs output = out.generate(num_axes=num_output_axes) @@ -740,7 +755,7 @@ def __repr__(self) -> str: def inverse(self, name: str | None = None) -> SequenceEdge | None: inverted: list[BaseTransfEdge] = [] - for t in self.transformations: + for t in reversed(self.transformations): inv = t.inverse() if inv is None: return None diff --git a/tests/transformations/edge/__init__.py b/tests/transformations/edge/__init__.py new file mode 100644 index 000000000..e69de29bb diff --git a/tests/transformations/edge/conftest.py b/tests/transformations/edge/conftest.py new file mode 100644 index 000000000..75de85b49 --- /dev/null +++ b/tests/transformations/edge/conftest.py @@ -0,0 +1,19 @@ +from __future__ import annotations + +from spatialdata.transformations.graph.vert import Axis, CoordSystem + +x_axis = Axis(name="x", type="space", unit="micrometer") +y_axis = Axis(name="y", type="space", unit="micrometer") +z_axis = Axis(name="z", type="space", unit="micrometer") +c_axis = Axis(name="c", type="channel") + +x_cs = CoordSystem(name="x", axes=[x_axis]) +y_cs = CoordSystem(name="y", axes=[y_axis]) +z_cs = CoordSystem(name="z", axes=[z_axis]) +c_cs = CoordSystem(name="c", axes=[c_axis]) +xy_cs = CoordSystem(name="xy", axes=[x_axis, y_axis]) +yx_cs = CoordSystem(name="yx", axes=[y_axis, x_axis]) +xyz_cs = CoordSystem(name="xyz", axes=[x_axis, y_axis, z_axis]) +zyx_cs = CoordSystem(name="zyx", axes=[z_axis, y_axis, x_axis]) +xyc_cs = CoordSystem(name="xyc", axes=[x_axis, y_axis, c_axis]) +cyx_cs = CoordSystem(name="cyx", axes=[c_axis, y_axis, x_axis]) diff --git a/tests/transformations/edge/test_edge_transformations.py b/tests/transformations/edge/test_edge_transformations.py new file mode 100644 index 000000000..8283369d3 --- /dev/null +++ b/tests/transformations/edge/test_edge_transformations.py @@ -0,0 +1,365 @@ +# pyright: strict + +from __future__ import annotations + +import numpy as np +import pytest + +from spatialdata._core.transformation_manager.exceptions import ( + AxisRedefinitionError, + EmptyTransformSequenceError, + IncompatibleCoordSystemsError, + MissingAxisError, + NotUnimodularError, + UnexpectedShapeError, + UnmappedAxisError, +) +from spatialdata.transformations.graph.edge import ( + AffineEdge, + ByDimensionEdge, + IdentityEdge, + MapAxisEdge, + ProjectAxisEdge, + RotationEdge, + ScaleEdge, + SequenceEdge, + TranslationEdge, +) +from tests.transformations.edge.conftest import ( + x_cs, + xy_cs, + xyc_cs, + xyz_cs, + y_axis, + y_cs, + yx_cs, + z_axis, + zyx_cs, +) + +POINTS_2D = np.array([[1.0, 2.0], [3.0, 4.0]]) +POINTS_3D = np.array([[1.0, 2.0, 3.0], [4.0, 5.0, 6.0]]) + + +class TestAffineEdge: + def test_constructor_rejects_wrong_linear_shape(self): + with pytest.raises(UnexpectedShapeError): + AffineEdge(linear=np.eye(3), input=xy_cs, output=xy_cs) + + def test_constructor_rejects_wrong_translation_shape(self): + with pytest.raises(UnexpectedShapeError): + AffineEdge(linear=np.eye(2), translation=np.zeros(3), input=xy_cs, output=xy_cs) + + def test_from_affine_matrix(self): + # fmt: off + affine_matrix = np.array([ + [2.0, 0.0, 1.0], + [0.0, 3.0, 5.0], + [0.0, 0.0, 1.0] + ]) + linear = [ + [2.0, 0.0], + [0.0, 3.0], + ] + translation = [ + 1.0, + 5.0 + ] + # fmt: on + edge = AffineEdge.from_affine_matrix(name=None, affine_matrix=affine_matrix, input=xy_cs, output=xy_cs) + np.testing.assert_equal(edge.linear, linear) + np.testing.assert_equal(edge.translation, translation) + np.testing.assert_allclose(edge.transform_points(POINTS_2D), np.array([[3.0, 11.0], [7.0, 17.0]])) + + def test_mapping_classmethod_builds_permutation_matrix(self): + edge = AffineEdge.mapping(input=xy_cs, output=yx_cs) + np.testing.assert_allclose(edge.linear, np.array([[0.0, 1.0], [1.0, 0.0]])) + np.testing.assert_allclose(edge.transform_points(POINTS_2D), np.array([[2.0, 1.0], [4.0, 3.0]])) + + def test_transform_points_scale_and_translate(self): + edge = AffineEdge( + linear=np.array([[2.0, 0.0], [0.0, 3.0]]), translation=np.array([1.0, 1.0]), input=xy_cs, output=xy_cs + ) + np.testing.assert_allclose(edge.transform_points(POINTS_2D), np.array([[3.0, 7.0], [7.0, 13.0]])) + + def test_transform_points_2d_to_3d(self): + linear = np.array([[1.0, 0.0], [0.0, 1.0], [1.0, 1.0]]) + edge = AffineEdge(linear=linear, input=xy_cs, output=xyz_cs) + np.testing.assert_allclose(edge.transform_points(POINTS_2D), np.array([[1.0, 2.0, 3.0], [3.0, 4.0, 7.0]])) + + def test_transform_points_rejects_wrong_shape(self): + edge = AffineEdge(linear=np.eye(2), input=xy_cs, output=xy_cs) + with pytest.raises(UnexpectedShapeError): + # expecting (n, 2), gets (2,) + edge.transform_points(np.asarray([1.0, 2.0])) + with pytest.raises(UnexpectedShapeError): + # expecting (n, 2), gets (n, 3) + edge.transform_points(POINTS_3D) + + def test_inverse_roundtrips(self): + edge = AffineEdge( + linear=np.array([[2.0, 0.0], [0.0, 4.0]]), translation=np.array([1.0, -1.0]), input=xy_cs, output=xy_cs + ) + inv = edge.inverse() + assert inv is not None + np.testing.assert_allclose(inv.transform_points(edge.transform_points(POINTS_2D)), POINTS_2D) + + def test_inverse_returns_none_for_singular_matrix(self): + edge = AffineEdge(linear=np.zeros((2, 2)), input=xy_cs, output=xy_cs) + assert edge.inverse() is None + + +class TestIdentityEdge: + def test_constructor_rejects_mismatched_number_of_axes(self): + with pytest.raises(IncompatibleCoordSystemsError): + IdentityEdge(None, input=xy_cs, output=xyz_cs) + + def test_transform_points_is_a_noop(self): + edge = IdentityEdge(None, input=xy_cs, output=xy_cs) + np.testing.assert_allclose(edge.transform_points(POINTS_2D), POINTS_2D) + + def test_to_affine_matches_transform_points(self): + edge = IdentityEdge(None, input=xy_cs, output=xy_cs) + affine = edge.to_affine() + np.testing.assert_allclose(affine.linear, np.eye(2)) + np.testing.assert_allclose(affine.transform_points(POINTS_2D), POINTS_2D) + + def test_inverse_is_still_the_identity(self): + edge = IdentityEdge(None, input=xy_cs, output=xy_cs) + inv = edge.inverse() + np.testing.assert_allclose(inv.transform_points(POINTS_2D), POINTS_2D) + + +class TestMapAxisEdge: + def test_constructor_rejects_different_sets_of_axes(self): + with pytest.raises(IncompatibleCoordSystemsError): + MapAxisEdge(input=xy_cs, output=xyz_cs) + + def test_transform_points_swaps_axes(self): + edge = MapAxisEdge(input=xy_cs, output=yx_cs) + np.testing.assert_allclose(edge.transform_points(POINTS_2D), np.array([[2.0, 1.0], [4.0, 3.0]])) + + def test_transform_points_permutation_of_three_axes(self): + edge = MapAxisEdge(input=xyz_cs, output=zyx_cs) + np.testing.assert_allclose(edge.transform_points(POINTS_3D), np.array([[3.0, 2.0, 1.0], [6.0, 5.0, 4.0]])) + + def test_to_affine_matches_transform_points(self): + edge = MapAxisEdge(input=xy_cs, output=yx_cs) + affine = edge.to_affine() + np.testing.assert_allclose(affine.transform_points(POINTS_2D), edge.transform_points(POINTS_2D)) + + def test_inverse_roundtrips(self): + edge = MapAxisEdge(input=xyz_cs, output=zyx_cs) + inv = edge.inverse() + np.testing.assert_allclose(inv.transform_points(edge.transform_points(POINTS_3D)), POINTS_3D) + + +class TestProjectAxisEdge: + def test_constructor_rejects_dropped_axis_missing_from_input(self): + with pytest.raises(MissingAxisError): + ProjectAxisEdge(input=xy_cs, output=xy_cs, dropped_inputs={z_axis}, created_outputs=set()) + + def test_constructor_rejects_created_axis_missing_from_output(self): + with pytest.raises(MissingAxisError): + ProjectAxisEdge(input=xy_cs, output=xy_cs, dropped_inputs=set(), created_outputs={z_axis}) + + def test_general_coord_system_incompatibility(self): + with pytest.raises(IncompatibleCoordSystemsError): + ProjectAxisEdge(input=xy_cs, output=xyz_cs, dropped_inputs=set(), created_outputs=set()) + with pytest.raises(IncompatibleCoordSystemsError): + ProjectAxisEdge(input=xyz_cs, output=xy_cs, dropped_inputs=set(), created_outputs=set()) + + def test_to_affine_for_identity_case(self): + edge = ProjectAxisEdge(input=xy_cs, output=xy_cs, dropped_inputs=set(), created_outputs=set()) + affine = edge.to_affine() + np.testing.assert_allclose(affine.linear, np.eye(2)) + np.testing.assert_allclose(edge.transform_points(POINTS_2D), POINTS_2D) + + def test_dropping_one_axis(self): + edge = ProjectAxisEdge(input=xyz_cs, output=xy_cs, dropped_inputs={z_axis}, created_outputs=set()) + affine = edge.to_affine() + assert affine is not None + np.testing.assert_allclose(edge.transform_points(POINTS_3D), np.array([[1.0, 2.0], [4.0, 5.0]])) + + def test_creating_one_axis(self): + edge = ProjectAxisEdge(input=xy_cs, output=xyz_cs, dropped_inputs=set(), created_outputs={z_axis}) + affine = edge.to_affine() + assert affine is not None + np.testing.assert_allclose(edge.transform_points(POINTS_2D), np.array([[1.0, 2.0, 0.0], [3.0, 4.0, 0.0]])) + + def test_inverse_none_when_axes_are_dropped_or_created(self): + dropped = ProjectAxisEdge(input=xyz_cs, output=xy_cs, dropped_inputs={z_axis}, created_outputs=set()) + assert dropped.inverse() is None + created = ProjectAxisEdge(input=xy_cs, output=xyz_cs, dropped_inputs=set(), created_outputs={z_axis}) + assert created.inverse() is None + + def test_inverse_roundtrip(self): + edge = ProjectAxisEdge(input=xy_cs, output=yx_cs, dropped_inputs=set(), created_outputs=set()) + inv = edge.inverse() + assert inv is not None + np.testing.assert_allclose(inv.transform_points(edge.transform_points(POINTS_2D)), POINTS_2D) + + +class TestTranslationEdge: + def test_constructor_rejects_mismatched_number_of_axes(self): + with pytest.raises(IncompatibleCoordSystemsError): + TranslationEdge(translation=np.array([1.0, 2.0]), input=xy_cs, output=xyz_cs) + + def test_transform_points_adds_translation(self): + edge = TranslationEdge(translation=np.array([10.0, 20.0]), input=xy_cs, output=xy_cs) + np.testing.assert_allclose(edge.transform_points(POINTS_2D), np.array([[11.0, 22.0], [13.0, 24.0]])) + + def test_to_affine_matches_transform_points(self): + edge = TranslationEdge(translation=np.array([10.0, 20.0]), input=xy_cs, output=xy_cs) + affine = edge.to_affine() + assert affine is not None + np.testing.assert_allclose(affine.transform_points(POINTS_2D), edge.transform_points(POINTS_2D)) + + def test_inverse_roundtrips(self): + edge = TranslationEdge(translation=np.array([10.0, 20.0]), input=xy_cs, output=xy_cs) + inv = edge.inverse() + np.testing.assert_allclose(inv.transform_points(edge.transform_points(POINTS_2D)), POINTS_2D) + + +class TestScaleEdge: + def test_constructor_rejects_wrong_scale_shape(self): + with pytest.raises(UnexpectedShapeError): + ScaleEdge(scale=np.array([1.0, 2.0, 3.0]), input=xy_cs, output=xy_cs) + + def test_constructor_rejects_mismatched_number_of_axes(self): + with pytest.raises(IncompatibleCoordSystemsError): + ScaleEdge(scale=np.array([1.0, 2.0]), input=xy_cs, output=xyz_cs) + + def test_transform_points(self): + edge = ScaleEdge(scale=np.array([2.0, 4.0]), input=xy_cs, output=xy_cs) + np.testing.assert_allclose(edge.transform_points(POINTS_2D), np.array([[2.0, 8.0], [6.0, 16.0]])) + + def test_to_affine_matches_transform_points(self): + edge = ScaleEdge(scale=np.array([2.0, 4.0]), input=xy_cs, output=xy_cs) + affine = edge.to_affine() + np.testing.assert_allclose(affine.transform_points(POINTS_2D), edge.transform_points(POINTS_2D)) + + def test_inverse_roundtrips(self): + edge = ScaleEdge(scale=np.array([2.0, 4.0]), input=xy_cs, output=xy_cs) + inv = edge.inverse() + assert inv is not None + np.testing.assert_allclose(inv.transform_points(edge.transform_points(POINTS_2D)), POINTS_2D) + + def test_inverse_returns_none_when_scale_is_zero(self): + edge = ScaleEdge(scale=np.array([0.0, 4.0]), input=xy_cs, output=xy_cs) + assert edge.inverse() is None + + +class TestRotationEdge: + def test_constructor_rejects_mismatched_number_of_axes(self): + with pytest.raises(IncompatibleCoordSystemsError): + RotationEdge(linear_matrix=np.eye(2), input=xy_cs, output=xyz_cs) + + def test_constructor_rejects_wrong_shape(self): + with pytest.raises(UnexpectedShapeError): + RotationEdge(linear_matrix=np.eye(3), input=xy_cs, output=xy_cs) + + def test_constructor_rejects_non_unimodular_matrix(self): + with pytest.raises(NotUnimodularError): + RotationEdge(linear_matrix=np.array([[1.0, 0.0], [0.0, -1.0]]), input=xy_cs, output=xy_cs) + + def test_transform_points_rotates_90_degrees(self): + edge = RotationEdge(linear_matrix=np.array([[0.0, -1.0], [1.0, 0.0]]), input=xy_cs, output=xy_cs) + np.testing.assert_allclose(edge.transform_points(POINTS_2D), np.array([[-2.0, 1.0], [-4.0, 3.0]])) + + def test_to_affine_matches_transform_points(self): + edge = RotationEdge(linear_matrix=np.array([[0.0, -1.0], [1.0, 0.0]]), input=xy_cs, output=xy_cs) + affine = edge.to_affine() + np.testing.assert_allclose(affine.transform_points(POINTS_2D), edge.transform_points(POINTS_2D)) + + def test_inverse_roundtrips(self): + edge = RotationEdge(linear_matrix=np.array([[0.0, -1.0], [1.0, 0.0]]), input=xy_cs, output=xy_cs) + inv = edge.inverse() + np.testing.assert_allclose(inv.transform_points(edge.transform_points(POINTS_2D)), POINTS_2D) + + +class TestSequenceEdge: + def test_constructor_rejects_empty_sequence(self): + with pytest.raises(EmptyTransformSequenceError): + SequenceEdge(transformations=[]) + + def test_constructor_rejects_incompatible_neighbors(self): + first = TranslationEdge(translation=np.array([1.0, 2.0]), input=xy_cs, output=xy_cs) + second = TranslationEdge(translation=np.array([1.0, 2.0, 3.0]), input=xyz_cs, output=xyz_cs) + with pytest.raises(IncompatibleCoordSystemsError): + SequenceEdge(transformations=[first, second]) + + def test_transform_points_composes_in_order(self): + translate = TranslationEdge(translation=np.array([1.0, 2.0]), input=xy_cs, output=xy_cs) + scale = ScaleEdge(scale=np.array([3.0, 4.0]), input=xy_cs, output=xy_cs) + edge = SequenceEdge(transformations=[translate, scale]) + expected = (POINTS_2D + np.array([1.0, 2.0])) * np.array([3.0, 4.0]) + np.testing.assert_allclose(edge.transform_points(POINTS_2D), expected) + + def test_to_affine_matches_transform_points(self): + translate = TranslationEdge(translation=np.array([1.0, 2.0]), input=xy_cs, output=xy_cs) + scale = ScaleEdge(scale=np.array([3.0, 4.0]), input=xy_cs, output=xy_cs) + edge = SequenceEdge(transformations=[translate, scale]) + affine = edge.to_affine() + np.testing.assert_allclose(affine.transform_points(POINTS_2D), edge.transform_points(POINTS_2D)) + + def test_inverse_roundtrips(self): + translate = TranslationEdge(translation=np.array([1.0, 2.0]), input=xy_cs, output=xy_cs) + scale = ScaleEdge(scale=np.array([3.0, 4.0]), input=xy_cs, output=xy_cs) + edge = SequenceEdge(transformations=[translate, scale]) + inv = edge.inverse() + assert inv is not None + np.testing.assert_allclose(inv.transform_points(edge.transform_points(POINTS_2D)), POINTS_2D) + + def test_inverse_returns_none_if_any_component_is_not_invertible(self): + scale = ScaleEdge(scale=np.array([0.0, 4.0]), input=xy_cs, output=xy_cs) + translate = TranslationEdge(translation=np.array([1.0, 2.0]), input=xy_cs, output=xy_cs) + edge = SequenceEdge(transformations=[scale, translate]) + assert edge.inverse() is None + + +class TestByDimensionEdge: + def test_constructor_rejects_input_axis_missing_from_overall_input(self): + sub = IdentityEdge(None, input=xyc_cs, output=xyc_cs) + with pytest.raises(MissingAxisError): + ByDimensionEdge(transformations=[sub], input=xy_cs, output=xyc_cs) + + def test_constructor_rejects_output_axis_missing_from_overall_output(self): + sub = IdentityEdge(None, input=xy_cs, output=xy_cs) + with pytest.raises(MissingAxisError): + ByDimensionEdge(transformations=[sub], input=xy_cs, output=x_cs) + + def test_constructor_rejects_output_axis_defined_more_than_once(self): + first = IdentityEdge(None, input=x_cs, output=x_cs) + second = IdentityEdge(None, input=x_cs, output=x_cs) + with pytest.raises(AxisRedefinitionError): + ByDimensionEdge(transformations=[first, second], input=x_cs, output=x_cs) + + def test_constructor_rejects_unmapped_output_axis(self): + sub = IdentityEdge(None, input=x_cs, output=x_cs) + try: + ByDimensionEdge(transformations=[sub], input=xy_cs, output=xy_cs) + raise AssertionError(f"Expected {UnmappedAxisError.__name__} to be raised") + except UnmappedAxisError as e: + assert e.axis == y_axis + + def test_transform_points_applies_each_transformation_per_axis(self): + scale_x = ScaleEdge(scale=np.array([2.0]), input=x_cs, output=x_cs) + translate_y = TranslationEdge(translation=np.array([5.0]), input=y_cs, output=y_cs) + edge = ByDimensionEdge(transformations=[scale_x, translate_y], input=xy_cs, output=xy_cs) + np.testing.assert_allclose(edge.transform_points(POINTS_2D), np.array([[2.0, 7.0], [6.0, 9.0]])) + + def test_to_affine_matches_transform_points(self): + scale_x = ScaleEdge(scale=np.array([2.0]), input=x_cs, output=x_cs) + translate_y = TranslationEdge(translation=np.array([5.0]), input=y_cs, output=y_cs) + edge = ByDimensionEdge(transformations=[scale_x, translate_y], input=xy_cs, output=xy_cs) + affine = edge.to_affine() + np.testing.assert_allclose(affine.transform_points(POINTS_2D), edge.transform_points(POINTS_2D)) + + def test_inverse_roundtrips(self): + scale_x = ScaleEdge(scale=np.array([2.0]), input=x_cs, output=x_cs) + translate_y = TranslationEdge(translation=np.array([5.0]), input=y_cs, output=y_cs) + edge = ByDimensionEdge(transformations=[scale_x, translate_y], input=xy_cs, output=xy_cs) + inv = edge.inverse() + assert inv is not None + np.testing.assert_allclose(inv.transform_points(edge.transform_points(POINTS_2D)), POINTS_2D) From 8ee025ea4050773fc5edcfb74264660dfb5ee982 Mon Sep 17 00:00:00 2001 From: Tomaz Vieira Date: Wed, 9 Sep 2026 12:08:42 +0000 Subject: [PATCH 11/26] Adds basic .sel test to parsed multiscale --- tests/io/test_ngff_06.py | 3 +++ 1 file changed, 3 insertions(+) diff --git a/tests/io/test_ngff_06.py b/tests/io/test_ngff_06.py index 46c2aea75..7e4288023 100644 --- a/tests/io/test_ngff_06.py +++ b/tests/io/test_ngff_06.py @@ -42,4 +42,7 @@ def test_parse_multiscale(): expected = np.asarray(scale.scale) * indices + translate.translation assert np.allclose(expected, point) + sliced_data_tree = data_tree.sel(x=slice(0, 256), y=slice(0, 256), method="nearest") + assert sliced_data_tree.equals(data_tree) + return data_tree, ms From ad1d8091f0173ef04001e4816ed012456badb69d Mon Sep 17 00:00:00 2001 From: Tomaz Vieira Date: Wed, 9 Sep 2026 12:08:57 +0000 Subject: [PATCH 12/26] Removes unused class --- src/spatialdata/transformations/graph/vert.py | 3 --- 1 file changed, 3 deletions(-) diff --git a/src/spatialdata/transformations/graph/vert.py b/src/spatialdata/transformations/graph/vert.py index 128c50706..c3401b76c 100644 --- a/src/spatialdata/transformations/graph/vert.py +++ b/src/spatialdata/transformations/graph/vert.py @@ -108,9 +108,6 @@ class CoordSystem: """A virtual coordinate system exists as an intermediate step between non-virtual coordinate systems and is usually ignored during serialization""" - class LegacyAxes: - pass - def __init__(self, name: str, axes: Sequence[Axis], virtual: bool = False): self.name = name self.axes = tuple(axes) From 99ecb22eb10dfbf32f5e4e87649e96ac657f0b91 Mon Sep 17 00:00:00 2001 From: Tomaz Vieira Date: Wed, 9 Sep 2026 12:26:17 +0000 Subject: [PATCH 13/26] Adds more slicing tests to parsed multiscales --- tests/io/test_ngff_06.py | 12 +++++++++--- 1 file changed, 9 insertions(+), 3 deletions(-) diff --git a/tests/io/test_ngff_06.py b/tests/io/test_ngff_06.py index 7e4288023..222bae81e 100644 --- a/tests/io/test_ngff_06.py +++ b/tests/io/test_ngff_06.py @@ -11,9 +11,11 @@ def test_parse_multiscale(): - data = np.random.default_rng(seed=SEED).random((256, 256)) + full_shape = (256, 256) + data = np.random.default_rng(seed=SEED).random(full_shape) omz_img = OMEZarrImage(data=data, axes="yx") - ms = OMEZarrMultiscale(image=omz_img, scale_factors=(2, 4, 8, 16)) + scale_factors = (2, 4, 8, 16) + ms = OMEZarrMultiscale(image=omz_img, scale_factors=scale_factors) data_tree, _transforms = try_parse_ngff06_multiscale(ms) for xr_scale_node, ngff_scale, ngff_meta in zip( @@ -45,4 +47,8 @@ def test_parse_multiscale(): sliced_data_tree = data_tree.sel(x=slice(0, 256), y=slice(0, 256), method="nearest") assert sliced_data_tree.equals(data_tree) - return data_tree, ms + quarter_tree = data_tree.sel(x=slice(0, 128), y=slice(0, 128), method="nearest") + for xr_scale_node, scale_factor in zip(quarter_tree.children.values(), [1, *scale_factors], strict=True): + xr_scale = xr_scale_node["image"] + expected_shape = np.asarray(full_shape) / 2 / scale_factor + np.testing.assert_equal(xr_scale.shape, expected_shape) From 7361863d8586cd2f2f37dec3df7f51a83014f760 Mon Sep 17 00:00:00 2001 From: Tomaz Vieira Date: Wed, 9 Sep 2026 13:07:48 +0000 Subject: [PATCH 14/26] Adds/fixes comments in try_parse_ngff06_multiscales --- src/spatialdata/_io/io_raster.py | 7 ++++++- 1 file changed, 6 insertions(+), 1 deletion(-) diff --git a/src/spatialdata/_io/io_raster.py b/src/spatialdata/_io/io_raster.py index 884a39e0a..655c1da9a 100644 --- a/src/spatialdata/_io/io_raster.py +++ b/src/spatialdata/_io/io_raster.py @@ -189,9 +189,9 @@ def try_parse_ngff06_multiscale(multiscale: oz.OMEZarrMultiscale) -> tuple[DataT assert in_cs_id is not None assert out_cs_ref is not None + # FIXME: not handling references into labels yet, which use name and path in_cs_name = in_cs_id.name out_cs_name = out_cs_ref.name - # FIXME: not handling references into labels yet assert in_cs_name is not None assert out_cs_name is not None @@ -211,6 +211,11 @@ def try_parse_ngff06_multiscale(multiscale: oz.OMEZarrMultiscale) -> tuple[DataT intrinsic_cs = name_to_cs[multiscale.metadata.intrinsic_coordinate_system.name] assert transf.input is not None assert transf.input.path is not None + + # This coord system doesn't exist explicitly in the NGFF file, nor will + # it exist in our graph of transformations; It is only created here + # for the sake of creating the transformations that will be expressed + # in levels of a xr.DataTree pixel_cs = CoordSystem( name=transf.input.path, axes=[Axis(name=ax.name, type=ax.type) for ax in intrinsic_cs.axes], From 0718fb72a6cbba13e82da6d97ce7835319d079f1 Mon Sep 17 00:00:00 2001 From: Tomaz Vieira Date: Wed, 9 Sep 2026 13:14:24 +0000 Subject: [PATCH 15/26] Removes to_model methods for now --- src/spatialdata/transformations/graph/edge.py | 85 ------------------- src/spatialdata/transformations/graph/vert.py | 23 ----- 2 files changed, 108 deletions(-) diff --git a/src/spatialdata/transformations/graph/edge.py b/src/spatialdata/transformations/graph/edge.py index 35370cc28..854fe18f5 100644 --- a/src/spatialdata/transformations/graph/edge.py +++ b/src/spatialdata/transformations/graph/edge.py @@ -103,10 +103,6 @@ def compose_with(self, transformation: BaseTransfEdge, name: str | None) -> Base """ return SequenceEdge(transformations=[self, transformation], name=name) - @abstractmethod - def to_model(self) -> ozm06trans.AnyTransform: - pass - class AffineEdge(BaseTransfEdge): """The Affine transformation from the NGFF specification.""" @@ -243,14 +239,6 @@ def to_affine(self, name: str | None = None) -> AffineEdge: input=self.input, output=self.output, linear=self.linear, translation=self.translation, name=name ) - def to_model(self) -> ozm06trans.Affine: - return ozm06trans.Affine( - name=self.name, - affine=tuple(tuple(row) for row in self.affine[:-1, :]), - input=self.input.to_model_cs_ident(), - output=self.output.to_model_cs_ident(), - ) - class IdentityEdge(BaseTransfEdge): """The Identity transformation from the NGFF specification.""" @@ -293,13 +281,6 @@ def to_affine(self, name: str | None = None) -> AffineEdge: name=name, ) - def to_model(self) -> ozm06trans.Identity: - return ozm06trans.Identity( - name=self.name, - input=self.input.to_model_cs_ident(), - output=self.output.to_model_cs_ident(), - ) - class MapAxisEdge(BaseTransfEdge): """The MapAxis transformation from the NGFF specification.""" @@ -353,15 +334,6 @@ def transform_points(self, points: ArrayLike) -> ArrayLike: def to_affine(self, name: str | None = None) -> AffineEdge: return AffineEdge.mapping(input=self.input, output=self.output, name=name) - def to_model(self) -> ozm06trans.MapAxis: - mapAxis: list[int] = [self.input.axes.index(out_ax) for out_ax in self.output.axes] - return ozm06trans.MapAxis( - name=self.name, - mapAxis=tuple(mapAxis), - input=self.input.to_model_cs_ident(), - output=self.output.to_model_cs_ident(), - ) - class ProjectAxisEdge(BaseTransfEdge): dropped_inputs: Final[set[Axis]] @@ -447,12 +419,6 @@ def inverse(self, name: str | None = None) -> BaseTransfEdge | None: name=name, ) - def to_model(self) -> ozm06trans.ProjectAxis: - return ozm06trans.ProjectAxis( - createdOutputs=tuple(self.output.axes.index(co) for co in self.created_outputs) or None, - droppedInputs=tuple(self.input.axes.index(di) for di in self.dropped_inputs) or None, - ) - def parse_project_axis( model: ozm06trans.ProjectAxis, @@ -536,14 +502,6 @@ def to_affine(self, name: str | None = None) -> AffineEdge: name=name, ) - def to_model(self) -> ozm06trans.Translation: - return ozm06trans.Translation( - name=self.name, - input=self.input.to_model_cs_ident(), - output=self.output.to_model_cs_ident(), - translation=tuple(self.translation), - ) - class ScaleEdge(BaseTransfEdge): """The Scale transformation from the NGFF specification.""" @@ -609,14 +567,6 @@ def to_affine(self, name: str | None = None) -> AffineEdge: name=name, ) - def to_model(self) -> ozm06trans.Scale: - return ozm06trans.Scale( - name=self.name, - input=self.input.to_model_cs_ident(), - output=self.output.to_model_cs_ident(), - scale=tuple(self.scale), - ) - class RotationEdge(BaseTransfEdge): """The Rotation transformation from the NGFF specification.""" @@ -691,14 +641,6 @@ def to_affine(self, name: str | None = None) -> AffineEdge: name=name, ) - def to_model(self) -> ozm06trans.Rotation: - return ozm06trans.Rotation( - name=self.name, - input=self.input.to_model_cs_ident(), - output=self.output.to_model_cs_ident(), - rotation=tuple(tuple(row) for row in self.rotation), - ) - class SequenceEdge(BaseTransfEdge): """The Sequence transformation from the NGFF specification.""" @@ -777,14 +719,6 @@ def to_affine(self, name: str | None = None) -> AffineEdge: def transform_points(self, points: ArrayLike) -> ArrayLike: return self.to_affine().transform_points(points) # FIXME - def to_model(self) -> ozm06trans.Sequence: - return ozm06trans.Sequence( - name=self.name, - input=self.input.to_model_cs_ident(), - output=self.output.to_model_cs_ident(), - transformations=tuple(t.to_model() for t in self.transformations), - ) - class ByDimensionEdge(BaseTransfEdge): """The ByDimension transformation from the NGFF specification.""" @@ -896,25 +830,6 @@ def to_affine(self, name: str | None = None) -> AffineEdge: name=name, ) - def to_model(self) -> ozm06trans.ByDimension: - by_dim_transfs: list[ozm06trans.ByDimensionTransform] = [] - for t in self.transformations: - input_axes = tuple(self.input.axes_names.index(ax_name) for ax_name in t.input.axes_names) - output_axes = tuple(self.output.axes_names.index(ax_name) for ax_name in t.output.axes_names) - by_dim_transfs.append( - ozm06trans.ByDimensionTransform( - input_axes=input_axes, - output_axes=output_axes, - transformation=t.to_model(), - ) - ) - return ozm06trans.ByDimension( - name=self.name, - input=self.input.to_model_cs_ident(), - output=self.output.to_model_cs_ident(), - transformations=tuple(by_dim_transfs), - ) - class CsGen: """A coordinate system generator diff --git a/src/spatialdata/transformations/graph/vert.py b/src/spatialdata/transformations/graph/vert.py index c3401b76c..1f9886b7d 100644 --- a/src/spatialdata/transformations/graph/vert.py +++ b/src/spatialdata/transformations/graph/vert.py @@ -75,15 +75,6 @@ def try_from_model(cls, model: ozm06ct.Axis) -> Axis: long_name=model.longName, ) - def to_model(self) -> ozm06ct.Axis: - return ozm06ct.Axis( - discrete=False, - longName=self.long_name, - name=self.name, - type=self.type, - unit=self.unit, - ) - class CoordSystemParsingException(Exception): pass @@ -139,20 +130,6 @@ def try_from_model_or_default[T](cls, model: ozi.CoordinateSystem | None, *, def return CoordSystem.try_from_model(model) return default - def to_model(self) -> ozi.CoordinateSystem | None: - if self.virtual: - return None - return ozi.CoordinateSystem( - name=self.name, - axes=tuple(ax.to_model() for ax in self.axes), - ) - - def to_model_cs_ident(self) -> ozi.CoordinateSystemIdentifier | None: - input = self.to_model() - if input is None: - return None - return ozi.CoordinateSystemIdentifier(name=input.name) - @property def num_axes(self) -> int: return len(self.axes) From df63609fd1780326b104f329fb0b07d231829d12 Mon Sep 17 00:00:00 2001 From: Tomaz Vieira Date: Wed, 9 Sep 2026 13:15:19 +0000 Subject: [PATCH 16/26] Moves AxisParsingException to exceptions.py --- src/spatialdata/_core/transformation_manager/exceptions.py | 4 ++++ src/spatialdata/transformations/graph/vert.py | 4 +--- 2 files changed, 5 insertions(+), 3 deletions(-) diff --git a/src/spatialdata/_core/transformation_manager/exceptions.py b/src/spatialdata/_core/transformation_manager/exceptions.py index a9dec744b..8224bd60d 100644 --- a/src/spatialdata/_core/transformation_manager/exceptions.py +++ b/src/spatialdata/_core/transformation_manager/exceptions.py @@ -285,3 +285,7 @@ def __init__(self, axis: Axis, cs: CoordSystem) -> None: self.axis = axis self.cs = cs super().__init__(f"Axis {axis.name} from coordinate system {cs.name} is not mapped to anything") + + +class AxisParsingException(Exception): + pass diff --git a/src/spatialdata/transformations/graph/vert.py b/src/spatialdata/transformations/graph/vert.py index 1f9886b7d..a6231813e 100644 --- a/src/spatialdata/transformations/graph/vert.py +++ b/src/spatialdata/transformations/graph/vert.py @@ -6,9 +6,7 @@ import ome_zarr.classes.image as ozi import ome_zarr_models.v06.coordinate_transforms as ozm06ct - -class AxisParsingException(Exception): - pass +from spatialdata._core.transformation_manager.exceptions import AxisParsingException class Axis: From 76f7c6b65ddd8010a1809ce0ee05707873562a50 Mon Sep 17 00:00:00 2001 From: Tomaz Vieira Date: Wed, 9 Sep 2026 13:51:41 +0000 Subject: [PATCH 17/26] Moves AxisParsingException back into vert.py to prevent circular imports --- src/spatialdata/_core/transformation_manager/exceptions.py | 4 ---- src/spatialdata/transformations/graph/vert.py | 4 +++- 2 files changed, 3 insertions(+), 5 deletions(-) diff --git a/src/spatialdata/_core/transformation_manager/exceptions.py b/src/spatialdata/_core/transformation_manager/exceptions.py index 8224bd60d..a9dec744b 100644 --- a/src/spatialdata/_core/transformation_manager/exceptions.py +++ b/src/spatialdata/_core/transformation_manager/exceptions.py @@ -285,7 +285,3 @@ def __init__(self, axis: Axis, cs: CoordSystem) -> None: self.axis = axis self.cs = cs super().__init__(f"Axis {axis.name} from coordinate system {cs.name} is not mapped to anything") - - -class AxisParsingException(Exception): - pass diff --git a/src/spatialdata/transformations/graph/vert.py b/src/spatialdata/transformations/graph/vert.py index a6231813e..1f9886b7d 100644 --- a/src/spatialdata/transformations/graph/vert.py +++ b/src/spatialdata/transformations/graph/vert.py @@ -6,7 +6,9 @@ import ome_zarr.classes.image as ozi import ome_zarr_models.v06.coordinate_transforms as ozm06ct -from spatialdata._core.transformation_manager.exceptions import AxisParsingException + +class AxisParsingException(Exception): + pass class Axis: From 4f725b1d30194c292d8589e4b156b94d7aae1514 Mon Sep 17 00:00:00 2001 From: Tomaz Vieira Date: Wed, 9 Sep 2026 15:54:26 +0000 Subject: [PATCH 18/26] Fixes docstrings, list exceptions in "Raises" Co-authored-by: Claude Sonnet 5 --- src/spatialdata/transformations/graph/edge.py | 126 +++++++++++++++++- src/spatialdata/transformations/graph/vert.py | 71 +++++++--- 2 files changed, 177 insertions(+), 20 deletions(-) diff --git a/src/spatialdata/transformations/graph/edge.py b/src/spatialdata/transformations/graph/edge.py index 854fe18f5..143469112 100644 --- a/src/spatialdata/transformations/graph/edge.py +++ b/src/spatialdata/transformations/graph/edge.py @@ -55,10 +55,10 @@ def transform_points(self, points: ArrayLike) -> ArrayLike: """ Transform points (coordinates). - Notes - ------- - This function will check if the dimensionality of the input and output coordinate systems of the - transformation are compatible with the given points. + Raises + ------ + UnexpectedShapeError + if `points`'s shape is incompatible with this transformation's input shape """ @abstractmethod @@ -100,6 +100,11 @@ def compose_with(self, transformation: BaseTransfEdge, name: str | None) -> Base Notes ------- Self is applied first, then the transformation passed as argument. + + Raises + ------ + IncompatibleCoordSystemsError + if this transformation's output coordinate system doesn't match `transformation`'s input """ return SequenceEdge(transformations=[self, transformation], name=name) @@ -134,6 +139,12 @@ def __init__( Input coordinate system of the transformation. output Output coordinate system of the transformation. + + Raises + ------ + UnexpectedShapeError + if `linear`'s shape isn't (output.num_axes, input.num_axes) or + if `translation`'s shape isn't (output.num_axes,) """ num_inputs = input.num_axes num_outputs = output.num_axes @@ -189,6 +200,11 @@ def from_affine_matrix( Input coordinate system of the transformation. output Output coordinate system of the transformation. + + Raises + ------ + UnexpectedShapeError + if `affine_matrix`'s shape isn't (output.num_axes + 1, input.num_axes + 1) """ return AffineEdge( linear=affine_matrix[:-1, :-1], @@ -200,6 +216,7 @@ def from_affine_matrix( @classmethod def mapping(cls, input: CoordSystem, output: CoordSystem, name: str | None = None) -> AffineEdge: + """Create an AffineEdge that maps input axes to output axes of the same name.""" linear: ArrayLike = np.zeros((output.num_axes, input.num_axes), dtype=float) for i, des_axis in enumerate(output.axes): for j, src_axis in enumerate(input.axes): @@ -227,6 +244,14 @@ def inverse(self, name: str | None = None) -> BaseTransfEdge | None: ) def transform_points(self, points: ArrayLike) -> ArrayLike: + """ + Transform points (coordinates). + + Raises + ------ + UnexpectedShapeError + if `points`'s shape is incompatible with this transformation's input shape + """ self._validate_transform_points_shapes(points) p = np.vstack([points.T, np.ones(points.shape[0])]) q = self.affine @ p @@ -259,6 +284,11 @@ def __init__( Input coordinate system of the transformation. output Output coordinate system of the transformation. + + Raises + ------ + IncompatibleCoordSystemsError + if `input` and `output` don't have the same number of axes """ if input.num_axes != output.num_axes: raise IncompatibleCoordSystemsError( @@ -270,6 +300,14 @@ def inverse(self, name: str | None = None) -> BaseTransfEdge: return IdentityEdge(input=self.output, output=self.input, name=name) def transform_points(self, points: ArrayLike) -> ArrayLike: + """ + Transform points (coordinates). + + Raises + ------ + UnexpectedShapeError + if `points`'s shape is incompatible with this transformation's input shape + """ self._validate_transform_points_shapes(points) return points @@ -302,6 +340,11 @@ def __init__( output Output coordinate system of the transformation, whose axes must be a shuffling of `input` + + Raises + ------ + IncompatibleCoordSystemsError + if `input` and `output` don't have the same set of axes """ if set(input.axes) != set(output.axes): @@ -325,6 +368,14 @@ def inverse(self, name: str | None = None) -> BaseTransfEdge: ) def transform_points(self, points: ArrayLike) -> ArrayLike: + """ + Transform points (coordinates). + + Raises + ------ + UnexpectedShapeError + if `points`'s shape is incompatible with this transformation's input shape + """ self._validate_transform_points_shapes(points) new_indices = [self.input.axes.index(out_ax) for out_ax in self.output.axes] mapped = points[:, new_indices] @@ -401,6 +452,14 @@ def to_affine(self, name: str | None = None) -> AffineEdge: return AffineEdge(name=name, input=self.input, output=self.output, linear=linear) def transform_points(self, points: ArrayLike) -> ArrayLike: + """ + Transform points (coordinates). + + Raises + ------ + UnexpectedShapeError + if `points`'s shape is incompatible with this transformation's input shape + """ return self.to_affine().transform_points(points) def inverse(self, name: str | None = None) -> BaseTransfEdge | None: @@ -426,6 +485,7 @@ def parse_project_axis( input: CoordSystem, out: CoordSystem | CsGen, ) -> ProjectAxisEdge: + """Parse a `ProjectAxis` NGFF transformation model into a `ProjectAxisEdge`.""" if isinstance(out, CoordSystem): output = out else: @@ -470,11 +530,20 @@ def __init__( ------ IncompatibleCoordSystemsError If the input and output have different number of dimensions + UnexpectedShapeError + if `translation`'s shape doesn't match `input` """ if input.num_axes != output.num_axes: raise IncompatibleCoordSystemsError( input=input, output=output, message="Number of input and output axes must be the same" ) + expected_translation_shape = (input.num_axes,) + if translation.shape != expected_translation_shape: + raise UnexpectedShapeError( + array_name="translation", + array_shape=translation.shape, + expected_shape=expected_translation_shape, + ) self.translation = translation super().__init__(input=input, output=output, name=name) @@ -490,6 +559,14 @@ def inverse(self, name: str | None = None) -> BaseTransfEdge: ) def transform_points(self, points: ArrayLike) -> ArrayLike: + """ + Transform points (coordinates). + + Raises + ------ + UnexpectedShapeError + if `points`'s shape is incompatible with this transformation's input shape + """ self._validate_transform_points_shapes(points) return points + self.translation @@ -556,6 +633,14 @@ def inverse(self, name: str | None = None) -> ScaleEdge | None: ) def transform_points(self, points: ArrayLike) -> ArrayLike: + """ + Transform points (coordinates). + + Raises + ------ + UnexpectedShapeError + if `points`'s shape is incompatible with this transformation's input shape + """ self._validate_transform_points_shapes(points) return points * self.scale @@ -628,6 +713,14 @@ def inverse(self, name: str | None = None) -> BaseTransfEdge: ) def transform_points(self, points: ArrayLike) -> ArrayLike: + """ + Transform points (coordinates). + + Raises + ------ + UnexpectedShapeError + if `points`'s shape is incompatible with this transformation's input shape + """ self._validate_transform_points_shapes(points) res = (self.rotation @ points.T).T assert isinstance(res, np.ndarray) @@ -717,6 +810,14 @@ def to_affine(self, name: str | None = None) -> AffineEdge: ) def transform_points(self, points: ArrayLike) -> ArrayLike: + """ + Transform points (coordinates). + + Raises + ------ + UnexpectedShapeError + if `points`'s shape is incompatible with this transformation's input shape + """ return self.to_affine().transform_points(points) # FIXME @@ -799,6 +900,14 @@ def inverse(self, name: str | None = None) -> BaseTransfEdge | None: ) def transform_points(self, points: ArrayLike) -> ArrayLike: + """ + Transform points (coordinates). + + Raises + ------ + UnexpectedShapeError + if `points`'s shape is incompatible with this transformation's input shape + """ input_axes = self.input.axes_names output_axes = self.output.axes_names self._validate_transform_points_shapes(points) @@ -882,6 +991,7 @@ def parse_identity( input: CoordSystem, out: CoordSystem | CsGen, ) -> IdentityEdge: + """Parse an `Identity` NGFF transformation model into an `IdentityEdge`.""" output = out.generate_like(input) if isinstance(out, CsGen) else out return IdentityEdge(name=model.name, input=input, output=output) @@ -892,6 +1002,7 @@ def parse_translation( input: CoordSystem, out: CoordSystem | CsGen, ) -> TranslationEdge: + """Parse a `Translation` NGFF transformation model into a `TranslationEdge`.""" output = out.generate_like(input) if isinstance(out, CsGen) else out return TranslationEdge( translation=np.asarray(model.translation, dtype=float), @@ -907,6 +1018,7 @@ def parse_scale( input: CoordSystem, out: CoordSystem | CsGen, ) -> ScaleEdge: + """Parse a `Scale` NGFF transformation model into a `ScaleEdge`.""" output = out.generate_like(input) if isinstance(out, CsGen) else out return ScaleEdge( scale=np.asarray(model.scale, dtype=float), @@ -922,6 +1034,7 @@ def parse_map_axis( input: CoordSystem, out: CoordSystem | CsGen, ) -> MapAxisEdge: + """Parse a `MapAxis` NGFF transformation model into a `MapAxisEdge`.""" if isinstance(out, CoordSystem): output = out else: @@ -944,6 +1057,7 @@ def parse_affine( input: CoordSystem, output: CoordSystem | CsGen, ) -> AffineEdge: + """Parse an `Affine` NGFF transformation model into an `AffineEdge`.""" num_output_axes = len(model.affine_matrix) # spec doesn't save last row output = output.generate(num_axes=num_output_axes) if isinstance(output, CsGen) else output affine_array = np.asarray(model.affine_matrix, dtype=float) @@ -958,6 +1072,7 @@ def parse_rotation( input: CoordSystem, out: CoordSystem | CsGen, ) -> RotationEdge: + """Parse a `Rotation` NGFF transformation model into a `RotationEdge`.""" num_output_axes = len(model.rotation_matrix) output = out.generate(num_axes=num_output_axes) if isinstance(out, CsGen) else out return RotationEdge( @@ -974,6 +1089,7 @@ def parse_sequence( input: CoordSystem, output: CoordSystem | CsGen, ) -> SequenceEdge: + """Parse a `Sequence` NGFF transformation model into a `SequenceEdge`.""" parsed_inners: list[BaseTransfEdge] = [] base_name = "intermediate" + ("" if not model.name else f"_for_{model.name}") @@ -1002,6 +1118,7 @@ def parse_by_dimension( input: CoordSystem, output: CoordSystem | CsGen, ) -> ByDimensionEdge: + """Parse a `ByDimension` NGFF transformation model into a `ByDimensionEdge`.""" if not isinstance(output, CoordSystem): max_out_idx = max(ax_idx for t in model.transformations for ax_idx in t.output_axes) output = output.generate(num_axes=max_out_idx + 1) @@ -1038,6 +1155,7 @@ def parse_ngff_transf( model: ozm06trans.AnyTransform, output: CoordSystem | CsGen, ) -> BaseTransfEdge: + """Parse an NGFF coordinate transformation model into a `BaseTransfEdge`""" if isinstance(model, ozm06trans.Identity): return parse_identity(model, input=input, out=output) elif isinstance(model, ozm06trans.Translation): diff --git a/src/spatialdata/transformations/graph/vert.py b/src/spatialdata/transformations/graph/vert.py index 1f9886b7d..177a8d8a9 100644 --- a/src/spatialdata/transformations/graph/vert.py +++ b/src/spatialdata/transformations/graph/vert.py @@ -61,6 +61,15 @@ def __eq__(self, value: object, /) -> bool: @classmethod def try_from_model(cls, model: ozm06ct.Axis) -> Axis: + """ + Parse an `Axis` from an ome-zarr-models axis model. + + Raises + ------ + AxisParsingException + if `model` doesn't have a name, has a type other than "channel" or "space", or has a unit + that isn't a string or None + """ name = model.name if name is None: raise AxisParsingException("Axis doesn't have a name") @@ -80,16 +89,15 @@ class CoordSystemParsingException(Exception): pass +class DuplicateAxisNameError(Exception): + def __init__(self, *, axis_name: str) -> None: + self.axis_name = axis_name + super().__init__(f"Axis name '{axis_name}' is used more than once") + + class CoordSystem: """ - Representation of a coordinate system, following the NGFF specification. - - Parameters - ---------- - name - name of the coordinate system - axes - names of the axes of the coordinate system + Representation of a coordinate system. """ name: Final[str] @@ -100,11 +108,29 @@ class CoordSystem: non-virtual coordinate systems and is usually ignored during serialization""" def __init__(self, name: str, axes: Sequence[Axis], virtual: bool = False): + """ + Parameters + ---------- + name + name of the coordinate system + axes + axes of the coordinate system + virtual + vitual coordinate systems don't serialize to NGFF + + Raises + ------ + DuplicateAxisNameError + if `axes` contains axes with duplicate names + """ self.name = name self.axes = tuple(axes) self.virtual = virtual - if len(self.axes) != len({axis.name for axis in self.axes}): - raise ValueError("Axes names must be unique") + seen_names: set[str] = set() + for axis in self.axes: + if axis.name in seen_names: + raise DuplicateAxisNameError(axis_name=axis.name) + seen_names.add(axis.name) def __repr__(self) -> str: return f"{self.__class__.__name__}({self.name!r}, {self.axes})" @@ -114,18 +140,31 @@ def __hash__(self) -> int: @classmethod def try_from_model(cls, model: ozi.CoordinateSystem) -> CoordSystem: - axes: list[Axis] = [] - for axis in model.axes: - if isinstance(parsed := Axis.try_from_model(axis), Exception): - raise CoordSystemParsingException(parsed) # FIXME - axes.append(parsed) + """ + Parse a `CoordSystem` from an ome-zarr-models coordinate system model. + + Raises + ------ + AxisParsingException + if any axis in `model.axes` fails to parse + DuplicateAxisNameError + if `model.axes` contains axes with duplicate names + """ return CoordSystem( name=model.name, - axes=axes, + axes=[Axis.try_from_model(axis) for axis in model.axes], ) @classmethod def try_from_model_or_default[T](cls, model: ozi.CoordinateSystem | None, *, default: T) -> CoordSystem | T: + """ + Parse a `CoordSystem` from `model`, or return `default` if `model` is None. + + Raises + ------ + AxisParsingException + if `model` is not None and any of its axes fails to parse + """ if model is not None: return CoordSystem.try_from_model(model) return default From 9591fa2cbd132315b002b977e0481cd7f85997aa Mon Sep 17 00:00:00 2001 From: Tomaz Vieira Date: Wed, 9 Sep 2026 16:03:39 +0000 Subject: [PATCH 19/26] Incorporates parse_project_axis --- src/spatialdata/transformations/graph/edge.py | 48 +++++++++---------- 1 file changed, 24 insertions(+), 24 deletions(-) diff --git a/src/spatialdata/transformations/graph/edge.py b/src/spatialdata/transformations/graph/edge.py index 143469112..d21fbef7d 100644 --- a/src/spatialdata/transformations/graph/edge.py +++ b/src/spatialdata/transformations/graph/edge.py @@ -479,30 +479,6 @@ def inverse(self, name: str | None = None) -> BaseTransfEdge | None: ) -def parse_project_axis( - model: ozm06trans.ProjectAxis, - *, - input: CoordSystem, - out: CoordSystem | CsGen, -) -> ProjectAxisEdge: - """Parse a `ProjectAxis` NGFF transformation model into a `ProjectAxisEdge`.""" - if isinstance(out, CoordSystem): - output = out - else: - num_dropped_inputs = len(model.droppedInputs or ()) - num_created_outputs = len(model.createdOutputs or ()) - num_output_axes = input.num_axes - num_dropped_inputs + num_created_outputs - output = out.generate(num_axes=num_output_axes) - - return ProjectAxisEdge( - created_outputs={output.axes[i] for i in model.createdOutputs or ()}, - dropped_inputs={input.axes[i] for i in model.droppedInputs or ()}, - input=input, - output=output, - name=model.name, - ) - - class TranslationEdge(BaseTransfEdge): """The Translation transformation from the NGFF specification.""" @@ -1012,6 +988,28 @@ def parse_translation( ) +def parse_project_axis( + model: ozm06trans.ProjectAxis, + *, + input: CoordSystem, + output: CoordSystem | CsGen, +) -> ProjectAxisEdge: + """Parse a `ProjectAxis` NGFF transformation model into a `ProjectAxisEdge`.""" + if isinstance(output, CsGen): + num_dropped_inputs = len(model.droppedInputs or ()) + num_created_outputs = len(model.createdOutputs or ()) + num_output_axes = input.num_axes - num_dropped_inputs + num_created_outputs + output = output.generate(num_axes=num_output_axes) + + return ProjectAxisEdge( + created_outputs={output.axes[i] for i in model.createdOutputs or ()}, + dropped_inputs={input.axes[i] for i in model.droppedInputs or ()}, + input=input, + output=output, + name=model.name, + ) + + def parse_scale( model: ozm06trans.Scale, *, @@ -1172,5 +1170,7 @@ def parse_ngff_transf( return parse_sequence(model, input=input, output=output) elif isinstance(model, ozm06trans.ByDimension): return parse_by_dimension(model, input=input, output=output) + elif isinstance(model, ozm06trans.ProjectAxis): + return parse_project_axis(model, input=input, output=output) else: raise NotImplementedError(f"Unsupported transformation: {model.type}") From 17eae96ad808beb4d7791aecfdd0e840ab504830 Mon Sep 17 00:00:00 2001 From: Tomaz Vieira Date: Wed, 9 Sep 2026 16:18:46 +0000 Subject: [PATCH 20/26] Fixes parsing translatiion model getting bad name Co-authored-by: Claude Sonnet 5 --- src/spatialdata/transformations/graph/edge.py | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/src/spatialdata/transformations/graph/edge.py b/src/spatialdata/transformations/graph/edge.py index d21fbef7d..74a277aca 100644 --- a/src/spatialdata/transformations/graph/edge.py +++ b/src/spatialdata/transformations/graph/edge.py @@ -984,7 +984,7 @@ def parse_translation( translation=np.asarray(model.translation, dtype=float), input=input, output=output, - name=input.name, + name=model.name, ) From 3a42bef804ee1e1a43628170c72c49eb8994de97 Mon Sep 17 00:00:00 2001 From: Tomaz Vieira Date: Wed, 9 Sep 2026 16:29:37 +0000 Subject: [PATCH 21/26] Adds test for parsing *Edge transforms from ngff Co-authored-by: Claude Sonnet 5 --- .../transformations/edge/test_edge_parsing.py | 226 ++++++++++++++++++ 1 file changed, 226 insertions(+) create mode 100644 tests/transformations/edge/test_edge_parsing.py diff --git a/tests/transformations/edge/test_edge_parsing.py b/tests/transformations/edge/test_edge_parsing.py new file mode 100644 index 000000000..5f95a4d61 --- /dev/null +++ b/tests/transformations/edge/test_edge_parsing.py @@ -0,0 +1,226 @@ +# pyright: strict + +from __future__ import annotations + +import numpy as np +import ome_zarr_models.v06.coordinate_transforms as ozm06trans + +from spatialdata.transformations.graph.edge import ( + AffineEdge, + ByDimensionEdge, + CsGen, + MapAxisEdge, + RotationEdge, + ScaleEdge, + SequenceEdge, + TranslationEdge, + parse_affine, + parse_by_dimension, + parse_identity, + parse_map_axis, + parse_project_axis, + parse_rotation, + parse_scale, + parse_sequence, + parse_translation, +) +from tests.transformations.edge.conftest import ( + xy_cs, + xyz_cs, + yx_cs, + z_axis, +) + +# The ome-zarr-models transform models used below all leave `input`/`output` at their +# default of None: those fields identify coordinate systems by name/path in the NGFF +# metadata and aren't consumed by the parse_* functions, which instead take the already +# resolved `CoordSystem` (or a `CsGen`) as separate arguments. Only the transform-specific +# fields (e.g. `translation`, `scale`, `name`) are exercised here. + +POINTS_2D = np.array([[1.0, 2.0], [3.0, 4.0]]) + + +class TestParseIdentity: + def test_parses_name(self): + model = ozm06trans.Identity(name="my-identity") + edge = parse_identity(model, input=xy_cs, out=xy_cs) + assert edge.name == "my-identity" + assert edge.input is xy_cs + assert edge.output is xy_cs + + def test_csgen_output_mirrors_input(self): + model = ozm06trans.Identity(name="my-identity") + edge = parse_identity(model, input=xy_cs, out=CsGen(base_name="gen")) + assert edge.output.axes_names == xy_cs.axes_names + assert edge.output.virtual is True + + +class TestParseTranslation: + def test_parses_translation_vector(self): + model = ozm06trans.Translation(name="my-translation", translation=(1.0, 2.0)) + edge = parse_translation(model, input=xy_cs, out=xy_cs) + np.testing.assert_allclose(edge.translation, [1.0, 2.0]) + assert edge.name == model.name + assert edge.input is xy_cs + assert edge.output is xy_cs + + def test_csgen_output_mirrors_input(self): + model = ozm06trans.Translation(translation=(1.0, 2.0)) + edge = parse_translation(model, input=xy_cs, out=CsGen(base_name="gen")) + assert edge.output.axes_names == xy_cs.axes_names + assert edge.output.virtual is True + + +class TestParseProjectAxis: + def test_parses_dropped_inputs(self): + model = ozm06trans.ProjectAxis(name="drop-z", droppedInputs=(2,)) + edge = parse_project_axis(model, input=xyz_cs, output=xy_cs) + assert edge.dropped_inputs == {z_axis} + assert edge.created_outputs == set() + assert edge.name == "drop-z" + + def test_parses_created_outputs(self): + model = ozm06trans.ProjectAxis(name="create-z", createdOutputs=(2,)) + edge = parse_project_axis(model, input=xy_cs, output=xyz_cs) + assert edge.created_outputs == {z_axis} + assert edge.dropped_inputs == set() + + def test_csgen_output_computes_dimensionality(self): + model = ozm06trans.ProjectAxis(droppedInputs=(2,)) + edge = parse_project_axis(model, input=xyz_cs, output=CsGen(base_name="gen")) + assert edge.output.num_axes == 2 + assert edge.dropped_inputs == {z_axis} + + +class TestParseScale: + def test_parses_scale_vector_and_name(self): + model = ozm06trans.Scale(name="my-scale", scale=(2.0, 3.0)) + edge = parse_scale(model, input=xy_cs, out=xy_cs) + assert isinstance(edge, ScaleEdge) + np.testing.assert_allclose(edge.scale, [2.0, 3.0]) + assert edge.name == "my-scale" + + def test_csgen_output_mirrors_input(self): + model = ozm06trans.Scale(scale=(2.0, 3.0)) + edge = parse_scale(model, input=xy_cs, out=CsGen(base_name="gen")) + assert edge.output.axes_names == xy_cs.axes_names + assert edge.output.virtual is True + + +class TestParseMapAxis: + def test_parses_with_explicit_output(self): + model = ozm06trans.MapAxis(name="my-map-axis", mapAxis=(1, 0)) + edge = parse_map_axis(model, input=xy_cs, out=yx_cs) + assert isinstance(edge, MapAxisEdge) + assert edge.name == "my-map-axis" + assert edge.output is yx_cs + np.testing.assert_allclose(edge.transform_points(POINTS_2D), np.array([[2.0, 1.0], [4.0, 3.0]])) + + def test_csgen_output_is_permuted_input_axes(self): + model = ozm06trans.MapAxis(mapAxis=(1, 0)) + edge = parse_map_axis(model, input=xy_cs, out=CsGen(base_name="gen")) + assert edge.output.axes_names == yx_cs.axes_names + assert edge.output.virtual is True + + +class TestParseAffine: + def test_parses_linear_translation_and_name(self): + model = ozm06trans.Affine(name="my-affine", affine=((2.0, 0.0, 1.0), (0.0, 3.0, 5.0))) + edge = parse_affine(model, input=xy_cs, output=xy_cs) + assert isinstance(edge, AffineEdge) + np.testing.assert_allclose(edge.linear, [[2.0, 0.0], [0.0, 3.0]]) + np.testing.assert_allclose(edge.translation, [1.0, 5.0]) + assert edge.name == "my-affine" + + def test_csgen_output_derived_from_matrix_row_count(self): + model = ozm06trans.Affine(affine=((2.0, 0.0, 1.0), (0.0, 3.0, 5.0))) + edge = parse_affine(model, input=xy_cs, output=CsGen(base_name="gen")) + assert edge.output.num_axes == 2 + assert edge.output.virtual is True + + +class TestParseRotation: + def test_parses_rotation_matrix_and_name(self): + model = ozm06trans.Rotation(name="my-rotation", rotation=((0.0, -1.0), (1.0, 0.0))) + edge = parse_rotation(model, input=xy_cs, out=xy_cs) + assert isinstance(edge, RotationEdge) + np.testing.assert_allclose(edge.rotation, [[0.0, -1.0], [1.0, 0.0]]) + assert edge.name == "my-rotation" + + def test_csgen_output_derived_from_matrix_row_count(self): + model = ozm06trans.Rotation(rotation=((0.0, -1.0), (1.0, 0.0))) + edge = parse_rotation(model, input=xy_cs, out=CsGen(base_name="gen")) + assert edge.output.num_axes == 2 + assert edge.output.virtual is True + + +class TestParseSequence: + def test_parses_name_and_composes_inner_transformations_in_order(self): + inner_translation = ozm06trans.Translation(translation=(1.0, 2.0)) + inner_scale = ozm06trans.Scale(scale=(2.0, 3.0)) + model = ozm06trans.Sequence(name="my-sequence", transformations=(inner_translation, inner_scale)) + edge = parse_sequence(model, input=xy_cs, output=xy_cs) + assert isinstance(edge, SequenceEdge) + assert edge.name == "my-sequence" + assert edge.input is xy_cs + assert edge.output is xy_cs + assert [type(t) for t in edge.transformations] == [TranslationEdge, ScaleEdge] + expected = (POINTS_2D + np.array([1.0, 2.0])) * np.array([2.0, 3.0]) + np.testing.assert_allclose(edge.transform_points(POINTS_2D), expected) + + def test_single_inner_transformation_uses_given_output_directly(self): + inner_translation = ozm06trans.Translation(translation=(1.0, 2.0)) + model = ozm06trans.Sequence(transformations=(inner_translation,)) + edge = parse_sequence(model, input=xy_cs, output=xy_cs) + assert len(edge.transformations) == 1 + assert edge.transformations[0].output is xy_cs + + def test_multiple_inner_transformations_use_virtual_intermediate_coord_systems(self): + inner_translation = ozm06trans.Translation(translation=(1.0, 2.0)) + inner_scale = ozm06trans.Scale(scale=(2.0, 3.0)) + model = ozm06trans.Sequence(transformations=(inner_translation, inner_scale)) + edge = parse_sequence(model, input=xy_cs, output=xy_cs) + intermediate = edge.transformations[0].output + assert intermediate is not xy_cs + assert intermediate.virtual is True + assert intermediate.axes_names == xy_cs.axes_names + assert edge.transformations[1].output is xy_cs + + def test_csgen_output_mirrors_input(self): + inner_translation = ozm06trans.Translation(translation=(1.0, 2.0)) + inner_scale = ozm06trans.Scale(scale=(2.0, 3.0)) + model = ozm06trans.Sequence(transformations=(inner_translation, inner_scale)) + edge = parse_sequence(model, input=xy_cs, output=CsGen(base_name="gen")) + assert edge.output.axes_names == xy_cs.axes_names + assert edge.output.virtual is True + + +class TestParseByDimension: + def test_parses_name_and_splits_axes_between_inner_transformations(self): + scale_x = ozm06trans.Scale(scale=(2.0,)) + translate_y = ozm06trans.Translation(translation=(5.0,)) + model = ozm06trans.ByDimension( + name="my-by-dim", + transformations=( + ozm06trans.ByDimensionTransform(transformation=scale_x, input_axes=(0,), output_axes=(0,)), + ozm06trans.ByDimensionTransform(transformation=translate_y, input_axes=(1,), output_axes=(1,)), + ), + ) + edge = parse_by_dimension(model, input=xy_cs, output=xy_cs) + assert isinstance(edge, ByDimensionEdge) + assert edge.name == "my-by-dim" + assert [type(t) for t in edge.transformations] == [ScaleEdge, TranslationEdge] + np.testing.assert_allclose(edge.transform_points(POINTS_2D), np.array([[2.0, 7.0], [6.0, 9.0]])) + + def test_csgen_output_derived_from_max_output_axis_index(self): + scale_x = ozm06trans.Scale(scale=(2.0,)) + translate_y = ozm06trans.Translation(translation=(5.0,)) + model = ozm06trans.ByDimension( + transformations=( + ozm06trans.ByDimensionTransform(transformation=scale_x, input_axes=(0,), output_axes=(0,)), + ozm06trans.ByDimensionTransform(transformation=translate_y, input_axes=(1,), output_axes=(1,)), + ), + ) + edge = parse_by_dimension(model, input=xy_cs, output=CsGen(base_name="gen")) + assert edge.output.num_axes == 2 + assert edge.output.virtual is True From 2849a29672459e7dd630a1d7cf09a483313b9764 Mon Sep 17 00:00:00 2001 From: Tomaz Vieira Date: Fri, 11 Sep 2026 08:17:09 +0000 Subject: [PATCH 22/26] Fixes dependencies on ome-zarr and ome-zarr-models --- pyproject.toml | 3 ++- 1 file changed, 2 insertions(+), 1 deletion(-) diff --git a/pyproject.toml b/pyproject.toml index 400d884b1..04947e348 100644 --- a/pyproject.toml +++ b/pyproject.toml @@ -39,7 +39,8 @@ dependencies = [ "networkx", "numba>=0.55", "numpy", - "ome-zarr>=0.18.0", + "ome-zarr>=0.19.2", + "ome-zarr-models>=1.8", "pandas", "pooch", "pyarrow", From 4c6c8d90993d25c78c15e0d7a05af5b9457ab07a Mon Sep 17 00:00:00 2001 From: Tomaz Vieira Date: Fri, 11 Sep 2026 09:57:01 +0000 Subject: [PATCH 23/26] Cleans up transform_points for rotation and affine --- src/spatialdata/_io/io_raster.py | 16 +- src/spatialdata/transformations/graph/edge.py | 147 ++++++++++-------- src/spatialdata/transformations/graph/vert.py | 7 +- 3 files changed, 97 insertions(+), 73 deletions(-) diff --git a/src/spatialdata/_io/io_raster.py b/src/spatialdata/_io/io_raster.py index 655c1da9a..0b336999d 100644 --- a/src/spatialdata/_io/io_raster.py +++ b/src/spatialdata/_io/io_raster.py @@ -34,7 +34,7 @@ RasterFormatType, get_ome_zarr_format, ) -from spatialdata._types import ELEMENT_TYPE, ELEMENT_TYPE_RASTER, GROUP_NAME +from spatialdata._types import ELEMENT_TYPE, ELEMENT_TYPE_RASTER from spatialdata._utils import get_pyramid_levels from spatialdata.models.models import ATTRS_KEY from spatialdata.models.pyramids_utils import dask_arrays_to_datatree @@ -44,7 +44,7 @@ _set_transformations, compute_coordinates, ) -from spatialdata.transformations.graph.edge import BaseTransfEdge, parse_ngff_transf +from spatialdata.transformations.graph.edge import BaseTransformationEdge, parse_ngff_transf from spatialdata.transformations.graph.vert import Axis, CoordSystem @@ -168,20 +168,20 @@ def _prepare_storage_options( return prepared_options -def try_read_ngff06_multiscale(store: Path) -> tuple[DataTree, Sequence[BaseTransfEdge]]: +def try_read_ngff06_multiscale(store: Path) -> tuple[DataTree, Sequence[BaseTransformationEdge]]: multiscale = oz.OMEZarrMultiscale.from_ome_zarr(str(store)) assert isinstance(multiscale, oz.OMEZarrMultiscale) # disambiguate from OMEZarrLabel return try_parse_ngff06_multiscale(multiscale) -def try_parse_ngff06_multiscale(multiscale: oz.OMEZarrMultiscale) -> tuple[DataTree, Sequence[BaseTransfEdge]]: +def try_parse_ngff06_multiscale(multiscale: oz.OMEZarrMultiscale) -> tuple[DataTree, Sequence[BaseTransformationEdge]]: """Parse an OMEZarMultiscale into a DataTree and collects Multiscale-level transforms.""" name_to_cs: dict[str, CoordSystem] = {} for cs in multiscale.metadata.coordinateSystems or (): parsed_cs = CoordSystem.try_from_model(cs) name_to_cs[cs.name] = parsed_cs - parsed_transfs: list[BaseTransfEdge] = [] + parsed_transfs: list[BaseTransformationEdge] = [] for transf in multiscale.metadata.coordinateTransformations or (): in_cs_id = transf.input out_cs_ref = transf.output @@ -364,9 +364,7 @@ def _get_multiscale_nodes(image_nodes: list[Node], nodes: list[Node]) -> list[No return nodes -def _get_raster_element_group( - raster_type: ELEMENT_TYPE_RASTER, group: zarr.Group, element_name: str -) -> zarr.Group: +def _get_raster_element_group(raster_type: ELEMENT_TYPE_RASTER, group: zarr.Group, element_name: str) -> zarr.Group: """Get the Zarr group holding a raster element that has just been written. Labels are nested one level deeper than images: ome-zarr writes them inside a "labels" group, so for them the @@ -385,7 +383,7 @@ def _get_raster_element_group( ------- The Zarr group of the raster element. """ - if raster_type != "labels": + if raster_type != ELEMENT_TYPE.LABELS: return group labels_group = group["labels"] if not isinstance(labels_group, zarr.Group): diff --git a/src/spatialdata/transformations/graph/edge.py b/src/spatialdata/transformations/graph/edge.py index 74a277aca..2201b257b 100644 --- a/src/spatialdata/transformations/graph/edge.py +++ b/src/spatialdata/transformations/graph/edge.py @@ -3,7 +3,7 @@ from __future__ import annotations from abc import ABC, abstractmethod -from collections.abc import Sequence +from collections.abc import Mapping, Sequence from typing import Final import numpy as np @@ -22,7 +22,7 @@ from spatialdata.transformations.graph.vert import Axis, CoordSystem -class BaseTransfEdge(ABC): +class BaseTransformationEdge(ABC): """Base class for all the transformations defined by the NGFF specification.""" input: Final[CoordSystem] @@ -47,7 +47,7 @@ def __repr__(self) -> str: return f"{type(self).__name__} ({domain} -> {codomain})" @abstractmethod - def inverse(self, name: str | None = None) -> BaseTransfEdge | None: + def inverse(self, name: str | None = None) -> BaseTransformationEdge | None: """Return the inverse of the transformation if it exists""" @abstractmethod @@ -84,7 +84,7 @@ def _validate_transform_points_shapes(self, points: ArrayLike) -> None: ) # order of the composition: self is applied first, then the transformation passed as argument - def compose_with(self, transformation: BaseTransfEdge, name: str | None) -> BaseTransfEdge: + def compose_with(self, transformation: BaseTransformationEdge, name: str | None) -> BaseTransformationEdge: """ Compose the transfomation object with another transformation @@ -109,7 +109,7 @@ def compose_with(self, transformation: BaseTransfEdge, name: str | None) -> Base return SequenceEdge(transformations=[self, transformation], name=name) -class AffineEdge(BaseTransfEdge): +class AffineEdge(BaseTransformationEdge): """The Affine transformation from the NGFF specification.""" linear: Final[ArrayLike] @@ -229,7 +229,7 @@ def mapping(cls, input: CoordSystem, output: CoordSystem, name: str | None = Non name=name, ) - def inverse(self, name: str | None = None) -> BaseTransfEdge | None: + def inverse(self, name: str | None = None) -> BaseTransformationEdge | None: try: # FIXME: I think there are more efficient/precise ways to invert a matrix inv = np.linalg.inv(self.affine) @@ -253,11 +253,7 @@ def transform_points(self, points: ArrayLike) -> ArrayLike: if `points`'s shape is incompatible with this transformation's input shape """ self._validate_transform_points_shapes(points) - p = np.vstack([points.T, np.ones(points.shape[0])]) - q = self.affine @ p - res = q[: self.output.num_axes, :].T - assert isinstance(res, np.ndarray) - return res + return (points @ self.linear.T) + self.translation def to_affine(self, name: str | None = None) -> AffineEdge: return AffineEdge( @@ -265,7 +261,7 @@ def to_affine(self, name: str | None = None) -> AffineEdge: ) -class IdentityEdge(BaseTransfEdge): +class IdentityEdge(BaseTransformationEdge): """The Identity transformation from the NGFF specification.""" def __init__( @@ -296,7 +292,7 @@ def __init__( ) super().__init__(input=input, output=output, name=name) - def inverse(self, name: str | None = None) -> BaseTransfEdge: + def inverse(self, name: str | None = None) -> BaseTransformationEdge: return IdentityEdge(input=self.output, output=self.input, name=name) def transform_points(self, points: ArrayLike) -> ArrayLike: @@ -313,14 +309,14 @@ def transform_points(self, points: ArrayLike) -> ArrayLike: def to_affine(self, name: str | None = None) -> AffineEdge: return AffineEdge( - linear=np.eye(self.input.num_axes), + linear=np.identity(self.input.num_axes), input=self.input, output=self.output, name=name, ) -class MapAxisEdge(BaseTransfEdge): +class MapAxisEdge(BaseTransformationEdge): """The MapAxis transformation from the NGFF specification.""" def __init__( @@ -329,6 +325,7 @@ def __init__( name: str | None = None, input: CoordSystem, output: CoordSystem, + input_to_output: Mapping[Axis, Axis], ) -> None: """ Parameters @@ -338,32 +335,52 @@ def __init__( input Input coordinate system of the transformation. output - Output coordinate system of the transformation, whose axes - must be a shuffling of `input` + Output coordinate system of the transformation. Must + have the same number of axes as `input` + input_to_output + A mapping from Raises ------ IncompatibleCoordSystemsError - if `input` and `output` don't have the same set of axes + if `input` and `output` don't have the number of axes """ - if set(input.axes) != set(output.axes): + if input.num_axes != output.num_axes: raise IncompatibleCoordSystemsError( - input=input, output=output, message="Input and output must have the same axes" + input=input, output=output, message="Input and output must have the same number of axes" ) + + unmapped_inputs = set(input.axes) + unmapped_outputs = set(output.axes) + for inp_ax, out_ax in input_to_output.items(): + try: + unmapped_inputs.remove(inp_ax) + except KeyError as e: + raise MissingAxisError(axis=inp_ax, cs=self.input) from e + try: + unmapped_outputs.remove(out_ax) + except KeyError as e: + raise MissingAxisError(axis=out_ax, cs=self.output) from e + + for ax in unmapped_inputs: + raise UnmappedAxisError(axis=ax, cs=self.input) + for ax in unmapped_outputs: + raise UnmappedAxisError(axis=ax, cs=self.output) + + self.input_to_output = input_to_output super().__init__(input=input, output=output, name=name) def __repr__(self) -> str: s = super().__repr__() + "\n" - s += "\n".join( - f" {out.name} <- {inp.name}\n" for out, inp in zip(self.output.axes, self.input.axes, strict=True) - ) + s += "\n".join(f" {out.name} <- {inp.name}\n" for out, inp in self.input_to_output.items()) return s - def inverse(self, name: str | None = None) -> BaseTransfEdge: + def inverse(self, name: str | None = None) -> BaseTransformationEdge: return MapAxisEdge( input=self.output, output=self.input, + input_to_output={out: inp for inp, out in self.input_to_output.items()}, name=name, ) @@ -377,16 +394,23 @@ def transform_points(self, points: ArrayLike) -> ArrayLike: if `points`'s shape is incompatible with this transformation's input shape """ self._validate_transform_points_shapes(points) - new_indices = [self.input.axes.index(out_ax) for out_ax in self.output.axes] + output_to_input = {out_ax: in_ax for in_ax, out_ax in self.input_to_output.items()} + new_indices = [self.input.axes.index(output_to_input[out_ax]) for out_ax in self.output.axes] + assert len(new_indices) == self.output.num_axes mapped = points[:, new_indices] assert isinstance(mapped, np.ndarray) return mapped def to_affine(self, name: str | None = None) -> AffineEdge: - return AffineEdge.mapping(input=self.input, output=self.output, name=name) + linear = np.zeros((self.output.num_axes, self.input.num_axes)) + for inp_ax, out_ax in self.input_to_output.items(): + inp_idx = self.input.axes.index(inp_ax) + out_idx = self.output.axes.index(out_ax) + linear[out_idx, inp_idx] = 1 + return AffineEdge(linear=linear, input=self.input, output=self.output) -class ProjectAxisEdge(BaseTransfEdge): +class ProjectAxisEdge(BaseTransformationEdge): dropped_inputs: Final[set[Axis]] created_outputs: Final[set[Axis]] @@ -462,7 +486,7 @@ def transform_points(self, points: ArrayLike) -> ArrayLike: """ return self.to_affine().transform_points(points) - def inverse(self, name: str | None = None) -> BaseTransfEdge | None: + def inverse(self, name: str | None = None) -> BaseTransformationEdge | None: # FIXME: there may be other cases where this is invertible if self.input.num_axes != self.output.num_axes: return None @@ -479,7 +503,7 @@ def inverse(self, name: str | None = None) -> BaseTransfEdge | None: ) -class TranslationEdge(BaseTransfEdge): +class TranslationEdge(BaseTransformationEdge): """The Translation transformation from the NGFF specification.""" def __init__( @@ -526,7 +550,7 @@ def __init__( def __repr__(self) -> str: return super().__repr__() + str(self.translation) - def inverse(self, name: str | None = None) -> BaseTransfEdge: + def inverse(self, name: str | None = None) -> BaseTransformationEdge: return TranslationEdge( translation=-self.translation, input=self.output, @@ -556,7 +580,7 @@ def to_affine(self, name: str | None = None) -> AffineEdge: ) -class ScaleEdge(BaseTransfEdge): +class ScaleEdge(BaseTransformationEdge): """The Scale transformation from the NGFF specification.""" def __init__( @@ -629,7 +653,7 @@ def to_affine(self, name: str | None = None) -> AffineEdge: ) -class RotationEdge(BaseTransfEdge): +class RotationEdge(BaseTransformationEdge): """The Rotation transformation from the NGFF specification.""" rotation: Final[ArrayLike] @@ -646,7 +670,7 @@ def __init__( Parameters ---------- linear_matrix - an array of shape (output.num_axes, input.num_axes) representing the rotation + And orthonormal matrix of shape (output.num_axes, input.num_axes) input Input coordinate system of the transformation. output @@ -657,8 +681,10 @@ def __init__( if linear_matrix's shape isn't (output.num_axes, input.num_axes) IncompatibleCoordSystemsError if input and output don't have the same number of axes - NotUnimodularError + DeterminantDifferentFromOne if linear_matrix doesn't have determinant ~= 1 + NotOrthonormal + if linear_matrix is not orthonormal """ if input.num_axes != output.num_axes: raise IncompatibleCoordSystemsError( @@ -680,7 +706,7 @@ def __repr__(self) -> str: s += "\n".join(str(row) for row in self.rotation) return s - def inverse(self, name: str | None = None) -> BaseTransfEdge: + def inverse(self, name: str | None = None) -> BaseTransformationEdge: return RotationEdge( linear_matrix=self.rotation.T, input=self.output, @@ -698,9 +724,7 @@ def transform_points(self, points: ArrayLike) -> ArrayLike: if `points`'s shape is incompatible with this transformation's input shape """ self._validate_transform_points_shapes(points) - res = (self.rotation @ points.T).T - assert isinstance(res, np.ndarray) - return res + return points @ self.rotation.T def to_affine(self, name: str | None = None) -> AffineEdge: return AffineEdge( @@ -711,14 +735,14 @@ def to_affine(self, name: str | None = None) -> AffineEdge: ) -class SequenceEdge(BaseTransfEdge): +class SequenceEdge(BaseTransformationEdge): """The Sequence transformation from the NGFF specification.""" def __init__( self, *, name: str | None = None, - transformations: Sequence[BaseTransfEdge], + transformations: Sequence[BaseTransformationEdge], ) -> None: """ Init the NgffSequence object. @@ -765,7 +789,7 @@ def __repr__(self) -> str: return out def inverse(self, name: str | None = None) -> SequenceEdge | None: - inverted: list[BaseTransfEdge] = [] + inverted: list[BaseTransformationEdge] = [] for t in reversed(self.transformations): inv = t.inverse() if inv is None: @@ -797,16 +821,16 @@ def transform_points(self, points: ArrayLike) -> ArrayLike: return self.to_affine().transform_points(points) # FIXME -class ByDimensionEdge(BaseTransfEdge): +class ByDimensionEdge(BaseTransformationEdge): """The ByDimension transformation from the NGFF specification.""" - transformations: Final[Sequence[BaseTransfEdge]] + transformations: Final[Sequence[BaseTransformationEdge]] def __init__( self, *, name: str | None = None, - transformations: Sequence[BaseTransfEdge], + transformations: Sequence[BaseTransformationEdge], input: CoordSystem, output: CoordSystem, ) -> None: @@ -861,8 +885,8 @@ def __repr__(self) -> str: out += "]" return out - def inverse(self, name: str | None = None) -> BaseTransfEdge | None: - inverse_transformations: list[BaseTransfEdge] = [] + def inverse(self, name: str | None = None) -> BaseTransformationEdge | None: + inverse_transformations: list[BaseTransformationEdge] = [] for t in self.transformations: inv = t.inverse() if inv is None: @@ -931,13 +955,13 @@ def __init__(self, base_name: str): def generate(self, *, num_axes: int) -> CoordSystem: out = CoordSystem( name=f"{self._base_name}{self._cs_count}", - axes=[ + axes=tuple( Axis( name=f"axis_{ax_idx}", type="space", # FIXME ) for ax_idx in range(num_axes) - ], + ), virtual=True, ) self._cs_count += 1 @@ -946,7 +970,7 @@ def generate(self, *, num_axes: int) -> CoordSystem: def generate_like(self, other: CoordSystem) -> CoordSystem: out = CoordSystem( name=f"{self._base_name}{self._cs_count}", - axes=[ + axes=tuple( Axis( name=axis.name, type=axis.type, @@ -954,7 +978,7 @@ def generate_like(self, other: CoordSystem) -> CoordSystem: long_name=axis.long_name, ) for axis in other.axes - ], + ), virtual=True, ) self._cs_count += 1 @@ -1039,14 +1063,13 @@ def parse_map_axis( dummy_cs = out.generate(num_axes=len(model.mapAxis)) output = CoordSystem( name=dummy_cs.name, - axes=[input.axes[i] for i in model.mapAxis], + axes=tuple(input.axes[i] for i in model.mapAxis), virtual=True, ) - return MapAxisEdge( - input=input, - output=output, - name=model.name, - ) + input_to_output = { + input.axes[inp_idx]: output_axis for inp_idx, output_axis in zip(model.mapAxis, output.axes, strict=True) + } + return MapAxisEdge(input=input, output=output, name=model.name, input_to_output=input_to_output) def parse_affine( @@ -1088,7 +1111,7 @@ def parse_sequence( output: CoordSystem | CsGen, ) -> SequenceEdge: """Parse a `Sequence` NGFF transformation model into a `SequenceEdge`.""" - parsed_inners: list[BaseTransfEdge] = [] + parsed_inners: list[BaseTransformationEdge] = [] base_name = "intermediate" + ("" if not model.name else f"_for_{model.name}") cs_gen: CsGen = output if isinstance(output, CsGen) else CsGen(base_name=base_name) @@ -1118,19 +1141,19 @@ def parse_by_dimension( ) -> ByDimensionEdge: """Parse a `ByDimension` NGFF transformation model into a `ByDimensionEdge`.""" if not isinstance(output, CoordSystem): - max_out_idx = max(ax_idx for t in model.transformations for ax_idx in t.output_axes) + max_out_idx = max(ax_idx for t in model.transformations for ax_idx in t.outputAxes) output = output.generate(num_axes=max_out_idx + 1) - piecewise_transforms: list[BaseTransfEdge] = [] + piecewise_transforms: list[BaseTransformationEdge] = [] for t in model.transformations: - inp_axes = [input.axes[i] for i in t.input_axes] + inp_axes = tuple(input.axes[i] for i in t.inputAxes) partial_input = CoordSystem( axes=inp_axes, name=f"{input.name}_{','.join(ax.name for ax in inp_axes)}", virtual=True, ) - out_axes = [output.axes[i] for i in t.output_axes] + out_axes = tuple(output.axes[i] for i in t.outputAxes) partial_out = CoordSystem( axes=out_axes, name=f"{output.name}_{','.join(ax.name for ax in inp_axes)}", @@ -1152,7 +1175,7 @@ def parse_ngff_transf( input: CoordSystem, model: ozm06trans.AnyTransform, output: CoordSystem | CsGen, -) -> BaseTransfEdge: +) -> BaseTransformationEdge: """Parse an NGFF coordinate transformation model into a `BaseTransfEdge`""" if isinstance(model, ozm06trans.Identity): return parse_identity(model, input=input, out=output) diff --git a/src/spatialdata/transformations/graph/vert.py b/src/spatialdata/transformations/graph/vert.py index 177a8d8a9..a78046b90 100644 --- a/src/spatialdata/transformations/graph/vert.py +++ b/src/spatialdata/transformations/graph/vert.py @@ -70,16 +70,19 @@ def try_from_model(cls, model: ozm06ct.Axis) -> Axis: if `model` doesn't have a name, has a type other than "channel" or "space", or has a unit that isn't a string or None """ + name = model.name if name is None: raise AxisParsingException("Axis doesn't have a name") - if model.type != "channel" and model.type != "space": + if model.type is None: + raise AxisParsingException("Axis doesn't have a type") + if model.type not in ("channel", "space"): raise AxisParsingException(f"Can't handle axis of type {model.type}") if not isinstance(model.unit, (str, type(None))): raise AxisParsingException("Can't handle axis unit") return Axis( name=name, - type=model.type, + type=model.type, # type: ignore[arg-type] unit=model.unit, long_name=model.longName, ) From 0c6d7b12b76eaae238d36e1a113dbf03659480a7 Mon Sep 17 00:00:00 2001 From: Tomaz Vieira Date: Fri, 11 Sep 2026 10:15:45 +0000 Subject: [PATCH 24/26] Fixes renamed fields in ome-zarr-models --- tests/transformations/edge/test_edge_parsing.py | 8 ++++---- 1 file changed, 4 insertions(+), 4 deletions(-) diff --git a/tests/transformations/edge/test_edge_parsing.py b/tests/transformations/edge/test_edge_parsing.py index 5f95a4d61..a20acbaa3 100644 --- a/tests/transformations/edge/test_edge_parsing.py +++ b/tests/transformations/edge/test_edge_parsing.py @@ -202,8 +202,8 @@ def test_parses_name_and_splits_axes_between_inner_transformations(self): model = ozm06trans.ByDimension( name="my-by-dim", transformations=( - ozm06trans.ByDimensionTransform(transformation=scale_x, input_axes=(0,), output_axes=(0,)), - ozm06trans.ByDimensionTransform(transformation=translate_y, input_axes=(1,), output_axes=(1,)), + ozm06trans.ByDimensionTransform(transformation=scale_x, inputAxes=(0,), outputAxes=(0,)), + ozm06trans.ByDimensionTransform(transformation=translate_y, inputAxes=(1,), outputAxes=(1,)), ), ) edge = parse_by_dimension(model, input=xy_cs, output=xy_cs) @@ -217,8 +217,8 @@ def test_csgen_output_derived_from_max_output_axis_index(self): translate_y = ozm06trans.Translation(translation=(5.0,)) model = ozm06trans.ByDimension( transformations=( - ozm06trans.ByDimensionTransform(transformation=scale_x, input_axes=(0,), output_axes=(0,)), - ozm06trans.ByDimensionTransform(transformation=translate_y, input_axes=(1,), output_axes=(1,)), + ozm06trans.ByDimensionTransform(transformation=scale_x, inputAxes=(0,), outputAxes=(0,)), + ozm06trans.ByDimensionTransform(transformation=translate_y, inputAxes=(1,), outputAxes=(1,)), ), ) edge = parse_by_dimension(model, input=xy_cs, output=CsGen(base_name="gen")) From 113cf7d6e679a8b5c7b083e54337de13e7586526 Mon Sep 17 00:00:00 2001 From: Tomaz Vieira Date: Thu, 10 Sep 2026 12:35:59 +0000 Subject: [PATCH 25/26] Use frozen dataclasses for CoordSystems and Axis. Fix MapAxisEdge --- src/spatialdata/_io/io_raster.py | 2 +- src/spatialdata/transformations/graph/edge.py | 4 +- src/spatialdata/transformations/graph/vert.py | 57 ++++--------- tests/transformations/edge/conftest.py | 82 ++++++++++++++++--- .../edge/test_edge_transformations.py | 32 ++++---- 5 files changed, 106 insertions(+), 71 deletions(-) diff --git a/src/spatialdata/_io/io_raster.py b/src/spatialdata/_io/io_raster.py index 0b336999d..1cfc1a942 100644 --- a/src/spatialdata/_io/io_raster.py +++ b/src/spatialdata/_io/io_raster.py @@ -218,7 +218,7 @@ def try_parse_ngff06_multiscale(multiscale: oz.OMEZarrMultiscale) -> tuple[DataT # in levels of a xr.DataTree pixel_cs = CoordSystem( name=transf.input.path, - axes=[Axis(name=ax.name, type=ax.type) for ax in intrinsic_cs.axes], + axes=tuple(Axis(name=ax.name, type=ax.type) for ax in intrinsic_cs.axes), virtual=True, ) diff --git a/src/spatialdata/transformations/graph/edge.py b/src/spatialdata/transformations/graph/edge.py index 2201b257b..7bea0b82e 100644 --- a/src/spatialdata/transformations/graph/edge.py +++ b/src/spatialdata/transformations/graph/edge.py @@ -681,10 +681,8 @@ def __init__( if linear_matrix's shape isn't (output.num_axes, input.num_axes) IncompatibleCoordSystemsError if input and output don't have the same number of axes - DeterminantDifferentFromOne + NotUnimodularError if linear_matrix doesn't have determinant ~= 1 - NotOrthonormal - if linear_matrix is not orthonormal """ if input.num_axes != output.num_axes: raise IncompatibleCoordSystemsError( diff --git a/src/spatialdata/transformations/graph/vert.py b/src/spatialdata/transformations/graph/vert.py index a78046b90..62b118b16 100644 --- a/src/spatialdata/transformations/graph/vert.py +++ b/src/spatialdata/transformations/graph/vert.py @@ -1,9 +1,9 @@ from __future__ import annotations from collections.abc import Sequence +from dataclasses import dataclass from typing import Final, Literal -import ome_zarr.classes.image as ozi import ome_zarr_models.v06.coordinate_transforms as ozm06ct @@ -11,34 +11,16 @@ class AxisParsingException(Exception): pass +@dataclass(frozen=True) class Axis: - """ - Representation of an axis, following the NGFF specification. - - Attributes - ---------- - name - name of the axis. - type - type of the axis. Should be in ["channel", "space"]. - unit - unit of the axis. For a set of valid options see https://ngff.openmicroscopy.org/ - long_name: - a longer, human-friendly name for this axis - """ + """Representation of a coordinate system axis""" name: Final[str] type: Final[Literal["space", "channel"]] - unit: Final[str | None] - long_name: Final[str | None] - - def __init__( - self, *, name: str, type: Literal["space", "channel"], unit: str | None = None, long_name: str | None = None - ): - self.name = name - self.type = type - self.unit = unit - self.long_name = long_name + unit: Final[str | None] = None + "unit of the axis. For a set of valid options see https://ngff.openmicroscopy.org/" + long_name: Final[str | None] = None + "a longer, human-friendly name for this axis" def cloned_with(self, *, unit: str | None) -> Axis: return Axis(name=self.name, type=self.type, unit=unit or self.unit, long_name=self.long_name) @@ -47,7 +29,7 @@ def __hash__(self) -> int: return hash((self.name, self.type, self.unit, self.long_name)) def __repr__(self) -> str: - return f"NgffAxis(name={self.name}, type={self.type})" + return f"Axis(name={self.name}, type={self.type})" def __eq__(self, value: object, /) -> bool: if not isinstance(value, Axis): @@ -98,6 +80,7 @@ def __init__(self, *, axis_name: str) -> None: super().__init__(f"Axis name '{axis_name}' is used more than once") +@dataclass(frozen=True) class CoordSystem: """ Representation of a coordinate system. @@ -105,30 +88,17 @@ class CoordSystem: name: Final[str] axes: Final[tuple[Axis, ...]] - virtual: Final[bool] """A virtual coordinate system exists as an intermediate step between non-virtual coordinate systems and is usually ignored during serialization""" - def __init__(self, name: str, axes: Sequence[Axis], virtual: bool = False): + def __post_init__(self) -> None: """ - Parameters - ---------- - name - name of the coordinate system - axes - axes of the coordinate system - virtual - vitual coordinate systems don't serialize to NGFF - Raises ------ DuplicateAxisNameError if `axes` contains axes with duplicate names """ - self.name = name - self.axes = tuple(axes) - self.virtual = virtual seen_names: set[str] = set() for axis in self.axes: if axis.name in seen_names: @@ -142,7 +112,7 @@ def __hash__(self) -> int: return hash((self.name, self.axes, self.virtual)) @classmethod - def try_from_model(cls, model: ozi.CoordinateSystem) -> CoordSystem: + def try_from_model(cls, model: ozm06ct.CoordinateSystem) -> CoordSystem: """ Parse a `CoordSystem` from an ome-zarr-models coordinate system model. @@ -155,11 +125,12 @@ def try_from_model(cls, model: ozi.CoordinateSystem) -> CoordSystem: """ return CoordSystem( name=model.name, - axes=[Axis.try_from_model(axis) for axis in model.axes], + axes=tuple(Axis.try_from_model(axis) for axis in model.axes), + virtual=False, ) @classmethod - def try_from_model_or_default[T](cls, model: ozi.CoordinateSystem | None, *, default: T) -> CoordSystem | T: + def try_from_model_or_default[T](cls, model: ozm06ct.CoordinateSystem | None, *, default: T) -> CoordSystem | T: """ Parse a `CoordSystem` from `model`, or return `default` if `model` is None. diff --git a/tests/transformations/edge/conftest.py b/tests/transformations/edge/conftest.py index 75de85b49..5ff56a219 100644 --- a/tests/transformations/edge/conftest.py +++ b/tests/transformations/edge/conftest.py @@ -7,13 +7,75 @@ z_axis = Axis(name="z", type="space", unit="micrometer") c_axis = Axis(name="c", type="channel") -x_cs = CoordSystem(name="x", axes=[x_axis]) -y_cs = CoordSystem(name="y", axes=[y_axis]) -z_cs = CoordSystem(name="z", axes=[z_axis]) -c_cs = CoordSystem(name="c", axes=[c_axis]) -xy_cs = CoordSystem(name="xy", axes=[x_axis, y_axis]) -yx_cs = CoordSystem(name="yx", axes=[y_axis, x_axis]) -xyz_cs = CoordSystem(name="xyz", axes=[x_axis, y_axis, z_axis]) -zyx_cs = CoordSystem(name="zyx", axes=[z_axis, y_axis, x_axis]) -xyc_cs = CoordSystem(name="xyc", axes=[x_axis, y_axis, c_axis]) -cyx_cs = CoordSystem(name="cyx", axes=[c_axis, y_axis, x_axis]) +a_axis = Axis(name="a", type="space", unit="micrometer") +b_axis = Axis(name="b", type="space", unit="micrometer") +c_axis = Axis(name="c", type="space", unit="micrometer") + +x_cs = CoordSystem(name="x", axes=(x_axis,), virtual=False) +y_cs = CoordSystem(name="y", axes=(y_axis,), virtual=False) +z_cs = CoordSystem(name="z", axes=(z_axis,), virtual=False) +c_cs = CoordSystem(name="c", axes=(c_axis,), virtual=False) +xy_cs = CoordSystem( + name="xy", + axes=( + x_axis, + y_axis, + ), + virtual=False, +) +yx_cs = CoordSystem( + name="yx", + axes=( + y_axis, + x_axis, + ), + virtual=False, +) +xyz_cs = CoordSystem( + name="xyz", + axes=( + x_axis, + y_axis, + z_axis, + ), + virtual=False, +) +zyx_cs = CoordSystem( + name="zyx", + axes=( + z_axis, + y_axis, + x_axis, + ), + virtual=False, +) +xyc_cs = CoordSystem( + name="xyc", + axes=( + x_axis, + y_axis, + c_axis, + ), + virtual=False, +) +cyx_cs = CoordSystem( + name="cyx", + axes=( + c_axis, + y_axis, + x_axis, + ), + virtual=False, +) + +abc_cs = CoordSystem(name="abc", axes=(a_axis, b_axis, c_axis), virtual=False) +cba_cs = CoordSystem( + name="cba", + axes=( + c_axis, + b_axis, + a_axis, + ), + virtual=False, +) +cba_cs = CoordSystem(name="cba", axes=(c_axis, b_axis, a_axis), virtual=False) diff --git a/tests/transformations/edge/test_edge_transformations.py b/tests/transformations/edge/test_edge_transformations.py index 8283369d3..22f6f15d1 100644 --- a/tests/transformations/edge/test_edge_transformations.py +++ b/tests/transformations/edge/test_edge_transformations.py @@ -26,6 +26,11 @@ TranslationEdge, ) from tests.transformations.edge.conftest import ( + a_axis, + abc_cs, + b_axis, + c_axis, + x_axis, x_cs, xy_cs, xyc_cs, @@ -71,11 +76,6 @@ def test_from_affine_matrix(self): np.testing.assert_equal(edge.translation, translation) np.testing.assert_allclose(edge.transform_points(POINTS_2D), np.array([[3.0, 11.0], [7.0, 17.0]])) - def test_mapping_classmethod_builds_permutation_matrix(self): - edge = AffineEdge.mapping(input=xy_cs, output=yx_cs) - np.testing.assert_allclose(edge.linear, np.array([[0.0, 1.0], [1.0, 0.0]])) - np.testing.assert_allclose(edge.transform_points(POINTS_2D), np.array([[2.0, 1.0], [4.0, 3.0]])) - def test_transform_points_scale_and_translate(self): edge = AffineEdge( linear=np.array([[2.0, 0.0], [0.0, 3.0]]), translation=np.array([1.0, 1.0]), input=xy_cs, output=xy_cs @@ -131,25 +131,29 @@ def test_inverse_is_still_the_identity(self): class TestMapAxisEdge: - def test_constructor_rejects_different_sets_of_axes(self): - with pytest.raises(IncompatibleCoordSystemsError): - MapAxisEdge(input=xy_cs, output=xyz_cs) - def test_transform_points_swaps_axes(self): - edge = MapAxisEdge(input=xy_cs, output=yx_cs) - np.testing.assert_allclose(edge.transform_points(POINTS_2D), np.array([[2.0, 1.0], [4.0, 3.0]])) + edge = MapAxisEdge( + input=xyz_cs, + output=abc_cs, + input_to_output={x_axis: b_axis, y_axis: c_axis, z_axis: a_axis}, + ) + np.testing.assert_equal(edge.transform_points(POINTS_3D), np.array([[3.0, 1.0, 2.0], [6.0, 4.0, 5.0]])) def test_transform_points_permutation_of_three_axes(self): - edge = MapAxisEdge(input=xyz_cs, output=zyx_cs) + edge = MapAxisEdge( + input=xyz_cs, output=abc_cs, input_to_output={x_axis: c_axis, y_axis: b_axis, z_axis: a_axis} + ) np.testing.assert_allclose(edge.transform_points(POINTS_3D), np.array([[3.0, 2.0, 1.0], [6.0, 5.0, 4.0]])) def test_to_affine_matches_transform_points(self): - edge = MapAxisEdge(input=xy_cs, output=yx_cs) + edge = MapAxisEdge(input=xy_cs, output=xy_cs, input_to_output={x_axis: y_axis, y_axis: x_axis}) affine = edge.to_affine() np.testing.assert_allclose(affine.transform_points(POINTS_2D), edge.transform_points(POINTS_2D)) def test_inverse_roundtrips(self): - edge = MapAxisEdge(input=xyz_cs, output=zyx_cs) + edge = MapAxisEdge( + input=xyz_cs, output=zyx_cs, input_to_output={x_axis: z_axis, y_axis: y_axis, z_axis: x_axis} + ) inv = edge.inverse() np.testing.assert_allclose(inv.transform_points(edge.transform_points(POINTS_3D)), POINTS_3D) From 876c176b0870a0775368d7fd37ca7f68a9b9fa8f Mon Sep 17 00:00:00 2001 From: Tomaz Vieira Date: Thu, 10 Sep 2026 13:49:16 +0000 Subject: [PATCH 26/26] Enforces Rotation matrix is orthonormal with det==1 --- .../_core/transformation_manager/exceptions.py | 8 +++++++- src/spatialdata/transformations/graph/edge.py | 11 ++++++++--- .../transformations/edge/test_edge_transformations.py | 11 ++++++++--- 3 files changed, 23 insertions(+), 7 deletions(-) diff --git a/src/spatialdata/_core/transformation_manager/exceptions.py b/src/spatialdata/_core/transformation_manager/exceptions.py index a9dec744b..06fa5e96b 100644 --- a/src/spatialdata/_core/transformation_manager/exceptions.py +++ b/src/spatialdata/_core/transformation_manager/exceptions.py @@ -264,12 +264,18 @@ def __init__( super().__init__(message) -class NotUnimodularError(Exception): +class DeterminantDifferentFromOne(Exception): def __init__(self, matrix: ArrayLike) -> None: self.matrix = matrix super().__init__("Matrix does not have det(M) == 1") +class NotOrthonormalError(Exception): + def __init__(self, matrix: ArrayLike) -> None: + self.matrix = matrix + super().__init__("Matrix is not orthonormal") + + class EmptyTransformSequenceError(Exception): def __init__(self) -> None: super().__init__("Empty sequence of transformations") diff --git a/src/spatialdata/transformations/graph/edge.py b/src/spatialdata/transformations/graph/edge.py index 7bea0b82e..9f20151e6 100644 --- a/src/spatialdata/transformations/graph/edge.py +++ b/src/spatialdata/transformations/graph/edge.py @@ -11,10 +11,11 @@ from spatialdata._core.transformation_manager.exceptions import ( AxisRedefinitionError, + DeterminantDifferentFromOne, EmptyTransformSequenceError, IncompatibleCoordSystemsError, MissingAxisError, - NotUnimodularError, + NotOrthonormalError, UnexpectedShapeError, UnmappedAxisError, ) @@ -681,8 +682,10 @@ def __init__( if linear_matrix's shape isn't (output.num_axes, input.num_axes) IncompatibleCoordSystemsError if input and output don't have the same number of axes - NotUnimodularError + DeterminantDifferentFromOne if linear_matrix doesn't have determinant ~= 1 + NotOrthonormal + if linear_matrix is not orthonormal """ if input.num_axes != output.num_axes: raise IncompatibleCoordSystemsError( @@ -693,8 +696,10 @@ def __init__( raise UnexpectedShapeError( array_name="linear_matrix", array_shape=linear_matrix.shape, expected_shape=expected_shape ) + if not np.allclose(linear_matrix.T @ linear_matrix, np.identity(input.num_axes)): + raise NotOrthonormalError(matrix=linear_matrix) if not np.isclose(np.linalg.det(linear_matrix), 1.0): - raise NotUnimodularError(matrix=linear_matrix) + raise DeterminantDifferentFromOne(matrix=linear_matrix) linear_matrix.flags.writeable = False self.rotation = linear_matrix super().__init__(input=input, output=output, name=name) diff --git a/tests/transformations/edge/test_edge_transformations.py b/tests/transformations/edge/test_edge_transformations.py index 22f6f15d1..2427da620 100644 --- a/tests/transformations/edge/test_edge_transformations.py +++ b/tests/transformations/edge/test_edge_transformations.py @@ -7,10 +7,11 @@ from spatialdata._core.transformation_manager.exceptions import ( AxisRedefinitionError, + DeterminantDifferentFromOne, EmptyTransformSequenceError, IncompatibleCoordSystemsError, MissingAxisError, - NotUnimodularError, + NotOrthonormalError, UnexpectedShapeError, UnmappedAxisError, ) @@ -263,10 +264,14 @@ def test_constructor_rejects_wrong_shape(self): with pytest.raises(UnexpectedShapeError): RotationEdge(linear_matrix=np.eye(3), input=xy_cs, output=xy_cs) - def test_constructor_rejects_non_unimodular_matrix(self): - with pytest.raises(NotUnimodularError): + def test_constructor_rejects_matrix_with_det_diff_from_1(self): + with pytest.raises(DeterminantDifferentFromOne): RotationEdge(linear_matrix=np.array([[1.0, 0.0], [0.0, -1.0]]), input=xy_cs, output=xy_cs) + def test_constructor_rejects_non_orthonormal_matrix(self): + with pytest.raises(NotOrthonormalError): + RotationEdge(linear_matrix=np.array([[1.0, 0.0], [0.0, 2.0]]), input=xy_cs, output=xy_cs) + def test_transform_points_rotates_90_degrees(self): edge = RotationEdge(linear_matrix=np.array([[0.0, -1.0], [1.0, 0.0]]), input=xy_cs, output=xy_cs) np.testing.assert_allclose(edge.transform_points(POINTS_2D), np.array([[-2.0, 1.0], [-4.0, 3.0]]))