ITK  6.0.0
Insight Toolkit
itkFastMarchingImageFilter.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 itkFastMarchingImageFilter_h
19 #define itkFastMarchingImageFilter_h
20 
21 #include "itkImageToImageFilter.h"
23 #include "itkLevelSet.h"
24 #include "itkMath.h"
25 #include "ITKFastMarchingExport.h"
26 
27 #include <functional>
28 #include <queue>
29 
30 namespace itk
31 {
37 {
38 public:
46  enum class Label : uint8_t
47  {
48  FarPoint = 0,
49  AlivePoint,
50  TrialPoint,
53  };
54 };
55 // Define how to print enumeration
56 extern ITKFastMarching_EXPORT std::ostream &
57  operator<<(std::ostream & out, const FastMarchingImageFilterEnums::Label value);
58 
134 template <typename TLevelSet, typename TSpeedImage = Image<float, TLevelSet::ImageDimension>>
135 class ITK_TEMPLATE_EXPORT FastMarchingImageFilter : public ImageToImageFilter<TSpeedImage, TLevelSet>
136 {
137 public:
138  ITK_DISALLOW_COPY_AND_MOVE(FastMarchingImageFilter);
139 
145 
147  itkNewMacro(Self);
148 
150  itkOverrideGetNameOfClassMacro(FastMarchingImageFilter);
151 
163  using OutputSpacingType = typename LevelSetImageType::SpacingType;
166 
167  class AxisNodeType : public NodeType
168  {
169  public:
170  AxisNodeType() = default;
171  int
172  GetAxis() const
173  {
174  return m_Axis;
175  }
176  void
177  SetAxis(int axis)
178  {
179  m_Axis = axis;
180  }
181  AxisNodeType &
182  operator=(const NodeType & node)
183  {
184  this->NodeType::operator=(node);
185  return *this;
186  }
187 
188  private:
189  int m_Axis{ 0 };
190  };
191 
193  using SpeedImageType = TSpeedImage;
194 
198 
200  static constexpr unsigned int SetDimension = LevelSetType::SetDimension;
201  static constexpr unsigned int SpeedImageDimension = SpeedImageType::ImageDimension;
202 
205 
207 #if !defined(ITK_LEGACY_REMOVE)
208 
209  static constexpr LabelEnum FarPoint = LabelEnum::FarPoint;
210  static constexpr LabelEnum AlivePoint = LabelEnum::AlivePoint;
211  static constexpr LabelEnum TrialPoint = LabelEnum::TrialPoint;
212  static constexpr LabelEnum InitialTrialPoint = LabelEnum::InitialTrialPoint;
213  static constexpr LabelEnum OutsidePoint = LabelEnum::OutsidePoint;
214 #endif
215 
218 
221 
222  template <typename TPixel>
223  void
225  {
226  using InternalImageType = Image<TPixel, SetDimension>;
227  using InternalRegionIterator = ImageRegionConstIteratorWithIndex<InternalImageType>;
228  InternalRegionIterator b_it(iImage, iImage->GetLargestPossibleRegion());
229  b_it.GoToBegin();
230 
231  TPixel zero_value{};
232  typename NodeContainer::ElementIdentifier NumberOfPoints = 0;
233 
234  NodeType node;
235  node.SetValue(0.);
236 
237  while (!b_it.IsAtEnd())
238  {
239  if (Math::ExactlyEquals(b_it.Get(), zero_value))
240  {
241  if (NumberOfPoints == 0)
242  {
243  m_OutsidePoints = NodeContainer::New();
244  }
245  node.SetIndex(b_it.GetIndex());
246  m_OutsidePoints->InsertElement(NumberOfPoints++, node);
247  }
248  ++b_it;
249  }
250  this->Modified();
251  }
252 
254  void
256  {
257  m_OutsidePoints = points;
258  this->Modified();
259  }
264  void
266  {
267  m_AlivePoints = points;
268  this->Modified();
269  }
273  NodeContainerPointer
275  {
276  return m_AlivePoints;
277  }
278 
281  void
283  {
284  m_TrialPoints = points;
285  this->Modified();
286  }
290  NodeContainerPointer
292  {
293  return m_TrialPoints;
294  }
295 
297  LabelImagePointer
299  {
300  return m_LabelImage;
301  }
302 
306  void
307  SetSpeedConstant(double value)
308  {
309  m_SpeedConstant = value;
310  m_InverseSpeed = -1.0 * itk::Math::sqr(1.0 / m_SpeedConstant);
311  this->Modified();
312  }
316  itkGetConstReferenceMacro(SpeedConstant, double);
317 
322  itkSetMacro(NormalizationFactor, double);
323  itkGetConstMacro(NormalizationFactor, double);
329  itkSetMacro(StoppingValue, double);
330 
332  itkGetConstReferenceMacro(StoppingValue, double);
333 
338  itkSetMacro(CollectPoints, bool);
339 
341  itkGetConstReferenceMacro(CollectPoints, bool);
342  itkBooleanMacro(CollectPoints);
349  NodeContainerPointer
351  {
352  return m_ProcessedPoints;
353  }
354 
361  virtual void
363  {
364  m_OutputRegion = size;
365  }
366  virtual OutputSizeType
368  {
369  return m_OutputRegion.GetSize();
370  }
371  itkSetMacro(OutputRegion, OutputRegionType);
372  itkGetConstReferenceMacro(OutputRegion, OutputRegionType);
373  itkSetMacro(OutputSpacing, OutputSpacingType);
374  itkGetConstReferenceMacro(OutputSpacing, OutputSpacingType);
375  itkSetMacro(OutputDirection, OutputDirectionType);
376  itkGetConstReferenceMacro(OutputDirection, OutputDirectionType);
377  itkSetMacro(OutputOrigin, OutputPointType);
378  itkGetConstReferenceMacro(OutputOrigin, OutputPointType);
379  itkSetMacro(OverrideOutputInformation, bool);
380  itkGetConstReferenceMacro(OverrideOutputInformation, bool);
381  itkBooleanMacro(OverrideOutputInformation);
384 #ifdef ITK_USE_CONCEPT_CHECKING
385  // Begin concept checking
388  itkConceptMacro(DoubleConvertibleToLevelSetCheck, (Concept::Convertible<double, PixelType>));
389  itkConceptMacro(LevelSetOStreamWritableCheck, (Concept::OStreamWritable<PixelType>));
390  // End concept checking
391 #endif
392 
393 protected:
395  ~FastMarchingImageFilter() override = default;
396  void
397  PrintSelf(std::ostream & os, Indent indent) const override;
398 
399  virtual void
400  Initialize(LevelSetImageType *);
401 
402  virtual void
403  UpdateNeighbors(const IndexType & index, const SpeedImageType *, LevelSetImageType *);
404 
405  virtual double
406  UpdateValue(const IndexType & index, const SpeedImageType *, LevelSetImageType *);
407 
408  const AxisNodeType &
409  GetNodeUsedInCalculation(unsigned int idx) const
410  {
411  return m_NodesUsed[idx];
412  }
413 
414  void
415  GenerateData() override;
416 
418  void
419  GenerateOutputInformation() override;
420 
421  void
422  EnlargeOutputRequestedRegion(DataObject * output) override;
423 
428  itkGetConstReferenceMacro(LargeValue, PixelType);
429 
430  OutputRegionType m_BufferedRegion{};
432  LevelSetIndexType m_StartIndex{};
433  LevelSetIndexType m_LastIndex{};
434 
435  itkGetConstReferenceMacro(StartIndex, LevelSetIndexType);
436  itkGetConstReferenceMacro(LastIndex, LevelSetIndexType);
437 
438 private:
439  NodeContainerPointer m_AlivePoints{};
440  NodeContainerPointer m_TrialPoints{};
441  NodeContainerPointer m_OutsidePoints{};
442 
443  LabelImagePointer m_LabelImage{};
444 
445  double m_SpeedConstant{};
446  double m_InverseSpeed{};
447  double m_StoppingValue{};
448 
449  bool m_CollectPoints{};
450  NodeContainerPointer m_ProcessedPoints{};
451 
452  OutputRegionType m_OutputRegion{};
453  OutputPointType m_OutputOrigin{};
454  OutputSpacingType m_OutputSpacing{};
455  OutputDirectionType m_OutputDirection{};
456  bool m_OverrideOutputInformation{};
457 
458  typename LevelSetImageType::PixelType m_LargeValue{};
459  AxisNodeType m_NodesUsed[SetDimension]{};
460 
464  using HeapContainer = std::vector<AxisNodeType>;
465  using NodeComparer = std::greater<AxisNodeType>;
466  using HeapType = std::priority_queue<AxisNodeType, HeapContainer, NodeComparer>;
467 
468  HeapType m_TrialHeap{};
469 
470  double m_NormalizationFactor{};
471 };
472 } // namespace itk
473 
474 #ifndef ITK_MANUAL_INSTANTIATION
475 # include "itkFastMarchingImageFilter.hxx"
476 #endif
477 
478 #endif
Pointer
SmartPointer< Self > Pointer
Definition: itkAddImageFilter.h:93
itk::FastMarchingImageFilter::PixelType
typename LevelSetType::PixelType PixelType
Definition: itkFastMarchingImageFilter.h:156
itk::FastMarchingImageFilter::HeapType
std::priority_queue< AxisNodeType, HeapContainer, NodeComparer > HeapType
Definition: itkFastMarchingImageFilter.h:466
ConstPointer
SmartPointer< const Self > ConstPointer
Definition: itkAddImageFilter.h:94
itk::LevelSetTypeDefault
Level set type information.
Definition: itkLevelSet.h:40
itk::Index
Represent a n-dimensional index in a n-dimensional image.
Definition: itkIndex.h:68
itk::Concept::OStreamWritable
Definition: itkConceptChecking.h:638
itk::GTest::TypedefsAndConstructors::Dimension2::DirectionType
ImageBaseType::DirectionType DirectionType
Definition: itkGTestTypedefsAndConstructors.h:52
itk::FastMarchingImageFilter::LabelImagePointer
typename LabelImageType::Pointer LabelImagePointer
Definition: itkFastMarchingImageFilter.h:220
itk::detail::VectorContainer
Define a front-end to the STL "vector" container that conforms to the IndexedContainerInterface.
Definition: itkVectorContainer.h:51
itkImageRegionConstIteratorWithIndex.h
itk::ImageRegionConstIteratorWithIndex
A multi-dimensional iterator templated over image type that walks an image region and is specialized ...
Definition: itkImageRegionConstIteratorWithIndex.h:130
itk::FastMarchingImageFilter::SetTrialPoints
void SetTrialPoints(NodeContainer *points)
Definition: itkFastMarchingImageFilter.h:282
itk::operator<<
std::ostream & operator<<(std::ostream &os, const Array< TValue > &arr)
Definition: itkArray.h:216
itk::FastMarchingImageFilter::SetSpeedConstant
void SetSpeedConstant(double value)
Definition: itkFastMarchingImageFilter.h:307
itk::FastMarchingImageFilter::LevelSetPointer
typename LevelSetType::LevelSetPointer LevelSetPointer
Definition: itkFastMarchingImageFilter.h:155
itk::FastMarchingImageFilter::AxisNodeType::operator=
AxisNodeType & operator=(const NodeType &node)
Definition: itkFastMarchingImageFilter.h:182
itk::FastMarchingImageFilter::SetOutsidePoints
void SetOutsidePoints(NodeContainer *points)
Definition: itkFastMarchingImageFilter.h:255
itk::FastMarchingImageFilterEnums::Label::FarPoint
itk::FastMarchingImageFilter::SpeedImagePointer
typename SpeedImageType::Pointer SpeedImagePointer
Definition: itkFastMarchingImageFilter.h:196
itk::GTest::TypedefsAndConstructors::Dimension2::PointType
ImageBaseType::PointType PointType
Definition: itkGTestTypedefsAndConstructors.h:51
itk::Math::ExactlyEquals
bool ExactlyEquals(const TInput1 &x1, const TInput2 &x2)
Return the result of an exact comparison between two scalar values of potentially different types.
Definition: itkMath.h:726
itk::FastMarchingImageFilter::OutputSpacingType
typename LevelSetImageType::SpacingType OutputSpacingType
Definition: itkFastMarchingImageFilter.h:163
itk::FastMarchingImageFilter::GetTrialPoints
NodeContainerPointer GetTrialPoints()
Definition: itkFastMarchingImageFilter.h:291
itk::GTest::TypedefsAndConstructors::Dimension2::SizeType
ImageBaseType::SizeType SizeType
Definition: itkGTestTypedefsAndConstructors.h:49
itk::SmartPointer< Self >
itk::Indent
Control indentation during Print() invocation.
Definition: itkIndent.h:49
itk::FastMarchingImageFilter::NodeComparer
std::greater< AxisNodeType > NodeComparer
Definition: itkFastMarchingImageFilter.h:465
itk::LevelSetNode
Represent a node in a level set.
Definition: itkLevelSetNode.h:45
itk::FastMarchingImageFilter::GetAlivePoints
NodeContainerPointer GetAlivePoints()
Definition: itkFastMarchingImageFilter.h:274
itk::FastMarchingImageFilter::LevelSetIndexType
typename LevelSetImageType::IndexType LevelSetIndexType
Definition: itkFastMarchingImageFilter.h:431
itk::FastMarchingImageFilter::NodeIndexType
typename NodeType::IndexType NodeIndexType
Definition: itkFastMarchingImageFilter.h:158
itk::FastMarchingImageFilter::SetBinaryMask
void SetBinaryMask(Image< TPixel, SetDimension > *iImage)
Definition: itkFastMarchingImageFilter.h:224
itk::Concept::SameDimension
Definition: itkConceptChecking.h:696
itk::FastMarchingImageFilter::NodeType
typename LevelSetType::NodeType NodeType
Definition: itkFastMarchingImageFilter.h:157
itk::GTest::TypedefsAndConstructors::Dimension2::IndexType
ImageBaseType::IndexType IndexType
Definition: itkGTestTypedefsAndConstructors.h:50
itk::ImageToImageFilter
Base class for filters that take an image as input and produce an image as output.
Definition: itkImageToImageFilter.h:108
itk::FastMarchingImageFilterEnums::Label::InitialTrialPoint
itk::FastMarchingImageFilterEnums
Contains all enum classes used by the FastMarchingImageFilter class.
Definition: itkFastMarchingImageFilter.h:36
itk::LevelSetTypeDefault::LevelSetPointer
typename TLevelSet::Pointer LevelSetPointer
Definition: itkLevelSet.h:51
itk::LightObject
Light weight base class for most itk classes.
Definition: itkLightObject.h:55
itk::FastMarchingImageFilterEnums::Label
Label
Definition: itkFastMarchingImageFilter.h:46
itk::FastMarchingImageFilter::SetOutputSize
virtual void SetOutputSize(const OutputSizeType &size)
Definition: itkFastMarchingImageFilter.h:362
itk::FastMarchingImageFilter
Solve an Eikonal equation using Fast Marching.
Definition: itkFastMarchingImageFilter.h:135
itk::FastMarchingImageFilter::SpeedImageType
TSpeedImage SpeedImageType
Definition: itkFastMarchingImageFilter.h:193
itk::LevelSetTypeDefault::NodeContainerPointer
typename NodeContainer::Pointer NodeContainerPointer
Definition: itkLevelSet.h:64
itk::GTest::TypedefsAndConstructors::Dimension2::RegionType
ImageBaseType::RegionType RegionType
Definition: itkGTestTypedefsAndConstructors.h:54
itk::FastMarchingImageFilter::OutputDirectionType
typename LevelSetImageType::DirectionType OutputDirectionType
Definition: itkFastMarchingImageFilter.h:164
itk::FastMarchingImageFilter::AxisNodeType
Definition: itkFastMarchingImageFilter.h:167
itk::FastMarchingImageFilter::AxisNodeType::SetAxis
void SetAxis(int axis)
Definition: itkFastMarchingImageFilter.h:177
itkImageToImageFilter.h
itk::FastMarchingImageFilterEnums::Label::AlivePoint
itk::FastMarchingImageFilter::OutputRegionType
typename LevelSetImageType::RegionType OutputRegionType
Definition: itkFastMarchingImageFilter.h:162
itk::FastMarchingImageFilter::GetOutputSize
virtual OutputSizeType GetOutputSize() const
Definition: itkFastMarchingImageFilter.h:367
itk::ImageBase::GetLargestPossibleRegion
virtual const RegionType & GetLargestPossibleRegion() const
Definition: itkImageBase.h:279
itk::FastMarchingImageFilterEnums::Label::TrialPoint
itk::ImageConstIteratorWithIndex::GoToBegin
void GoToBegin()
itk::LevelSetTypeDefault::LevelSetImageType
TLevelSet LevelSetImageType
Definition: itkLevelSet.h:45
itk::FastMarchingImageFilter::HeapContainer
std::vector< AxisNodeType > HeapContainer
Definition: itkFastMarchingImageFilter.h:464
itkConceptMacro
#define itkConceptMacro(name, concept)
Definition: itkConceptChecking.h:65
itk
The "itk" namespace contains all Insight Segmentation and Registration Toolkit (ITK) classes....
Definition: itkAnnulusOperator.h:24
itk::FastMarchingImageFilter::OutputPointType
typename LevelSetImageType::PointType OutputPointType
Definition: itkFastMarchingImageFilter.h:165
Label
itkLevelSet.h
itk::FastMarchingImageFilter::NodeContainer
typename LevelSetType::NodeContainer NodeContainer
Definition: itkFastMarchingImageFilter.h:159
itk::FastMarchingImageFilter::OutputSizeType
typename LevelSetImageType::SizeType OutputSizeType
Definition: itkFastMarchingImageFilter.h:161
itk::FastMarchingImageFilter::NodeContainerPointer
typename LevelSetType::NodeContainerPointer NodeContainerPointer
Definition: itkFastMarchingImageFilter.h:160
itk::Concept::Convertible
Definition: itkConceptChecking.h:216
itk::FastMarchingImageFilter::GetLabelImage
LabelImagePointer GetLabelImage() const
Definition: itkFastMarchingImageFilter.h:298
itk::FastMarchingImageFilter::LevelSetImageType
typename LevelSetType::LevelSetImageType LevelSetImageType
Definition: itkFastMarchingImageFilter.h:154
itk::FastMarchingImageFilterEnums::Label::OutsidePoint
itk::Image
Templated n-dimensional image class.
Definition: itkImage.h:88
New
static Pointer New()
itk::FastMarchingImageFilter::AxisNodeType::GetAxis
int GetAxis() const
Definition: itkFastMarchingImageFilter.h:172
itk::FastMarchingImageFilter::SpeedImageConstPointer
typename SpeedImageType::ConstPointer SpeedImageConstPointer
Definition: itkFastMarchingImageFilter.h:197
itk::LevelSetTypeDefault::PixelType
typename TLevelSet::PixelType PixelType
Definition: itkLevelSet.h:55
itk::FastMarchingImageFilter::GetNodeUsedInCalculation
const AxisNodeType & GetNodeUsedInCalculation(unsigned int idx) const
Definition: itkFastMarchingImageFilter.h:409
itk::FastMarchingImageFilter::SetAlivePoints
void SetAlivePoints(NodeContainer *points)
Definition: itkFastMarchingImageFilter.h:265
itkMath.h
itk::DataObject
Base class for all data objects in ITK.
Definition: itkDataObject.h:293
itk::FastMarchingImageFilter::GetProcessedPoints
NodeContainerPointer GetProcessedPoints() const
Definition: itkFastMarchingImageFilter.h:350