GRD2GST_BASE - Base function for GRD2GST.

Contents

Syntax

   [gx2, gy2, gxy] = GRD2GST_BASE(gx, gy, rho, int, hsize, ...
                                  samp, tn, thez, sigt, eign);
   [gx2, gy2, gxy, mag, or] = GRD2GST_BASE(gx, gy, rho, int, hsize, ...
                                         samp, tn, thez, sigt, eign);

See also

Ressembles: GRD2GST. Requires: GSTDECOMP, GSTFEATURE_BASE, SMOOTHFILT_BASE.

Function implementation

function [gx2, gy2, gxy, varargout] = ...
    grd2gst_base(gx, gy, rho, int, hsize, samp, tn, thez, sigt, eign)

[X Y C] = size(gx);                                                    %#ok

compute the component of J(grad u_sigma) while taking into account the multichannel information by averaging over the different channels

gx2 = mean(gx.^2, 3);
gxy = mean(gx.*gy, 3);
gy2 = mean(gy.^2, 3);

following the approach of [Scheun03], average the gradient over the channels, which is nothing else than the gradient of the average of the channels

if tn
    sigma0 = 16; % initial isotropy
    eta0 = 12; % intial anisotropy
    [l1,l2,e1,e2] = gstdecomp(gx2, gy2, gxy);
    l1 = l1*0 + sqrt(sigma0*eta0);
    l2 = l2*0 + sqrt(sigma0/eta0);
    [gx2, gy2, gxy] = gstdecomp(l1,l2,e1,e2);
end

compute the direction of smoothing (along edges) in the anisotropic case

if ischar(int) && any(strcmp(int,{'ani','hourglass'}))
    theta = mod(atan2(mean(gy,3),mean(gx,3)),pi);
else
    theta =[];
end

perform tensor smoothing: isotropic (classical) or anisotropic (improved) spatial averaging of the tensor and of the directional derivatives

if rho~=0 && ~(islogical(int) && ~int)
    gx2 = smoothfilt_base( gx2, rho, int, hsize, samp, thez, sigt, theta );
    gy2 = smoothfilt_base( gy2, rho, int, hsize, samp, thez, sigt, theta );
    gxy = smoothfilt_base( gxy, rho, int, hsize, samp, thez, sigt, theta );
end

if nargout >= 4
    % norm based on the eigenvalues
    varargout{1} = gstfeature_base(gx2, gy2, gxy, 'norm', eign, [], []);
    %  double angle of the gradient tensor
    if nargout==5
        varargout{2} = gstfeature_base(gx2, gy2, gxy, 'orient', [], [], []);
    end
end

end % end of grd2gst_base