[Insight-users] I need your help about registration optimizer algorithm

[于荷峰]

Hello everyone,

I do a registration example in ITK,it is “ImageRegistration1.cxx”.But I
want to use “itkParticleSwarmOptimizer.h” instead of
"itkRegularStepGradientDescentOptimizer.h".I have modified this
programming,but it has been an error.



The code is as follows:



#include "itkImageRegistrationMethod.h"

#include "itkTranslationTransform.h"

#include "itkMeanSquaresImageToImageMetric.h"

#include "itkParticleSwarmOptimizer.h"

// Software Guide : EndCodeSnippet



#include "itkImageFileReader.h"

#include "itkImageFileWriter.h"


#include "itkResampleImageFilter.h"

#include "itkCastImageFilter.h"

#include "itkRescaleIntensityImageFilter.h"

#include "itkSubtractImageFilter.h"




class CommandIterationUpdate : public itk::Command

{

public:

  typedef  CommandIterationUpdate   Self;

  typedef  itk::Command             Superclass;

  typedef itk::SmartPointer<Self>  Pointer;

  itkNewMacro( Self );


protected:

  CommandIterationUpdate() {};


public:


  typedef itk::ParticleSwarmOptimizer     OptimizerType;

  typedef const OptimizerType                         *OptimizerPointer;


  void Execute(itk::Object *caller, const itk::EventObject & event)

  {

    Execute( (const itk::Object *)caller, event);

  }


  void Execute(const itk::Object * object, const itk::EventObject & event)

  {

    OptimizerPointer optimizer =

                         dynamic_cast< OptimizerPointer >( object );


    if( ! itk::IterationEvent().CheckEvent( &event ) )

      {

      return;

      }


    //std::cout << optimizer->GetMaximalNumberOfIterations () << " = ";

    std::cout << optimizer->GetValue() << " : ";

    std::cout << optimizer->GetCurrentPosition() << std::endl;

  }


};



int main( int argc, char *argv[] )

{



  const    unsigned int    Dimension = 2;

  typedef  float           PixelType;



  typedef itk::Image< PixelType, Dimension >  FixedImageType;

  typedef itk::Image< PixelType, Dimension >  MovingImageType;



  typedef itk::TranslationTransform< double, Dimension > TransformType;



  typedef itk::ParticleSwarmOptimizer       OptimizerType;

  static OptimizerType::RandomVariateGeneratorType::IntegerType seedOffset
= 0;


  typedef itk::MeanSquaresImageToImageMetric<

                                    FixedImageType,

                                    MovingImageType >    MetricType;



  typedef itk:: LinearInterpolateImageFunction<

                                    MovingImageType,

                                    double          >    InterpolatorType;



  typedef itk::ImageRegistrationMethod<

                                    FixedImageType,

                                    MovingImageType >    RegistrationType;



  MetricType::Pointer         metric        = MetricType::New();

  TransformType::Pointer      transform     = TransformType::New();

  OptimizerType::Pointer      optimizer     = OptimizerType::New();

  InterpolatorType::Pointer   interpolator  = InterpolatorType::New();

  RegistrationType::Pointer   registration  = RegistrationType::New();



  registration->SetMetric(        metric        );

  registration->SetOptimizer(     optimizer     );

  registration->SetTransform(     transform     );

  registration->SetInterpolator(  interpolator  );

  // Software Guide : EndCodeSnippet



  typedef itk::ImageFileReader< FixedImageType  > FixedImageReaderType;

  typedef itk::ImageFileReader< MovingImageType > MovingImageReaderType;

  FixedImageReaderType::Pointer  fixedImageReader  =
FixedImageReaderType::New();

  MovingImageReaderType::Pointer movingImageReader =
MovingImageReaderType::New();


  fixedImageReader->SetFileName(  "BrainProtonDensitySliceBorder20.png" );

  movingImageReader->SetFileName(
"BrainProtonDensitySliceShifted13x17y.png" );





  registration->SetFixedImage(    fixedImageReader->GetOutput()    );

  registration->SetMovingImage(   movingImageReader->GetOutput()   );



  fixedImageReader->Update();

  registration->SetFixedImageRegion(

                    fixedImageReader->GetOutput()->GetBufferedRegion() );



  typedef RegistrationType::ParametersType ParametersType;

  ParametersType initialParameters( transform->GetNumberOfParameters() );


  initialParameters[0] = 0.0;  // Initial offset in mm along X

  initialParameters[1] = 0.0;  // Initial offset in mm along Y


  registration->SetInitialTransformParameters( initialParameters );





  optimizer->UseSeedOn();

  optimizer->SetSeed(8775070 + seedOffset++);


  OptimizerType::ParameterBoundsType bounds;

  bounds.push_back( std::make_pair( 0, 1 ) );

  optimizer->SetParameterBounds( bounds );

  optimizer->SetNumberOfParticles( 10 );//numberOfParticles

  optimizer->SetMaximalNumberOfIterations( 100 );//maxIterations

  optimizer->SetParametersConvergenceTolerance(
0.1,transform->GetNumberOfParameters() );//1、xTolerance

  optimizer->SetFunctionConvergenceTolerance( 0.001 );//fTolerance




  CommandIterationUpdate::Pointer observer = CommandIterationUpdate::New();

  optimizer->AddObserver( itk::IterationEvent(), observer );





  try

    {

    registration->Update();

std::cout << "Optimizer stop condition: "

              << registration->GetOptimizer()->GetStopConditionDescription()

              << std::endl;

    }

  catch( itk::ExceptionObject & err )

    {

    std::cerr << "ExceptionObject caught !" << std::endl;

    std::cerr << err << std::endl;

    return EXIT_FAILURE;

    }



  ParametersType finalParameters =
registration->GetLastTransformParameters();



  const double TranslationAlongX = finalParameters[0];

  const double TranslationAlongY = finalParameters[1];





  const double bestValue = optimizer->GetValue();



  std::cout << "Result = " << std::endl;

  std::cout << " Translation X = " << TranslationAlongX  << std::endl;

  std::cout << " Translation Y = " << TranslationAlongY  << std::endl;

  //std::cout << " Iterations    = " << numberOfIterations << std::endl;

  std::cout << " Metric value  = " << bestValue          << std::endl;





  typedef itk::ResampleImageFilter<

                            MovingImageType,

                            FixedImageType >    ResampleFilterType;



  ResampleFilterType::Pointer resampler = ResampleFilterType::New();

  resampler->SetInput( movingImageReader->GetOutput() );



  resampler->SetTransform( registration->GetOutput()->Get() );



  FixedImageType::Pointer fixedImage = fixedImageReader->GetOutput();

  resampler->SetSize( fixedImage->GetLargestPossibleRegion().GetSize() );

  resampler->SetOutputOrigin(  fixedImage->GetOrigin() );

  resampler->SetOutputSpacing( fixedImage->GetSpacing() );

  resampler->SetOutputDirection( fixedImage->GetDirection() );

  resampler->SetDefaultPixelValue( 100 );



  typedef unsigned char OutputPixelType;

  typedef itk::Image< OutputPixelType, Dimension > OutputImageType;

  typedef itk::CastImageFilter<

                        FixedImageType,

                        OutputImageType > CastFilterType;

  typedef itk::ImageFileWriter< OutputImageType >  WriterType;



  WriterType::Pointer      writer =  WriterType::New();

  CastFilterType::Pointer  caster =  CastFilterType::New();




  writer->SetFileName( "output.png" );





  caster->SetInput( resampler->GetOutput() );

  writer->SetInput( caster->GetOutput()   );

  writer->Update();



  typedef itk::SubtractImageFilter<

                                  FixedImageType,

                                  FixedImageType,

                                  FixedImageType > DifferenceFilterType;


  DifferenceFilterType::Pointer difference = DifferenceFilterType::New();


  difference->SetInput1( fixedImageReader->GetOutput() );

  difference->SetInput2( resampler->GetOutput() );



  typedef itk::RescaleIntensityImageFilter<

                                  FixedImageType,

                                  OutputImageType >   RescalerType;


  RescalerType::Pointer intensityRescaler = RescalerType::New();


  intensityRescaler->SetInput( difference->GetOutput() );

  intensityRescaler->SetOutputMinimum(   0 );

  intensityRescaler->SetOutputMaximum( 255 );


  resampler->SetDefaultPixelValue( 1 );



  WriterType::Pointer writer2 = WriterType::New();

  writer2->SetInput( intensityRescaler->GetOutput() );

  // Software Guide : EndCodeSnippet



  /*if( argc > 4 )

    {*/

    writer2->SetFileName( "differenceAfter.png" );

    writer2->Update();

    /*}*/






  TransformType::Pointer identityTransform = TransformType::New();

  identityTransform->SetIdentity();

  resampler->SetTransform( identityTransform );




  /*if( argc > 5 )

    {*/

    writer2->SetFileName( "differenceBefor.png" );

    writer2->Update();

    /*}*/








  return EXIT_SUCCESS;

}



But an error been caught ,when you run it.The error information is as
follow:

Location:“void _thiscall
itk::ParticleSwarmOptimizerBase::validateSettings<vod>”

File:D:\..\itkParticleSwarmOptimizerBase.cxx

Line:360

Description:itk::ERROR:ParticleSwarmOptimizer<02A022C8>:cost function and
parameter bounds dimensions mismatch



Please give some help............................
-------------- next part --------------
An HTML attachment was scrubbed...
URL: <http://www.itk.org/pipermail/insight-users/attachments/20130510/6a0bfb05/attachment.htm>