Registration with a Versor with Python and without SimpleITK

Beginner Questions
,
1

Continuing the discussion from 3D image registration with Python and without SimpleITK:

Hello!

Instead of the translation registration, I am now trying a registration with a versor. How can I do that? This is my code (MWE):

# Copyright 2022 edgar
# License:
#
#    This program is free software: you can redistribute it
#    and/or modify it under the terms of the GNU General
#    Public License as published by the Free Software
#    Foundation, version 3 of the License.
#
#    This program is distributed in the hope that it will be
#    useful, but WITHOUT ANY WARRANTY; without even the
#    implied warranty of MERCHANTABILITY or FITNESS FOR A
#    PARTICULAR PURPOSE. See the GNU General Public License
#    for more details.
#
#    You should have received a copy of the GNU General
#    Public License along with this program. If not, see
#    <https://www.gnu.org/licenses/>.
#
# Copyright Insight Software Consortium
#
# Licensed under the Apache License, Version 2.0 (the
# "License"); you may not use this file except in compliance
# with the License.  You may obtain a copy of the License at
#
#        http://www.apache.org/licenses/LICENSE-2.0.txt
#
# Unless required by applicable law or agreed to in writing,
# software distributed under the License is distributed on
# an "AS IS" BASIS, WITHOUT WARRANTIES OR CONDITIONS OF ANY
# KIND, either express or implied.  See the License for the
# specific language governing permissions and limitations
# under the License.

import sys
import itk
import numpy as np
from distutils.version import StrictVersion as VS
# from benchmark3d.primitives import displ_series_obj3d
from skimage.morphology import disk
from scipy.ndimage import shift


# Draw a 3D disk
side_size = 10
obj2d = disk(side_size / 2, dtype=np.uint8)
obj3d = np.asarray([disk(side_size / 2, dtype=np.uint8)
                    for idx in range(5)])
# Create 4 frames
series_length = 4
holder = np.asarray([obj3d for idx in range(series_length)])
# Move the last frame
holder[-1] = shift(holder[-1],
                   (series_length - 1, 0, 0))

# Set the right data type
pixel_type = itk.ctype("float")
dimension = 3
image_type = itk.Image[pixel_type, dimension]

# Load first and last frames as ITK data (images);
# swap axes ITK: i, j, k, Numpy: k, i, j
# (https://discourse.itk.org/t/
#   importing-image-from-array-and-axis-reorder/1192/)
fixed_image_arr = holder[0].transpose((2, 1, 0))
moving_image_arr = holder[-1].transpose((2, 1, 0))
fixed_image = itk.image_from_array(
    fixed_image_arr.astype(np.float32))
moving_image = itk.image_from_array(
    moving_image_arr.astype(np.float32))

dimension = fixed_image.GetImageDimension()
# fixed_image_type = itk.Image[pixel_type, dimension]
# moving_image_type = itk.Image[pixel_type, dimension]
fixed_image_type = type(fixed_image)
moving_image_type = type(moving_image)

versor_transform = itk.VersorTransform[itk.D]
initial_transform = versor_transform.New()

versor_optimizer = itk.VersorTransformOptimizer.New(
    # # Not for this optimiser
    # LearningRate=4,
    MinimumStepLength=0.001,
    RelaxationFactor=0.5,
    NumberOfIterations=200)

mse_img2img_metric = itk.MeanSquaresImageToImageMetricv4
mse_metric = mse_img2img_metric[fixed_image_type,
                                moving_image_type].New()

registration_method = itk.ImageRegistrationMethodv4[
    fixed_image_type, moving_image_type]
registration = registration_method.New(
        FixedImage=fixed_image,
        MovingImage=moving_image,
        Metric=mse_metric,
        Optimizer=versor_optimizer,
        InitialTransform=initial_transform)

moving_initial_transform = versor_transform.New()
initial_parameters = moving_initial_transform.GetParameters()

, but I get this error:

  File "/usr/lib/python3.10/itk/itkImageRegistrationMethodv4Python.py", line 363, in New
    template_class.New(obj, *args, **kargs)
  File "/usr/lib/python3.10/itk/support/template_class.py", line 800, in New
    itk.set_inputs(self, args, kargs)
  File "/usr/lib/python3.10/itk/support/extras.py", line 1252, in set_inputs
    attrib(itk.output(value))
  TypeError: in method 'itkImageRegistrationMethodv4REGv4F3F3_SetOptimizer', argument 2 of type 'itkObjectToObjectOptimizerBaseTemplateD *'

If I try with double:

  # Copyright 2022 edgar
  # License:
  #
  #    This program is free software: you can redistribute it
  #    and/or modify it under the terms of the GNU General
  #    Public License as published by the Free Software
  #    Foundation, version 3 of the License.
  #
  #    This program is distributed in the hope that it will be
  #    useful, but WITHOUT ANY WARRANTY; without even the
  #    implied warranty of MERCHANTABILITY or FITNESS FOR A
  #    PARTICULAR PURPOSE. See the GNU General Public License
  #    for more details.
  #
  #    You should have received a copy of the GNU General
  #    Public License along with this program. If not, see
  #    <https://www.gnu.org/licenses/>.
  #
  # Copyright Insight Software Consortium
  #
  # Licensed under the Apache License, Version 2.0 (the
  # "License"); you may not use this file except in compliance
  # with the License.  You may obtain a copy of the License at
  #
  #        http://www.apache.org/licenses/LICENSE-2.0.txt
  #
  # Unless required by applicable law or agreed to in writing,
  # software distributed under the License is distributed on
  # an "AS IS" BASIS, WITHOUT WARRANTIES OR CONDITIONS OF ANY
  # KIND, either express or implied.  See the License for the
  # specific language governing permissions and limitations
  # under the License.

  pixel_type = itk.ctype("double")
  dimension = 3
  image_type = itk.Image[pixel_type, dimension]

  # Load first and last frames as ITK data (images);
  # swap axes ITK: i, j, k, Numpy: k, i, j
  # (https://discourse.itk.org/t/
  #   importing-image-from-array-and-axis-reorder/1192/)
  fixed_image_arr = holder[0].transpose((2, 1, 0))
  moving_image_arr = holder[-1].transpose((2, 1, 0))
  fixed_image = itk.image_from_array(
      fixed_image_arr.astype(np.double))
  moving_image = itk.image_from_array(
      moving_image_arr.astype(np.double))

  dimension = fixed_image.GetImageDimension()
  # fixed_image_type = itk.Image[pixel_type, dimension]
  # moving_image_type = itk.Image[pixel_type, dimension]
  fixed_image_type = type(fixed_image)
  moving_image_type = type(moving_image)

  versor_transform = itk.VersorTransform[itk.D]
  initial_transform = versor_transform.New()

  versor_optimizer = itk.VersorTransformOptimizer.New(
      # # Not for this optimiser
      # LearningRate=4,
      MinimumStepLength=0.001,
      RelaxationFactor=0.5,
      NumberOfIterations=200)

  mse_img2img_metric = itk.MeanSquaresImageToImageMetricv4
  mse_metric = mse_img2img_metric[fixed_image_type,
                                  moving_image_type].New()

  registration_method = itk.ImageRegistrationMethodv4[
      fixed_image_type, moving_image_type]
  registration = registration_method.New(
          FixedImage=fixed_image,
          MovingImage=moving_image,
          Metric=mse_metric,
          Optimizer=versor_optimizer,
          InitialTransform=initial_transform)

  moving_initial_transform = versor_transform.New()
  initial_parameters = moving_initial_transform.GetParameters()

I get this

Traceback (most recent call last):
  File "/usr/lib/python3.10/itk/support/template_class.py", line 525, in __getitem__
    this_item = self.__template__[key]
KeyError: (<class 'itk.itkImagePython.itkImageD3'>, <class 'itk.itkImagePython.itkImageD3'>)

During handling of the above exception, another exception occurred:

Traceback (most recent call last):
  File "<string>", line 17, in __PYTHON_EL_eval
  File "<string>", line 1, in <module>
  File "/usr/lib/python3.10/itk/support/template_class.py", line 529, in __getitem__
    raise itk.TemplateTypeError(self, key)
itk.support.extras.TemplateTypeError: itk.MeanSquaresImageToImageMetricv4 is not wrapped for input type `itk.Image[itk.D,3], itk.Image[itk.D,3]`.

To limit the size of the package, only a limited number of
types are available in ITK Python. To print the supported
types, run the following command in your python environment:

    itk.MeanSquaresImageToImageMetricv4.GetTypes()

Possible solutions:
,* If you are an application user:
,** Convert your input image into a supported format (see below).
,** Contact developer to report the issue.
,* If you are an application developer, force input images to be
loaded in a supported pixel type.

    e.g.: instance = itk.MeanSquaresImageToImageMetricv4[itk.Image[itk.F,2], itk.Image[itk.F,2]].New(my_input)

,* (Advanced) If you are an application developer, build ITK Python yourself and
turned to `ON` the corresponding CMake option to wrap the pixel type or image
dimension you need. When configuring ITK with CMake, you can set
`ITK_WRAP_${type}` (replace ${type} with appropriate pixel type such as
`double`). If you need to support images with 4 or 5 dimensions, you can add
these dimensions to the list of dimensions in the CMake variable
`ITK_WRAP_IMAGE_DIMS`.

Supported input types:

itk.Image[itk.F,2]
itk.Image[itk.F,3]

I can compile again with -DITK_WRAP_double:BOOL=ON, but may be float32 is a bit less memory expensive. How can a registration be run with versors and np.float32? Thanks!

2

What happens if you use versor_transform = itk.VersorTransform[itk.F] (transform which uses floats internally)?

3

Hi! Thanks!

If everything stays as float (pixel_type, fixed_image, moving_image), the error is similar, but now, at the end:

...
  File "/usr/lib/python3.10/itk/support/template_class.py", line 529, in __getitem__
    raise itk.TemplateTypeError(self, key)
itk.support.extras.TemplateTypeError: itk.VersorTransform is not wrapped for input type `itk.F`.
...
    e.g.: instance = itk.VersorTransform[itk.D].New(my_input)
...
Supported input types:

itk.D

I can see that my CMakeCache.txt has ITK_WRAP_float:BOOL=ON and -DITK_WRAP_IMAGE_DIMS:STRING="2;3;4"

4

Transforms are normally used with double parameters, so this is not strange. It was worth a shot.

@brad-t-moore, do you have any advice about TypeError: in method 'itkImageRegistrationMethodv4REGv4F3F3_SetOptimizer', argument 2 of type 'itkObjectToObjectOptimizerBaseTemplateD *'?

1 Like
5

As a side note, I think you do not want to transpose the dimensions when creating fixed_image_arr and moving_image_arr. Although ITK conventions have the indices in reverse order relative to numpy, that applies when using itk::Index values and itk::Size values. I suspect that holder[0].shape and fixed_image.shape do not match in the present code, but they should. (That fixed_image.GetLargestPossibleRegion().GetSize() returns the indices in reverse order to these shapes is okay.)

6

Yes, you may be right. I will check that. Thanks.

7

The VersorTransformOptimizer is not a v4 optimizer is my guess.

8

Thanks.

  1. What would be the appropriate way to run a registration with a versor? (what modifications are needed in the code that I posted?).
  2. If a versor is not a feasible option, what other type of transform is available for rotation as a v4 transform?
9

A complete list is here. Try RegularStepGradientDescentOptimizerv4, LBFGSOptimizerv4, PowellOptimizerv4 or AmoebaOptimizerv4.

1 Like
10

Thanks, I will try those :slight_smile: .

11

So the VersorTransform can still be used but you need to use a v4 optimizer. The ITK software guide has a section that looks relevant. Pg. 223 of Book II (second half of the PDF)

https://itk.org/ItkSoftwareGuide.pdf

1 Like
12

Thank you, I will look into that. I can already see that:

13

For those interested, it can also be found here:
https://itk.org/ITKSoftwareGuide/html/Book2/ITKSoftwareGuide-Book2ch3.html#x26-1050003.6