1 #ifndef CLITKIMAGERESAMPLEGENERICFILTER2_CXX
2 #define CLITKIMAGERESAMPLEGENERICFILTER2_CXX
5 -------------------------------------------------------------------
6 * @file clitkImageResampleGenericFilter.cxx
7 * @author David Sarrut <David.Sarrut@creatis.insa-lyon.fr>
8 * @date 23 Feb 2008 08:37:53
12 -------------------------------------------------------------------*/
14 #include "clitkImageResampleGenericFilter.h"
16 //--------------------------------------------------------------------
17 clitk::ImageResampleGenericFilter::ImageResampleGenericFilter():
18 ImageToImageGenericFilter<Self>("ImageResample") {
19 mApplyGaussianFilterBefore = false;
20 mDefaultPixelValue = 0.0;
21 mInterpolatorName = "NN";
23 InitializeImageTypeWithDim<2>();
24 InitializeImageTypeWithDim<3>();
25 InitializeImageTypeWithDim<4>();
27 //--------------------------------------------------------------------
30 //--------------------------------------------------------------------
31 template<unsigned int Dim>
32 void clitk::ImageResampleGenericFilter::InitializeImageTypeWithDim() {
33 ADD_IMAGE_TYPE(Dim, char);
34 ADD_IMAGE_TYPE(Dim, short);
35 ADD_IMAGE_TYPE(Dim, int);
36 ADD_IMAGE_TYPE(Dim, float);
38 //--------------------------------------------------------------------
41 //--------------------------------------------------------------------
42 template<class InputImageType>
43 void clitk::ImageResampleGenericFilter::UpdateWithInputImageType() {
46 typedef typename InputImageType::SizeType SizeType;
47 typedef typename InputImageType::SpacingType SpacingType;
48 typedef typename InputImageType::PointType PointType;
49 typedef typename InputImageType::PixelType PixelType;
50 static unsigned int dim = InputImageType::ImageDimension;
53 typename InputImageType::Pointer input = this->GetInput<InputImageType>(0);
56 if (!std::numeric_limits<PixelType>::is_signed) {
57 if ((mInterpolatorName == "bspline") || (mInterpolatorName == "blut")) {
58 std::cerr << "Warning : input pixel type is not signed, use bspline interpolation at your own risk ..." << std::endl;
63 if (mOutputSize.size() != dim) {
64 std::cerr << "Please set size with " << dim << " dimensions." << std::endl;
67 if (mOutputSpacing.size() != dim) {
68 std::cerr << "Please set spacing with " << dim << " dimensions." << std::endl;
71 mOutputOrigin.resize(dim);
73 if (mApplyGaussianFilterBefore && mSigma.size() != dim) {
74 std::cerr << "Please set sigma with " << dim << " dimensions." << std::endl;
78 // Create Image Filter
79 typedef itk::ResampleImageFilter<InputImageType,InputImageType> FilterType;
80 typename FilterType::Pointer filter = FilterType::New();
82 // Instance of the transform object to be passed to the resample
83 // filter. By default, identity transform is applied
84 typedef itk::AffineTransform<double, InputImageType::ImageDimension> TransformType;
85 typename TransformType::Pointer transform = TransformType::New();
86 filter->SetTransform(transform);
88 // Set filter's parameters
90 SpacingType outputSpacing;
91 PointType outputOrigin;
92 for(unsigned int i=0; i<InputImageType::ImageDimension; i++) {
93 outputSize[i] = mOutputSize[i];
94 outputSpacing[i] = mOutputSpacing[i];
95 outputOrigin[i] = input->GetOrigin()[i];
98 filter->SetSize(outputSize);
99 filter->SetOutputSpacing(outputSpacing);
100 filter->SetOutputOrigin(outputOrigin);
101 filter->SetDefaultPixelValue(static_cast<PixelType>(mDefaultPixelValue));//DS TODO//JV comme ça?
103 // Select interpolator
104 if (mInterpolatorName == "nn") {
105 typedef itk::NearestNeighborInterpolateImageFunction<InputImageType, double> InterpolatorType;
106 typename InterpolatorType::Pointer interpolator = InterpolatorType::New();
107 filter->SetInterpolator(interpolator);
110 if (mInterpolatorName == "linear") {
111 typedef itk::LinearInterpolateImageFunction<InputImageType, double> InterpolatorType;
112 typename InterpolatorType::Pointer interpolator = InterpolatorType::New();
113 filter->SetInterpolator(interpolator);
116 if (mInterpolatorName == "bspline") {
117 typedef itk::BSplineInterpolateImageFunction<InputImageType, double> InterpolatorType;
118 typename InterpolatorType::Pointer interpolator = InterpolatorType::New();
119 interpolator->SetSplineOrder(mBSplineOrder);
120 filter->SetInterpolator(interpolator);
123 if (mInterpolatorName == "blut") {
124 typedef itk::BSplineInterpolateImageFunctionWithLUT<InputImageType, double> InterpolatorType;
125 typename InterpolatorType::Pointer interpolator = InterpolatorType::New();
126 interpolator->SetSplineOrder(mBSplineOrder);
127 interpolator->SetLUTSamplingFactor(mSamplingFactors[0]);
128 filter->SetInterpolator(interpolator);
131 std::cerr << "Sorry, I do not know the interpolator '" << mInterpolatorName
132 << "'. Known interpolators are : nn, linear, bspline, blut" << std::endl;
139 // Build initial Gaussian bluring (if needed)
140 typedef itk::RecursiveGaussianImageFilter<InputImageType, InputImageType> GaussianFilterType;
141 std::vector<typename GaussianFilterType::Pointer> gaussianFilters;
142 if (mApplyGaussianFilterBefore) {
143 for(unsigned int i=0; i<InputImageType::ImageDimension; i++) {
145 gaussianFilters.push_back(GaussianFilterType::New());
147 gaussianFilters[i]->SetDirection(i);
148 gaussianFilters[i]->SetOrder(GaussianFilterType::ZeroOrder);
149 gaussianFilters[i]->SetNormalizeAcrossScale(false);
150 gaussianFilters[i]->SetSigma(mSigma[i]); // in millimeter !
152 if (i==0) gaussianFilters[i]->SetInput(input);
153 else gaussianFilters[i]->SetInput(gaussianFilters[i-1]->GetOutput());
155 filter->SetInput(gaussianFilters[InputImageType::ImageDimension-1]->GetOutput());
158 filter->SetInput(input);
165 catch( itk::ExceptionObject & err ) {
166 std::cerr << "Error while filtering " << mInputFilenames[0].c_str()
167 << " " << err << std::endl;
172 typename InputImageType::Pointer outputImage = filter->GetOutput();
174 // Write/save results
175 this->SetNextOutput<InputImageType>(outputImage);
178 //--------------------------------------------------------------------
181 //--------------------------------------------------------------------
182 void clitk::ImageResampleGenericFilter::SetOutputSize(const std::vector<int> & size) {
183 mOutputSize.resize(size.size());
184 std::copy(size.begin(), size.end(), mOutputSize.begin());
186 //--------------------------------------------------------------------
188 //--------------------------------------------------------------------
189 void clitk::ImageResampleGenericFilter::SetOutputSpacing(const std::vector<double> & spacing) {
190 mOutputSpacing.resize(spacing.size());
191 std::copy(spacing.begin(), spacing.end(), mOutputSpacing.begin());
193 //--------------------------------------------------------------------
195 //--------------------------------------------------------------------
196 void clitk::ImageResampleGenericFilter::SetInterpolationName(const std::string & inter) {
197 mInterpolatorName = inter;
199 //--------------------------------------------------------------------
201 //--------------------------------------------------------------------
202 void clitk::ImageResampleGenericFilter::SetGaussianSigma(const std::vector<double> & sigma) {
203 mApplyGaussianFilterBefore = true;
204 mSigma.resize(sigma.size());
205 std::copy(sigma.begin(), sigma.end(), mSigma.begin());
207 //--------------------------------------------------------------------