# Struggling to do simple 3D rigid registration

**URL:** https://discourse.itk.org/t/struggling-to-do-simple-3d-rigid-registration/3217
**Category:** Beginner Questions
**Tags:** simpleitk
**Created:** [June 22, 2020, 5:57pm UTC](https://discourse.itk.org/t/struggling-to-do-simple-3d-rigid-registration/3217 "2020-06-22T17:57:16Z")
**Posts on this page:** 5
**Page:** 1

<div class="post-metadata">

### Author: ![msb](https://discourse.itk.org/user_avatar/discourse.itk.org/msb/32/1588_2.png) [@msb](https://discourse.itk.org/u/msb)
#### Post date: [June 22, 2020, 5:57pm UTC](https://discourse.itk.org/t/struggling-to-do-simple-3d-rigid-registration/3217/1 "2020-06-22T17:57:16Z")

</div>

Hello,

I’m am struggling to perform a simple 3D rigid registration between 2 images. I have started with [one of the tutorials](https://simpleitk.org/doxygen/latest/html/ImageRegistrationMethod1_2ImageRegistrationMethod1_8py-example.html), and tried to adapt it to my situation. Particularly I want to adjust it to do a rigid registration (i.e. including rotations), but I can’t get it to work.

I have taken a PET brain image, manually transformed it by [15, 15, -15] mm (no rotations for now), and then I’m trying to register it to the original image. If I set the SetInitialTransform to be TranslationTransform then it works as expected. However if I set it to be Euler3DTransform then it doesn’t seem to work and I don’t understand the results.

Here is my code:

```
#!/usr/bin/env python
from __future__ import print_function
import SimpleITK as sitk
import sys
import os

def command_iteration(method):
    print("{0:3} = {1:10.5g} : {2}".format(method.GetOptimizerIteration(),
                                           method.GetMetricValue(),
                                           method.GetOptimizerPosition()))

if len(sys.argv) < 4:
    print("Usage:", sys.argv[0], "<fixedImageFilter> <movingImageFile>",
          "<outputTransformFile>")
    sys.exit(1)

fixed = sitk.ReadImage(sys.argv[1], sitk.sitkFloat32)
moving = sitk.ReadImage(sys.argv[2], sitk.sitkFloat32)

# image size: [75, 75, 73] pixels
# pixel size: [4, 4, 2.78] mm
# image extent: [300, 300, 202.94] mm

fixed.SetSpacing((4.0, 4.0, 2.78)); # pixel size in mm
moving.SetSpacing((4.0, 4.0, 2.78));

# fixed.SetOrigin((150, 150, 101.47));
# moving.SetOrigin((150, 150, 101.47));

R = sitk.ImageRegistrationMethod()
R.SetMetricAsMeanSquares()
R.SetOptimizerAsRegularStepGradientDescent(4.0, 0.001, 50, 0.5, 1e-7)

tx = sitk.CenteredTransformInitializer(fixed, moving, sitk.Euler3DTransform(),
    sitk.CenteredTransformInitializerFilter.GEOMETRY)
R.SetInitialTransform(tx)
# R.SetInitialTransform(sitk.TranslationTransform(fixed.GetDimension()))
R.SetInterpolator(sitk.sitkLinear)

R.AddCommand(sitk.sitkIterationEvent, lambda: command_iteration(R))
outTx = R.Execute(fixed, moving)

print("-------")
print(outTx)
print("Optimizer stop condition: {0}"
      .format(R.GetOptimizerStopConditionDescription()))
print(" Iteration: {0}".format(R.GetOptimizerIteration()))
print(" Metric value: {0}".format(R.GetMetricValue()))

sitk.WriteTransform(outTx, sys.argv[3])

```

And the output is:

> ```
> itk::simple::Transform
> Euler3DTransform (0x7fc363631460)
> RTTI typeinfo: itk::Euler3DTransform<double>
> Reference Count: 3
> Modified Time: 2401
> Debug: Off
> Object Name: 
> Observers: 
> none
> Matrix: 
> 0.968088 -0.235134 -0.0867007 
> 0.203394 0.939297 -0.276318 
> 0.14641 0.249866 0.957148 
> Offset: [49.8245, 7.73086, -55.7376]
> Center: [148, 148, 100.08]
> Translation: [1.62477, 1.19532, -1.37756]
> Inverse: 
> 0.968088 0.203394 0.14641 
> -0.235134 0.939297 0.249866 
> -0.0867007 -0.276318 0.957148 
> Singular: 0
> Euler's angles: AngleX=2.88905 AngleY=-3.29338 AngleZ=3.38688
> m_ComputeZYX = 0
> 
> Optimizer stop condition: RegularStepGradientDescentOptimizerv4: Maximum number of iterations (50) exceeded.
> Iteration: 50
> Metric value: 3.75348540451e-05
> 
> ```

Here are my questions:

1. Is my code correct to do a rigid registration with possible translations and/or rotations, with a given pixel size?
2. I expected “Translation” or “Offset” to be [15, 15, -15]. I actually don’t understand the difference between these two.
3. Should I set SetOrigin((150, 150, 101.47)) for the images before the registration? (I would like to do this because in a real case I will usually have an offset center)
4. What is “Center” in the output? I imagine it’s the estimated image center, but if I SetOrigin as above then it just gets added to this “Center”.
5. Are the “Euler’s angles” the rotational part of the registration parameters? If so I expect them to be near 0 or 2\*pi in this case, but they are near pi.

I’m pretty sure that I am just missing something trivial but I could not figure it out after searching the examples, tutorials, documentation, and this forum.

Many thanks for your help.

---

<div class="post-metadata">

### Author: ![zivy](https://discourse.itk.org/user_avatar/discourse.itk.org/zivy/32/1726_2.png) [@zivy](https://discourse.itk.org/u/zivy)
#### Post date: [June 22, 2020, 6:25pm UTC](https://discourse.itk.org/t/struggling-to-do-simple-3d-rigid-registration/3217/2 "2020-06-22T18:25:43Z")

</div>

Hello @msb,

Before you go any farther, please go over the following resources (possibly you missed some of the relevant information):

1. [Fundamental concepts](https://simpleitk.readthedocs.io/en/master/fundamentalConcepts.html).
2. [Registration overview](https://simpleitk.readthedocs.io/en/master/registrationOverview.html).
3. Set of registration notebooks (6\* series) on the [SimpleITK Jupyter notebook repository](https://github.com/InsightSoftwareConsortium/SimpleITK-Notebooks).

With respect to your registration:

1. Why are you manually modifying the image spacing after reading? This is very unusual (this information is always obtained directly from the image when it is read).
2. You have premature termination (stopping condition is “Maximum number of iterations (50) exceeded.”).
3. The nomenclature of “Center”, “Offset”, “Translation” that appear in the printed transformation are described in the fundamental concepts page referenced above and are specific to ITK/SimpleITK.

---

<div class="post-metadata">

### Author: ![msb](https://discourse.itk.org/user_avatar/discourse.itk.org/msb/32/1588_2.png) [@msb](https://discourse.itk.org/u/msb)
#### Post date: [June 22, 2020, 11:19pm UTC](https://discourse.itk.org/t/struggling-to-do-simple-3d-rigid-registration/3217/3 "2020-06-22T23:19:21Z")

</div>

Hello @zivy,

Thank you for your reply. I reviewed those documents and links again in more detail and they did indeed help clarify some things that I missed when I went through them initially.

To answer your first question: I am manually modifying the image parameters because those parameters are not present in the header information. I saved this images from Matlab as nifti images and creating a nifti header from scratch is not trivial in Matlab, so I just edited the relevant values in Python.

I now know how the image origin should be set. I set the image parameters like this:

> ```
> # image size: [75, 75, 73] pixels
> # pixel size: [4, 4, 2.78] mm
> # image extent: [300, 300, 202.94] mm
> fixed.SetSpacing((4.0, 4.0, 2.78)); # pixel size in mm
> moving.SetSpacing((4.0, 4.0, 2.78));
> fixed.SetOrigin((-148, -148, -100.08));
> moving.SetOrigin((-148, -148, -100.08));
> 
> ```

Setting the origin to be half a pixel in from the image extent. The image direction vectors are the identity, which should be correct.

I know that I have premature termination. Below I ran it again until convergence (at 750 iterations), but the result is still incorrect. Being quite a simple registration problem I believe that I must be doing something else wrong.

I set up my registration like this:

> ```
> R = sitk.ImageRegistrationMethod()
> R.SetMetricAsMeanSquares()
> R.SetOptimizerAsRegularStepGradientDescent(4.0, 0.01, 2000, 0.5, 1e-6)
> R.SetInitialTransform(sitk.CenteredTransformInitializer(fixed, moving,
> sitk.Euler3DTransform(), sitk.CenteredTransformInitializerFilter.GEOMETRY))
> R.SetInterpolator(sitk.sitkLinear)
> 
> ```

And the result is

> ```
> itk::simple::Transform
> Euler3DTransform (0x7f838ab99f50)
> RTTI typeinfo: itk::Euler3DTransform<double>
> Reference Count: 2
> Modified Time: 5226
> Debug: Off
> Object Name: 
> Observers: 
> none
> Matrix: 
> 0.994857 -0.0991812 -0.0205532 
> 0.0967002 0.990413 -0.0986452 
> 0.0301399 0.0961504 0.99491 
> Offset: [9.50359, 7.69879, -9.5203]
> Center: [0, 0, 0]
> Translation: [9.50359, 7.69879, -9.5203]
> Inverse: 
> 0.994857 0.0967002 0.0301399 
> -0.0991812 0.990413 0.0961504 
> -0.0205532 -0.0986452 0.99491 
> Singular: 0
> Euler's angles: AngleX=3.04529 AngleY=-3.17188 AngleZ=3.2414
> m_ComputeZYX = 0
> 
> Optimizer stop condition: RegularStepGradientDescentOptimizerv4: Step too small after 751 iterations. Current step (0.0078125) is less than minimum step (0.01).
> Iteration: 752
> Metric value: 8.91753355411e-06
> 
> ```

I believe that using sitk.CenteredTransformInitializerFilter.GEOMETRY is the right choice since this is a PET image being registered to a slightly shifted version of itself. I expect the registration to return a transformation of [15, 15, -15] and angles of near [0, 0, 0].

Can you possibly see if I am doing something wrong? Or can you suggest a simple case that I can test with?

Many thanks!

---

<div class="post-metadata">

### Author: ![zivy](https://discourse.itk.org/user_avatar/discourse.itk.org/zivy/32/1726_2.png) [@zivy](https://discourse.itk.org/u/zivy)
#### Post date: [June 24, 2020, 1:19pm UTC](https://discourse.itk.org/t/struggling-to-do-simple-3d-rigid-registration/3217/4 "2020-06-24T13:19:01Z")

</div>

Hello @msb,

First step in debugging a registration is to visualize the initial condition after you used the `CenteredTransformInitializer`.

Second look at your code, it seems you are not calling `R.SetOptimizerScalesFromPhysicalShift()` this is a critical step if your transformation is rigid as the units of the rotation parameters and translation parameters are not commensurate (a unit change in rotation is not equivalent to a unit change in translation). If you are only doing translation then you don’t need to call the method because all optimized parameters are the same.

In your case you were describing the problem as a translation only, but in your code you use `Euler3DTransform`.

---

<div class="post-metadata">

### Author: ![msb](https://discourse.itk.org/user_avatar/discourse.itk.org/msb/32/1588_2.png) [@msb](https://discourse.itk.org/u/msb)
#### Post date: [June 24, 2020, 11:54pm UTC](https://discourse.itk.org/t/struggling-to-do-simple-3d-rigid-registration/3217/5 "2020-06-24T23:54:14Z")

</div>

Ha! It was the call to SetOptimizerScalesFromPhysicalShift()! That has fixed it it seems.

Thanks very much for your help @zivy.