X-Git-Url: https://git.creatis.insa-lyon.fr/pubgit/?a=blobdiff_plain;f=common%2FclitkCoeffsToDVF.h;h=911e8a9f7b7dccb8aa41746ea34617a338be2c5e;hb=1d70c4971d5dd20b399252a4ea339b7c878d9b74;hp=0fe30bb8c2d85127afc837227e0b1987329e7347;hpb=5be83a695bafa815cca85e180c0070e30020fa4c;p=clitk.git

diff --git a/common/clitkCoeffsToDVF.h b/common/clitkCoeffsToDVF.h
index 0fe30bb..911e8a9 100644
--- a/common/clitkCoeffsToDVF.h
+++ b/common/clitkCoeffsToDVF.h
@@ -18,69 +18,119 @@
 #ifndef clitkCoeffsToDVF_h
 #define clitkCoeffsToDVF_h
 
+#include "itkImageFileReader.h"
+#include "itkImageIOBase.h"
+
 #include "clitkBSplineDeformableTransform.h"
+#include "clitkResampleBSplineDeformableTransformImageFilter.h"
 #if ITK_VERSION_MAJOR >= 4
 #include "itkTransformToDisplacementFieldSource.h"
 #else
 #include "itkTransformToDeformationFieldSource.h"
 #endif
+#include "itkBSplineDeformableTransform.h"
 
-//-------------------------------------------------------------------
-// Convert Coefficient image to DVF
-//-------------------------------------------------------------------
-template<class DisplacementFieldType>
-typename DisplacementFieldType::Pointer
-CoeffsToDVF(std::string fileName, std::string likeFileName, bool verbose = false)
+namespace clitk
 {
-  typedef clitk::BSplineDeformableTransform<double, DisplacementFieldType::ImageDimension, DisplacementFieldType::ImageDimension> TransformType;
-  typedef typename TransformType::CoefficientImageType CoefficientImageType;
+  //-------------------------------------------------------------------
+  // Initialize transform from coefficient images
+  //-------------------------------------------------------------------
+  template <class TransformType>
+  void 
+  SetInitialTransformParameters(typename TransformType::Pointer transform, const typename TransformType::CoefficientImageType::Pointer coefficientImage, typename TransformType::CoefficientImageType::SpacingType outputSpacing) 
+  {
+    unsigned int dim = TransformType::CoefficientImageType::ImageDimension;
+    transform->SetSplineOrder(3);
+    transform->SetGridRegion( coefficientImage->GetLargestPossibleRegion() );
+    transform->SetGridOrigin( coefficientImage->GetOrigin() );
+    transform->SetGridSpacing( coefficientImage->GetSpacing() );
+    transform->SetGridDirection( coefficientImage->GetDirection() );
+    typename TransformType::RegionType::SizeType samplingFactors;
+    for (unsigned int i=0; i< dim; i++) {
+        samplingFactors[i]= (int) ( coefficientImage->GetSpacing()[i]/ outputSpacing[i]);
+    }
+    transform->SetLUTSamplingFactors(samplingFactors);
+    
+    typedef typename TransformType::ParametersType     ParametersType;
+    const unsigned int numberOfParameters = transform->GetNumberOfParameters();
+    ParametersType params(numberOfParameters);
+    params.Fill( 0.0 );
 
-  typedef itk::ImageFileReader<CoefficientImageType> CoeffReaderType;
-  typename CoeffReaderType::Pointer reader = CoeffReaderType::New();
-  reader->SetFileName(fileName);
-  reader->Update();
+    typedef itk::ImageRegionConstIterator<typename TransformType::CoefficientImageType> Iterator;
+    Iterator it (coefficientImage, coefficientImage->GetLargestPossibleRegion() );
+    it.GoToBegin();
+    unsigned int k = 0;
+    while (! it.IsAtEnd()) {
+        for (unsigned int j = 0; j < dim; j++)
+            params[k+j]=it.Get()[j];
 
-  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
+        ++it;
+        k += dim;
+    }
 
-  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];
+    transform->SetParameters(params);
+    transform->SetBulkTransform(NULL);
   }
-  
-  if (verbose) {
-    std::cout << "Interpolating coefficients with grid:" << std::endl;
-    std::cout << output_size << output_spacing << std::endl;
+
+  //-------------------------------------------------------------------
+  // Convert Coefficient image to DVF
+  //-------------------------------------------------------------------
+  template<class DisplacementFieldType>
+  typename DisplacementFieldType::Pointer
+  BLUTCoeffsToDVF(std::string fileName, std::string likeFileName, bool verbose = false)
+  {
+    const unsigned int dim = DisplacementFieldType::ImageDimension;
+    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();
+
+  #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 < dim; 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];
+    }
+
+    typename CoefficientImageType::Pointer coeffs = reader->GetOutput();
+    typename TransformType::Pointer transform = TransformType::New();
+    SetInitialTransformParameters<TransformType>(transform, coeffs, output_spacing);
+
+    if (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();
   }
-  
-  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();
 }
 
 #endif