[Insight-users] Flaw in ITK? : orientation of images
Koen Van Leemput
koen.vanleemput@hus.fi
Fri, 25 Oct 2002 11:37:32 +0300
Hi,
In ITK, itk::Image comes with an IndexToPhysicalTransform and its inverse=
=20
PhysicalToIndexTransform that defines how coordinates in the index space =
map=20
into coordinates in the "real-world" physical space and vice versa. In=20
itk::Image, these transforms are AffineTransforms that take into account =
the=20
spacing of the voxels (scaling) and the position of the origin (translati=
on).
In many medical image processing applications, however, there is a need f=
or=20
additional positioning information, defining how the patient is oriented =
with=20
respect to the image coordinate system. Examples include flipping and=20
swapping of image axes, gantry tilt in CT images, interfacing with other=20
image processing tools (SPM, 3D Slicer, ...), and proper handling of vari=
ous=20
medical image formats (DICOM, SPM, ...).=20
Because of these reasons, I decided to go for it and derive a class from=20
itk::Image that builds m_IndexToPhysicalTransform (and=20
m_PhysicalToIndexTransform by) by composing the standard scaling and=20
translation with an additional affine transform that defines orientation.=
=20
Apart from some convenience methods to manipulate this orientation, I=20
overloaded itk::Image::RebuildTransforms() and itk::Image::CopyInformatio=
n(),=20
and voila, it works: the orientation information travels nicely through t=
he=20
ITK pipeline.
BUT... things are messed up by what seems to me a design flaw in=20
itk::ImageFunction. Indeed, rather than using=20
itk::Image::TransformPhysicalPointToContinuousIndex() to convert a point =
in=20
physical space into an index in the image grid, itk::ImageFunction for so=
me=20
reason provides its own implementation=20
itk::ImageFunction::ConvertPointToContinuousIndex that noisily ignores=20
whatever PhysicalToIndexTransform the image at hand actually uses. This=20
messes up all the registraton metrics and itk::ResampleImageFilter, where=
the=20
following construction is typically used:
m_FixedImage->TransformIndexToPhysicalPoint( inputIndex, inputPoint )=
;
transformedPoint =3D m_Transform->TransformPoint( inputPoint );
if( m_Interpolator->IsInsideBuffer( transformedPoint ) )
{
movingValue =3D m_Interpolator->Evaluate( transformedPoint );
[clip clip] =09=09=09
}
(Note that also itk::ImageFunction::IsInsideBuffer() doesn't handle thing=
s=20
properly here neither)
As I see it, itk::ImageFunction::IsInsideBuffer() and=20
itk::ImageFunction::ConvertPointToContinuousIndex() should use the=20
implementation provided by itk::Image. To avoid inefficiencies induced by=
=20
converting from physical space into index space twice, the following=20
construction should probably be considered:
m_FixedImage->TransformIndexToPhysicalPoint( inputIndex, inputPoint )=
;
transformedPoint =3D m_Transform->TransformPoint( inputPoint );
m_MovingImage->TransformPhysicalPointToContinuousIndex(
transformedPoint, transformedIndex );
if( m_Interpolator->IsInsideBuffer( transformedIndex ) )
{
movingValue =3D m_Interpolator->Evaluate( transformedIndex );
[clip clip] =09=09=09
}
More in general, is support in ITK for orientation planned any time in th=
e=20
future? After all, the fiddling around with itk::PermuteAxesImageFilters=20
itk::FlipImageFilters in the registration examples seems quite cumbersome=
=2E..
Cheers,
Koen