Canny算子理解,及Matlab实现
2016-04-10 23:53
447 查看
JohnCanny于1986年提出Canny算子,它与Marr(LoG)边缘检测方法类似,也属于是先平滑后求导数的方法。本文对根据上述的边缘检测过程对Canny检测算法的原理进行介绍。 并结合实验,对canny算法进行理解,指出参数选取中容易产生的问题。
canny边缘检测一共四个部分:
1. 用高斯滤波器平滑图像;(图像去噪)
2. 用一阶偏导有限差分计算梯度幅值和方向;(特征增强)
3. 对梯度幅值进行非极大值抑制 ;(边缘检测)
4. 用双阈值算法检测和连接边缘。(形态学处理)
值得注意的是,在选取滤波窗口w时,根据我的实验(无噪声情况实验,若有噪声,也会影响后续的非极大值抑制),滤波窗口w边长为偶数时,将导致得到的梯度只有单像素最大幅值。而对边缘分析可知,应该产生双像素的梯度幅值才符合实际情况,所以滤波窗口w的选取一定要为基数边长。
实例如下,当选偶数边长时,幅值结果如下
当选基数边长时,幅值结果如下
值得一提的是,强边缘的阈值设置后并不需要包括所有边缘,只需要在弱边缘上有“点状”的分布即可,后续会根据与弱边缘的8连通域情况,得到效果很好的边缘。 另外,对得到的强边缘进行开操作和闭操作可以得到很好的去噪效果,而同样的操作却不能对弱边缘使用。
弱边缘图像:
强边缘图像:
代码如下:
1) 0:左边 和 右边
2)45:右上 和 左下
3)90: 上边 和 下边
4)135: 左上 和 右下
这样做的好处是简单, 但是这种简化的方法无法达到最好的效果, 因为,自然图像中的边缘梯度方向不一定是沿着这四个方向的。因此,就有很大的必要进行插值,找出在一个像素点上最能吻合其所在梯度方向的两侧的像素值。
非极大值抑制是canny算法的精髓,对此理解不深的同学,可以参考(/article/1497953.html) 这里的解释,解释得非常到位。
值得一提的是,为了在噪声中得到真实的双像素边缘,可以在非极大值抑制时不选取最大值为1,只要大小为周围最大值的80%即可。
代码如下:
原图:
canny算法检测到的边界:
代码如下
[rstrong,cstrong] = find(strong_edge);
edge = bwselect(weak_edge, cstrong, rstrong, 8);
在Matlab “edge()”函数中,最后还做了一个形态学瘦化处理,edge= bwmorph(edge, ‘thin’, 1); 但是我觉得得到的边缘已经足够细了,双像素边缘在之后对边缘区域的操作中也会更加有用。
canny边缘检测一共四个部分:
1. 用高斯滤波器平滑图像;(图像去噪)
2. 用一阶偏导有限差分计算梯度幅值和方向;(特征增强)
3. 对梯度幅值进行非极大值抑制 ;(边缘检测)
4. 用双阈值算法检测和连接边缘。(形态学处理)
平滑去噪
canny边缘检测的前两步相对不复杂,所以我就直接调用系统函数了。IM=imread('2.tif'); IM=imnoise(IM); %imnoise对uint8类型加噪为0-255,e对double类型加噪为0-1, [m,n]=size(IM); IM=double(IM); %高斯滤波 w=fspecial('gaussian',[9 9]); img=imfilter(IM,w,'replicate');
值得注意的是,在选取滤波窗口w时,根据我的实验(无噪声情况实验,若有噪声,也会影响后续的非极大值抑制),滤波窗口w边长为偶数时,将导致得到的梯度只有单像素最大幅值。而对边缘分析可知,应该产生双像素的梯度幅值才符合实际情况,所以滤波窗口w的选取一定要为基数边长。
实例如下,当选偶数边长时,幅值结果如下
当选基数边长时,幅值结果如下
计算幅值
%%sobel边缘检测 w_h=[1,2,1;0,0,0;-1,-2,-1]; %Sobel算子,梯度方向朝上为正 img_h=imfilter(img,w_h,'replicate'); %梯度是竖着的边缘,即横边缘 w_w=[-1,0,1;-2,0,2;-1,0,1]; %Sobel算子,梯度方向朝右为正 img_w=imfilter(img,w_w,'replicate'); %梯度是横着的边缘,,即竖边缘 img_grad=sqrt(img_w.^2+img_h.^2); %梯度的绝对值 grad_max=max(max(img_grad)); img_grad=img_grad/grad_max; %归一化
阈值选取
阈值用于后续的边缘连接处理,高阈值代表强边缘,低阈值代表弱边缘,在OpenCV中,阈值需要人工选取。在Matlab中,阈值根据梯度幅值的直方图选取,默认选择幅值大小在前70%的最后一个幅值作为高阈值,低阈值与高阈值的比例为4:10.在我的程序中,由于需要检测很强的边缘,故直方图比例设置得很高。值得一提的是,强边缘的阈值设置后并不需要包括所有边缘,只需要在弱边缘上有“点状”的分布即可,后续会根据与弱边缘的8连通域情况,得到效果很好的边缘。 另外,对得到的强边缘进行开操作和闭操作可以得到很好的去噪效果,而同样的操作却不能对弱边缘使用。
弱边缘图像:
强边缘图像:
代码如下:
k=1.2; PercentOfPixelsNotEdges = 1-k*(m+n)/(m*n);%0.995; %强边缘的比例 ThresholdRatio = 0.52; %强弱边缘的比例 thresh=[ ]; if isempty(thresh)%通过直方图自动计算高低阈值大小 counts=imhist(img_grad, 256); highThresh = find(cumsum(counts) > PercentOfPixelsNotEdges*m*n,1,'first') / 256; %PercentOfPixelsNotEdges=0.8,即不是边界的比例 lowThresh = ThresholdRatio*highThresh; thresh = [lowThresh highThresh]; elseif length(thresh)==1 highThresh = thresh; if thresh>=1 error(message('images:edge:thresholdMustBeLessThanOne')) end lowThresh = ThresholdRatio*thresh; thresh = [lowThresh highThresh]; elseif length(thresh)==2 lowThresh = thresh(1); highThresh = thresh(2); if (lowThresh >= highThresh) || (highThresh >= 1) error(message('images:edge:thresholdOutOfRange')) end end
非极大值抑制
在John Canny提出的Canny算子的论文中,非最大值抑制就只是在0、90、45、135四个梯度方向上进行的,每个像素点梯度方向按照相近程度用这四个方向来代替。这种情况下,非最大值抑制所比较的相邻两个像素就是:1) 0:左边 和 右边
2)45:右上 和 左下
3)90: 上边 和 下边
4)135: 左上 和 右下
这样做的好处是简单, 但是这种简化的方法无法达到最好的效果, 因为,自然图像中的边缘梯度方向不一定是沿着这四个方向的。因此,就有很大的必要进行插值,找出在一个像素点上最能吻合其所在梯度方向的两侧的像素值。
非极大值抑制是canny算法的精髓,对此理解不深的同学,可以参考(/article/1497953.html) 这里的解释,解释得非常到位。
值得一提的是,为了在噪声中得到真实的双像素边缘,可以在非极大值抑制时不选取最大值为1,只要大小为周围最大值的80%即可。
代码如下:
%%下面是非极大值抑制 weak_edge2=zeros(m+2,n+2); img_grad2 = padarray(img_grad,[1 1],'replicate'); %为了计算边界上的点,要扩充一下 img_h2 = padarray(img_h,[1 1],'replicate'); img_w2 = padarray(img_w,[1 1],'replicate'); for i=2:m+1 % 图像边缘的像素不能计算 for j=2:n+1 if img_grad2(i,j)<=lowThresh M1=2; M2=2;%这时检测点肯定是不大于1的,故isbigger=0 else %img_grad(i,j)~=0 Mx=img_w2(i,j);%梯度是横着的 My=img_h2(i,j);%梯度是竖着的 if abs(My)>=abs(Mx) %y方向梯度比x方向梯度大,梯度方向“竖着”的情况 g1=img_grad2(i-1,j);%g1是检测点正上方一行 g3=img_grad2(i+1,j);%g3是检测点正下方一行 weight=(abs(Mx)/abs(My)); if Mx*My>=0 %x,y方向梯度同号,梯度线分布在第1、3象限,当Mx*My==0时,Mx=0,weight=0 g2=img_grad2(i-1,j+1);%g2是检测点上右方 g4=img_grad2(i+1,j-1);%g4是检测点下左方 end if Mx*My<0 %x,y方向梯度异号,梯度线分布在第2、4象限 g2=img_grad2(i-1,j-1);%g2是检测点上左方 g4=img_grad2(i+1,j+1);%g4是检测点下右方 end M1=g1*(1-weight)+g2*weight;%M1是上方的插值 M2=g3*(1-weight)+g4*weight;%M2是下方的插值 end if abs(My)<abs(Mx)%x方向梯度比y方向梯度大,梯度方向“横着”的情况 g1=img_grad2(i,j+1);%g1是检测点正右方 g3=img_grad2(i,j-1);%g3是检测点正左方 weight=(abs(My)/abs(Mx)); if Mx*My>=0 %x,y方向梯度同号,梯度线分布在第1、3象限 g2=img_grad2(i-1,j+1);%g2是检测点上右方 g4=img_grad2(i+1,j-1);%g4是检测点下左方 end if Mx*My<0 %x,y方向梯度异号,梯度线分布在第2、4象限 g2=img_grad2(i+1,j+1);%g2是检测点下右方 g4=img_grad2(i-1,j-1);%g4是检测点上左方 end M1=g1*(1-weight)+g2*weight;%M1是上方的插值 M2=g3*(1-weight)+g4*weight;%M2是下方的插值 end end %此处选择与经典canny不同,考虑到噪声干扰,为了找到“双线边缘”,假定误差不超过0.9的都同为最大值 %也可根据加入噪声水平确定该值 isbigger=(img_grad2(i,j)>=(0.8*M1))&&(img_grad2(i,j)>=(0.8*M2)); %如果当前点比两边点都大,则置为白色 if isbigger weak_edge2(i,j)=255; end end end weak_edge=weak_edge2(2:m+1,2:n+1); figure(1); imshow(weak_edge) [rstrong,cstrong] = find(img_grad>highThresh & weak_edge); strong_edge=zeros(m,n); for i=1:length(rstrong) r=rstrong(i); c=cstrong(i); strong_edge(r,c)=weak_edge(r,c); end se=strel('square',2); strong_edge=imopen(strong_edge,se); se=strel('square',3); strong_edge=imclose(strong_edge,se); figure(2); imshow(strong_edge)
形态学处理
最后,以强边缘为种子,在弱边缘中寻找其8连通域,作为边缘即可。原图:
canny算法检测到的边界:
代码如下
[rstrong,cstrong] = find(strong_edge);
edge = bwselect(weak_edge, cstrong, rstrong, 8);
在Matlab “edge()”函数中,最后还做了一个形态学瘦化处理,edge= bwmorph(edge, ‘thin’, 1); 但是我觉得得到的边缘已经足够细了,双像素边缘在之后对边缘区域的操作中也会更加有用。
相关文章推荐
- Machine Learning Logistic Regression and Newton's Method Andrew Ng 课程练习 Matlab Script 详细解析
- 给win7、win8安装matlab 2014
- MATLAB链接MinGW编译器
- 【matlab】meshgrid生成网格原理1
- ubuntu 15.10 安装matlab2014b
- matlab mean函数
- 【Matlab】matlab与matplotlib作图比较
- matlab 终止正在运行的程序
- 如何降低自己的gcc版本_caffe_matlabconfigure
- matlab中save,load使用方法
- 读取指定文件夹中所有文件名以及文件路径,并读到matlab
- Matlab GUI图像图像基础
- matlab安装遇到问题/install/Matlab/bin/util/oscheck.sh: /lib64/libc.so.6: not found
- Matlab并行编程方法
- matlab实现分水岭算法处理图像分割
- Matlab size函数
- matlab reshape函数
- matlab中实现Gabor滤波器
- MATLAB曲线绘制
- Matlab中image、imagesc和imshow函数用法解析