From 75d074a4caf86e1c036bc7879557a1600ad0ceac Mon Sep 17 00:00:00 2001 From: Romulo Pinho <romulo.pinho@lyon.unicancer.fr> Date: Thu, 19 Jan 2012 17:39:42 +0100 Subject: [PATCH] image warping using B-Spline coefficients - uses the resolution of the input image for the output --- tools/clitkWarpImage.cxx | 7 +++ tools/clitkWarpImage.ggo | 5 +- tools/clitkWarpImageGenericFilter.h | 2 +- tools/clitkWarpImageGenericFilter.txx | 89 ++++++++++++++++++++++++--- 4 files changed, 91 insertions(+), 12 deletions(-) diff --git a/tools/clitkWarpImage.cxx b/tools/clitkWarpImage.cxx index 3360661..673a583 100644 --- a/tools/clitkWarpImage.cxx +++ b/tools/clitkWarpImage.cxx @@ -40,6 +40,13 @@ int main(int argc, char * argv[]) GGO(clitkWarpImage, args_info); CLITK_INIT; + if (args_info.coeff_given) { + if (args_info.verbose_flag) + std::cout<< "Using coefficient images, so forcing sapcing = 1.\n" << std::endl; + + args_info.spacing_arg = 1; + } + // Filter clitk::WarpImageGenericFilter::Pointer genericFilter=clitk::WarpImageGenericFilter::New(); diff --git a/tools/clitkWarpImage.ggo b/tools/clitkWarpImage.ggo index 8a572c4..12aa18e 100644 --- a/tools/clitkWarpImage.ggo +++ b/tools/clitkWarpImage.ggo @@ -10,7 +10,10 @@ section "I/O" option "input" i "Input image filename" string yes option "output" o "Output image filename" string yes -option "vf" - "Vector field filename" string yes + +defgroup "DVFoption" groupdesc="an option of this group is required" required +groupoption "vf" - "Vector field filename" string yes group="DVFoption" +groupoption "coeff" - "B-Spline coefficients filename" string yes group="DVFoption" section "Options" diff --git a/tools/clitkWarpImageGenericFilter.h b/tools/clitkWarpImageGenericFilter.h index 93eda73..9394d8d 100644 --- a/tools/clitkWarpImageGenericFilter.h +++ b/tools/clitkWarpImageGenericFilter.h @@ -98,7 +98,7 @@ namespace clitk //---------------------------------------- template <unsigned int Dimension> void UpdateWithDim(std::string PixelType); template <unsigned int Dimension, class PixelType> void UpdateWithDimAndPixelType(); - + template<class DisplacementFieldType> typename DisplacementFieldType::Pointer CoeffsToDVF(std::string fileName, std::string likeFileName); //---------------------------------------- // Data members diff --git a/tools/clitkWarpImageGenericFilter.txx b/tools/clitkWarpImageGenericFilter.txx index 8544df2..ce979aa 100644 --- a/tools/clitkWarpImageGenericFilter.txx +++ b/tools/clitkWarpImageGenericFilter.txx @@ -27,6 +27,13 @@ * ===================================================*/ +#include "itkVectorResampleImageFilter.h" +#include "clitkBSplineDeformableTransform.h" +#if ITK_VERSION_MAJOR >= 4 +#include "itkTransformToDisplacementFieldSource.h" +#else +#include "itkTransformToDeformationFieldSource.h" +#endif namespace clitk { @@ -64,6 +71,63 @@ WarpImageGenericFilter::UpdateWithDim(std::string PixelType) } } +//------------------------------------------------------------------- +// Convert Coefficient image to DVF +//------------------------------------------------------------------- +template<class DisplacementFieldType> +typename DisplacementFieldType::Pointer +WarpImageGenericFilter::CoeffsToDVF(std::string fileName, std::string likeFileName) +{ + typedef clitk::BSplineDeformableTransform<double, DisplacementFieldType::ImageDimension, DisplacementFieldType::ImageDimension> TransformType; + typedef typename TransformType::CoefficientImageType CoefficientImageType; + + typedef itk::ImageFileReader<CoefficientImageType> CoeffReaderType; + typename CoeffReaderType::Pointer reader = CoeffReaderType::New(); + reader->SetFileName(fileName); + reader->Update(); + + typename TransformType::Pointer transform = TransformType::New(); + transform->SetCoefficientImage(reader->GetOutput()); + +#if ITK_VERSION_MAJOR >= 4 + typedef itk::TransformToDisplacementFieldSource<DisplacementFieldType, double> ConvertorType; +#else + typedef itk::TransformToDeformationFieldSource<DisplacementFieldType, double> ConvertorType; +#endif + + typedef itk::ImageIOBase ImageIOType; + typename ImageIOType::Pointer imageIO = itk::ImageIOFactory::CreateImageIO(likeFileName.c_str(), itk::ImageIOFactory::ReadMode); + imageIO->SetFileName(likeFileName); + imageIO->ReadImageInformation(); + + typename ConvertorType::Pointer convertor= ConvertorType::New(); + typename ConvertorType::SizeType output_size; + typename ConvertorType::SpacingType output_spacing; + typename ConvertorType::OriginType output_origin; + typename ConvertorType::DirectionType output_direction; + for (unsigned int i = 0; i < DisplacementFieldType::ImageDimension; i++) { + output_size[i] = imageIO->GetDimensions(i); + output_spacing[i] = imageIO->GetSpacing(i); + output_origin[i] = imageIO->GetOrigin(i); + for (unsigned int j = 0; j < DisplacementFieldType::ImageDimension; j++) + output_direction[i][j] = imageIO->GetDirection(i)[j]; + } + + if (m_Verbose) { + std::cout << "Interpolating coefficients with grid:" << std::endl; + std::cout << output_size << output_spacing << std::endl; + } + + convertor->SetNumberOfThreads(1); + convertor->SetTransform(transform); + convertor->SetOutputOrigin(output_origin); + convertor->SetOutputSpacing(output_spacing); + convertor->SetOutputSize(output_size); + convertor->SetOutputDirection(output_direction); + convertor->Update(); + + return convertor->GetOutput(); +} //------------------------------------------------------------------- // Update with the number of dimensions and the pixeltype @@ -78,21 +142,26 @@ WarpImageGenericFilter::UpdateWithDimAndPixelType() typedef itk::Image<PixelType, Dimension> OutputImageType; typedef itk::Vector<float, Dimension> DisplacementType; typedef itk::Image<DisplacementType, Dimension> DeformationFieldType; - - + // Read the input typedef itk::ImageFileReader<InputImageType> InputReaderType; typename InputReaderType::Pointer reader = InputReaderType::New(); reader->SetFileName( m_InputFileName); reader->Update(); typename InputImageType::Pointer input= reader->GetOutput(); - - //Read the deformation field - typedef itk::ImageFileReader<DeformationFieldType> DeformationFieldReaderType; - typename DeformationFieldReaderType::Pointer deformationFieldReader= DeformationFieldReaderType::New(); - deformationFieldReader->SetFileName(m_ArgsInfo.vf_arg); - deformationFieldReader->Update(); - typename DeformationFieldType::Pointer deformationField =deformationFieldReader->GetOutput(); + + typename DeformationFieldType::Pointer deformationField; + if (m_ArgsInfo.coeff_given) { + deformationField = CoeffsToDVF<DeformationFieldType>(m_ArgsInfo.coeff_arg, m_InputFileName); + } + else { + //Read the deformation field + typedef itk::ImageFileReader<DeformationFieldType> DeformationFieldReaderType; + typename DeformationFieldReaderType::Pointer deformationFieldReader= DeformationFieldReaderType::New(); + deformationFieldReader->SetFileName(m_ArgsInfo.vf_arg); + deformationFieldReader->Update(); + deformationField =deformationFieldReader->GetOutput(); + } // Intensity interpolator typedef clitk::GenericVectorInterpolator<args_info_clitkWarpImage, DeformationFieldType, double> GenericInterpolatorType; @@ -139,7 +208,7 @@ WarpImageGenericFilter::UpdateWithDimAndPixelType() // ------------------------------------------- // Spacing like input // ------------------------------------------- - else { + else if (!m_ArgsInfo.coeff_given) { // Get size typename DeformationFieldType::SizeType newSize; for (unsigned int i=0 ; i <Dimension; i++) -- 2.47.1