Units and the new transforms classes

[LucaMarconato]

@kevinyamauchi @giovp @ivirshup

I plan to re-add support for storing units, here is how.

At the moment we store the transformations inside .attrs[TRANFORM_KEY] (or .uns[TRANSFOM_KEY] for AnnData objects for storing shapes), in the form of a dict[str, BaseTransformation]. The key is the name of the coordinate system the object aligns to.

I propose to add a new slot .attrs[COORDINATE_SYSTEMS_KEY] of type dict[str, NgffCoordinateSystem] to store unit and other Axis information.

This is similar to the old ways of storing coordinate systems, but with one difference. The transformation make zero use of this information, they work even without. In particular the order of the axes is irrelevant, and the fact that a coordinate system has an axis or not (say the z axes), doesn't influence the transformation, which can be applied anyway.

For instance consider the following:

my_image: SpatialImage = sdata.image['my_image']
assert get_dims(my_image) == ('c', 'y', 'x')
scale = Scale([2., 2.], axes=('x', 'y'))
sdata.set_transformation(image, scale, 'global')
my_cs = NgffCoordinateSystem(
    name='global', axes=[
        Axis(name='x', type='space', unit='meters'), 
        Axis(name='y', type='space', unit='meters'), 
        Axis(name='c', type='channel')
    ]
)
# the name 'global' is redundant, I'll maybe omit this from the function call
sdata.set_coordinate_system(image, 'global', my_cs)  
# in this case, the following line is equivalent to t = scale
t = sdata.map_coordinate_systems(my_image, 'global')
my_image_scaled = t.transform(my_image)
assert get_dims(my_image_scaled) == ('c', 'y', 'x')

What we have is that the information of the order of the axes of 'global' is not used. If we were doing the same for transforming points ('x', 'y', 'z') to 'global', we would obtain points ('x', 'y', 'z'), with x and y coordinates scaled.

Still when visualizing the data (napari/static plotting), or doing spatial queries, we will be able to make use of the knowledge that the 'x' axis is in micrometers. If an axis is specified by the coordinate system (like nothing is said about 'z'), then simply no information regarding 'z' would be retrieved.

Finally, when saving things to NGFF I will write something that conforms to the specs.

[giovp]

thanks @LucaMarconato for the clear explanation.

I think this is an option but, wouldn't it be better to than just have a CoordinateSystem object even as the key of the transform?

e.g.

sdata: SpatialData ...

sdata.images["image0"].attrs[CoordinateSystem(name="czyx"):Affine(....)]

This is similar to the old ways of storing coordinate systems, but with one difference. The transformation make zero use of this information, they work even without. In particular the order of the axes is irrelevant, and the fact that a coordinate system has an axis or not (say the z axes), doesn't influence the transformation, which can be applied anyway.

This I totally agree

[LucaMarconato]

I thought about it but I think that the ergonomics would suffer a lot, especially for all those cases in which the coordinate system could be anything and the user cares just about the name. Also, from a code perspective, having an extra key .attrs[COORDINATE_SYSTEMS_KEY] would be quicker to implement.

[giovp]

I see the point. I think we'd have then to check that all keys in .attrs[COORDINATE_SYSTEMS_KEY] are also in .attrs[TRANSFORM_KEY]

[LucaMarconato]

I would not require that, but I would raise a warning when a coordinate system is specified without a transformation mapping to it.

Also I would raise a warning or an exception in the following case (this applies also to the case in which we use NgffCoordinateSystem as keys of attrs[TRANSFORM_KEY]: different elements contain coordinate systems with the same name but different axes/units.