[Insight-users] Plotting the metric

Luis Ibanez luis.ibanez at kitware.com
Wed Jun 18 09:45:32 EDT 2008


Hi Giorgos,

Thanks for the detailed description of the problem.

  It seems that you are not setting the center of rotation
  of the Euler2DTransform. Therefore the rotation will happen
  around the (0,0) point, which may be located in the corner
  of your images, or even outside of the images, depending of
  what the image origin settings are in your Fixed and Moving
  image.

  You should consider introducing the CenteredTransformInitializer
  to initialize the parameters of the Transform, in particular its
  center of rotation.


  Note that, if you do this, then you should remove the lines

 >           parameters[1] = 0;
 >           parameters[2] = 0;

  from the body of the for-loop, and only keep the line:

 >           parameters[0] = i;




   Regards,


      Luis



---------------------
Giorgos Pichis wrote:
> Dear all,
> 
> I am trying to check the metric (Mattes MI) values between two images, 
> while rotating the one over the other.
> Therefore I use Euler2DTransform.
> The results though of the plot aren't the expected.
> I am comparing two test images where the one is rotated 10 degrees over 
> the other,
> and the metric minimum is appeared in the value 1.
> Do you know what might be the problem, or what should I change in the 
> following code?
> 
> Thanks a lot
> Best Regards
> Giorgos
> ........................................................................
> #include "itkImage.h"
> #include "itkImageFileReader.h"
> #include "itkImageFileWriter.h"
> #include "itkMattesMutualInformationImageToImageMetric.h"
> #include "itkEuler2DTransform.h"
> #include "itkLinearInterpolateImageFunction.h"
> 
> int main( int argc, char * argv[] )
> {
>   if( argc < 3 )
>     {
>     std::cerr << "Usage: " << std::endl;
>     std::cerr << argv[0] << "  fixedImage  movingImage" << std::endl;
>     return 1;
>     }
> 
>   const     unsigned int   Dimension = 2;
>   typedef   unsigned char  PixelType;
>   typedef itk::Image< PixelType, Dimension >   ImageType;
>   typedef itk::Image< PixelType, Dimension >   ImageType;
>   typedef itk::ImageFileReader< ImageType >  ReaderType;
> 
>   ReaderType::Pointer fixedReader  = ReaderType::New();
>   ReaderType::Pointer movingReader = ReaderType::New();
> 
>   fixedReader->SetFileName(  argv[ 1 ] );
>   movingReader->SetFileName( argv[ 2 ] );
> 
>   try
>     {
>     fixedReader->Update();
>     movingReader->Update();
>     }
>   catch( itk::ExceptionObject & excep )
>     {
>     std::cerr << "Exception catched !" << std::endl;
>     std::cerr << excep << std::endl;
>     }
> 
>   typedef itk::MattesMutualInformationImageToImageMetric< ImageType, 
> ImageType >  MetricType;
>   MetricType::Pointer metric = MetricType::New();
> 
>   typedef itk::Euler2DTransform< double >  TransformType;
>   TransformType::Pointer transform = TransformType::New();
> 
>   typedef itk::LinearInterpolateImageFunction<
>                                  ImageType, double >  InterpolatorType;
>   InterpolatorType::Pointer interpolator = InterpolatorType::New();
> 
>   metric->SetInterpolator( interpolator );
>   metric->SetTransform( transform );
> 
>   metric->SetNumberOfHistogramBins( 20 );
>   metric->SetNumberOfSpatialSamples( 10000 );
> 
>  transform->SetIdentity();
> 
>   ImageType::ConstPointer fixedImage  = fixedReader->GetOutput();
>   ImageType::ConstPointer movingImage = movingReader->GetOutput();
> 
>   metric->SetFixedImage(  fixedImage  );
>   metric->SetMovingImage( movingImage );
> 
>   metric->SetFixedImageRegion(  fixedImage->GetBufferedRegion()  );
> 
>   try
>     {
>     metric->Initialize();
>     }
>   catch( itk::ExceptionObject & excep )
>     {
>     std::cerr << "Exception catched !" << std::endl;
>     std::cerr << excep << std::endl;
>     return -1;
>     }
> 
>       MetricType::TransformParametersType parameters( 
> transform->GetNumberOfParameters() );
>  
>       int RangeOfDegrees=50;
>       double RangeInRadians=RangeOfDegrees * atan(1.0)/45.0;
>      
>       for (double i=-RangeInRadians; i<=RangeInRadians; i+=atan(1.0)/45.0)
>       {
>           parameters[0] = i;
>           parameters[1] = 0;
>           parameters[2] = 0;
>          
>           const double value = metric->GetValue( parameters );
>           double finalAngleInDegrees = i * 45.0 / atan(1.0);
>           std::cout << finalAngleInDegrees << "         "  << value << 
> std::endl;
>       }
>   std::cout << std::endl;
>   return 0;
> }
> 
> 
> ------------------------------------------------------------------------
> 
> _______________________________________________
> Insight-users mailing list
> Insight-users at itk.org
> http://www.itk.org/mailman/listinfo/insight-users


More information about the Insight-users mailing list