function varargout = EMP8_erode(s)
% EMP8_erode       Erode small components and/or border region on label image 
%                  EMP8: Energy Minimization using Potts model and 8 neighbors
%
% This function may use:
%   EMP8_erode_mex(..)    to do the actual work (if s.numPix2 > 0)
%   EMP8_erode255_mex(..) to erode along borders (if s.borderWidth > 0)
%   EMP8_cost_mex(..)     to find cost for a labeling
%
% To do only erode along borders set s.numPix2 = 0  (upper limit for segments to
% erode in EMP8_erode_mex) and s.borderWidth > 0.
%
% Use of function:
% ----------------
%   L = EMP8_erode(s)
%   [L,c] = EMP8_erode(s)
% 
%     L   is the resulting new labeling after erode
%     c   is a struct with cost results in each field
%     s   is a struct as returned by EMP8_argin, see help for that file
% 
% Example:  (see also EMP8_ex1)
% -------------------------------------------------------
%  clear all; tim = 4;  sigma = 5;  
%  s = EMP8_argin(tim, sigma, tim, 'v',2);   % use Ltrue as well
%  [L,c] = EMP8_erode(s);
%  s.borderWidth = 5; [L,c] = EMP8_erode(s);
%  figure(1); clf;  % and plot results in a figure
%  figTitle = ['Labeling of test image ',int2str(tim),' after using EMP8-erode.'];
%  LIS_plotcol(reshape(L,s.I,s.J), figTitle, c.errorRate, c.erodeNofComp, c.Etot, c.Edata, c.E4+c.E8);

%----------------------------------------------------------------------
% Copyright (c) 2013.  Karl Skretting.  All rights reserved.
% University of Stavanger, Signal Processing Group
% Mail:  karl.skretting@uis.no   Homepage:  http://www.ux.uis.no/~karlsk/
% 
% HISTORY:  dd.mm.yyyy
% Ver. 1.0  20.06.2013  Make function based on LIS_erode
% Ver. 1.1  11.12.2013  Some changes, and use s (struct) as input
%----------------------------------------------------------------------

% some identities that may be helpful
% i = rem(pix-1,I)+1;  j = floor((pix-1)/I)+1;  pix = (j-1)*I+i;

mfile = 'EMP8_erode';

if ( ~exist('EMP8_erode_mex','file') || ...
     ~exist('EMP8_cost_mex','file') || ...     
     ~exist('EMP8_erode255_mex','file') )
    c = [mfile,': one or more needed mex-files are not availble.'];
    disp(c);
    if (nargout >= 1); varargout{1} = []; end;
    if (nargout >= 2); varargout{2} = c; end;
    return;
end

%% do the calulations
param = [s.M,s.I,s.J,-1, s.bF4,s.bF8,s.swUpq,0, s.lambda, s.dR, s.minsumRi, s.maxsumRi, ...
         s.numPix1, s.numPix2, s.shapeFactor, -1];   % 16 parameters

if ((s.bF4 < 0) || (s.bF8 < 0))
    % this should normally not happen
    disp([mfile,': Values are not assigned to bF4 or bF8 yet, do that.']);
    disp([' bF4=',int2str(s.bF4),', bF8=',int2str(s.bF8),', lambda=',num2str(s.lambda), ...
         ', bF8percent=',int2str(s.bF8percent)  ]);
    if (s.bF4 < 0) 
        param(6) = s.bF8percent; 
    else   % (s.bF8 < 0)  % and (bF4 >= 0)
        param(6) = floor(s.bF4*(s.bF8percent/100)+0.5);
    end
    if (s.swUpq == 1)
        Ci = EMP8_cost_mex(s.Ri, s.Li, param);
    elseif (s.swUpq == 2)
        Ci = EMP8_cost_mex(s.Ri, s.Li, s.Upq, param);
    elseif (s.swUpq == 3)
        Ci = EMP8_cost_mex(s.Ri, s.Li, s.Uim, s.Utab, param);
    end
    if (s.bF4 < 0)
        s.bF4 = Ci(5);
        s.bF8 = Ci(6);
    end
    if (s.dR < 0)
        s.dR = Ci(14); 
        s.minsumRi = Ci(15); 
        s.maxsumRi = Ci(16);
    end
    param = [s.M,s.I,s.J,-1, s.bF4,s.bF8,s.swUpq,0, s.lambda, s.dR, s.minsumRi, s.maxsumRi, ...
             s.numPix1, s.numPix2, s.shapeFactor, -1];   % 16 parameters
    if (s.verbose >= 1)
        disp([mfile,': The values assigned to bF4 or bF8 are']);
        disp([' bF4=',int2str(s.bF4),', bF8=',int2str(s.bF8),', lambda=',num2str(s.lambda), ...
            ', bF8percent=',int2str(s.bF8percent)  ]);
    end
elseif (s.verbose >= 1)
    disp(' ');
    disp([mfile,': some values are bF4 = ',int2str(s.bF4),', bF8 = ',int2str(s.bF8),...
          ', lambda = ',num2str(s.lambda),', bF8percent = ',int2str(s.bF8percent), ...
          ', swUpq = ',int2str(s.swUpq)]);
end

% do erode in mex-file, EMP8_erode_mex
tic;
if s.numPix2 > 0
    if (s.swUpq == 1)
        [L, c] = EMP8_erode_mex(s.Ri, s.Li, param);
    elseif (s.swUpq == 2)
        [L, c] = EMP8_erode_mex(s.Ri, s.Li, s.Upq, param);
    elseif (s.swUpq == 3)
        [L, c] = EMP8_erode_mex(s.Ri, s.Li, s.Uim, s.Utab, param);
    end
    erodeUpdated = c(17);
    erodeNofComp = c(18);
    erodeCompEdges = c(19);
    if (s.verbose >= 1)  % very verbose
        if (numel(s.Lt) == (s.I*s.J))
            er = nnz(L ~= s.Lt)/(s.I*s.J);
        else
            er = -1;
        end
        erodeCost = sprintf(s.formatCost, 'Lerode', c([10,13,11,12,9]), 100*er);
        if (er < 0) % it is nicer to print nothing rather the negative error rate
            erodeCost = [erodeCost(1:72),'.'];
        end
    end
else
    disp([mfile,': numPix2 = ',num2str(s.numPix2),' -> do not call EMP8_erode_mex(..)']);
    L = s.Li;
end
erodeTime = toc;
if ((s.verbose >= 2) && exist('erodeUpdated','var'))
    disp([' Part EMP8_erode_mex(..) updated ',int2str(erodeUpdated),...
          ' pixels using ',num2str(erodeTime),' seconds.']);
    if (s.verbose >= 2); disp(erodeCost); end;
end

% do erode in mex-file, EMP8_erode255_mex
if s.borderWidth > 0
    tic;
    I = s.I;  J = s.J;
    L = reshape(L,I,J);
    mask = [false(I,1), (L(:,1:(J-1)) ~= L(:,2:J))] | ...
        [(L(:,1:(J-1)) ~= L(:,2:J)), false(I,1)] | ...
        [false(1,J); (L(1:(I-1),:) ~= L(2:I,:))] | ...
        [(L(1:(I-1),:) ~= L(2:I,:)); false(1,J)] ;
    mask = imdilate(mask, strel('disk',s.borderWidth));
    if nnz(mask) == (I*J)
        disp([mfile,': borderWidth = ',num2str(s.borderWidth),...
            ' -> all pixels are on border. Can not erode.']);
        L = L(:);  % set this as column vector again.
    else
        L(mask) = 255;
        if (s.swUpq == 1)
            [L, c] = EMP8_erode255_mex(s.Ri, L(:), param);
        elseif (s.swUpq == 2)
            [L, c] = EMP8_erode255_mex(s.Ri, L(:), s.Upq, param);
        elseif (s.swUpq == 3)
            [L, c] = EMP8_erode255_mex(s.Ri, L(:), s.Uim, s.Utab, param);
        end
        erode255Processed = c(17);
        if (s.verbose >= 1)  % verbose
            if (numel(s.Lt) == (s.I*s.J))
                er = nnz(L ~= s.Lt)/(s.I*s.J);
            else
                er = -1;
            end
            erode255Cost = sprintf(s.formatCost, 'Le255', c([10,13,11,12,9]), 100*er);
            if (er < 0) % it is nicer to print nothing rather the negative error rate
                erode255Cost = [erode255Cost(1:72),'.'];
            end
        end
    end
    erode255Time = toc;
    if ((s.verbose >= 2) && exist('erode255Processed','var'))
        disp([' Part EMP8_erode255_mex(..) processed ',num2str(erode255Processed),...
              ' pixels (borderWidth = ',num2str(s.borderWidth),...
              ') using ',num2str(erode255Time),' seconds.']);
        if (s.verbose >= 2); disp(erode255Cost); end;
    end
else
    erode255Time = 0;
end

% if useICM
%     L = reshape(L,I,J);
%     mask = [false(I,1), (L(:,1:(J-1)) ~= L(:,2:J))] | ...
%         [(L(:,1:(J-1)) ~= L(:,2:J)), false(I,1)] | ...
%         [false(1,J); (L(1:(I-1),:) ~= L(2:I,:))] | ...
%         [(L(1:(I-1),:) ~= L(2:I,:)); false(1,J)] ;
%     param(4) = nnz(mask);
%     if (swUpq == 1)
%         [L, C] = EMP8_icm_mex(Ri, L(:), int32(find(mask)-1), param);
%     elseif (swUpq == 2)
%         [L, C] = EMP8_icm_mex(Ri, L(:), int32(find(mask)-1), Upq, param);
%     elseif (swUpq == 3)
%         [L, C] = EMP8_icm_mex(Ri, L(:), int32(find(mask)-1), Uim, Utab, param);
%     end
%     if (verbose >= 2)
%         disp([' ICM finally relabeled ',int2str(C(17)),' pixels.']);
%     end
% end

%% return results
if (nargout >= 1)
    varargout{1} = L;  
end
if (nargout >= 2)
    if (numel(s.Lt) == (s.I*s.J))
        er = nnz(L ~= s.Lt)/(s.I*s.J);
    else
        er = -1;
    end
    updated = nnz(L ~= s.Li);
    if exist('c','var')   % and isnumeric(c)
        % den normale situasjonen
        c(17:20) = [updated; erodeTime; erode255Time; er];
        c = costTilStruct(c);
    else
        disp([mfile,': probably nothing was done.']);
        if isfield(s,'Ci')
            c = s.Ci;
        else
            if (s.swUpq == 1)
                c = EMP8_cost_mex(s.Ri, L, param);
            elseif (s.swUpq == 2)
                c = EMP8_cost_mex(s.Ri, L, s.Upq, param);
            elseif (s.swUpq == 3)
                c = EMP8_cost_mex(s.Ri, L, s.Uim, s.Utab, param);
            end
            c(17:20) = [updated; erodeTime; erode255Time; er];
            c = costTilStruct(c);
        end
    end
    if exist('erodeUpdated','var'); c.erodeUpdated = erodeUpdated; end;
    if exist('erodeNofComp','var'); c.erodeNofComp = erodeNofComp; end;
    if exist('erodeCompEdges','var'); c.erodeCompEdges = erodeCompEdges; end;
    if exist('erode255Processed','var'); c.erode255Processed = erode255Processed; end;
    if exist('erodeCost','var'); c.erodeCost = erodeCost; end;
    if exist('erode255Cost','var'); c.erode255Cost = erode255Cost; end;
    if (s.verbose >= 1)
        ttt = 0;
        if exist('erodeTime','var'); ttt = ttt + erodeTime; end;
        if exist('erode255Time','var'); ttt = ttt + erode255Time; end;
        disp([mfile,': done. Updated ',int2str(updated),' pixels using ',...
              num2str(ttt),' seconds.']);
    end

    varargout{2} = c;
end

return



function s = costTilStruct(c)
if (numel(c)<20)
    c = [c(:);-ones(20-numel(c),1)];
end
s = struct('Ctot',c(1),'Cdata',c(2),'C4',c(3),'C8',c(4),...
           'bF4',c(5),'bF8',c(6),'n4',c(7),'n8',c(8),...
           'Etot',c(9),'Edata',c(10),'E4',c(11),'E8',c(12),...
           'lambda',c(13),'dR',c(14),'minsumRi',c(15),'maxsumRi',c(16),...
           'updated',c(17),'time',c(18)+c(19),...
           'erodeTime',c(18),'erode255Time',c(19),'errorRate',c(20));
return;
       
