itkCuberilleImageToMeshFilter.h

Go to the documentation of this file.

/*=========================================================================
 *
 *  Copyright Insight Software Consortium
 *
 *  Licensed under the Apache License, Version 2.0 (the "License");
 *  you may not use this file except in compliance with the License.
 *  You may obtain a copy of the License at
 *
 *         http://www.apache.org/licenses/LICENSE-2.0.txt
 *
 *  Unless required by applicable law or agreed to in writing, software
 *  distributed under the License is distributed on an "AS IS" BASIS,
 *  WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
 *  See the License for the specific language governing permissions and
 *  limitations under the License.
 *
 *=========================================================================*/
#ifndef itkCuberilleImageToMeshFilter_h
#define itkCuberilleImageToMeshFilter_h

#define DEBUG_PRINT 0
#define USE_GRADIENT_RECURSIVE_GAUSSIAN 0
#define USE_ADVANCED_PROJECTION 0
#define USE_LINESEARCH_PROJECTION 0

#include "itkMacro.h"
#include "itkMesh.h"
#include "itkImageToMeshFilter.h"
#include "itkCellInterface.h"
#include "itkTriangleCell.h"
#include "itkQuadrilateralCell.h"
#include "itkDefaultStaticMeshTraits.h"
#include "itkConstShapedNeighborhoodIterator.h"
#include "itkLinearInterpolateImageFunction.h"
#include "itkGradientImageFilter.h"
#include "itkGradientRecursiveGaussianImageFilter.h"
#include "itkVectorLinearInterpolateImageFunction.h"

namespace itk
{

template < class TInputImage, class TOutputMesh, class TInterpolator=itk::LinearInterpolateImageFunction<TInputImage> >
class CuberilleImageToMeshFilter : public ImageToMeshFilter< TInputImage, TOutputMesh >
{
public:
  ITK_DISALLOW_COPY_AND_ASSIGN(CuberilleImageToMeshFilter);

  using Self = CuberilleImageToMeshFilter;
  using Superclass = ImageToMeshFilter< TInputImage, TOutputMesh >;
  using Pointer = SmartPointer<Self>;
  using ConstPointer = SmartPointer<const Self>;

  itkNewMacro(Self);

  itkTypeMacro(CuberilleImageToMeshFilter, ImageToMeshFilter);

  using OutputMeshType = TOutputMesh;
  using OutputMeshPointer = typename OutputMeshType::Pointer;
  using OutputMeshTraits = typename OutputMeshType::MeshTraits;
  using OutputPointType = typename OutputMeshType::PointType;
  using OutputPixelType = typename OutputMeshTraits::PixelType;
  using CellTraits = typename OutputMeshType::CellTraits;
  using PointsContainerPointer = typename OutputMeshType::PointsContainerPointer;
  using PointsContainer = typename OutputMeshType::PointsContainer;
  using CellsContainerPointer = typename OutputMeshType::CellsContainerPointer;
  using CellsContainer = typename OutputMeshType::CellsContainer;
  using PointIdentifier = typename OutputMeshType::PointIdentifier;
  using CellIdentifier = typename OutputMeshType::CellIdentifier;
  using CellInterfaceType = CellInterface<OutputPixelType, CellTraits>;
  using TriangleCellType = TriangleCell<CellInterfaceType>;
  using TriangleAutoPointer = typename TriangleCellType::SelfAutoPointer;
  using TriangleCellAutoPointer = typename TriangleCellType::CellAutoPointer;
  using QuadrilateralCellType = QuadrilateralCell<CellInterfaceType>;
  using QuadrilateralAutoPointer = typename QuadrilateralCellType::SelfAutoPointer;
  using QuadrilateralCellAutoPointer = typename QuadrilateralCellType::CellAutoPointer;

  using InputImageType = TInputImage;
  using InputImagePointer = typename InputImageType::Pointer;
  using InputImageConstPointer = typename InputImageType::ConstPointer;
  using InputPixelType = typename InputImageType::PixelType;
  using SizeType = typename InputImageType::SizeType;
  using SpacingType = typename InputImageType::SpacingType;
  using SpacingValueType = typename InputImageType::SpacingValueType;
  using IndexType = typename InputImageType::IndexType;
  using PointType = typename OutputMeshType::PointType;

  using InterpolatorType = TInterpolator;
  using InterpolatorPointer = typename InterpolatorType::Pointer;
  using InterpolatorOutputType = typename InterpolatorType::OutputType;

  using InputImageIteratorType = ConstShapedNeighborhoodIterator< InputImageType >;
#if USE_GRADIENT_RECURSIVE_GAUSSIAN
  using GradientFilterType = GradientRecursiveGaussianImageFilter< InputImageType >;
#else
  using GradientFilterType = GradientImageFilter< InputImageType >;
#endif
  using GradientFilterPointer = typename GradientFilterType::Pointer;
  using GradientImageType = typename GradientFilterType::OutputImageType;
  using GradientImagePointer = typename GradientImageType::Pointer;
  using GradientPixelType = typename GradientFilterType::OutputPixelType;
  using GradientInterpolatorType = itk::VectorLinearInterpolateImageFunction< GradientImageType >;
  using GradientInterpolatorPointer = typename GradientInterpolatorType::Pointer;

  itkGetMacro( IsoSurfaceValue, InputPixelType );
  itkSetMacro( IsoSurfaceValue, InputPixelType );

  using Superclass::SetInput;
  virtual void SetInput( const InputImageType * inputImage );

  itkGetConstObjectMacro( Interpolator, InterpolatorType );
  itkSetObjectMacro( Interpolator, InterpolatorType );

  itkGetMacro( GenerateTriangleFaces, bool );
  itkSetMacro( GenerateTriangleFaces, bool );
  itkBooleanMacro( GenerateTriangleFaces );

  itkGetMacro( ProjectVerticesToIsoSurface, bool );
  itkSetMacro( ProjectVerticesToIsoSurface, bool );
  itkBooleanMacro( ProjectVerticesToIsoSurface );

  itkGetMacro( ProjectVertexSurfaceDistanceThreshold, double );
  itkSetClampMacro( ProjectVertexSurfaceDistanceThreshold, double, 0.0, NumericTraits<InputPixelType>::max() );

  itkGetMacro( ProjectVertexStepLength, double );
  itkSetClampMacro( ProjectVertexStepLength, double, 0.0, 100000.0 );

  itkGetMacro( ProjectVertexStepLengthRelaxationFactor, double );
  itkSetClampMacro( ProjectVertexStepLengthRelaxationFactor, double, 0.0, 1.0 );

  itkGetMacro( ProjectVertexMaximumNumberOfSteps, unsigned int );
  itkSetMacro( ProjectVertexMaximumNumberOfSteps, unsigned int );

protected:
  CuberilleImageToMeshFilter();
  ~CuberilleImageToMeshFilter() override;
  void PrintSelf(std::ostream& os, Indent indent) const override;

  void GenerateData() override;
  void GenerateOutputInformation() override{ }; // do nothing

private:
  class VertexLookupNode
    {
    public:
      using Self = VertexLookupNode;

      VertexLookupNode() : m_X( 0 ), m_Y( 0 ) { }
      VertexLookupNode(unsigned long x, unsigned long y) : m_X( x ), m_Y( y ) { }

      unsigned long GetX() { return m_X; }
      unsigned long GetY() { return m_Y; }

      bool operator> ( const Self& node ) const
        { return (m_Y > node.m_Y) || ( (m_Y == node.m_Y) && (m_X > node.m_X) ); }
      bool operator>= ( const Self& node ) const
        { return (m_Y >= node.m_Y) || ( (m_Y == node.m_Y) && (m_X >= node.m_X) ); }
      bool operator< ( const Self& node ) const
        { return (m_Y < node.m_Y) || ( (m_Y == node.m_Y) && (m_X < node.m_X) ); }
      bool operator<= ( const Self& node ) const
        { return (m_Y <= node.m_Y) || ( (m_Y == node.m_Y) && (m_X <= node.m_X) ); }

    private:
      unsigned long m_X;
      unsigned long m_Y;
  };

  template< class TMeshType >
  class VertexLookupMap
    {
    public:
      using Self = VertexLookupMap;
      using PointIdentifier = typename TMeshType::PointIdentifier;
      using MapType = std::map< VertexLookupNode, PointIdentifier >;

      VertexLookupMap() { }

      void Clear() { m_Map.clear(); }

      void AddVertex( unsigned int x, unsigned int y, PointIdentifier id )
        {
          VertexLookupNode node( x, y );
          m_Map.insert( typename MapType::value_type(node, id) );
        }

      bool GetVertex( unsigned int x, unsigned int y, PointIdentifier &id )
        {
        bool result = false;
        VertexLookupNode node( x, y );
        auto it = m_Map.find( node );
        if ( it != m_Map.end() )
          {
          result = true;
          id = it->second;
          }
        return result;
        }

    private:
      MapType m_Map;
  };

  using VertexLookupMapType = VertexLookupMap< OutputMeshType >;

  inline void ComputeGradientImage();

  inline void SetVerticesFromFace( unsigned int face, bool *vertexHasQuad );

  inline IndexType GetVertexLookupIndex( unsigned int vertex, IndexType index );

  inline void ProjectVertexToIsoSurface( PointType &vertex );

  inline void AddVertex( PointIdentifier &id, IndexType index, const InputImageType* image, OutputMeshType* mesh );

  inline void AddQuadFace( CellIdentifier &id, PointIdentifier f[4], OutputMeshType* mesh );

  InputPixelType              m_IsoSurfaceValue;
  InterpolatorPointer         m_Interpolator;
  GradientInterpolatorPointer m_GradientInterpolator;
  SpacingValueType            m_MaxSpacing;
  bool                        m_GenerateTriangleFaces;
  bool                        m_ProjectVerticesToIsoSurface;
  double                      m_ProjectVertexSurfaceDistanceThreshold;
  double                      m_ProjectVertexStepLength;
  double                      m_ProjectVertexStepLengthRelaxationFactor;
  unsigned int                m_ProjectVertexMaximumNumberOfSteps;
#if DEBUG_PRINT
  unsigned int m_ProjectVertexTerminate[3];
#endif
};

} // end namespace itk

#ifndef ITK_MANUAL_INSTANTIATION
#include "itkCuberilleImageToMeshFilter.hxx"
#endif

#endif