Code covered by the BSD License

Hessian based Frangi Vesselness filter

Dirk-Jan Kroon (view profile)

11 Jun 2009 (Updated )

Enhancement of Vessel/ridge like structures in 2D/3D image using hessian eigen values

Hessian3D(Volume,Sigma)
```function [Dxx, Dyy, Dzz, Dxy, Dxz, Dyz] = Hessian3D(Volume,Sigma)
%  This function Hessian3D filters the image with an Gaussian kernel
%  followed by calculation of 2nd order gradients, which aprroximates the
%  2nd order derivatives of the image.
%
% [Dxx, Dyy, Dzz, Dxy, Dxz, Dyz] = Hessian3D(I,Sigma)
%
% inputs,
%   I : The image volume, class preferable double or single
%   Sigma : The sigma of the gaussian kernel used. If sigma is zero
%           no gaussian filtering.
%
% outputs,
%   Dxx, Dyy, Dzz, Dxy, Dxz, Dyz: The 2nd derivatives
%
% Function is written by D.Kroon University of Twente (June 2009)
% defaults
if nargin < 2, Sigma = 1; end

if(Sigma>0)
F=imgaussian(Volume,Sigma);
else
F=Volume;
end

% Create first and second order diferentiations
clear Dz;

clear Dy;

clear Dx;

% This function does the same as the default matlab "gradient" function
% but with one direction at the time, less cpu and less memory usage.
%
% Example:
%

[k,l,m] = size(F);
D  = zeros(size(F),class(F));

switch lower(option)
case 'x'
% Take forward differences on left and right edges
D(1,:,:) = (F(2,:,:) - F(1,:,:));
D(k,:,:) = (F(k,:,:) - F(k-1,:,:));
% Take centered differences on interior points
D(2:k-1,:,:) = (F(3:k,:,:)-F(1:k-2,:,:))/2;
case 'y'
D(:,1,:) = (F(:,2,:) - F(:,1,:));
D(:,l,:) = (F(:,l,:) - F(:,l-1,:));
D(:,2:l-1,:) = (F(:,3:l,:)-F(:,1:l-2,:))/2;
case 'z'
D(:,:,1) = (F(:,:,2) - F(:,:,1));
D(:,:,m) = (F(:,:,m) - F(:,:,m-1));
D(:,:,2:m-1) = (F(:,:,3:m)-F(:,:,1:m-2))/2;
otherwise
disp('Unknown option')
end
```