Slicer  4.10
Slicer is a multi-platform, free and open source software package for visualization and medical image computing
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 
21 #include "itkBinaryThresholdImageFilter.h"
22 #include "itkConnectedComponentImageFilter.h"
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:
74  typedef ImageToImageFilter< TImage, TImage > Superclass;
75  typedef SmartPointer< Self > Pointer;
76  typedef Image< typename TImage::PixelType, TImage::ImageDimension - 1 > SliceType;
77 
79  itkNewMacro( Self );
80 
82  itkTypeMacro( MorphologicalContourInterpolator, ImageToImageFilter );
83 
85  itkSetMacro( Label, typename TImage::PixelType );
86 
88  itkGetMacro( Label, typename TImage::PixelType );
89 
91  itkGetConstMacro( Label, typename TImage::PixelType );
92 
94  itkSetMacro( Axis, int );
95 
97  itkGetMacro( Axis, int );
98 
100  itkGetConstMacro( Axis, int );
101 
104  itkSetMacro( HeuristicAlignment, bool );
105 
108  itkGetMacro( HeuristicAlignment, bool );
109 
112  itkGetConstMacro( HeuristicAlignment, bool );
113 
117  itkSetMacro( UseDistanceTransform, bool );
118 
122  itkGetMacro( UseDistanceTransform, bool );
123 
127  itkGetConstMacro( UseDistanceTransform, bool );
128 
131  itkSetMacro( UseCustomSlicePositions, bool );
132 
135  itkGetMacro( UseCustomSlicePositions, bool );
136 
139  itkGetConstMacro( UseCustomSlicePositions, bool );
140 
142  void
144  {
145  if ( useBall != m_UseBallStructuringElement )
146  {
147  m_UseBallStructuringElement = useBall;
148  m_ConnectedComponents->SetFullyConnected( useBall );
149  this->Modified();
150  }
151  }
152 
154  itkGetMacro( UseBallStructuringElement, bool );
155 
157  itkGetConstMacro( UseBallStructuringElement, bool );
158 
163  void
165 
167  typedef std::set< typename TImage::IndexValueType > SliceSetType;
168 
170  void
172  {
173  m_LabeledSlices.clear();
174  m_LabeledSlices.resize( TImage::ImageDimension );
175  this->Modified();
176  }
177 
180  void
181  SetLabeledSliceIndices( unsigned int axis,
182  typename TImage::PixelType label,
183  const std::vector< typename TImage::IndexValueType >& indices )
184  {
185  m_LabeledSlices[axis][label] = SliceSetType().insert( indices.begin(), indices.end() );
186  this->Modified();
187  }
188 
191  void
192  SetLabeledSliceIndices( unsigned int axis, typename TImage::PixelType label, const SliceSetType& indices )
193  {
194  m_LabeledSlices[axis][label] = indices;
195  this->Modified();
196  }
197 
200  GetLabeledSliceIndices( unsigned int axis, typename TImage::PixelType label )
201  {
202  return m_LabeledSlices[axis][label];
203  }
204 
205  // each label gets a set of slices in which it is present
206  typedef itksys::hash_map< typename TImage::PixelType, SliceSetType > LabeledSlicesType;
207  typedef std::vector< LabeledSlicesType > SliceIndicesType;
208 
212  {
213  return m_LabeledSlices;
214  }
215 
216 protected:
219  typename TImage::PixelType m_Label;
220  int m_Axis;
225  IdentifierType m_MinAlignIters; // minimum number of iterations in align method
226  IdentifierType m_MaxAlignIters; // maximum number of iterations in align method
227  IdentifierType m_ThreadCount; // for thread local instances
228  SliceIndicesType m_LabeledSlices; // one for each axis
229 
231  typedef Image< bool, TImage::ImageDimension > BoolImageType;
232  typedef Image< float, TImage::ImageDimension - 1 > FloatSliceType;
233  typedef Image< bool, TImage::ImageDimension - 1 > BoolSliceType;
234 
236  bool
237  ImagesEqual( typename BoolSliceType::Pointer& a, typename BoolSliceType::Pointer& b );
238 
240  virtual void
241  GenerateData() ITK_OVERRIDE;
242 
244  void
245  InterpolateBetweenTwo( int axis,
246  TImage* out,
247  typename TImage::PixelType label,
248  typename TImage::IndexValueType i,
249  typename TImage::IndexValueType j,
250  typename SliceType::Pointer& iconn,
251  typename SliceType::Pointer& jconn,
252  ThreadIdType threadId );
253 
259  void
260  InterpolateAlong( int axis, TImage* out );
261 
263  void
264  Extrapolate( int axis,
265  TImage* out,
266  typename TImage::PixelType label,
267  typename TImage::IndexValueType i,
268  typename TImage::IndexValueType j,
269  typename SliceType::Pointer& iConn,
270  typename TImage::PixelType iRegionId,
271  ThreadIdType threadId );
272 
274  typename FloatSliceType::Pointer
275  MaurerDM( typename BoolSliceType::Pointer& inImage, ThreadIdType threadId );
276 
279  std::vector< typename BoolSliceType::Pointer >
281  typename BoolSliceType::Pointer& end,
282  ThreadIdType threadId );
283 
285  typename BoolSliceType::Pointer
286  FindMedianImageDilations( typename BoolSliceType::Pointer& intersection,
287  typename BoolSliceType::Pointer& iMask,
288  typename BoolSliceType::Pointer& jMask,
289  ThreadIdType threadId );
290 
292  typename BoolSliceType::Pointer
293  FindMedianImageDistances( typename BoolSliceType::Pointer& intersection,
294  typename BoolSliceType::Pointer& iMask,
295  typename BoolSliceType::Pointer& jMask,
296  ThreadIdType threadId );
297 
299  void
300  Interpolate1to1( int axis,
301  TImage* out,
302  typename TImage::PixelType label,
303  typename TImage::IndexValueType i,
304  typename TImage::IndexValueType j,
305  typename SliceType::Pointer& iConn,
306  typename TImage::PixelType iRegionId,
307  typename SliceType::Pointer& jConn,
308  typename TImage::PixelType jRegionId,
309  const typename SliceType::IndexType& translation,
310  bool recursive,
311  ThreadIdType threadId );
312 
313  typedef std::vector< typename TImage::PixelType > PixelList;
314 
316  void
317  Interpolate1toN( int axis,
318  TImage* out,
319  typename TImage::PixelType label,
320  typename TImage::IndexValueType i,
321  typename TImage::IndexValueType j,
322  typename SliceType::Pointer& iConn,
323  typename TImage::PixelType iRegionId,
324  typename SliceType::Pointer& jConn,
325  const PixelList& jRegionIds,
326  const typename SliceType::IndexType& translation,
327  ThreadIdType threadId );
328 
330  typename SliceType::Pointer
331  TranslateImage( typename SliceType::Pointer& image,
332  const typename SliceType::IndexType& translation,
333  typename SliceType::RegionType newRegion );
334 
337  IdentifierType
338  CardSymDifference( typename BoolSliceType::Pointer& shape1, typename BoolSliceType::Pointer& shape2 );
339 
341  virtual void
342  AllocateOutputs() ITK_OVERRIDE;
343 
345  typename SliceType::IndexType
346  Centroid( typename SliceType::Pointer& conn, const PixelList& regionIds );
347 
350  void
351  IntersectionRegions( const typename SliceType::IndexType& translation,
352  typename SliceType::RegionType& iRegion,
353  typename SliceType::RegionType& jRegion );
354 
356  IdentifierType
357  Intersection( typename SliceType::Pointer& iConn,
358  typename TImage::PixelType iRegionId,
359  typename SliceType::Pointer& jConn,
360  const PixelList& jRegionIds,
361  const typename SliceType::IndexType& translation );
362 
364  typename SliceType::IndexType
365  Align( typename SliceType::Pointer& iConn,
366  typename TImage::PixelType iRegionId,
367  typename SliceType::Pointer& jConn,
368  const PixelList& jRegionIds );
369 
370  typedef itksys::hash_map< typename TImage::PixelType, typename TImage::RegionType > BoundingBoxesType;
371  BoundingBoxesType m_BoundingBoxes; // bounding box for each label
372 
374  typename SliceType::RegionType
375  BoundingBox( itk::SmartPointer< SliceType > image );
376 
380  template< typename T2 >
381  void
382  ExpandRegion( typename T2::RegionType& region, const typename T2::IndexType& index );
383 
385  typename SliceType::Pointer
386  RegionedConnectedComponents( const typename TImage::RegionType& region,
387  typename TImage::PixelType label,
388  IdentifierType& objectCount );
389 
391  typename BoolSliceType::Pointer
392  Dilate1( typename BoolSliceType::Pointer& seed, typename BoolSliceType::Pointer& mask, ThreadIdType threadId );
393 
394  typedef ExtractImageFilter< TImage, SliceType > RoiType;
395  typename RoiType::Pointer m_RoI;
396 
397  typedef BinaryThresholdImageFilter< SliceType, BoolSliceType > BinarizerType;
399 
400  typedef ConnectedComponentImageFilter< BoolSliceType, SliceType > ConnectedComponentsType;
402 
403 private:
404  MorphologicalContourInterpolator( const Self & ) ITK_DELETED_FUNCTION;
405  void
406  operator=( const Self & ) ITK_DELETED_FUNCTION;
407 };
408 } // namespace itk
409 
410 #ifndef ITK_MANUAL_INSTANTIATION
411 #include "itkMorphologicalContourInterpolator.hxx"
412 #endif
413 
414 #endif // itkMorphologicalContourInterpolator_h
IdentifierType CardSymDifference(typename BoolSliceType::Pointer &shape1, typename BoolSliceType::Pointer &shape2)
SliceType::Pointer RegionedConnectedComponents(const typename TImage::RegionType &region, typename TImage::PixelType label, IdentifierType &objectCount)
Image< float, TImage::ImageDimension - 1 > FloatSliceType
void IntersectionRegions(const typename SliceType::IndexType &translation, typename SliceType::RegionType &iRegion, typename SliceType::RegionType &jRegion)
std::vector< typename TImage::PixelType > PixelList
itksys::hash_map< typename TImage::PixelType, typename TImage::RegionType > BoundingBoxesType
BoolSliceType::Pointer FindMedianImageDistances(typename BoolSliceType::Pointer &intersection, typename BoolSliceType::Pointer &iMask, typename BoolSliceType::Pointer &jMask, ThreadIdType threadId)
SliceSetType GetLabeledSliceIndices(unsigned int axis, typename TImage::PixelType label)
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, ThreadIdType threadId)
itkSetMacro(Label, typename TImage::PixelType)
LRU Cache.
Simplified inverse ITK transforms.
itkTypeMacro(MorphologicalContourInterpolator, ImageToImageFilter)
BoolSliceType::Pointer Dilate1(typename BoolSliceType::Pointer &seed, typename BoolSliceType::Pointer &mask, ThreadIdType threadId)
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, ThreadIdType threadId)
Image< bool, TImage::ImageDimension - 1 > BoolSliceType
void SetLabeledSliceIndices(unsigned int axis, typename TImage::PixelType label, const std::vector< typename TImage::IndexValueType > &indices)
ExtractImageFilter< TImage, SliceType > RoiType
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, ThreadIdType threadId)
std::vector< typename BoolSliceType::Pointer > GenerateDilationSequence(typename BoolSliceType::Pointer &begin, typename BoolSliceType::Pointer &end, ThreadIdType threadId)
std::set< typename TImage::IndexValueType > SliceSetType
virtual void GenerateData() ITK_OVERRIDE
itkGetConstMacro(Label, typename TImage::PixelType)
bool ImagesEqual(typename BoolSliceType::Pointer &a, typename BoolSliceType::Pointer &b)
Image< typename TImage::PixelType, TImage::ImageDimension - 1 > SliceType
BoolSliceType::Pointer FindMedianImageDilations(typename BoolSliceType::Pointer &intersection, typename BoolSliceType::Pointer &iMask, typename BoolSliceType::Pointer &jMask, ThreadIdType threadId)
SliceType::IndexType Centroid(typename SliceType::Pointer &conn, const PixelList &regionIds)
BinaryThresholdImageFilter< SliceType, BoolSliceType > BinarizerType
SliceType::Pointer TranslateImage(typename SliceType::Pointer &image, const typename SliceType::IndexType &translation, typename SliceType::RegionType newRegion)
FloatSliceType::Pointer MaurerDM(typename BoolSliceType::Pointer &inImage, ThreadIdType threadId)
Interpolates contours between slices. Based on a paper by Albu et al.
ConnectedComponentImageFilter< BoolSliceType, SliceType > ConnectedComponentsType
virtual void AllocateOutputs() ITK_OVERRIDE
Image< bool, TImage::ImageDimension > BoolImageType
itksys::hash_map< typename TImage::PixelType, SliceSetType > LabeledSlicesType
void ExpandRegion(typename T2::RegionType &region, const typename T2::IndexType &index)
itkGetMacro(Label, typename TImage::PixelType)
void SetLabeledSliceIndices(unsigned int axis, typename TImage::PixelType label, const SliceSetType &indices)
SliceType::IndexType Align(typename SliceType::Pointer &iConn, typename TImage::PixelType iRegionId, typename SliceType::Pointer &jConn, const PixelList &jRegionIds)
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, ThreadIdType threadId)
SliceType::RegionType BoundingBox(itk::SmartPointer< SliceType > image)
IdentifierType Intersection(typename SliceType::Pointer &iConn, typename TImage::PixelType iRegionId, typename SliceType::Pointer &jConn, const PixelList &jRegionIds, const typename SliceType::IndexType &translation)
void InterpolateAlong(int axis, TImage *out)