File: ResampleImageFilter.cxx

package info (click to toggle)
insighttoolkit4 4.13.3withdata-dfsg1-4
  • links: PTS, VCS
  • area: main
  • in suites: bullseye
  • size: 489,260 kB
  • sloc: cpp: 557,342; ansic: 146,850; fortran: 34,788; python: 16,572; sh: 2,187; lisp: 2,070; tcl: 993; java: 362; perl: 200; makefile: 129; csh: 81; pascal: 69; xml: 19; ruby: 10
file content (451 lines) | stat: -rw-r--r-- 15,686 bytes parent folder | download | duplicates (5)
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
193
194
195
196
197
198
199
200
201
202
203
204
205
206
207
208
209
210
211
212
213
214
215
216
217
218
219
220
221
222
223
224
225
226
227
228
229
230
231
232
233
234
235
236
237
238
239
240
241
242
243
244
245
246
247
248
249
250
251
252
253
254
255
256
257
258
259
260
261
262
263
264
265
266
267
268
269
270
271
272
273
274
275
276
277
278
279
280
281
282
283
284
285
286
287
288
289
290
291
292
293
294
295
296
297
298
299
300
301
302
303
304
305
306
307
308
309
310
311
312
313
314
315
316
317
318
319
320
321
322
323
324
325
326
327
328
329
330
331
332
333
334
335
336
337
338
339
340
341
342
343
344
345
346
347
348
349
350
351
352
353
354
355
356
357
358
359
360
361
362
363
364
365
366
367
368
369
370
371
372
373
374
375
376
377
378
379
380
381
382
383
384
385
386
387
388
389
390
391
392
393
394
395
396
397
398
399
400
401
402
403
404
405
406
407
408
409
410
411
412
413
414
415
416
417
418
419
420
421
422
423
424
425
426
427
428
429
430
431
432
433
434
435
436
437
438
439
440
441
442
443
444
445
446
447
448
449
450
451
/*=========================================================================
 *
 *  Copyright Insight Software Consortium
 *
 *  Licensed under the Apache License, Version 2.0 (the "License");
 *  you may not use this file except in compliance with the License.
 *  You may obtain a copy of the License at
 *
 *         http://www.apache.org/licenses/LICENSE-2.0.txt
 *
 *  Unless required by applicable law or agreed to in writing, software
 *  distributed under the License is distributed on an "AS IS" BASIS,
 *  WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
 *  See the License for the specific language governing permissions and
 *  limitations under the License.
 *
 *=========================================================================*/

//  Software Guide : BeginCommandLineArgs
//    INPUTS:  {BrainProtonDensitySlice.png}
//    OUTPUTS: {ResampleImageFilterOutput1.png}
//    ARGUMENTS:    0
//  Software Guide : EndCommandLineArgs
//  Software Guide : BeginCommandLineArgs
//    INPUTS:  {BrainProtonDensitySlice.png}
//    OUTPUTS: {ResampleImageFilterOutput2.png}
//    ARGUMENTS:    1
//  Software Guide : EndCommandLineArgs
//  Software Guide : BeginCommandLineArgs
//    INPUTS:  {BrainProtonDensitySlice.png}
//    OUTPUTS: {ResampleImageFilterOutput3.png}
//    ARGUMENTS:    2
//  Software Guide : EndCommandLineArgs
//  Software Guide : BeginCommandLineArgs
//    INPUTS:  {BrainProtonDensitySlice.png}
//    OUTPUTS: {ResampleImageFilterOutput4.png}
//    ARGUMENTS:    3
//  Software Guide : EndCommandLineArgs
//
//  Software Guide : BeginLatex
//
//  Resampling an image is a very important task in image analysis. It is
//  especially important in the frame of image registration. The
//  \doxygen{ResampleImageFilter} implements image resampling through the use
//  of \doxygen{Transform}s. The inputs expected by this filter are an image,
//  a transform and an interpolator. The space coordinates of the image are
//  mapped through the transform in order to generate a new image. The extent
//  and spacing of the resulting image are selected by the user. Resampling
//  is performed in space coordinates, not pixel/grid coordinates. It is
//  quite important to ensure that image spacing is properly set on the
//  images involved. The interpolator is required since the mapping from one
//  space to the other will often require evaluation of the intensity of the
//  image at non-grid positions.
//
//  \index{itk::ResampleImageFilter}
//
//  Software Guide : EndLatex


#include "itkImage.h"
#include "itkImageFileReader.h"
#include "itkImageFileWriter.h"


//  Software Guide : BeginLatex
//
//  The header file corresponding to this filter should be included first.
//
//  \index{itk::ResampleImageFilter!header}
//
//  Software Guide : EndLatex


// Software Guide : BeginCodeSnippet
#include "itkResampleImageFilter.h"
// Software Guide : EndCodeSnippet


//  Software Guide : BeginLatex
//
//  The header files corresponding to the transform and interpolator must
//  also be included.
//
//  \index{itk::AffineTransform!header}
//  \index{itk::Nearest\-Neighbor\-Interpolate\-Image\-Function!header}
//
//  Software Guide : EndLatex


// Software Guide : BeginCodeSnippet
#include "itkAffineTransform.h"
#include "itkNearestNeighborInterpolateImageFunction.h"
// Software Guide : EndCodeSnippet


int main( int argc, char * argv[] )
{
  if( argc < 4 )
    {
    std::cerr << "Usage: " << std::endl;
    std::cerr << argv[0] << "  inputImageFile  outputImageFile";
    std::cerr << "  [exampleAction={0,1,2,3}]" << std::endl;
    return EXIT_FAILURE;
    }

  int exampleAction = 0;

  if( argc >= 4 )
    {
    exampleAction = atoi( argv[3] );
    }

  //  Software Guide : BeginLatex
  //
  //  The dimension and pixel types for input and output image must be
  //  defined and with them the image types can be instantiated.
  //
  //  Software Guide : EndLatex

  // Software Guide : BeginCodeSnippet
  const unsigned int                               Dimension = 2;
  typedef unsigned char                            InputPixelType;
  typedef unsigned char                            OutputPixelType;
  typedef itk::Image< InputPixelType,  Dimension > InputImageType;
  typedef itk::Image< OutputPixelType, Dimension > OutputImageType;
  // Software Guide : EndCodeSnippet


  typedef itk::ImageFileReader< InputImageType  >  ReaderType;
  typedef itk::ImageFileWriter< OutputImageType >  WriterType;

  ReaderType::Pointer reader = ReaderType::New();
  WriterType::Pointer writer = WriterType::New();

  reader->SetFileName( argv[1] );
  writer->SetFileName( argv[2] );


  //  Software Guide : BeginLatex
  //
  //  Using the image and transform types it is now possible to instantiate
  //  the filter type and create the filter object.
  //
  //  \index{itk::ResampleImageFilter!instantiation}
  //  \index{itk::ResampleImageFilter!New()}
  //  \index{itk::ResampleImageFilter!Pointer}
  //
  //  Software Guide : EndLatex

  // Software Guide : BeginCodeSnippet
  typedef itk::ResampleImageFilter<InputImageType,OutputImageType> FilterType;
  FilterType::Pointer filter = FilterType::New();
  // Software Guide : EndCodeSnippet


  //  Software Guide : BeginLatex
  //
  //  The transform type is typically defined using the image dimension
  //  and the type used for representing space coordinates.
  //
  //  \index{itk::AffineTransform!instantiation}
  //
  //  Software Guide : EndLatex


  // Software Guide : BeginCodeSnippet
  typedef itk::AffineTransform< double, Dimension >  TransformType;
  // Software Guide : EndCodeSnippet


  //  Software Guide : BeginLatex
  //
  //  An instance of the transform object is instantiated and passed to the
  //  resample filter. By default, the parameters of the transform are set to
  //  represent the identity transform.
  //
  //  \index{itk::ResampleImageFilter!SetTransform()}
  //
  //  Software Guide : EndLatex


  // Software Guide : BeginCodeSnippet
  TransformType::Pointer transform = TransformType::New();
  filter->SetTransform( transform );
  // Software Guide : EndCodeSnippet


  //  Software Guide : BeginLatex
  //
  //  The interpolator type is defined using the full image type and the type
  //  used for representing space coordinates.
  //
  //  \index{itk::Nearest\-Neighbor\-Interpolate\-Image\-Function!instantiation}
  //
  //  Software Guide : EndLatex

  // Software Guide : BeginCodeSnippet
  typedef itk::NearestNeighborInterpolateImageFunction<
                       InputImageType, double >  InterpolatorType;
  // Software Guide : EndCodeSnippet


  //  Software Guide : BeginLatex
  //
  //  An instance of the interpolator object is instantiated and passed to
  //  the resample filter.
  //
  //  \index{itk::ResampleImageFilter!SetInterpolator()}
  //
  //  Software Guide : EndLatex


  // Software Guide : BeginCodeSnippet
  InterpolatorType::Pointer interpolator = InterpolatorType::New();
  filter->SetInterpolator( interpolator );
  // Software Guide : EndCodeSnippet


  //  Software Guide : BeginLatex
  //
  //  Given that some pixels of the output image may end up being mapped
  //  outside the extent of the input image it is necessary to decide what
  //  values to assign to them. This is done by invoking the
  //  \code{SetDefaultPixelValue()} method.
  //
  //  \index{itk::ResampleImageFilter!SetDefaultPixelValue()}
  //
  //  Software Guide : EndLatex

  // Software Guide : BeginCodeSnippet
  filter->SetDefaultPixelValue( 0 );
  // Software Guide : EndCodeSnippet


  //  Software Guide : BeginLatex
  //
  //  The sampling grid of the output space is specified with the spacing along
  //  each dimension and the origin.
  //
  //  \index{itk::ResampleImageFilter!SetOutputOrigin()}
  //  \index{itk::ResampleImageFilter!SetOutputSpacing()}
  //
  //  Software Guide : EndLatex

  // Software Guide : BeginCodeSnippet
  // pixel spacing in millimeters along X and Y
  const double spacing[ Dimension ] = { 1.0, 1.0 };
  filter->SetOutputSpacing( spacing );

  // Physical space coordinate of origin for X and Y
  const double origin[ Dimension ] = { 0.0, 0.0 };
  filter->SetOutputOrigin( origin );
  // Software Guide : EndCodeSnippet


  // Software Guide : BeginCodeSnippet
  InputImageType::DirectionType direction;
  direction.SetIdentity();
  filter->SetOutputDirection( direction );
  // Software Guide : EndCodeSnippet


  //  Software Guide : BeginLatex
  //
  //  The extent of the sampling grid on the output image is defined by a
  //  \code{SizeType} and is set using the \code{SetSize()} method.
  //
  //  \index{itk::ResampleImageFilter!SetSize()}
  //
  //  Software Guide : EndLatex

  // Software Guide : BeginCodeSnippet
  InputImageType::SizeType   size;

  size[0] = 300;  // number of pixels along X
  size[1] = 300;  // number of pixels along Y

  filter->SetSize( size );
  // Software Guide : EndCodeSnippet


  //  Software Guide : BeginLatex
  //
  //  The input to the filter can be taken from any other filter, for example
  //  a reader. The output can be passed down the pipeline to other filters,
  //  for example a writer. An update call on any downstream filter will
  //  trigger the execution of the resampling filter.
  //
  //  \index{itk::ResampleImageFilter!SetInput()}
  //  \index{itk::ResampleImageFilter!GetOutput()}
  //
  //  Software Guide : EndLatex


  // Software Guide : BeginCodeSnippet
  filter->SetInput( reader->GetOutput() );
  writer->SetInput( filter->GetOutput() );
  writer->Update();
  // Software Guide : EndCodeSnippet


  //  Software Guide : BeginLatex
  //
  // \begin{figure}
  // \center
  // \includegraphics[width=0.44\textwidth]{BrainProtonDensitySlice}
  // \includegraphics[width=0.44\textwidth]{ResampleImageFilterOutput1}
  // \itkcaption[Effect of the Resample filter]{Effect of the resample filter.}
  // \label{fig:ResampleImageFilterOutput1}
  // \end{figure}
  //
  // \begin{figure}
  // \center
  // \includegraphics[width=\textwidth]{ResampleImageFilterOutput1Analysis}
  // \itkcaption[Analysis of resampling in common coordinate system]{Analysis of
  // the resample image done in a common coordinate system.}
  // \label{fig:ResampleImageFilterOutput1Analysis}
  // \end{figure}
  //
  //  Figure \ref{fig:ResampleImageFilterOutput1} illustrates the effect of
  //  this filter on a slice of MRI brain image using an affine transform
  //  containing an identity transform. Note that any analysis of the
  //  behavior of this filter must be done on the space coordinate system in
  //  millimeters, not with respect to the sampling grid in pixels. The
  //  figure shows the resulting image in the lower left quarter of the
  //  extent. This may seem odd if analyzed in terms of the image grid but is
  //  quite clear when seen with respect to space coordinates.  Figure
  //  \ref{fig:ResampleImageFilterOutput1} is particularly misleading
  //  because the images are rescaled to fit nicely on the text of this book.
  //  Figure \ref{fig:ResampleImageFilterOutput1Analysis} clarifies the
  //  situation. It shows the two same images placed on an equally-scaled
  //  coordinate system.  It becomes clear here that an identity transform is
  //  being used to map the image data, and that simply, we have requested to
  //  resample additional empty space around the image. The input image is
  //  $181 \times 217$ pixels in size and we have requested an output of $300
  //  \times 300$ pixels. In this case, the input and output images both have
  //  spacing of $1mm \times 1mm$ and origin of $(0.0,0.0)$.
  //
  //  Software Guide : EndLatex


  //  Software Guide : BeginLatex
  //
  //  Let's now set values on the transform. Note that the supplied transform
  //  represents the mapping of points from the output space to the input
  //  space.  The following code sets up a translation.
  //
  //  \index{itk::AffineTransform!Translate()}
  //
  //  Software Guide : EndLatex

  // Software Guide : BeginCodeSnippet
  TransformType::OutputVectorType translation;
  translation[0] = -30;  // X translation in millimeters
  translation[1] = -50;  // Y translation in millimeters
  transform->Translate( translation );
  // Software Guide : EndCodeSnippet


  if( exampleAction == 1 )
    {
    writer->Update();
    }


  //  Software Guide : BeginLatex
  //
  // \begin{figure}
  // \center
  // \includegraphics[width=0.44\textwidth]{BrainProtonDensitySlice}
  // \includegraphics[width=0.44\textwidth]{ResampleImageFilterOutput2}
  // \itkcaption[ResampleImageFilter with a translation by
  // $(-30,-50)$]{ResampleImageFilter with a translation by $(-30,-50)$.}
  // \label{fig:ResampleImageFilterOutput2}
  // \end{figure}
  //
  // \begin{figure}
  // \center
  // \includegraphics[width=\textwidth]{ResampleImageFilterOutput2Analysis}
  // \itkcaption[ResampleImageFilter. Analysis of a translation by
  // $(-30,-50)$]{ResampleImageFilter. Analysis of a translation by
  // $(-30,-50)$.}
  // \label{fig:ResampleImageFilterOutput2Analysis}
  // \end{figure}
  //
  // The output image resulting from the translation can be seen in Figure
  // \ref{fig:ResampleImageFilterOutput2}. Again, it is better to interpret
  // the result in a common coordinate system as illustrated in Figure
  // \ref{fig:ResampleImageFilterOutput2Analysis}.
  //
  // Probably the most important thing to keep in mind when resampling images
  // is that the transform is used to map points from the \textbf{output}
  // image space into the \textbf{input} image space. In this case, Figure
  // \ref{fig:ResampleImageFilterOutput2Analysis} shows that the translation
  // is applied to every point of the output image and the resulting position
  // is used to read the intensity from the input image. In this way, the
  // gray level of the point $P$ in the output image is taken from the point
  // $T(P)$ in the input image. Where $T$ is the transformation. In the
  // specific case of the Figure
  // \ref{fig:ResampleImageFilterOutput2Analysis}, the value of point
  // $(105,188)$ in the output image is taken from the point $(75,138)$ of
  // the input image because the transformation applied was a translation of
  // $(-30,-50)$.
  //
  //  Software Guide : EndLatex


  //  Software Guide : BeginLatex
  //
  //  It is sometimes useful to intentionally set the default output value to
  //  a distinct gray value in order to highlight the mapping of the image
  //  borders. For example, the following code sets the default external
  //  value of $100$.  The result is shown in the right side of Figure
  //  \ref{fig:ResampleImageFilterOutput3Analysis}.
  //
  //  \index{itk::ResampleImageFilter!SetDefaultPixelValue()}
  //
  //  Software Guide : EndLatex

  // Software Guide : BeginCodeSnippet
  filter->SetDefaultPixelValue( 100 );
  // Software Guide : EndCodeSnippet


  //  Software Guide : BeginLatex
  //
  // \begin{figure}
  // \center
  // \includegraphics[width=\textwidth]{ResampleImageFilterOutput3Analysis}
  // \itkcaption[ResampleImageFilter highlighting image
  // borders]{ResampleImageFilter highlighting image borders with
  // SetDefaultPixelValue().}
  // \label{fig:ResampleImageFilterOutput3Analysis}
  // \end{figure}
  //
  //  With this change we can better appreciate the effect of the previous
  //  translation transform on the image resampling. Figure
  //  \ref{fig:ResampleImageFilterOutput3Analysis} illustrates how the point
  //  $(30,50)$ of the output image gets its gray value from the point $(0,0)$
  //  of the input image.
  //
  //  Software Guide : EndLatex

  if( exampleAction == 2 )
    {
    writer->Update();
    }


  return EXIT_SUCCESS;
}