ITK  5.0.0
Insight Segmentation and Registration Toolkit
itkMorphologicalContourInterpolator.h
Go to the documentation of this file.
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 #ifndef itkMorphologicalContourInterpolator_h
19 #define itkMorphologicalContourInterpolator_h
20 
23 #include "itkExtractImageFilter.h"
24 #include "itkImageToImageFilter.h"
25 #include "itksys/hash_map.hxx"
26 
27 namespace itk
28 {
64 template< typename TImage >
66  public ImageToImageFilter< TImage, TImage >
67 {
68  template< typename T >
70 
71 public:
72  ITK_DISALLOW_COPY_AND_ASSIGN(MorphologicalContourInterpolator);
73 
78  using SliceType = Image< typename TImage::PixelType, TImage::ImageDimension - 1 >;
79 
81  itkNewMacro( Self );
82 
85 
87  itkSetMacro( Label, typename TImage::PixelType );
88 
90  itkGetMacro( Label, typename TImage::PixelType );
91 
93  itkGetConstMacro( Label, typename TImage::PixelType );
94 
96  itkSetMacro( Axis, int );
97 
99  itkGetMacro( Axis, int );
100 
102  itkGetConstMacro( Axis, int );
103 
106  itkSetMacro( HeuristicAlignment, bool );
107 
110  itkGetMacro( HeuristicAlignment, bool );
111 
114  itkGetConstMacro( HeuristicAlignment, bool );
115 
119  itkSetMacro( UseDistanceTransform, bool );
120 
124  itkGetMacro( UseDistanceTransform, bool );
125 
129  itkGetConstMacro( UseDistanceTransform, bool );
130 
133  itkSetMacro( UseCustomSlicePositions, bool );
134 
137  itkGetMacro( UseCustomSlicePositions, bool );
138 
141  itkGetConstMacro( UseCustomSlicePositions, bool );
142 
144  void
146  {
147  if ( useBall != m_UseBallStructuringElement )
148  {
149  m_UseBallStructuringElement = useBall;
150  m_ConnectedComponents->SetFullyConnected( useBall );
151  this->Modified();
152  }
153  }
155 
157  itkGetMacro( UseBallStructuringElement, bool );
158 
160  itkGetConstMacro( UseBallStructuringElement, bool );
161 
166  void
168 
170  using SliceSetType = std::set< typename TImage::IndexValueType >;
171 
173  void
175  {
176  m_LabeledSlices.clear();
177  m_LabeledSlices.resize( TImage::ImageDimension );
178  this->Modified();
179  }
181 
184  void
185  SetLabeledSliceIndices( unsigned int axis,
186  typename TImage::PixelType label,
187  const std::vector< typename TImage::IndexValueType >& indices )
188  {
189  SliceSetType sliceSet;
190  sliceSet.insert( indices.begin(), indices.end() );
191  m_LabeledSlices[axis][label] = sliceSet;
192  this->Modified();
193  }
195 
198  void
199  SetLabeledSliceIndices( unsigned int axis, typename TImage::PixelType label, const SliceSetType& indices )
200  {
201  m_LabeledSlices[axis][label] = indices;
202  this->Modified();
203  }
205 
208  GetLabeledSliceIndices( unsigned int axis, typename TImage::PixelType label )
209  {
210  return m_LabeledSlices[axis][label];
211  }
212 
213  // each label gets a set of slices in which it is present
214  using LabeledSlicesType = itksys::hash_map< typename TImage::PixelType, SliceSetType >;
215  using SliceIndicesType = std::vector< LabeledSlicesType >;
216 
220  {
221  return m_LabeledSlices;
222  }
223 
224 protected:
227  typename TImage::PixelType m_Label;
228  int m_Axis;
233  IdentifierType m_MinAlignIters; // minimum number of iterations in align method
234  IdentifierType m_MaxAlignIters; // maximum number of iterations in align method
235  IdentifierType m_ThreadCount; // for thread local instances
236  SliceIndicesType m_LabeledSlices; // one for each axis
237 
240  using FloatSliceType = Image< float, TImage::ImageDimension - 1 >;
241  using BoolSliceType = Image< bool, TImage::ImageDimension - 1 >;
242 
244  bool
245  ImagesEqual( typename BoolSliceType::Pointer& a, typename BoolSliceType::Pointer& b );
246 
248  void
249  GenerateData() override;
250 
252  void
253  InterpolateBetweenTwo( int axis,
254  TImage* out,
255  typename TImage::PixelType label,
256  typename TImage::IndexValueType i,
257  typename TImage::IndexValueType j,
258  typename SliceType::Pointer& iconn,
259  typename SliceType::Pointer& jconn );
260 
266  void
267  InterpolateAlong( int axis, TImage* out, float startProgress, float endProgress);
268 
270  void
271  Extrapolate( int axis,
272  TImage* out,
273  typename TImage::PixelType label,
274  typename TImage::IndexValueType i,
275  typename TImage::IndexValueType j,
276  typename SliceType::Pointer& iConn,
277  typename TImage::PixelType iRegionId );
278 
280  typename FloatSliceType::Pointer
281  MaurerDM( typename BoolSliceType::Pointer& inImage );
282 
285  std::vector< typename BoolSliceType::Pointer >
287  typename BoolSliceType::Pointer& end );
288 
290  typename BoolSliceType::Pointer
291  FindMedianImageDilations( typename BoolSliceType::Pointer& intersection,
292  typename BoolSliceType::Pointer& iMask,
293  typename BoolSliceType::Pointer& jMask );
294 
296  typename BoolSliceType::Pointer
297  FindMedianImageDistances( typename BoolSliceType::Pointer& intersection,
298  typename BoolSliceType::Pointer& iMask,
299  typename BoolSliceType::Pointer& jMask );
300 
302  void
303  Interpolate1to1( int axis,
304  TImage* out,
305  typename TImage::PixelType label,
306  typename TImage::IndexValueType i,
307  typename TImage::IndexValueType j,
308  typename SliceType::Pointer& iConn,
309  typename TImage::PixelType iRegionId,
310  typename SliceType::Pointer& jConn,
311  typename TImage::PixelType jRegionId,
312  const typename SliceType::IndexType& translation,
313  bool recursive );
314 
315  using PixelList = std::vector< typename TImage::PixelType >;
316 
318  void
319  Interpolate1toN( int axis,
320  TImage* out,
321  typename TImage::PixelType label,
322  typename TImage::IndexValueType i,
323  typename TImage::IndexValueType j,
324  typename SliceType::Pointer& iConn,
325  typename TImage::PixelType iRegionId,
326  typename SliceType::Pointer& jConn,
327  const PixelList& jRegionIds,
328  const typename SliceType::IndexType& translation );
329 
331  typename SliceType::Pointer
332  TranslateImage( typename SliceType::Pointer& image,
333  const typename SliceType::IndexType& translation,
334  typename SliceType::RegionType newRegion );
335 
339  CardinalSymmetricDifference( typename BoolSliceType::Pointer& shape1, typename BoolSliceType::Pointer& shape2 );
340 
342  void
343  AllocateOutputs() override;
344 
346  typename SliceType::IndexType
347  Centroid( typename SliceType::Pointer& conn, const PixelList& regionIds );
348 
351  void
352  IntersectionRegions( const typename SliceType::IndexType& translation,
353  typename SliceType::RegionType& iRegion,
354  typename SliceType::RegionType& jRegion );
355 
358  Intersection( typename SliceType::Pointer& iConn,
359  typename TImage::PixelType iRegionId,
360  typename SliceType::Pointer& jConn,
361  const PixelList& jRegionIds,
362  const typename SliceType::IndexType& translation );
363 
365  typename SliceType::IndexType
366  Align( typename SliceType::Pointer& iConn,
367  typename TImage::PixelType iRegionId,
368  typename SliceType::Pointer& jConn,
369  const PixelList& jRegionIds );
370 
371  using BoundingBoxesType = itksys::hash_map< typename TImage::PixelType, typename TImage::RegionType >;
372  BoundingBoxesType m_BoundingBoxes; // bounding box for each label
373 
375  typename SliceType::RegionType
377 
381  template< typename T2 >
382  void
383  ExpandRegion( typename T2::RegionType& region, const typename T2::IndexType& index );
384 
386  typename SliceType::Pointer
387  RegionedConnectedComponents( const typename TImage::RegionType& region,
388  typename TImage::PixelType label,
389  IdentifierType& objectCount );
390 
392  typename BoolSliceType::Pointer
393  Dilate1( typename BoolSliceType::Pointer& seed, typename BoolSliceType::Pointer& mask );
394 
397 
400 
403 };
404 } // namespace itk
405 
406 #ifndef ITK_MANUAL_INSTANTIATION
407 #include "itkMorphologicalContourInterpolator.hxx"
408 #endif
409 
410 #endif // itkMorphologicalContourInterpolator_h
SliceType::Pointer RegionedConnectedComponents(const typename TImage::RegionType &region, typename TImage::PixelType label, IdentifierType &objectCount)
Light weight base class for most itk classes.
Decrease the image size by cropping the image to the selected region bounds.
void IntersectionRegions(const typename SliceType::IndexType &translation, typename SliceType::RegionType &iRegion, typename SliceType::RegionType &jRegion)
void InterpolateBetweenTwo(int axis, TImage *out, typename TImage::PixelType label, typename TImage::IndexValueType i, typename TImage::IndexValueType j, typename SliceType::Pointer &iconn, typename SliceType::Pointer &jconn)
SliceSetType GetLabeledSliceIndices(unsigned int axis, typename TImage::PixelType label)
std::set< typename TImage::IndexValueType > SliceSetType
void SetLabeledSliceIndices(unsigned int axis, typename TImage::PixelType label, const std::vector< typename TImage::IndexValueType > &indices)
void Interpolate1toN(int axis, TImage *out, typename TImage::PixelType label, typename TImage::IndexValueType i, typename TImage::IndexValueType j, typename SliceType::Pointer &iConn, typename TImage::PixelType iRegionId, typename SliceType::Pointer &jConn, const PixelList &jRegionIds, const typename SliceType::IndexType &translation)
std::vector< typename TImage::PixelType > PixelList
Binarize an input image by thresholding.
bool ImagesEqual(typename BoolSliceType::Pointer &a, typename BoolSliceType::Pointer &b)
std::vector< typename BoolSliceType::Pointer > GenerateDilationSequence(typename BoolSliceType::Pointer &begin, typename BoolSliceType::Pointer &end)
SliceType::IndexType Centroid(typename SliceType::Pointer &conn, const PixelList &regionIds)
SizeValueType IdentifierType
Definition: itkIntTypes.h:87
SliceType::Pointer TranslateImage(typename SliceType::Pointer &image, const typename SliceType::IndexType &translation, typename SliceType::RegionType newRegion)
itksys::hash_map< typename TImage::PixelType, typename TImage::RegionType > BoundingBoxesType
signed long IndexValueType
Definition: itkIntTypes.h:90
BoolSliceType::Pointer FindMedianImageDilations(typename BoolSliceType::Pointer &intersection, typename BoolSliceType::Pointer &iMask, typename BoolSliceType::Pointer &jMask)
typename Superclass::RegionType RegionType
Definition: itkImage.h:139
itksys::hash_map< typename TImage::PixelType, SliceSetType > LabeledSlicesType
Interpolates contours between slices. Based on a paper by Albu et al.
BoolSliceType::Pointer Dilate1(typename BoolSliceType::Pointer &seed, typename BoolSliceType::Pointer &mask)
virtual void Modified() const
void Extrapolate(int axis, TImage *out, typename TImage::PixelType label, typename TImage::IndexValueType i, typename TImage::IndexValueType j, typename SliceType::Pointer &iConn, typename TImage::PixelType iRegionId)
BoolSliceType::Pointer FindMedianImageDistances(typename BoolSliceType::Pointer &intersection, typename BoolSliceType::Pointer &iMask, typename BoolSliceType::Pointer &jMask)
typename Superclass::IndexType IndexType
Definition: itkImage.h:121
void ExpandRegion(typename T2::RegionType &region, const typename T2::IndexType &index)
Base class for filters that take an image as input and produce an image as output.
void SetLabeledSliceIndices(unsigned int axis, typename TImage::PixelType label, const SliceSetType &indices)
IdentifierType CardinalSymmetricDifference(typename BoolSliceType::Pointer &shape1, typename BoolSliceType::Pointer &shape2)
void Interpolate1to1(int axis, TImage *out, typename TImage::PixelType label, typename TImage::IndexValueType i, typename TImage::IndexValueType j, typename SliceType::Pointer &iConn, typename TImage::PixelType iRegionId, typename SliceType::Pointer &jConn, typename TImage::PixelType jRegionId, const typename SliceType::IndexType &translation, bool recursive)
FloatSliceType::Pointer MaurerDM(typename BoolSliceType::Pointer &inImage)
SliceType::IndexType Align(typename SliceType::Pointer &iConn, typename TImage::PixelType iRegionId, typename SliceType::Pointer &jConn, const PixelList &jRegionIds)
void InterpolateAlong(int axis, TImage *out, float startProgress, float endProgress)
SliceType::RegionType BoundingBox(itk::SmartPointer< SliceType > image)
Templated n-dimensional image class.
Definition: itkImage.h:75
IdentifierType Intersection(typename SliceType::Pointer &iConn, typename TImage::PixelType iRegionId, typename SliceType::Pointer &jConn, const PixelList &jRegionIds, const typename SliceType::IndexType &translation)
Label the objects in a binary image.