function [e0,e1,theta] = det2MM(im,sigmaI,sigmaO) % function [e0,e1,theta] = det2MM(im,sigmaI,sigmaO) % % Compute the eigenspectrum of the spatially averaged second moment % matrix. % % INPUT % im Image. % sigmaI Inner scale (sigma for image derivatives). % sigmaO Outer scale (sigma for spatial averaging). % % OUTPUT % e0,e1 Smaller,larger eigenvalues. % theta Orientation of fist eigenvector + pi/2 % (i.e. orientation of possible edge). % % David R. Martin % March 2003 if ndims(im)>2, im = rgb2gray(im); end; idiag = norm(size(im)); if nargin<2, sigmaI=2; end if nargin<3, sigmaO=sigmaI; end % Scott Konishi says this is close % to optimal. sigmaI = max(0.5,sigmaI); sigmaO = max(0.5,sigmaO); % compute x and y image derivatives at inner scale fb = cell(2,1); fb{1} = oeFilter(sigmaI,3,pi/2,1); fb{2} = fb{1}'; fim = fbRun(fb,im); dx = fim{1}; dy = fim{2}; % compute smoothed squared image derivatives at outer scale f = oeFilter(sigmaO,3); dx2 = applyFilter(f,dx.^2); dy2 = applyFilter(f,dy.^2); dxy = applyFilter(f,dx.*dy); % compute eigenvalues of the spatially averaged 2nd moment matrix % and the orientations of the eigenvectors k = sqrt( (dx2-dy2).^2 + 4.*dxy.^2 ); eig0 = (dx2 + dy2 - k) / 2; eig1 = (dx2 + dy2 + k) / 2; t0 = atan2( dx2-eig0, -dxy ); t1 = atan2( dx2-eig1, -dxy ); % order eigenvalues by their absolute value, so e0<=e1, and pick % out the orientation corresponding to the largest eigenvalue x = (abs(eig1) > abs(eig0)); e0 = abs(eig0.*x + eig1.*~x); e1 = abs(eig1.*x + eig0.*~x); theta = t1.*x + t0.*~x; theta = mod(theta+pi/2,pi); % check postconditions if any(e0>e1), error('e0>e1'); end