Transformation

Transformation associated algorithms.

Transformation objects can be averaged, composed, created or inverted. They are meant to be applied to image, intensity or labelled.

timagetk.algorithms.trsf.NULL_TRSF = array([[0., 0., 0., 0.],        [0., 0., 0., 0.],        [0., 0., 0., 0.],        [0., 0., 0., 1.]])

The null transformation

timagetk.algorithms.trsf.apply_trsf(image, trsf=None, template_img=None, interpolation=None, **kwargs)[source]

Apply a transformation to an image.

Parameters:
  • image (timagetk.SpatialImage or timagetk.LabelledImage) – Image to transform.

  • trsf (Trsf, optional) – Input transformation, default is identity, goes from the result image towards the input image.

  • template_img (timagetk.SpatialImage, optional) – Template image used for the dimensions & voxelsize of the output transformation vectorfield.

  • interpolation ({"nearest", "linear", "cspline", "cellbased"}) – Interpolation method to use. See the “Notes” section for a detailled explanations.

  • isotropic_voxel (float) – If set, ignore given trsf & template_img and performs isometric resampling to given voxelsize.

  • cell_based_sigma (float) – Required when using “cellbased” interpolation method, sigma for cell-based interpolation. In voxel units!

  • params (str, optional) – CLI parameter string used by vt.apply_trsf method.

  • single_channel (str) – Use this with a MultiChannelImage to apply the transformation only on the selected channel. If the template_img is a MultiChannelImage, this selected the reference channel for the template.

Returns:

timagetk.SpatialImage or timagetk.LabelledImage or timagetk.MultiChannelImage or timagetk.TissueImage2D or timagetk.TissueImage3D – The transformed image.

Raises:

TypeError – If image is not a SpatialImage. If trsf is not a Trsf, if not None. If template_img is not a SpatialImage, if not None.

Notes

To apply a transformation to a segmented image, use ‘-cellbased’ in params instead of default DEF_APPLY_TRSF.

Interpolation methods definitions:

  • nearest: nearest neighbor interpolation (for binary or label images)

  • linear: bi- or tri-linear interpolation (for intensity images)

  • cspline: cubic spline (for intensity images)

  • cellbased: cell based interpolation (for label images)

Example

>>> from timagetk.algorithms.trsf import apply_trsf
>>> from timagetk.algorithms.trsf import create_trsf
>>> from timagetk.io import imread
>>> from timagetk.io.util import shared_dataset
>>> from timagetk.synthetic_data.labelled_image import example_layered_sphere_labelled_image
>>> from timagetk.visu.mplt import grayscale_imshow
>>> img_path = shared_dataset("p58")[0]
>>> image = imread(img_path)
>>> print(image)
>>> image_iso = apply_trsf(image, isotropic_voxel=0.3)
>>> print(image_iso)
>>> trsf = create_trsf('random', template_img=image)
>>> trsf.print()
>>> out_image = apply_trsf(image, trsf, interpolation='linear')
>>> grayscale_imshow([image, out_image], 5, title=['Original', 'Transformed'], suptitle="Arbitrary transformation", val_range="auto")
>>> from timagetk import LabelledImage
>>> image = example_layered_sphere_labelled_image()
>>> out_nearest = apply_trsf(image, interpolation='nearest', quiet=True)
>>> out_nearest = apply_trsf(image, trsf, interpolation='nearest')
>>> sigma=0.3
>>> out_cellbased = apply_trsf(image, trsf, interpolation='cellbased', cell_based_sigma=sigma)
>>> grayscale_imshow([image, out_nearest, out_cellbased], 5, title=['Original', 'Neareest interpolation', f'Cell based interpolation (sigma={sigma})'], suptitle="Arbitrary transformation to labelled image", cmap='viridis', val_range="auto")
>>> trsf = create_trsf('sinus', template_img=image, trsf_type='vectorfield')
>>> out_cellbased = apply_trsf(image, trsf, interpolation='cellbased', cell_based_sigma=sigma)
>>> out_cellbased = apply_trsf(image, trsf, template_img=image, interpolation='cellbased', cell_based_sigma=sigma)
timagetk.algorithms.trsf.compose_trsf(list_trsf, template_img=None, **kwargs)[source]

Composition of transformations.

Parameters:
  • list_trsf (list of Trsf) – List of transformations to compose.

  • template_img (timagetk.SpatialImage or timagetk.LabelledImage, optional) – The template to use, required for vectorfields composition.

  • params (str) – CLI parameter string used by vt.compose_trsf method.

Returns:

Trsf – Composition of given transformations.

Raises:

TypeError – If list_trsf is not a list of Trsf. If template_img is not a SpatialImage, when specified.

Examples

>>> from timagetk.algorithms.trsf import compose_trsf
>>> from timagetk.algorithms.trsf import create_trsf
>>> trsf_1 = create_trsf('random', trsf_type='rigid')
>>> trsf_2 = create_trsf('random', trsf_type='affine')
>>> comp_trsf = compose_trsf([trsf_1, trsf_2])
>>> comp_trsf.print()
>>> from timagetk.array_util import random_spatial_image
>>> img = random_spatial_image([15, 40, 40], voxelsize=[0.5, 0.21, 0.21])
>>> trsf_3 = create_trsf('sinus', template_img=img)
>>> comp_trsf = compose_trsf([trsf_1, trsf_2, trsf_3])
>>> comp_trsf.print()
timagetk.algorithms.trsf.create_trsf(trsf='identity', template_img=None, trsf_type='affine', **kwargs)[source]

Create a transformation, can be identity, random or sinus.

Parameters:
  • trsf (str in {"identity", "random", "sinus"}, optional) – Name of the transformation to create, default is identity. Use random with linear trnasformations. Use sinus with non-linear transformation

  • template_img (timagetk.SpatialImage or timagetk.LabelledImage, optional) – Image used to defines the template shape and voxelsize, required with non-linear transformations (vectorfield).

  • trsf_type (str in {"rigid", "affine", "vectorfield"}, optional) – Type of transformation to create, default is DEF_TRSF_TYPE.

  • seed (int) – The seed to use for the random generator of linear transformations.

  • angle_range ([float, float]) – Angle range of the linear tranformation, expressed in radian. Limited to “rigid” & “affine” tranformations.

  • translation_range ([float, float]) – Angle range of the linear tranformation, expressed in ??? FIXME. Limited to “rigid” & “affine” tranformations.

  • shear_range ([float, float]) – Shear range of the linear tranformation, expressed in ??? FIXME. Limited to “affine” tranformations.

  • scale_range ([float, float]) – Scale range of the linear tranformation, expressed in ??? FIXME. Limited to “affine” tranformations.

  • sinusoid_amplitude ([float, [float, [float]]]) – Sinusoid amplitude of the non-linear tranformation, expressed in radian. Limited to “vectorfield” tranformations.

  • sinusoid_period ([float, [float, [float]]]) – Sinusoid period of the non-linear tranformation, expressed in radian. Limited to “vectorfield” tranformations.

  • params (str, optional) – CLI parameter string used by vt.create_trsf method.

Returns:

Trsf – The created vt transformation object.

Raises:
  • TypeError – If trsf is not a str. If trsf_type is not a str. If template_img is not a SpatialImage, if not None.

  • ValueError – If trsf is not in CREATE_TRSF. If trsf_type is not in TRSF_TYPE.

Notes

To create “random” linear transformations, use ‘random’.

To create “random” non-linear transformations, use ‘sinus’ and provide a template image.

Examples

>>> from timagetk.algorithms.trsf import create_trsf
>>> # Example #1 - Create identity matrix, *i.e.* no-deformation:
>>> identity_trsf = create_trsf()
>>> identity_trsf.print()
transformation type is AFFINE_3D
transformation unit is REAL_UNIT
   1.000000000000000    0.000000000000000    0.000000000000000    0.000000000000000
   0.000000000000000    1.000000000000000    0.000000000000000    0.000000000000000
   0.000000000000000    0.000000000000000    1.000000000000000    0.000000000000000
   0.000000000000000    0.000000000000000    0.000000000000000    1.000000000000000
>>> # Example #2 - Create random AFFINE transformation matrix:
>>> identity_trsf = create_trsf('random', trsf_type='affine', seed=1)
>>> identity_trsf.print()
transformation type is AFFINE_3D
transformation unit is REAL_UNIT
   0.798876800000000   -0.721781810000000    0.074926290000000    0.268018200000000
   0.315911350000000    0.308944260000000   -0.604184840000000    9.044594500000001
   0.702414430000000    0.746748760000000    0.664103850000000    8.323901360000001
   0.000000000000000    0.000000000000000    0.000000000000000    1.000000000000000
>>> identity_trsf.is_rigid()
False
>>> identity_trsf.is_affine()
True
>>> # Example #3 - Create random RIGID transformation matrix:
>>> trsf = create_trsf('random', trsf_type='rigid', angle_range=[0., 0.25], translation_range=[0., 1.5], seed=1)
>>> trsf.print()
transformation type is AFFINE_3D
transformation unit is REAL_UNIT
   0.991243410000000   -0.126818630000000   -0.036790380000000    1.367471040000000
   0.119926610000000    0.981213850000000   -0.151119090000000    0.296327050000000
   0.055263950000000    0.145383660000000    0.987830700000000    0.502834130000000
   0.000000000000000    0.000000000000000    0.000000000000000    1.000000000000000
>>> # Example #4 - Create random VECTORFIELD transformation:
>>> from timagetk.array_util import random_spatial_image
>>> img = random_spatial_image([15, 40, 40], voxelsize=[0.5, 0.21, 0.21])
>>> trsf = create_trsf('sinus', template_img=img, trsf_type='vectorfield')
>>> trsf.print()
timagetk.algorithms.trsf.inv_trsf(trsf, template_img=None, **kwargs)[source]

Inversion of a transformation matrix.

Parameters:
  • trsf (Trsf) – Transformation object to invert.

  • template_img (timagetk.SpatialImage, optional.) – Control output transformation geometry, used only with vectorfield transformations.

  • error (float) – Absolute error (in real world unit) to determine convergence.

  • iteration (int) – Maximal number of iterations to reach convergence.

  • params (str, optional) – C style command line parameters used by vt.inv_trsf method, DEF_INV_TRSF by default.

Returns:

Trsf – The inverted vt transformation object.

Raises:

TypeError – If trsf is not a Trsf. If template_img is not a SpatialImage, if not None.

Example

>>> from timagetk.algorithms.trsf import create_trsf
>>> from timagetk.algorithms.trsf import inv_trsf
>>> trsf_in = create_trsf('random', trsf_type='rigid')
>>> trsf_in.print()
>>> trsf_out = inv_trsf(trsf_in)
>>> type(trsf_out)
>>> trsf_out.print()
>>> trsf_in = create_trsf('random')
>>> trsf_out = inv_trsf(trsf_in)
>>> type(trsf_out)
>>> trsf_out.print()
>>> from timagetk.array_util import random_spatial_image
>>> img = random_spatial_image([15, 40, 40], voxelsize=[0.5, 0.21, 0.21])
>>> trsf = create_trsf('sinus', template_img=img, trsf_type='vectorfield')
>>> trsf_out = inv_trsf(trsf)
timagetk.algorithms.trsf.trsfs_averaging(list_trsf, method='mean', **kwargs)[source]

Transformations averaging algorithm.

Parameters:
  • list_trsf (list) – List of vt transformations.

  • method ({'mean', 'robust-mean'}) – Averaging method to use.

  • trsf_type ({"rigid", "affine", "vectorfield"}) – Type of transformation to create.

  • fluid_sigma (float or list of float) – Sigma for fluid regularization, i.e. field interpolation and regularization for pairings (only for vector field)

  • estimator_type ({"wlts", "lts", "wls", "ls"}) – Transformation estimator. See the “Notes” section for a detailled explanations.

  • lts_fraction (float) – If defined, set the fraction of pairs that are kept, used with trimmed estimations.

  • lts_deviation (float) – If defined, set the threshold to discard pairings (see notes), used with trimmed estimations.

  • lts_iterations (int) – If defined, set the maximal number of iterations, used with trimmed estimations.

Returns:

Trsf – Averaged transformation.

See also

timagetk.third_party.vt_parser.trsfs_averaging_kwargs, timagetk.third_party.vt_parser.general_kwargs, timagetk.third_party.vt_parser.parallel_kwargs

TODO

create and use temporary filename with out_file_path kwarg when using list of str as input then load it to return Trsf

Example

>>> from timagetk.algorithms.trsf import trsfs_averaging
>>> from timagetk.algorithms.trsf import create_trsf
>>> trsf_1 = create_trsf('random', trsf_type='affine', seed=1)
>>> trsf_2 = create_trsf('random', trsf_type='affine', seed=2)
>>> trsf_out = trsfs_averaging([trsf_1, trsf_2])
>>> trsf_out.print()
>>> trsf_1 = create_trsf('random', trsf_type='rigid', seed=1)
>>> trsf_2 = create_trsf('random', trsf_type='affine', seed=2)
>>> trsf_out = trsfs_averaging([trsf_1, trsf_2])
>>> trsf_out.print()
>>> from timagetk.array_util import random_spatial_image
>>> img = random_spatial_image([15, 40, 40], voxelsize=[0.5, 0.21, 0.21])
>>> trsf_1 = create_trsf('sinus', trsf_type='vectorfield', template_img=img)
>>> trsf_2 = create_trsf('sinus', trsf_type='vectorfield', template_img=img)
>>> trsf_out = trsfs_averaging([trsf_1, trsf_2], trsf_type='vectorfield', fluid_sigma=1.5)
>>> # Use files instead of Trsf:
>>> from tempfile import NamedTemporaryFile as Tmp
>>> trsf_fname_1 = Tmp().name + '.trsf'
>>> trsf_fname_2 = Tmp().name + '.trsf'
>>> from timagetk.algorithms.trsf import create_trsf
>>> create_trsf('random', trsf_type='affine', seed=1).write(trsf_fname_1)
>>> create_trsf('random', trsf_type='affine', seed=2).write(trsf_fname_2)
>>> from vt import mean_trsfs
>>> from tempfile import gettempdir
>>> from os.path import join
>>> mean_fname = join(gettempdir(), 'mean.trsf')
>>> mean_trsfs([trsf_fname_1, trsf_fname_2], out_file_path=mean_fname, params="-trsf-type affine -mean")
>>> from timagetk import Trsf
>>> mean_trsf = Trsf(mean_fname)
>>> print(mean_trsf.get_array())
>>> # Use files instead of Trsf:
>>> from tempfile import NamedTemporaryFile as Tmp
>>> trsf_fname_1 = Tmp().name + '.trsf'
>>> trsf_fname_2 = Tmp().name + '.trsf'
>>> from timagetk.algorithms.trsf import create_trsf
>>> from timagetk.array_util import random_spatial_image
>>> img = random_spatial_image([15, 40, 40], voxelsize=[0.5, 0.21, 0.21])
>>> trsf_1 = create_trsf('sinus', template_img=img, trsf_type='vectorfield')
>>> trsf_1.write(trsf_fname_1)
>>> trsf_2 = create_trsf('sinus', template_img=img, trsf_type='vectorfield')
>>> trsf_2.write(trsf_fname_2)
>>> from vt import mean_trsfs
>>> from tempfile import gettempdir
>>> from os.path import join
>>> mean_fname = join(gettempdir(), 'mean.trsf')
>>> mean_trsfs([trsf_fname_1, trsf_fname_2], out_file_path=mean_fname, params="-trsf-type vectorfield -mean")
>>> from timagetk import Trsf
>>> mean_trsf = Trsf(mean_fname)
>>> mean_trsf.print()
>>> trsf_out = trsfs_averaging([trsf_fname_1, trsf_fname_2])
>>> trsf_out.print()