% CCAREA   Compute the area enclosed by a chain code
%
% CCAREA(CC) computes the area of the object represented by CC.
%
% For example:
%    a = label(~threshold(readim('cermet')))==33;
%    cc = dip_imagechaincode(img,2,33);
%    ccarea(cc.chain)
%    sum(a)
%
% The method is based on the method by P. Zamperoni (Signal
% Processing 3(3):267-271, 1981). However, because the method was
% not well explained in the original article, I had trouble
% reproducing it with the increments given, and had to derive a
% slightly different set of rules.

% (c)2010, Cris Luengo, Centre for Image Analysis, Uppsala, Sweden.

function A = ccarea(chaincode)

M = ['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'];

if isempty(chaincode)
   A = 1;
   return; % if the object has 1 pixel, the chain code is empty.
end

B=1000; A=0;
chaincode = [chaincode(end),chaincode];
for ii=2:length(chaincode)
   switch M(chaincode(ii-1)+1,chaincode(ii)+1)
      case 'X'
         error('This is not a valid chain code!?')
      case 'A'
         A = A-B+1;
      case 'B'
         A = A-B;
      case 'C'
         A = A+B;
      case 'D'
         A = A+B-1;
   end
   switch chaincode(ii)
      case 0
         A = A+B;
      case 1
         B = B+1;
         A = A+B;
      case 2
         B = B+1;
      case 3
         B = B+1;
         A = A-B+1;
      case 4
         A = A-B+1;
      case 5
         B = B-1;
         A = A-B+1;
      case 6
         B = B-1;
      case 7
         B = B-1;
         A = A+B;
   end
end   
