[Insight-users] NonRigidRegistration using MutualInformation
Luis Ibanez
luis.ibanez at kitware.com
Fri Nov 9 16:42:46 EST 2007
Hi Matsuo,
1) Are you building your application for "Release" ?
(make sure that it is not build for "Debug".
2) How much time in seconds is taking per iteration ?
3) Is the metric value improving between iterations ?
It is normal for the combination that you are using
to be computationally demanding. However, depending
to your answer to question (2), there may be things
that you may setup better.
Please let us know,
Thanks
Luis
----------------
松尾薫 wrote:
> Dear Insight user,
>
> I am trying to 3D NonRigidRegistration(CT-CT).
> The components i'm using are
> BSplineDeformableTransform,
> LBFGSBOptimizer,
> NormalizedMutualInformationHistogramImageToImageMetric
> (don't use MultiResolution)
>
> I have a 3D images 280*280*201 and spacing 1*1*1.
> The time it takes to calculate each itertion is veeeeery
> slow.(Processing doesn't end......)
>
> How can I solve that????
>
> Regards.
> Matsuo
>
> Code----------------------------------------------------------------------------------
> const unsigned int ImageDimension = 3;
> typedef signed short PixelType;
>
> typedef itk::Image< PixelType, ImageDimension > FixedImageType;
> typedef itk::Image< PixelType, ImageDimension > MovingImageType;
>
> const unsigned int SpaceDimension = ImageDimension;
> const unsigned int SplineOrder = 3;
> typedef double CoordinateRepType;
> typedef itk::BSplineDeformableTransform< CoordinateRepType,
> SpaceDimension, SplineOrder > TransformType;
> typedef itk::LBFGSBOptimizer OptimizerType;
> typedef itk::NormalizedMutualInformationHistogramImageToImageMetric<
> FixedImageType, MovingImageType > MetricType;
> typedef itk::ImageRegistrationMethod< FixedImageType,
> MovingImageType > RegistrationType;
> typedef itk::LinearInterpolateImageFunction< MovingImageType, double
> > InterpolatorType;
>
> MetricType::Pointer metric = MetricType::New();
> OptimizerType::Pointer optimizer = OptimizerType::New();
> InterpolatorType::Pointer interpolator = InterpolatorType::New();
> RegistrationType::Pointer registration = RegistrationType::New();
> TransformType::Pointer transform = TransformType::New();
>
> registration->SetMetric( metric );
> registration->SetOptimizer( optimizer );
> registration->SetInterpolator( interpolator );
> registration->SetTransform( transform );
>
> typedef itk::ImageFileReader< FixedImageType > FixedImageReaderType;
> typedef itk::ImageFileReader< MovingImageType > MovingImageReaderType;
>
> FixedImageReaderType::Pointer fixedImageReader =
> FixedImageReaderType::New();
> MovingImageReaderType::Pointer movingImageReader =
> MovingImageReaderType::New();
>
> fixedImageReader->SetFileName( argv[1] );
> movingImageReader->SetFileName( argv[2] );
> FixedImageType::ConstPointer fixedImage = fixedImageReader->GetOutput();
>
> registration->SetFixedImage( fixedImage );
> registration->SetMovingImage( movingImageReader->GetOutput() );
>
> fixedImageReader->Update();
>
> FixedImageType::RegionType fixedRegion =
> fixedImage->GetBufferedRegion();
> registration->SetFixedImageRegion( fixedRegion );
>
> typedef TransformType::RegionType RegionType;
> RegionType bsplineRegion;
> RegionType::SizeType gridSizeOnImage;
> RegionType::SizeType gridBorderSize;
> RegionType::SizeType totalGridSize;
>
> gridSizeOnImage.Fill( 8 );
> gridBorderSize.Fill( 3 );
> totalGridSize = gridSizeOnImage + gridBorderSize;
>
> bsplineRegion.SetSize( totalGridSize );
>
> typedef TransformType::SpacingType SpacingType;
> SpacingType spacing = fixedImage->GetSpacing();
>
> typedef TransformType::OriginType OriginType;
> OriginType origin = fixedImage->GetOrigin();;
>
> FixedImageType::SizeType fixedImageSize = fixedRegion.GetSize();
>
> for(unsigned int r=0; r<ImageDimension; r++)
> {
> spacing[r] *= floor( static_cast<double>(fixedImageSize[r] - 1) /
> static_cast<double>(gridSizeOnImage[r] - 1) );
> origin[r] -= spacing[r];
> }
>
> transform->SetGridSpacing( spacing );
> transform->SetGridOrigin( origin );
> transform->SetGridRegion( bsplineRegion );
>
>
> typedef TransformType::ParametersType ParametersType;
> const unsigned int numberOfParameters =
> transform->GetNumberOfParameters();
> ParametersType parameters( numberOfParameters );
> parameters.Fill( 0.0 );
> transform->SetParameters( parameters );
> registration->SetInitialTransformParameters(
> transform->GetParameters() );
>
> unsigned int numberOfHistogramBins = 256;
> MetricType::HistogramType::SizeType histogramSize;
> histogramSize[0] = numberOfHistogramBins;
> histogramSize[1] = numberOfHistogramBins;
> metric->SetHistogramSize( histogramSize );
>
> typedef MetricType::ScalesType ScalesType;
> ScalesType scales( numberOfParameters );
>
> scales.Fill ( 1.0 );
>
> metric->SetDerivativeStepLengthScales(scales);
>
> OptimizerType::BoundSelectionType boundSelect(
> transform->GetNumberOfParameters() );
> OptimizerType::BoundValueType upperBound(
> transform->GetNumberOfParameters() );
> OptimizerType::BoundValueType lowerBound(
> transform->GetNumberOfParameters() );
>
> boundSelect.Fill( 0 );
> upperBound.Fill( 0.0 );
> lowerBound.Fill( 0.0 );
>
> optimizer->SetBoundSelection( boundSelect );
> optimizer->SetUpperBound( upperBound );
> optimizer->SetLowerBound( lowerBound );
>
> optimizer->SetCostFunctionConvergenceFactor( 1e+7);
> optimizer->SetProjectedGradientTolerance(1e-4 );
> optimizer->SetMaximumNumberOfIterations( 500 );
> optimizer->SetMaximumNumberOfEvaluations( 500 );
> optimizer->SetMaximumNumberOfCorrections( 12 );
>
> CommandIterationUpdate::Pointer observer =
> CommandIterationUpdate::New();
> optimizer->AddObserver( itk::IterationEvent(), observer );
>
> try
> {
> registration->StartRegistration();
> }
> catch( itk::ExceptionObject & err )
> {
> std::cerr << "ExceptionObject caught !" << std::endl;
> std::cerr << err << std::endl;
> return -1;
> }
> ::
> ::
> --------------------------------------------------------------------------------------
>
>
>
>
>
>
>
>
>
>
>
>
>
>
>
>
>
>
>
>
> ------------------------------------------------------------------------
>
> _______________________________________________
> 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