lib/cpExtensions/Algorithms/RasterContourFilter.hxx

   1 // -------------------------------------------------------------------------
   2 // @author Leonardo Florez-Valencia (florez-l@javeriana.edu.co)
   3 // ------------------------------------------------------------------------
   4
   5 // Inclusion test taken from:
   6 // https://www.ecse.rpi.edu/Homepages/wrf/Research/Short_Notes/pnpoly.html
   7
   8 #ifndef __cpExtensions__Algorithms__RasterContourFilter__hxx__
   9 #define __cpExtensions__Algorithms__RasterContourFilter__hxx__
  10
  11 #include <itkImageRegionIteratorWithIndex.h>
  12 #include <vtkPolygon.h>
  13
  14 // -------------------------------------------------------------------------
  15 template< class _TImage >
  16 void cpExtensions::Algorithms::RasterContourFilter< _TImage >::
  17 AddPoint( double x, double y )
  18 {
  19   this->m_Polygon->GetPoints( )->InsertNextPoint( x, y, 0 );
  20   this->Modified( );
  21 }
  22
  23 // -------------------------------------------------------------------------
  24 template< class _TImage >
  25 void cpExtensions::Algorithms::RasterContourFilter< _TImage >::
  26 AddPoint( double p[ 2 ] )
  27 {
  28   this->m_Polygon->GetPoints( )->InsertNextPoint( p[ 0 ], p[ 1 ], 0 );
  29   this->Modified( );
  30 }
  31
  32 // -------------------------------------------------------------------------
  33 template< class _TImage >
  34 void cpExtensions::Algorithms::RasterContourFilter< _TImage >::
  35 AddPoints( vtkPoints* points )
  36 {
  37   double p[ 3 ];
  38   for( unsigned long i = 0; i < points->GetNumberOfPoints( ); ++i )
  39   {
  40     points->GetPoint( i, p );
  41     this->m_Polygon->GetPoints( )->InsertNextPoint( p[ 0 ], p[ 1 ], 0 );
  42
  43   } // rof
  44   this->Modified( );
  45 }
  46
  47 // -------------------------------------------------------------------------
  48 template< class _TImage >
  49 void cpExtensions::Algorithms::RasterContourFilter< _TImage >::
  50 ClearPoints( )
  51 {
  52   this->m_Polygon = vtkSmartPointer< vtkPolygon >::New( );
  53 }
  54
  55 // -------------------------------------------------------------------------
  56 template< class _TImage >
  57 cpExtensions::Algorithms::RasterContourFilter< _TImage >::
  58 RasterContourFilter( )
  59   : Superclass( ),
  60     m_InsideValue( TPixel( 1 ) ),
  61     m_OutsideValue( TPixel( 0 ) )
  62 {
  63   this->ClearPoints( );
  64 }
  65
  66 // -------------------------------------------------------------------------
  67 template< class _TImage >
  68 cpExtensions::Algorithms::RasterContourFilter< _TImage >::
  69 ~RasterContourFilter( )
  70 {
  71 }
  72
  73 // -------------------------------------------------------------------------
  74 template< class _TImage >
  75 void cpExtensions::Algorithms::RasterContourFilter< _TImage >::
  76 AllocateOutputs( )
  77 {
  78   _TImage* out = this->GetOutput( 0 );
  79   out->SetSpacing( this->m_Template->GetSpacing( ) );
  80   out->SetRegions( this->m_Template->GetRequestedRegion( ) );
  81   out->SetOrigin( this->m_Template->GetOrigin( ) );
  82   out->SetDirection( this->m_Template->GetDirection( ) );
  83   out->Allocate( );
  84 }
  85
  86 // -------------------------------------------------------------------------
  87 template< class _TImage >
  88 void cpExtensions::Algorithms::RasterContourFilter< _TImage >::
  89 BeforeThreadedGenerateData( )
  90 {
  91 }
  92
  93 // -------------------------------------------------------------------------
  94 template< class _TImage >
  95 void cpExtensions::Algorithms::RasterContourFilter< _TImage >::
  96 AfterThreadedGenerateData( )
  97 {
  98 }
  99
 100 // -------------------------------------------------------------------------
 101 template< class _TImage >
 102 void cpExtensions::Algorithms::RasterContourFilter< _TImage >::
 103 ThreadedGenerateData( const TRegion& region, itk::ThreadIdType id )
 104 {
 105   vtkPoints* points = this->m_Polygon->GetPoints( );
 106   unsigned long nPoints = points->GetNumberOfPoints( );
 107   double* arr =
 108     static_cast< double* >( points->GetData()->GetVoidPointer( 0 ) );
 109   double p[ 3 ], n[ 3 ], b[ 6 ];
 110   this->m_Polygon->GetPoints( )->GetBounds( b );
 111   n[ 0 ] = n[ 1 ] = double( 0 );
 112   n[ 2 ] = double( 1 );
 113
 114   _TImage* out = this->GetOutput( );
 115   itk::ImageRegionIteratorWithIndex< _TImage > iIt( out, region );
 116   for( iIt.GoToBegin( ); !iIt.IsAtEnd( ); ++iIt )
 117   {
 118     TIndex idx = iIt.GetIndex( );
 119     TPoint pnt;
 120     out->TransformIndexToPhysicalPoint( idx, pnt );
 121     p[ 0 ] = pnt[ 0 ];
 122     p[ 1 ] = pnt[ 1 ];
 123     p[ 2 ] = double( 0 );
 124
 125     int i, j, c = 0;
 126     int nvert = nPoints;
 127     for( i = 0, j = nvert - 1; i < nvert; j = i++ )
 128     {
 129       double pi[ 3 ], pj[ 3 ];
 130       points->GetPoint( i, pi );
 131       points->GetPoint( j, pj );
 132
 133       if(
 134         ( ( pi[ 1 ] > p[ 1 ] ) != ( pj[ 1 ] > p[ 1 ] ) ) &&
 135         ( p[ 0 ] < ( pj[ 0 ] - pi[ 0 ] ) * ( p[ 1 ] - pi[ 1 ] ) / ( pj[ 1 ] - pi[ 1 ] ) + pi[ 0 ] )
 136         )
 137         c = !c;
 138
 139     } // rof
 140
 141     if( c != 0 )
 142       iIt.Set( this->m_InsideValue );
 143     else
 144       iIt.Set( this->m_OutsideValue );
 145
 146   } // rof
 147 }
 148
 149 #endif // __cpExtensions__Algorithms__RasterContourFilter__hxx__
 150
 151 // eof - $RCSfile$