1 /*=========================================================================
2 *
3 * Copyright Insight Software Consortium
4 *
5 * Licensed under the Apache License, Version 2.0 (the "License");
6 * you may not use this file except in compliance with the License.
7 * You may obtain a copy of the License at
8 *
9 * http://www.apache.org/licenses/LICENSE-2.0.txt
10 *
11 * Unless required by applicable law or agreed to in writing, software
12 * distributed under the License is distributed on an "AS IS" BASIS,
13 * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
14 * See the License for the specific language governing permissions and
15 * limitations under the License.
16 *
17 *=========================================================================*/
18
19 // Software Guide : BeginCommandLineArgs
20 // INPUTS: {VisibleWomanEyeSlice.png}
21 // OUTPUTS: {ImageRegionIteratorWithIndexOutput.png}
22 // Software Guide : EndCommandLineArgs
23
24 // Software Guide : BeginLatex
25 //
26 // \index{Iterators!speed}
27 // The ``WithIndex'' family of iterators was designed for algorithms that
28 // use both the value and the location of image pixels in calculations. Unlike
29 // \doxygen{ImageRegionIterator}, which calculates an index only when
30 // asked for, \doxygen{ImageRegionIteratorWithIndex} maintains its
31 // index location as a member variable that is updated during the increment or
32 // decrement process. Iteration speed is penalized, but the index queries are
33 // more efficient.
34 //
35 // \index{itk::ImageRegionIteratorWithIndex!example of using|(}
36 //
37 // The following example illustrates the use of
38 // ImageRegionIteratorWithIndex. The algorithm mirrors
39 // a 2D image across its $x$-axis (see \doxygen{FlipImageFilter} for an ND
40 // version). The algorithm makes extensive use of the \code{GetIndex()}
41 // method.
42 //
43 // We start by including the proper header file.
44 //
45 // Software Guide : EndLatex
46
47 #include "itkImage.h"
48 #include "itkRGBPixel.h"
49 // Software Guide : BeginCodeSnippet
50 #include "itkImageRegionIteratorWithIndex.h"
51 // Software Guide : EndCodeSnippet
52 #include "itkImageFileReader.h"
53 #include "itkImageFileWriter.h"
54
main(int argc,char * argv[])55 int main( int argc, char *argv[] )
56 {
57 // Verify the number of parameters on the command line.
58 if ( argc < 3 )
59 {
60 std::cerr << "Missing parameters. " << std::endl;
61 std::cerr << "Usage: " << std::endl;
62 std::cerr << argv[0]
63 << " inputImageFile outputImageFile"
64 << std::endl;
65 return EXIT_FAILURE;
66 }
67
68 // Software Guide : BeginLatex
69 //
70 // For this example, we will use an RGB pixel type so that we can process color
71 // images. Like most other ITK image iterator,
72 // ImageRegionIteratorWithIndex class expects the image type as its
73 // single template parameter.
74 //
75 // Software Guide : EndLatex
76
77 // Software Guide : BeginCodeSnippet
78 constexpr unsigned int Dimension = 2;
79
80 using RGBPixelType = itk::RGBPixel< unsigned char >;
81 using ImageType = itk::Image< RGBPixelType, Dimension >;
82
83 using IteratorType = itk::ImageRegionIteratorWithIndex< ImageType >;
84 // Software Guide : EndCodeSnippet
85
86 using ReaderType = itk::ImageFileReader< ImageType >;
87 using WriterType = itk::ImageFileWriter< ImageType >;
88
89 ImageType::ConstPointer inputImage;
90 ReaderType::Pointer reader = ReaderType::New();
91 reader->SetFileName( argv[1] );
92 try
93 {
94 reader->Update();
95 inputImage = reader->GetOutput();
96 }
97 catch ( itk::ExceptionObject &err)
98 {
99 std::cerr << "ExceptionObject caught !" << std::endl;
100 std::cerr << err << std::endl;
101 return EXIT_FAILURE;
102 }
103
104 // Software Guide : BeginLatex
105 //
106 // An \code{ImageType} smart pointer called \code{inputImage} points to the
107 // output of the image reader. After updating the image reader, we can
108 // allocate an output image of the same size, spacing, and origin as the
109 // input image.
110 //
111 // Software Guide : EndLatex
112
113 // Software Guide : BeginCodeSnippet
114 ImageType::Pointer outputImage = ImageType::New();
115 outputImage->SetRegions( inputImage->GetRequestedRegion() );
116 outputImage->CopyInformation( inputImage );
117 outputImage->Allocate();
118 // Software Guide : EndCodeSnippet
119
120 // Software Guide : BeginLatex
121 //
122 // Next we create the iterator that walks the output image. This algorithm
123 // requires no iterator for the input image.
124 //
125 // Software Guide : EndLatex
126
127 // Software Guide : BeginCodeSnippet
128 IteratorType outputIt( outputImage, outputImage->GetRequestedRegion() );
129 // Software Guide : EndCodeSnippet
130
131 // Software Guide: BeginLatex
132 //
133 // This axis flipping algorithm works by iterating through the output image,
134 // querying the iterator for its index, and copying the value from the input
135 // at an index mirrored across the $x$-axis.
136 //
137 // Software Guide : EndLatex
138
139 // Software Guide : BeginCodeSnippet
140 ImageType::IndexType requestedIndex =
141 outputImage->GetRequestedRegion().GetIndex();
142 ImageType::SizeType requestedSize =
143 outputImage->GetRequestedRegion().GetSize();
144
145 for ( outputIt.GoToBegin(); !outputIt.IsAtEnd(); ++outputIt)
146 {
147 ImageType::IndexType idx = outputIt.GetIndex();
148 idx[0] = requestedIndex[0] + requestedSize[0] - 1 - idx[0];
149 outputIt.Set( inputImage->GetPixel(idx) );
150 }
151 // Software Guide : EndCodeSnippet
152
153 WriterType::Pointer writer = WriterType::New();
154 writer->SetFileName( argv[2] );
155 writer->SetInput(outputImage);
156 try
157 {
158 writer->Update();
159 }
160 catch ( itk::ExceptionObject &err)
161 {
162 std::cerr << "ExceptionObject caught !" << std::endl;
163 std::cerr << err << std::endl;
164 return EXIT_FAILURE;
165 }
166
167 // Software Guide : BeginLatex
168 //
169 // Let's run this example on the image \code{VisibleWomanEyeSlice.png} found in
170 // the \code{Examples/Data} directory.
171 // Figure~\ref{fig:ImageRegionIteratorWithIndexExample} shows how the original
172 // image has been mirrored across its $x$-axis in the output.
173 //
174 // \begin{figure} \center
175 // \includegraphics[width=0.44\textwidth]{VisibleWomanEyeSlice}
176 // \includegraphics[width=0.44\textwidth]{ImageRegionIteratorWithIndexOutput}
177 // \itkcaption[Using the ImageRegionIteratorWithIndex]{Results of using
178 // ImageRegionIteratorWithIndex to mirror an image across an axis. The original
179 // image is shown at left. The mirrored output is shown at right.}
180 // \label{fig:ImageRegionIteratorWithIndexExample}
181 // \end{figure}
182 //
183 // \index{itk::ImageRegionIteratorWithIndex!example of using|)}
184 //
185 // Software Guide : EndLatex
186
187 return EXIT_SUCCESS;
188 }
189