Affine transformations are written to disk different from how they are stored in memory #1284

Activity

1. MeyerBender

   added the

   Bug

   issue type on Oct 23, 2026
2. ajkswamy commented on Nov 2, 2026

   

   ajkswamy

   on Nov 2, 2026

   Collaborator

   More actions

   Hi @MeyerBender Thanks for the detailed bug report.

   an extra dimension gets added implicitly.

   This cannot be avoided when writing to NGFF. It will also be the case once we implement support for NGFF 0.6 and above. In this case, the implicit behaviour of spatialdata --- that the data along the c dimension remains unchanged during transformations, becomes explicity defined when the 3x3 affine matrix gets converted to 4x4.

   It might be better to handle this explicitly, or at least warn the user that the transform is being changed when writing to disk
   The only information being lost during the round trip in your example is that the user did not define a transform along "c". Could you please elaborate more on the case where this information would be needed after a round trip?
3. MeyerBender commented on Nov 2, 2026

   

   MeyerBender

   on Nov 2, 2026

   ContributorAuthor

   More actions

   There's two ways I could imagine dealing with this:

   1. You could simply raise a warning when the spatialdata object is written to disk, just to make users aware that something is being changed in the background.
   2. In the long term, it might be best for spatialdata to only support transforms that are also supported by NGFF. If users try to set an incompatible transform, you could raise an error. For compatibility with older versions, this would likely require a make_transforms_compatible() function or something along those lines
4. ajkswamy commented on Nov 3, 2026

   

   ajkswamy

   on Nov 3, 2026

   Collaborator

   More actions

   Thank you @MeyerBender for possible solutions. Before we look at solutions, I would like to understand the issue better. Sorry that I made a mistake in formatting my previous comment and thereby failed to emphasize my main question, which I repeat here:

   The only information being lost during the round trip in your example is that the user did not define a transform along "c". Could you please elaborate more on the case where this information would be needed after a round trip?
5. MeyerBender commented on Nov 3, 2026

   

   MeyerBender

   on Nov 3, 2026

   Last edited by MeyerBender

   ContributorAuthor

   More actions

   The issue arises when users try to use the transfomration matrix for custom operations. For example, consider the example below, where I get the transformation matrix and apply it to a point. Depending on if I do this before or after writing, I get completely different results. I would argue that the object should not change this significantly just because of a write/read operation.

   import numpy as np
   import spatialdata as sd
   from spatialdata.datasets import blobs
   from spatialdata.transformations import Affine, get_transformation, set_transformation
   
   sdata = blobs()
   img = sdata["blobs_image"]  # dims (c, y, x)
   
   # A plain 2D translation: x += 100, y += 200
   matrix = np.array([[1, 0, 100], [0, 1, 200], [0, 0, 1]], dtype=float)
   t = Affine(matrix, input_axes=("x", "y"), output_axes=("x", "y"))
   set_transformation(img, t, to_coordinate_system="global")
   
   p = np.array([10.0, 20.0, 1.0])  # point (x, y) in homogeneous coordinates
   
   before = get_transformation(img)
   print(before.matrix @ p)  # [110. 220. 1.]  correct
   
   sdata.write("blobs_tmp.zarr", overwrite=True)
   after = get_transformation(sd.read_zarr("blobs_tmp.zarr")["blobs_image"])
   print(after.input_axes, "->", after.output_axes)  # ('c', 'y', 'x') -> ('c', 'x', 'y')
   print(after.matrix @ p)  # ValueError: matmul: ... (size 3 is different from 4)
   

   Another note to make here: simply dropping the c coordinate would silently lead to wrong results, since the x and y coordinate are swapped after the write/read operation.
6. ajkswamy commented on Nov 3, 2026

   

   ajkswamy

   on Nov 3, 2026

   Collaborator

   More actions

   Thank you @MeyerBender for the use case. Now I understand the issue better. This issue will be solved when we provide transform_points method for transformation objects and remove the need for users to directly use the affine matrix to transform point data represented in numpy arrays. Such methods will be added as part of the upcoming refactor of the transformation system. See here for an early draft.

   Another solution would be to use the current PointsModel.parse and transform APIs. The results are dask DataFrame objects, which have very similar API to pandas DataFrame objects. The input to PointsModel.parse can also be DataFrame objects, which would make the code more readable and avoid the implicit interpretation of axis labels and order.

   import numpy as np
   import spatialdata as sd
   from spatialdata import transform
   from spatialdata.datasets import blobs
   from spatialdata.models import PointsModel
   from spatialdata.transformations import Affine, get_transformation, set_transformation
   
   sdata = blobs()
   img = sdata["blobs_image"]  # dims (c, y, x)
   
   # A plain 2D translation: x += 100, y += 200
   matrix = np.array([[1, 0, 100], [0, 1, 200], [0, 0, 1]], dtype=float)
   t = Affine(matrix, input_axes=("x", "y"), output_axes=("x", "y"))
   set_transformation(img, t, to_coordinate_system="global")
   
   # p = np.array([10.0, 20.0, 1.0])  # point (x, y) in homogeneous coordinates
   
   # before = get_transformation(img)
   # print(before.matrix @ p)  # [110. 220. 1.]  correct
   
   points = PointsModel.parse(np.array([[10.0, 20.0]]), transformations={"cs1": t})
   # each row of numpy array is interpreted as a point in xy order
   # see https://spatialdata.scverse.org/en/stable/tutorials/notebooks/notebooks/examples/models2.html#points
   print(points.compute())
   # dask DataFrame
   #       x     y
   # 0  10.0  20.0
   transformed_points_before = transform(points, to_coordinate_system="cs1")
   print(transformed_points_before.compute())
   # dask DataFrame
   #        x      y
   # 0  110.0  220.0
   
   sdata.write("blobs_tmp.zarr", overwrite=True)
   after = get_transformation(sd.read_zarr("blobs_tmp.zarr")["blobs_image"])
   print(after.input_axes, "->", after.output_axes)  # ('c', 'y', 'x') -> ('c', 'x', 'y')
   # print(after.matrix @ p)  # ValueError: matmul: ... (size 3 is different from 4)
   
   set_transformation(points, transformation=after, to_coordinate_system='cs2')
   transformed_points_after = transform(points, to_coordinate_system="cs2")
   print(transformed_points_after.compute())
   # dask DataFrame
   #        x      y
   # 0  110.0  220.0
7. ajkswamy

   changed the issue type from

   Bug

   to

   Enhancement

   on Nov 3, 2026