function [edgeOut, g1Out, g2Out] = hdledge(frameIn, method, threshold)

%  HDLEDGE Find edges of objects in an image
%  [edgeOut,G1,G2] = hdledge(frmInput) finds edges in the input image
%  by using the Sobel method with a threshold of 20. The input image
%  data type can be any numeric or logical type. The function returns
%  a binary image of edges detected and orthogonal gradient images
%  in the full-precision data type.
%
%  [edgeOut,G1,G2] = hdledge(frmInput,method,threshold)
%  finds edges in the input image by using the method and gradient
%  threshold that you specify. Method must be "Sobel",
%  "Roberts", or "Prewitt".
%
%   Example:
%   Apply edge detector to an image.
%   frmInput = imread('rice.png');
%   [edgeOut,G1,G2] = hdledge(frmInput,'Sobel',20);
%   figure, imshow(frmInput), figure, imshow(edgeOut);
%
%
%   See also visionhdl.EdgeDetector, edge
%
%   Copyright 2023 The MathWorks, Inc.

%#codegen

%--- Input Validation -----------------------------------------------------

% Validate method and boundary method by enforcing a string set via
% musbeValidMethod function

arguments
    frameIn
    method(1,:) char {mustBeValidMethod} = 'Sobel';
    threshold = 20;
end

funcName = 'hdledge';
paddingMethod = 'Replicate';


validateattributes(frameIn, {'half', 'single', 'double', 'int8', 'uint8', 'int16', 'uint16', 'int32', 'uint32', 'embedded.fi', 'logical'}, ... % Class
    {'2d', 'nonsparse', 'real', 'nonnan', 'finite'}, ... % Attributes
    funcName, 'frameIn', 1);  %funcName, argName, argIndx

validateattributes(threshold, {'half', 'single', 'double', 'int8', 'uint8', 'int16', 'uint16', 'int32', 'uint32', 'embedded.fi'}, ... % Class
    {'scalar','nonnan','finite','real'}, ... % Attributes
    funcName, 'threshold', 2);  %funcName, argName, argIndx

% No doubles for HDL
if coder.target('hdl')
    coder.internal.assert(~isa(frameIn, 'double'),'visionhdl:hdledge:DoubleType', funcName);
end

% Fixed point checks
if(isfi(frameIn) && (frameIn.WordLength > 32))
    coder.internal.error('visionhdl:hdledge:InputWordLength', funcName, '32');
end

% --- Setup --------------------------------------------------

% Declare compile time constants for the kernel core based on the input method.
switch lower(method)
    case 'sobel'
        kx = 3;
        ky = 3;
        scale = 4;
        extraBitsWord = 3;
        extraBitsFrac = 3;

    case 'prewitt'
        kx = 3;
        ky = 3;
        scale = 4;
        extraBitsWord = 19;
        extraBitsFrac = 18;

    case 'roberts'
        kx = 2;
        ky = 2;
        scale = 6;
        extraBitsWord = 1;
        extraBitsFrac = 1;

    otherwise
        kx = 3;
        ky = 3;
        scale = 4;
        extraBitsWord = 3;
        extraBitsFrac = 3;


end

% Input Pre-processing
[frameInCasted, thresholdSqCasted] = parse_inputs(frameIn, threshold, extraBitsWord, extraBitsFrac, funcName);

% --- Filter --------------------------------------------------
[edgeImage,g1, g2] = EdgeDetectorCore(frameInCasted, method, kx, ky, paddingMethod);

% --- Threshold Computation for logical input -----------------------------------------

if(islogical(frameIn))
    edgeImageSum = accumulator(edgeImage);
    frameHeight = size(edgeImage,1);
    frameWidth = size(edgeImage,2);
    cutoff = single(scale * edgeImageSum) * single(1/(frameHeight*frameWidth));
    thresholdCasted = fi(sqrt(cutoff), 0, 16, 16, 'RoundingMethod','Floor','OverflowAction', 'Wrap');
    thresholdSqCasted = fi(thresholdCasted*thresholdCasted, issigned(thresholdSqCasted), thresholdSqCasted.WordLength, thresholdSqCasted.FractionLength);
    if(coder.target('MATLAB'))
        disp(['Threshold used for the image:',num2str(thresholdCasted)]);
    end

end

%--- Output Post processing----------------------------------------

if(islogical(frameIn))
    g1Out = fi(g1, 1, 17, 16);
    g2Out = fi(g2, 1, 17, 16);
else
    g1Out = g1;
    g2Out = g2;
end

% Threshold the edge output based on the cutoff value
edgeOut = edgeImage > thresholdSqCasted;

end

%--- Pre-processing Core Function ------------------------------------------
function [frameInCasted, thresholdSq] = parse_inputs(frameIn, threshold,extraBitsWord, extraBitsFrac, funcName)
% Casts input frame and threshold to the corresponding fixed data type.

inputDataType = class(frameIn);

if(isfloat(frameIn) || isfi(frameIn))
    frameInCasted = frameIn;
else
    switch inputDataType

        case 'int8'
            signed = 1;
            wordLength = 8;
            fractionLength = 0;

        case 'uint8'
            signed = 0;
            wordLength = 8;
            fractionLength = 0;

        case 'int16'
            signed = 1;
            wordLength = 16;
            fractionLength = 0;

        case 'uint16'
            signed = 0;
            wordLength = 16;
            fractionLength = 0;

        case 'int32'
            signed = 1;
            wordLength = 32;
            fractionLength = 0;

        case 'uint32'
            signed = 0;
            wordLength = 32;
            fractionLength = 0;

        case 'logical'
            signed = 0;
            wordLength = 1;
            fractionLength = 0;

        otherwise
            signed = 0;
            wordLength = 8;
            fractionLength = 0;  

    end

    frameInCasted = fi(frameIn, signed, wordLength, fractionLength, hdlfimath);

end


if(isfloat(frameIn))
    gradientOutput = cast(0, 'like', frameIn);
else
    gradientOutput = fi(0, 1, frameInCasted.WordLength + extraBitsWord, ...
        frameInCasted.FractionLength+extraBitsFrac, 'RoundingMethod','Floor','OverflowAction', 'Wrap');
end

thresholdD = sum(gradientOutput*gradientOutput,1,'native');

% Limit wordLength for fixed point threshold to 128.
if isfi(thresholdD)
    thresholdDBK = fi(0,...
        issigned(thresholdD),...
        min(128,(thresholdD.WordLength)*2),...
        (thresholdD.FractionLength),...
        'RoundingMethod','Floor',...
        'OverflowAction', 'Wrap');
else
    thresholdDBK = thresholdD;
end

thresholdSquare = threshold*threshold;
thresholdSq = cast(thresholdSquare, 'like', thresholdDBK);

if (isfi(thresholdDBK) && ~islogical(frameIn))
    if single(realmax(thresholdDBK)) < single((threshold)*(threshold))
        coder.internal.warning('visionhdl:hdledge:ThresholdOverflow', funcName);
    elseif (abs(single(thresholdSq)-single((threshold)*(threshold))) > 0) && ...
            (abs(single(thresholdSq)-single((threshold)*(threshold))) < single(eps(thresholdDBK)))
        coder.internal.warning('visionhdl:hdledge:ThresholdPrecision', funcName);
    end
end

end

%--- NPUFun ----------------------------------------------------------------
function [edgeImage, g1Out, g2Out] = EdgeDetectorCore(frameIn, method, kx, ky, boundaryMethod)
% Configure the npu core based on the input method.

% Declare compile time kernel constants based on the input data type.
if(isfloat(frameIn))
    sobelCoeffX = cast([1,0,-1;2,0,-2;1,0,-1]./8, 'like', frameIn);
    prewittCoeffX = cast(1/6, 'like', frameIn);
    robertsCoeff45 = cast([1,0;0,-1]/2, 'like', frameIn);
else
    if(isfi(frameIn))
        frameInfimath = fimath('RoundingMethod', frameIn.RoundingMethod,'OverflowAction', frameIn.OverflowAction);
    else
        frameInfimath = hdlfimath;
    end
    sobelCoeffX = fi([1,0,-1;2,0,-2;1,0,-1]./8, 1, 4, 3, frameInfimath);
    prewittCoeffX = fi(1/6, 0, 16, 18, frameInfimath);
    robertsCoeff45 = fi([1,0;0,-1]/2, 1, 2, 1, frameInfimath);
end


switch lower(method)

    case 'sobel'
        [edgeImage,g1Out, g2Out] =  hdl.npufun(@sobelCore, [kx,  ky], frameIn, 'KernelArg', sobelCoeffX, 'KernelArg', sobelCoeffX', 'BoundaryMethod', boundaryMethod);

    case 'prewitt'
        [edgeImage, g1Out, g2Out] = hdl.npufun(@prewittCore, [kx,ky], frameIn, 'KernelArg', prewittCoeffX, 'BoundaryMethod', boundaryMethod);

    case 'roberts'
        [edgeImage,g1Out, g2Out] =  hdl.npufun(@robertsCore, [kx,  ky], frameIn, 'KernelArg', robertsCoeff45, 'KernelArg', fliplr(robertsCoeff45), 'BoundaryMethod', boundaryMethod);

    otherwise
        [edgeImage,g1Out, g2Out] =  hdl.npufun(@sobelCore, [kx,  ky], frameIn, 'KernelArg', sobelCoeffX, 'KernelArg', sobelCoeffX', 'BoundaryMethod', boundaryMethod);

end

end


function [edgeImage, g1Out, g2Out] = sobelCore(frameIn, coeff1, coeff2)

g1Out = coder.hdl.pipeline(sum(frameIn(:).*coeff1(:), 'native'));
g2Out = coder.hdl.pipeline(sum(frameIn(:).*coeff2(:), 'native'));

edgeImage = coder.hdl.pipeline(g1Out*g1Out + g2Out*g2Out);


end

function [edgeImage, g1Out, g2Out] = prewittCore(Img, prewittCoeff)

% Update frameIn to break the signed range during computation of the output
if ~isfloat(Img)
    frameIn = fi(Img, 1,Img.WordLength+1, Img.FractionLength);
else
    frameIn = Img;
end

% Pipeline Stage 1 for g1Out
g1Stage11 = sum(coder.hdl.pipeline(frameIn(1)+frameIn(2)),'native');
g1Stage12 = sum(coder.hdl.pipeline(frameIn(3)-frameIn(7)),'native');
g1Stage13 = sum(coder.hdl.pipeline(frameIn(8)+frameIn(9)),'native');

% Pipeline Stage 2 for g1Out
g1Stage2 = sum(coder.hdl.pipeline(g1Stage11+coder.hdl.pipeline(g1Stage12-g1Stage13)),'native');

g1Out = sum(coder.hdl.pipeline(g1Stage2*prewittCoeff),'native');

% Pipeline Stage 1 for g2Out
g2Stage11 = sum(coder.hdl.pipeline(frameIn(1)+frameIn(4)),'native');
g2Stage12 = sum(coder.hdl.pipeline(frameIn(7)-frameIn(3)),'native');
g2Stage13 = sum(coder.hdl.pipeline(frameIn(6)+frameIn(9)),'native');

% Pipeline Stage 2 for g2Out
g2Stage2 = sum(coder.hdl.pipeline(g2Stage11+coder.hdl.pipeline(g2Stage12-g2Stage13)),'native');

g2Out = sum(coder.hdl.pipeline(g2Stage2*prewittCoeff),'native');

edgeImage = g1Out*g1Out + g2Out*g2Out;

end

function [edgeImage, g1Out, g2Out] = robertsCore(frameIn, coeff1, coeff2)

g1Out = coder.hdl.pipeline(sum(frameIn(:).*coeff1(:), 'native'));
g2Out = coder.hdl.pipeline(sum(frameIn(:).*coeff2(:), 'native'));

edgeImage = coder.hdl.pipeline(g1Out*g1Out + g2Out*g2Out);

end

function  edgeImageSum = accumulator(edgeImage)

out = fi(0,0,32,16);
edgeImageSum =  hdl.iteratorfun(@accumulatorCore, edgeImage, out);

end


function edgeImageSum =  accumulatorCore(edgeImage, edgeImageSum, idx)
edgeImageSum = edgeImageSum + edgeImage;
end



% --- Method Validation ---------------------------------------
function mustBeValidMethod(method)

coder.internal.assert(strcmpi(method,'Sobel') || ...
    strcmpi(method,'Prewitt') || ...
    strcmpi(method,'Roberts'), ...
    'visionhdl:hdledge:InvalidMethod', 'hdledge');
end
