ITK  5.0.0
Insight Segmentation and Registration Toolkit
itkFastMarchingImageFilter.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 itkFastMarchingImageFilter_h
19 #define itkFastMarchingImageFilter_h
20 
21 #include "itkImageToImageFilter.h"
23 #include "itkLevelSet.h"
24 #include "itkMath.h"
25 
26 #include <functional>
27 #include <queue>
28 #include "itkMath.h"
29 
30 namespace itk
31 {
106 template<
107  typename TLevelSet,
108  typename TSpeedImage = Image< float, TLevelSet ::ImageDimension > >
109 class ITK_TEMPLATE_EXPORT FastMarchingImageFilter:
110  public ImageToImageFilter< TSpeedImage, TLevelSet >
111 {
112 public:
113  ITK_DISALLOW_COPY_AND_ASSIGN(FastMarchingImageFilter);
114 
120 
122  itkNewMacro(Self);
123 
126 
138  using OutputSpacingType = typename LevelSetImageType::SpacingType;
141 
142  class AxisNodeType:public NodeType
143  {
144 public:
146  int GetAxis() const { return m_Axis; }
147  void SetAxis(int axis) { m_Axis = axis; }
148  const AxisNodeType & operator=(const NodeType & node)
149  { this->NodeType::operator=(node); return *this; }
150 
151 private:
152  int m_Axis{0};
153  };
154 
156  using SpeedImageType = TSpeedImage;
157 
159  using SpeedImagePointer = typename SpeedImageType::Pointer;
160  using SpeedImageConstPointer = typename SpeedImageType::ConstPointer;
161 
163  static constexpr unsigned int SetDimension = LevelSetType::SetDimension;
164  static constexpr unsigned int SpeedImageDimension = SpeedImageType::ImageDimension;
165 
168 
173  enum LabelType { FarPoint = 0, AlivePoint,
174  TrialPoint, InitialTrialPoint, OutsidePoint };
175 
178 
181 
182  template< typename TPixel >
184  {
185  using InternalImageType = Image< TPixel, SetDimension >;
186  using InternalRegionIterator = ImageRegionConstIteratorWithIndex<InternalImageType>;
187  InternalRegionIterator b_it( iImage, iImage->GetLargestPossibleRegion() );
188  b_it.GoToBegin();
189 
190  TPixel zero_value = NumericTraits< TPixel >::ZeroValue();
191  typename NodeContainer::ElementIdentifier NumberOfPoints = 0;
192 
193  NodeType node;
194  node.SetValue( 0. );
195 
196  while( !b_it.IsAtEnd() )
197  {
198  if( Math::ExactlyEquals(b_it.Get(), zero_value) )
199  {
200  if( NumberOfPoints == 0 )
201  {
202  m_OutsidePoints = NodeContainer::New();
203  }
204  node.SetIndex( b_it.GetIndex() );
205  m_OutsidePoints->InsertElement( NumberOfPoints++, node );
206 
207  }
208  ++b_it;
209  }
210  this->Modified();
211  }
212 
215  {
216  m_OutsidePoints = points;
217  this->Modified();
218  }
220 
224  {
225  m_AlivePoints = points;
226  this->Modified();
227  }
229 
232  {
233  return m_AlivePoints;
234  }
235 
239  {
240  m_TrialPoints = points;
241  this->Modified();
242  }
244 
247  {
248  return m_TrialPoints;
249  }
250 
253  {
254  return m_LabelImage;
255  }
256 
260  void SetSpeedConstant(double value)
261  {
262  m_SpeedConstant = value;
263  m_InverseSpeed = -1.0 * itk::Math::sqr(1.0 / m_SpeedConstant);
264  this->Modified();
265  }
267 
269  itkGetConstReferenceMacro(SpeedConstant, double);
270 
275  itkSetMacro(NormalizationFactor, double);
276  itkGetConstMacro(NormalizationFactor, double);
278 
282  itkSetMacro(StoppingValue, double);
283 
285  itkGetConstReferenceMacro(StoppingValue, double);
286 
291  itkSetMacro(CollectPoints, bool);
292 
294  itkGetConstReferenceMacro(CollectPoints, bool);
295  itkBooleanMacro(CollectPoints);
297 
303  {
304  return m_ProcessedPoints;
305  }
306 
313  virtual void SetOutputSize(const OutputSizeType & size)
314  { m_OutputRegion = size; }
316  { return m_OutputRegion.GetSize(); }
317  itkSetMacro(OutputRegion, OutputRegionType);
318  itkGetConstReferenceMacro(OutputRegion, OutputRegionType);
319  itkSetMacro(OutputSpacing, OutputSpacingType);
320  itkGetConstReferenceMacro(OutputSpacing, OutputSpacingType);
321  itkSetMacro(OutputDirection, OutputDirectionType);
322  itkGetConstReferenceMacro(OutputDirection, OutputDirectionType);
323  itkSetMacro(OutputOrigin, OutputPointType);
324  itkGetConstReferenceMacro(OutputOrigin, OutputPointType);
325  itkSetMacro(OverrideOutputInformation, bool);
326  itkGetConstReferenceMacro(OverrideOutputInformation, bool);
327  itkBooleanMacro(OverrideOutputInformation);
329 
330 #ifdef ITK_USE_CONCEPT_CHECKING
331  // Begin concept checking
332  itkConceptMacro( SameDimensionCheck,
334  itkConceptMacro( SpeedConvertibleToDoubleCheck,
336  itkConceptMacro( DoubleConvertibleToLevelSetCheck,
338  itkConceptMacro( LevelSetOStreamWritableCheck,
340  // End concept checking
341 #endif
342 
343 protected:
345  ~FastMarchingImageFilter() override = default;
346  void PrintSelf(std::ostream & os, Indent indent) const override;
347 
348  virtual void Initialize(LevelSetImageType *);
349 
350  virtual void UpdateNeighbors(const IndexType & index,
351  const SpeedImageType *, LevelSetImageType *);
352 
353  virtual double UpdateValue(const IndexType & index,
354  const SpeedImageType *, LevelSetImageType *);
355 
356  const AxisNodeType & GetNodeUsedInCalculation(unsigned int idx) const
357  { return m_NodesUsed[idx]; }
358 
359  void GenerateData() override;
360 
362  void GenerateOutputInformation() override;
363 
364  void EnlargeOutputRequestedRegion(DataObject *output) override;
365 
370  itkGetConstReferenceMacro(LargeValue, PixelType);
371 
376 
377  itkGetConstReferenceMacro(StartIndex, LevelSetIndexType);
378  itkGetConstReferenceMacro(LastIndex, LevelSetIndexType);
379 
380 private:
384 
386 
390 
393 
399 
400  typename LevelSetImageType::PixelType m_LargeValue;
401  AxisNodeType m_NodesUsed[SetDimension];
402 
406  using HeapContainer = std::vector< AxisNodeType >;
407  using NodeComparer = std::greater< AxisNodeType >;
408  using HeapType =
409  std::priority_queue< AxisNodeType, HeapContainer, NodeComparer >;
410 
412 
414 };
415 } // namespace itk
416 
417 #ifndef ITK_MANUAL_INSTANTIATION
418 #include "itkFastMarchingImageFilter.hxx"
419 #endif
420 
421 #endif
typename LevelSetImageType::DirectionType OutputDirectionType
void SetBinaryMask(Image< TPixel, SetDimension > *iImage)
typename LevelSetImageType::SizeType OutputSizeType
typename LevelSetType::NodeType NodeType
typename NodeContainer::Pointer NodeContainerPointer
Definition: itkLevelSet.h:65
Light weight base class for most itk classes.
typename SpeedImageType::Pointer SpeedImagePointer
typename LevelSetType::LevelSetPointer LevelSetPointer
LevelSetImageType::PixelType m_LargeValue
typename LevelSetType::LevelSetImageType LevelSetImageType
std::priority_queue< AxisNodeType, HeapContainer, NodeComparer > HeapType
typename LevelSetImageType::RegionType OutputRegionType
Define numeric traits for std::vector.
Represent a n-dimensional index in a n-dimensional image.
Definition: itkIndex.h:66
virtual const RegionType & GetLargestPossibleRegion() const
Definition: itkImageBase.h:252
virtual void SetOutputSize(const OutputSizeType &size)
void SetOutsidePoints(NodeContainer *points)
bool ExactlyEquals(const TInput1 &x1, const TInput2 &x2)
Return the result of an exact comparison between two scalar values of potetially different types...
Definition: itkMath.h:707
typename TLevelSet::Pointer LevelSetPointer
Definition: itkLevelSet.h:51
void SetTrialPoints(NodeContainer *points)
typename LevelSetType::NodeContainerPointer NodeContainerPointer
std::greater< AxisNodeType > NodeComparer
A multi-dimensional iterator templated over image type that walks an image region and is specialized ...
typename TLevelSet::PixelType PixelType
Definition: itkLevelSet.h:55
virtual OutputSizeType GetOutputSize() const
std::vector< AxisNodeType > HeapContainer
Represent a node in a level set.
typename LabelImageType::Pointer LabelImagePointer
const AxisNodeType & GetNodeUsedInCalculation(unsigned int idx) const
const AxisNodeType & operator=(const NodeType &node)
typename LevelSetType::NodeContainer NodeContainer
LabelImagePointer GetLabelImage() const
Base class for filters that take an image as input and produce an image as output.
Define a front-end to the STL &quot;vector&quot; container that conforms to the IndexedContainerInterface.
typename LevelSetType::PixelType PixelType
typename SpeedImageType::ConstPointer SpeedImageConstPointer
Control indentation during Print() invocation.
Definition: itkIndent.h:49
typename LevelSetImageType::SpacingType OutputSpacingType
typename NodeType::IndexType NodeIndexType
Solve an Eikonal equation using Fast Marching.
typename LevelSetImageType::IndexType LevelSetIndexType
#define itkConceptMacro(name, concept)
Level set type information.
Definition: itkLevelSet.h:40
Base class for all data objects in ITK.
Templated n-dimensional image class.
Definition: itkImage.h:75
void SetAlivePoints(NodeContainer *points)
typename LevelSetImageType::PointType OutputPointType
NodeContainerPointer GetProcessedPoints() const