ITK  6.0.0
Insight Toolkit
Examples/RegistrationITKv4/DeformableRegistration10.cxx
/*=========================================================================
*
* Copyright NumFOCUS
*
* Licensed under the Apache License, Version 2.0 (the "License");
* you may not use this file except in compliance with the License.
* You may obtain a copy of the License at
*
* https://www.apache.org/licenses/LICENSE-2.0.txt
*
* Unless required by applicable law or agreed to in writing, software
* distributed under the License is distributed on an "AS IS" BASIS,
* WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
* See the License for the specific language governing permissions and
* limitations under the License.
*
*=========================================================================*/
#include "itkConfigure.h"
#ifndef ITK_USE_FFTWD
# error "This program needs FFTWD to work."
#endif
constexpr unsigned int Dimension = 3;
// The following section of code implements a Command observer
// that will monitor the evolution of the registration process.
//
class CommandIterationUpdate : public itk::Command
{
public:
using Self = CommandIterationUpdate;
itkNewMacro(CommandIterationUpdate);
protected:
CommandIterationUpdate() = default;
using InternalImageType = itk::Image<float, Dimension>;
using VectorPixelType = itk::Vector<float, Dimension>;
using DisplacementFieldType = itk::Image<VectorPixelType, Dimension>;
using RegistrationFilterType = itk::CurvatureRegistrationFilter<
InternalImageType,
InternalImageType,
DisplacementFieldType,
InternalImageType,
InternalImageType,
DisplacementFieldType>>;
public:
void
Execute(itk::Object * caller, const itk::EventObject & event) override
{
Execute((const itk::Object *)caller, event);
}
void
Execute(const itk::Object * object, const itk::EventObject & event) override
{
const auto * filter = static_cast<const RegistrationFilterType *>(object);
if (!(itk::IterationEvent().CheckEvent(&event)))
{
return;
}
std::cout << filter->GetMetric() << std::endl;
}
};
int
main(int argc, char * argv[])
{
if (argc < 4)
{
std::cerr << "Missing Parameters " << std::endl;
std::cerr << "Usage: " << argv[0];
std::cerr << " fixedImageFile movingImageFile ";
std::cerr << " outputImageFile " << std::endl;
return EXIT_FAILURE;
}
using PixelType = short;
using FixedImageType = itk::Image<PixelType, Dimension>;
using MovingImageType = itk::Image<PixelType, Dimension>;
using FixedImageReaderType = itk::ImageFileReader<FixedImageType>;
using MovingImageReaderType = itk::ImageFileReader<MovingImageType>;
auto fixedImageReader = FixedImageReaderType::New();
auto movingImageReader = MovingImageReaderType::New();
fixedImageReader->SetFileName(argv[1]);
movingImageReader->SetFileName(argv[2]);
using InternalPixelType = float;
using InternalImageType = itk::Image<InternalPixelType, Dimension>;
using FixedImageCasterType =
using MovingImageCasterType =
auto fixedImageCaster = FixedImageCasterType::New();
auto movingImageCaster = MovingImageCasterType::New();
fixedImageCaster->SetInput(fixedImageReader->GetOutput());
movingImageCaster->SetInput(movingImageReader->GetOutput());
using MatchingFilterType =
auto matcher = MatchingFilterType::New();
matcher->SetInput(movingImageCaster->GetOutput());
matcher->SetReferenceImage(fixedImageCaster->GetOutput());
matcher->SetNumberOfHistogramLevels(1024);
matcher->SetNumberOfMatchPoints(7);
matcher->ThresholdAtMeanIntensityOn();
using VectorPixelType = itk::Vector<float, Dimension>;
using DisplacementFieldType = itk::Image<VectorPixelType, Dimension>;
using RegistrationFilterType = itk::CurvatureRegistrationFilter<
InternalImageType,
InternalImageType,
DisplacementFieldType,
InternalImageType,
InternalImageType,
DisplacementFieldType>>;
auto filter = RegistrationFilterType::New();
filter->SetTimeStep(1);
filter->SetConstraintWeight(0.1);
auto observer = CommandIterationUpdate::New();
filter->AddObserver(itk::IterationEvent(), observer);
using MultiResRegistrationFilterType =
InternalImageType,
DisplacementFieldType>;
multires->SetRegistrationFilter(filter);
multires->SetNumberOfLevels(3);
multires->SetFixedImage(fixedImageCaster->GetOutput());
multires->SetMovingImage(matcher->GetOutput());
unsigned int nIterations[4] = { 10, 20, 50, 50 };
multires->SetNumberOfIterations(nIterations);
multires->Update();
using DisplacementFieldTransformType =
auto displacementTransform = DisplacementFieldTransformType::New();
displacementTransform->SetDisplacementField(multires->GetOutput());
using WarperType = itk::ResampleImageFilter<MovingImageType,
MovingImageType,
InternalPixelType,
InternalPixelType>;
using InterpolatorType =
auto warper = WarperType::New();
auto interpolator = InterpolatorType::New();
FixedImageType::Pointer fixedImage = fixedImageReader->GetOutput();
warper->SetInput(movingImageReader->GetOutput());
warper->SetInterpolator(interpolator);
warper->SetOutputSpacing(fixedImage->GetSpacing());
warper->SetOutputOrigin(fixedImage->GetOrigin());
warper->SetOutputDirection(fixedImage->GetDirection());
warper->SetTransform(displacementTransform);
// Write warped image out to file
using OutputPixelType = unsigned short;
using OutputImageType = itk::Image<OutputPixelType, Dimension>;
using CastFilterType =
auto writer = WriterType::New();
auto caster = CastFilterType::New();
writer->SetFileName(argv[3]);
caster->SetInput(warper->GetOutput());
writer->SetInput(caster->GetOutput());
writer->Update();
if (argc > 4) // if a fourth line argument has been provided...
{
auto fieldWriter = FieldWriterType::New();
fieldWriter->SetFileName(argv[4]);
fieldWriter->SetInput(multires->GetOutput());
try
{
fieldWriter->Update();
}
catch (const itk::ExceptionObject & e)
{
e.Print(std::cerr);
}
}
return EXIT_SUCCESS;
}
itkFastSymmetricForcesDemonsRegistrationFunction.h
Pointer
SmartPointer< Self > Pointer
Definition: itkAddImageFilter.h:93
itk::CastImageFilter
Casts input pixels to output pixel type.
Definition: itkCastImageFilter.h:100
itkDisplacementFieldTransform.h
itkLinearInterpolateImageFunction.h
itk::HistogramMatchingImageFilter
Normalize the grayscale values for a source image by matching the shape of the source image histogram...
Definition: itkHistogramMatchingImageFilter.h:75
itkHistogramMatchingImageFilter.h
itkMultiResolutionPDEDeformableRegistration.h
itk::Vector
A templated class holding a n-Dimensional vector.
Definition: itkVector.h:62
itkImageFileReader.h
itk::FastSymmetricForcesDemonsRegistrationFunction
Definition: itkFastSymmetricForcesDemonsRegistrationFunction.h:47
itk::SmartPointer
Implements transparent reference counting.
Definition: itkSmartPointer.h:51
itkCastImageFilter.h
itk::DisplacementFieldTransform
Provides local/dense/high-dimensionality transformation via a a displacement field.
Definition: itkDisplacementFieldTransform.h:87
itk::ImageFileReader
Data source that reads image data from a single file.
Definition: itkImageFileReader.h:75
itk::LinearInterpolateImageFunction
Linearly interpolate an image at specified positions.
Definition: itkLinearInterpolateImageFunction.h:51
itk::Command
Superclass for callback/observer methods.
Definition: itkCommand.h:45
itk::ImageFileWriter
Writes image data to a single file.
Definition: itkImageFileWriter.h:90
itk::CurvatureRegistrationFilter
Deformably register two images using the fast curvature algorithm.
Definition: itkCurvatureRegistrationFilter.h:101
itk::Command
class ITK_FORWARD_EXPORT Command
Definition: itkObject.h:42
itk::Command::Execute
virtual void Execute(Object *caller, const EventObject &event)=0
itkImageFileWriter.h
itkCurvatureRegistrationFilter.h
itk::MultiResolutionPDEDeformableRegistration
Framework for performing multi-resolution PDE deformable registration.
Definition: itkMultiResolutionPDEDeformableRegistration.h:88
itk::ResampleImageFilter
Resample an image via a coordinate transform.
Definition: itkResampleImageFilter.h:90
itk::Object
Base class for most ITK classes.
Definition: itkObject.h:61
itk::Math::e
static constexpr double e
Definition: itkMath.h:56
itk::Image
Templated n-dimensional image class.
Definition: itkImage.h:88
itk::EventObject
Abstraction of the Events used to communicating among filters and with GUIs.
Definition: itkEventObject.h:58
New
static Pointer New()
AddImageFilter
Definition: itkAddImageFilter.h:81
itk::GTest::TypedefsAndConstructors::Dimension2::Dimension
constexpr unsigned int Dimension
Definition: itkGTestTypedefsAndConstructors.h:44
Superclass
BinaryGeneratorImageFilter< TInputImage1, TInputImage2, TOutputImage > Superclass
Definition: itkAddImageFilter.h:90