lib/ivq/ITK/Simple3DCurve.cxx

   1 /* =======================================================================
   2  * @author: Leonardo Florez-Valencia
   3  * @email: florez-l@javeriana.edu.co
   4  * =======================================================================
   5  */
   6
   7 #include <ivq/ITK/Simple3DCurve.h>
   8
   9 // -------------------------------------------------------------------------
  10 template< class _TScalar >
  11 void ivq::ITK::Simple3DCurve< _TScalar >::
  12 Modified( ) const
  13 {
  14   this->Superclass::Modified( );
  15   this->m_FramesUpdated = false;
  16 }
  17
  18 // -------------------------------------------------------------------------
  19 template< class _TScalar >
  20 void ivq::ITK::Simple3DCurve< _TScalar >::
  21 Clear( )
  22 {
  23   this->m_Points.clear( );
  24   this->m_Frames.clear( );
  25   this->Modified( );
  26 }
  27
  28 // -------------------------------------------------------------------------
  29 template< class _TScalar >
  30 void ivq::ITK::Simple3DCurve< _TScalar >::
  31 Smooth( const unsigned int& kernel_size )
  32 {
  33   std::vector< TPoint > new_points;
  34
  35   TPoint p;
  36   int k = int( kernel_size );
  37   unsigned int nPoints = this->m_Points.size( );
  38   for( unsigned int i = 0; i < nPoints; ++i )
  39   {
  40     typename TPoint::VectorType v( TScalar( 0 ) );
  41     for( int j = -k; j <= k; ++j )
  42     {
  43       int l = int( i ) + j;
  44       if     ( l < 0        ) p = this->m_Points[ 0 ];
  45       else if( l >= nPoints ) p = this->m_Points[ nPoints - 1 ];
  46       else                    p = this->m_Points[ l ];
  47       v += p.GetVectorFromOrigin( );
  48
  49     } // rof
  50     v /= _TScalar( ( k << 1 ) + 1 );
  51
  52     p.Fill( 0 );
  53     p += v;
  54     new_points.push_back( p );
  55
  56   } // rof
  57   this->m_Points = new_points;
  58   this->m_Frames.clear( );
  59   this->Modified( );
  60 }
  61
  62 // -------------------------------------------------------------------------
  63 template< class _TScalar >
  64 unsigned long ivq::ITK::Simple3DCurve< _TScalar >::
  65 GetNumberOfPoints( ) const
  66 {
  67   return( this->m_Points.size( ) );
  68 }
  69
  70 // -------------------------------------------------------------------------
  71 template< class _TScalar >
  72 const typename ivq::ITK::Simple3DCurve< _TScalar >::
  73 TMatrix& ivq::ITK::Simple3DCurve< _TScalar >::
  74 GetFrame( unsigned int id ) const
  75 {
  76   if( !( this->m_FramesUpdated ) )
  77   {
  78     unsigned long N = this->m_Points.size( );
  79
  80     std::vector< TVector > tg( N );
  81     tg[ 0 ] = this->m_Points[ 1 ] - this->m_Points[ 0 ];
  82     tg[ N - 1 ] = this->m_Points[ N - 1 ] - this->m_Points[ N - 2 ];
  83     for( unsigned int i = 1; i < N - 1; ++i )
  84       tg[ i ] = this->m_Points[ i + 1 ] - this->m_Points[ i - 1 ];
  85
  86     std::vector< TVector > no( N );
  87     std::vector< TVector > bn( N );
  88     TVector prev_tg( _TScalar( 0 ) ), prev_no( _TScalar( 0 ) );
  89     prev_tg[ 0 ] = _TScalar( 1 );
  90     prev_no[ 1 ] = _TScalar( 1 );
  91     for( unsigned int i = 0; i < N; ++i )
  92     {
  93       auto ntg = tg[ i ].GetNorm( );
  94       if( double( ntg ) > double( 0 ) )
  95         tg[ i ] /= ntg;
  96
  97       _TScalar ct = prev_tg * tg[ i ];
  98       TVector a = itk::CrossProduct( prev_tg, tg[ i ] );
  99       _TScalar st = a.GetNorm( );
 100       if( st > _TScalar( 0 ) )
 101         a /= st;
 102
 103       no[ i ] =
 104         ( prev_no * ct ) +
 105         ( itk::CrossProduct( a, prev_no ) * st ) +
 106         ( a * ( ( a * prev_no ) * ( _TScalar( 1 ) - ct ) ) );
 107       auto nno = no[ i ].GetNorm( );
 108       if( double( nno ) > double( 0 ) )
 109         no[ i ] /= nno;
 110
 111       bn[ i ] = itk::CrossProduct( tg[ i ], no[ i ] );
 112
 113       prev_tg = tg[ i ];
 114       prev_no = no[ i ];
 115
 116     } // rof
 117
 118     this->m_Frames.clear( );
 119     for( unsigned int i = 0; i < N; ++i )
 120     {
 121       TMatrix m;
 122       for( unsigned int d = 0; d < 3; ++d )
 123       {
 124         m[ d ][ 0 ] = tg[ i ][ d ];
 125         m[ d ][ 1 ] = no[ i ][ d ];
 126         m[ d ][ 2 ] = bn[ i ][ d ];
 127
 128       } // rof
 129       this->m_Frames.push_back( m );
 130
 131     } // rof
 132
 133     this->m_FramesUpdated = true;
 134
 135   } // fi
 136   return( this->m_Frames[ id ] );
 137 }
 138
 139 // -------------------------------------------------------------------------
 140 template< class _TScalar >
 141 const typename ivq::ITK::Simple3DCurve< _TScalar >::
 142 TPoint& ivq::ITK::Simple3DCurve< _TScalar >::
 143 GetPoint( unsigned int id ) const
 144 {
 145   return( this->m_Points[ id ] );
 146 }
 147
 148 // -------------------------------------------------------------------------
 149 template< class _TScalar >
 150 typename ivq::ITK::Simple3DCurve< _TScalar >::
 151 TVector ivq::ITK::Simple3DCurve< _TScalar >::
 152 GetVector( unsigned int id ) const
 153 {
 154   return( this->m_Points[ id ].GetVectorFromOrigin( ) );
 155 }
 156
 157 // -------------------------------------------------------------------------
 158 template< class _TScalar >
 159 ivq::ITK::Simple3DCurve< _TScalar >::
 160 Simple3DCurve( )
 161   : Superclass( ),
 162     m_FramesUpdated( false )
 163 {
 164 }
 165
 166 // -------------------------------------------------------------------------
 167 template< class _TScalar >
 168 ivq::ITK::Simple3DCurve< _TScalar >::
 169 ~Simple3DCurve( )
 170 {
 171 }
 172
 173 // -------------------------------------------------------------------------
 174 template class ivq::ITK::Simple3DCurve< float >;
 175 template class ivq::ITK::Simple3DCurve< double >;
 176
 177 // eof - $RCSfile$