【图像检测】基于mom方法结合Hessian和曲线拟合的方法实现血管的直径并输出测量图像附matlab代码

matlab科研助手

已于 2022-08-03 15:58:25 修改

阅读量916

点赞数

CC 4.0 BY-SA版权

分类专栏：图像处理文章标签： matlab 开发语言

于 2022-07-09 21:49:24 首次发布

本文链接：https://blog.youkuaiyun.com/matlab_dingdang/article/details/125699363

图像处理专栏收录该内容

872 篇文章

订阅专栏

1 内容介绍

在骨架化提取出血管中心线的基础上,提出一种基于MoM评价模型的冠脉血管直径的跟踪测量方法.该算法利用血管的两条边缘线相对于中心线的对称性和最优化评价思想,通过对实际造影图像的量化测量并将其与实际直径值进行对比,发现使用提出的算法能够准确地测量血管的直径,并且可以有效地解决血管分支和重叠处直径测量问题.

2 仿真代码

% 交互选择 clc;clear all;close all %% hf = figure(1); A = imread('C:\Users\Administrator\Desktop\22.bmp'); A = rgb2gray(A); %-此句注释后下面的处理就相当于对真彩色图像进行 imshow(A) rect = getrect(hf); w = round(rect(3)); h = round(rect(4)); if w >= 1 && h >= 1 r = [rect(1),rect(1)+w,rect(1)+w,rect(1);rect(2),rect(2),rect(2)+h,rect(2)+h]; rectangle('Position',[rect(1),rect(2),w,h], 'edgecolor','red'); figure(2); bw = roipoly(A,r(1,:),r(2,:)); AA = reshape(A, [(size(A,1)*size(A,2)), size(A,3)]); BB = AA(bw, :); B = reshape(BB, [h, w, size(A,3)]); imshow(B); end %%

function [outIm,whatScale,Direction] = FrangiFilter2D(I, options) defaultoptions = struct('FrangiScaleRange', [1 10], 'FrangiScaleRatio', 2, 'FrangiBetaOne', 0.5, 'FrangiBetaTwo', 15, 'verbose',true,'BlackWhite',true);

if(~exist('options','var')), options=defaultoptions; else tags = fieldnames(defaultoptions); for i=1:length(tags) if(~isfield(options,tags{i})), options.(tags{i})=defaultoptions.(tags{i}); end end if(length(tags)~=length(fieldnames(options))), warning('FrangiFilter2D:unknownoption','unknown options found'); end end

sigmas=options.FrangiScaleRange(1):options.FrangiScaleRatio:options.FrangiScaleRange(2); sigmas = sort(sigmas, 'ascend'); beta = 2*options.FrangiBetaOne^2; c = 2*options.FrangiBetaTwo^2; ALLfiltered=zeros([size(I) length(sigmas)]); ALLangles=zeros([size(I) length(sigmas)]);

for i = 1:length(sigmas), if(options.verbose) %disp(['当前尺度空间的sigma值为: ' num2str(sigmas(i)) ]); end [Dxx,Dxy,Dyy] = Hessian2D(I,sigmas(i)); Dxx = (sigmas(i)^2)*Dxx; Dxy = (sigmas(i)^2)*Dxy; Dyy = (sigmas(i)^2)*Dyy; [Lambda2,Lambda1,Ix,Iy]=eig2image(Dxx,Dxy,Dyy); angles = atan2(Ix,Iy); Lambda1(Lambda1==0) = eps; Rb = (Lambda2./Lambda1).^2; S2 = Lambda1.^2 + Lambda2.^2; Ifiltered = exp(-Rb/beta) .*(ones(size(I))-exp(-S2/c)); if(options.BlackWhite) Ifiltered(Lambda1<0)=0; else Ifiltered(Lambda1>0)=0; end ALLfiltered(:,:,i) = Ifiltered; ALLangles(:,:,i) = angles; end

if length(sigmas) > 1, [outIm,whatScale] = max(ALLfiltered,[],3); outIm = reshape(outIm,size(I)); if(nargout>1) whatScale = reshape(whatScale,size(I)); end if(nargout>2) Direction = reshape(ALLangles((1:numel(I))'+(whatScale(:)-1)*numel(I)),size(I)); end else outIm = reshape(ALLfiltered,size(I)); if(nargout>1) whatScale = ones(size(I)); end if(nargout>2) Direction = reshape(ALLangles,size(I)); end end

function [k_s1,k_s2,new_xW,new_yW] = Seg_all(xW, yW)
sortflag = 1; % 递增标志1，递减标志0
jj = 2; % 同方向的不同数
kk = 1; % 换向换行
new_yW{1,1} = yW(1);
new_xW{1,1} = xW(1);
for i = 2:length(yW)

if yW(i) >= yW(i-1) && sortflag == 1
% if yW(i) == yW(i-1)
% new_yW{kk,jj} = yW(i)+0.001*jj;
% new_xW{kk,jj} = xW(i);
% else
new_yW{kk,jj} = yW(i);
new_xW{kk,jj} = xW(i);
% end
% ak(kk) = 1; % ak(kk) = 1;正向
jj = jj + 1;
sortflag = 1;
else
if sortflag == 1
kk = kk+1;
jj = 1;
end
sortflag = 0;
if yW(i) <= yW(i-1) && sortflag == 0
% if yW(i) == yW(i-1)
% new_yW{kk,jj} = yW(i)-0.001*jj;
% new_xW{kk,jj} = xW(i);
% else
new_yW{kk,jj} = yW(i);
new_xW{kk,jj} = xW(i);
% end
% ak(kk) = 0; % ak(kk) = 0;反向
jj = jj + 1;
sortflag = 0;
else
kk = kk+1;
jj = 1;
sortflag = 1;
% 由反转为正则+1
new_yW{kk,jj} = yW(i); % 转向补植
new_xW{kk,jj} = xW(i);
jj = jj+1;
end
end
end

for j = 1 : size(new_yW,1)
newy = cell2mat(new_yW(j,:));
newx = cell2mat(new_xW(j,:));
nny = categorical(newy);
samenum_y = tabulate(nny); % 判断同一值出现频率
hold on
if ~isempty(find(cell2mat(samenum_y(:,3)) > 50)) % 若同一x值出现次数过多
% sp=spap2(1,2,newy,newx);
p=polyfit(newy,newx,2);
% sp=csapi(newy,newx);
y = linspace(min(newy),max(newy),length(newy));
for i = 1:length(y)
k_s2{j,i} = 2*p(1)*y(i)+p(2); % 血管流向斜率
k_s1{j,i} = -1./k_s2{j,i}; % 法一：拟合曲线求法线斜率
end

else
% sp=spap2(1,4,newy,newx); % 拟合
sp=csapi(newy,newx);
fnplt(sp,'b-');
y = linspace(min(newy),max(newy),length(newy));
for i = 1:length(y)
k_s1{j,i} = fnval(fnder(sp), y(i)); % 血管流向斜率
k_s2{j,i} = -1./k_s1{j,i}; % 法一：拟合曲线求法线斜率
end
end

% sp=spap2(1,4,newx,newy); % 拟合
% plot(newy,newx,'r.')
hold off
end

%
% % % figure;
% for j = 1 : size(new_yW,1)
% newy = cell2mat(new_yW(j,:));
% newx = cell2mat(new_xW(j,:));
% if length(newy)<2
% k_s1 = 1;
% k_s2 = -1;
% else
% k_s1{j,1} = (newy(2)-newy(1))/(newx(2)-newx(1));
% k_s2{j,1} = -1/k_s1{j,1};
% for k = 2:length(newy)
% k_s1{j,k} = (newy(k)-newy(k-1))/(newx(k)-newx(k-1));
% k_s2{j,k}= -1/k_s1{j,k};
% end
% % hold on
% % plot(newy, newx, 'r.')
% end
% end

% 图像掩模 clc;clear all;close all; %% 读图 Image = imread('..\test_imgs\IMG_0062.tif'); im = imread('..\seg_result\IMG_0062_out.png'); se = strel('disk',1); im = imerode(im,se); IM(:,:,1) = im2double(im) .* double(Image(:,:,1)); IM(:,:,2) = im2double(im) .* double(Image(:,:,2)); IM(:,:,3) = im2double(im) .* double(Image(:,:,3)); % figure; % subplot(1,2,1);imshow(Image) % subplot(1,2,2);imshow(IM) %% 做标记 E = double(Image(:,:,1))./double(max(Image(:))); I = double(IM(:,:,1))./double(max(IM(:))); % Here is overlaying using AlphaData (Opacity):