#include "clitkVFResampleGenericFilter.h"
//--------------------------------------------------------------------
-clitk::VFResampleGenericFilter::VFResampleGenericFilter() {
+clitk::VFResampleGenericFilter::VFResampleGenericFilter():
+ clitk::ImageToImageGenericFilter<Self>("VFResample") {
+ InitializeImageType<2>();
+ InitializeImageType<3>();
+ InitializeImageType<4>();
mApplyGaussianFilterBefore = false;
mDefaultPixelValue = 0.0;
mInterpolatorName = "NN";
}
//--------------------------------------------------------------------
+
+//--------------------------------------------------------------------
+template<unsigned int Dim>
+void clitk::VFResampleGenericFilter::InitializeImageType() {
+ ADD_IMAGE_TYPE(Dim, float);
+}
+//--------------------------------------------------------------------
+
+
+//--------------------------------------------------------------------
+template<class ImageType>
+void clitk::VFResampleGenericFilter::UpdateWithInputImageType() {
+
+ if (mNbOfComponents == 1) {
+ std::cerr << "Error, only one components ? Use clitkImageResample instead." << std::endl;
+ exit(0);
+ }
+ typedef typename ImageType::PixelType PixelType;
+ if (mNbOfComponents == 2) Update_WithDimAndPixelTypeAndComponent<ImageType::ImageDimension,PixelType,2>();
+ if (mNbOfComponents == 3) Update_WithDimAndPixelTypeAndComponent<ImageType::ImageDimension,PixelType,3>();
+ if (mNbOfComponents == 4) Update_WithDimAndPixelTypeAndComponent<ImageType::ImageDimension,PixelType,4>();
+}
+//--------------------------------------------------------------------
+
+
+//--------------------------------------------------------------------
+template<unsigned int Dim, class PixelType, unsigned int DimCompo>
+void clitk::VFResampleGenericFilter::Update_WithDimAndPixelTypeAndComponent() {
+ // Reading input
+ typedef itk::Vector<PixelType, DimCompo> DisplacementType;
+ typedef itk::Image< DisplacementType, Dim > ImageType;
+
+ typename ImageType::Pointer input = clitk::readImage<ImageType>(mInputFilenames, mIOVerbose);
+
+ // Main filter
+ typename ImageType::Pointer outputImage = ComputeImage<ImageType>(input);
+
+ // Write results
+ SetNextOutput<ImageType>(outputImage);
+}
+//--------------------------------------------------------------------
+
//--------------------------------------------------------------------
-void clitk::VFResampleGenericFilter::Update() {
- // Load image header
- itk::ImageIOBase::Pointer header = clitk::readImageHeader(mInputFilenames[0]);
+template<class ImageType>
+typename ImageType::Pointer
+clitk::VFResampleGenericFilter::ComputeImage(typename ImageType::Pointer inputImage) {
+
+ // Check options
+ static unsigned int dim = ImageType::ImageDimension;
+ if (mOutputSize.size() != dim) {
+ std::cerr << "Please set size with " << dim << " dimensions." << std::endl;
+ return NULL;
+ }
+ if (mOutputSpacing.size() != dim) {
+ std::cerr << "Please set spacing with " << dim << " dimensions." << std::endl;
+ return NULL;
+ }
+ mOutputOrigin.resize(dim);
+
+ if (mApplyGaussianFilterBefore && mSigma.size() != dim) {
+ std::cerr << "Please set sigma with " << dim << " dimensions." << std::endl;
+ return NULL;
+ }
+
+ // Some typedefs
+ typedef typename ImageType::SizeType SizeType;
+ typedef typename ImageType::SpacingType SpacingType;
+ typedef typename ImageType::PointType PointType;
+
+ // Create Image Filter
+ typedef itk::VectorResampleImageFilter<ImageType,ImageType> FilterType;
+ typename FilterType::Pointer filter = FilterType::New();
+
+ // Instance of the transform object to be passed to the resample
+ // filter. By default, identity transform is applied
+ typedef itk::AffineTransform<double, ImageType::ImageDimension> TransformType;
+ typename TransformType::Pointer transform = TransformType::New();
+ filter->SetTransform(transform);
+
+ // Set filter's parameters
+ SizeType outputSize;
+ SpacingType outputSpacing;
+ PointType outputOrigin;
+ for(unsigned int i=0; i<ImageType::ImageDimension; i++) {
+ outputSize[i] = mOutputSize[i];
+ outputSpacing[i] = mOutputSpacing[i];
+ outputOrigin[i] = inputImage->GetOrigin()[i];
+ }
+
+ filter->SetSize(outputSize);
+ filter->SetOutputSpacing(outputSpacing);
+ filter->SetOutputOrigin(outputOrigin);
+ filter->SetDefaultPixelValue(static_cast<typename ImageType::PixelType>(mDefaultPixelValue));
+
+ // Select interpolator
+ if (mInterpolatorName == "nn") {
+ typedef itk::VectorNearestNeighborInterpolateImageFunction<ImageType, double> InterpolatorType;
+ typename InterpolatorType::Pointer interpolator = InterpolatorType::New();
+ filter->SetInterpolator(interpolator);
+ }
+ else {
+ if (mInterpolatorName == "linear") {
+ typedef itk::VectorLinearInterpolateImageFunction<ImageType, double> InterpolatorType;
+ typename InterpolatorType::Pointer interpolator = InterpolatorType::New();
+ filter->SetInterpolator(interpolator);
+ }
+ else {
+ std::cerr << "Sorry, I do not know the interpolator (for vector field) '" << mInterpolatorName
+ << "'. Known interpolators are : nn, linear" << std::endl;
+ exit(0);
+ }
+ }
+
+ // Build initial Gaussian bluring (if needed)
+ typedef itk::RecursiveGaussianImageFilter<ImageType, ImageType> GaussianFilterType;
+ std::vector<typename GaussianFilterType::Pointer> gaussianFilters;
+ if (mApplyGaussianFilterBefore) {
+ for(unsigned int i=0; i<ImageType::ImageDimension; i++) {
+ // Create filter
+ gaussianFilters.push_back(GaussianFilterType::New());
+ // Set options
+ gaussianFilters[i]->SetDirection(i);
+ gaussianFilters[i]->SetOrder(GaussianFilterType::ZeroOrder);
+ gaussianFilters[i]->SetNormalizeAcrossScale(false);
+ gaussianFilters[i]->SetSigma(mSigma[i]); // in millimeter !
+ // Set input
+ if (i==0) gaussianFilters[i]->SetInput(inputImage);
+ else gaussianFilters[i]->SetInput(gaussianFilters[i-1]->GetOutput());
+ }
+ filter->SetInput(gaussianFilters[ImageType::ImageDimension-1]->GetOutput());
+ }
+ else {
+ filter->SetInput(inputImage);
+ }
+
+ // Go !
+ try {
+ filter->Update();
+ }
+ catch( itk::ExceptionObject & err ) {
+ std::cerr << "Error while filtering " << mInputFilenames[0].c_str()
+ << " " << err << std::endl;
+ exit(0);
+ }
+
+ // Return result
+ return filter->GetOutput();
- // Determine dim, pixel type, number of components
- mDim = header->GetNumberOfDimensions();
- mPixelTypeName = header->GetComponentTypeAsString(header->GetComponentType());
- mNbOfComponents = header->GetNumberOfComponents();
-
- // Switch by dimension
- if (mDim == 2) { Update_WithDim<2>(); return; }
- if (mDim == 3) { Update_WithDim<3>(); return; }
- if (mDim == 4) { Update_WithDim<4>(); return; }
-
- std::cerr << "Error, dimension of input image is " << mDim << ", but I only work with 2,3,4." << std::endl;
- exit(0);
}
//--------------------------------------------------------------------
+
//--------------------------------------------------------------------
void clitk::VFResampleGenericFilter::SetOutputSize(const std::vector<int> & size) {
mOutputSize.resize(size.size());
git://git.creatis.insa-lyon.fr / clitk.git / blobdiff