[ITK] [ITK-users] 2D rigid transformation

Bill Lorensen bill.lorensen at gmail.com
Tue May 19 11:20:06 EDT 2015


It would be useful if you included a compilable, minimal example.

I include your code snippet in a small program that reads two files
and performs the optimization. I ran it with the fixed and moving set
to the same image. Here are my results:

  InitialPosition = [-6.283185307179586, 0, 0]
  CurrentValue = 3528.15
  NumberOfSteps = [6, 0, 0]
  Stop = 1
  StepLength = 1
  CurrentIndex = [0, 0, 0]
  MaximumMetricValue = 3528.15
  MinimumMetricValue = 4.50585e-27
  MinimumMetricValuePosition = [-6.283185307179586, 0, 0]
  MaximumMetricValuePosition = [-3.141592653589793, 0, 0]

I suspect that your images are filled with zeroes. You do not show how
you populate the images.

Here is my code:
#include "itkImageFileReader.h"
#include "itkEuler2DTransform.h"
#include "itkExhaustiveOptimizerv4.h"
#include "itkMeanSquaresImageToImageMetricv4.h"
#include "itkCenteredTransformInitializer.h"
#include "itkImageRegistrationMethodv4.h"
#include "itkImage.h"

int main (int argc, char *argv[])
{
  if (argc < 3)
    {
    std::cout << "Usage: " << argv[0] << " fixedImage movingImage" << std::endl;
    return EXIT_FAILURE;
    }
  typedef itk::Image<double, 2>                 FixedImageType;
  typedef itk::Image<double, 2>                 MovingImageType;
  typedef itk::ImageFileReader<FixedImageType>  FixedImageReaderType;
  typedef itk::ImageFileReader<MovingImageType> MovingImageReaderType;
  typedef itk::Euler2DTransform< double >       TransformType;
  typedef itk::ExhaustiveOptimizerv4< double >  OptimizerType;
  typedef itk::MeanSquaresImageToImageMetricv4< FixedImageType,
MovingImageType >
                                                MetricType;
  typedef itk::CenteredTransformInitializer< TransformType,
FixedImageType,  MovingImageType >
                                                TransformInitializerType;
  typedef itk::ImageRegistrationMethodv4< FixedImageType,
MovingImageType, TransformType >
                                                RegistrationType;

  FixedImageReaderType::Pointer     fixedImageReader    =
FixedImageReaderType::New();
  MovingImageReaderType::Pointer    movingImageReader    =
MovingImageReaderType::New();
  FixedImageType::Pointer           fixedImage    = FixedImageType::New();
  MovingImageType::Pointer          movingImage    = MovingImageType::New();
  TransformType::Pointer            transform    = TransformType::New();
  MetricType::Pointer               metric       = MetricType::New();
  OptimizerType::Pointer            optimizer    = OptimizerType::New();
  RegistrationType::Pointer         registration = RegistrationType::New();
  TransformInitializerType::Pointer initializer  =
TransformInitializerType::New();

  fixedImageReader->SetFileName (argv[1]);
  fixedImageReader->Update();
  fixedImage = fixedImageReader->GetOutput();

  movingImageReader->SetFileName (argv[2]);
  movingImageReader->Update();
  movingImage = movingImageReader->GetOutput();

  unsigned int angles = 12;
  OptimizerType::StepsType steps( transform->GetNumberOfParameters() );
  steps[0] = int(angles/2);
  steps[1] = 0;
  steps[2] = 0;
  optimizer->SetNumberOfSteps( steps );

  OptimizerType::ScalesType scales( transform->GetNumberOfParameters() );
  scales[0] = 2.0 * vnl_math::pi / angles;
  scales[1] = 1.0;
  scales[2] = 1.0;

  optimizer->SetScales( scales );

  initializer->SetTransform(   transform );
  initializer->SetFixedImage(  fixedImage );
  initializer->SetMovingImage( movingImage );
  initializer->InitializeTransform();

// Initialize registration
  registration->SetMetric(           metric );
  registration->SetOptimizer(        optimizer );
  registration->SetFixedImage(       fixedImage );
  registration->SetMovingImage(      movingImage );
  registration->SetInitialTransform( transform );

  try
    {
    registration->Update();
    optimizer->Print(std::cout);
    }
  catch( itk::ExceptionObject & err )
    {
    std::cerr << "ExceptionObject caught !" << std::endl;
    std::cerr << err << std::endl;
    return EXIT_FAILURE;
    }
  return EXIT_SUCCESS;
}


On Tue, May 19, 2015 at 5:41 AM, Pietro Nardelli
<p.nardelli at umail.ucc.ie> wrote:
> Hi guys,
>
> thank you very much for your help. I was trying to follow the examples you
> sent me (converting it in native C++ ITK code), but it seems that the
> optimizer does not work properly. In particular, no matter the number of
> angles and the step length I choose, it always gives me the same result. It
> seems to me that the optimizer is never updating the angles, therefore it is
> not following the grid I create. Also the metric is always 0. There is
> definitely something wrong! Here is part of the script I am using and the
> result I always get. Am I overlooking anything?
>
>   typedef itk::Euler2DTransform< double >
> TransformType;
>   typedef itk::ExhaustiveOptimizerv4< double >
> OptimizerType;
>   typedef itk::MeanSquaresImageToImageMetricv4< FixedInputImageType,
> MovingInputImageType >               MetricType;
>   typedef itk::CenteredTransformInitializer< TransformType,
> FixedInputImageType,  MovingInputImageType > TransformInitializerType;
>   typedef itk::ImageRegistrationMethodv4< FixedInputImageType,
> MovingInputImageType, TransformType > RegistrationType;
>
>   typename TransformType::Pointer               transform    =
> TransformType::New();
>   typename MetricType::Pointer                      metric       =
> MetricType::New();
>   typename OptimizerType::Pointer                optimizer    =
> OptimizerType::New();
>   typename RegistrationType::Pointer            registration =
> RegistrationType::New();
>   typename TransformInitializerType::Pointer initializer  =
> TransformInitializerType::New();
>
>   unsigned int angles = 12;
>   OptimizerType::StepsType steps( transform->GetNumberOfParameters() );
>   steps[0] = int(angles/2);
>   steps[1] = 0;
>   steps[2] = 0;
>   optimizer->SetNumberOfSteps( steps );
>
>   OptimizerType::ScalesType scales( transform->GetNumberOfParameters() );
>   scales[0] = 2.0 * vnl_math::pi / angles;
>   scales[1] = 1.0;
>   scales[2] = 1.0;
>
>   optimizer->SetScales( scales );
>
>   initializer->SetTransform(   transform );
>   initializer->SetFixedImage(  fixImage );
>   initializer->SetMovingImage( movImage );
>   //initializer->GeometryOn();
>   initializer->InitializeTransform();
>
>   // Initialize registration
>   registration->SetMetric(                  metric    );
>   registration->SetOptimizer(            optimizer );
>   registration->SetFixedImage(         fixedImage  );
>   registration->SetMovingImage(      movingImage  );
>   registration->SetInitialTransform(   transform );
>
>   try
>   {
>     registration->Update();
>   }
>   catch( itk::ExceptionObject & err )
>   {
>     std::cerr << "ExceptionObject caught !" << std::endl;
>     std::cerr << err << std::endl;
>     return EXIT_FAILURE;
>   }
>
> Result:
>
> Final parameters: [-1.8849555921538759, -40, 40]
> Result =
> Metric value  = 0
> Angle (rad) =  -1.88496
> Angle (degrees) =  -108
> Iterations    = 13
> Rotation Center = [225.5, 185.5]
> ExhaustiveOptimizerv4: Completed sampling of parametric space of size 3
>
>
> Thank you very much,
> Pietro
>
> Pietro Nardelli, MEngSc, BE
>
> Ph.D Candidate
>
> School of Engineering
>
> Electrical & Electronic Engineering
>
> University College Cork
>
> College Road
>
> Cork, Ireland
>
>
> On Fri, May 15, 2015 at 7:21 PM, Yaniv, Ziv Rafael (NIH/NLM/LHC) [C]
> <zivrafael.yaniv at nih.gov> wrote:
>>
>> Hello Pietro,
>>
>> You should use the ExhaustiveOptimizerv4
>> (http://www.itk.org/Doxygen/html/classitk_1_1ExhaustiveOptimizerv4.html)
>> which allows you to set a grid on which the similarity metric is evaluated.
>>
>> If you are familiar with python, then the following SimpleITK notebook may
>> be of use to you (see last section):
>> https://github.com/zivy/SimpleITK-Notebook-Staging/blob/master/registration3.ipynb
>>
>>      regards
>>               Ziv
>>
>> From: Pietro Nardelli
>> <p.nardelli at umail.ucc.ie<mailto:p.nardelli at umail.ucc.ie>>
>> Date: Friday, May 15, 2015 at 2:13 PM
>> To: "insight-users at itk.org<mailto:insight-users at itk.org>"
>> <insight-users at itk.org<mailto:insight-users at itk.org>>
>> Subject: [ITK-users] 2D rigid transformation
>>
>> Hello guys,
>>
>> is there a way to have a 2D rigid registration that uses a specific number
>> of rotations and chooses the best one? I have two images that are simply
>> rotated with respect to each other, and I would like to register them using
>> for example 36 rotations (therefore computing the mean squared error every
>> 10 degrees). At the moment I am using the 2DRigidTransform with a specified
>> center, with a regular step descent optimizer and the
>> ImageRegistrationMethodv4. I saw that the transform has the function
>> SetFixedParameters() but I am not really sure whether I understand correctly
>> that that would tell the optimizer the angles (and translations) to use at
>> every iteration. Could anyone please clarify this?
>>
>> Thank you very much,
>> Pietro
>>
>>
>
>
> _____________________________________
> Powered by www.kitware.com
>
> Visit other Kitware open-source projects at
> http://www.kitware.com/opensource/opensource.html
>
> Kitware offers ITK Training Courses, for more information visit:
> http://www.kitware.com/products/protraining.php
>
> Please keep messages on-topic and check the ITK FAQ at:
> http://www.itk.org/Wiki/ITK_FAQ
>
> Follow this link to subscribe/unsubscribe:
> http://public.kitware.com/mailman/listinfo/insight-users
>



-- 
Unpaid intern in BillsBasement at noware dot com
_____________________________________
Powered by www.kitware.com

Visit other Kitware open-source projects at
http://www.kitware.com/opensource/opensource.html

Kitware offers ITK Training Courses, for more information visit:
http://www.kitware.com/products/protraining.php

Please keep messages on-topic and check the ITK FAQ at:
http://www.itk.org/Wiki/ITK_FAQ

Follow this link to subscribe/unsubscribe:
http://public.kitware.com/mailman/listinfo/insight-users


More information about the Community mailing list