FreemanChain.ih

// IMPLEMENTATION of inline methods.




// Iterator on points


template <typename TInteger>
inline
  DGtal::FreemanChain<TInteger>::ConstIterator::ConstIterator( 
      const FreemanChain & aChain, unsigned int n)
  : myFc( &aChain ), myPos( 0 )
{
  if ( n < myFc->chain.size() ) 
  {
    myXY.at(0)=aChain.x0;
    myXY.at(1)=aChain.y0;
    while ( myPos < n ) this->next();
  } 
  else 
  {// iterator end() 
    myXY.at(0)=aChain.xn;
    myXY.at(1)=aChain.yn;
    myPos = (typename DGtal::FreemanChain<TInteger>::Index)(myFc->chain.size()+1);
  }
}

template <typename TInteger>
inline
typename DGtal::FreemanChain<TInteger>::ConstIterator & 
DGtal::FreemanChain<TInteger>::ConstIterator::operator= ( const ConstIterator & other )
{       
  if ( this != &other )
  {
    myFc = other.myFc;
    myPos = other.myPos;
    myXY = other.myXY;
  }
  return *this;
}


template <typename TInteger>
inline
void DGtal::FreemanChain<TInteger>::ConstIterator::next()
{
  if ( myPos < myFc->chain.size() )
  {
    switch ( myFc->code( myPos ) )
    {
      case '0': (myXY.at(0))++; break;
      case '1': (myXY.at(1))++; break;
      case '2': (myXY.at(0))--; break;
      case '3': (myXY.at(1))--; break;
    }
    ++myPos;
  } else ++myPos;
}




template <typename TInteger>
inline
void DGtal::FreemanChain<TInteger>::ConstIterator::nextInLoop()
{
  if ( myPos < myFc->chain.size() )
  {
    switch ( myFc->code( myPos ) )
    {
      case '0': (myXY.at(0))++; break;
      case '1': (myXY.at(1))++; break;
      case '2': (myXY.at(0))--; break;
      case '3': (myXY.at(1))--; break;
    }
    myPos = ( myPos + 1 ) % myFc->chain.size();
  }
}


template <typename TInteger>
inline
void DGtal::FreemanChain<TInteger>::ConstIterator::previous()
{
  if ( myPos == myFc->chain.size() + 1 ) 
  {
    myXY = myFc->lastPoint();
    --myPos;
  } 
  else 
  {
    if ( myPos >= 1 )
    --myPos;
    if ( myPos < myFc->chain.size() )  
    {
      switch ( myFc->code( myPos ) ) 
      {
        case '0': (myXY.at(0))--; break;
        case '1': (myXY.at(1))--; break;
        case '2': (myXY.at(0))++; break;
        case '3': (myXY.at(1))++; break;
      }
    }
  }
}


template <typename TInteger>
inline
void DGtal::FreemanChain<TInteger>::ConstIterator::previousInLoop()
{
  if ( myPos == 0 ) myPos = (typename DGtal::FreemanChain<TInteger>::Index)(myFc->chain.size() - 1);
  else --myPos;
  switch ( myFc->code( myPos ) )
  {
    case '0': (myXY.at(0))--; break;
    case '1': (myXY.at(1))--; break;
    case '2': (myXY.at(0))++; break;
    case '3': (myXY.at(1))++; break;
  }
}



// ----------------------- Standard services ------------------------------



template <typename TInteger>
DGtal::FreemanChain<TInteger>::FreemanChain( const std::string & s, TInteger x, TInteger y )
  : chain( s ), x0( x ), y0( y ), xn( x ), yn( y )
{
  DGtal::FreemanChain<TInteger>::computeLastPoint();
}


template <typename TInteger>
DGtal::FreemanChain<TInteger>::
FreemanChain( const std::vector<Point>& vectPoints )
{
  readFromPointsRange(vectPoints.begin(), vectPoints.end(), *this); 
}



template <typename TInteger>
DGtal::FreemanChain<TInteger>::FreemanChain(std::istream & in ){
  read(in, *this);
  computeLastPoint();
}


template <typename TInteger>
DGtal::FreemanChain<TInteger>::FreemanChain( const FreemanChain<TInteger> & aOther )
  : chain( aOther.chain ), x0( aOther.x0 ), y0( aOther.y0 ), 
    xn( aOther.xn), yn( aOther.yn)
{ }




template <typename TInteger>
typename DGtal::FreemanChain<TInteger> & 
DGtal::FreemanChain<TInteger>::operator=( const FreemanChain<TInteger> & aOther )
{
  if ( this != &aOther )
    {
      chain = aOther.chain;
      x0 = aOther.x0;
      y0 = aOther.y0;
      xn = aOther.xn;
      yn = aOther.yn;
    }
  return *this;
}





// Iterators services


template <typename TInteger>
typename DGtal::FreemanChain<TInteger>::ConstIterator 
DGtal::FreemanChain<TInteger>::begin() const
{
  return ConstIterator( *this, 0 );
}


template <typename TInteger>
typename DGtal::FreemanChain<TInteger>::ConstIterator 
DGtal::FreemanChain<TInteger>::end() const
{
  Point p(0,0); // *(end()) is invalid
  return ConstIterator( *this, (typename DGtal::FreemanChain<TInteger>::Index)(chain.size() + 1), p );
}


template <typename TInteger>
typename DGtal::FreemanChain<TInteger>::Index
DGtal::FreemanChain<TInteger>::next( Index aPos ) const
{
  ++aPos;
  if ( aPos >= chain.size() )
    aPos = 0;
  return aPos;
}

template <typename TInteger>
typename DGtal::FreemanChain<TInteger>::Index
DGtal::FreemanChain<TInteger>::previous( Index aPos ) const
{
  if ( aPos == 0 ) aPos = chain.size() - 1;
  else --aPos;
  return aPos;
}

template <typename TInteger>
char
DGtal::FreemanChain<TInteger>::code( Index aPos ) const
{
  ASSERT( aPos < chain.size() );
  //return chain[ aPos ] - '0';
  return chain[ aPos ];
}


template <typename TInteger>
inline
typename DGtal::FreemanChain<TInteger>::Index
DGtal::FreemanChain<TInteger>::size() const
{
  return (typename DGtal::FreemanChain<TInteger>::Size)(chain.size());
}


template <typename TInteger>
inline
void DGtal::FreemanChain<TInteger>::computeBoundingBox( TInteger & min_x,
    TInteger& min_y, TInteger& max_x, TInteger& max_y ) const
{
  min_x = max_x = x0;
  min_y = max_y = y0;
  for ( ConstIterator it = begin();
      it != end();
      ++it )
  {
    Point p( *it );
    if ( p.at(0) < min_x )
      min_x = p.at(0);
    else
      if ( p.at(0) > max_x )
        max_x = p.at(0);
    if ( p.at(1) < min_y )
      min_y = p.at(1);
    else
      if ( p.at(1) > max_y )
        max_y = p.at(1);
  }
}


template <typename TInteger>
inline
typename FreemanChain<TInteger>::ConstIterator
DGtal::FreemanChain<TInteger>::findQuadrantChange( OrderedAlphabet & A ) const
{
  ConstIterator it = begin();
  ConstIterator it_end = end();
  // find first letters a and b.
  char code1 = it.getCode();
  it.next();
  while ( ( it != it_end ) && ( it.getCode() == code1 ) )
    it.next();
  ASSERT( ( it != it_end )
      && "[DGtal::FreemanChain::findQuadrantChange( OrderedAlphabet & A ) const] 1-letter freeman chain." );
  char  code2 = it.getCode();
  // find third letter c.
  while ( ( it != it_end ) && ( ( it.getCode() == code1 )
        || ( it.getCode() == code2 ) ) )
    it.next();
  ASSERT( ( it != it_end )
      && "[DGtal::FreemanChain::findQuadrantChange( OrderedAlphabet & A ) const] 2-letters Freeman chain." );
  // reorder a and b.
  it.previous();
  if ( it.getCode() != code2 )
    swap( code1, code2 );
  // find first a.
  do
  {
    it.previous();
  }
  while ( it.getCode() == code2 );
  char a_char = chain[ it.getPosition() ];
  // the next is the first b.
  it.next();
  char b_char = chain[ it.getPosition() ];
  // Reorder the alphabet to match the quadrant change.
  while ( A.order( b_char ) != 1 )
    A.shiftLeft();
  if ( A.order( a_char ) == 0 )
  {
    A.reverse();
    while ( A.order( b_char ) != 1 )
      A.shiftLeft();
  }
  ASSERT( ( A.order( b_char ) == 1 )
      && ( A.order( a_char ) == 2 )
      && "[DGtal::FreemanChain::findQuadrantChange( OrderedAlphabet & A ) const] Internal error: invalid Quadrant change found." );
  return it;
}


template <typename TInteger>
inline
typename FreemanChain<TInteger>::ConstIterator
DGtal::FreemanChain<TInteger>::findQuadrantChange4( OrderedAlphabet & A ) const
{
  ConstIterator it = begin();
  ConstIterator it_end = end();
  // find first letters a and b.
  uint8_t code1 = it.getCode();
  it.next();
  while ( ( it != it_end ) && ( it.getCode() == code1 ) )
    it.next();
  ASSERT( ( it != it_end )
      && "[DGtal::FreemanChain::findQuadrantChange( OrderedAlphabet & A ) const] 1-letter freeman chain." );
  uint8_t code2 = it.getCode();
  // find third letter c.
  while ( ( it != it_end ) && ( ( it.getCode() == code1 )
        || ( it.getCode() == code2 ) ) )
    it.next();
  ASSERT( ( it != it_end )
      && "[DGtal::FreemanChain::findQuadrantChange( OrderedAlphabet & A ) const] 2-letters Freeman chain." );
  uint8_t code3 = it.getCode();
  // find fourth letter d.
  while ( ( it != it_end ) && ( ( it.getCode() == code1 )
        || ( it.getCode() == code2 )
        || ( it.getCode() == code3 ) ) )
    it.next();
  ASSERT( ( it != it_end )
      && "[DGtal::FreemanChain::findQuadrantChange( OrderedAlphabet & A ) const] 3-letters Freeman chain." );
  // define true c.
  it.previous();
  code3 = it.getCode();
  // find first b.
  do
  {
    it.previous();
  }
  while ( it.getCode() == code3 );
  char a_char = chain[ it.getPosition() ];
  // the next is the first c.
  it.next();
  char b_char = chain[ it.getPosition() ];
  // Reorder the alphabet to match the quadrant change.
  while ( A.order( b_char ) != 1 )
    A.shiftLeft();
  if ( A.order( a_char ) == 0 )
  {
    A.reverse();
    while ( A.order( b_char ) != 1 )
      A.shiftLeft();
  }
  ASSERT( ( A.order( b_char ) == 1 )
      && ( A.order( a_char ) == 2 )
      && "[DGtal::FreemanChain::findQuadrantChange( OrderedAlphabet & A ) const] Internal error: invalid Quadrant change found." );
  return it;
}


template <typename TInteger>
inline
int DGtal::FreemanChain<TInteger>::isClosed() const
{
  return (x0 == xn) && (y0 == yn);
}


template <typename TInteger>
inline
int DGtal::FreemanChain<TInteger>::ccwLoops() const
{
  ConstIterator it = this->begin();
  ConstIterator it_end = this->end();
  --it_end;
  ConstIterator it_suiv = it;
  Point spos = *it;
  int nb_ccw_turns = 0;
  while ( it != it_end )
  {
    int code1 = it.getCode();
    it_suiv.nextInLoop();
    int code2 = it_suiv.getCode();
    char diff = ( code2 - code1 + 4 ) % 4;
    if ( diff == 1 )
      ++nb_ccw_turns;
    else
      if ( diff == 3 )
        --nb_ccw_turns;
      else
        if ( diff == 2 )
          return 0;
    ++it;
  }
  if ( spos == *it_suiv )
    return nb_ccw_turns / 4;
  else
    return 0;
}


template <typename TInteger>
inline
typename DGtal::FreemanChain<TInteger>
DGtal::FreemanChain<TInteger>::subChain(Index pos, Size n) const
{
  Self newChain;
  Point startPoint = getPoint(pos);
  Point endPoint = getPoint(pos+n);
  newChain.chain = chain.substr(pos,n);
  newChain.x0 = startPoint.at(0);
  newChain.y0 = startPoint.at(1);
  newChain.xn = endPoint.at(0);
  newChain.yn = endPoint.at(1);
  return newChain;
}


template <typename TInteger>
inline
typename DGtal::FreemanChain<TInteger>
DGtal::FreemanChain<TInteger>::operator+(const Self& aOther) const
{
  Self newChain;
  newChain.chain = chain + aOther.chain;
  newChain.x0 = x0;
  newChain.y0 = y0;
  Point newLastPoint = lastPoint() + ( aOther.lastPoint() - aOther.firstPoint() );
  newChain.xn = newLastPoint.at(0);
  newChain.yn = newLastPoint.at(1);
  return newChain;
}


template <typename TInteger>
inline
typename DGtal::FreemanChain<TInteger>&
DGtal::FreemanChain<TInteger>::operator+=(const Self& aOther) 
{
  chain += aOther.chain;
  Point newLastPoint = lastPoint() + ( aOther.lastPoint() - aOther.firstPoint() );
  xn = newLastPoint.at(0);
  yn = newLastPoint.at(1);
  return *this;
}


template <typename TInteger>
inline
typename DGtal::FreemanChain<TInteger>::Point 
DGtal::FreemanChain<TInteger>::getPoint(Index aPos) const
{
  ConstIterator it;
  Size n = this->size();
  if ( aPos < n / 2 ) {
    it = begin();
    for (unsigned int i=0; i<aPos; ++i) 
      it++;
  }
  else 
  {
    it = end();
    it--;
    n = n-aPos;
    for (unsigned int i=0; i<n; ++i) {
      it--;
    }
  }
  return *it;
}



template <typename TInteger>
inline
typename DGtal::FreemanChain<TInteger> &
DGtal::FreemanChain<TInteger>::extend(char aCode)
{
  chain += std::string( &aCode, 1);
  int dx=0, dy=0;
  //displacement( dx, dy, aCode - '0' );
  displacement( dx, dy, aCode );
  xn += dx;
  yn += dy;
  return *this;
}





template <typename TInteger>
inline
typename DGtal::FreemanChain<TInteger> &
DGtal::FreemanChain<TInteger>::retract(Size n)
{
  ConstIterator it = end();
  for (unsigned int i=0; i<=n; ++i)
    it--;
  xn = (*it).at(0);
  yn = (*it).at(1);
  Size mySize = size();
  ASSERT( (n <= mySize) && "Tried to shorten a FreemanChain by more then its length");
  chain.resize( mySize - n );
  return *this;
}




// ----------------------- Interface --------------------------------------


template <typename TInteger>
inline
void DGtal::FreemanChain<TInteger>::selfDisplay ( std::ostream & out ) const
{
  out << x0 << " " << y0 << " "  << chain;
}


template <typename TInteger>
inline
bool DGtal::FreemanChain<TInteger>::isValid () const
{
  return true;
}









// Static services


template <typename TInteger>
inline
void DGtal::FreemanChain<TInteger>::read( std::istream & in, FreemanChain & c ) 
{
  string str;
  while ( true )
  {
    getline( in, str );
    if ( ! in.good() )
      return;
    if ( ( str.size() > 0 ) && ( str[ 0 ] != '#' ) )
    {
      istringstream str_in( str );
      str_in >> c.x0 >> c.y0 >> c.chain;
      return;
    }
  }
}

template <typename TInteger>
template <typename TConstIterator>
inline
void DGtal::FreemanChain<TInteger>::readFromPointsRange( const TConstIterator& itBegin, const TConstIterator& itEnd, FreemanChain & c ) 
{
  TConstIterator it( itBegin ); 

  if (it != itEnd) 
  { //if the range is not empty
    Point startingPt( *it ); 
    ++it; 
    
    if (it != itEnd)
    { //if there is more than one element
      ostringstream oss;
      Point pt( startingPt );
      do 
      {
        
        Point ptSuiv( *it ); 
        Integer dx = ptSuiv[0] - pt[0]; 
        Integer dy = ptSuiv[1] - pt[1]; 
        short number = ( dx != 0 ? (1 - dx) : (2 - dy) );
        if ( (number < 0) || (number > 3) )
        {
          std::cerr << "not connected points (method readFromPointsRange of FreemanChain)" << std::endl;
          throw ConnectivityException(); 
        }
        oss << number;
        pt = ptSuiv;
        ++it;         
        
      } while ( it != itEnd ); 
      
      c.xn=pt[0];
      c.yn=pt[1];
      c.chain = oss.str(); 
    }
    
    c.x0=startingPt[0];
    c.y0=startingPt[1];  
  }

}

template <typename TInteger>
template <typename TRange>
inline
void DGtal::FreemanChain<TInteger>::readFromPointsRange( const TRange& aRange, FreemanChain & c ) 
{
  BOOST_CONCEPT_ASSERT(( CConstSinglePassRange<TRange> ));
  readFromPointsRange( aRange.begin(), aRange.end(), c ); 
}

template <typename TInteger>
inline
void DGtal::FreemanChain<TInteger>::getContourPoints(
    const FreemanChain & fc, std::vector<Point> & aVContour) 
{
  aVContour.clear();
  for ( ConstIterator it = fc.begin();
      it != fc.end();
      ++it )
  {
    aVContour.push_back(*it);
  }
}



template <typename TInteger>
void 
DGtal::FreemanChain<TInteger>::movePointFromFC(Point & aPoint, char aCode ){
  switch ( aCode )
    {
    case '0': aPoint.at(0)++; break;
    case '1': aPoint.at(1)++; break;
    case '2': aPoint.at(0)--; break;
    case '3': aPoint.at(1)--; break;
    }
}

// Depreciated
//
//template <typename TInteger>
//void 
//DGtal::FreemanChain<TInteger>::alphabet( char & aZero, char & aOne, char aQuadrant )
//{
//  switch ( aQuadrant )
//    {
//      case '0':
//        aZero = '0';
//        aOne = '1';
//        break;
//      case '1':
//        aZero = '1';
//        aOne = '2';
//        break;
//      case '2':
//        aZero = '2';
//        aOne = '3';
//        break;
//      case '3':
//        aZero = '3';
//        aOne = '0';
//        break;
//    }
//
//};



template <typename TInteger>
inline
char DGtal::FreemanChain<TInteger>::movement( char aCode1,
    char aCode2, bool ccw )
{
  unsigned int cfg = ( ccw ? 0 : 16 ) + ( (aCode1 - '0') << 2 ) + (aCode2 - '0');
  static const char tbl[ 32 ] =
  {
    '2', '1', '0', '3', '3', '2', '1', '0',
    '0', '3', '2', '1', '1', '0', '3', '2',
    '2', '3', '0', '1', '1', '2', '3', '0',
    '0', '1', '2', '3', '3', '0', '1', '2'
  };
  return tbl[ cfg ];
}



template <typename TInteger>
inline
char DGtal::FreemanChain<TInteger>::addToCode( char code, int n)
{
  return '0' + ( ( (code - '0' ) + n ) & 3 );
}





template <typename TInteger>
inline
void
DGtal::FreemanChain<TInteger>::displacement( int & dx, int & dy, char aCode )
{
  switch ( aCode ) 
    {
    case '0': dx = 1; dy = 0; break;
    case '1': dx = 0; dy = 1; break;
    case '2': dx = -1; dy = 0; break;
    case '3': dx = 0; dy = -1; break;
    }
}


template <typename TInteger>
inline
typename DGtal::PointVector<2,TInteger>
DGtal::FreemanChain<TInteger>::displacement( char aCode )
{
#ifdef CPP11_INITIALIZER_LIST  
  switch ( aCode ) 
  {
    case '0':  return Point({1,0}); 
    case '1':  return Point({0,1}); 
    case '2':  return Point({-1,0});
    case '3':  return Point({0,-1});
  }
  return  Point({0,0}); 
#else
  int t[2] = {0,0};

  switch ( aCode ) 
  {
    case '0':  t[0]++; return Point( t ); 
    case '1':  t[1]++; return Point( t ); 
    case '2':  t[0]--; return Point( t );
    case '3':  t[1]--; return Point( t );
  }
  return  Point( t ); 
#endif //CPP11_INITIALIZER_LIST
}



template <typename TInteger>
inline
char
DGtal::FreemanChain<TInteger>::turnedCode( char aCode, bool ccw )
{
  if ( ccw ) return ( aCode == '3' ) ? '0' : ( aCode + 1 );
  else       return ( aCode == '0' ) ? '3' : ( aCode - 1 );
}
//DGtal::FreemanChain<TInteger>::turnedCode( unsigned int aCode, bool ccw )
//{
//  if ( ccw ) return ( aCode + 1 ) & 0x3;
//  else       return ( aCode - 1 ) & 0x3;
//}







template <typename TInteger>
void DGtal::FreemanChain<TInteger>::innerContour( 
    FreemanChain & aInnerChain,
    std::vector<unsigned int> & aOuter2inner,
    std::vector<unsigned int> & aInner2outer,
    const FreemanChain & aOuterChain,
    bool ccw )
{
  unsigned int nb = (unsigned int)aOuterChain.chain.size();
  unsigned int j = 0;
  aOuter2inner.clear();
  aOuter2inner.reserve( nb );
  // aInnerChain.chain.reserve( nb + 4 );
  aInnerChain.chain = "";
  aInner2outer.clear();
  aInner2outer.reserve( nb + ( ccw ? 4 : -4 ) );
  int dx0=0, dy0=0;
  int dx1=0, dy1=0;
  FreemanChain<TInteger>::displacement( dx0, dy0, aOuterChain.code( 0 ) );
  int turn = ccw ? 1 : 3;
  //FreemanChain<TInteger>::displacement( dx1, dy1, ( aOuterChain.code( 0 ) + turn ) % 4 );
  FreemanChain<TInteger>::displacement( dx1, dy1, addToCode( aOuterChain.code( 0 ) , turn ) );
  dx0 += dx1;
  dy0 += dy1;
  aInnerChain.x0 = dx0 > 0 ? aOuterChain.x0 : aOuterChain.x0 - 1;
  aInnerChain.y0 = dy0 > 0 ? aOuterChain.y0 : aOuterChain.y0 - 1;

  ConstIterator it_begin = aOuterChain.begin();
  ConstIterator it = it_begin;
  it.next();
  for ( unsigned int i = 0; i < nb; ++i )
  {
    // Check if contour is open.
    // cerr << "i=" << i << " code=" << aOuterChain.code( i ) << endl;
    switch ( movement( aOuterChain.code( i ), aOuterChain.code( ( i + 1 ) % nb ), ccw ) )
    {
      case '0':
        // contour going in then out.
        aInnerChain.chain += aOuterChain.chain[ i ];
        //aInnerChain.chain += ( ( ( (unsigned int) ( aOuterChain.chain[ i ] - '0' )
        //        + ( ccw ? 3 : 1 ) ) )
        //    % 4 ) + '0';
        aInnerChain.chain += addToCode ( aOuterChain.chain[ i ],  (ccw) ? 3 : 1  );
        aInnerChain.chain += aOuterChain.chain[ ( i + 1 ) % nb ];
        aOuter2inner.push_back( j );
        aInner2outer.push_back( i );
        aInner2outer.push_back( i + 1 );
        aInner2outer.push_back( i + 1 );
        j += 3;
        break;

      case '1':
        // contour turning toward its inside.
        aOuter2inner.push_back( j );
        break;

      case '2':
        // contour going straight ahead
        aInnerChain.chain += aOuterChain.chain[ i ];
        aOuter2inner.push_back( j );
        aInner2outer.push_back( i );
        ++j;
        break;

      case '3':
        // contour turning toward its outside.
        aInnerChain.chain += aOuterChain.chain[ i ];
        aInnerChain.chain += aOuterChain.chain[ ( i + 1 ) % nb ];
        aOuter2inner.push_back( j );
        aInner2outer.push_back( i );
        aInner2outer.push_back( i + 1 );
        j += 2;
        break;
    }

    // Advances along contour and check if it is a closed contour.
    it.next();
    if ( ( i == nb - 1 )
        && ( *it_begin != *it ) )
      // freeman chain is *not* a closed loop.
    {
      aInnerChain.chain += aOuterChain.chain[ i ];
      aOuter2inner.push_back( j );
      aInner2outer.push_back( i );
      ++i;
      ++j;
      break;
    }
  }
}




template <typename TInteger>
bool DGtal::FreemanChain<TInteger>::cleanOuterSpikes( 
    FreemanChain & aCleanC,
    std::vector<unsigned int> & aC2clean,
    std::vector<unsigned int> & aClean2c,
    const FreemanChain & c,
    bool ccw )
{
  unsigned int nb = (unsigned int)c.chain.size();
  if ( nb == 0 )
  {
    cerr << "[DGtal::FreemanChain::cleanOuterSpikes]"
      << " cleanOuterSpikes: Empty input chain"
      << endl;
    return false;
  }

  ModuloComputer< DGtal::int32_t > mc( nb );
  ModuloComputer< DGtal::int32_t >::UnsignedInteger i = 0;
  ModuloComputer< DGtal::int32_t >::UnsignedInteger j = 0;
  vector<unsigned int> c2cleanTMP;
  aCleanC.chain.reserve( nb );
  aCleanC.chain = "";
  aC2clean.clear();
  aClean2c.clear();
  aC2clean.reserve( nb );
  aClean2c.reserve( nb );
  c2cleanTMP.reserve( nb );
  ConstIterator it = c.begin();
  ConstIterator itn = c.begin();
  itn.nextInLoop();
  // Find a consistent starting point.
  unsigned int n;
  unsigned int size_spike = 0;
  for ( n = 0; n < nb; ++n )
  {
    size_spike = 0;
    while ( movement( it.getCode(), itn.getCode(), ccw ) == '0' )
    {
      it.previousInLoop();
      itn.nextInLoop();
      mc.increment( i );
      size_spike += 2;
      if ( size_spike >= nb )
      {
        cerr << "[DGtal::FreemanChain::cleanOuterSpikes]"
          << " Spike is longer than contour !"
          << " size_spike=" << size_spike
          << " nb=" << nb
          << endl;
        return false;
      }
    }
    mc.increment( i );
    it = itn;
    itn.nextInLoop();
    if ( size_spike > 0 )
      break;
  }
  unsigned int start_idx = it.getPosition();
  i = start_idx;
  // JOL : 2009/07/07, added starting coordinates
  // XP  : 2011/09/06, added ending coordinates
  Point P = *it;
  aCleanC.x0 = P.at(0);
  aCleanC.y0 = P.at(1);
  aCleanC.xn = P.at(0);
  aCleanC.yn = P.at(1);

  // cerr << "Starting point is " << i << endl;
  ASSERT( ( n < nb ) || ( i == 0 ) );
  if (  n == nb  )
  { // do nothing
    aCleanC.chain = c.chain;
    for ( unsigned int ni = 0; ni < nb; ++ni )
    {
      aC2clean.push_back( ni );
      aClean2c.push_back( ni );
    }
    if ( size_spike != 0 )
      cerr << "[DGtal::FreemanChain::cleanOuterSpikes]"
        << "No starting point found (only spikes !)" << endl;

    return size_spike == 0;
  }
  // Loops over all letters.
  ConstIterator it_begin = it;
  deque<unsigned int> clean_code;
  deque<unsigned int> clean_idx;
  vector<unsigned int> begin_outer_spike;
  vector<unsigned int> end_outer_spike;
  // i follows iterator it.
  do
  {
    clean_code.push_back( it.getCode() );
    clean_idx.push_back( i );
    itn = it;
    it.nextInLoop();
    mc.increment( i );
    // cerr << "- i=" << i << " (" << clean_code.back()
    //     << it.getCode() << ") ";
    size_spike = 0;
    unsigned int last_spike_idx = end_outer_spike.empty() ?
      start_idx :
      end_outer_spike.back();
    j = i;
    while ( ( ! clean_code.empty() )
        && ( j != last_spike_idx )
        && ( movement( clean_code.back(), it.getCode(), ccw ) == '0' )
        && ( it != it_begin ) )
    {
      clean_code.pop_back();
      clean_idx.pop_back();
      mc.increment( i );
      mc.decrement( j );
      it.nextInLoop();
      itn.previousInLoop();
      size_spike += 2;
    }
    // cerr << "i=" << i << " size_spike=" << size_spike
    //     << " last_spike_idx=" << last_spike_idx
    //     << endl;
    if ( size_spike != 0 )
    {
      // There is a spike. Is it an outer one ?
      unsigned int previous_code = itn.getCode();
      unsigned int previous_idx = itn.getPosition();
      // JOL : do not
      // consider any more "cleaned contour" since we cannot go
      // further than the last spike.
      // unsigned int previous_code =
      //   clean_code.empty() ? itn.getCode() : clean_code.back();
      // unsigned int previous_idx =
      //   clean_code.empty() ? itn.getPosition() : clean_idx.back();
      itn = it;
      itn.previousInLoop();
      unsigned int move1 = movement( previous_code, addToCode ( itn.getCode() , 2 ), ccw );
      unsigned int move2 = movement( itn.getCode(), it.getCode() , ccw );
      bool return_spike = ( move1 == '0' ) || ( move2 == '0' );
      bool outer_spike  = ( move1 == '3' ) || ( move2 == '3' );
      //          if ( return_spike )
      //            cerr << "[DGtal::FreemanChain::cleanOuterSpikes] return spike."
      //                 << endl;
      //          if ( ! ( ( outer_spike && ( move1 != 1 ) && ( move2 != 1 ) )
      //                   || ( ! outer_spike && ( move1 != 3 ) && ( move2 != 3 ) ) ) )
      //            cerr << "[DGtal::FreemanChain::cleanOuterSpikes] "
      //                 << "Weird spike. Invalid contour (expected 3 3) ?"
      //                 << " move1=" << move1
      //                 << " move2=" << move2
      //                 << " ccw=" << ccw
      //                 << " start_idx=" << start_idx
      //                 << " size_spike=" << size_spike
      //                 << " it=" << it.getPosition()
      //                 << " itp=" << previous_idx
      //                 << endl
      //                 << c.chain << endl;
      // Process waiting steps.
      if ( outer_spike || return_spike )
      {
        begin_outer_spike.push_back( mc.next( previous_idx ) );
        end_outer_spike.push_back( i );
        // cout << " outer spike [" << begin_outer_spike.back()
        //         << "," << end_outer_spike.back() << "[  " << endl;
      }
    }
  }
  while ( it != it_begin );

  // Once outer spikes are known, we can create the new contour.
  aC2clean.resize( nb );
  i = start_idx % nb;
  j = 0;
  unsigned int nb_spikes = (unsigned int)begin_outer_spike.size();
  unsigned int k = 0;
  n = 0;
  while ( n < nb )
  {
    if ( ( k == nb_spikes ) || ( i != begin_outer_spike[ k ] ) )
    {
      aCleanC.chain.push_back( c.chain[ i ] );
      aC2clean[ i ] = j;
      aClean2c.push_back( i );
      mc.increment( i );
      ++j;
      ++n;
    }
    else
    {
      while ( i != end_outer_spike[ k ] )
      {
        aC2clean[ i ] = j;
        mc.increment( i );
        ++n;
      }
      ++k;
    }
  }
  for ( unsigned int ii = 0; ii < nb; ++ii )
    if ( aC2clean[ ii ] >= aCleanC.chain.size() )
    { 
      if ( aC2clean[ ii ] == aCleanC.chain.size() )
        aC2clean[ ii ] = 0;
      else
      {
        cerr << "[DGtal::FreemanChain::cleanOuterSpikes]"
          << "Bad correspondence for aC2clean[" << ii << "]"
          << " = " << aC2clean[ ii ] << " >= " << aCleanC.chain.size()
          << endl;
        aC2clean[ ii ] = aC2clean[ ii ] % aCleanC.chain.size();
      }
    }

  for ( unsigned int jj = 0; j < aCleanC.chain.size(); ++jj )
    if ( aClean2c[ jj ] >= nb )
    {
      cerr << "[DGtal::FreemanChain::cleanOuterSpikes]"
        << "Bad correspondence for aClean2c[" << jj << "]"
        << " = " << aClean2c[ jj ] << " >= " << nb
        << endl;
      aClean2c[ jj ] = aClean2c[ jj ] % nb;
    }
  return true;
}









//     /**
//      * Given a Freeman chain [c] coding a 4-connected pixel loop, computes
//      * its subsampling by the transformation:
//      * X = ( x - x0 ) div h,
//      * Y = ( y - y0 ) div v.
//      *
//      * @param aSubc (output) the subsampled Freeman chain code (may
//      * contain spikes)
//      *
//      * @param aC2subc (output) the mapping associating an element to
//      * its subsampled element.
//      *
//      * @param aSubc2c (output) the inverse mapping associating a
//      * subsampled element to its element. More precisely, subc2c[ j ]
//      * is the last pointel to be in j.
//      *
//      * @param c the input chain code.
//      *
//      * @param h the subsampling along x
//      * @param v the subsampling along y
//      * @param x0 the x-origin of the frame (X,Y) in (x,y)
//      * @param y0 the y-origin of the frame (X,Y) in (x,y)
//      *
//      * @return 'false' if initial contour was empty or if [subc] is empty,
//      * 'true' otherwise.
//      */
//     static bool subsample( FreemanChain & aSubc,
//          std::vector<unsigned int> & aC2subc,
//          std::vector<unsigned int> & aSubc2c,
//          const FreemanChain & c,
//          unsigned int h, unsigned int v,
//          int x0, int y0 ){
//       if ( ( h == 0 ) || ( v == 0 ) ) return false;
//       FreemanChain<TInteger>::ConstIterator it = c.begin();
//   unsigned int j = 0;
//   unsigned int nb = c.chain.size();
//   if ( nb == 0 ) return false;
//   Point fxy( it.get() );
//   Point fXY;
//   fXY.at(0)= ( fxy.at(0) - x0 ) / h;
//   fXY.at(1)= ( fxy.at(1) - y0 ) / v;
//   aSubc.x0 = fXY.at(0);
//   aSubc.y0 = fXY.at(1);
//   aC2subc.clear();
//   aC2subc.reserve( nb );
//   aSubc2c.clear();
//   aSubc2c.reserve( nb );
//   for ( unsigned int i = 0; i < nb; ++i )
//     {
//       aC2subc.push_back( j );
//       it.nextInLoop();
//       Point nxy( it.get() );
//       Point nXY;
//       nXY.at(0)= ( nxy.at(0) - x0 ) / h;
//       nXY.at(1)= ( nxy.at(1) - y0 ) / v;
//       if ( nXY != fXY )
//   {
//     aSubc2c.push_back( i );
//     char code;
//     if ( nXY.at(0) > fXY.at(0) )       code = '0';
//     else if ( nXY.at(0) < fXY.at(0) )  code = '2';
//     else if ( nXY.at(1) > fXY.at(1) )  code = '1';
//     else                           code = '3';
//     aSubc.chain += code;
//     ++j;
//     fXY = nXY;
//   }
//     }
// //   aC2subc.push_back( j );
// //   it.nextInLoop();
// //   for ( unsigned int i = 1; i <= nb; ++i )
// //     {
// //       // JOL
// //       //aC2subc.push_back( j );
// //       Vector2i nxy( it.get() );
// //       Vector2i nXY( ( nxy.x() - x0 ) / h, ( nxy.y() - y0 ) / v );
// //       if ( nXY != fXY )
// //   {
// //     char code;
// //     if ( nXY.x() > fXY.x() )       code = '0';
// //     else if ( nXY.x() < fXY.x() )  code = '2';
// //     else if ( nXY.y() > fXY.y() )  code = '1';
// //     else                           code = '3';
// //     aSubc.chain += code;
// //     aSubc2c.push_back( i - 1 );
// //     ++j;
// //     fXY = nXY;
// //   }
// //       if ( i != nb ) aC2subc.push_back( j );
// //       it.nextInLoop();
// //     }
// //   // TODO : enhance this.
//   unsigned int nbsub =  aSubc.chain.size();
//   // Last correspondence may be in fact to 0 instead of nbsub.
//   if ( nbsub != 0 )
//     for ( unsigned int i = 0; i < nb; ++i )
//       if ( aC2subc[ i ] >= nbsub ) aC2subc[ i ] -= nbsub;
//   ASSERT( aC2subc.size() == nb );
//   return nbsub != 0;
//     }





// Drawing services                 

template <typename TInteger>
inline
std::string
DGtal::FreemanChain<TInteger>::className() const
{
  return "FreemanChain";
}

template <typename TInteger>
inline
void 
DGtal::FreemanChain<TInteger>::computeLastPoint() 
{
  ConstIterator it = this->begin();
  for (unsigned int i=0; i < chain.size(); ++i)
    it++;
  xn = (*it).at(0);
  yn = (*it).at(1);
}






// Implementation of inline functions and external operators                 //

template <typename TInteger>
inline
  std::ostream&
DGtal::operator<< ( std::ostream & out,
    const FreemanChain<TInteger> & aObject )
{
  aObject.selfDisplay ( out );
  return out;
}

//                                                                           //