function Ires = oprec(input, se)

% OPREC Perform opening by reconstruction on bulk images
% Syntax: Ires = oprec(input, se, n);
% Input:
% - input : input image
% - se : strel object
% Output:
% - Ires : matrix of the image after reconstruction
%
% Matlab v6.5 compatible (at least op. by rec. on binary images)

error(nargchk(2, 2, nargin));
Ires = opclrec(+uint8(input), se, 2, 'o');


%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% OPCLREC : opening or closing by reconstruction
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
function Ires = opclrec(I, struct_el, color, opcl)

str = struct_el(1);
if length(struct_el) == 2
    str_el.se = struct_el(2);
else
    str_el.se = 0;
end % IF
switch opcl
case {'c', 'C', '0'}
    if color == 2
        I = +uint8(imcomplement(logical(double(I))));
        str_el.marker = imopen(I, str);
        Ires = +uint8(imcomplement(logical(double(conddil(I, str_el, 0, color)))));
    else
        I = imcomplement(I);
        str_el.marker = imopen(I, str);
        Ires = imcomplement(conddil(I, str_el, 0, color));
    end % IF
otherwise
    str_el.marker = double(imopen(I, str));
    Ires = double(conddil(I, str_el, 0, color));
end % SWITCH


%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% CONDDIL : perform a conditional dilation on a binary image (image reconstruction)
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
function [Ires, comment] = conddil(I, str_el, n, color)

se = str_el.se;
Ires = str_el.marker;
if color == 2
    try
        se == 0; % Use no structuring element, but 8-neighbour connectivity
        Ires = double(imreconstruct(+uint8(Ires), I));
    catch
        if n == 0
            iseq = 0;
            while iseq == 0
                Ires0 = Ires;
                Ires = double(imdilate(Ires, se) & I);
                iseq = isequal(Ires, Ires0);
            end % WHILE
        else
            for i=1:n
                Ires = double(imdilate(Ires, se) & I);
            end % FOR
        end % IF
    end % TRY
else % Grayscale images
    try
        se == 0; % Use no structuring element, but 8-neighbour connectivity
        Ires = imreconstruct(+uint8(Ires), I);
    catch
        if n == 0
            iseq = 0;
            while iseq == 0
                Ires0 = Ires;
                Ires = min(imdilate(Ires, se), I);
                iseq = isequal(Ires, Ires0);
            end % WHILE
        else
            for i=1:n
                Ires = min(imdilate(Ires, se), I);
            end % FOR
        end % IF
    end % TRY
end % IF