【图像增强】基于 hessian特征和Frangi滤波实现血管图像增强附matlab代码

本文涉及的产品
视觉智能开放平台,分割抠图1万点
视觉智能开放平台,图像资源包5000点
视觉智能开放平台,视频资源包5000点
简介: 【图像增强】基于 hessian特征和Frangi滤波实现血管图像增强附matlab代码

✅作者简介:热爱科研的Matlab仿真开发者,修心和技术同步精进,matlab项目合作可私信。

🍎个人主页:Matlab科研工作室

🍊个人信条:格物致知。

更多Matlab仿真内容点击👇

智能优化算法  神经网络预测雷达通信 无线传感器

信号处理图像处理路径规划元胞自动机无人机

⛄ 内容介绍

In this paper we extend the Frangi filter[1] to recognize edges and do not enhance them. We give a theoretical framework for optimal scale selection and choice of the free parameters. We discuss discretization details concerning especially the discrete kernel used for building the scale-space and the choice of discrete scales. We present several experiments on phantom data to objectively and quantitatively compare and judge the filters. Experiments on real coronary angiograms enhance the improvement reached by the integration of the edge indicator.

⛄ 部分代码

function blobness = blobness3D(I, sigmas, spacing, tau, whiteondark)

% calculates blobness probability map (local sphericity) of a 3D input

% image

%

% blobness = blobness3D(V, sigmas, spacing, tau, whiteondark)

%

% inputs,

%   I : 3D image

%   sigmas : vector of scales on which the vesselness is computed

%   spacing : input image spacing resolution - during hessian matrix

%       computation, the gaussian filter kernel size in each dimension can

%       be adjusted to account for different image spacing for different

%       dimensions

%   tau : (between 0.5 and 1) : parameter that controls response uniformity

%       - lower tau -> more intense output response            

%   whiteondark: (true/false) : are vessels (tubular structures) bright on

%       dark background or dark on bright (default for 3D is true)

%

% outputs,

%   blobness: maximum blobness response over scales sigmas

%

% example:

%   B = blobness3D(I, 1:5, [1;1;1], 1, true);

%

% Function was written by T. Jerman, University of Ljubljana (October 2014)

% Based on code by D. Kroon, University of Twente (May 2009)


verbose = 1;


if nargin<5

   whiteondark = true; % default

end


I(~isfinite(I)) = 0;

I = single(I);


for j = 1:length(sigmas)

   

   if verbose  

       disp(['Current Filter Sigma: ' num2str(sigmas(j)) ]);

   end


   [Lambda1, Lambda2, Lambda3] = volumeEigenvalues(I,sigmas(j),spacing,whiteondark);

   

   % proposed filter

   Lambda3M = Lambda3;

   Lambda3M(Lambda3>=min(Lambda3(:)).*tau) = min(Lambda3(:)).*tau;            

   response = ((Lambda1.^2).*Lambda3M.* 27) ./ ((2*Lambda1+Lambda3M).^3);

   response(abs(Lambda1)>abs(Lambda3M)) = 1;


   response(Lambda1>=0) = 0;

   response(Lambda2>=0) = 0;

   response(Lambda3>=0) = 0;        

   response(~isfinite(response)) = 0;


   %keep max response

   if(j==1)

       blobness = response;

   else        

       blobness = max(blobness,response);

   end

       

   clear response Lambda1 Lambda2 Lambda3 Lambda3M    

   

end


blobness = blobness ./ max(blobness(:));

blobness(blobness < 1e-2) = 0;



function [Lambda1, Lambda2, Lambda3] = volumeEigenvalues(V,sigma,spacing,whiteondark)

% calculates the three eigenvalues for each voxel in a volume


% Calculate 3D hessian

[Hxx, Hyy, Hzz, Hxy, Hxz, Hyz] = Hessian3D(V,sigma,spacing);


% Correct for scaling

c=sigma.^2;

Hxx = c*Hxx; Hxy = c*Hxy;

Hxz = c*Hxz; Hyy = c*Hyy;

Hyz = c*Hyz; Hzz = c*Hzz;


if whiteondark == false

   c=-1;

   Hxx = c*Hxx; Hxy = c*Hxy;

   Hxz = c*Hxz; Hyy = c*Hyy;

   Hyz = c*Hyz; Hzz = c*Hzz;    

end


% reduce computation by computing vesselness only where needed

% S.-F. Yang and C.-H. Cheng, 揊ast computation of Hessian-based

% enhancement filters for medical images,?Comput. Meth. Prog. Bio., vol.

% 116, no. 3, pp. 215?25, 2014.

B1 = - (Hxx + Hyy + Hzz);

B2 = Hxx .* Hyy + Hxx .* Hzz + Hyy .* Hzz - Hxy .* Hxy - Hxz .* Hxz - Hyz .* Hyz;

B3 = Hxx .* Hyz .* Hyz + Hxy .* Hxy .* Hzz + Hxz .* Hyy .* Hxz - Hxx .* Hyy .* Hzz - Hxy .* Hyz .* Hxz - Hxz .* Hxy .* Hyz;


T = ones(size(B1));

T(B1<=0) = 0;

T(B2<=0) = 0;

T(B3<=0) = 0;

T(B1 .* B2 <= B3) = 0;


clear B1 B2 B3;


indeces = find(T==1);


Hxx = Hxx(indeces);

Hyy = Hyy(indeces);

Hzz = Hzz(indeces);

Hxz = Hxz(indeces);

Hyz = Hyz(indeces);

Hxy = Hxy(indeces);


% Calculate eigen values

[Lambda1i,Lambda2i,Lambda3i]=eig3volume(Hxx,Hxy,Hxz,Hyy,Hyz,Hzz);


% Free memory

clear Hxx Hyy Hzz Hxy Hxz Hyz;


Lambda1 = zeros(size(T));

Lambda2 = zeros(size(T));

Lambda3 = zeros(size(T));


Lambda1(indeces) = Lambda1i;

Lambda2(indeces) = Lambda2i;

Lambda3(indeces) = Lambda3i;


% some noise removal

Lambda1(~isfinite(Lambda1)) = 0;

Lambda2(~isfinite(Lambda2)) = 0;

Lambda3(~isfinite(Lambda3)) = 0;


Lambda1(abs(Lambda1) < 1e-4) = 0;

Lambda2(abs(Lambda2) < 1e-4) = 0;

Lambda3(abs(Lambda3) < 1e-4) = 0;



function [Dxx, Dyy, Dzz, Dxy, Dxz, Dyz] = Hessian3D(Volume,Sigma,spacing)

%  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(Volume,Sigma,spacing)

%

% 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.

%   spacing : input image spacing

%

% 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=imbigaussian(Volume,Sigma,0.5);

   F=imgaussian(Volume,Sigma,spacing);

else

   F=Volume;

end


% Create first and second order diferentiations

Dz=gradient3(F,'z');

Dzz=(gradient3(Dz,'z'));

clear Dz;


Dy=gradient3(F,'y');

Dyy=(gradient3(Dy,'y'));

Dyz=(gradient3(Dy,'z'));

clear Dy;


Dx=gradient3(F,'x');

Dxx=(gradient3(Dx,'x'));

Dxy=(gradient3(Dx,'y'));

Dxz=(gradient3(Dx,'z'));

clear Dx;


function D = gradient3(F,option)

% 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:

%

% Fx = gradient3(F,'x');


[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

       

function I=imgaussian(I,sigma,spacing,siz)

% IMGAUSSIAN filters an 1D, 2D color/greyscale or 3D image with an

% Gaussian filter. This function uses for filtering IMFILTER or if

% compiled the fast  mex code imgaussian.c . Instead of using a

% multidimensional gaussian kernel, it uses the fact that a Gaussian

% filter can be separated in 1D gaussian kernels.

%

% J=IMGAUSSIAN(I,SIGMA,SIZE)

%

% inputs,

%   I: The 1D, 2D greyscale/color, or 3D input image with

%           data type Single or Double

%   SIGMA: The sigma used for the Gaussian kernel

%   SIZE: Kernel size (single value) (default: sigma*6)

%

% outputs,

%   J: The gaussian filtered image

%

% note, compile the code with: mex imgaussian.c -v

%

% example,

%   I = im2double(imread('peppers.png'));

%   figure, imshow(imgaussian(I,10));

%

% Function is written by D.Kroon University of Twente (September 2009)


if(~exist('siz','var')), siz=sigma*6; end


if(sigma>0)


   % Filter each dimension with the 1D Gaussian kernels\

   x=-ceil(siz/spacing(1)/2):ceil(siz/spacing(1)/2);

   H = exp(-(x.^2/(2*(sigma/spacing(1))^2)));

   H = H/sum(H(:));    

   Hx=reshape(H,[length(H) 1 1]);

   

   x=-ceil(siz/spacing(2)/2):ceil(siz/spacing(2)/2);

   H = exp(-(x.^2/(2*(sigma/spacing(2))^2)));

   H = H/sum(H(:));    

   Hy=reshape(H,[1 length(H) 1]);


   x=-ceil(siz/spacing(3)/2):ceil(siz/spacing(3)/2);

   H = exp(-(x.^2/(2*(sigma/spacing(3))^2)));

   H = H/sum(H(:));    

   Hz=reshape(H,[1 1 length(H)]);

   

   I=imfilter(imfilter(imfilter(I,Hx, 'same' ,'replicate'),Hy, 'same' ,'replicate'),Hz, 'same' ,'replicate');

end

⛄ 运行结果

⛄ 参考文献

1. [T. Jerman, F. Pernus, B. Likar, Z. Spiclin, "*Enhancement of Vascular Structures in 3D and 2D Angiographic Images*", IEEE Transactions on Medical Imaging, 35(9), p. 2107-2118 (2016)

2. [T. Jerman, F. Pernus, B. Likar, Z. Spiclin, "*Blob Enhancement and Visualization for Improved Intracranial Aneurysm Detection*", IEEE Transactions on Visualization and Computer Graphics, 22(6), p. 1705-1717 (2016)

3. [T. Jerman, F. Pernus, B. Likar, Z. Spiclin, "*Beyond Frangi: an improved multiscale vesselness filter*", Proc. SPIE 9413, Medical Imaging 2015: Image Processing, 94132A (2015)

❤️ 关注我领取海量matlab电子书和数学建模资料
❤️部分理论引用网络文献,若有侵权联系博主删除


相关文章
|
4天前
|
算法 数据安全/隐私保护 计算机视觉
基于FPGA的图像双线性插值算法verilog实现,包括tb测试文件和MATLAB辅助验证
本项目展示了256×256图像通过双线性插值放大至512×512的效果,无水印展示。使用Matlab 2022a和Vivado 2019.2开发,提供完整代码及详细中文注释、操作视频。核心程序实现图像缩放,并在Matlab中验证效果。双线性插值算法通过FPGA高效实现图像缩放,确保质量。
|
1月前
|
算法 数据安全/隐私保护 计算机视觉
基于Retinex算法的图像去雾matlab仿真
本项目展示了基于Retinex算法的图像去雾技术。完整程序运行效果无水印,使用Matlab2022a开发。核心代码包含详细中文注释和操作步骤视频。Retinex理论由Edwin Land提出,旨在分离图像的光照和反射分量,增强图像对比度、颜色和细节,尤其在雾天条件下表现优异,有效解决图像去雾问题。
|
1月前
|
算法 人机交互 数据安全/隐私保护
基于图像形态学处理和凸包分析法的指尖检测matlab仿真
本项目基于Matlab2022a实现手势识别中的指尖检测算法。测试样本展示无水印运行效果,完整代码含中文注释及操作视频。算法通过图像形态学处理和凸包检测(如Graham扫描法)来确定指尖位置,但对背景复杂度敏感,需调整参数PARA1和PARA2以优化不同手型的检测精度。
|
6月前
|
安全
【2023高教社杯】D题 圈养湖羊的空间利用率 问题分析、数学模型及MATLAB代码
本文介绍了2023年高教社杯数学建模竞赛D题的圈养湖羊空间利用率问题,包括问题分析、数学模型建立和MATLAB代码实现,旨在优化养殖场的生产计划和空间利用效率。
272 6
【2023高教社杯】D题 圈养湖羊的空间利用率 问题分析、数学模型及MATLAB代码
|
4月前
|
算法 数据安全/隐私保护
织物图像的配准和拼接算法的MATLAB仿真,对比SIFT,SURF以及KAZE
本项目展示了织物瑕疵检测中的图像拼接技术,使用SIFT、SURF和KAZE三种算法。通过MATLAB2022a实现图像匹配、配准和拼接,最终检测并分类织物瑕疵。SIFT算法在不同尺度和旋转下保持不变性;SURF算法提高速度并保持鲁棒性;KAZE算法使用非线性扩散滤波器构建尺度空间,提供更先进的特征描述。展示视频无水印,代码含注释及操作步骤。
|
5月前
|
算法 数据可视化 数据安全/隐私保护
基于LK光流提取算法的图像序列晃动程度计算matlab仿真
该算法基于Lucas-Kanade光流方法,用于计算图像序列的晃动程度。通过计算相邻帧间的光流场并定义晃动程度指标(如RMS),可量化图像晃动。此版本适用于Matlab 2022a,提供详细中文注释与操作视频。完整代码无水印。
|
6月前
|
数据采集 存储 移动开发
【2023五一杯数学建模】 B题 快递需求分析问题 建模方案及MATLAB实现代码
本文介绍了2023年五一杯数学建模竞赛B题的解题方法,详细阐述了如何通过数学建模和MATLAB编程来分析快递需求、预测运输数量、优化运输成本,并估计固定和非固定需求,提供了完整的建模方案和代码实现。
138 0
【2023五一杯数学建模】 B题 快递需求分析问题 建模方案及MATLAB实现代码
|
6月前
|
存储 算法 搜索推荐
【2022年华为杯数学建模】B题 方形件组批优化问题 方案及MATLAB代码实现
本文提供了2022年华为杯数学建模竞赛B题的详细方案和MATLAB代码实现,包括方形件组批优化问题和排样优化问题,以及相关数学模型的建立和求解方法。
162 3
【2022年华为杯数学建模】B题 方形件组批优化问题 方案及MATLAB代码实现
|
9月前
|
数据安全/隐私保护
耐震时程曲线,matlab代码,自定义反应谱与地震波,优化源代码,地震波耐震时程曲线
地震波格式转换、时程转换、峰值调整、规范反应谱、计算反应谱、计算持时、生成人工波、时频域转换、数据滤波、基线校正、Arias截波、傅里叶变换、耐震时程曲线、脉冲波合成与提取、三联反应谱、地震动参数、延性反应谱、地震波缩尺、功率谱密度
基于混合整数规划的微网储能电池容量规划(matlab代码)
基于混合整数规划的微网储能电池容量规划(matlab代码)