forked from InsightSoftwareConsortium/ITK
-
Notifications
You must be signed in to change notification settings - Fork 0
/
NeighborhoodConnectedImageFilter.cxx
332 lines (264 loc) · 11.5 KB
/
NeighborhoodConnectedImageFilter.cxx
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
/*=========================================================================
Program: Insight Segmentation & Registration Toolkit
Module: NeighborhoodConnectedImageFilter.cxx
Language: C++
Date: $Date$
Version: $Revision$
Copyright (c) Insight Software Consortium. All rights reserved.
See ITKCopyright.txt or http://www.itk.org/HTML/Copyright.htm for details.
This software is distributed WITHOUT ANY WARRANTY; without even
the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR
PURPOSE. See the above copyright notices for more information.
=========================================================================*/
#if defined(_MSC_VER)
#pragma warning ( disable : 4786 )
#endif
#ifdef __BORLANDC__
#define ITK_LEAN_AND_MEAN
#endif
// Software Guide : BeginLatex
//
// The following example illustrates the use of the
// \doxygen{NeighborhoodConnectedImageFilter}. This filter is a close variant
// of the \doxygen{ConnectedThresholdImageFilter}. On one hand, the
// ConnectedThresholdImageFilter accepts a pixel in the region if its intensity
// is in the interval defined by two user-provided threshold values. The
// NeighborhoodConnectedImageFilter, on the other hand, will only accept a
// pixel if \textbf{all} its neighbors have intensities that fit in the
// interval. The size of the neighborhood to be considered around each pixel is
// defined by a user-provided integer radius.
//
// The reason for considering the neighborhood intensities instead of only the
// current pixel intensity is that small structures are less likely to be
// accepted in the region. The operation of this filter is equivalent to
// applying the ConnectedThresholdImageFilter followed by mathematical
// morphology erosion using a structuring element of the same shape as
// the neighborhood provided to the NeighborhoodConnectedImageFilter.
//
// Software Guide : EndLatex
// Software Guide : BeginCodeSnippet
#include "itkNeighborhoodConnectedImageFilter.h"
// Software Guide : EndCodeSnippet
#include "itkImage.h"
#include "itkCastImageFilter.h"
// Software Guide : BeginLatex
//
// The \doxygen{CurvatureFlowImageFilter} is used here to smooth the image
// while preserving edges.
//
// Software Guide : EndLatex
// Software Guide : BeginCodeSnippet
#include "itkCurvatureFlowImageFilter.h"
// Software Guide : EndCodeSnippet
#include "itkImageFileReader.h"
#include "itkImageFileWriter.h"
int main( int argc, char *argv[] )
{
if( argc < 7 )
{
std::cerr << "Missing Parameters " << std::endl;
std::cerr << "Usage: " << argv[0];
std::cerr << " inputImage outputImage seedX seedY lowerThreshold upperThreshold" << std::endl;
return 1;
}
// Software Guide : BeginLatex
//
// We now define the image type using a particular pixel type and image
// dimension. In this case the \code{float} type is used for the pixels due
// to the requirements of the smoothing filter.
//
// Software Guide : EndLatex
// Software Guide : BeginCodeSnippet
typedef float InternalPixelType;
const unsigned int Dimension = 2;
typedef itk::Image< InternalPixelType, Dimension > InternalImageType;
// Software Guide : EndCodeSnippet
typedef unsigned char OutputPixelType;
typedef itk::Image< OutputPixelType, Dimension > OutputImageType;
typedef itk::CastImageFilter< InternalImageType, OutputImageType >
CastingFilterType;
CastingFilterType::Pointer caster = CastingFilterType::New();
// We instantiate reader and writer types
//
typedef itk::ImageFileReader< InternalImageType > 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
//
// The smoothing filter type is instantiated using the image type as
// a template parameter.
//
// Software Guide : EndLatex
// Software Guide : BeginCodeSnippet
typedef itk::CurvatureFlowImageFilter<InternalImageType, InternalImageType>
CurvatureFlowImageFilterType;
// Software Guide : EndCodeSnippet
// Software Guide : BeginLatex
//
// Then, the filter is created by invoking the \code{New()} method and
// assigning the result to a \doxygen{SmartPointer}.
//
// Software Guide : EndLatex
// Software Guide : BeginCodeSnippet
CurvatureFlowImageFilterType::Pointer smoothing =
CurvatureFlowImageFilterType::New();
// Software Guide : EndCodeSnippet
// Software Guide : BeginLatex
//
// We now declare the type of the region growing filter. In this case it is
// the NeighborhoodConnectedImageFilter.
//
// Software Guide : EndLatex
// Software Guide : BeginCodeSnippet
typedef itk::NeighborhoodConnectedImageFilter<InternalImageType,
InternalImageType > ConnectedFilterType;
// Software Guide : EndCodeSnippet
// Software Guide : BeginLatex
//
// One filter of this class is constructed using the \code{New()} method.
//
// Software Guide : EndLatex
// Software Guide : BeginCodeSnippet
ConnectedFilterType::Pointer neighborhoodConnected = ConnectedFilterType::New();
// Software Guide : EndCodeSnippet
// Software Guide : BeginLatex
//
// Now it is time to create a simple, linear data processing pipeline. A
// file reader is added at the beginning of the pipeline and a cast
// filter and writer are added at the end. The cast filter is required
// to convert \code{float} pixel types to integer types since only a
// few image file formats support \code{float} types.
//
// Software Guide : EndLatex
// Software Guide : BeginCodeSnippet
smoothing->SetInput( reader->GetOutput() );
neighborhoodConnected->SetInput( smoothing->GetOutput() );
caster->SetInput( neighborhoodConnected->GetOutput() );
writer->SetInput( caster->GetOutput() );
// Software Guide : EndCodeSnippet
// Software Guide : BeginLatex
//
// The CurvatureFlowImageFilter requires a couple of parameters to
// be defined. The following are typical values for $2D$ images. However
// they may have to be adjusted depending on the amount of noise present in
// the input image.
//
// Software Guide : EndLatex
// Software Guide : BeginCodeSnippet
smoothing->SetNumberOfIterations( 5 );
smoothing->SetTimeStep( 0.125 );
// Software Guide : EndCodeSnippet
// Software Guide : BeginLatex
//
// The NeighborhoodConnectedImageFilter requires that two main parameters
// are specified. They are the lower and upper thresholds of the interval
// in which intensity values must fall to be included in the
// region. Setting these two values too close will not allow enough
// flexibility for the region to grow. Setting them too far apart will
// result in a region that engulfs the image.
//
// \index{itk::NeighborhoodConnectedImageFilter!SetLower()}
// \index{itk::NeighborhoodConnectedImageFilter!SetUppder()}
//
// Software Guide : EndLatex
const InternalPixelType lowerThreshold = atof( argv[5] );
const InternalPixelType upperThreshold = atof( argv[6] );
// Software Guide : BeginCodeSnippet
neighborhoodConnected->SetLower( lowerThreshold );
neighborhoodConnected->SetUpper( upperThreshold );
// Software Guide : EndCodeSnippet
// Software Guide : BeginLatex
//
// Here, we add the crucial parameter that defines the neighborhood size
// used to determine whether a pixel lies in the region. The larger the
// neighborhood, the more stable this filter will be against noise in the
// input image, but also the longer the computing time will be. Here we
// select a filter of radius $2$ along each dimension. This results in a
// neighborhood of $5 \times 5$ pixels.
//
// Software Guide : EndLatex
// Software Guide : BeginCodeSnippet
InternalImageType::SizeType radius;
radius[0] = 2; // two pixels along X
radius[1] = 2; // two pixels along Y
neighborhoodConnected->SetRadius( radius );
// Software Guide : EndCodeSnippet
// Software Guide : BeginLatex
//
// As in the ConnectedThresholdImageFilter we must now provide the
// intensity value to be used for the output pixels accepted in the region
// and at least one seed point to define the initial region.
//
// \index{itk::NeighborhoodConnectedImageFilter!SetSeed()}
// \index{itk::NeighborhoodConnectedImageFilter!SetReplaceValue()}
//
// Software Guide : EndLatex
InternalImageType::IndexType index;
index[0] = atoi( argv[3] );
index[1] = atoi( argv[4] );
// Software Guide : BeginCodeSnippet
neighborhoodConnected->SetSeed( index );
neighborhoodConnected->SetReplaceValue( 255 );
// Software Guide : EndCodeSnippet
// Software Guide : BeginLatex
//
// The invocation of the \code{Update()} method on the writer triggers the
// execution of the pipeline. It is usually wise to put update calls in a
// \code{try/catch} block in case errors occur and exceptions are thrown.
//
// Software Guide : EndLatex
// Software Guide : BeginCodeSnippet
try
{
writer->Update();
}
catch( itk::ExceptionObject & excep )
{
std::cerr << "Exception caught !" << std::endl;
std::cerr << excep << std::endl;
}
// Software Guide : EndCodeSnippet
// Software Guide : BeginLatex
//
// Now we'll run this example using the image
// \code{BrainProtonDensitySlice.png} as input available from the
// directory \code{Examples/Data}. We can easily segment the major
// anatomical structures by providing seeds in the appropriate locations
// and defining values for the lower and upper thresholds. For example
//
// \begin{center}
// \begin{tabular}{|l|c|c|c|c|}
// \hline
// Structure & Seed Index & Lower & Upper & Output Image \\ \hline
// White matter & $(60,116)$ & 150 & 180 & Second from left in Figure \ref{fig:NeighborhoodConnectedImageFilterOutput} \\ \hline
// Ventricle & $(81,112)$ & 210 & 250 & Third from left in Figure \ref{fig:NeighborhoodConnectedImageFilterOutput} \\ \hline
// Gray matter & $(107,69)$ & 180 & 210 & Fourth from left in Figure \ref{fig:NeighborhoodConnectedImageFilterOutput} \\ \hline
// \end{tabular}
// \end{center}
//
// \begin{figure} \center
// \includegraphics[width=0.24\textwidth]{BrainProtonDensitySlice.eps}
// \includegraphics[width=0.24\textwidth]{NeighborhoodConnectedImageFilterOutput1.eps}
// \includegraphics[width=0.24\textwidth]{NeighborhoodConnectedImageFilterOutput2.eps}
// \includegraphics[width=0.24\textwidth]{NeighborhoodConnectedImageFilterOutput3.eps}
// \itkcaption[NeighborhoodConnected segmentation results ]{Segmentation results
// of the NeighborhoodConnectedImageFilter for various seed points.}
// \label{fig:NeighborhoodConnectedImageFilterOutput}
// \end{figure}
//
// As with the ConnectedThresholdImageFilter, several seeds could
// be provided to the filter by using the \code{AddSeed()} method.
// Compare the output of Figure
// \ref{fig:NeighborhoodConnectedImageFilterOutput} with those of Figure
// \ref{fig:ConnectedThresholdOutput} produced by the
// ConnectedThresholdImageFilter. You may want to play with the
// value of the neighborhood radius and see how it affect the smoothness of
// the segmented object borders, the size of the segmented region and how
// much that costs in computing time.
//
// Software Guide : EndLatex
return 0;
}