[Image processing] Image skeleton extraction based on central axis transformation with matlab code

1 内容介绍

介绍了一种中轴变换算法.该算法修正了现存算法的缺点,可以快速去除二值图像中的冗余信息,准确抽取图像的单像素骨架,保留基本特征,特别是斜线和圆形焊盘特征.大量实验表明,该方法降低图像的复杂度,仅保留图像形状的拓扑信息,为后续处理提供理想的图像,提高整个系统的性能.​

2 完整代码

% 中轴算法 Medial Axis Transfrom%% by HPC_ZY 20190802clear; close all; clc​load data​%% matlab% 求骨架ticskeleton = bwmorph(im,'skeleton',Inf);t1 = toc;% 去毛刺ticskeleton1 = bwmorph(skeleton,'spur',16); t2 = toc;​figureshowres(im,skeleton,t1,1,2,1,'matlab库')showres(im,skeleton1,t2,1,2,2,'matlab库（消除毛刺）')​%% 迭代删除% 迭代删除法ticskeleton1 = method1(im);t1 = toc;​% 迭代删除法（查找表）ticskeleton2 = method2(im);t2 = toc;​% 迭代删除法（简化）ticskeleton3 = method3(im);t3 = toc;​figureshowres(im,skeleton1,t1,1,3,1,'迭代删除法')showres(im,skeleton2,t2,1,3,2,'迭代删除法（查找表）')showres(im,skeleton3,t3,1,3,3,'迭代删除法（简化）')​​%% ------------function--------------​% 迭代删除法function skeleton = method1(bw)​bw = double(bw); % 转为非逻辑类型,方便后面求和[M,N] = size(bw);skeleton = bw;deleCount = 1; % 用于记录一次迭代删除的像素个数while deleCount>0    % 删（右,下）    deleCount = 0;    label = zeros(M,N,'logical');    for i = 2:M-1        for j = 2:N-1            if skeleton(i,j) % 如果该点为1,判断周围像素                p = [skeleton(i-1,j),skeleton(i-1,j+1),skeleton(i,j+1),skeleton(i+1,j+1),...                    skeleton(i+1,j),skeleton(i+1,j-1),skeleton(i,j-1),skeleton(i-1,j-1)];                Np = sum(p);                Tp = length(find(([p(2:end),p(1)]-p) == 1));                if Np>1 && Np<7 && Tp==1 && p(1)*p(3)*p(5)==0 && p(7)*p(3)*p(5)==0                    deleCount = deleCount+1;                    label(i,j) = 1;                end            end        end    end    skeleton(label) = 0;    % 删（左,上）    label = zeros(M,N,'logical');    for i = 2:M-1        for j = 2:N-1            if skeleton(i,j) % 如果该点为1,判断周围像素                p = [skeleton(i-1,j),skeleton(i-1,j+1),skeleton(i,j+1),skeleton(i+1,j+1),...                    skeleton(i+1,j),skeleton(i+1,j-1),skeleton(i,j-1),skeleton(i-1,j-1)];                Np = sum(p);                Tp = length(find(([p(2:end),p(1)]-p) == 1));                if Np>1 && Np<7 && Tp==1 && p(1)*p(3)*p(7)==0 && p(1)*p(5)*p(7)==0                    deleCount = deleCount+1;                    label(i,j) = 1;                end            end        end    end    skeleton(label) = 0;   end​​end​% 迭代删除法（查找表）function skeleton = method2(bw)​bw = double(bw);[M,N] = size(bw);skeleton = bw;​% 生成查找表[List1,List2] = getList1();​[x,y] = find(bw);a = min(x);b = max(x);c = min(y);d = max(y);​deleCount = 1;mat = [1,2,4,8,16,32,64,128];while deleCount>0        deleCount = 0;    % 删（右,下）    label1 = zeros(M,N,'logical');    for i = a:b        for j = c:d            if skeleton(i,j)                p = [skeleton(i-1,j),skeleton(i-1,j+1),skeleton(i,j+1),skeleton(i+1,j+1),...                    skeleton(i+1,j),skeleton(i+1,j-1),skeleton(i,j-1),skeleton(i-1,j-1)];                idx = sum(p.*mat)+1;                if List1(idx)                    label1(i,j) = true;                    deleCount = deleCount+1;                end            end        end    end    skeleton(label1) = 0;        % 删（左,上）    label2 = zeros(M,N,'logical');    for i = a:b        for j = c:d            if skeleton(i,j)                p = [skeleton(i-1,j),skeleton(i-1,j+1),skeleton(i,j+1),skeleton(i+1,j+1),...                    skeleton(i+1,j),skeleton(i+1,j-1),skeleton(i,j-1),skeleton(i-1,j-1)];                idx = sum(p.*mat)+1;                if List2(idx)                    label2(i,j) = true;                    deleCount = deleCount+1;                end            end        end    end        skeleton(label2) = 0;end​end​% 迭代删除法（简化）function skeleton = method3(bw)​bw = double(bw);[M,N] = size(bw);skeleton = bw;% 生成查找表List = getList2();% 收缩搜索范围[x,y] = find(bw);a = min(x);b = max(x);c = min(y);d = max(y);% 开始迭代deleCount = 1;mat = [1,2,4,8,16,32,64,128];while deleCount>0     deleCount = 0;    label = zeros(M,N,'logical');    for i = a:b        for j = c:d            if skeleton(i,j)                p = [skeleton(i-1,j),skeleton(i-1,j+1),skeleton(i,j+1),skeleton(i+1,j+1),...                    skeleton(i+1,j),skeleton(i+1,j-1),skeleton(i,j-1),skeleton(i-1,j-1)];                idx = sum(p.*mat)+1;                if List(idx)                    label(i,j) = true;                    deleCount = deleCount+1;                end            end        end    end    skeleton(label) = 0;end​​end​% 生成查找表function [List1,List2] = getList1()​List1 = zeros(256,1,'logical');List2 = zeros(256,1,'logical');for n = 0:255    p = bitget(n,1:8);    Np = sum(p);    Tp = length(find(([p(2:end),p(1)]-p) == 1));    if Np>1 && Np<7 && Tp==1        % 右,下        if p(1)*p(3)*p(5)==0 && p(7)*p(3)*p(5)==0            List1(n+1) = true;        end        % 左,上        if p(1)*p(3)*p(7)==0 && p(1)*p(5)*p(7)==0            List2(n+1) = true;        end    end end​end​% 生成查找表（简化）function List = getList2()​List = zeros(256,1,'logical');for n = 0:255    p = bitget(n,1:8);    Np = sum(p);    Tp = length(find(([p(2:end),p(1)]-p) == 1));    if Np>1 && Np<7 && Tp==1        if (p(1)*p(3)*p(5)==0 && p(7)*p(3)*p(5)==0)||...                (p(1)*p(3)*p(7)==0 && p(1)*p(5)*p(7)==0)            List(n+1) = true;        end    end end​end​% 显示结果function showres(im,skeleton,t,i,j,k,titlename)​subplot(i,j,k),imshow(im)hold on[x,y] = find(skeleton);plot(y,x,'r.')​title(['\fontsize{16}',titlename])xlabel(['\fontsize{16}耗时(秒)：',num2str(t)])end​

3 运行结果

4 参考文献

[1]史聪伟, 赵杰煜, 常俊生. 基于中轴变换的骨架特征提取算法[J]. 计算机工程, 2019, 45(7):9.

博主简介：擅长智能优化算法、神经网络预测、信号处理、元胞自动机、图像处理、路径规划、无人机等多种领域的Matlab仿真,相关matlab代码问题可私信交流.

部分理论引用网络文献,若有侵权联系博主删除.