insightsoftwareconsortium/itkelastix

An ITK Python interface to elastix, a toolbox for rigid and nonrigid registration of images

ITKElastix

Overview

Provides an ITK Python, JavaScript, and WebAssembly interface to elastix, a toolbox for rigid and nonrigid registration of images.

elastix is open source software, based on the well-known Insight Toolkit (ITK). The software consists of a collection of algorithms that are commonly used to solve (medical) image registration problems. The modular design of elastix allows the user to quickly configure, test, and compare different registration methods for a specific application.

👨‍💻 Live JavaScript API Demo ✨

Installation

Install cross-platform native binary Python packages with pip:

pip install itk-elastix

Experimental WebAssembly Python packages can be installed across platforms with:

pip install itkwasm-elastix

Note The API for the WebAssembly and native binary packages are slightly different. For the WebAssembly interface, see the package's Sphinx documentation. For the interface to the native binaries, see the examples/ directory in this repository.

JavaScript / TypesScript packages cas be installed with

npm install @itk-wasm/elastix

Usage

To register two images with the native Python binaries, traditionally called the fixed image and the moving image:

import itk

fixed_image = itk.imread('path/to/fixed_image.mha')
moving_image = itk.imread('path/to/moving_image.mha')

registered_image, params = itk.elastix_registration_method(fixed_image, moving_image)

Interactive examples and tutorial material can be found in the examples directory. Run the examples in free cloud compute containers on MyBinder or clone the repository and run the notebooks locally in Jupyter Notebook or Jupyter Lab. Try out the experimental GPU packages on Paperspace Gradient.

ITKElastix can be used with both the procedural and the object oriented method, as shown in the example notebooks. The procedural method is shorter, less explicit and currently slightly less functional than the object oriented method, however the execution time and output do not differ apart from possible differences due to the stochastic nature of the Elastix algorithm.

To find parameters that work well with specific datasets, see the elastix Model Zoo.

For a graphical user interface in a desktop application, see the napari plugin.

Additional documentation is available for the WebAssembly JavaScript bindings and WebAssembly Python bindings.

Acknowledgements

ITKElastix was developed in part with support from:

• NIH NIMH BRAIN Initiative under award 1RF1MH126732.
• Chan Zuckerberg Initiative (CZI) Essential Open Source Software for Science award for Open Source Image Registration: The elastix Toolbox.

The lead developers of elastix are Stefan Klein and Marius Staring.

This software was initially developed at the Image Sciences Institute, under supervision of Josien P.W. Pluim. Today, many have contributed to elastix.

If you use this software anywhere we would appreciate if you cite the following articles:

• K. Ntatsis, N. Dekker, V. Valk, T. Birdsong, D. Zukić, S. Klein, M Staring, M McCormick, "itk-elastix: Medical image registration in Python", Proceedings of the 22nd Python in Science Conference, pp. 101 - 105, 2023, https://doi.org/10.25080/gerudo-f2bc6f59-00d.

• D.P. Shamonin, E.E. Bron, B.P.F. Lelieveldt, M. Smits, S. Klein and M. Staring, "Fast Parallel Image Registration on CPU and GPU for Diagnostic Classification of Alzheimer's Disease", Frontiers in Neuroinformatics, vol. 7, no. 50, pp. 1-15, January 2014.

This ITK module is based on SimpleElastix, created by Kasper Marstal. For more information, see:

• Kasper Marstal, Floris Berendsen, Marius Staring and Stefajkn Klein, "SimpleElastix: A user-friendly, multi-lingual library for medical image registration", International Workshop on Biomedical Image Registration (WBIR), Las Vegas, Nevada, USA, 2016

itkelastix's Issues

Acknowledgement link

Hi, I just noticed that you acknowledge SimpleElastix in the README. Thanks for this, much appreciated!!

Is it possible to have the README link to my Github profile github.com/kaspermarstal the same way that you link to Stefan's github.com/stefanklein and Marius' github.com/mstaring?

No option to load parameter file from memory

It is not possible to load a transform parameter file directly from memory, the file has to be save to a .txt file first before it can be read and given to elastix. See example notebook 11.

Dependency propblem:

Hello.

I seem to have a problem running ITKElastix both with my own registration script and with itkElastixRegistrationMethodTest.py from master.

$ python3 itkElastixRegistrationMethodTest.py global_0001.png global_0002.png  test.png
Traceback (most recent call last):
  File "itkElastixRegistrationMethodTest.py", line 61, in <module>
    registered, transform = itk.elastix_registration_method(fixed, moving)
  File "/home/lobacheo/.local/lib/python3.7/site-packages/itkLazy.py", line 52, in __getattribute__
    itkBase.LoadModule(module, namespace)
  File "/home/lobacheo/.local/lib/python3.7/site-packages/itkBase.py", line 82, in LoadModule
    LoadModule(dep, namespace)
  File "/home/lobacheo/.local/lib/python3.7/site-packages/itkBase.py", line 78, in LoadModule
    data = module_data[name]
KeyError: 'ITKSmoothing'

I have Python bindings to ITK from the latest pip.

$ pip3 install --upgrade --user pip
Requirement already up-to-date: pip in /home/lobacheo/.local/lib/python3.7/site-packages (20.0.2)
$ pip3 list | grep itk | awk '{print $1}' | xargs pip3 install --user --upgrade

$ pip3 list | grep itk
itk                          5.1rc1.post1       
itk-anisotropicdiffusionlbr  1.0.2              
itk-core                     5.1rc1.post1       
itk-elastix                  0.5.1              
itk-elastix-opencl           0.5.1              
itk-filtering                5.1rc1.post1       
itk-io                       5.1rc1.post1       
itk-meshtopolydata           0.5.1              
itk-numerics                 5.1rc1.post1       
itk-registration             5.1rc1.post1       
itk-segmentation             5.1rc1.post1       
itkwidgets                   0.25.2             

$ uname -a
Linux iana1 4.19.0-6-amd64 #1 SMP Debian 4.19.67-2+deb10u1 (2019-09-20) x86_64 GNU/Linux
$ cat /etc/debian_version 
10.2
$ python3 --version
Python 3.7.3

python -v dump is attached: py-v_dump.txt.

I can provide the images at the request, but it lloks to be as the poblem is not the input. It is also not the file format. It is around the parameter file object, as far as I can see.

Best,
Oleg

Turn KWStyle back on

It was disabled by #7.

Deformation Field calculation not possible with procedural method of Transformix

module 'itk' has no attribute 'imread'

hi, I just installed the itk_elastix_opencl as the way you said :
pip install itk-elastix-opencl
and I got Attribute Error:module 'itk' has no attribute 'imread'
when i do fixed_image = itk.imread('path/to/fixed_image.mha'),why i got this issue?

Problem when opening a big TIF file

Hello,

I am doing mouse brains alignement to the Allen Brain Atlas. Samples are optically cleared and images are obtained on a light sheet microscope.

It works nicely on small files (25x25x25 um voxel : 100 Mo in size).

However, when I am trying to apply the transformation to a high resolution data set of 70 Go, I have a problem: I cannot open the file at all:

• when using the view module from ITKWIDGETS, the big image looks empty.
• when trying to change the voxel spacing, it works on small size file but not on the big one.

If I crop the big one, it is opened correctly... That makes me think that the data type is good, but that it is some sort of memory problem.

I am working on Windows 10, with 512 Go Ram.

Thank you very much for your help!!

Best,

Jérémie

Initial transform fails for transformix

A similar segfault (#68) occurs when loading an initial transform in transformix.

registration with mask does not work

I tried to use a binary mask as the fixed mask to guide the registration, but the log said fMask was unspecified. What Should I do to make the mask work?

How to get a metric for the quality of the registration?

Is there a way to output the value of the metric that was used to achieve a registration?
In the notebooks are some examples on how to set which metric should be used, but I didn't find anything about getting the information after the registration.

After pip install: ImportError: DLL load failed while importing _ITKCommonPython

Could you please help us out? We (@ViktorvdValk @ntatsisk, me) tried to pip install itk-elastix on Windows. It did work well with itk-elastix 0.12.0, but no longer with the latest version, 0.14.1 or 0.14.0. It might be a Windows issue...? See output from the Python terminal, below.

Python 3.9.12 (main, Apr  4 2022, 05:22:27) [MSC v.1916 64 bit (AMD64)] :: Anaconda, Inc. on win32
Type "help", "copyright", "credits" or "license" for more information.
>>> import itk
>>> itk.elastix_registration_method
Traceback (most recent call last):
  File "<stdin>", line 1, in <module>
  File "C:\ProgramData\Anaconda3\envs\ITKElastixPy39_14\lib\site-packages\itk\support\lazy.py", line 138, in __getattribute__
    base.itk_load_swig_module(module, namespace)
  File "C:\ProgramData\Anaconda3\envs\ITKElastixPy39_14\lib\site-packages\itk\support\base.py", line 96, in itk_load_swig_module
    itk_load_swig_module(dep, namespace)
  File "C:\ProgramData\Anaconda3\envs\ITKElastixPy39_14\lib\site-packages\itk\support\base.py", line 96, in itk_load_swig_module
    itk_load_swig_module(dep, namespace)
  File "C:\ProgramData\Anaconda3\envs\ITKElastixPy39_14\lib\site-packages\itk\support\base.py", line 96, in itk_load_swig_module
    itk_load_swig_module(dep, namespace)
  File "C:\ProgramData\Anaconda3\envs\ITKElastixPy39_14\lib\site-packages\itk\support\base.py", line 132, in itk_load_swig_module
    l_module = loader.load(swig_module_name)
  File "C:\ProgramData\Anaconda3\envs\ITKElastixPy39_14\lib\site-packages\itk\support\base.py", line 291, in load
    l_spec.loader.exec_module(l_module)  # pytype: disable=attribute-error
  File "<frozen importlib._bootstrap_external>", line 850, in exec_module
  File "<frozen importlib._bootstrap>", line 228, in _call_with_frames_removed
  File "C:\ProgramData\Anaconda3\envs\ITKElastixPy39_14\lib\site-packages\itk\support\..\ITKCommonPython.py", line 13, in <module>
    from . import _ITKCommonPython
ImportError: DLL load failed while importing _ITKCommonPython: The specified module could not be found.
>>> itk.elastix_registration_method

ITKElastix segfaults if an InitialTransformParameterFile itself contains `InitialTransformParametersFileName`

Same issue as in SimpleITK: SuperElastix/SimpleElastix#121

Python wrappings don't work in RelWithDebInfo mode

import itk
parameter_object = itk.ParameterObject.New()

produces:

Traceback (most recent call last):
  File "C:\Dev\ITK-py\Wrapping\Generators\Python\itk\support\lazy.py", line 137, in __getattribute__
    base.itk_load_swig_module(module, namespace)
  File "C:\Dev\ITK-py\Wrapping\Generators\Python\itk\support\base.py", line 110, in itk_load_swig_module
    l_module = loader.load(swig_module_name)
  File "C:\Dev\ITK-py\Wrapping\Generators\Python\itk\support\base.py", line 259, in load
    l_spec.loader.exec_module(l_module)  # pytype: disable=attribute-error
  File "<frozen importlib._bootstrap_external>", line 728, in exec_module
  File "<frozen importlib._bootstrap>", line 219, in _call_with_frames_removed
  File "C:\Dev\ITK-py\Wrapping\Generators\Python\itk\support\..\ElastixPython.py", line 13, in <module>
    from . import _ElastixPython
ImportError: cannot import name '_
' from 'itk' (C:\Dev\ITK-py\Wrapping\Generators\Python\itk\__init__.py)
python-BaseException

Process finished with exit code -1

_ElastixPython.pyd exists on disk.

If Release mode is used, everything works as expected. If this is ever fixed, a humane warning in #142 should be removed.

Preprocessing example notebooks

A few notebooks that demonstrate how preprocessing steps can be used to improve registration. Apply itk denoising filters, etc.

Work with nii files

Hi!
I would like to know if it's possible to load .nii files directly.

For example: fixed_image = itk.imread('image.nii')

Thank you very much in advance!

Transforming the Transform Parameters to a homogeneous transformation matrix

I would like a way of outputting the transformation parameters as a single transformation matrix.
As an intermediate solution, a clear explanation on how to use the information given in the Transform Parameters to achieve the same result as the output of the registration.

Currently I am using the Center of Rotation (z, y, x) and the 6 Transform values (z_rot, y_rot, x_rot, z_trans, y_trans, x_trans)
to create this matrix myself, but there are slight differences in the result of the image transformed with my matrix and the direct result of the registration. What is the order of the rotations? Is the description of the 6 Transform values even correct?

---

Edit:

After a week of trying and just 4 hours after typing this issue I finally figured out what elastix was doing:

CenterOfRotation (z, y, x)
TransformValues (- z_rot, - y_rot, - x_rot, z_trans, y_trans, x_trans)

To get the Homogeneous Transformation Matrix you need to shift the volume to the center (0,0,0) using the
CenterOfRotation.

CoR_mat = np.array([[1,0,0, -x ],
         [0,1,0, -y ],
         [0,0,1, -z ], 
         [0,0,0,1]])

(np being numpy)
Then rotate the volume and shift it back to its original position. Remember that you need the value from the parameter file multiplied with -1.

x_rot = np.array([[1, 0, 0, 0],
                      [0, math.cos(x_rot), -math.sin(x_rot), 0],
                      [0, math.sin(x_rot), math.cos(x_rot), 0], 
                      [0, 0, 0, 1]])
    
y_rot = np.array([[math.cos(y_rot), 0, math.sin(y_rot), 0],
                      [0, 1, 0, 0],
                      [-math.sin(y_rot), 0, math.cos(y_rot), 0],
                      [0, 0, 0, 1]])
    
z_rot = np.array([[math.cos(z_rot), -math.sin(z_rot), 0, 0],
                      [math.sin(z_rot), math.cos(z_rot), 0, 0],
                      [0, 0, 1, 0],
                      [0, 0, 0, 1]])
    
rotation_mat =  z_rot @ y_rot @ x_rot

neg_CoR_mat = np.array([[1,0,0, x * scale],
                                 [0,1,0, y * scale],
                                 [0,0,1, z * scale], 
                                 [0,0,0,1]])

At last just apply the translation from the parameter file and multiply it all together:

translation_m = np.array([[1,0,0, x_trans],
                                 [0,1,0, y_trans],
                                 [0,0,1, z_trans], 
                                 [0,0,0,1]])

HMT = translation_m @ neg_CoR_mat  @ rotation_mat @ CoR_mat 

BUT!: This is a Pull transformation, things like scipy take this matrix to apply it to a volume, mathematically you would need to
inverse it to apply it to a vector to get the transformed vector you wanted.

This gave me a headache...

Connect to project MONAI

Investigate the connection with project MONAI

Procedural interface Transformix not working in 3D

The procedural interface of transformix is not working in 3D, see example notebook 9.
For transformix to work in 3D in the object oriented interface an explicit call to the 3D transformix class is necessary

Error when using GradientDifference or NormalizedGradientCorrelation Metric

I am trying to use the GradientDifference and NormalizedGradientCorrelation metrics but could not quite figure out how to prepare the input. I am following the ITK_Example07_MultimetricMultiImageRegistration example and use a parameters_Bspline_Multimetric.txt where I have only changed the metric.
When I am using the NormalizedGradientCorrelation metric alone with 2D float32 images in the same way as in the example, I get the error

File: _deps/elx-src/Components/Metrics/NormalizedGradientCorrelation/elxNormalizedGradientCorrelationMetric.hxx
Line: 55
Description: ITK ERROR: NormalizedGradientCorrelationMetric(0x555fb32f4fd0): FixedImage must be 3D

So I figured that it probably needs the two gradient channels as input. When I run the following code

fixed_grad = itk.GradientImageFilter(fixed)
moving_grad = itk.GradientImageFilter(moving)

elastix_object = itk.ElastixRegistrationMethod.New(fixed_grad, moving_grad,
                                                   log_to_console=True,
                                                   number_of_threads=12)

it throws itk.support.extras.TemplateTypeError: itk.ElastixRegistrationMethod is not wrapped for input type None.
The images from the gradient filter have vector pixel type so, combining the two error messages, the next thing I have tried is to convert the gradient filter to 3D float images. But then I end up with the following error:

File: _deps/elx-src/Components/Metrics/NormalizedGradientCorrelation/elxNormalizedGradientCorrelationMetric.hxx
Line: 61
Description: ITK ERROR: NormalizedGradientCorrelationMetric(0x5649875b3f90): Metric can only be used for 2D-3D registration. FixedImageSize[2] must be 1

Do you have any idea what I am missing or could point me to resources that I could check out? Ideally, could you consider adding an example for the gradient-based metrics?
Thanks and best wishes,
Philipp

[Progress Bars] Including them to show status of iterations for each parameter

Issue

• Working with large 3D volumes leads to a time consuming registration process.

Proposed Solution

• Having a progress bar (like tqdm) would be a useful indicator of how long the registration shall be running for
• The MaximumNumberOfIterations flag of the optimizer in the parameter files could be used to indicate the current status of each parameter.

Example

• The 3D CT, monomodal (for head and neck) preset contains two parameter maps
  • Rigid [with MaximumNumberOfIterations=250]
  • BSpline [with MaximumNumberOfIterations=1000]

Risks

• I am not sure if the iterations for the optimization process actually reaches the MaximumNumberOfIterations every single time.
• Though I have seen the following message in the logs

  • Time spent in resolution 1 (ITK initialization and iterating): 77.464 s.
    Stopping condition: Maximum number of iterations has been reached.

Reduce wheel size below 100MB

With v0.14.0,itk_elastix-0.14.0-cp310-cp310-win_amd64.whl is 119MB and therefore over the PyPI limit. I yanked the release as only 2 wheels were uploaded.

The result is the same with and without pointset

Hello ITKElastix team,

Thank you very much for your useful library. I want to try the 3D mask registration based on 4 points. My code is as follows. But when I try with and without a pointset, the output is the same.

for m in a:
    for n in b:
        fixed_image = itk.imread(f'mask_prepration_registration_corect/image{m}.nii.gz', itk.F)
        moving_image = itk.imread(f'mask_prepration_registration_corect/image{n}.nii.gz', itk.F)

        moving_mask = itk.imread(f'mask_prepration_registration_corect/mask{n}.nii.gz', itk.F)
        parameter_object = itk.ParameterObject.New()
        parameter_map_rigid = parameter_object.GetDefaultParameterMap('rigid')
        parameter_object.AddParameterMap(parameter_map_rigid)
        parameter_map_affine = parameter_object.GetDefaultParameterMap('affine')
        parameter_object.AddParameterMap(parameter_map_affine)
        parameter_map_bspline = parameter_object.GetDefaultParameterMap('bspline')
        parameter_object.AddParameterMap(parameter_map_bspline)

        result_image, result_transform_parameters = itk.elastix_registration_method(
            fixed_image, moving_image,
            fixed_point_set_file_name=f'data/landmarks/targets{m}.txt',
            moving_point_set_file_name=f'data/landmarks/targets{n}.txt',
            parameter_object=parameter_object,
            log_to_console=False)

        result_transform_parameters.SetParameter('FinalBSplineInterpolationOrder', '0')
        result_moving_mask = itk.transformix_filter(moving_mask, result_transform_parameters)

        itk.imwrite(result_moving_mask, f'masks_registration_result/{m}/reg_{m}_{n}.nii.gz')

Where is my solution wrong?

ITK's function image_view_from_array produces 4D image with array of >10*9 inputs.

This is the cause for the memory error in PR #99

The following script returns an image with DimSize = 10 10 3 2074148458.
While it should return an image with DimSize = 10 10 10 3.

def image_generator_3D(x1, x2, y1, y2, z1, z2):
    image = np.zeros([10, 10, 10], np.float32)
    image[z1:z2, y1:y2, x1:x2] = 1
    return image

vector_of_images = np.zeros([3, 10, 10, 10], np.float32)
i = 0
for x in range(0, 5, 2):
    image = image_generator_3D(x, x+5, x, x+5, x, x+5)
    vector_of_images[i] = image
    i += 1

vector_itk = itk.image_view_from_array(vector_of_images.T)

Or did I make a mistake here?
Note: Making a 3D image from 3 2D images does work this way.

Zarr spatial transformation serialization notebook.

Depends on #79

Problem using ITKElastix in another ITK remote module

I am trying to use ITKElastix from HASI. To accomplish this, I am adding ITKElastix as a remote module to ITK. ITK builds fine, but HASI runs into link errors (cannot find elastix_lib.lib). In an attempt to fix this, I tried adding CMake export code to ITKElastix, but this causes circular dependency of transformix_lib and elastix_lib when trying to generate ITK project from CMake:

Selecting Windows SDK version 10.0.18362.0 to target Windows 10.0.18363.
Running module dependency checks
Enabled BoneEnhancement, requested by Module_BoneEnhancement.
Enabled Cuberille, requested by Module_Cuberille.
Enabled Elastix, requested by Module_Elastix.
Enabled ITKAnisotropicSmoothing, requested by ITK_BUILD_DEFAULT_MODULES.
Enabled ITKAntiAlias, requested by ITK_BUILD_DEFAULT_MODULES.
...
Enabled ITKVtkGlue, requested by Module_ITKVtkGlue.
Enabled ITKWatersheds, requested by ITK_BUILD_DEFAULT_MODULES.
Enabled ITKZLIB, requested by ITK_BUILD_DEFAULT_MODULES.
Elastix_LIBRARIES: elastix_lib;transformix_lib
elx_BINARY_DIR: C:/Dev/ITK-git-2015/_deps/elx-build
Checking for InitOnceExecuteOnce:
Performing Test InitOnceExecuteOnce - Success
Warnings Configuration:
Filter ZLIB is ON
Could NOT find Perl (missing: PERL_EXECUTABLE) 
Configuring done
CMake Error at C:/Dev/ITK-git-2015/_deps/elx-src/Core/CMakeLists.txt:184 (add_library):
  Target "transformix_lib" links to itself.

CMake Error at C:/Dev/ITK-git-2015/_deps/elx-src/Core/CMakeLists.txt:184 (add_library):
  Target "transformix_lib" links to itself.

CMake Error at C:/Dev/ITK-git-2015/_deps/elx-src/Core/CMakeLists.txt:154 (add_library):
  Target "elastix_lib" links to itself.

CMake Error at C:/Dev/ITK-git-2015/_deps/elx-src/Core/CMakeLists.txt:154 (add_library):
  Target "elastix_lib" links to itself.

CMake Error at C:/Dev/ITK-git-2015/_deps/elx-src/Core/CMakeLists.txt:154 (add_library):
  Target "elastix_lib" links to itself.

CMake Error at C:/Dev/ITK-git-2015/_deps/elx-src/Core/CMakeLists.txt:184 (add_library):
  Target "transformix_lib" links to itself.

Generating done

enumeration support for the wrapping of ParameterObject

It would be helpful to have a way of enumerating a list/dict of a ParameterObject like
for i, parameterMap in enumerate(parameterObject):.
Currently, this needs to be done as:

    for i in range(parameterObject.GetNumberOfParameterMaps()):
       parameterMap = parameterObject.GetParameterMap(i)

SetFixedImage doesn't work for 3D images

I was having a look at the notebooks (good way to start learning on how to use this package) and trying to run some of them.

I tried the 3D registration example with some data I have and it returns the following error at itk.ElastixRegistrationMethod.New().SetFixedImage: Expecting argument of type itkImageF2 or itkImageSourceIF2.
I tried this after loading a 3D nifti image using imread - itk.itkImagePython.itkImageF3.

Does this mean it's not possible to run registration in 3D images using the elastix filter?

[Elastix Parameters --> ITK Transforms] Provide function to convert Elastix's TransformParameter.txt to ITK transforms

Context

• 3D Visualization software like 3D Slicer provide a larger number of visualization modules for different purposes.
• One of the modules, the Transform module allows loading of .tfm files. These files can be created by using for e.g. sitk.WriteTransform(sitk.BSplineTransform(3), 'BSpline.tfm'). This Transform module is quite handy as it allows for easy visualization of the deformed grids of BSpline.

Is there a roadmap to include these kind of Elastix -> ITK Transform helper functions in ITKElastix?

• Example 1: Convert Affine + BSpline into a Composite Transform

Example notebook for reading / writing point sets from standard file formats

Depends on #75. Use itk.meshread, itk.meshwrite to read write itk.PointSet's for registration.

Postprocessing example notebook

A few post-processing example notebooks to demonstrate how a resulting registration can be used for analysis. Quantification in regions from a labeled atlas. Multiple features for segmentation. Use of the transformation to align related data structures. Etc.

elastix resampled image from transformix and ItkTransform resampled image produce different results for same transformations

Hi all, specifically @ViktorvdValk @N-Dekker

I've been monitoring a bit the GitHub issues and upcoming PRs for ITKElastix and elastix itself and thought it was timely to bring this issue up considering the direction things are going for interoperability between elastix and ITK transforms. I wasn't sure exactly where to file this so I will put it here and have tagged who I think would be interested.

Related issues:
#79
elastix:
SuperElastix/elastix#387

Background

I am registering some 2D microscopy images and having great results with elastix, however, I'd like to use a bit more of ITK for managing point set transformations and composing long and complex sequences of transforms. However, in my initial tests I've discovered that the resampling in ITK and elastix can produce different images with spatial shifts that could affect point set and image transformations.

Reproducible example & data

Here is a full example with data (data link: https://drive.google.com/file/d/1-bUfi69YmBVGFPzX4XRTPPCkIja9Qur8/view?usp=sharing, sorry it is microscopy so rather large 1.1 GB, includes both the moving image, and ITK and elastix transformed images) using SimpleITK but I think it would apply to ITK or elastix implementation one uses. Apologies for the length, I wanted everything to be as clear as possible.

import numpy as np
import SimpleITK as sitk
import napari


def elx_to_itk_euler2d(transform_data):
    t_params = np.asarray(transform_data["TransformParameters"]).astype(np.float)
    t_center = np.asarray(transform_data["CenterOfRotationPoint"]).astype(np.float)

    theta = t_params[0]
    c, s = np.cos(theta), np.sin(theta)
    rotmat = np.array(((c, -s), (s, c)))

    affine = sitk.AffineTransform(2)
    affine.SetTranslation(t_params[-2:].tolist())
    affine.SetMatrix(rotmat.flatten().tolist())
    affine.SetCenter(t_center.tolist())
    return affine


def elx_to_itk_affine2d(transform_data):
    t_params = np.asarray(transform_data["TransformParameters"]).astype(np.float)
    t_center = np.asarray(transform_data["CenterOfRotationPoint"]).astype(np.float)

    affine = sitk.AffineTransform(2)
    affine.SetTranslation(t_params[-2:].tolist())
    affine.SetMatrix(t_params[:4].flatten().tolist())
    affine.SetCenter(t_center.tolist())
    return affine


# read elastix parameter files
euler_pmap = sitk.ReadParameterFile("./example-data/TransformParameters.0.txt")
affine_pmap = sitk.ReadParameterFile("./example-data/TransformParameters.1.txt")

# read moving image and set spacing
moving_im = sitk.ReadImage("./example-data/moving_image.tif")
moving_im.SetSpacing((0.65, 0.65))

# set up fixed image from last transform in sequence
size_out = np.asarray(affine_pmap["Size"]).astype(np.int)
spacing_out = np.asarray(affine_pmap["Spacing"]).astype(np.float)
# direction and origin are default
fixed_im = sitk.Image(size_out.tolist(), sitk.sitkUInt16)
fixed_im.SetSpacing(spacing_out.tolist())

# convert elastix transform parameter maps to ITK representations
itk_euler = elx_to_itk_euler2d(euler_pmap)
itk_affine = elx_to_itk_affine2d(affine_pmap)

# compose sequential transforms
compos_tform = sitk.Transform(2, sitk.sitkComposite)
compos_tform.AddTransform(itk_euler)
compos_tform.AddTransform(itk_affine)

# apply transform
itk_resampled_im = sitk.Resample(
    moving_im, fixed_im, compos_tform, sitk.sitkNearestNeighbor, 0
)

# load elastix resampled image
elx_resampled_im = sitk.ReadImage("./example-data/elx_resampled_image.tif")

# optionally read the ITK resampled image from disk rather than computing it
# itk_resampled_im = sitk.ReadImage("./example-data/itk_resampled_image.tif")

# visualize with napari
with napari.gui_qt():
    viewer = napari.Viewer()
    viewer.add_image(
        sitk.GetArrayFromImage(itk_resampled_im),
        name="itk resampled",
        contrast_limits=[0, 5000],
        colormap="red",
        blending="additive"

    )
    viewer.add_image(
        sitk.GetArrayFromImage(elx_resampled_im),
        name="elx resampled",
        contrast_limits=[0, 5000],
        colormap="green",
        blending="additive"
    )

Results

Here are the results with red-green (elastix, itk) additive color blending of the different resamplings:
low-res is ok

high-res shows some issues, a small positive shift of itk image

I've also looking at corresponding areas between the original moving image and the different resamplings and elastix's resampled image has higher fidelity to the original pixel intensities and shapes and is a closer alignment with the fixed image.

Not sure what is the cause of this. Other things I've explored producing the same result:

• downsampling of the image
• precision of numbers in TransformParameters
• sequence transformation rather than ITK composition transform, i.e., transform with euler -> interpolate -> affine -> interpolate

Sorry for the long issue.. hope it is clear.

itk.imread() gives non-fatal dynamic cast errors in the lldb python debugger.

The lldb warnings/errors thrown by itk.imread:

Process 6134 launched: '/Users/viktorvandervalk/anaconda3/envs/ITKDebugNielsFinal/bin/python' (x86_64)
2021-02-18 17:09:41.070552+0100 python[6134:292558] dynamic_cast error 1: Both of the following type_info's should have public visibility. At least one of them is hidden. N3itk10DataObjectE, N3itk5ImageIfLj2EEE.
2021-02-18 17:09:41.070616+0100 python[6134:292558] dynamic_cast error 2: One or more of the following type_info's has hidden visibility or is defined in more than one translation unit. They should all have public visibility. N3itk10DataObjectE, N3itk5ImageIfLj2EEE, N3itk9ImageBaseILj2EEE.
2021-02-18 17:09:41.078145+0100 python[6134:292558] dynamic_cast error 1: Both of the following type_info's should have public visibility. At least one of them is hidden. N3itk10DataObjectE, N3itk5ImageIfLj2EEE.
2021-02-18 17:09:41.078190+0100 python[6134:292558] dynamic_cast error 2: One or more of the following type_info's has hidden visibility or is defined in more than one translation unit. They should all have public visibility. N3itk10DataObjectE, N3itk5ImageIfLj2EEE, N3itk9ImageBaseILj2EEE.

Unable to run example notebooks

Attempting to run the example notebooks causes the kernel to die when executing the registration step.

After running into issues attempting my own registrations, I gave the examples a try without any luck. Not sure that the issues are necessarily related, but figured I'd start here with the examples.

I'm using the python package for itk: itk 5.1.0.post3 (<-- perhaps this is the issue???) and the most recent ITKElastix (0.6.2).

Perhaps I'm missing something simple, but any help, tips, or guidance is appreciated! Thanks!

object-oriented interface similar to SimpleElastix

Is there the possibility to make use of ITKElastix with an object-oriented interface similar to SimpleElastix (as described here: https://simpleelastix.readthedocs.io/HelloWorld.html)?
For example to transfer an SE script to ITKElastix, I am looking for the correct way to translate e.g. these lines:
https://github.com/romangrothausmann/elastix_scripts/blob/ba86940716bacfed96641a48a39053014519cb54/recRegStack.py#L153-L157

"HasKey" or similar method for the wrapped ParameterObject

For pythonic way of checking if a key is in a list/dict like
if 'InitialTransformParametersFileName' in parameterObject:
it would be helpful to have a "HasKey" or similar method for the wrapped ParameterObject.

Installation not working on M1 Macs

I have tried installing using pip install itk-elastix on my M1 Macbook, but I get the error:

ERROR: Could not find a version that satisfies the requirement itk-elastix (from versions: none) ERROR: No matching distribution found for itk-elastix

I have tried with Python 3.7, 3.8 and 3.9 with no success. I have also tried compiling from source with no success.
My experience with compiling these sort of packages is minimal? Could you provide some help?

Enable all components of elastix

Add ELX_USE_ALL_COMPONENTS flag to elastix build.

Is the log file really closed?

Hi,
I got the impression that the log file generated from running (Windows10, jupyter notebook)

result_image, result_transform_parameters = itk.elastix_registration_method(
fixed, moving,
...
log_to_console=False,
log_to_file=True,
...)

is not really closed after finishing registration. Is this a known issue or did I forget to add some code?

Compile with different pixel types?

Is it possible to compile the wrapped elastix to support pixel types other than float32? I do a lot of registration with very large microscopy planes and having just an unsigned short pixel type would save a lot of memory. I'm happy to try things out on my end if anyone could point me to the relevant code..

Thanks!

Get Affine Transformation Matrix from TransformParameters

I'm using ITK to perform a affine coregistration on 3D MRI scans.
Now similar to what is suggested in https://github.com/InsightSoftwareConsortium/ITKElastix/blob/master/examples/ITK_Example08_SimpleTransformix.ipynb I would like to perform the affine transformation on another image but we are using another library for that.

So I would need a full affine matrix for the transformation that describes the affine coregistration.

I noticed that the registration returns the transform parameters:

result_image, result_transform_parameters = itk.elastix_registration_method( ...
parameter_map = transform_parameters.GetParameterMap(0)
transform_parameters = np.array(parameter_map['TransformParameters'], dtype=float)

From other GitHub issues I gathered that the last 3 values of the 12 values in transform_parameters is for the translation and the first 9 are for the rotation. Not sure if that is correct.

Based on this I started playing around and got close to the correct affine transformation matrix (the resulting images are similar but not quite right).

rotation = transform_parameters[:9].reshape(3, 3)
translation = transform_parameters[-3:][..., np.newaxis]
translation = -np.flip(translation)
reg_affine: np.ndarray = np.append(rotation, translation, axis=1) 
reg_affine = np.append(reg_affine, [[0,0 ,0,1]], axis=0) 

The translation might be right but the rotation is slightly off. What am I missing / what is wrong? Also why do I need to flip and negate the translation? Do I need to do something similar with the rotation matrix?

Would really appreciate some insight. The ITK methods are a blackbox to me.

---

Not sure what example data I can share. Here some example values for the above variables:

transform_parameters = [ 9.98573286e-01, -3.92533565e-03, -7.45617495e-03,  6.11828178e-04,
                        9.98081930e-01,  7.21948277e-03,  2.74329272e-03, -1.22292851e-02,
                        1.00763651e+00, -1.52612188e-01,  8.89476641e-01,  2.55258072e+00]

with the resulting affine matrix:

reg_affine = [[ 1.00763651e+00, -1.22292851e-02,  2.74329272e-03, -2.55258072e+00],
             [ 7.21948277e-03,  9.98081930e-01,  6.11828178e-04, -8.89476641e-01],
             [-7.45617495e-03, -3.92533565e-03,  9.98573286e-01, 1.52612188e-01],
             [ 0.00000000e+00,  0.00000000e+00,  0.00000000e+00, 1.00000000e+00]

Getting the result image in the object oriented approach fails

GetOutput() is shadowed from ImageToImageFilter.

Transfer PointSet from Moving Image to Fixed Image

Hi @thewtex @blezek @bradking @cpatrick @arnaudgelas @dzenanz @N-Dekker

In the example, ITK_Example09_PointSetAndMaskTransformation.ipynb, it is mentioned ...

When transforming a point set, this happens the other way around: from the fixed image domain to the moving image domain. This is useful for example when you want to know where a point maps to. Point sets therefor need not be transformed with the backwards mapping protocol, but can instead be transformed with regular forward transformation (fixed -> moving). Transforming point sets can be used to get the regions of interest (ROI) in the moving image, based on ROI of the fixed image.

Can we transfer points in moving image to the fixed image domain ?
How to do the backwards mapping in this case ?
What change should be made here ?

# Procedural interface of transformix filter
result_point_set = itk.transformix_pointset(
    moving_image, result_transform_parameters,
    fixed_point_set_file_name='data/CT_3D_lung_fixed_point_set_corrected.txt',
    output_directory = './exampleoutput')

Thanks

Need Help: ITK ERROR: ElastixRegistrationMethod(00000187F232D140)

Dear Developers:
I am a Ph.D. candidate from Fudan University, China.
Thanks for your kind providing us with the good package: ITK-Elastix, which is very helpful for medical registration. However, I encountered some questions after the installation.
My programing environment is Anaconda+Python 3.7.12. After successful installation via the command:

"pip install itk-elastix", I met some errors when I tried to register the images via the following codes:

# Import Packages
import numpy as np
import matplotlib.pyplot as plt
import itk
b0_path = r"C:\Users\yunfe\Desktop\Image Registration\b0.mha"
b800_path = r"C:\Users\yunfe\Desktop\Image Registration\b800.mha"

# Load images with itk floats (itk.F). Necessary for elastix
fixed_image = itk.imread(b0_path, itk.F)
moving_image = itk.imread(b800_path, itk.F)
parameter_object = itk.ParameterObject.New()
default_rigid_parameter_map = parameter_object.GetDefaultParameterMap('rigid')
parameter_object.AddParameterMap(default_rigid_parameter_map)
# Call registration function
result_image, result_transform_parameters = itk.elastix_registration_method(
    fixed_image, moving_image,
    parameter_object=parameter_object,
    log_to_console=False)_

The errors are as follows:

When I tried with other codes:

# Load Elastix Image Filter Object
elastix_object = itk.ElastixRegistrationMethod.New(fixed_image, moving_image)
# elastix_object.SetFixedImage(fixed_image)
# elastix_object.SetMovingImage(moving_image)
elastix_object.SetParameterObject(parameter_object)

# Set additional options
elastix_object.SetLogToConsole(False)

# Update filter object (required)
elastix_object.UpdateLargestPossibleRegion()

# Results of Registration
result_image = elastix_object.GetOutput()
result_transform_parameters = elastix_object.GetTransformParameterObject()

Similar errors also appeared:

4D groupwise registration not supported in currrent ITKElastix

ITKElastix should be rebuild with support for 4D images to use the 4D groupwise image registration of Elastix

transformix_jacobian tries to use files

Also seems to be not working in itkTransformFilterTest.py, commented in #109

25:   File "itkTransformixFilterTest.py", line 51, in <module>
25:     jacobian = itk.transformix_jacobian(moving, parameters, log_to_console=True)
25:   File "/home/matt/bin/ITK-Wrap-Release/Wrapping/Generators/Python/itk/support/../ElastixPython.py", line 99, in transformix_jacobian
25:     det_spatial_jacobian = itk.imread('spatialJacobian.nii', itk.F)
25:   File "/home/matt/bin/ITK-Wrap-Release/Wrapping/Generators/Python/itk/support/itkExtras.py", line 692, in imread
25:     raise RuntimeError("No ImageIO is registered to handle the given file.")
25: RuntimeError: No ImageIO is registered to handle the given file.
25: itkTestDriver: Process exited with return value: 1

Example notebook demonstrating integration with (py)clesperanto.

Error on grid spacing

Hello,
We are trying to register 2 3D images with ITK-Elastix with the following code

# Load images in itk workspace
fixed_image = itk.imread('0N_ROI.tif')
moving_image = itk.imread('1450N_ROI.tif')

# Cast images to an elastix compatible one
fixed_image = fixed_image.astype(itk.F)
moving_image = moving_image.astype(itk.F)

#Load and set parameter file
parameter_object = itk.ParameterObject.New()
parameter_map_rigid = parameter_object.GetDefaultParameterMap('affine',1)
parameter_object.AddParameterMap(parameter_map_rigid)
parameter_map_Bspline = parameter_object.GetDefaultParameterMap('bspline',5)
parameter_object.AddParameterMap(parameter_map_Bspline)
parameter_object.SetParameter("FinalGridSpacingInPhysicalUnits", "16")
parameter_object.SetParameter("Metric1Weight","10000")
parameter_object.SetParameter("MaximumNumberOfIterations","1500")

# check the parameter object
print(parameter_object)


# Call registration function
result_image, result_transform_parameters = itk.elastix_registration_method(
    fixed_image, moving_image,
    parameter_object=parameter_object,
    log_to_console=True,
    output_directory='Res/')I`

I get an strange error :
Description: ITK ERROR: ParameterMapInterface(0x6018de0): ERROR: Casting entry number 0 for the parameter "GridSpacingSchedule" failed!
You tried to cast "3,952542" from std::string to d

Do you have any idea of what is happening ?

How can I get the BSspline transform and apply to other images?

Hi, I want to use BSpline DIR, and get the transform. I can get the transform parameters in examples/2_RegistrationParameters.ipynb, but how can I apply it to other images in the same coordinates.
Like in elastix, transformix cmd can finish this. Is there any pure python code?

No option to output transformed pointset file to memory in transformix.

See example notebook 9.
Pointset file is saved as .txt file but no output in memory is available for further calculations.

2D similarity transform parameter map doesn't keep "CenterOfRotationPoint" parameter

Hi all, thanks for putting together this package, I'm really enjoying the ease of installation.

It seems the similarity transform as wrapped in ITKElastix doesn't return the CenterOfRotationPoint transform parameter (although it is written to file) when doing 2D to 2D registration.

I encountered this issue using SimpleElastix as well so it may be in elastix codebase itself.

import numpy as np
import itk

zero_image = np.zeros([1024,1024])

similarity_parameter_map = parameter_object.GetDefaultParameterMap("rigid")

# set method to similarity & 2D
similarity_parameter_map["Transform"] = ["SimilarityTransform"]
similarity_parameter_map["MovingImageDimension"] = ["2"]
similarity_parameter_map["FixedImageDimension"] = ["2"]

parameter_object = itk.ParameterObject.New()
parameter_object.AddParameterMap(similarity_parameter_map)

sim_result_image, sim_params = itk.elastix_registration_method(
    zero_image.astype(np.float32),
    zero_image.astype(np.float32),
    parameter_object=parameter_object,
    log_to_console=True,
    output_directory="./transforms_similarity",
)

transform_cent_rot = np.array(
    sim_params.GetParameter(0, "CenterOfRotationPoint")
).astype(np.float32)[::-1]

assert len(transform_cent_rot) == 2

print(sim_params)

ParameterObject (0000028751D7CE00)
  RTTI typeinfo:   class elastix::ParameterObject
  Reference Count: 1
  Modified Time: 2007229
  Debug: Off
  Object Name: 
  Observers: 
    none
ParameterMap 0: 
  (CenterOfRotationPoint)
  (CompressResultImage "false")
  (DefaultPixelValue 0)
  (Direction 1 0 0 1)
  (FinalBSplineInterpolationOrder 3)
  (FixedImageDimension 2)
  (FixedInternalImagePixelType "float")
  (HowToCombineTransforms "Compose")
  (Index 0 0)
  (InitialTransformParametersFileName "NoInitialTransform")
  (MovingImageDimension 2)
  (MovingInternalImagePixelType "float")
  (NumberOfParameters 4)
  (Origin 0 0)
  (ResampleInterpolator "FinalBSplineInterpolator")
  (Resampler "DefaultResampler")
  (ResultImageFormat "nii")
  (ResultImagePixelType "float")
  (Size 1024 1024)
  (Spacing 1 1)
  (Transform "SimilarityTransform")
  (TransformParameters 1 0 0 0)
  (UseDirectionCosines "true")

file on disk

(Transform "SimilarityTransform")
(NumberOfParameters 4)
(TransformParameters 1.000000 0.000000 0.000000 0.000000)
(InitialTransformParametersFileName "NoInitialTransform")
(UseBinaryFormatForTransformationParameters "false")
(HowToCombineTransforms "Compose")

// Image specific
(FixedImageDimension 2)
(MovingImageDimension 2)
(FixedInternalImagePixelType "float")
(MovingInternalImagePixelType "float")
(Size 1024 1024)
(Index 0 0)
(Spacing 1.0000000000 1.0000000000)
(Origin 0.0000000000 0.0000000000)
(Direction 1.0000000000 0.0000000000 0.0000000000 1.0000000000)
(UseDirectionCosines "true")

// SimilarityTransform specific
(CenterOfRotationPoint 511.5000000000 511.5000000000)

// ResampleInterpolator specific
(ResampleInterpolator "FinalBSplineInterpolator")
(FinalBSplineInterpolationOrder 3)

// Resampler specific
(Resampler "DefaultResampler")
(DefaultPixelValue 0.000000)
(ResultImageFormat "nii")
(ResultImagePixelType "float")
(CompressResultImage "false")