Hello,
I was trying to implement deform registration from example7 : https://itk.org/Doxygen/html/Examples_2RegistrationITKv4_2DeformableRegistration7_8cxx-example.html
//Deformable registration using BSpline Transform
#include "itkImage.h"
#include "itkImageFileReader.h"
#include "itkImageFileWriter.h"
#include "itkImageRegistrationMethodv4.h"
#include "itkMeanSquaresImageToImageMetricv4.h"
#include "itkBSplineTransform.h"
#include "itkLBFGSBOptimizerv4.h"
#include "itkResampleImageFilter.h"
#include "itkCastImageFilter.h"
#include "itkSquaredDifferenceImageFilter.h"
#include "itkIdentityTransform.h"
#include "itkBSplineTransformInitializer.h"
#include "itkTransformToDisplacementFieldFilter.h"
//A command observer used to monitor the evolution of the registration process.
#include "itkCommand.h"
class CommandIterationUpdate : public itk::Command
{
public:
using Self = CommandIterationUpdate;
using Superclass = itk::Command;
using Pointer = itk::SmartPointer<Self>;
itkNewMacro(Self);
protected:
CommandIterationUpdate() = default;
public:
using OptimizerType = itk::LBFGSBOptimizerv4;
using OptimizerPointer = const OptimizerType *;
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
{
auto optimizer = static_cast<OptimizerPointer>(object);
if (!(itk::IterationEvent().CheckEvent(&event)))
{
return;
}
std::cout << optimizer->GetCurrentIteration() << " ";
std::cout << optimizer->GetCurrentMetricValue() << " ";
std::cout << optimizer->GetInfinityNormOfProjectedGradient() << std::endl;
}
};
int main(int argc,char * argv[]){
const unsigned int nDims = 3 ;
const unsigned int SpaceDimension = nDims;
const unsigned int SplineOrder = 3;
unsigned int numberOfGridNodesInOneDimension = 5; //ex uses 8
typedef itk::Image<double, nDims > ImageType;
//using PixelType = float;
typedef itk::ImageFileReader<ImageType>readerType; //fixed
readerType::Pointer reader1 = readerType::New();
reader1->SetFileName("../AffineFiles/AffineTemplate21.nii");
reader1->Update();
ImageType::Pointer fixedImage = reader1->GetOutput();
typedef itk::ImageFileReader<ImageType>readerType;
readerType::Pointer reader2 = readerType::New();
reader2->SetFileName("../AffineFiles/KKI2009-01-MPRAGE.nii"); //moving
reader2->Update();
ImageType::Pointer movingImage = reader2->GetOutput();
//typedef itk::CoordinateRep<double>CoordinateRepType;
using CoordinateRepType = double;
typedef itk::BSplineTransform<double, 3, 3> TransformType;
typedef itk::LBFGSBOptimizerv4 OptimizerType;
typedef itk::MeanSquaresImageToImageMetricv4<ImageType, ImageType> MetricType;
typedef itk::ImageRegistrationMethodv4<ImageType,ImageType> RegistrationType;
typedef itk::BSplineTransformInitializer<TransformType,ImageType> InitializerType;
typedef itk::ResampleImageFilter<ImageType,ImageType> ResampleFilterType;
TransformType::Pointer outputBSplineTransform = TransformType::New();
OptimizerType::Pointer optimizer = OptimizerType::New();
MetricType::Pointer metric = MetricType::New();
RegistrationType::Pointer registration = RegistrationType::New();
InitializerType::Pointer initializer = InitializerType::New(); //transformInitializer
ResampleFilterType::Pointer resample = ResampleFilterType::New();
registration->SetMetric(metric);
registration->SetOptimizer(optimizer);
registration->SetMovingImage(movingImage);
registration->SetFixedImage(fixedImage);
TransformType::MeshSizeType meshSize;
meshSize.Fill(numberOfGridNodesInOneDimension - SplineOrder);
initializer->SetTransform(outputBSplineTransform);
initializer->SetImage(fixedImage);
initializer->SetTransformDomainMeshSize(meshSize);
initializer->InitializeTransform();
//Set transform to identity
using ParametersType = TransformType::ParametersType;
const unsigned int numberOfParameters = outputBSplineTransform->GetNumberOfParameters();
ParametersType parameters(numberOfParameters);
parameters.Fill(0.0);
outputBSplineTransform->SetParameters(parameters);
registration->SetInitialTransform(outputBSplineTransform);
registration->InPlaceOn();
// A single level registration process is run using
// the shrink factor 1 and smoothing sigma 0.
const unsigned int numberOfLevels = 1;
RegistrationType::ShrinkFactorsArrayType shrinkFactorsPerLevel;
shrinkFactorsPerLevel.SetSize(numberOfLevels);
shrinkFactorsPerLevel[0] = 1;
RegistrationType::SmoothingSigmasArrayType smoothingSigmasPerLevel;
smoothingSigmasPerLevel.SetSize(numberOfLevels);
smoothingSigmasPerLevel[0] = 0;
registration->SetNumberOfLevels(numberOfLevels);
registration->SetSmoothingSigmasPerLevel(smoothingSigmasPerLevel);
registration->SetShrinkFactorsPerLevel(shrinkFactorsPerLevel);
const unsigned int numParameters = outputBSplineTransform->GetNumberOfParameters();
OptimizerType::BoundSelectionType boundSelect(numParameters);
OptimizerType::BoundValueType upperBound(numParameters);
OptimizerType::BoundValueType lowerBound(numParameters);
boundSelect.Fill(OptimizerType::UNBOUNDED);
upperBound.Fill(0.0);
lowerBound.Fill(0.0);
optimizer->SetBoundSelection(boundSelect);
optimizer->SetUpperBound(upperBound);
optimizer->SetLowerBound(lowerBound);
optimizer->SetCostFunctionConvergenceFactor(1e+12);
optimizer->SetGradientConvergenceTolerance(1.0e-35);
optimizer->SetNumberOfIterations(500);
optimizer->SetMaximumNumberOfFunctionEvaluations(500);
optimizer->SetMaximumNumberOfCorrections(5);
//Create the Command observer and register it with the optimizer.
CommandIterationUpdate::Pointer observer = CommandIterationUpdate::New();
optimizer->AddObserver(itk::IterationEvent(), observer);
std::cout << "Starting Registration " << std::endl;
//Run the registration
try
{
registration->Update();
std::cout << "Optimizer stop condition = "
<< registration->GetOptimizer()->GetStopConditionDescription()
<< std::endl;
}
catch (itk::ExceptionObject & err)
{
std::cerr << "ExceptionObject caught !" << std::endl;
std::cerr << err << std::endl;
return EXIT_FAILURE;
}
OptimizerType::ParametersType finalParameters = outputBSplineTransform->GetParameters();
std::cout << "Last Transform Parameters" << std::endl;
std::cout << finalParameters << std::endl;
// Finally we use the last transform in order to resample the image.
resample->SetTransform(outputBSplineTransform);
resample->SetInput(movingImage);
resample->SetSize(fixedImage->GetLargestPossibleRegion().GetSize());
resample->SetOutputOrigin(fixedImage->GetOrigin());
resample->SetOutputSpacing(fixedImage->GetSpacing());
resample->SetOutputDirection(fixedImage->GetDirection());
resample->SetDefaultPixelValue(100); // 0 or 100
//using OutputPixelType = unsigned char;
//using OutputImageType = itk::Image<OutputPixelType, ImageDimension>;
typedef itk::CastImageFilter<ImageType,ImageType>CastFilterType;
typedef itk::ImageFileWriter<ImageType>writerType;
CastFilterType::Pointer caster = CastFilterType::New();
writerType::Pointer writer = writerType::New();
caster->SetInput(resample->GetOutput());
writer->SetInput(caster->GetOutput());
writer->SetFileName("../DeformReg.nii");
// writer->Update();
try
{
writer->Update();
}
catch (itk::ExceptionObject & err)
{
std::cerr << "ExceptionObject caught !" << std::endl;
std::cerr << err << std::endl;
return EXIT_FAILURE;
}
typedef itk::SquaredDifferenceImageFilter<ImageType,ImageType,ImageType>SquaredDifferenceImageFilterType;
SquaredDifferenceImageFilterType::Pointer difference = SquaredDifferenceImageFilterType::New();
writerType::Pointer writer2 = writerType::New();
writer2->SetInput(difference->GetOutput());
// Compute the difference image between the
// fixed and resampled moving image.
if (argc >= 5)
{
difference->SetInput1(fixedImage->GetOutput());
difference->SetInput2(resample->GetOutput());
writer2->SetFileName(argv[4]);
try
{
writer2->Update();
}
catch (itk::ExceptionObject & err)
{
std::cerr << "ExceptionObject caught !" << std::endl;
std::cerr << err << std::endl;
return EXIT_FAILURE;
}
}
// Compute the difference image between the
// fixed and moving image before registration.
if (argc >= 6)
{
writer2->SetFileName(argv[5]);
difference->SetInput1(fixedImage->GetOutput());
difference->SetInput2(movingImage->GetOutput());
try
{
writer2->Update();
}
catch (itk::ExceptionObject & err)
{
std::cerr << "ExceptionObject caught !" << std::endl;
std::cerr << err << std::endl;
return EXIT_FAILURE;
}
}
// Generate the explicit deformation field resulting from
// the registration.
using VectorPixelType = itk::Vector<float, nDims>;
using DisplacementFieldImageType = itk::Image<VectorPixelType, nDims>;
// using DisplacementFieldGeneratorType = itk::TransformToDisplacementFieldFilter<DisplacementFieldImageType, \
// CoordinateRepType>;
typedef itk::TransformToDisplacementFieldFilter<ImageType,double>DisplacementFieldGeneratorType;
DisplacementFieldGeneratorType::Pointer dispfieldGenerator = DisplacementFieldGeneratorType::New();
dispfieldGenerator->UseReferenceImageOn();
dispfieldGenerator->SetReferenceImage(fixedImage);
dispfieldGenerator->SetTransform(outputBSplineTransform);
try
{
dispfieldGenerator->Update();
}
catch (itk::ExceptionObject & err)
{
std::cerr << "Exception detected while generating deformation field";
std::cerr << " : " << err << std::endl;
return EXIT_FAILURE;
}
using FieldWriterType = itk::ImageFileWriter<DisplacementFieldImageType>;
FieldWriterType::Pointer fieldWriter = FieldWriterType::New();
fieldWriter->SetInput(dispfieldGenerator->GetOutput());
if (argc >= 7)
{
fieldWriter->SetFileName(argv[6]);
try
{
fieldWriter->Update();
}
catch (itk::ExceptionObject & excp)
{
std::cerr << "Exception thrown " << std::endl;
std::cerr << excp << std::endl;
return EXIT_FAILURE;
}
}
return EXIT_SUCCESS;
//Done
// return 0;
}
I am getting the following errors:
[banikr@gw345 DeformReg-build]$ make
Scanning dependencies of target DeformReg
[ 50%] Building CXX object CMakeFiles/DeformReg.dir/DeformReg.cxx.o
/home/banikr/DeformReg/DeformReg.cxx: In function ‘int main(int, char**)’:
/home/banikr/DeformReg/DeformReg.cxx:216:41: error: ‘itk::SmartPointer<itk::Image<double, 3u> >::ObjectType {aka class itk::Image<double, 3u>}’ has no member named ‘GetOutput’
difference->SetInput1(fixedImage->GetOutput());
^
/home/banikr/DeformReg/DeformReg.cxx:237:41: error: ‘itk::SmartPointer<itk::Image<double, 3u> >::ObjectType {aka class itk::Image<double, 3u>}’ has no member named ‘GetOutput’
difference->SetInput1(fixedImage->GetOutput());
^
/home/banikr/DeformReg/DeformReg.cxx:238:42: error: ‘itk::SmartPointer<itk::Image<double, 3u> >::ObjectType {aka class itk::Image<double, 3u>}’ has no member named ‘GetOutput’
difference->SetInput2(movingImage->GetOutput());
^
In file included from /home/banikr/DeformReg/DeformReg.cxx:14:0:
/accre/arch/easybuild/software/Compiler/GCC/5.4.0-2.26/ITK/4.12.2-Python-2.7.12/include/ITK-4.12/itkTransformToDisplacementFieldFilter.h: In instantiation of ‘class itk::TransformToDisplacementFieldFilter<itk::Image<double, 3u>, double>’:
/home/banikr/DeformReg/DeformReg.cxx:260:33: required from here
/accre/arch/easybuild/software/Compiler/GCC/5.4.0-2.26/ITK/4.12.2-Python-2.7.12/include/ITK-4.12/itkTransformToDisplacementFieldFilter.h:84:51: error: ‘itk::TransformToDisplacementFieldFilter<itk::Image<double, 3u>, double>::PixelType {aka double}’ is not a class, struct, or union type
typedef typename PixelType::ValueType PixelValueType;
^
In file included from /accre/arch/easybuild/software/Compiler/GCC/5.4.0-2.26/ITK/4.12.2-Python-2.7.12/include/ITK-4.12/itkLightObject.h:21:0,
from /accre/arch/easybuild/software/Compiler/GCC/5.4.0-2.26/ITK/4.12.2-Python-2.7.12/include/ITK-4.12/itkObject.h:31,
from /accre/arch/easybuild/software/Compiler/GCC/5.4.0-2.26/ITK/4.12.2-Python-2.7.12/include/ITK-4.12/itkRegion.h:31,
from /accre/arch/easybuild/software/Compiler/GCC/5.4.0-2.26/ITK/4.12.2-Python-2.7.12/include/ITK-4.12/itkImageRegion.h:31,
from /accre/arch/easybuild/software/Compiler/GCC/5.4.0-2.26/ITK/4.12.2-Python-2.7.12/include/ITK-4.12/itkImage.h:21,
from /home/banikr/DeformReg/DeformReg.cxx:2:
/accre/arch/easybuild/software/Compiler/GCC/5.4.0-2.26/ITK/4.12.2-Python-2.7.12/include/ITK-4.12/itkTransformToDisplacementFieldFilter.h: In instantiation of ‘constexpr const unsigned int itk::TransformToDisplacementFieldFilter<itk::Image<double, 3u>, double>::PixelDimension’:
/accre/arch/easybuild/software/Compiler/GCC/5.4.0-2.26/ITK/4.12.2-Python-2.7.12/include/ITK-4.12/itkTransformToDisplacementFieldFilter.h:153:3: required from ‘class itk::TransformToDisplacementFieldFilter<itk::Image<double, 3u>, double>’
/home/banikr/DeformReg/DeformReg.cxx:260:33: required from here
/accre/arch/easybuild/software/Compiler/GCC/5.4.0-2.26/ITK/4.12.2-Python-2.7.12/include/ITK-4.12/itkTransformToDisplacementFieldFilter.h:151:3: error: ‘Dimension’ is not a member of ‘itk::TransformToDisplacementFieldFilter<itk::Image<double, 3u>, double>::PixelType {aka double}’
itkStaticConstMacro(PixelDimension, unsigned int,
^
/home/banikr/DeformReg/DeformReg.cxx:278:56: error: no matching function for call to ‘itk::ImageFileWriter<itk::Image<itk::Vector<float, 3u>, 3u> >::SetInput(itk::ImageSource<itk::Image<double, 3u> >::OutputImageType*)’
fieldWriter->SetInput(dispfieldGenerator->GetOutput());
^
In file included from /accre/arch/easybuild/software/Compiler/GCC/5.4.0-2.26/ITK/4.12.2-Python-2.7.12/include/ITK-4.12/itkImageFileWriter.h:226:0,
from /home/banikr/DeformReg/DeformReg.cxx:4:
/accre/arch/easybuild/software/Compiler/GCC/5.4.0-2.26/ITK/4.12.2-Python-2.7.12/include/ITK-4.12/itkImageFileWriter.hxx:62:1: note: candidate: void itk::ImageFileWriter<TInputImage>::SetInput(const InputImageType*) [with TInputImage = itk::Image<itk::Vector<float, 3u>, 3u>; itk::ImageFileWriter<TInputImage>::InputImageType = itk::Image<itk::Vector<float, 3u>, 3u>]
ImageFileWriter< TInputImage >
^
/accre/arch/easybuild/software/Compiler/GCC/5.4.0-2.26/ITK/4.12.2-Python-2.7.12/include/ITK-4.12/itkImageFileWriter.hxx:62:1: note: no known conversion for argument 1 from ‘itk::ImageSource<itk::Image<double, 3u> >::OutputImageType* {aka itk::Image<double, 3u>*}’ to ‘const InputImageType* {aka const itk::Image<itk::Vector<float, 3u>, 3u>*}’
In file included from /accre/arch/easybuild/software/Compiler/GCC/5.4.0-2.26/ITK/4.12.2-Python-2.7.12/include/ITK-4.12/itkImageBase.hxx:34:0,
from /accre/arch/easybuild/software/Compiler/GCC/5.4.0-2.26/ITK/4.12.2-Python-2.7.12/include/ITK-4.12/itkImageBase.h:786,
from /accre/arch/easybuild/software/Compiler/GCC/5.4.0-2.26/ITK/4.12.2-Python-2.7.12/include/ITK-4.12/itkNeighborhoodAccessorFunctor.h:22,
from /accre/arch/easybuild/software/Compiler/GCC/5.4.0-2.26/ITK/4.12.2-Python-2.7.12/include/ITK-4.12/itkImage.h:28,
from /home/banikr/DeformReg/DeformReg.cxx:2:
/accre/arch/easybuild/software/Compiler/GCC/5.4.0-2.26/ITK/4.12.2-Python-2.7.12/include/ITK-4.12/itkProcessObject.h:521:16: note: candidate: virtual void itk::ProcessObject::SetInput(const DataObjectIdentifierType&, itk::DataObject*)
virtual void SetInput(const DataObjectIdentifierType & key, DataObject *input);
^
/accre/arch/easybuild/software/Compiler/GCC/5.4.0-2.26/ITK/4.12.2-Python-2.7.12/include/ITK-4.12/itkProcessObject.h:521:16: note: candidate expects 2 arguments, 1 provided
In file included from /accre/arch/easybuild/software/Compiler/GCC/5.4.0-2.26/ITK/4.12.2-Python-2.7.12/include/ITK-4.12/itkTransformToDisplacementFieldFilter.h:197:0,
from /home/banikr/DeformReg/DeformReg.cxx:14:
/accre/arch/easybuild/software/Compiler/GCC/5.4.0-2.26/ITK/4.12.2-Python-2.7.12/include/ITK-4.12/itkTransformToDisplacementFieldFilter.hxx: In instantiation of ‘void itk::TransformToDisplacementFieldFilter<TOutputImage, TParametersValueType>::LinearThreadedGenerateData(const OutputImageRegionType&, itk::ThreadIdType) [with TOutputImage = itk::Image<double, 3u>; TParametersValueType = double; itk::TransformToDisplacementFieldFilter<TOutputImage, TParametersValueType>::OutputImageRegionType = itk::ImageRegion<3u>; itk::ThreadIdType = unsigned int]’:
/accre/arch/easybuild/software/Compiler/GCC/5.4.0-2.26/ITK/4.12.2-Python-2.7.12/include/ITK-4.12/itkTransformToDisplacementFieldFilter.hxx:171:5: required from ‘void itk::TransformToDisplacementFieldFilter<TOutputImage, TParametersValueType>::ThreadedGenerateData(const OutputImageRegionType&, itk::ThreadIdType) [with TOutputImage = itk::Image<double, 3u>; TParametersValueType = double; itk::TransformToDisplacementFieldFilter<TOutputImage, TParametersValueType>::OutputImageRegionType = itk::ImageRegion<3u>; itk::ThreadIdType = unsigned int]’
/home/banikr/DeformReg/DeformReg.cxx:301:1: required from here
/accre/arch/easybuild/software/Compiler/GCC/5.4.0-2.26/ITK/4.12.2-Python-2.7.12/include/ITK-4.12/itkTransformToDisplacementFieldFilter.hxx:283:20: error: cannot convert ‘itk::Vector<double, 3u>’ to ‘itk::TransformToDisplacementFieldFilter<itk::Image<double, 3u>, double>::PixelType {aka double}’ in assignment
displacement = transformedPoint - outputPoint;
^
/accre/arch/easybuild/software/Compiler/GCC/5.4.0-2.26/ITK/4.12.2-Python-2.7.12/include/ITK-4.12/itkTransformToDisplacementFieldFilter.hxx: In instantiation of ‘void itk::TransformToDisplacementFieldFilter<TOutputImage, TParametersValueType>::NonlinearThreadedGenerateData(const OutputImageRegionType&, itk::ThreadIdType) [with TOutputImage = itk::Image<double, 3u>; TParametersValueType = double; itk::TransformToDisplacementFieldFilter<TOutputImage, TParametersValueType>::OutputImageRegionType = itk::ImageRegion<3u>; itk::ThreadIdType = unsigned int]’:
/accre/arch/easybuild/software/Compiler/GCC/5.4.0-2.26/ITK/4.12.2-Python-2.7.12/include/ITK-4.12/itkTransformToDisplacementFieldFilter.hxx:177:3: required from ‘void itk::TransformToDisplacementFieldFilter<TOutputImage, TParametersValueType>::ThreadedGenerateData(const OutputImageRegionType&, itk::ThreadIdType) [with TOutputImage = itk::Image<double, 3u>; TParametersValueType = double; itk::TransformToDisplacementFieldFilter<TOutputImage, TParametersValueType>::OutputImageRegionType = itk::ImageRegion<3u>; itk::ThreadIdType = unsigned int]’
/home/banikr/DeformReg/DeformReg.cxx:301:1: required from here
/accre/arch/easybuild/software/Compiler/GCC/5.4.0-2.26/ITK/4.12.2-Python-2.7.12/include/ITK-4.12/itkTransformToDisplacementFieldFilter.hxx:213:18: error: cannot convert ‘itk::Vector<double, 3u>’ to ‘itk::TransformToDisplacementFieldFilter<itk::Image<double, 3u>, double>::PixelType {aka double}’ in assignment
displacement = transformedPoint - outputPoint;
^
make[2]: *** [CMakeFiles/DeformReg.dir/DeformReg.cxx.o] Error 1
make[1]: *** [CMakeFiles/DeformReg.dir/all] Error 2
make: *** [all] Error 2
[banikr@gw345 DeformReg-build]$
The problem I am trying to address: deform register 21 images to make a deform registration atlas.
I am using both typedef and using since I am not familiar with using syntax but the code in the example follows that.