影象平滑與邊緣檢測
阿新 • • 發佈:2018-12-20
一、概要:
使用Canny邊緣檢測演算法作為例子,介紹影象的平滑方法和邊緣檢測。
Canny邊緣檢測演算法分為四步:
step1:用高斯濾波器平滑影象; step2:用一階偏導的有限差分來計算梯度的幅值和方向;(在橫豎兩個方向上計算邊緣,再求平方和的開方) step3:對梯度幅值進行非極大值抑制; step4:用雙閾值演算法檢測和連線邊緣。
demo&效果:
原圖(lenna.bmp):
高斯濾波後的影象:
初步求邊緣後的影象:
非極大值抑制後的影象:
雙閾值檢測後的影象:
以下對這四步進行詳細介紹:
二、影象平滑與高斯濾波
影象平滑的目的是消除或儘量減少噪聲的影響,改善影象的質量。在假定加性噪聲是隨機獨立分佈的條件下,利用鄰域的平均或加權平均可以有效的抑制噪聲干擾。
方法是做鄰域計算,比如:
高斯濾波:
採用高斯函式作為加權函式。 原因一:二維高斯函式具有旋轉對稱性,保證濾波時各方向平滑程度相同; 原因二:離中心點越遠權值越小。確保邊緣細節不被模糊。
設定σ2和n,確定高斯模板權值。如σ2 =2和n=5,整數化和歸一化後得:
程式碼:
1.直接用系統函式的方法
- img=imread('lenna.bmp');
- imshow(img);
- [m n]=size(img);
- img=double(img);
- %高斯濾波
- w=fspecial('gaussian',[5 5]); %[ 5 5 ]是模板尺寸,預設是[ 3 3 ] 模板即上文所提的模板
- img=imfilter(img,w,'replicate');
- figure;
- imshow(uint8(img))
2.自己實現的高斯濾波(使用上面的模板)
- img=imread('lenna.bmp');
- img=im2double(img);
- [m n]=size(img);
- h=zeros(m,n);
- for i=1:m
- for j=1:n
- if i<3 || i>(m-3) || j<3 || j>(n-3)
- h(i,j)=img(i,j);
- else %程式碼長度原因,分行相加。
- h(i,j)=h(i,j)+img(i-2
,j-2)+img(i-2,j+2)+img(i+2,j-2)+img(i+2,j+2); - h(i,j)=h(i,j)+(img(i-1,j-2)+img(i+1,j-2)+img(i-2,j-1)+img(i+2,j-1)+img(i-2,j+1)+img(i+2,j+1)+img(i-1,j+2)+img(i+1,j+2))*2;
- h(i,j)=h(i,j)+3*(img(i,j-2)+img(i,j+2)+img(i+2,j)+img(i-2,j));
- h(i,j)=h(i,j)+4*(img(i-1,j-1)+img(i-1,j+1)+img(i+1,j-1)+img(i+1,j+1));
- h(i,j)=h(i,j)+6*(img(i,j-1)+img(i,j+1)+img(i+1,j)+img(i-1,j));
- h(i,j)=h(i,j)+7*img(i,j);
- h(i,j)=h(i,j)/79;
- end
- end
- end
- imshow(h,[]);
三、初步邊緣檢測
對高斯平滑後的影象進行sobel邊緣檢測。這裡需要求橫的和豎的還有聯合的,所以一共三個需要sobel邊緣檢測影象。
程式碼:
- img=imread('lenna.bmp');
- [m n]=size(img);
- img=double(img);
- w=fspecial('gaussian',[5 5]);
- img=imfilter(img,w,'replicate');
- figure;
- %%sobel邊緣檢測
- w=fspecial('sobel');
- img_w=imfilter(img,w,'replicate'); %求橫邊緣
- w=w'; %轉置矩陣
- img_h=imfilter(img,w,'replicate'); %求豎邊緣
- img=sqrt(img_w.^2+img_h.^2); %平方和再開方 .^表示每一位都和自己乘,不清楚的自己百度
- figure;
- imshow(uint8(img))
四、非極大值抑制以及雙閾值檢測邊緣
什麼是非極大值抑制?
非極大值抑制是在梯度方向上的極大值。
程式碼:(覺得看看註釋就差不多明白思路了)
- img=imread('lenna.bmp');
- imshow(img);
- [m n]=size(img);
- img=double(img);
- %%canny邊緣檢測的前兩步相對不復雜,所以我就直接呼叫系統函數了
- %%高斯濾波
- w=fspecial('gaussian',[5 5]);
- img=imfilter(img,w,'replicate');
- figure;
- imshow(uint8(img))
- %%sobel邊緣檢測
- w=fspecial('sobel');
- img_w=imfilter(img,w,'replicate'); %求橫邊緣
- w=w';
- img_h=imfilter(img,w,'replicate'); %求豎邊緣
- img=sqrt(img_w.^2+img_h.^2); %注意這裡不是簡單的求平均,而是平方和在開方。我曾經好長一段時間都搞錯了
- figure;
- imshow(uint8(img))
- %%下面是非極大抑制
- new_edge=zeros(m,n);
- for i=2:m-1
- for j=2:n-1
- Mx=img_w(i,j);
- My=img_h(i,j);
- if My~=0
- o=atan(Mx/My); %邊緣的法線弧度
- elseif My==0 && Mx>0
- o=pi/2;
- else
- o=-pi/2;
- end
- %Mx處用My和img進行插值
- adds=get_coords(o);%邊緣畫素法線一側求得的兩點座標,插值需要
- M1=My*img(i+adds(2),j+adds(1))+(Mx-My)*img(i+adds(4),j+adds(3)); %插值後得到的畫素,用此畫素和當前畫素比較
- adds=get_coords(o+pi);%邊緣法線另一側求得的兩點座標,插值需要
- M2=My*img(i+adds(2),j+adds(1))+(Mx-My)*img(i+adds(4),j+adds(3));%另一側插值得到的畫素,同樣和當前畫素比較
- isbigger=(Mx*img(i,j)>M1)*(Mx*img(i,j)>=M2)+(Mx*img(i,j)<M1)*(Mx*img(i,j)<=M2); %如果當前點比兩邊點都大置1
- if isbigger
- new_edge(i,j)=img(i,j);
- end
- end
- end
- figure;
- imshow(uint8(new_edge))
- %%下面是滯後閾值處理
- up=120; %上閾值
- low=100; %下閾值
- set(0,'RecursionLimit',10000); %設定最大遞迴深度
- for i=1:m
- for j=1:n
- if new_edge(i,j)>up &&new_edge(i,j)~=255 %判斷上閾值
- new_edge(i,j)=255;
- new_edge=connect(new_edge,i,j,low);
- end
- end
- end
- figure;
- imshow(new_edge==255)
get_coords.m:- function re=get_coords(angle) %angle是邊緣法線角度,返回法線前後兩點
- sigma=0.000000001;
- x1=ceil(cos(angle+pi/8)*sqrt(2)-0.5-sigma);
- y1=ceil(-sin(angle-pi/8)*sqrt(2)-0.5-sigma);
- x2=ceil(cos(angle-pi/8)*sqrt(2)-0.5-sigma);
- y2=ceil(-sin(angle-pi/8)*sqrt(2)-0.5-sigma);
- re=[x1 y1 x2 y2];
- end
connect.m:
- functionnedge=connect(nedge,y,x,low) %種子定位後的連通分析
- neighbour=[-1 -1;-1 0;-1 1;0 -1;0 1;1 -1;1 0;1 1]; %八連通搜尋
- [m n]=size(nedge);
- for k=1:8
- yy=y+neighbour(k,1);
- xx=x+neighbour(k,2);
- if yy>=1 &&yy<=m &&xx>=1 && xx<=n
- if nedge(yy,xx)>=low && nedge(yy,xx)~=255 %判斷下閾值
- nedge(yy,xx)=255;
- nedge=connect(nedge,yy,xx,low);
- end
- end
- end
- end
參考資料:
1.http://www.cnblogs.com/tiandsp/archive/2012/12/13/2817240.html
2.南京大學軟體學院數字影象處理課程PPT
一些更詳細的原理性的東西可以參考:
http://blog.csdn.net/likezhaobin/article/details/6892176