essay on programming languages, computer science, information techonlogies and all.

Wednesday, July 24, 2013

Properties of chain code i.e. perimeter, area and center of mass

The chain code is compact and efficient to store and transfer. Though it is not that straight forward to calculate some properties of shape .i.e area, perimeter and center of mass.

Perimeter is relatively easy and from Bernd Jähne book ( Digital Image Processing 5th: 508-511 ), it can be calculated as below code.
  double CalculatePerimeter( const ChainCode& chains )
  {
    double p = 0;
    double sqrt2 = sqrt( 2.0 );

    for( size_t i=0; i < chains.size(); ++i )
    {
      if( chains[i]%2 == 0 ) 
      {
        p += 1.0;
      }
      else
      {
        p += sqrt2;
      }
    }

    return p;
  }
Bernd Jähne also showed how to get the area using a chain code but it doesn't work well when the chain code is 1 pixel width or height e.g. ( 0 0 4 4 ) will return 2 instead of 3. Refer below image how the area end up 1 pixel short.


Cris Luengo ( Cris's Image Analysis Blog : More chain code measures ) shows a matlab code that works well with 1 pixel wide chain code. In above case of 3 pixels, it adds two more addition when the direction changes which add up 1 more pixel. Refer below image.
Here is the C/C++ code snippet that is converted from the above Cris Luengo's matlab code.

// making a circular chain code by inserting the last code to the front 
ChainCode cc;
cc.push_back( chains[ chains.size()-1 ] );
for( size_t i=0; i < chains.size(); ++i )
{
  cc.push_back( chains[i] );
}
m_Area = CalculateArea( cc );


int CalculateArea( const ChainCode& cc ) // a circular chain
{
  int (&M)[8][8]( TableAreaIncrement );
  int B = 10, A = 0;

  for( size_t i=1; i < cc.size(); ++i )
  {
    uint8_t c_i = cc[i];
    uint8_t c_im1 = cc[i-1];

    switch( M[c_im1][c_i] )
    {
    case 'X' :
      throw Exception( StringFormat::As( _T("Invalid chain code found as %d->%d"), 
        c_im1, c_i ));
    case 'A': A += -B+1;  break;
    case 'B': A += -B;    break;
    case 'C': A += B;    break;
    case 'D': A += B-1;    break;
    }

    switch( c_i )
    {
    case 0: A+=B;      break;
    case 1: A+=(++B);    break;
    case 2: B++;      break;
    case 3: A+=-(++B)+1;  break;
    case 4: A+=-B+1;    break;
    case 5: A+=-(--B)+1;  break;
    case 6: B--;      break;
    case 7: A+=(--B);    break;
    }
  }

  return A;
}
  
int ChainCodeProperties::TableAreaIncrement[8][8] = 
{
  {'0','1','B','X','A','A','6','7'},
  {'0','1','B','B','X','A','6','7'},
  {'C','C','2','3','4','X','C','C'},
  {'C','C','2','3','4','5','X','C'},
  {'C','C','2','3','4','5','D','X'},
  {'X','C','2','3','4','5','D','D'},
  {'0','X','A','A','A','A','6','7'},
  {'0','1','X','A','A','A','6','7'},
};

The center of mass can be calculated by summing up x and y coordinates. But just simply adding up does not cover 1 pixel wide chain code as below image. The (1,1) pixel added two times instead one time. The compensation can be done by adding more pixels when a direction changes in opposite way as the idea of Cris Luengo's area calculation.
Below is the C/C++ code snippet that calculate the centroid.
void CalculateCentroid( const ChainCode& cc, double* cx, double *cy )
{
  int (&O)[8][2]( TableOffsetIncrement );
  int vpc = 0;      // visitied pixel count
  int sx = 0, sy = 0;    // sum of x and y
  int px = 0, py = 0;
  
  for( size_t i=1; i < cc.size(); ++i )
  {
    uint8_t c_i = cc[i];
    uint8_t c_im1 = cc[i-1];

    if( c_i != c_im1 && (c_i & 0x3 ) == ( c_im1 & 0x3 ) )
    {
      // when direction changes in opposite direction 
      // e.g. 0 -> 4, 4 -> 0, 1 -> 5, ...
      vpc++; sx += px;  sy+=py;
    }

    px += O[ c_i ][0];
    py += O[ c_i ][1];

    vpc++; sx += px;  sy+=py;
  }

  *cx = (double)sx / (double)vpc;
  *cy = (double)sy / (double)vpc;
}

int ChainCodeProperties::TableOffsetIncrement[8][2] =
{
  { 1,  0  },
  { 1,  -1  },
  { 0,  -1  },
  { -1,   -1  },
  { -1,  0  },
  { -1,  1 },
  { 0,  1 },
  { 1,  1 },
};  
It works well most of the case but it fails if there is a pixel which is visited more than 2 times as below image. It has (1,1) pixel counted 3 times but it should be counted 2 times.


No comments: