角点过滤：一则简单地利用骨架局部线性回归的matlab实例

由于精确计数需要，需要对二值图腐蚀后提取骨架、去分叉后，调用matlab自带函数corner进行角点提取。涉及到的操作均为自带函数：

提取骨架函数：J=bwmorph(I,'skel',refine_time); refine_time为细化次数，取值为Inf时细化到只剩一个像素宽的骨架，输入和输出图象I和J均为二值图像，下同；
去骨架分叉函数：J=bwmorph(I,'spur',spur_time);spur为去分叉次数。
角点检测函数：C=corner(I,corner_num); corner_num为要提取的角点数，默认为200.输出为corner_num*2大小矩阵，每一行代表一个角点坐标。

二值图、骨架图、去分叉图、角点检测图依次如下：

图中红色星号为角点。角点显示代码如下：

    C=corner(Ispur,corner_num);
    figure
    imshow(Ispur)
    hold on
    plot(C(:,1),C(:,2),'r*');

由于输出的角点总数与输入的预期角点总数完全一致（除非预期角点总数超过所有角点数），所以需要对角点进行进一步的过滤：过滤掉分布在直线骨架附近的角点，并保留位于拐弯处、交叉处（分叉需要提前处理干净，不然影响结果）的角点，也是本文章的核心。算法介绍和示意图：

假设上图为输入的二值图像，方块表示一个像素，红色方块表示为一个角点像素，白色方块为构成骨架的一个像素，该算法选取以目标角点为中心，一定范围内（图中蓝色框，大小为2R*2R）的所有骨架上的像素（图中被标记为1,2,...6的白色像素），对这些像素点（包括角点像素在内）的坐标（xi,yi）进行线性回归，得到线性函数 $\hat{y}=\hat{a}x+\dot{\hat{b}}$ 以及与样本点间的回归系数、区间估计、残差、置信区间、判定系数、F统计量观测值、p值、误差方差的估计和显著性水平等。根据得到的某一个或多个值的线性组合进行判断是否大于或小于一个训练过或根据经验预设的阈值thre，本算法以判定系数R^2为例，设置阈值thre,如果在以该角点为中心一定范围内的骨架像素点线性回归的判定系数大于这个阈值thre，则认为该角点位于直线或类直线上，否则认定为位于拐弯处或交叉处。

扫描二维码关注公众号，回复： 13016146 查看本文章

PS：线性回归可通过自写函数实现，或者直接调用函数regress，调用方式为：

[b,bint,r,rint,stats]=regress(Y,X)

其中Y为因变量的列向量，X为增补后的自变量的列向量（即第一列为全1，第二列为x）。输出b为回归系数，bint为回归系数的区间估计，r为残差，rint为置信区间，stats为四4个值的矩阵，用于检验回归模型，stats的四个值分别解释为判定系数、F统计量观测值，检验的p的值，误差方差的估计。

测试图像截取为去分叉后的图像的局部：

设几个角点坐标（147,124）（163,146）（181,142)）（139,111）【下图红色圆圈为角点位置】

代码实现：

function JudgeList=LinearJudge(I,C,R)
    [m,n]=size(I);
    JudgeList=zeros(1,length(C));
    
    for c=1:length(C)
        pn=0;
        ci=C(c,2);
        cj=C(c,1);
        PointList=[];
        for i =ci-R:ci+R
            for j=cj-R:cj+R
                if i>=1 && i<=m && j>=1 && j<=n
                    if I(i,j)==255
                        pn=pn+1;
                        PointList(pn,:)=[i,j];
                    end
                end
            end
        end
        [b,bint,r,rint,stats]=regress(PointList(:,2),[ones(pn,1),PointList(:,1)]);
        JudgeList(c)=stats(1);
    end
end

输入为二值图像I（在这转为只有0和255的图像，如果为0和1的二值图像，则将255改为1即可），角点坐标矩阵C和范围半径R。PointList为每个角点范围内的骨架像素点坐标，输出为所有角点对应的值，在这取每次的判定系数stats(1)。

测试实例代码：

I=rgb2gray(imread('spur_01.jpg'));
I=uint8(im2bw(I)*255);%转为只有0和255的二值图像
C=[[147,124];[163,146];[181,142];[139,111]];
imshow(I)
hold on
plot(C(:,1),C(:,2),'ro')
R=ceil(sum(sum(I==255))/20);
result=LinearJudge(I,C,R);
thre=0.75;
JudgeList=LinearJudge(I,C,R);%调用算法
C2=C(find(JudgeList<=thre),:);
figure,
imshow(I)
hold on
plot(C2(:,1),C2(:,2),'go')

其中R取的是骨架长度除以20的向上取整。C2为满足条件（输出大于thre）的角点坐标，满足条件的角点显示为绿色圆圈：

应用于整个图像：(R=5,thre=0.72)

PS：仅为测试，具体两个参数取值可通过调试等方式获得。算法输出属性值可以选择残差等。

角点过滤：一则简单地利用骨架局部线性回归的matlab实例

猜你喜欢