[Insight-users] Model-Based-Segmentation/Registration (2nd)

Gunther Sudra Gunther.Sudra at web.de
Thu, 12 Feb 2004 17:12:11 +0100


Hi,

this time I included my code. My Question is now how to tell the optimizer to use the radius of 
the ellipse additionally for registration?

Background: 
Goal is to register the ellipse in the synthetic image with the model-ellipse.

/////////////////////////////////////////////////////////////////////////////////////////////////

template < class TOptimizer >
class IterationCallback : public itk::Command 
{

... as in example model-based-registraton

}

template <typename TFixedImage, typename TMovingSpatialObject>
class SimpleImageToSpatialObjectMetric : 
    public itk::ImageToSpatialObjectMetric<TFixedImage,TMovingSpatialObject>
{

... as in example model-based-registraton

}


my function:
{

                // synthetic image:
	// one ellipse at (64,64) with radius (10,20)
	// in one group
	
	typedef   float           InternalPixelType;
			
	typedef itk::Image< InternalPixelType, 2 >  InternalImageType; 
	
	typedef itk::EllipseSpatialObject< 2 >   EllipseType;

	EllipseType::Pointer ellipse_toFind = EllipseType::New();		
	
	double radius_toFind[2] = {10.0, 20.0};
	 
	ellipse_toFind->SetRadius( radius_toFind );
	
	EllipseType::TransformType::OffsetType offset;

	offset[ 0 ] = 64;
	offset[ 1 ] = 64;
	ellipse_toFind->GetTransform()->SetOffset(offset);
	ellipse_toFind->ComputeGlobalTransform();
	
	typedef itk::GroupSpatialObject< 2 >     GroupType;
	
	GroupType::Pointer group_toFind = GroupType::New();

	group_toFind->AddSpatialObject( ellipse_toFind );



	typedef itk::SpatialObjectToImageFilter<  
                                   GroupType,
                                   InternalImageType >   SpatialObjectToImageFilterType;
  

	SpatialObjectToImageFilterType::Pointer imageFilter = 
                                     SpatialObjectToImageFilterType::New();
  
	imageFilter->SetInput(  group_toFind  );
  
  
	InternalImageType::SizeType size;
	size[ 0 ] = 128;
	size[ 1 ] = 128;
 
	imageFilter->SetSize( size );
  
	imageFilter->Update();
  

	//model which should be registered with the synthetic image
	//one ellipse at (64,64) with radius (10,30)
	//goal is to reduce radius to (10,20)


	EllipseType::Pointer ellipse_result = EllipseType::New();

	double radius_result[2] = {10.0, 30.0};
	 
	ellipse_result->SetRadius( radius_result );
	
	
	offset[ 0 ] = 64;
	offset[ 1 ] = 64;
	ellipse_result->GetTransform()->SetOffset(offset);
	ellipse_result->ComputeGlobalTransform();

                GroupType::Pointer group_result = GroupType::New();

	group_result->AddSpatialObject( ellipse_result );


	//registration method
	typedef itk::ImageToSpatialObjectRegistrationMethod<
                                      InternalImageType,
                                      GroupType  >  RegistrationType;

	RegistrationType::Pointer registration = RegistrationType::New();
  
	//metric from softwareguide example model-based registration
	typedef SimpleImageToSpatialObjectMetric<  InternalImageType,
                                             GroupType   > MetricType;

	MetricType::Pointer metric = MetricType::New();
  
	//linear interpolator
	typedef itk::LinearInterpolateImageFunction< 
                                         InternalImageType,
                                         double     >  InterpolatorType;

	InterpolatorType::Pointer interpolator = InterpolatorType::New();
  
	//optimizer
	typedef itk::OnePlusOneEvolutionaryOptimizer  OptimizerType;

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


	
	typedef itk::Euler2DTransform<> TransformType;

	TransformType::Pointer transform = TransformType::New();
  
	itk::Statistics::NormalVariateGenerator::Pointer generator 
                      = itk::Statistics::NormalVariateGenerator::New();
  
	generator->Initialize(12345);
  
	optimizer->SetNormalVariateGenerator( generator );
	
	optimizer->Initialize( 20 );
	optimizer->SetMaximumIteration( 400 );
  

	TransformType::ParametersType parametersScale;
	parametersScale.resize(3);

	parametersScale[0] = 1000; // angle scale

	for( unsigned int i=1; i<3; i++ )
	    {
		parametersScale[i] = 2; // offset scale
		}
	optimizer->SetScales( parametersScale );
  
	typedef IterationCallback< OptimizerType >   IterationCallbackType;

	IterationCallbackType::Pointer callback = IterationCallbackType::New();

	callback->SetOptimizer( optimizer );
  
	  
	registration->SetFixedImage( imageFilter->GetOutput() );
	registration->SetMovingSpatialObject( group_result );
	registration->SetTransform( transform );
	registration->SetInterpolator( interpolator );
	registration->SetOptimizer( optimizer );
	registration->SetMetric( metric );
  

	TransformType::ParametersType initialParameters( 
                      transform->GetNumberOfParameters() );
  
	initialParameters[0] = 0.0;     // Angle
	initialParameters[1] = 0.0;     // Offset X
	initialParameters[2] = 0.0;     // Offset Y 

	registration->SetInitialTransformParameters(initialParameters);
  

	std::cout << "Initial Parameters  : " << initialParameters << std::endl;

	optimizer->MaximizeOn();
  
	try {
	    registration->StartRegistration();
	     }
	catch( itk::ExceptionObject & exp ) {
	    std::cerr << "Exception caught ! " << std::endl;
	    std::cerr << exp << std::endl;
                    }
	RegistrationType::ParametersType finalParameters 
                         = registration->GetLastTransformParameters();

	std::cout << "Final Solution is : " << finalParameters << std::endl;
}
______________________________________________________________________________
Nachrichten, Musik und Spiele schnell und einfach per Quickstart im 
WEB.DE Screensaver - Gratis downloaden: http://screensaver.web.de/?mc=021110