Add affine transforms

README.md

@@ -113,6 +113,14 @@ imsave("result.tif", cle.pull(labeled))
 
 </td></tr><tr><td>
 
+<img src="https://github.com/clEsperanto/pyclesperanto_prototype/raw/master/docs/images/screenshot_affine_transforms.png" width="300"/>
+
+</td><td>
+
+[Rotation, scaling, translation, affine transforms](https://github.com/clEsperanto/pyclesperanto_prototype/tree/master/demo/transforms/affine_transforms.ipynb)
+
+</td></tr><tr><td>
+
 <img src="https://github.com/clEsperanto/pyclesperanto_prototype/raw/master/docs/images/screenshot_multiply_vectors_and_matrices.png" width="300"/>
 
 </td><td>

docs/release_notes.md

@@ -39,6 +39,11 @@ results of the GPU-accelerated statistics are dictionaries which contain the sam
 * `mode_of_proximal_neighbors_map` and alias `mode_of_distal_neighbors_map`
 * `standard_deviation_of_n_nearest_neighbors_map`
 * `standard_deviation_of_proximal_neighbors_map` and alias `standard_deviation_of_distal_neighbors_map`
+* `affine_transform` (yet without shearing)
+* `translate`
+* `rotate`
+* `scale`
+* `rigid_transform`
 
 ### Backwards compatibility breaking changes
 * `statistics_of_labelled_pixels` and `statistics_of_background_and_labelled_pixels` produce different output now. 
@@ -49,6 +54,9 @@ dictionaries,e.g. `stats['area']` instead of `[s.area for s in stats]`. If not p
 * `n_closest_points` ignores the background position and the distance to labels themselves per default now. 
   Consider passing `ignore_background=False` and `ignore_self=False` to go back to former functionality.
 
+### Deprecation warnings
+* `resample` has been deprecated. Use `scale` instead.
+
 ### Bug fixes
 * `imshow` in 3D caused an error
 * `push_regionprops_column` didn't actually push

pyclesperanto_prototype/__init__.py

@@ -4,6 +4,7 @@
 from ._tier3 import *
 from ._tier4 import *
 from ._tier5 import *
+from ._tier8 import *
 from ._tier9 import *
 
 __version__ = "0.6.0"

pyclesperanto_prototype/_tier1/_resample.py

@@ -30,6 +30,12 @@ def resample(source : Image, destination : Image = None, factor_x : float = 1, f
     ----------
     .. [1] https://clij.github.io/clij2-docs/reference_resample
     """
+    import warnings
+    warnings.warn(
+        "Deprecated: pyclesperanto_prototype.resample() is deprecated. Use scale() instead",
+        DeprecationWarning
+    )
+
     import numpy as np
 
     source_dimensions = source.shape

pyclesperanto_prototype/_tier8/_AffineTransform3D.py

@@ -0,0 +1,179 @@
+import numpy as np
+import transforms3d
+
+class AffineTransform3D:
+    """This class is a convenience class to setup affine transform matrices.
+    When initialized, this object corresponds to a null transform. Afterwards,
+    you can append transforms, e.g. by calling `transform.translate(10, 20)`.
+
+    The API aims to be compatible to Imglib2 AffineTransform3D.
+
+    Inspired by: https://github.com/imglib/imglib2-realtransform/blob/master/src/main/java/net/imglib2/realtransform/AffineTransform3D.java
+
+    """
+
+    def __init__(self):
+        self._matrix = transforms3d.zooms.zfdir2aff(1)
+
+
+    def scale(self, scale_x: float = None, scale_y: float = None, scale_z: float = None):
+        """
+        Scaling the current affine transform matrix.
+
+        Parameters
+        ----------
+        scale_x : float
+            scaling along x axis
+        scale_y : float
+            scaling along y axis
+        scale_z : float
+            scaling along z axis
+
+        Returns
+        -------
+        self
+
+        """
+        if scale_x is not None:
+            self._concatenate(transforms3d.zooms.zfdir2aff(scale_x, direction=(1, 0, 0), origin=(0, 0, 0)))
+        if scale_y is not None:
+            self._concatenate(transforms3d.zooms.zfdir2aff(scale_y, direction=(0, 1, 0), origin=(0, 0, 0)))
+        if scale_z is not None:
+            self._concatenate(transforms3d.zooms.zfdir2aff(scale_z, direction=(0, 0, 1), origin=(0, 0, 0)))
+
+        return self
+
+    def rotate(self, axis : int = 2, angle_in_degrees : float = 0):
+        """
+        Rotation around a given axis (default: z-axis, meaning rotation in x-y-plane)
+
+        Parameters
+        ----------
+        axis : int
+            axis to rotate around (0=x, 1=y, 2=z)
+        angle_in_degrees : int
+            angle in degrees. To convert radians to degrees use this formula:
+            angle_in_deg = angle_in_rad / numpy.pi * 180.0
+
+        Returns
+        -------
+        self
+        """
+        angle_in_rad = angle_in_degrees * np.pi / 180.0
+
+        if axis == 0:
+            self._concatenate(self._3x3_to_4x4(transforms3d.euler.euler2mat(angle_in_rad, 0, 0)))
+        if axis == 1:
+            self._concatenate(self._3x3_to_4x4(transforms3d.euler.euler2mat(0, angle_in_rad, 0)))
+        if axis == 2:
+            self._concatenate(self._3x3_to_4x4(transforms3d.euler.euler2mat(0, 0, angle_in_rad)))
+
+        return self
+
+    def translate(self, translate_x: float = 0, translate_y: float = 0, translate_z : float = 0):
+        """Translation along axes.
+
+        Parameters
+        ----------
+        translate_x : float
+            translation along x-axis
+        translate_y : float
+            translation along y-axis
+        translate_z : float
+            translation along z-axis
+
+        Returns
+        -------
+        self
+
+        """
+        self._concatenate(np.asarray([
+            [1, 0, 0, translate_x],
+            [0, 1, 0, translate_y],
+            [0, 0, 1, translate_z],
+            [0, 0, 0, 1],
+        ]))
+        return self
+
+    def center(self, shape, undo : bool = False):
+        """Change the center of the image to the root of the coordinate system
+
+        Parameters
+        ----------
+        shape : iterable
+            shape of the image which should be centered
+        undo : bool, optional
+            if False (default), the image is moved so that the center of the image is in the root of the coordinate system
+            if True, it is translated in the opposite direction
+
+        Returns
+        -------
+        self
+
+        """
+
+        presign = 1
+        if not undo:
+            presign = -1
+
+        if len(shape) == 2:
+            self.translate(
+                translate_x = presign * shape[1] / 2,
+                translate_y = presign * shape[0] / 2
+            )
+        else: # 3 dimensional image
+            self.translate(
+                translate_x = presign * shape[2] / 2,
+                translate_y = presign * shape[1] / 2,
+                translate_z = presign * shape[0] / 2
+            )
+
+        return self
+
+    def shear(self):
+        raise NotImplementedError("Shearing has not been implemented yet. \n"
+                                  "See https://github.com/clEsperanto/pyclesperanto_prototype/issues/90")
+
+    def _3x3_to_4x4(self, matrix):
+        # I bet there is an easier way to do this.
+        # But I don't know what to google for :D
+        #                            haesleinhuepf
+        return np.asarray([
+            [matrix[0,0], matrix[0,1], matrix[0,2], 0],
+            [matrix[1,0], matrix[1,1], matrix[1,2], 0],
+            [matrix[2,0], matrix[2,1], matrix[2,2], 0],
+            [0, 0, 0, 1]
+        ])
+
+    def _concatenate(self, matrix):
+        self._matrix = np.matmul(matrix, self._matrix)
+
+    def inverse(self):
+        """Computes the inverse of the transformation.
+
+        This can be useful, e.g. when you want to know the transformation
+        from image B to image A but you only know the transformation from
+        A to B.
+
+        Returns
+        -------
+        self
+
+        """
+        self._matrix = np.linalg.inv(self._matrix)
+        return self
+
+    def copy(self):
+        """Makes a copy of the current transform which can then be
+        manipulated without changing the source.
+
+        Returns
+        -------
+        copy of the current transform
+        """
+        a_copy = AffineTransform3D()
+        a_copy._matrix = np.copy(self._matrix)
+        return a_copy
+
+    def __array__(self):
+        return self._matrix

pyclesperanto_prototype/_tier8/__init__.py

@@ -0,0 +1,8 @@
+from ._affine_transform import affine_transform
+from ._AffineTransform3D import AffineTransform3D
+from ._rigid_transform import rigid_transform
+from ._rotate import rotate
+from ._scale import scale
+from ._translate import translate
+
+

pyclesperanto_prototype/_tier8/_affine_transform.py

@@ -0,0 +1,90 @@
+from typing import Union
+
+from .._tier0 import plugin_function
+from .._tier0 import Image
+from .._tier0 import push_zyx
+from ._AffineTransform3D import AffineTransform3D
+from skimage.transform import AffineTransform
+import numpy as np
+
+@plugin_function
+def affine_transform(source : Image, destination : Image = None, transform : Union[np.ndarray, AffineTransform3D, AffineTransform] = None, linear_interpolation : bool = False):
+    """
+    Applies an affine transform to an image.
+
+    Parameters
+    ----------
+    source : Image
+        image to be transformed
+    destination : Image, optional
+        image where the transformed image should be written to
+    transform : 4x4 numpy array or AffineTransform3D object or skimage.transform.AffineTransform object
+        transform matrix or object describing the transformation
+    linear_interpolation: bool
+        not implemented yet
+
+    Returns
+    -------
+    destination
+
+    """
+    import numpy as np
+    from .._tier0 import empty_image_like
+    from .._tier0 import execute
+    from .._tier1 import copy
+    from .._tier0 import create
+    from .._tier1 import copy_slice
+
+    # deal with 2D input images
+    if len(source.shape) == 2:
+        source_3d = create([1, source.shape[0], source.shape[1]])
+        copy_slice(source, source_3d, 0)
+        source = source_3d
+
+    # deal with 2D output images
+    original_destination = destination
+    copy_back_after_transforming = False
+    if len(destination.shape) == 2:
+        destination = create([1, destination.shape[0], destination.shape[1]])
+        copy_slice(original_destination, destination, 0)
+        copy_back_after_transforming = True
+
+        # we invert the transform because we go from the target image to the source image to read pixels
+    if isinstance(transform, AffineTransform3D):
+        transform_matrix = np.asarray(transform.copy().inverse())
+    elif isinstance(transform, AffineTransform):
+        matrix = np.asarray(transform.params)
+        matrix = np.asarray([
+            [matrix[0,0], matrix[0,1], 0, matrix[0,2]],
+            [matrix[1,0], matrix[1,1], 0, matrix[1,2]],
+            [0, 0, 1, 0],
+            [matrix[2,0], matrix[2,1], 0, matrix[2,2]]
+        ])
+        transform_matrix = np.linalg.inv(matrix)
+    else:
+        transform_matrix = np.linalg.inv(transform)
+
+    gpu_transform_matrix = push_zyx(transform_matrix)
+
+    kernel_suffix = ''
+    if linear_interpolation:
+        image = empty_image_like(source)
+        copy(source, image)
+        source = image
+        kernel_suffix = '_interpolate'
+
+
+    parameters = {
+        "input": source,
+        "output": destination,
+        "mat": gpu_transform_matrix
+    }
+
+    execute(__file__, '../clij-opencl-kernels/kernels/affine_transform_' + str(len(destination.shape)) + 'd' + kernel_suffix + '_x.cl',
+            'affine_transform_' + str(len(destination.shape)) + 'd' + kernel_suffix, destination.shape, parameters)
+
+    # deal with 2D output images
+    if copy_back_after_transforming:
+        copy_slice(destination, original_destination, 0)
+
+    return original_destination

pyclesperanto_prototype/_tier8/_rigid_transform.py

@@ -0,0 +1,58 @@
+from .._tier0 import plugin_function
+from .._tier0 import Image
+
+@plugin_function
+def rigid_transform(source : Image, destination : Image = None, translate_x : float = 0, translate_y : float = 0, translate_z : float = 0, angle_around_x_in_degrees : float = 0, angle_around_y_in_degrees : float = 0, angle_around_z_in_degrees : float = 0, rotate_around_center=True):
+    """Translate the image by a given vector and rotate it by given angles.
+
+    Angles are given in degrees. To convert radians to degrees, use this formula:
+
+    angle_in_degrees = angle_in_radians / numpy.pi * 180.0
+
+    Parameters
+    ----------
+    source : Image
+        image to be transformed
+    destination : Image, optional
+        target image
+    translate_x : float
+        translation along x axis in pixels
+    translate_y : float
+        translation along y axis in pixels
+    translate_z : float
+        translation along z axis in pixels
+    angle_around_x_in_degrees : float
+        rotation around x axis in radians
+    angle_around_y_in_degrees : float
+        rotation around y axis in radians
+    angle_around_z_in_degrees : float
+        rotation around z axis in radians
+    rotate_around_center : boolean
+        if True: rotate image around center
+        if False: rotate image around origin
+
+    Returns
+    -------
+    destination
+
+    """
+    from ._AffineTransform3D import AffineTransform3D
+    from ._affine_transform import affine_transform
+
+    transform = AffineTransform3D()
+    if rotate_around_center:
+        transform.center(source.shape)
+
+    if angle_around_x_in_degrees != 0:
+        transform.rotate(0, angle_around_x_in_degrees)
+    if angle_around_y_in_degrees != 0:
+        transform.rotate(1, angle_around_y_in_degrees)
+    if angle_around_z_in_degrees != 0:
+        transform.rotate(2, angle_around_z_in_degrees)
+
+    if rotate_around_center:
+        transform.center(source.shape, undo=True)
+
+    transform.translate(translate_x, translate_y, translate_z)
+
+    return affine_transform(source, destination, transform)