CSci 6971: Image Registration Lecture 10: Registration Components February 13, 2004 Prof. Chuck Stewart, RPI Dr.

Download Report

Transcript CSci 6971: Image Registration Lecture 10: Registration Components February 13, 2004 Prof. Chuck Stewart, RPI Dr.

CSci 6971: Image Registration
Lecture 10: Registration Components
February 13, 2004
Prof. Chuck Stewart, RPI
Dr. Luis Ibanez, Kitware
Registration Components
Basic
Registration Framework
Image Registration
Lecture 10
2
Image Registration Framework
Fixed
Image
Metric
Moving
Image
Image Registration
Interpolator
Optimizer
Transform
Parameters
Lecture 10
3
Image Registration
#include ”itkImageRegistrationMethod.h”
#include ”itkTranslationTransform.h”
#include ”itkMeanSquaresImageToImageMetric.h”
#include ”itkLinearInterpolateImageFunction.h”
#include ”itkRegularStepGradientDescentOptimizer.h”
#include ”itkImage.h”
#include ”itkImageFileReader.h”
#include ”itkImageFileWriter.h”
#include ”itkResampleImageFilter.h”
Image Registration
Lecture 10
4
Image Registration
const unsigned int Dimension = 2;
typedef unsigned char PixelType;
typedef itk::Image< PixelType , Dimension >
typedef itk::Image< PixelType , Dimension >
FixedImageType;
MovingImageType;
typedef itk::TranslationTransform< double, Dimension >
typedef itk::RegularStepGradientDescentOptimizer
typedef itk::LinearInterpolateImageFunction<
MovingImageType , double >
TransformType;
OptimizerType;
InterpolatorType;
typedef itk::MeanSquaresImageToImageMetric<
FixedImageType , MovingImageType >
MetricType;
typedef itk::ImageRegistrationMethod<
FixedImageType , MovingImageType >
RegistrationType;
Image Registration
Lecture 10
5
Image Registration
TransformType::Pointer
OptimizerType::Pointer
InterpolatorType::Pointer
MetricType::Pointer
RegistrationType::Pointer
transform
optimizer
interpolator
metric
registrator
registrator->SetTransform(
registrator->SetOptimizer(
registrator->SetInterpolator(
registrator->SetMetric(
= TransformType::New();
= OptimizerType::New();
= InterpolatorType::New();
= MetricType::New();
= RegistrationType::New();
transform );
optimizer );
interpolator );
metric
);
registrator->SetFixedImage( fixedImageReader->GetOutput() );
registrator->SetMovingImage( movingImageReader->GetOutput() );
Image Registration
Lecture 10
6
Image Registration
registrator->SetFixedImageRegion(
fixedImageReader->GetOutput()->GetBufferedRegion() );
typedef RegistrationType::ParametersType
ParametersType;
transform->SetIdentity();
registrator->SetInitialTransformParameters(
transform->GetParameters() );
optimizer->SetMaximumStepLength( 4.00 );
optimizer->SetMinimumStepLength( 0.01 );
optimizer->SetNumberOfIterations( 100 );
optimizer->MaximizeOff();
Image Registration
Lecture 10
7
Image Registration
try
{
registrator->StartRegistration ();
}
catch( itk::ExceptionObject & excp )
{
std::cerr << “Error in registration” << std::endl;
std::cerr << excp << std::endl;
}
transform->SetParameters(
registrator->GetLastTransformParameters() );
Image Registration
Lecture 10
8
Image Registration Framework
Fixed
Image
Metric
Optimizer
Moving
Image
Interpolator
Parameters
Transform
Resampler
Moving
Registered
Image Registration
Lecture 10
9
Image Registration
typedef itk::ResampleImageFilter<
FixedImageType , MovingImageType >
ResamplerType;
ResamplerType ::Pointer resampler = ResamplerType::New();
resampler->SetTransform ( transform );
resampler->SetInput( movingImageReader->GetOutput() );
FixedImageType ::Pointer fixedImage = fixedImageReader->GetOutput();
resampler->SetOrigin( fixedImage->GetOrigin() );
resampler->SetSpacing( fixedImage->GetSpacing() );
resampler->SetSize(
fixedImage->GetLargestPossibleRegion()->GetSize() );
resampler->Update();
Image Registration
Lecture 10
10
Image Registration
Tracking the
Registration Process
Image Registration
Lecture 10
11
Observing Registration
#include ”itkCommand.h”
class CommandIteration : public itk::Command {
public:
typedef CommandIteration Self;
typedef itk::Command
SuperClass;
typedef itk::SmartPointer< Self > Pointer;
itkNewMacro( Self );
protected:
CommandIteration() {};
public:
typedef itk::RegularStepGradientDescentOptimizer OptimizerType;
typedef
const OptimizerType
* OptimizerPointer;
Image Registration
Lecture 10
12
Observing Registration
void Execute( itk::Object * caller, const itk::EventObject & event )
{
this-> Execute( (const itk::Object *) caller, event );
}
void Execute( const itk::Object * caller, const itk::EventObject & event )
{
OptimizerPointer optimizer =
dynamic_cast< OptimizerPointer >( caller );
if( typeid( event ) == typeid( itk::IterationEvent ) )
{
std::cout << optimizer->GetCurrentIteration() << “ : “;
std::cout << optimizer->GetValue() << “ : “;
std::cout << optimizer->GetCurrentPosition() << std::endl;
}
}
Image Registration
Lecture 10
13
Observing Registration
CommandIteration::Pointer observer = CommandIteration::New();
optimizer->AddObserver( itk::IterationEvent(), observer )
try
{
registrator->StartRegistration ();
}
catch( itk::ExceptionObject & excp )
{
std::cerr << “Error in registration” << std::endl;
std::cerr << excp << std::endl;
}
Image Registration
Lecture 10
14
Registration dirty little secrets
Rotation – Translation
Parameter Scaling
Image Registration
Lecture 10
15
Registration dirty little secrets
Fixed Image
Image Registration
Moving Image
Lecture 10
Registered
Moving Image
16
Gradient Descent Optimizer
f( x , y )
Image Registration
∆
G( x , y ) =
S = L ∙ G( x , y )
f( x , y )
Lecture 10
17
Gradient Descent Optimizer
f( x , y )
Image Registration
∆
G( x , y ) =
S = L ∙ G( x , y )
f( x , y )
Lecture 10
18
Gradient Descent Optimizer
f( x , y )
Image Registration
∆
G( x , y ) =
S = L ∙ G( x , y )
f( x , y )
Lecture 10
19
Gradient Descent Optimizer
f( x , y )
Image Registration
∆
G( x , y ) =
S = L ∙ G( x , y )
f( x , y )
Lecture 10
20
Gradient Descent Optimizer
f( x , y )
Image Registration
∆
G( x , y ) =
S = L ∙ G( x , y )
f( x , y )
Lecture 10
21
Gradient Descent Optimizer
f( x , y )
How about a factor
Image Registration
∆
G( x , y ) =
1:100
?
S = L ∙ G( x , y )
f( x , y )
Lecture 10
22
Radians & Millimeters
y
1.57 rad
0.0 rad
x
180 mm
Image Registration
Lecture 10
23
Radians & Millimeters
1.57 rad
y
1:100
180 mm
1.57 rad
0.0 rad
x
180 mm
Image Registration
Lecture 10
24
Other Transforms
Centered Rigid 2D
Image Registration
Lecture 10
25
Rotation around the Center
y
y
Transform
10o
x
x
Fixed
Image Registration
Moving
Lecture 10
26
Centered Rigid 2D Transform
#include ”itkImageRegistrationMethod.h”
#include ”itkCenteredRigid2DTransform.h”
#include ”itkMeanSquaresImageToImageMetric.h”
#include ”itkLinearInterpolateImageFunction.h”
#include ”itkRegularStepGradientDescentOptimizer.h”
#include ”itkCenteredTransformInitializer.h”
#include ”itkImage.h”
#include ”itkImageFileReader.h”
#include ”itkImageFileWriter.h”
#include ”itkResampleImageFilter.h”
Image Registration
Lecture 10
27
Centered Rigid 2D Transform
const unsigned int Dimension = 2;
typedef unsigned char PixelType;
typedef itk::Image< PixelType , Dimension >
typedef itk::Image< PixelType , Dimension >
FixedImageType;
MovingImageType;
typedef itk::CenteredRigid2DTransform< double >
TransformType;
typedef itk:: CenteredTransformInitializer<
TransformType ,
FixedImageType ,
MovingImageType
>
InitializerType;
Image Registration
Lecture 10
28
Centered Rigid 2D Transform
TransformType::Pointer
InitializerType::Pointer
OptimizerType::Pointer
InterpolatorType::Pointer
MetricType::Pointer
RegistrationType::Pointer
transform
initializer
optimizer
interpolator
metric
registrator
registrator->SetTransform(
registrator->SetOptimizer(
registrator->SetInterpolator(
registrator->SetMetric(
= TransformType::New();
= InitializerType::New();
= OptimizerType::New();
= InterpolatorType::New();
= MetricType::New();
= RegistrationType::New();
transform );
optimizer );
interpolator );
metric
);
registrator->SetFixedImage( fixedImageReader->GetOutput() );
registrator->SetMovingImage( movingImageReader->GetOutput() );
Image Registration
Lecture 10
29
Centered Transform Initializer
Geometry On
y
y
Transform
x
x
Fixed Image
Image Registration
Moving Image
Lecture 10
30
Centered Transform Initializer
Moments On
y
y
Transform
x
x
Fixed Image
Moving Image
Xm = Sum( x ∙ I ) / Sum( I )
Ym = Sum( y ∙ I ) / Sum( I )
Image Registration
Lecture 10
31
Centered Rigid 2D Transform
registrator->SetFixedImageRegion(
fixedImageReader->GetOutput()->GetBufferedRegion() );
initializer->SetTransform ( transform );
initializer->SetFixedImage( fixedImageReader->GetOutput() );
initializer->SetMovingImage( movingImageReader->GetOutput() );
initializer->MomentsOn();
initializer->InitializeTransform();
registrator->SetInitialTransformParameters(
transform->GetParameters() );
Image Registration
Lecture 10
32
Centered Rigid 2D Transform
Array of
y
Parameters
x
Image Registration
Lecture 10
P[ 0 ]
Angle
P[ 1 ]
Center X
P[ 2 ]
Center Y
P[ 3 ]
Trans X
P[ 4 ]
Trans Y
33
Centered Rigid 2D Transform
typedef OptimizerType::ScaleType
OptimizerScalesType;
OptimizerScalesType optimizerScales(
optimizer->SetMaximumStepLength() );
const double translationScale = 1.0 / 1000.0 ;
optimizerScales[ 0 ] = 1.0;
optimizerScales[ 1 ] = translationScale;
optimizerScales[ 2 ] = translationScale;
optimizerScales[ 3 ] = translationScale;
optimizerScales[ 4 ] = translationScale;
optimizer->SetScales( optimizerScales );
Image Registration
Lecture 10
34
Centered Rigid 2D Transform
try
{
registrator->StartRegistration ();
}
catch( itk::ExceptionObject & excp )
{
std::cerr << “Error in registration” << std::endl;
std::cerr << excp << std::endl;
}
transform->SetParameters(
registrator->GetLastTransformParameters() );
Image Registration
Lecture 10
35
Image Registration Framework
Fixed
Image
Metric
Optimizer
Moving
Image
Interpolator
Parameters
Transform
Resampler
Moving
Registered
Image Registration
Lecture 10
36
Centered Rigid 2D Transform
typedef itk::ResampleImageFilter<
FixedImageType , MovingImageType >
ResamplerType;
ResamplerType ::Pointer resampler = ResamplerType::New();
resampler->SetTransform ( transform );
resampler->SetInput( movingImageReader->GetOutput() );
FixedImageType ::Pointer fixedImage = fixedImageReader->GetOutput();
resampler->SetOrigin( fixedImage->GetOrigin() );
resampler->SetSpacing( fixedImage->GetSpacing() );
resampler->SetSize(
fixedImage->GetLargestPossibleRegion()->GetSize() );
resampler->Update();
Image Registration
Lecture 10
37
Other Transforms
Affine 2D
Image Registration
Lecture 10
38
Affine Transforms
Image Registration
M11
M12
T1
M21
M22
T2
Lecture 10
39
Affine Transforms
Qx
M11
M12
∙
=
Qy
Image Registration
Px
M21
M22
Lecture 10
T1
+
Py
T2
40
Affine Transform – Shearing
Qx
1
a
∙
=
Qy
Image Registration
Px
0
1
Lecture 10
T1
+
Py
T2
41
Shearing Component
y
a
( x + ay , y )
1
(x,y)
Image Registration
Lecture 10
x
42
Shearing Component
y’
y
x’
x
Image Registration
Lecture 10
43
Affine Transform – Shearing
Qx
1
0
∙
=
Qy
Image Registration
Px
b
1
Lecture 10
T1
+
Py
T2
44
Shearing Component
y
( x , y + bx )
b
(x,y)
1
x
Image Registration
Lecture 10
45
Shearing Component
y
y’
x’
x
Image Registration
Lecture 10
46
Coefficients - Orders of Magnitude
Qx
1.0 0.1
∙
=
Qy
Image Registration
Px
0.1
1.0
Lecture 10
100
+
Py
100
47
Affine Transform – Rotation
Qx
cosθ
-sinθ
∙
=
Qy
Image Registration
Px
sinθ
cosθ
Lecture 10
T1
+
Py
T2
48
Affine Transform – Scaling
Qx
Sx
0
∙
=
Qy
Image Registration
Px
0
Sy
Lecture 10
T1
+
Py
T2
49
Affine Transform - Parameterization
[
Image Registration
P0
M11
M12
T1
M21
M22
T2
P1
P2
P3
Lecture 10
P4
P5
]
50
Parameter Scaling
[
P0
M11
M12
T1
M21
M22
T2
P1
P2
P3
P4
P5
]
range [
1.0 0.1 0.1 1.0
100 100
scale [
1.0 10
0.01 0.01 ]
Image Registration
10 1.0
Lecture 10
]
51
Centered Affine Transform
#include ”itkImageRegistrationMethod.h”
#include ”itkCenteredAffineTransform.h”
#include ”itkMeanSquaresImageToImageMetric.h”
#include ”itkLinearInterpolateImageFunction.h”
#include ”itkRegularStepGradientDescentOptimizer.h”
#include ”itkCenteredTransformInitializer.h”
#include ”itkImage.h”
#include ”itkImageFileReader.h”
#include ”itkImageFileWriter.h”
#include ”itkResampleImageFilter.h”
Image Registration
Lecture 10
52
Centered Affine Transform
const unsigned int Dimension = 2;
typedef unsigned char PixelType;
typedef itk::Image< PixelType , Dimension >
typedef itk::Image< PixelType , Dimension >
typedef itk::CenteredAffineTransform< double >
FixedImageType;
MovingImageType;
TransformType;
typedef itk:: CenteredTransformInitializer<
TransformType ,
FixedImageType ,
MovingImageType
>
InitializerType;
Image Registration
Lecture 10
53
Centered Affine Transform
TransformType::Pointer
InitializerType::Pointer
OptimizerType::Pointer
InterpolatorType::Pointer
MetricType::Pointer
RegistrationType::Pointer
transform
initializer
optimizer
interpolator
metric
registrator
registrator->SetTransform(
registrator->SetOptimizer(
registrator->SetInterpolator(
registrator->SetMetric(
= TransformType::New();
= InitializerType::New();
= OptimizerType::New();
= InterpolatorType::New();
= MetricType::New();
= RegistrationType::New();
transform );
optimizer );
interpolator );
metric
);
registrator->SetFixedImage( fixedImageReader->GetOutput() );
registrator->SetMovingImage( movingImageReader->GetOutput() );
Image Registration
Lecture 10
54
Centered Transform Initializer
Geometry On
y
y
Transform
x
x
Fixed Image
Image Registration
Moving Image
Lecture 10
55
Centered Transform Initializer
Moments On
y
y
Transform
x
x
Fixed Image
Moving Image
Xm = Sum( x ∙ I ) / Sum( I )
Ym = Sum( y ∙ I ) / Sum( I )
Image Registration
Lecture 10
56
Centered Affine Transform
registrator->SetFixedImageRegion(
fixedImageReader->GetOutput()->GetBufferedRegion() );
initializer->SetTransform ( transform );
initializer->SetFixedImage( fixedImageReader->GetOutput() );
initializer->SetMovingImage( movingImageReader->GetOutput() );
initializer->MomentsOn();
initializer->InitializeTransform();
registrator->SetInitialTransformParameters(
transform->GetParameters() );
Image Registration
Lecture 10
57
Affine Transform - Parameterization
M11
M12
T1
M21
M22
T2
P[0] P[1] P[2] P[3]
Image Registration
Lecture 10
P[4]
P[5]
58
Centered Affine Transform
typedef OptimizerType::ScaleType
OptimizerScalesType;
OptimizerScalesType optimizerScales(
optimizer->GetNumberOfParameters() );
optimizerScales[ 0 ] =
optimizerScales[ 1 ] =
optimizerScales[ 2 ] =
optimizerScales[ 3 ] =
optimizerScales[ 4 ] =
optimizerScales[ 5 ] =
1.0;
10.0;
10.0;
1.0;
0.01;
0.01;
optimizer->SetScales( optimizerScales );
Image Registration
Lecture 10
59
Centered Rigid 2D Transform
try
{
registrator->StartRegistration ();
}
catch( itk::ExceptionObject & excp )
{
std::cerr << “Error in registration” << std::endl;
std::cerr << excp << std::endl;
}
transform->SetParameters(
registrator->GetLastTransformParameters() );
Image Registration
Lecture 10
60
Image Registration Framework
Fixed
Image
Metric
Optimizer
Moving
Image
Interpolator
Parameters
Transform
Resampler
Moving
Registered
Image Registration
Lecture 10
61
Final Resampling
typedef itk::ResampleImageFilter<
FixedImageType , MovingImageType >
ResamplerType;
ResamplerType ::Pointer resampler = ResamplerType::New();
resampler->SetTransform ( transform );
resampler->SetInput( movingImageReader->GetOutput() );
FixedImageType ::Pointer fixedImage = fixedImageReader->GetOutput();
resampler->SetOrigin( fixedImage->GetOrigin() );
resampler->SetSpacing( fixedImage->GetSpacing() );
resampler->SetSize(
fixedImage->GetLargestPossibleRegion()->GetSize() );
resampler->Update();
Image Registration
Lecture 10
62
End
Enjoy ITK !
Image Registration
Lecture 10
63