ITK  6.0.0
Insight Toolkit
itkContourExtractor2DImageFilter.h
Go to the documentation of this file.
1/*=========================================================================
2 *
3 * Copyright NumFOCUS
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 * https://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#ifndef itkContourExtractor2DImageFilter_h
19#define itkContourExtractor2DImageFilter_h
20
21#include <deque>
22#include <list>
23#include <unordered_map>
24#include "itkConceptChecking.h"
26#include "itkImageRegionRange.h"
28#include "itkIndexRange.h"
30
31namespace itk
32{
107template <typename TInputImage>
108class ITK_TEMPLATE_EXPORT ContourExtractor2DImageFilter
109 : public ImageToPathFilter<TInputImage, PolyLineParametricPath<2>>
110{
111public:
112 ITK_DISALLOW_COPY_AND_MOVE(ContourExtractor2DImageFilter);
113
115 static constexpr unsigned int InputImageDimension{ TInputImage::ImageDimension };
116
118 using InputImageType = TInputImage;
120
126
128 itkNewMacro(Self);
129
131 itkOverrideGetNameOfClassMacro(ContourExtractor2DImageFilter);
132
137 using InputOffsetType = typename InputImageType::OffsetType;
138 using InputPixelType = typename InputImageType::PixelType;
144 using VertexValueType = typename VertexType::ValueType;
145
148
155
158 itkSetMacro(ReverseContourOrientation, bool);
159 itkGetConstReferenceMacro(ReverseContourOrientation, bool);
160 itkBooleanMacro(ReverseContourOrientation);
166 itkSetMacro(VertexConnectHighPixels, bool);
167 itkGetConstReferenceMacro(VertexConnectHighPixels, bool);
168 itkBooleanMacro(VertexConnectHighPixels);
172 itkSetMacro(LabelContours, bool);
173 itkGetConstReferenceMacro(LabelContours, bool);
174 itkBooleanMacro(LabelContours);
179 void
181
182 itkGetConstReferenceMacro(RequestedRegion, InputRegionType);
183 void
185
188 itkSetMacro(ContourValue, InputRealType);
189 itkGetConstReferenceMacro(ContourValue, InputRealType);
192#ifdef ITK_USE_CONCEPT_CHECKING
193 // Begin concept checking
195 itkConceptMacro(InputPixelTypeComparable, (Concept::Comparable<InputPixelType>));
196 itkConceptMacro(InputHasPixelTraitsCheck, (Concept::HasPixelTraits<InputPixelType>));
197 itkConceptMacro(InputHasNumericTraitsCheck, (Concept::HasNumericTraits<InputPixelType>));
198 // End concept checking
199#endif
200
201protected:
203 ~ContourExtractor2DImageFilter() override = default;
204 void
205 PrintSelf(std::ostream & os, Indent indent) const override;
206
207 void
208 GenerateData() override;
209
213 void
215
216private:
217 using LabelsContainer = std::vector<InputPixelType>;
218 using LabelsConstIterator = typename LabelsContainer::const_iterator;
219 using LabelsIterator = typename LabelsContainer::iterator;
220
221 InputRealType m_ContourValue{};
222 bool m_ReverseContourOrientation{};
223 bool m_VertexConnectHighPixels{};
224 bool m_LabelContours{};
225 bool m_UseCustomRegion{};
226 InputRegionType m_RequestedRegion{};
227
228 // Represent each contour as deque of vertices to facilitate addition of
229 // nodes at beginning or end. At the end of the processing, we will copy
230 // the contour into a PolyLineParametricPath.
231 // We subclass the deque to store an additional bit of information: an
232 // identification number for each growing contour. We use this number so
233 // that when it becomes necessary to merge two growing contours, we can
234 // merge the newer one into the older one. This helps because then we can
235 // guarantee that the output contour list is ordered from left to right,
236 // top to bottom (in terms of the first pixel of the contour encountered
237 // by the marching squares). Currently we make no guarantees that this
238 // pixel is the first pixel in the contour list, just that the contours
239 // are so ordered in the output. Ensuring this latter condition (first
240 // pixel traversed = first pixel in contour) would be possible by either
241 // changing the merging rules, which would make the contouring operation
242 // slower, or by storing additional data as to which pixel was first.
243 class ContourType : public std::deque<VertexType>
244 {
245 public:
246 unsigned int m_ContourNumber;
247 };
248
249 // Store all the growing contours in a list. We may need to delete contours
250 // from anywhere in the sequence (when we merge them together), so we need to
251 // use a list instead of a vector or similar.
252 using ContourConstIterator = typename ContourType::const_iterator;
253 using ContourContainerType = std::list<ContourType>;
254 using ContourContainerIterator = typename ContourContainerType::iterator;
255 using ContourContainerConstIterator = typename ContourContainerType::const_iterator;
256
257 // We use a hash to associate the endpoints of each contour with the
258 // contour itself. This makes it easy to look up which contour we should add
259 // a new arc to.
260
261 // We can't store the contours themselves in the hashtable because we
262 // need to have two tables (one to hash from beginpoint -> contour and one
263 // for endpoint -> contour), and sometimes will remove a contour from the
264 // tables (if it has been closed or merged with another contour). Because
265 // sometimes we will need to merge contours, we need to be able to quickly
266 // remove contours from our list when they have been merged into
267 // another. Thus, we store an iterator pointing to the contour in the list.
268
270 {
271 using CoordinateType = typename VertexType::CoordRepType;
272 inline size_t
273 operator()(const VertexType & v) const noexcept
274 {
275 return std::hash<CoordinateType>{}(v[0]) ^ (std::hash<CoordinateType>{}(v[1]) << 1);
276 }
277 };
278 using VertexToContourContainerIteratorMap = std::unordered_map<VertexType, ContourContainerIterator, VertexHash>;
279 using VertexToContourContainerIteratorMapIterator = typename VertexToContourContainerIteratorMap::iterator;
280 using VertexToContourContainerIteratorMapConstIterator = typename VertexToContourContainerIteratorMap::const_iterator;
281 using VertexToContourContainerIteratorMapKeyValuePair = typename VertexToContourContainerIteratorMap::value_type;
282
283 // Subroutine to create contours for a single label.
284 void
286 const InputImageType * input,
287 const InputRegionType extendedRegion,
288 SizeValueType totalNumberOfPixels,
289 ContourContainerType & contoursOutput);
290
291 // Subroutine to handle the case that the supplied values are
292 // labels, which are *not* compared to a contour value.
293 void
295
298 InputPixelType toValue,
299 InputIndexType fromIndex,
300 InputOffsetType toOffset);
301
303 {
307 SizeValueType m_NumberOfContoursCreated = 0;
308 };
309
310 void
311 AddSegment(const VertexType from, const VertexType to, ContourData & contourData);
312
313 void
314 FillOutputs(const std::vector<InputPixelType> & allLabels,
315 std::unordered_map<InputPixelType, ContourContainerType> & labelsContoursOutput);
316
317 InputPixelType m_UnusedLabel{};
318};
319} // end namespace itk
320
321#ifndef ITK_MANUAL_INSTANTIATION
322# include "itkContourExtractor2DImageFilter.hxx"
323#endif
324
325#endif
A templated class holding a point in n-Dimensional image space.
Computes a list of PolyLineParametricPath objects from the contours in a 2D image.
typename InputImageType::IndexType InputIndexType
typename InputImageType::PixelType InputPixelType
typename ContourContainerType::const_iterator ContourContainerConstIterator
void SetRequestedRegion(const InputRegionType region)
~ContourExtractor2DImageFilter() override=default
typename VertexToContourContainerIteratorMap::value_type VertexToContourContainerIteratorMapKeyValuePair
typename VertexToContourContainerIteratorMap::iterator VertexToContourContainerIteratorMapIterator
void GenerateInputRequestedRegion() override
typename LabelsContainer::iterator LabelsIterator
VertexType InterpolateContourPosition(InputPixelType fromValue, InputPixelType toValue, InputIndexType fromIndex, InputOffsetType toOffset)
typename NumericTraits< InputPixelType >::RealType InputRealType
typename InputImageType::SizeType InputSizeType
typename OutputPathType::VertexListType VertexListType
void AddSegment(const VertexType from, const VertexType to, ContourData &contourData)
void PrintSelf(std::ostream &os, Indent indent) const override
typename OutputPathType::Pointer OutputPathPointer
typename InputImageType::Pointer InputImagePointer
typename LabelsContainer::const_iterator LabelsConstIterator
typename InputImageType::OffsetType InputOffsetType
typename VertexListType::ConstPointer VertexListConstPointer
typename ContourContainerType::iterator ContourContainerIterator
typename VertexToContourContainerIteratorMap::const_iterator VertexToContourContainerIteratorMapConstIterator
typename OutputPathType::VertexType VertexType
typename InputImageType::RegionType InputRegionType
void FillOutputs(const std::vector< InputPixelType > &allLabels, std::unordered_map< InputPixelType, ContourContainerType > &labelsContoursOutput)
std::unordered_map< VertexType, ContourContainerIterator, VertexHash > VertexToContourContainerIteratorMap
typename ContourType::const_iterator ContourConstIterator
void CreateSingleContour(InputPixelType label, const InputImageType *input, const InputRegionType extendedRegion, SizeValueType totalNumberOfPixels, ContourContainerType &contoursOutput)
A multi-dimensional iterator templated over image type that walks a region of pixels.
A multi-dimensional iterator templated over image type that walks a region of pixels.
Base class for filters that take an image as input and produce an path as output.
Control indentation during Print() invocation.
Definition: itkIndent.h:50
Light weight base class for most itk classes.
Define additional traits for native types such as int or float.
Represent a path of line segments through ND Space.
Define a front-end to the STL "vector" container that conforms to the IndexedContainerInterface.
SmartPointer< const Self > ConstPointer
SmartPointer< Self > Pointer
#define itkConceptMacro(name, concept)
The "itk" namespace contains all Insight Segmentation and Registration Toolkit (ITK) classes....
unsigned long SizeValueType
Definition: itkIntTypes.h:86
vcl_size_t operator()(const VertexType &v) const noexcept