lib/fpa/Image/ExtractEndPointsAndBifurcationsFromMinimumSpanningTree.hxx

   1 #ifndef __FPA__IMAGE__EXTRACTENDPOINTSANDBIFURCATIONSFROMMINIMUMSPANNINGTREE__HXX__
   2 #define __FPA__IMAGE__EXTRACTENDPOINTSANDBIFURCATIONSFROMMINIMUMSPANNINGTREE__HXX__
   3
   4 #include <cmath>
   5 #include <cpPlugins_ITKInstances/ImageIterators.h>
   6
   7 // -------------------------------------------------------------------------
   8 template< class _TImage, class _TMST >
   9 const typename
  10 fpa::Image::ExtractEndPointsAndBifurcationsFromMinimumSpanningTree< _TImage, _TMST >::
  11 TImage*
  12 fpa::Image::ExtractEndPointsAndBifurcationsFromMinimumSpanningTree< _TImage, _TMST >::
  13 GetCostsImage( )
  14 {
  15   return(
  16     dynamic_cast< const _TImage* >( this->itk::ProcessObject::GetInput( 1 ) )
  17     );
  18 }
  19
  20 // -------------------------------------------------------------------------
  21 template< class _TImage, class _TMST >
  22 const typename
  23 fpa::Image::ExtractEndPointsAndBifurcationsFromMinimumSpanningTree< _TImage, _TMST >::
  24 TImage*
  25 fpa::Image::ExtractEndPointsAndBifurcationsFromMinimumSpanningTree< _TImage, _TMST >::
  26 GetDistanceMap( )
  27 {
  28   return(
  29     dynamic_cast< const _TImage* >( this->itk::ProcessObject::GetInput( 2 ) )
  30     );
  31 }
  32
  33 // -------------------------------------------------------------------------
  34 template< class _TImage, class _TMST >
  35 void
  36 fpa::Image::ExtractEndPointsAndBifurcationsFromMinimumSpanningTree< _TImage, _TMST >::
  37 SetCostsImage( TImage* image )
  38 {
  39   this->itk::ProcessObject::SetNthInput( 1, const_cast< _TImage* >( image ) );
  40 }
  41
  42 // -------------------------------------------------------------------------
  43 template< class _TImage, class _TMST >
  44 void
  45 fpa::Image::ExtractEndPointsAndBifurcationsFromMinimumSpanningTree< _TImage, _TMST >::
  46 SetDistanceMap( TImage* image )
  47 {
  48   this->itk::ProcessObject::SetNthInput( 2, const_cast< _TImage* >( image ) );
  49 }
  50
  51 // -------------------------------------------------------------------------
  52 template< class _TImage, class _TMST >
  53 fpa::Image::ExtractEndPointsAndBifurcationsFromMinimumSpanningTree< _TImage, _TMST >::
  54 ExtractEndPointsAndBifurcationsFromMinimumSpanningTree( )
  55   : Superclass( ),
  56     m_SquaredDistanceMap( false )
  57 {
  58   this->SetNumberOfRequiredInputs( 3 );
  59 }
  60
  61 // -------------------------------------------------------------------------
  62 template< class _TImage, class _TMST >
  63 fpa::Image::ExtractEndPointsAndBifurcationsFromMinimumSpanningTree< _TImage, _TMST >::
  64 ~ExtractEndPointsAndBifurcationsFromMinimumSpanningTree( )
  65 {
  66 }
  67
  68 // -------------------------------------------------------------------------
  69 template< class _TImage, class _TMST >
  70 void
  71 fpa::Image::ExtractEndPointsAndBifurcationsFromMinimumSpanningTree< _TImage, _TMST >::
  72 GenerateData( )
  73 {
  74   auto image = this->GetCostsImage( );
  75
  76   this->m_MarkImage = TMarkImage::New( );
  77   this->m_MarkImage->SetLargestPossibleRegion( image->GetLargestPossibleRegion( ) );
  78   this->m_MarkImage->SetRequestedRegion( image->GetRequestedRegion( ) );
  79   this->m_MarkImage->SetBufferedRegion( image->GetBufferedRegion( ) );
  80   this->m_MarkImage->SetDirection( image->GetDirection( ) );
  81   this->m_MarkImage->SetOrigin( image->GetOrigin( ) );
  82   this->m_MarkImage->SetSpacing( image->GetSpacing( ) );
  83   this->m_MarkImage->Allocate( );
  84   this->m_MarkImage->FillBuffer( 0 );
  85
  86   this->m_SkeletonImage = TMarkImage::New( );
  87   this->m_SkeletonImage->SetLargestPossibleRegion( image->GetLargestPossibleRegion( ) );
  88   this->m_SkeletonImage->SetRequestedRegion( image->GetRequestedRegion( ) );
  89   this->m_SkeletonImage->SetBufferedRegion( image->GetBufferedRegion( ) );
  90   this->m_SkeletonImage->SetDirection( image->GetDirection( ) );
  91   this->m_SkeletonImage->SetOrigin( image->GetOrigin( ) );
  92   this->m_SkeletonImage->SetSpacing( image->GetSpacing( ) );
  93   this->m_SkeletonImage->Allocate( );
  94   this->m_SkeletonImage->FillBuffer( 0 );
  95
  96   this->Superclass::GenerateData( );
  97 }
  98
  99 // -------------------------------------------------------------------------
 100 template< class _TImage, class _TMST >
 101 void
 102 fpa::Image::ExtractEndPointsAndBifurcationsFromMinimumSpanningTree< _TImage, _TMST >::
 103 _MarkSkeleton( const TVertex& v, const unsigned long& l )
 104 {
 105   if( this->m_SkeletonImage->GetRequestedRegion( ).IsInside( v ) )
 106     this->m_SkeletonImage->SetPixel( v, l );
 107 }
 108
 109 // -------------------------------------------------------------------------
 110 template< class _TImage, class _TMST >
 111 void
 112 fpa::Image::ExtractEndPointsAndBifurcationsFromMinimumSpanningTree< _TImage, _TMST >::
 113 _MarkSphere(
 114   const TVertex& v,
 115   const double& r,
 116   const unsigned long& l
 117   )
 118 {
 119   // Get marking region
 120   auto rreg = this->m_MarkImage->GetRequestedRegion( );
 121   auto spac = this->m_MarkImage->GetSpacing( );
 122   TVertex ireg = rreg.GetIndex( );
 123   TVertex jreg = ireg + rreg.GetSize( );
 124   TVertex idx, jdx;
 125   typename _TImage::SizeType size;
 126   bool all = false;
 127   for( unsigned int d = 0; d < _TImage::ImageDimension; ++d )
 128   {
 129     unsigned long s = std::ceil( r / double( spac[ d ] ) );
 130     if( s < 3 )
 131     {
 132       s = 3;
 133       all = true;
 134
 135     } // fi
 136     s += 3;
 137     idx[ d ] = v[ d ] - s;
 138     jdx[ d ] = v[ d ] + s;
 139
 140     jreg[ d ]--;
 141     if( idx[ d ] < ireg[ d ] ) idx[ d ] = ireg[ d ];
 142     if( idx[ d ] > jreg[ d ] ) idx[ d ] = jreg[ d ];
 143     if( jdx[ d ] < ireg[ d ] ) jdx[ d ] = ireg[ d ];
 144     if( jdx[ d ] > jreg[ d ] ) jdx[ d ] = jreg[ d ];
 145     size[ d ] = jdx[ d ] - idx[ d ] + 1;
 146
 147   } // rof
 148   typename _TImage::RegionType region;
 149   region.SetIndex( idx );
 150   region.SetSize( size );
 151
 152   // Mark region
 153   typename _TImage::PointType c;
 154   this->m_MarkImage->TransformIndexToPhysicalPoint( v, c );
 155   itk::ImageRegionIteratorWithIndex< TMarkImage >
 156     spIt( this->m_MarkImage, region );
 157   for( spIt.GoToBegin( ); !spIt.IsAtEnd( ); ++spIt )
 158   {
 159     typename _TImage::PointType pnt;
 160     this->m_MarkImage->TransformIndexToPhysicalPoint( spIt.GetIndex( ), pnt );
 161     // TODO: if( double( pnt.EuclideanDistanceTo( c ) ) <= r || all )
 162     // TODO: if( spIt.Get( ) == 0 )
 163       spIt.Set( l );
 164
 165   } // rof
 166 }
 167
 168 // -------------------------------------------------------------------------
 169 template< class _TImage, class _TMST >
 170 unsigned long
 171 fpa::Image::ExtractEndPointsAndBifurcationsFromMinimumSpanningTree< _TImage, _TMST >::
 172 _Mark( const TVertex& v )
 173 {
 174   return( this->m_MarkImage->GetPixel( v ) );
 175 }
 176
 177 // -------------------------------------------------------------------------
 178 template< class _TImage, class _TMST >
 179 unsigned long fpa::Image::ExtractEndPointsAndBifurcationsFromMinimumSpanningTree< _TImage, _TMST >::
 180 _SkeletonMark( const TVertex& v )
 181 {
 182   return( this->m_SkeletonImage->GetPixel( v ) );
 183 }
 184
 185 // -------------------------------------------------------------------------
 186 template< class _TImage, class _TMST >
 187 double
 188 fpa::Image::ExtractEndPointsAndBifurcationsFromMinimumSpanningTree< _TImage, _TMST >::
 189 _Radius( const TVertex& v )
 190 {
 191   double d = double( this->GetDistanceMap( )->GetPixel( v ) );
 192   double s = std::fabs( d );
 193   if( s > double( 0 ) )
 194     s /= d;
 195   if( this->m_SquaredDistanceMap )
 196     return( s * std::sqrt( std::fabs( d ) ) );
 197   else
 198     return( d );
 199 }
 200
 201 #endif // __FPA__IMAGE__EXTRACTENDPOINTSANDBIFURCATIONSFROMMINIMUMSPANNINGTREE__HXX__
 202
 203 // eof - $RCSfile$