From 15d8f003bd99e86a15bf641e4374575547f0c26e Mon Sep 17 00:00:00 2001 From: Ajayrama Kumaraswamy Date: Mon, 10 Aug 2026 10:04:18 +0200 Subject: [PATCH 01/31] 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 | 15 + src/spatialdata/__init__.py | 1 - src/spatialdata/_core/spatialdata.py | 2 + .../_core/transformation_manager/__init__.py | 635 +++++++++++++++++ .../transformation_manager/exceptions.py | 229 +++++++ src/spatialdata/_io/io_raster.py | 40 +- src/spatialdata/_io/io_zarr.py | 35 +- src/spatialdata/_types.py | 35 +- tests/core/transformation_manager/conftest.py | 92 +++ .../test_transformation_manager.py | 643 ++++++++++++++++++ 10 files changed, 1693 insertions(+), 34 deletions(-) create mode 100644 src/spatialdata/_core/transformation_manager/__init__.py create mode 100644 src/spatialdata/_core/transformation_manager/exceptions.py create mode 100644 tests/core/transformation_manager/conftest.py create mode 100644 tests/core/transformation_manager/test_transformation_manager.py diff --git a/.gitignore b/.gitignore index d157ddaae..1c56e190e 100644 --- a/.gitignore +++ b/.gitignore @@ -54,3 +54,18 @@ _version.py # benchmarking and profiling .asv/ profile.speedscope.json + +# other +node_modules/ + +.mypy_cache +.ruff_cache + +# Kilo and plans +.kilo/ +plans/ + +# test coverage +.coverage +htmlcov/ + diff --git a/src/spatialdata/__init__.py b/src/spatialdata/__init__.py index 96e7032dd..1c04924a4 100644 --- a/src/spatialdata/__init__.py +++ b/src/spatialdata/__init__.py @@ -59,7 +59,6 @@ "polygon_query": "spatialdata._core.query.spatial_query", # _core.spatialdata "SpatialData": "spatialdata._core.spatialdata", - # _io._utils "get_dask_backing_files": "spatialdata._io._utils", # _io.format "SpatialDataFormatType": "spatialdata._io.format", diff --git a/src/spatialdata/_core/spatialdata.py b/src/spatialdata/_core/spatialdata.py index 9f5c91c7a..34e25aa64 100644 --- a/src/spatialdata/_core/spatialdata.py +++ b/src/spatialdata/_core/spatialdata.py @@ -22,6 +22,7 @@ from zarr.errors import GroupNotFoundError from spatialdata._core._elements import Images, Labels, Points, Shapes, Tables +from spatialdata._core.transformation_manager import TransformationManager from spatialdata._core.validation import ( check_all_keys_case_insensitively_unique, check_target_region_column_symmetry, @@ -132,6 +133,7 @@ def __init__( self._shapes: Shapes = Shapes(shared_keys=self._shared_keys) self._tables: Tables = Tables(shared_keys=self._shared_keys) self.attrs = attrs if attrs else {} # type: ignore[assignment] + self.transformation_manager = TransformationManager() # Initialize the TransformationManager element_names = list(chain.from_iterable([e.keys() for e in [images, labels, points, shapes] if e is not None])) diff --git a/src/spatialdata/_core/transformation_manager/__init__.py b/src/spatialdata/_core/transformation_manager/__init__.py new file mode 100644 index 000000000..079fe5cce --- /dev/null +++ b/src/spatialdata/_core/transformation_manager/__init__.py @@ -0,0 +1,635 @@ +from __future__ import annotations + +from collections.abc import Sequence + +import networkx as nx + +from spatialdata._core.transformation_manager.exceptions import ( + CannotRemoveCoordinateSystemError, + CoordinateSystemAlreadyExistsError, + CoordinateSystemHasElementsError, + CoordinateSystemHasTransformationsError, + CoordinateSystemNotFoundError, + ElementAlreadyExistsError, + ElementNotRegisteredToAnyCoordinateSystemError, + InvalidPathError, + TransformationNotFoundError, + TransformationPathAmbiguousError, + TransformationPathAmbiguousMultipleEdgeExpectedError, + TransformationPathAmbiguousNoEdgeExpectedError, + TransformationPathNotFoundError, + TransformationPathNotSimple, +) +from spatialdata.transformations import transformations as sd_transforms +from spatialdata.transformations.ngff.ngff_coordinate_system import NgffCoordinateSystem + +TRANSFORM_KEY = "transformation" +type EDGE_DEF = tuple[NgffCoordinateSystem, NgffCoordinateSystem, sd_transforms.BaseTransformation] + + +class TransformationManager: + def __init__(self) -> None: + """Initialize a TransformationManager with empty graph and mappings.""" + self._graph: nx.MultiDiGraph[NgffCoordinateSystem] = nx.MultiDiGraph() + # MultiDiGraph with NgffCoordinateSystem objects as nodes and transforms as edge attributes + self._element_to_cs_mapping: dict[str, NgffCoordinateSystem] = {} + # mapping element_name to the coordinate system to which the element belongs + + def assert_element_exists(self, element_name: str) -> None: + """ + Assert that an element exists in the transformation manager. + + Parameters + ---------- + element_name + The name of the element to check. + + Raises + ------ + ElementNotFoundError + If the element does not exist. + """ + if element_name not in self._element_to_cs_mapping: + raise ElementNotRegisteredToAnyCoordinateSystemError(element_name) + + def assert_element_does_not_exist(self, element_name: str) -> None: + """ + Assert than an element doesn't exist in the transformation manager. + + Parameters + ---------- + element_name + The name of the element to check. + + Raises + ------ + ElementAlreadyExistsError + If the element already exists. + """ + if element_name in self._element_to_cs_mapping: + raise ElementAlreadyExistsError(element_name) + + def assert_coordinate_system_exists(self, cs: NgffCoordinateSystem) -> None: + """ + Assert that coordinate system exists in the graph. + + Parameters + ---------- + cs + The coordinate system to check. + + Raises + ------ + CoordinateSystemNotFoundError + If the coordinate system does not exist. + """ + if cs not in self._graph: + raise CoordinateSystemNotFoundError(cs.name) + + def assert_coordinate_system_does_not_exist(self, cs: NgffCoordinateSystem) -> None: + """ + Assert that coordinate system doesn't exist in the graph. + + Parameters + ---------- + cs + The coordinate system to check. + + Raises + ------ + CoordinateSystemAlreadyExistsError + If the coordinate system already exists. + """ + if cs in self._graph: + raise CoordinateSystemAlreadyExistsError(cs.name) + + def assert_an_edge_exists_between_coordinate_systems( + self, source_cs: NgffCoordinateSystem, target_cs: NgffCoordinateSystem, edge_key: str | None = None + ) -> None: + """ + Assert that an edge exists between coordinate systems. + + Parameters + ---------- + source_cs + The input coordinate system. + target_cs + The output coordinate system. + edge_key + If specified, asserts that the specific edge exists between coordinate systems + + + Raises + ------ + CoordinateSystemNotFoundError + If either coordinate system does not exist. + TransformationNotFoundError + If the edge does not exist. + """ + self.assert_coordinate_system_exists(source_cs) + self.assert_coordinate_system_exists(target_cs) + if not self._graph.has_edge(source_cs, target_cs): + raise TransformationNotFoundError(source_cs.name, target_cs.name, edge_key) + + def assert_coordinate_system_has_no_transformations(self, cs: NgffCoordinateSystem) -> None: + """ + Assert that coordinate system has no associated transformations. + + Parameters + ---------- + cs + The coordinate system to check. + + Raises + ------ + CoordinateSystemNotFoundError + If either coordinate system does not exist. + CoordinateSystemHasTransformationsError + If the coordinate system has associated transformations. + + """ + self.assert_coordinate_system_exists(cs) + has_successors = next(self._graph.successors(cs), None) is not None + has_predecessors = next(self._graph.predecessors(cs), None) is not None + + if has_successors or has_predecessors: + raise CoordinateSystemHasTransformationsError(cs.name) + + def assert_coordinate_system_has_no_elements(self, cs: NgffCoordinateSystem) -> None: + """ + Assert that a coordinate system doesn't have elements belonging to it. + + Parameters + ---------- + cs + The coordinate system to check + + Raises + ------ + CoordinateSystemHasElementsError + if the coordinate system has elements belonging to it. + """ + elements = self._get_elements_belonging_to_cs(cs) + # also checks if cs exists + if elements: + raise CoordinateSystemHasElementsError(cs.name, elements) + + def add_coordinate_system(self, cs: NgffCoordinateSystem) -> None: + """ + Add a coordinate system to the transformation manager. + + Parameters + ---------- + cs + The coordinate system to add. + + Raises + ------ + CoordinateSystemAlreadyExistsError + If the coordinate system already exists. + """ + self.assert_coordinate_system_does_not_exist(cs) + self._graph.add_node(cs) + + def remove_coordinate_system(self, cs: NgffCoordinateSystem) -> None: + """ + Remove a coordinate system from the transformation manager. + + Parameters + ---------- + cs + The coordinate system to remove. + + Raises + ------ + CoordinateSystemNotFoundError + If the coordinate system is not found + CannotRemoveCoordinateSystemError + If the coordinate system cannot be removed. + Caused by + CoordinateSystemHasTransformationsError or + CoordinateSystemHasElementsError when the system has dependencies. + """ + try: + self.assert_coordinate_system_has_no_transformations(cs) + # also asserts that cs exists + self.assert_coordinate_system_has_no_elements(cs) + except (CoordinateSystemHasTransformationsError, CoordinateSystemHasElementsError) as err: + raise CannotRemoveCoordinateSystemError(cs.name) from err + + self._graph.remove_node(cs) + + def list_coordinate_systems(self) -> list[NgffCoordinateSystem]: + """ + List all registered coordinate systems. + + Returns + ------- + A list of coordinate system objects. + """ + return list(self._graph.nodes()) + + def add_element(self, element_name: str, coordinate_system: NgffCoordinateSystem) -> None: + """ + Register an element and associate it with a coordinate system. + + Parameters + ---------- + element_name + The name of the element. + coordinate_system + The coordinate system to which the element belongs. + + Raises + ------ + CoordinateSystemNotFoundError + If the coordinate system is not found. + """ + self.assert_element_does_not_exist(element_name) + self.assert_coordinate_system_exists(coordinate_system) + self._element_to_cs_mapping[element_name] = coordinate_system + + def get_element_coordinate_system(self, element_name: str) -> NgffCoordinateSystem: + """ + Get the coordinate system to which an element belongs. + + Parameters + ---------- + element_name + The name of the element. + + Returns + ------- + The coordinate system + + Raises + ------ + ElementNotFoundError + If the element does not exist. + """ + self.assert_element_exists(element_name) + return self._element_to_cs_mapping[element_name] + + def _unset_element(self, element_name: str) -> None: + """ + Unregister an element from the coordinate system to which it belongs. + + Parameters + ---------- + element_name + The name of the element. + + Raises + ------ + ElementNotFoundError + If the element has not been registered to any coordinate system. + """ + self.assert_element_exists(element_name) + del self._element_to_cs_mapping[element_name] + + @staticmethod + def _get_edge_key(edge_def: EDGE_DEF) -> str: + """ + Encode the definition of an edge into a string to be used as graph-wide edge identifier. + + Parameters + ---------- + edge_def + tuple, (source_cs, target_cs, transformation) + + Returns + ------- + edge_key + str, the encoded edge identifier + """ + source_cs, target_cs, transform = edge_def + return f"{source_cs.name}_{target_cs.name}_{id(transform)}" + + def add_transformation( + self, + source_cs: NgffCoordinateSystem, + target_cs: NgffCoordinateSystem, + transformation: sd_transforms.BaseTransformation, + ) -> None: + """ + Add a transformation between coordinate systems. + + Parameters + ---------- + source_cs + The input coordinate system. + target_cs + The output coordinate system. + transformation + The transformation to add. + + Raises + ------ + CoordinateSystemNotFoundError + If either coordinate system does not exist. + """ + self.assert_coordinate_system_exists(source_cs) + self.assert_coordinate_system_exists(target_cs) + + edge_key = self._get_edge_key((source_cs, target_cs, transformation)) + edge_attributes = {TRANSFORM_KEY: transformation} + self._graph.add_edge(source_cs, target_cs, key=edge_key, **edge_attributes) + + def get_direct_transformations( + self, source_cs: NgffCoordinateSystem, target_cs: NgffCoordinateSystem + ) -> list[sd_transforms.BaseTransformation]: + """ + Retrieve transformations directly defined between coordinate systems. + + Parameters + ---------- + source_cs + The input coordinate system. + target_cs + The output coordinate system. + + Returns + ------- + List of transformations (empty if there are no direct transformations). + + Raises + ------ + CoordinateSystemNotFoundError + If either coordinate system does not exist. + """ + self.assert_coordinate_system_exists(source_cs) + self.assert_coordinate_system_exists(target_cs) + + transforms = [] + if self._graph.has_edge(source_cs, target_cs): + for edge_data in self._graph[source_cs][target_cs].values(): + transform: sd_transforms.BaseTransformation = edge_data[TRANSFORM_KEY] + transforms.append(transform) + return transforms + + def remove_specific_transformation( + self, + source_cs: NgffCoordinateSystem, + target_cs: NgffCoordinateSystem, + transformation: sd_transforms.BaseTransformation, + ) -> None: + """ + Remove a specific transformation between coordinate systems. + + Parameters + ---------- + source_cs + The input coordinate system. + target_cs + The output coordinate system. + transformation + The transformation to remove. + (mainly useful for cases with multiple transformations between the same coordinate systems). + + Raises + ------ + CoordinateSystemNotFoundError + If either coordinate system does not exist. + TransformationNotFoundError + If the transformation does not exist. + """ + expected_edge_key = self._get_edge_key((source_cs, target_cs, transformation)) + self.assert_an_edge_exists_between_coordinate_systems(source_cs, target_cs, expected_edge_key) + # also checks if source_cs and target_cs exist + self._graph.remove_edge(source_cs, target_cs, key=expected_edge_key) + + def remove_all_transformations_between_coordinate_systems( + self, + source_cs: NgffCoordinateSystem, + target_cs: NgffCoordinateSystem, + ) -> None: + """ + Remove all transformation between coordinate systems. + + Parameters + ---------- + source_cs + The input coordinate system. + target_cs + The output coordinate system. + + Raises + ------ + CoordinateSystemNotFoundError + If either coordinate system does not exist. + TransformationNotFoundError + If no transformation exists between the coordiante systems + """ + self.assert_coordinate_system_exists(source_cs) + self.assert_coordinate_system_exists(target_cs) + if self._graph.has_edge(source_cs, target_cs): + for edge_key in list(self._graph[source_cs][target_cs].keys()): + # need to covert keys() to list to freeze it, else it will change during the following removal + self._graph.remove_edge(source_cs, target_cs, key=edge_key) + + def _get_transformation_sequences_from_simple_paths_after_disambiguation( + self, + paths: list[list[NgffCoordinateSystem]], + expected_intermediate_edges: Sequence[EDGE_DEF], + ) -> list[sd_transforms.Sequence]: + """ + Traverses paths to form sequence of Transformations. + + In case of ambiguity looks into `expected_intermediate_transformations` to disambiguate. + + Parameters + ---------- + paths: + sequence of list of nodes + expected_intermediate_edges: + list of edge definitions, for use when multiple edges are found between two coordinate systems in a path + An edge is defined as a tuple, (source_cs, target_cs, transform). + + Returns + ------- + list of transformations, each an instance of the transformation class Sequence + + Raises + ------ + TransformationPathNotSimple + if any path in `paths` is not simple, i.e., has recurring coordinate systems + """ + for path in paths: + if len(set(path)) != len(path): + raise TransformationPathNotSimple(path) + + expected_intermediate_transformation_edge_keys = set() + if expected_intermediate_edges is not None: + for edge_def in expected_intermediate_edges: + edge_key = self._get_edge_key(edge_def) + expected_intermediate_transformation_edge_keys |= {edge_key} + + all_sequences = [] + deduplicated_paths = list({repr(x): x for x in paths}.values()) + for path in deduplicated_paths: + if len(path) <= 1: + raise InvalidPathError(path) + + transformation_list = [] + for i in range(len(path) - 1): + source_cs_here, target_cs_here = path[i], path[i + 1] + edge_data = self._graph[source_cs_here][target_cs_here] + if len(edge_data) > 1: + # when there are multiple edges between a pair of coordinate systems in the path + # find the expected edge based on key + expected_intermediate_transformation_key_here = ( + expected_intermediate_transformation_edge_keys & set(edge_data.keys()) + ) + if len(expected_intermediate_transformation_key_here) == 0: + # transformation was not specified in `expected_intermediate_transformations` for disambiguation + raise TransformationPathAmbiguousNoEdgeExpectedError(source_cs_here.name, target_cs_here.name) + if len(expected_intermediate_transformation_key_here) > 1: + # multiple transformations were specified in `expected_intermediate_transformations` + # for disambiguation + raise TransformationPathAmbiguousMultipleEdgeExpectedError( + source_cs_here.name, target_cs_here.name, len(expected_intermediate_transformation_key_here) + ) + edge_key_to_use = list(expected_intermediate_transformation_key_here)[0] + transformation_list.append(edge_data[edge_key_to_use][TRANSFORM_KEY]) + else: + # Only one edge, no ambiguity + edge_key = next(iter(edge_data.keys())) + transformation_list.append(edge_data[edge_key][TRANSFORM_KEY]) + all_sequences.append(sd_transforms.Sequence(transformation_list)) + return all_sequences + + def get_all_shortest_transformation_sequences( + self, + source_cs: NgffCoordinateSystem, + target_cs: NgffCoordinateSystem, + expected_intermediate_edges: Sequence[EDGE_DEF] = (), + ) -> list[sd_transforms.Sequence]: + """ + Get all shortest sequences of transformations between two coordinate systems. + + Parameters + ---------- + source_cs + The source coordinate system. + target_cs + The target coordinate system. + expected_intermediate_edges + list of intermediate edges expected. + An edge is defined as a tuple, (source_cs, target_cs, transform). + Used to choose an edge when multiple edges are found between the same coordinate systems. + + Returns + ------- + list[Sequence] + All shortest sequences of transformations from source_cs to target_cs. + When multiple transformations are defined between the same coordinate systems, only those containing + a transformation among intermediate transformations are included. + + Raises + ------ + CoordinateSystemNotFoundError + If either coordinate system does not exist. + TransformationPathNotFoundError + If no path exists between the source and target coordinate systems. + TransformationPathAmbiguousError + When multiple transformations are defined between the same coordinate systems and transformations + are not specified in `expeceted_intermediate_transformations` for disambiguation. + """ + try: + paths = list(nx.all_shortest_paths(self._graph, source=source_cs, target=target_cs)) + except nx.NetworkXNoPath as nxe: + raise TransformationPathNotFoundError(source_cs.name, target_cs.name) from nxe + + try: + return self._get_transformation_sequences_from_simple_paths_after_disambiguation( + paths, expected_intermediate_edges + ) + except TransformationPathAmbiguousError as tpae: + raise TransformationPathAmbiguousError(source_cs.name, target_cs.name) from tpae + + def get_all_transformation_sequences( + self, + source_cs: NgffCoordinateSystem, + target_cs: NgffCoordinateSystem, + expected_intermediate_edges: Sequence[EDGE_DEF] = (), + ) -> list[sd_transforms.Sequence]: + """ + Get all existing sequences of transformations between two coordinate systems. + + Parameters + ---------- + source_cs + The source coordinate system. + target_cs + The target coordinate system. + expected_intermediate_edges + list of intermediate edges expected. + An edge is defined as a tuple, (source_cs, target_cs, transform). + Used to choose an edge when multiple edges are found between the same coordinate systems. + + Returns + ------- + list[Sequence] + All existing sequences of transformations from source_cs to target_cs. + When multiple transformations are defined between the same coordinate systems, only those containing + a transformation among intermediate transformations are included. + + Raises + ------ + CoordinateSystemNotFoundError + If either coordinate system does not exist. + TransformationPathNotFoundError + If no path exists between the source and target coordinate systems. + TransformationPathAmbiguousError + When multiple transformations are defined between the same coordinate systems and transformations + are not specified in `expected_intermediate_transformations` for disambiguation. + """ + try: + paths = list(nx.all_simple_paths(self._graph, source=source_cs, target=target_cs)) + except nx.NetworkXNoPath as nxe: + raise TransformationPathNotFoundError(source_cs.name, target_cs.name) from nxe + + try: + return self._get_transformation_sequences_from_simple_paths_after_disambiguation( + paths, expected_intermediate_edges + ) + except TransformationPathAmbiguousError as tpae: + raise tpae from TransformationPathAmbiguousError(source_cs.name, target_cs.name) + + def _get_elements_belonging_to_cs(self, cs: NgffCoordinateSystem) -> list[str]: + """ + Get all elements belonging to a coordinate system. + + Parameters + ---------- + cs + The coordinate system to check. + + Returns + ------- + List of element names belonging to the coordinate system. + + Raises + ------ + CoordinateSystemNotFoundError + If the coordinate system does not exist. + + """ + self.assert_coordinate_system_exists(cs) + + elements = [] + for element_name, element_cs in self._element_to_cs_mapping.items(): + if element_cs == cs: + elements.append(element_name) + + return elements + + def __repr__(self) -> str: + """Return a string representation of the TransformationManager.""" + return ( + f"TransformationManager(" + f" coordinate_systems={list(self._graph.nodes())}, " + f" coordinate_transforms={[x[TRANSFORM_KEY] for *_, x in self._graph.edges(data=True)]}, " + f" elements={list(self._element_to_cs_mapping.keys())}" + f")" + ) diff --git a/src/spatialdata/_core/transformation_manager/exceptions.py b/src/spatialdata/_core/transformation_manager/exceptions.py new file mode 100644 index 000000000..c8324032a --- /dev/null +++ b/src/spatialdata/_core/transformation_manager/exceptions.py @@ -0,0 +1,229 @@ +from __future__ import annotations + +from spatialdata.transformations.ngff.ngff_coordinate_system import NgffCoordinateSystem + + +class CoordinateSystemNotFoundError(ValueError): + """ + Exception raised when a coordinate system is not found in the transformation manager. + + Attributes + ---------- + name : str + The name of the coordinate system that was not found. + """ + + def __init__(self, name: str) -> None: + self.name = name + super().__init__(f"Coordinate system '{name}' not found in the transformation manager.") + + +class CoordinateSystemAlreadyExistsError(ValueError): + """ + Exception raised when a coordinate system already exists in the transformation manager. + + Attributes + ---------- + name : str + The name of the coordinate system to check. + """ + + def __init__(self, name: str) -> None: + self.name = name + super().__init__(f"Coordinate system '{name}' already exists in the transformation manager.") + + +class ElementNotRegisteredToAnyCoordinateSystemError(KeyError): + """ + Exception raised when an element is not found in the transformation manager. + + Attributes + ---------- + element_name : str + The name of the element that was not found. + """ + + def __init__(self, element_name: str) -> None: + self.element_name = element_name + super().__init__(f"Element '{element_name}' has not been registered to any coordinate system.") + + +class TransformationNotFoundError(KeyError): + """ + Exception raised when a transformation is not found between coordinate systems. + + Attributes + ---------- + input_cs_name : str + The name of the input coordinate system. + output_cs_name : str + The name of the output coordinate system. + edge_key: str or None + key used when adding transformation + """ + + def __init__(self, source_cs_name: str, target_cs_name: str, edge_key: str | None = None) -> None: + self.input_cs_name = source_cs_name + self.output_cs_name = target_cs_name + self.edge_key = edge_key + msg = f"Transformation from '{source_cs_name}' to '{target_cs_name}' not found" + if edge_key is not None: + msg += f" with key '{edge_key}'" + super().__init__(msg) + + +class ElementAlreadyExistsError(ValueError): + """ + Exception raised when trying to add an Element that already exists. + + Attributes + ---------- + name : str + The name of the Element that already exists. + """ + + def __init__(self, name: str) -> None: + self.name = name + super().__init__(f"Element '{name}' already exists in the transformation manager") + + +class InvalidPathError(ValueError): + """ + Exception raised when a path is defined with less than 2 nodes. + + Attributes + ---------- + invalid_path : list[NgffCoordinateSystem] + """ + + def __init__(self, invalid_path: list[NgffCoordinateSystem]) -> None: + super().__init__(f"Found an invalid path with less than 2 nodes: {invalid_path}") + + +class TransformationPathNotFoundError(ValueError): + """ + Exception raised when no transformation path exists between coordinate systems. + + Attributes + ---------- + source_cs_name : str + The name of the source coordinate system. + target_cs_name : str + The name of the target coordinate system. + """ + + def __init__(self, source_cs_name: str, target_cs_name: str) -> None: + self.source_cs_name = source_cs_name + self.target_cs_name = target_cs_name + super().__init__(f"No transformation path found from {source_cs_name} to {target_cs_name}") + + +class TransformationPathAmbiguousError(ValueError): + """ + Exception raised when multiple transformation path exists between coordinate systems. + + Attributes + ---------- + source_cs_name : str + The name of the source coordinate system. + target_cs_name : str + The name of the target coordinate system. + """ + + def __init__(self, source_cs_name: str, target_cs_name: str) -> None: + self.source_cs_name = source_cs_name + self.target_cs_name = target_cs_name + base_msg = f"Transformation Path ambiguous from {source_cs_name} to {target_cs_name}." + cause_of_confusion = self.cause_of_confusion() + msg = f"{base_msg} {cause_of_confusion}" if cause_of_confusion else base_msg + super().__init__(msg) + + def cause_of_confusion(self) -> str: + return "" + + +class TransformationPathAmbiguousNoEdgeExpectedError(TransformationPathAmbiguousError): + def __init__(self, source_cs_name: str, target_cs_name: str) -> None: + super().__init__(source_cs_name, target_cs_name) + + def cause_of_confusion(self) -> str: + + return "Multiple edges found. None of them were specified to be expected" + + +class TransformationPathAmbiguousMultipleEdgeExpectedError(TransformationPathAmbiguousError): + def __init__(self, source_cs_name: str, target_cs_name: str, number_of_edges_expected: int) -> None: + self.number_of_edges_expected = number_of_edges_expected + super().__init__(source_cs_name, target_cs_name) + + def cause_of_confusion(self) -> str: + return f"Multiple ({self.number_of_edges_expected}) edges were specified to be expected" + + +class TransformationPathNotSimple(ValueError): + """ + Exception raised when a path represented by a list of coordinate systems is not simple. + + A simple path is one in which each coordinate system appears only once + """ + + def __init__(self, path: list[NgffCoordinateSystem]) -> None: + self.path = path + css_formatted_one_per_line = "\n".join(repr(cs) for cs in path) + super().__init__( + f"Transformation Path not simple, i.e., some coordinate systems appear multiple times:\n" + f"{css_formatted_one_per_line}" + ) + + +class CannotRemoveCoordinateSystemError(ValueError): + """Exception raised when trying to remove a coordinate system. + + Attributes + ---------- + name: str + name of the coordinate system + """ + + def __init__(self, name: str) -> None: + self.name = name + super().__init__(f"Cannot remove coordinate system with name {name}.") + + +class CoordinateSystemHasTransformationsError(ValueError): + """ + Exception raised when trying to remove a coordinate system that has associated transformations. + + Attributes + ---------- + name : str + The name of the coordinate system that has transformations. + """ + + def __init__(self, name: str) -> None: + self.name = name + super().__init__(f"Coordinate System ('{name}') has transformations.") + + +class CoordinateSystemHasElementsError(ValueError): + """ + Exception raised when trying to remove a coordinate system that has associated elements. + + Attributes + ---------- + name : str + The name of the coordinate system that has elements. + associated_elements : list[str] + List of element names associated with the coordinate system. + """ + + def __init__(self, name: str, associated_elements: list[str]) -> None: + self.name = name + self.associated_elements = associated_elements + super().__init__(f"Coordinate system '{name}' has elements belonging to it: {associated_elements}") + + +class TransformationManagerWarning(UserWarning): + """Base warning category for TransformationManager.""" + + pass diff --git a/src/spatialdata/_io/io_raster.py b/src/spatialdata/_io/io_raster.py index 3806fabf7..d2cb9c70c 100644 --- a/src/spatialdata/_io/io_raster.py +++ b/src/spatialdata/_io/io_raster.py @@ -3,7 +3,7 @@ import warnings from collections.abc import Mapping, Sequence from pathlib import Path -from typing import Any, Literal, TypeGuard +from typing import Any, Literal, TypeGuard, get_args import dask.array as da import numpy as np @@ -30,6 +30,7 @@ RasterFormatType, get_ome_zarr_format, ) +from spatialdata._types import ELEMENT_TYPE, ELEMENT_TYPE_RASTER, GROUP_NAME from spatialdata._utils import get_pyramid_levels from spatialdata.models.models import ATTRS_KEY from spatialdata.models.pyramids_utils import dask_arrays_to_datatree @@ -162,10 +163,12 @@ def _prepare_storage_options( def _read_multiscale( - store: str | Path, raster_type: Literal["image", "labels"], reader_format: Format + store: str | Path, raster_type: ELEMENT_TYPE_RASTER, reader_format: Format ) -> DataArray | DataTree: assert isinstance(store, str | Path) - assert raster_type in ["image", "labels"] + assert raster_type in get_args(ELEMENT_TYPE_RASTER), ValueError( + f"{raster_type} is not a valid raster type. Must be one of {get_args(ELEMENT_TYPE_RASTER)}." + ) nodes: list[Node] = [] image_loc = ZarrLocation(store, fmt=reader_format) @@ -214,7 +217,7 @@ def _read_multiscale( transformations = _get_transformations_from_ngff_dict(encoded_ngff_transformations) # if image, read channels metadata channels: list[Any] | None = None - if raster_type == "image": + if raster_type == ELEMENT_TYPE.IMAGE: if legacy_channels_metadata is not None: channels = [d["label"] for d in legacy_channels_metadata["channels"]] if omero_metadata is not None: @@ -229,7 +232,7 @@ def _read_multiscale( data = loaded_node.array(resolution=datasets[0]) si = DataArray( data, - name="image", + name=ELEMENT_TYPE.IMAGE, dims=axes, coords={"c": channels} if channels is not None else {}, ) @@ -265,7 +268,7 @@ def _get_multiscale_nodes(image_nodes: list[Node], nodes: list[Node]) -> list[No def _get_raster_element_group( - raster_type: Literal["image", "labels"], group: zarr.Group, element_name: str + 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. @@ -299,7 +302,7 @@ def _get_raster_element_group( def _write_raster( - raster_type: Literal["image", "labels"], + raster_type: ELEMENT_TYPE_RASTER, raster_data: DataArray | DataTree, group: zarr.Group, name: str, @@ -332,12 +335,13 @@ def _write_raster( metadata Additional metadata for the raster element """ - if raster_type not in ["image", "labels"]: - raise ValueError(f"{raster_type} is not a valid raster type. Must be 'image' or 'labels'.") + assert raster_type in get_args(ELEMENT_TYPE_RASTER), ValueError( + f"{raster_type} is not a valid raster type. Must be one of {get_args(ELEMENT_TYPE_RASTER)}." + ) # "name" and "label_metadata" are only used for labels. "name" is written in write_multiscale_ngff() but ignored in # write_image_ngff() (possibly an ome-zarr-py bug). We only use "name" to ensure correct group access in the # ome-zarr API. - if raster_type == "labels": + if raster_type == ELEMENT_TYPE.LABELS: metadata["name"] = name metadata["label_metadata"] = label_metadata @@ -367,7 +371,7 @@ def _write_raster( raise ValueError("Not a valid labels object") group = _get_raster_element_group(raster_type, group, name) - if raster_type == "image": + if raster_type == ELEMENT_TYPE.IMAGE: # ome-zarr-py >= 0.18 no longer writes the omero channel metadata, so we write it ourselves. overwrite_channel_names(group, raster_data) if ATTRS_KEY not in group.attrs: @@ -463,7 +467,7 @@ def _update_dict_v3(d: dict[str, Any]) -> None: def _write_raster_dataarray( - raster_type: Literal["image", "labels"], + raster_type: ELEMENT_TYPE_RASTER, group: zarr.Group, element_name: str, raster_data: DataArray, @@ -493,7 +497,7 @@ def _write_raster_dataarray( metadata Additional metadata for the raster element """ - write_single_scale_ngff = write_image_ngff if raster_type == "image" else write_labels_ngff + write_single_scale_ngff = write_image_ngff if raster_type == ELEMENT_TYPE.IMAGE else write_labels_ngff data = raster_data.data transformations = _get_transformations(raster_data) @@ -532,7 +536,7 @@ def _write_raster_dataarray( def _write_raster_datatree( - raster_type: Literal["image", "labels"], + raster_type: ELEMENT_TYPE_RASTER, group: zarr.Group, element_name: str, raster_data: DataTree, @@ -562,7 +566,9 @@ def _write_raster_datatree( metadata Additional metadata for the raster element """ - write_multi_scale_ngff = write_multiscale_ngff if raster_type == "image" else write_multiscale_labels_ngff + write_multi_scale_ngff = ( + write_multiscale_ngff if raster_type == ELEMENT_TYPE.IMAGE else write_multiscale_labels_ngff + ) data = get_pyramid_levels(raster_data, attr="data") list_of_input_axes: list[Any] = get_pyramid_levels(raster_data, attr="dims") assert len(set(list_of_input_axes)) == 1 @@ -632,7 +638,7 @@ def write_image( ) _write_raster( - raster_type="image", + raster_type=ELEMENT_TYPE.IMAGE, raster_data=image, group=group, name=name, @@ -659,7 +665,7 @@ def write_labels( ) _write_raster( - raster_type="labels", + raster_type=ELEMENT_TYPE.LABELS, raster_data=labels, group=group, name=name, diff --git a/src/spatialdata/_io/io_zarr.py b/src/spatialdata/_io/io_zarr.py index 186d7c974..1301761d5 100644 --- a/src/spatialdata/_io/io_zarr.py +++ b/src/spatialdata/_io/io_zarr.py @@ -5,7 +5,7 @@ from collections.abc import Callable from json import JSONDecodeError from pathlib import Path -from typing import Any, Literal, cast +from typing import Any, Literal, cast, get_args import zarr.storage from anndata import AnnData @@ -28,7 +28,14 @@ from spatialdata._io.io_shapes import _read_shapes from spatialdata._io.io_table import _read_table from spatialdata._logging import logger -from spatialdata._types import JSONValue, Raster_T +from spatialdata._types import ( + ELEMENT_TYPE, + ELEMENT_TYPE_RASTER, + ELEMENT_TYPE_VECTOR, + GROUP_NAME, + JSONValue, + Raster_T, +) def _read_zarr_group_spatialdata_element( @@ -37,8 +44,8 @@ def _read_zarr_group_spatialdata_element( sdata_version: Literal["0.1", "0.2"], selector: set[str], read_func: Callable[..., Any], - group_name: Literal["images", "labels", "shapes", "points", "tables"], - element_type: Literal["image", "labels", "shapes", "points", "tables"], + group_name: GROUP_NAME, + element_type: ELEMENT_TYPE, element_container: (dict[str, Raster_T] | dict[str, DaskDataFrame] | dict[str, GeoDataFrame] | dict[str, AnnData]), on_bad_files: Literal[BadFileHandleMethod.ERROR, BadFileHandleMethod.WARN], ) -> None: @@ -76,14 +83,14 @@ def _read_zarr_group_spatialdata_element( ValueError, ), ): - if element_type in ["image", "labels"]: + if element_type in get_args(ELEMENT_TYPE_RASTER): reader_format = get_raster_format_for_read(elem_group, sdata_version) element = read_func( elem_group_path, element_type, reader_format, ) - elif element_type in ["shapes", "points", "tables"]: + elif element_type in get_args(ELEMENT_TYPE_VECTOR) + (ELEMENT_TYPE.TABLES,): element = read_func(elem_group_path) else: raise ValueError(f"Unknown element type {element_type}") @@ -204,19 +211,19 @@ def read_zarr( # we could make this more readable. One can get lost when looking at this dict and iteration over the items group_readers: dict[ - Literal["images", "labels", "shapes", "points", "tables"], + GROUP_NAME, tuple[ Callable[..., Any], - Literal["image", "labels", "shapes", "points", "tables"], + ELEMENT_TYPE, dict[str, Raster_T] | dict[str, DaskDataFrame] | dict[str, GeoDataFrame] | dict[str, AnnData], ], ] = { - # ome-zarr-py needs a kwargs that has "image" has key. So here we have "image" and not "images" - "images": (_read_multiscale, "image", images), - "labels": (_read_multiscale, "labels", labels), - "points": (_read_points, "points", points), - "shapes": (_read_shapes, "shapes", shapes), - "tables": (_read_table, "tables", tables), + # ome-zarr-py needs a kwargs that has "image" as key. So here we have "image" and not "images" + GROUP_NAME.IMAGES: (_read_multiscale, ELEMENT_TYPE.IMAGE, images), + GROUP_NAME.LABELS: (_read_multiscale, ELEMENT_TYPE.LABELS, labels), + GROUP_NAME.POINTS: (_read_points, ELEMENT_TYPE.POINTS, points), + GROUP_NAME.SHAPES: (_read_shapes, ELEMENT_TYPE.SHAPES, shapes), + GROUP_NAME.TABLES: (_read_table, ELEMENT_TYPE.TABLES, tables), } for group_name, ( read_func, diff --git a/src/spatialdata/_types.py b/src/spatialdata/_types.py index e41de9fe7..014587c9e 100644 --- a/src/spatialdata/_types.py +++ b/src/spatialdata/_types.py @@ -1,11 +1,22 @@ from __future__ import annotations -from typing import Any +from enum import StrEnum +from typing import Any, Literal import numpy as np from xarray import DataArray, DataTree -__all__ = ["ArrayLike", "ColorLike", "DTypeLike", "JSONValue", "Raster_T"] +__all__ = [ + "ArrayLike", + "ColorLike", + "DTypeLike", + "JSONValue", + "Raster_T", + "ELEMENT_TYPE", + "ELEMENT_TYPE_RASTER", + "ELEMENT_TYPE_VECTOR", + "GROUP_NAME", +] from numpy.typing import DTypeLike, NDArray @@ -26,3 +37,23 @@ # arrays, sets, DataFrames). Note that JSON has no tuples and only string keys, so a tuple is read back as a list and # a non-string key as a string. type JSONValue = dict[str, JSONValue] | list[JSONValue] | str | int | float | bool | None + + +class ELEMENT_TYPE(StrEnum): + IMAGE = "image" + LABELS = "labels" + SHAPES = "shapes" + POINTS = "points" + TABLES = "tables" + + +class GROUP_NAME(StrEnum): + IMAGES = "images" + LABELS = "labels" + SHAPES = "shapes" + POINTS = "points" + TABLES = "tables" + + +ELEMENT_TYPE_RASTER = Literal[ELEMENT_TYPE.IMAGE, ELEMENT_TYPE.LABELS] +ELEMENT_TYPE_VECTOR = Literal[ELEMENT_TYPE.POINTS, ELEMENT_TYPE.SHAPES] diff --git a/tests/core/transformation_manager/conftest.py b/tests/core/transformation_manager/conftest.py new file mode 100644 index 000000000..696935996 --- /dev/null +++ b/tests/core/transformation_manager/conftest.py @@ -0,0 +1,92 @@ +#!/usr/bin/env python3 + +""" +Fixtures for transformation manager tests. +""" + +from __future__ import annotations + +import pytest + +from spatialdata.transformations import Affine, Scale, Translation +from spatialdata.transformations.ngff.ngff_coordinate_system import NgffAxis, NgffCoordinateSystem + + +def get_ngff_coodinate_system(cs_name: str) -> NgffCoordinateSystem: + return NgffCoordinateSystem( + name=cs_name, + axes=[NgffAxis(name="x", type="space", unit="micrometer"), NgffAxis(name="y", type="space", unit="micrometer")], + ) + + +@pytest.fixture +def one_point_graph() -> tuple[list[NgffCoordinateSystem], list[Translation]]: + """Fixture providing a single point graph with one coordinate system and one transformation.""" + coordinate_systems = [get_ngff_coodinate_system("cs1")] + transformations = [Translation(translation=[1.0, 2.0], axes=("x", "y"))] + return coordinate_systems, transformations + + +@pytest.fixture +def fully_connected_two_point_graph() -> tuple[list[NgffCoordinateSystem], list[Translation]]: + """Fixture providing a fully connected two-point graph with two coordinate systems and transformations.""" + coordinate_systems = [get_ngff_coodinate_system("cs1"), get_ngff_coodinate_system("cs2")] + transformations = [ + Translation(translation=[-1.0, -2.0], axes=("x", "y")), + ] + return coordinate_systems, transformations + + +@pytest.fixture +def four_point_graph() -> tuple[list[NgffCoordinateSystem], list[Scale | Translation]]: + """Fixture providing a four-point graph with four coordinate systems and transformations.""" + coordinate_systems = [ + get_ngff_coodinate_system("cs1"), + get_ngff_coodinate_system("cs2"), + get_ngff_coodinate_system("cs3"), + get_ngff_coodinate_system("cs4"), + ] + transformations = [ + Translation(translation=[1.0, 2.0], axes=("x", "y")), # cs1 -> cs2 + Translation(translation=[3.0, 4.0], axes=("x", "y")), # cs2 -> cs3 + Scale( + scale=[ + 2, + 2, + ], + axes=("x", "y"), + ), # cs3 -> cs4; + Translation(translation=[4.0, 6.0], axes=("x", "y")), # cs1 -> cs3 (consistent with cs1->cs2 and cs2->cs3) + ] + return coordinate_systems, transformations + + +@pytest.fixture +def five_point_graph() -> tuple[list[NgffCoordinateSystem], list[Scale | Translation | Affine]]: + """Fixture providing a five-point graph with five coordinate systems and five transformations.""" + + coordinate_systems = [ + get_ngff_coodinate_system("cs1"), + get_ngff_coodinate_system("cs2"), + get_ngff_coodinate_system("cs3"), + get_ngff_coodinate_system("cs4"), + get_ngff_coodinate_system("cs5"), + ] + transformations = [ + Translation(translation=[1.0, 2.0], axes=("x", "y")), # cs1 -> cs2, cs4 -> cs3 + Translation(translation=[3.0, 4.0], axes=("x", "y")), # cs2 -> cs3, cs1 -> cs4 + Translation(translation=[4.0, 6.0], axes=("x", "y")), # cs1 -> cs3 (consistent with cs1->cs2 and cs2->cs3) + Scale( + scale=[ + 2, + 2, + ], + axes=("x", "y"), + ), # cs3 -> cs5; + Affine( + matrix=[[2, 0, 0], [0, 2, 0], [0, 0, 1]], + input_axes=("x", "y"), + output_axes=("x", "y"), + ), # cs3 -> cs5 + ] + return coordinate_systems, transformations diff --git a/tests/core/transformation_manager/test_transformation_manager.py b/tests/core/transformation_manager/test_transformation_manager.py new file mode 100644 index 000000000..c1f5d8300 --- /dev/null +++ b/tests/core/transformation_manager/test_transformation_manager.py @@ -0,0 +1,643 @@ +#!/usr/bin/env python3 + +""" +Unit tests for the TransformationManager class in _core/transformation_manager/__init__.py. + +This module tests all the functionality of the TransformationManager class to achieve 100% coverage. +""" + +from __future__ import annotations + +import pytest + +from spatialdata._core.transformation_manager import TRANSFORM_KEY, TransformationManager +from spatialdata._core.transformation_manager.exceptions import ( + CannotRemoveCoordinateSystemError, + CoordinateSystemAlreadyExistsError, + CoordinateSystemNotFoundError, + ElementAlreadyExistsError, + ElementNotRegisteredToAnyCoordinateSystemError, + InvalidPathError, + TransformationNotFoundError, + TransformationPathAmbiguousError, + TransformationPathAmbiguousMultipleEdgeExpectedError, + TransformationPathAmbiguousNoEdgeExpectedError, + TransformationPathNotFoundError, + TransformationPathNotSimple, +) +from spatialdata.transformations.transformations import Sequence + + +def test_initialization(): + """Test that TransformationManager initializes correctly.""" + + tm = TransformationManager() + assert len(tm._graph.nodes()) == 0 + assert len(tm._graph.edges()) == 0 + assert tm._element_to_cs_mapping == {} + + +def test_add_coordinate_system(one_point_graph): + """Test adding a coordinate system.""" + + tm = TransformationManager() + [cs1], _ = one_point_graph + + tm.add_coordinate_system(cs1) + assert cs1 in tm._graph.nodes() + + +def test_add_coordinate_system_duplicate(one_point_graph): + """Test adding an already existing coordinate system""" + tm = TransformationManager() + [cs1], _ = one_point_graph + + tm.add_coordinate_system(cs1) + with pytest.raises(CoordinateSystemAlreadyExistsError, match=f"Coordinate system '{cs1.name}' already exists"): + tm.add_coordinate_system(cs1) + + +def test_remove_coordinate_system(one_point_graph): + """Test removing a coordinate system.""" + tm = TransformationManager() + [cs1], _ = one_point_graph + + tm.add_coordinate_system(cs1) + tm.remove_coordinate_system(cs1) + + assert cs1 not in tm._graph.nodes() + + +def test_remove_coordinate_system_nonexistent(one_point_graph): + """Test that removing a non-existent coordinate system raises CoordinateSystemNotFoundError.""" + tm = TransformationManager() + [cs1], _ = one_point_graph + + # Add the coordinate system first + tm.add_coordinate_system(cs1) + + # Remove it + tm.remove_coordinate_system(cs1) + + # Try to remove it again - should raise CoordinateSystemNotFoundError + with pytest.raises(CoordinateSystemNotFoundError): + tm.remove_coordinate_system(cs1) + + +def test_remove_coordinate_system_with_associated_transformations(fully_connected_two_point_graph): + """ + Test that removing a coordinate system with associated transformations raises CannotRemoveCoordinateSystemError. + """ + tm = TransformationManager() + [cs1, cs2], [transformation] = fully_connected_two_point_graph + + tm.add_coordinate_system(cs1) + tm.add_coordinate_system(cs2) + + tm.add_transformation(cs1, cs2, transformation) + + # Should raise CannotRemoveCoordinateSystem when trying to remove a coordinate system with transformations + with pytest.raises(CannotRemoveCoordinateSystemError, match="Cannot remove coordinate system"): + tm.remove_coordinate_system(cs1) + + # Should raise CannotRemoveCoordinateSystem when trying to remove a coordinate system with transformations + with pytest.raises(CannotRemoveCoordinateSystemError, match="Cannot remove coordinate system"): + tm.remove_coordinate_system(cs2) + + +def test_remove_coordinate_system_with_belonging_elements(one_point_graph): + """Test that removing a coordinate system with element associations raises CannotRemoveCoordinateSystemError.""" + tm = TransformationManager() + [cs1], _ = one_point_graph + tm.add_coordinate_system(cs1) + tm.add_element("image1", cs1) + + with pytest.raises(CannotRemoveCoordinateSystemError, match="Cannot remove coordinate system"): + tm.remove_coordinate_system(cs1) + + +def test_list_coordinate_systems(fully_connected_two_point_graph): + """Test listing all coordinate systems.""" + tm = TransformationManager() + [cs1, cs2], _ = fully_connected_two_point_graph + tm.add_coordinate_system(cs1) + tm.add_coordinate_system(cs2) + + systems = tm.list_coordinate_systems() + assert len(systems) == 2 + assert cs1 in systems + assert cs2 in systems + + +def test_add_element(one_point_graph): + """Test adding an element with an existing coordinate system.""" + tm = TransformationManager() + [cs1], _ = one_point_graph + tm.add_coordinate_system(cs1) + tm.add_element("image1", cs1) + + mapping = tm._element_to_cs_mapping + assert "image1" in mapping + assert mapping["image1"] == cs1 + + +def test_add_element_duplicate(one_point_graph): + """Test that adding a duplicate element raises ElementAlreadyExistsError.""" + tm = TransformationManager() + [cs1], _ = one_point_graph + tm.add_coordinate_system(cs1) + element_name = "image1" + tm.add_element(element_name, cs1) + + # Try to add the same element again - should raise ElementAlreadyExistsError + with pytest.raises( + ElementAlreadyExistsError, match=f"Element '{element_name}' already exists in the transformation manager" + ): + tm.add_element(element_name, cs1) + + +def test_unset_element(one_point_graph): + """Test unsetting an element.""" + tm = TransformationManager() + [cs1], _ = one_point_graph + tm.add_coordinate_system(cs1) + tm.add_element("image1", cs1) + tm._unset_element("image1") + + assert "image1" not in tm._element_to_cs_mapping + + +def test_unset_element_nonexistent(one_point_graph): + """Test that unsetting a non-existent element raises ElementNotFoundError.""" + tm = TransformationManager() + [cs1], _ = one_point_graph + tm.add_coordinate_system(cs1) + element_name = "image1" + + # Try to unset non-existent element + with pytest.raises(ElementNotRegisteredToAnyCoordinateSystemError): + tm._unset_element(element_name) + + +def test_add_transformation(fully_connected_two_point_graph): + """Test adding a transformation between coordinate systems.""" + + tm = TransformationManager() + [cs1, cs2], [transform] = fully_connected_two_point_graph + + tm.add_coordinate_system(cs1) + tm.add_coordinate_system(cs2) + + tm.add_transformation(cs1, cs2, transform) + + assert tm._graph.has_edge(cs1, cs2) + # Get the first edge key and check the transformation + edge_key = tm._get_edge_key((cs1, cs2, transform)) + assert tm._graph[cs1][cs2][edge_key][TRANSFORM_KEY] == transform + + +def test_add_transformation_existing_edge(fully_connected_two_point_graph): + """Test adding a transformation between coordinate systems that already have an edge.""" + tm = TransformationManager() + [cs1, cs2], [transform] = fully_connected_two_point_graph + tm.add_coordinate_system(cs1) + tm.add_coordinate_system(cs2) + + tm.add_transformation(cs1, cs2, transform) + assert tm._graph.has_edge(cs1, cs2) + + # Add same transformation again between the same coordinate systems + tm.add_transformation(cs1, cs2, transform) + # Should have just rewritten the edge + assert tm._graph.has_edge(cs1, cs2) + assert len(tm._graph[cs1][cs2]) == 1 + + +def test_add_transformation_nonexistent_cs(fully_connected_two_point_graph): + """Test that adding a transformation with non-existent coordinate systems raises CoordinateSystemNotFoundError.""" + tm = TransformationManager() + [cs1, cs2], [transform] = fully_connected_two_point_graph + + expected_msg = "Coordinate system .* not found in the transformation manager." + with pytest.raises(CoordinateSystemNotFoundError, match=expected_msg): + tm.add_transformation(cs1, cs2, transform) + + # Add one coordinate system + tm.add_coordinate_system(cs1) + + with pytest.raises(CoordinateSystemNotFoundError, match=expected_msg): + tm.add_transformation(cs1, cs2, transform) + + +def test_get_direct_transformations(fully_connected_two_point_graph): + """Test getting an existing transformation.""" + tm = TransformationManager() + [cs1, cs2], [transform] = fully_connected_two_point_graph + tm.add_coordinate_system(cs1) + tm.add_coordinate_system(cs2) + + tm.add_transformation(cs1, cs2, transform) + + retrieved = tm.get_direct_transformations(cs1, cs2) + assert retrieved == [transform] + + +def test_get_direct_transformations_nonexistent(fully_connected_two_point_graph): + """Test getting a non-existent transformation.""" + tm = TransformationManager() + coordinate_systems, _transformations = fully_connected_two_point_graph + cs1, cs2 = coordinate_systems + tm.add_coordinate_system(cs1) + tm.add_coordinate_system(cs2) + + retrieved = tm.get_direct_transformations(cs1, cs2) + assert retrieved == [] + + +def test_remove_all_transformations_between_coordinate_systems(fully_connected_two_point_graph, mocker): + """Test removing all transformations between coordinate systems.""" + + tm = TransformationManager() + [cs1, cs2], [transform] = fully_connected_two_point_graph + tm.add_coordinate_system(cs1) + tm.add_coordinate_system(cs2) + + tm.add_transformation(cs1, cs2, transform) + + spy = mocker.spy(tm._graph, "remove_edge") + tm.remove_all_transformations_between_coordinate_systems(cs1, cs2) + assert spy.call_count == 1 + assert not tm._graph.has_edge(cs1, cs2) + + +def test_remove_all_transformation_nonexistent(fully_connected_two_point_graph, mocker): + """Test that removing non-existent transformations between coordinate systems raises TransformationNotFoundError.""" + + tm = TransformationManager() + [cs1, cs2], _ = fully_connected_two_point_graph + tm._graph.add_node(cs1) + tm._graph.add_node(cs2) + + spy = mocker.spy(tm._graph, "remove_edge") + tm.remove_all_transformations_between_coordinate_systems(cs1, cs2) + assert spy.call_count == 0 + + +def test_remove_specific_transformation_between_coordinate_systems(five_point_graph): + """Test removing specific transformation between coordinate systems.""" + + tm = TransformationManager() + [_cs1, _cs2, cs3, _cs4, cs5], [_transform1, _transform2, _transform3, transform4, transform5] = five_point_graph + tm.add_coordinate_system(cs3) + tm.add_coordinate_system(cs5) + + tm.add_transformation(cs3, cs5, transform4) + tm.add_transformation(cs3, cs5, transform5) + + tm.remove_specific_transformation(cs3, cs5, transform4) + edge_key = tm._get_edge_key((cs3, cs5, transform4)) + assert not tm._graph.has_edge(cs3, cs5, key=edge_key) + + +def test_remove_specific_transformation_between_coordinate_systems_non_existent(five_point_graph): + """ + Test that removing non-existent specific transformation between coordinate systems raises + TransformationNotFoundError. + """ + + tm = TransformationManager() + [_cs1, _cs2, cs3, _cs4, cs5], [_transform1, _transform2, _transform3, transform4, _transform5] = five_point_graph + tm.add_coordinate_system(cs3) + tm.add_coordinate_system(cs5) + + with pytest.raises( + TransformationNotFoundError, match=f"Transformation from '{cs3.name}' to '{cs5.name}' not found" + ): + tm.remove_specific_transformation(cs3, cs5, transform4) + + +def test_get_all_shortest_transformation_sequences_direct(four_point_graph): + """Test getting all shortest transformation sequences for a direct transformation path.""" + tm = TransformationManager() + [cs1, cs2, cs3, cs4], [transform1, transform2, _transform3, transform4] = four_point_graph + + tm.add_coordinate_system(cs1) + tm.add_coordinate_system(cs2) + tm.add_coordinate_system(cs3) + tm.add_coordinate_system(cs4) + + tm.add_transformation(cs1, cs2, transform1) + tm.add_transformation(cs2, cs3, transform2) + tm.add_transformation(cs1, cs3, transform4) + + sequences = tm.get_all_shortest_transformation_sequences(cs1, cs3) + assert len(sequences) == 1 + assert Sequence([transform4]) in sequences + + +def test_get_all_shortest_transformation_sequences_indirect_one_path(four_point_graph): + """Test getting all shortest transformation sequences for one indirect transformation path.""" + tm = TransformationManager() + [cs1, cs2, cs3, _cs4], [transform1, transform2, _transform3, _transform4] = four_point_graph + tm.add_coordinate_system(cs1) + tm.add_coordinate_system(cs2) + tm.add_coordinate_system(cs3) + + tm.add_transformation(cs1, cs2, transform1) + tm.add_transformation(cs2, cs3, transform2) + + sequences = tm.get_all_shortest_transformation_sequences(cs1, cs3) + assert len(sequences) == 1 + assert Sequence([transform1, transform2]) in sequences + + +def test_get_all_shortest_transformation_sequences_indirect_two_paths(five_point_graph): + """Test getting all shortest transformation sequences for two indirect transformation paths.""" + tm = TransformationManager() + ([cs1, cs2, cs3, cs4, _cs5], [transform1, transform2, _transform3, _transform4, _transform5]) = five_point_graph + tm.add_coordinate_system(cs1) + tm.add_coordinate_system(cs2) + tm.add_coordinate_system(cs3) + tm.add_coordinate_system(cs4) + + tm.add_transformation(cs1, cs2, transform1) + tm.add_transformation(cs2, cs3, transform2) + tm.add_transformation(cs1, cs4, transform2) + tm.add_transformation(cs4, cs3, transform1) + + sequences = tm.get_all_shortest_transformation_sequences(cs1, cs3) + assert len(sequences) == 2 + assert Sequence([transform1, transform2]) in sequences + assert Sequence([transform2, transform1]) in sequences + + +def test_get_all_shortest_transformation_sequences_no_path(four_point_graph): + """Test that getting all shortest transformation sequences with no path raises TransformationPathNotFoundError.""" + tm = TransformationManager() + [cs1, cs2, cs3, cs4], [transform1, transform2, _transform3, _transform4] = four_point_graph + tm.add_coordinate_system(cs1) + tm.add_coordinate_system(cs2) + tm.add_coordinate_system(cs3) + tm.add_coordinate_system(cs4) + + tm.add_transformation(cs1, cs2, transform1) + tm.add_transformation(cs2, cs3, transform2) + + expected_error_msg = f"No transformation path found from {cs1.name} to {cs4.name}" + with pytest.raises(TransformationPathNotFoundError, match=expected_error_msg): + tm.get_all_shortest_transformation_sequences(cs1, cs4) + + +def test_get_all_shortest_transformation_sequences_multiple_paths_multiple_edges_success(five_point_graph): + """Test getting all shortest transformation sequences, with multiple path and multiple edges between nodes.""" + tm = TransformationManager() + [cs1, cs2, cs3, cs4, cs5], [transform1, transform2, transform3, transform4, transform5] = five_point_graph + tm.add_coordinate_system(cs1) + tm.add_coordinate_system(cs2) + tm.add_coordinate_system(cs3) + tm.add_coordinate_system(cs4) + tm.add_coordinate_system(cs5) + + tm.add_transformation(cs1, cs2, transform1) + tm.add_transformation(cs2, cs3, transform2) + tm.add_transformation(cs1, cs4, transform2) + tm.add_transformation(cs4, cs3, transform1) + # tm.add_transformation(cs1, cs3, transform3) + tm.add_transformation(cs3, cs5, transform4) + tm.add_transformation(cs3, cs5, transform5) + + expected_intermediate_edges = [(cs3, cs5, transform4)] + + sequences = tm.get_all_shortest_transformation_sequences( + cs1, cs5, expected_intermediate_edges=expected_intermediate_edges + ) + assert len(sequences) == 2 + assert Sequence([transform1, transform2, transform4]) in sequences + assert Sequence([transform2, transform1, transform4]) in sequences + + +def test_get_all_shortest_transformation_sequences_multiple_paths_multiple_edges_failure(five_point_graph): + """Test getting all shortest transformation sequences, with multiple path and multiple edges between nodes.""" + tm = TransformationManager() + [cs1, cs2, cs3, cs4, cs5], [transform1, transform2, transform3, transform4, transform5] = five_point_graph + tm.add_coordinate_system(cs1) + tm.add_coordinate_system(cs2) + tm.add_coordinate_system(cs3) + tm.add_coordinate_system(cs4) + tm.add_coordinate_system(cs5) + + tm.add_transformation(cs1, cs2, transform1) + tm.add_transformation(cs2, cs3, transform2) + tm.add_transformation(cs1, cs4, transform2) + tm.add_transformation(cs4, cs3, transform1) + tm.add_transformation(cs1, cs3, transform3) + tm.add_transformation(cs3, cs5, transform4) + tm.add_transformation(cs3, cs5, transform5) + + expected_msg = f"Transformation Path ambiguous from {cs1.name} to {cs5.name}." + with pytest.raises(TransformationPathAmbiguousError, match=expected_msg): + tm.get_all_shortest_transformation_sequences(cs1, cs5) + + +def test_get_all_transformation_sequences(four_point_graph): + """Test getting all transformation sequences.""" + tm = TransformationManager() + [cs1, cs2, cs3, cs4], [transform1, transform2, transform3, transform4] = four_point_graph + tm.add_coordinate_system(cs1) + tm.add_coordinate_system(cs2) + tm.add_coordinate_system(cs3) + tm.add_coordinate_system(cs4) + + tm.add_transformation(cs1, cs2, transform1) + tm.add_transformation(cs2, cs3, transform2) + tm.add_transformation(cs3, cs4, transform3) + tm.add_transformation(cs1, cs3, transform4) + + sequences = tm.get_all_transformation_sequences(cs1, cs4) + assert len(sequences) == 2 + assert Sequence([transform1, transform2, transform3]) in sequences + assert Sequence([transform4, transform3]) in sequences + + +def test_get_all_transformation_sequences_multiple_paths(five_point_graph): + """Test getting all transformation sequences, with multiple path and multiple edges between nodes.""" + tm = TransformationManager() + [cs1, cs2, cs3, cs4, _cs5], [transform1, transform2, transform3, _transform4, _transform5] = five_point_graph + tm.add_coordinate_system(cs1) + tm.add_coordinate_system(cs2) + tm.add_coordinate_system(cs3) + tm.add_coordinate_system(cs4) + + tm.add_transformation(cs1, cs2, transform1) + tm.add_transformation(cs2, cs3, transform2) + tm.add_transformation(cs1, cs4, transform2) + tm.add_transformation(cs4, cs3, transform1) + tm.add_transformation(cs1, cs3, transform3) + + sequences = tm.get_all_transformation_sequences(cs1, cs3) + assert len(sequences) == 3 + assert Sequence([transform1, transform2]) in sequences + assert Sequence([transform2, transform1]) in sequences + assert Sequence([transform3]) in sequences + + +def test_get_all_transformation_sequences_multiple_paths_multiple_edges_success(five_point_graph): + """Test getting all transformation sequences, with multiple paths and multiple edges between nodes.""" + tm = TransformationManager() + [cs1, cs2, cs3, cs4, cs5], [transform1, transform2, transform3, transform4, transform5] = five_point_graph + tm.add_coordinate_system(cs1) + tm.add_coordinate_system(cs2) + tm.add_coordinate_system(cs3) + tm.add_coordinate_system(cs4) + tm.add_coordinate_system(cs5) + + tm.add_transformation(cs1, cs2, transform1) + tm.add_transformation(cs2, cs3, transform2) + tm.add_transformation(cs1, cs4, transform2) + tm.add_transformation(cs4, cs3, transform1) + tm.add_transformation(cs1, cs3, transform3) + tm.add_transformation(cs3, cs5, transform4) + tm.add_transformation(cs3, cs5, transform5) + + expected_intermediate_edges = [(cs3, cs5, transform4)] + sequences = tm.get_all_transformation_sequences(cs1, cs5, expected_intermediate_edges=expected_intermediate_edges) + assert len(sequences) == 3 + assert Sequence([transform1, transform2, transform4]) in sequences + assert Sequence([transform2, transform1, transform4]) in sequences + assert Sequence([transform3, transform4]) in sequences + + +def test_get_all_transformation_sequences_multiple_paths_multiple_edges_no_edge_expected(five_point_graph): + """Test getting all transformation sequences with multiple paths and multiple edges, no edge expected.""" + tm = TransformationManager() + [cs1, cs2, cs3, cs4, cs5], [transform1, transform2, transform3, transform4, transform5] = five_point_graph + tm.add_coordinate_system(cs1) + tm.add_coordinate_system(cs2) + tm.add_coordinate_system(cs3) + tm.add_coordinate_system(cs4) + tm.add_coordinate_system(cs5) + + tm.add_transformation(cs1, cs2, transform1) + tm.add_transformation(cs2, cs3, transform2) + tm.add_transformation(cs1, cs4, transform2) + tm.add_transformation(cs4, cs3, transform1) + tm.add_transformation(cs1, cs3, transform3) + tm.add_transformation(cs3, cs5, transform4) + tm.add_transformation(cs3, cs5, transform5) + + expected_msg = f"Transformation Path ambiguous from {cs3.name} to {cs5.name}." + with pytest.raises(TransformationPathAmbiguousNoEdgeExpectedError, match=expected_msg): + tm.get_all_transformation_sequences(cs1, cs5) + + +def test_get_all_transformation_sequences_multiple_paths_multiple_edges_multiple_expected(five_point_graph): + """Test getting all transformation sequences with multiple paths and multiple edges, multiple edges expected.""" + tm = TransformationManager() + [cs1, cs2, cs3, cs4, cs5], [transform1, transform2, transform3, transform4, transform5] = five_point_graph + tm.add_coordinate_system(cs1) + tm.add_coordinate_system(cs2) + tm.add_coordinate_system(cs3) + tm.add_coordinate_system(cs4) + tm.add_coordinate_system(cs5) + + tm.add_transformation(cs1, cs2, transform1) + tm.add_transformation(cs2, cs3, transform2) + tm.add_transformation(cs1, cs4, transform2) + tm.add_transformation(cs4, cs3, transform1) + tm.add_transformation(cs1, cs3, transform3) + tm.add_transformation(cs3, cs5, transform4) + tm.add_transformation(cs3, cs5, transform5) + + expected_intermediate_edges = [ + (cs3, cs5, transform4), + (cs3, cs5, transform5), + ] + + expected_msg = f"Transformation Path ambiguous from {cs3.name} to {cs5.name}." + with pytest.raises(TransformationPathAmbiguousMultipleEdgeExpectedError, match=expected_msg): + tm.get_all_transformation_sequences(cs1, cs5, expected_intermediate_edges=expected_intermediate_edges) + + +def test_get_transformation_sequences_from_simple_paths_after_disambiguation_with_non_simple_paths( + fully_connected_two_point_graph, +): + """Test _get_transformation_sequences_from_simple_paths_after_disambiguation with non-simple paths.""" + tm = TransformationManager() + [cs1, cs2], _ = fully_connected_two_point_graph + + tm.add_coordinate_system(cs1) + + paths = [[cs1, cs2, cs1, cs2]] + + with pytest.raises(TransformationPathNotSimple): + tm._get_transformation_sequences_from_simple_paths_after_disambiguation(paths, ()) + + +def test_get_transformation_sequences_from_simple_paths_after_disambiguation_with_single_node_path( + fully_connected_two_point_graph, +): + """Test _get_transformation_sequences_from_simple_paths_after_disambiguation with a path with only one node.""" + tm = TransformationManager() + [cs1, cs2], [transformation] = fully_connected_two_point_graph + + tm.add_coordinate_system(cs1) + + paths = [[cs1]] + with pytest.raises(InvalidPathError): + tm._get_transformation_sequences_from_simple_paths_after_disambiguation(paths, ()) + + +def test_get_all_transformation_sequences_no_path(four_point_graph): + """Test getting all transformation sequences.""" + tm = TransformationManager() + [cs1, cs2, cs3, cs4], [transform1, transform2, transform3, transform4] = four_point_graph + tm.add_coordinate_system(cs1) + tm.add_coordinate_system(cs2) + tm.add_coordinate_system(cs3) + tm.add_coordinate_system(cs4) + + expected_error_msg = f"No transformation path found from {cs1.name} to {cs4.name}" + with pytest.raises(TransformationPathNotFoundError, match=expected_error_msg): + tm.get_all_shortest_transformation_sequences(cs1, cs4) + + +def test_get_element_coordinate_system_success(one_point_graph): + """Test successfully getting an element's coordinate system.""" + tm = TransformationManager() + [cs1], _ = one_point_graph + tm.add_coordinate_system(cs1) + tm.add_element("image1", cs1) + + # This should successfully return the coordinate system + result = tm.get_element_coordinate_system("image1") + assert result == cs1 + + +def test_repr_method_comprehensive(fully_connected_two_point_graph): + """Test TransformationManager.__repr__ method comprehensively.""" + + # Test empty TransformationManager + tm_empty = TransformationManager() + repr_empty = repr(tm_empty) + assert "TransformationManager" in repr_empty + assert "coordinate_systems=[]" in repr_empty + assert "coordinate_transforms=[]" in repr_empty + assert "elements=[]" in repr_empty + + # Test TransformationManager with data + tm = TransformationManager() + [cs1, cs2], [transform] = fully_connected_two_point_graph + tm.add_coordinate_system(cs1) + tm.add_coordinate_system(cs2) + + tm.add_transformation(cs1, cs2, transform) + tm.add_element("image1", cs1) + tm.add_element("image2", cs2) + + repr_str = repr(tm) + assert "TransformationManager" in repr_str + assert "cs1" in repr_str + assert "cs2" in repr_str + assert repr_str.find(repr(transform)) + assert "image1" in repr_str + assert "image2" in repr_str From 89b8667b4fce2f28210e7921f47a88bbbc3a3159 Mon Sep 17 00:00:00 2001 From: Ajayrama Kumaraswamy Date: Mon, 10 Aug 2026 10:04:18 +0200 Subject: [PATCH 02/31] 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 db849b7b937e32cb5fca768a9423f77b368d827e Mon Sep 17 00:00:00 2001 From: Tomaz Vieira Date: Fri, 21 Aug 2026 11:25:40 +0000 Subject: [PATCH 03/31] 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 dccc3f6ea0bd5306b66f946bd1d77ed14e0a1832 Mon Sep 17 00:00:00 2001 From: Tomaz Vieira Date: Wed, 2 Sep 2026 16:11:38 +0000 Subject: [PATCH 04/31] 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 950ffd75946735e384a8d38312c33f804b812179 Mon Sep 17 00:00:00 2001 From: Tomaz Vieira Date: Thu, 3 Sep 2026 16:32:06 +0000 Subject: [PATCH 05/31] 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 0a2f374aa5dfd8bfe18710ca937ce603b9111dc9 Mon Sep 17 00:00:00 2001 From: Tomaz Vieira Date: Fri, 4 Sep 2026 09:16:37 +0000 Subject: [PATCH 06/31] 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 d29a2a485da7bace45b811c1aacac718d180bc1f Mon Sep 17 00:00:00 2001 From: Tomaz Vieira Date: Fri, 4 Sep 2026 11:13:27 +0000 Subject: [PATCH 07/31] 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 3618cdd7a2a8d5e799509fa9e72b690ffa2e1de4 Mon Sep 17 00:00:00 2001 From: Tomaz Vieira Date: Tue, 8 Sep 2026 08:59:35 +0000 Subject: [PATCH 08/31] 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 4313c563444599a39e82821a2db7235bc6d78b2b Mon Sep 17 00:00:00 2001 From: Tomaz Vieira Date: Tue, 8 Sep 2026 11:26:18 +0000 Subject: [PATCH 09/31] 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 9c2e1d3a0ddf0aabb52138ef41883538334a376c Mon Sep 17 00:00:00 2001 From: Tomaz Vieira Date: Tue, 8 Sep 2026 11:30:05 +0000 Subject: [PATCH 10/31] 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 f94498ed6460205d9099d259b1e8600b7eb044c2 Mon Sep 17 00:00:00 2001 From: Tomaz Vieira Date: Wed, 9 Sep 2026 11:00:58 +0000 Subject: [PATCH 11/31] 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 c0da25f6252b1b77b76ef25b337f97b4f606d06d Mon Sep 17 00:00:00 2001 From: Tomaz Vieira Date: Wed, 9 Sep 2026 12:08:42 +0000 Subject: [PATCH 12/31] 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 3a943827a48f67b8e644a3ca969510ea93d7357a Mon Sep 17 00:00:00 2001 From: Tomaz Vieira Date: Wed, 9 Sep 2026 12:08:57 +0000 Subject: [PATCH 13/31] 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 129fd2d51beb1f901443ccbbae36fe0b76b0d807 Mon Sep 17 00:00:00 2001 From: Tomaz Vieira Date: Wed, 9 Sep 2026 12:26:17 +0000 Subject: [PATCH 14/31] 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 e050c10bf840f86be06ac062c626b6bc5a4a523c Mon Sep 17 00:00:00 2001 From: Tomaz Vieira Date: Wed, 9 Sep 2026 13:07:48 +0000 Subject: [PATCH 15/31] 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 5836cd0860b94f204995e7692923af431da82d23 Mon Sep 17 00:00:00 2001 From: Tomaz Vieira Date: Wed, 9 Sep 2026 13:14:24 +0000 Subject: [PATCH 16/31] 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 de5718ce678137102b9581154ee79d2af09303e7 Mon Sep 17 00:00:00 2001 From: Tomaz Vieira Date: Wed, 9 Sep 2026 13:15:19 +0000 Subject: [PATCH 17/31] 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 d58038d77be7ee6e81ef260985b9c05616c03d33 Mon Sep 17 00:00:00 2001 From: Tomaz Vieira Date: Wed, 9 Sep 2026 13:51:41 +0000 Subject: [PATCH 18/31] 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 55e1ba587a5293bb8478bf5f300dfa9f29443d5a Mon Sep 17 00:00:00 2001 From: Tomaz Vieira Date: Wed, 9 Sep 2026 15:54:26 +0000 Subject: [PATCH 19/31] 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 d802aa83c66984a819c772ed84e760e09275b1fd Mon Sep 17 00:00:00 2001 From: Tomaz Vieira Date: Wed, 9 Sep 2026 16:03:39 +0000 Subject: [PATCH 20/31] 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 c18df09d0eebb782808e8e839312247e88c232b9 Mon Sep 17 00:00:00 2001 From: Tomaz Vieira Date: Wed, 9 Sep 2026 16:18:46 +0000 Subject: [PATCH 21/31] 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 c10c7583d6240516466c55e858a42e39149ebf25 Mon Sep 17 00:00:00 2001 From: Tomaz Vieira Date: Wed, 9 Sep 2026 16:29:37 +0000 Subject: [PATCH 22/31] 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 b6d4a66f0afb456ad727862a3a3ddb6816c6e414 Mon Sep 17 00:00:00 2001 From: Tomaz Vieira Date: Fri, 11 Sep 2026 08:17:09 +0000 Subject: [PATCH 23/31] 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 ec7f22bd0e9542ad56d5caed77c5533f050804e6 Mon Sep 17 00:00:00 2001 From: Tomaz Vieira Date: Fri, 11 Sep 2026 09:57:01 +0000 Subject: [PATCH 24/31] 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 05a28c7a6c530716121041cc27ce07f87f76b483 Mon Sep 17 00:00:00 2001 From: Tomaz Vieira Date: Fri, 11 Sep 2026 10:15:45 +0000 Subject: [PATCH 25/31] 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 98c1b37bb0c7b3c2c4bb072e88f17e0408927172 Mon Sep 17 00:00:00 2001 From: Tomaz Vieira Date: Thu, 10 Sep 2026 12:35:59 +0000 Subject: [PATCH 26/31] 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 834f921c68881fffbae387cc61e8af8b3a76cf04 Mon Sep 17 00:00:00 2001 From: Tomaz Vieira Date: Thu, 10 Sep 2026 13:49:16 +0000 Subject: [PATCH 27/31] 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]])) From 6ca7e2140c8b688b4b1510289b25f75bd32b2c3b Mon Sep 17 00:00:00 2001 From: Tomaz Vieira Date: Mon, 21 Sep 2026 13:14:34 +0000 Subject: [PATCH 28/31] Bump xarray so we can use RangeIndex --- pyproject.toml | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/pyproject.toml b/pyproject.toml index 04947e348..13ab89cc6 100644 --- a/pyproject.toml +++ b/pyproject.toml @@ -52,7 +52,7 @@ dependencies = [ "spatial-image>=1.2.3", "typing-extensions>=4.8", "universal-pathlib>=0.2.6", - "xarray>=2024.10", + "xarray>=2025.4", "xarray-spatial>=0.3.5", "zarr>=3", ] From 0c3ae6498b754045b81612dbc61967f495c85121 Mon Sep 17 00:00:00 2001 From: Tomaz Vieira Date: Mon, 21 Sep 2026 13:22:23 +0000 Subject: [PATCH 29/31] Removes try-from_model_or_defaul for now --- src/spatialdata/transformations/graph/vert.py | 14 -------------- 1 file changed, 14 deletions(-) diff --git a/src/spatialdata/transformations/graph/vert.py b/src/spatialdata/transformations/graph/vert.py index 62b118b16..a0c2b9b9b 100644 --- a/src/spatialdata/transformations/graph/vert.py +++ b/src/spatialdata/transformations/graph/vert.py @@ -129,20 +129,6 @@ def try_from_model(cls, model: ozm06ct.CoordinateSystem) -> CoordSystem: virtual=False, ) - @classmethod - 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. - - 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 - @property def num_axes(self) -> int: return len(self.axes) From 1951436ae22cb2d3bb1257c88a6600362063ea1c Mon Sep 17 00:00:00 2001 From: Tomaz Vieira Date: Mon, 21 Sep 2026 13:46:08 +0000 Subject: [PATCH 30/31] Emits warning on AffineEdge.inverse when condition is big --- src/spatialdata/transformations/graph/edge.py | 8 ++++++++ 1 file changed, 8 insertions(+) diff --git a/src/spatialdata/transformations/graph/edge.py b/src/spatialdata/transformations/graph/edge.py index 9f20151e6..822cd811b 100644 --- a/src/spatialdata/transformations/graph/edge.py +++ b/src/spatialdata/transformations/graph/edge.py @@ -5,6 +5,7 @@ from abc import ABC, abstractmethod from collections.abc import Mapping, Sequence from typing import Final +from warnings import warn import numpy as np import ome_zarr_models.v06.coordinate_transforms as ozm06trans @@ -231,6 +232,13 @@ def mapping(cls, input: CoordSystem, output: CoordSystem, name: str | None = Non ) def inverse(self, name: str | None = None) -> BaseTransformationEdge | None: + cond = np.linalg.cond(self.affine) + if cond > 1e10: + warn( + f"Inverting an affine matrix with large condition number ({cond:.2e})", + RuntimeWarning, + stacklevel=2, + ) try: # FIXME: I think there are more efficient/precise ways to invert a matrix inv = np.linalg.inv(self.affine) From 475d15b32d38bb87eaffa45f1f07cc5e85613b39 Mon Sep 17 00:00:00 2001 From: Tomaz Vieira Date: Mon, 21 Sep 2026 15:10:27 +0000 Subject: [PATCH 31/31] Uses matrix in exceptions --- src/spatialdata/_core/transformation_manager/exceptions.py | 4 ++-- 1 file changed, 2 insertions(+), 2 deletions(-) diff --git a/src/spatialdata/_core/transformation_manager/exceptions.py b/src/spatialdata/_core/transformation_manager/exceptions.py index 06fa5e96b..86adeb1f5 100644 --- a/src/spatialdata/_core/transformation_manager/exceptions.py +++ b/src/spatialdata/_core/transformation_manager/exceptions.py @@ -267,13 +267,13 @@ def __init__( class DeterminantDifferentFromOne(Exception): def __init__(self, matrix: ArrayLike) -> None: self.matrix = matrix - super().__init__("Matrix does not have det(M) == 1") + super().__init__(f"Matrix does not have det(M) == 1:\n{matrix}") class NotOrthonormalError(Exception): def __init__(self, matrix: ArrayLike) -> None: self.matrix = matrix - super().__init__("Matrix is not orthonormal") + super().__init__(f"Matrix is not orthonormal:\n{matrix}") class EmptyTransformSequenceError(Exception):