Hi, All :<br><br>I want to combine Affine transform and Mattes MI metric into my registration process.<br>The process is&nbsp; fine to&nbsp; work but the result of registration is incorrect.<br>I think this is caused by some coefficients setting that is not suitable for my application,
<br>but how to determine a applicable coefficients that is very suffering for me.<br>So, I wish any experienced friends can give me some indication that including how to set:<br>OptimizerScales (specially translation Scale) &lt;---How to set ???
<br>MaximumStepLength &lt;---How to set ???<br>
MinimumStepLength &lt;---How to set ???<br>NumberOfHistogramBins (now I set 20~50 in typically)<br>NumberOfSpatialSamples (now I set 20% of image pixels )<br><br><br>My code is posting as below:<br>Regards.<br><br>class CommandIterationUpdate : public itk::Command 
<br>{<br>public:<br>&nbsp; typedef&nbsp; CommandIterationUpdate&nbsp;&nbsp; Self;<br>&nbsp; typedef&nbsp; itk::Command&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp; Superclass;<br>&nbsp; typedef itk::SmartPointer&lt;Self&gt;&nbsp; Pointer;<br>&nbsp; itkNewMacro( Self );<br>protected:<br>&nbsp; CommandIterationUpdate() {};
<br>public:<br>&nbsp; typedef itk::RegularStepGradientDescentOptimizer&nbsp;&nbsp;&nbsp;&nbsp; OptimizerType;<br>&nbsp; typedef&nbsp;&nbsp; const OptimizerType&nbsp;&nbsp; *&nbsp;&nbsp;&nbsp; OptimizerPointer;<br><br>&nbsp; void Execute(itk::Object *caller, const itk::EventObject &amp; event)
<br>&nbsp; {<br>&nbsp;&nbsp;&nbsp; Execute( (const itk::Object *)caller, event);<br>&nbsp; }<br><br>&nbsp; void Execute(const itk::Object * object, const itk::EventObject &amp; event)<br>&nbsp; {<br>&nbsp;&nbsp;&nbsp; OptimizerPointer optimizer = <br>&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp; dynamic_cast&lt; OptimizerPointer &gt;( object );
<br>&nbsp;&nbsp;&nbsp; if( ! itk::IterationEvent().CheckEvent( &amp;event ) )<br>&nbsp;&nbsp;&nbsp;&nbsp;&nbsp; {<br>&nbsp;&nbsp;&nbsp;&nbsp;&nbsp; return;<br>&nbsp;&nbsp;&nbsp;&nbsp;&nbsp; }<br>&nbsp;&nbsp;&nbsp;&nbsp;&nbsp; std::cout &lt;&lt; optimizer-&gt;GetCurrentIteration() &lt;&lt; &quot;&nbsp;&nbsp; &quot;;<br>&nbsp;&nbsp;&nbsp;&nbsp;&nbsp; std::cout &lt;&lt; optimizer-&gt;GetValue() &lt;&lt; &quot;&nbsp;&nbsp; &quot;;
<br>&nbsp;&nbsp;&nbsp;&nbsp;&nbsp; std::cout &lt;&lt; optimizer-&gt;GetCurrentPosition();<br>&nbsp;&nbsp;&nbsp;&nbsp; <br>&nbsp;&nbsp;&nbsp;&nbsp;&nbsp; // Print the angle for the trace plot<br>&nbsp;&nbsp;&nbsp;&nbsp;&nbsp; vnl_matrix&lt;double&gt; p(2, 2);<br>&nbsp;&nbsp;&nbsp;&nbsp;&nbsp; p[0][0] = (double) optimizer-&gt;GetCurrentPosition()[0];
<br>&nbsp;&nbsp;&nbsp;&nbsp;&nbsp; p[0][1] = (double) optimizer-&gt;GetCurrentPosition()[1];<br>&nbsp;&nbsp;&nbsp;&nbsp;&nbsp; p[1][0] = (double) optimizer-&gt;GetCurrentPosition()[2];<br>&nbsp;&nbsp;&nbsp;&nbsp;&nbsp; p[1][1] = (double) optimizer-&gt;GetCurrentPosition()[3];<br>&nbsp;&nbsp;&nbsp;&nbsp;&nbsp; vnl_svd&lt;double&gt; svd(p);
<br>&nbsp;&nbsp;&nbsp;&nbsp;&nbsp; vnl_matrix&lt;double&gt; r(2, 2);<br>&nbsp;&nbsp;&nbsp;&nbsp;&nbsp; r = svd.U() * vnl_transpose(svd.V());<br>&nbsp;&nbsp;&nbsp;&nbsp;&nbsp; double angle = asin(r[1][0]);<br>&nbsp;&nbsp;&nbsp;&nbsp;&nbsp; std::cout &lt;&lt; &quot; AffineAngle: &quot; &lt;&lt; angle * 45.0 / atan(1.0) &lt;&lt; std::endl;
<br>&nbsp;&nbsp;&nbsp; }<br>};<br><br>int main( int argc, char *argv[] )<br>{<br>&nbsp; if( argc &lt; 4 )<br>&nbsp;&nbsp;&nbsp; {<br>&nbsp;&nbsp;&nbsp; std::cerr &lt;&lt; &quot;Missing Parameters &quot; &lt;&lt; std::endl;<br>&nbsp;&nbsp;&nbsp; std::cerr &lt;&lt; &quot;Usage: &quot; &lt;&lt; argv[0];
<br>&nbsp;&nbsp;&nbsp; std::cerr &lt;&lt; &quot;&nbsp;&nbsp; fixedImageFile&nbsp; movingImageFile outputImagefile&quot; &lt;&lt; std::endl;<br>&nbsp;&nbsp;&nbsp; return 1;<br>&nbsp;&nbsp;&nbsp; }<br><br>&nbsp; const&nbsp;&nbsp;&nbsp; unsigned int&nbsp;&nbsp;&nbsp; Dimension = 2;<br>&nbsp; typedef&nbsp; float&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp; PixelType;
<br><br>&nbsp; typedef itk::Image&lt; PixelType, Dimension &gt;&nbsp; FixedImageType;<br>&nbsp; typedef itk::Image&lt; PixelType, Dimension &gt;&nbsp; MovingImageType;<br><br>&nbsp; typedef itk::AffineTransform&lt; <br>&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp; double, 
<br>&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp; Dimension&nbsp; &gt;&nbsp;&nbsp;&nbsp;&nbsp; TransformType;<br><br>&nbsp; typedef itk::RegularStepGradientDescentOptimizer&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp; OptimizerType;<br><br>&nbsp; typedef itk::MattesMutualInformationImageToImageMetric&lt; <br>
&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp; FixedImageType, <br>&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp; MovingImageType &gt;&nbsp;&nbsp;&nbsp; MetricType;<br><br>&nbsp; typedef itk:: LinearInterpolateImageFunction&lt; <br>&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp; MovingImageType,
<br>&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp; double&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp; &gt;&nbsp;&nbsp;&nbsp; InterpolatorType;<br>&nbsp; typedef itk::ImageRegistrationMethod&lt; <br>&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp; FixedImageType, <br>&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp; MovingImageType &gt;&nbsp;&nbsp;&nbsp; RegistrationType;
<br>&nbsp; <br>&nbsp; MetricType::Pointer&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp; metric&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp; = MetricType::New();<br>&nbsp; OptimizerType::Pointer&nbsp;&nbsp;&nbsp;&nbsp;&nbsp; optimizer&nbsp;&nbsp;&nbsp;&nbsp; = OptimizerType::New();<br>&nbsp; InterpolatorType::Pointer&nbsp;&nbsp; interpolator&nbsp; = InterpolatorType::New();<br>
&nbsp; RegistrationType::Pointer&nbsp;&nbsp; registration&nbsp; = RegistrationType::New();<br><br>&nbsp; registration-&gt;SetMetric(&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp; metric&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp; );<br>&nbsp; registration-&gt;SetOptimizer(&nbsp;&nbsp;&nbsp;&nbsp; optimizer&nbsp;&nbsp;&nbsp;&nbsp; );<br>&nbsp; registration-&gt;SetInterpolator(&nbsp; interpolator&nbsp; );
<br><br>&nbsp; TransformType::Pointer&nbsp; transform = TransformType::New();<br>&nbsp; registration-&gt;SetTransform( transform );<br><br>&nbsp; typedef itk::ImageFileReader&lt; FixedImageType&nbsp; &gt; FixedImageReaderType;<br>&nbsp; typedef itk::ImageFileReader&lt; MovingImageType &gt; MovingImageReaderType;
<br>&nbsp; FixedImageReaderType::Pointer&nbsp; fixedImageReader&nbsp; = FixedImageReaderType::New();<br>&nbsp; MovingImageReaderType::Pointer movingImageReader = MovingImageReaderType::New();<br>&nbsp; fixedImageReader-&gt;SetFileName(&nbsp; argv[1] );
<br>&nbsp; movingImageReader-&gt;SetFileName( argv[2] );<br><br><br>&nbsp; registration-&gt;SetFixedImage(&nbsp;&nbsp;&nbsp; fixedImageReader-&gt;GetOutput()&nbsp;&nbsp;&nbsp; );<br>&nbsp; registration-&gt;SetMovingImage(&nbsp;&nbsp; movingImageReader-&gt;GetOutput()&nbsp;&nbsp; );<br>
&nbsp; fixedImageReader-&gt;Update();<br><br>&nbsp; registration-&gt;SetFixedImageRegion( <br>&nbsp;&nbsp;&nbsp;&nbsp; fixedImageReader-&gt;GetOutput()-&gt;GetBufferedRegion() );<br><br>&nbsp; typedef itk::CenteredTransformInitializer&lt; <br>&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp; TransformType, 
<br>&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp; FixedImageType, <br>&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp; MovingImageType &gt;&nbsp; TransformInitializerType;<br>&nbsp; TransformInitializerType::Pointer initializer = TransformInitializerType::New();
<br>&nbsp; initializer-&gt;SetTransform(&nbsp;&nbsp; transform );<br>&nbsp; initializer-&gt;SetFixedImage(&nbsp; fixedImageReader-&gt;GetOutput() );<br>&nbsp; initializer-&gt;SetMovingImage( movingImageReader-&gt;GetOutput() );<br>&nbsp; initializer-&gt;MomentsOn();
<br>&nbsp; initializer-&gt;InitializeTransform();<br>&nbsp; <br>&nbsp; registration-&gt;SetInitialTransformParameters( <br>&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp; transform-&gt;GetParameters() );<br><br>&nbsp; typedef OptimizerType::ScalesType&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp; OptimizerScalesType;
<br>&nbsp; OptimizerScalesType optimizerScales( transform-&gt;GetNumberOfParameters() );<br><br>//-----------------------------------------------------------------------------&nbsp;&nbsp;&nbsp; // How to determine these ?<br>&nbsp; metric-&gt;SetNumberOfHistogramBins( 50 );
<br>&nbsp; metric-&gt;SetNumberOfSpatialSamples( 10000 );<br>&nbsp; metric-&gt;ReinitializeSeed( 76926294 );<br><br>&nbsp; double translationScale = 1.0 / 1000.0;<br>&nbsp; optimizerScales[0] =&nbsp; 1.0;<br>&nbsp; optimizerScales[1] =&nbsp; 1.0;<br>&nbsp; optimizerScales[2] =&nbsp; 
1.0;<br>&nbsp; optimizerScales[3] =&nbsp; 1.0;<br>&nbsp; optimizerScales[4] =&nbsp; translationScale;<br>&nbsp; optimizerScales[5] =&nbsp; translationScale;<br>&nbsp; optimizer-&gt;SetScales( optimizerScales );<br><br>&nbsp; double steplength = 1;&nbsp;&nbsp;&nbsp; //0.1<br>&nbsp; unsigned int maxNumberOfIterations = 200;
<br>&nbsp; optimizer-&gt;SetMaximumStepLength( steplength ); <br>&nbsp; optimizer-&gt;SetMinimumStepLength( 0.0001 );&nbsp;&nbsp;&nbsp; //0.0001<br>&nbsp; optimizer-&gt;SetNumberOfIterations( maxNumberOfIterations );<br>//-----------------------------------------------------------------------------
<br><br><br>&nbsp; optimizer-&gt;MinimizeOn();<br><br>&nbsp; <br>&nbsp; CommandIterationUpdate::Pointer observer = CommandIterationUpdate::New();<br>&nbsp; optimizer-&gt;AddObserver( itk::IterationEvent(), observer );<br><br>&nbsp; try <br>&nbsp;&nbsp;&nbsp; { 
<br>&nbsp;&nbsp;&nbsp; registration-&gt;StartRegistration(); <br>&nbsp;&nbsp;&nbsp; } <br>&nbsp; catch( itk::ExceptionObject &amp; err ) <br>&nbsp;&nbsp;&nbsp; { <br>&nbsp;&nbsp;&nbsp; std::cerr &lt;&lt; &quot;ExceptionObject caught !&quot; &lt;&lt; std::endl; <br>&nbsp;&nbsp;&nbsp; std::cerr &lt;&lt; err &lt;&lt; std::endl; 
<br>&nbsp;&nbsp;&nbsp; return -1;<br>&nbsp;&nbsp;&nbsp; }&nbsp;&nbsp; <br><br><br><br>//-----------------------------------------------------------<br>&nbsp; OptimizerType::ParametersType finalParameters = <br>&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp; registration-&gt;GetLastTransformParameters();
<br><br>&nbsp; const double finalRotationCenterX = transform-&gt;GetCenter()[0];<br>&nbsp; const double finalRotationCenterY = transform-&gt;GetCenter()[1];<br>&nbsp; const double finalTranslationX&nbsp;&nbsp;&nbsp; = finalParameters[4];<br>&nbsp; const double finalTranslationY&nbsp;&nbsp;&nbsp; = finalParameters[5];
<br><br>&nbsp; const unsigned int numberOfIterations = optimizer-&gt;GetCurrentIteration();<br>&nbsp; <br>&nbsp; const double bestValue = optimizer-&gt;GetValue();<br>&nbsp; std::cout &lt;&lt; &quot;Result = &quot; &lt;&lt; std::endl;<br>&nbsp; std::cout &lt;&lt; &quot; Center X&nbsp;&nbsp;&nbsp;&nbsp;&nbsp; = &quot; &lt;&lt; finalRotationCenterX&nbsp; &lt;&lt; std::endl;
<br>&nbsp; std::cout &lt;&lt; &quot; Center Y&nbsp;&nbsp;&nbsp;&nbsp;&nbsp; = &quot; &lt;&lt; finalRotationCenterY&nbsp; &lt;&lt; std::endl;<br>&nbsp; std::cout &lt;&lt; &quot; Translation X = &quot; &lt;&lt; finalTranslationX&nbsp; &lt;&lt; std::endl;<br>&nbsp; std::cout &lt;&lt; &quot; Translation Y = &quot; &lt;&lt; finalTranslationY&nbsp; &lt;&lt; std::endl;
<br>&nbsp; std::cout &lt;&lt; &quot; Iterations&nbsp;&nbsp;&nbsp; = &quot; &lt;&lt; numberOfIterations &lt;&lt; std::endl;<br>&nbsp; std::cout &lt;&lt; &quot; Metric value&nbsp; = &quot; &lt;&lt; bestValue&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp; &lt;&lt; std::endl;<br>&nbsp; <br>