A function for computing Eigenvalues and the first Eigenvector of the structure tensor of an image.

structure_tensor2.m

%structure_tensor2
%
% Computes the eigenvalues and eigenvectors of the structure tensor of im.
%
function [e1 e2 ev1] = structure_tensor2(im, opt)

if ~exist('opt', 'var')
  opt = struct();
end

if ~isfield(opt, 'f_presmooth')
  opt.f_presmooth = fspecial('gaussian', 11,2);
end
if ~isfield(opt, 'f_tensorsmooth')
  opt.f_tensorsmooth = fspecial('gaussian', 11,2);
end

%
% build structure tensor for all image locations
%

% blur the image to get better gradients
im = imfilter(im, opt.f_presmooth, 'replicate', 'same');

[gx gy] = derivative5(im, 'x', 'y');

% compute structure tensor values
gx2 = gx.^2;
gxy = gx.*gy;
gy2 = gy.^2;


% handle RGB images
if ndims(gx) == 3
  gx2 = mean(gx2,3);
  gxy = mean(gxy,3);
  gy2 = mean(gy2,3);
end

% average structor tensor values
gx2 = imfilter(gx2, opt.f_tensorsmooth, 'same');
gxy = imfilter(gxy, opt.f_tensorsmooth, 'same');
gy2 = imfilter(gy2, opt.f_tensorsmooth, 'same');

%
% solve for eigenvalues.  Closed form for 2x2 matrix.
%

r1 = (gx2 + gy2)./2;
r2 = sqrt( 4.* gxy.^2 + (gx2 - gy2).^2)./2;

e1 = r1 + r2;
e2 = r1 - r2;

%
% solve for eigenvectors.
%

v1denom = gy2-e1;
v1ratio = -gxy./v1denom;

V1(:,:) = repmat(1, size(v1ratio));
V2(:,:) = v1ratio;

% hard code for left/right symmetry or no image texture
tmp = isinf(v1ratio) | isnan(v1ratio);
V1(tmp) = 0;
V2(tmp) = 1;

% make eigenvector unit length
ev1 = cat(3, V1, V2) ./ repmat(sqrt(V1.^2 + V2.^2),[1 1 2]);

Requires derivative5.m from https://www.peterkovesi.com/matlabfns/, but you could replace that with something else pretty easily.