（亲测可用）基于matlab的用自写函数来实现图像的灰度处理sobel canny算子边缘检测

I=imread(‘E:\matlab\toolbox\images\imdata\tupian4.jpg’); %图片的位置是自己设置的

subplot(2,5,1);

imshow(I);

title(‘原始图像’);

G=I(:,:,2);

subplot(2,5,2);

imshow(G);

title(‘灰度处理后图像’); %彩色图像需要灰度处理，黑白图像则不需要

c=avefilt(G,3); %对G用3x3的模板进行滤波

subplot(2,5,3);

imshow(uint8( c ));

title(‘均值滤波后的图像’);

hx=[-1 -2 -1;0 0 0 ;1 2 1]; %生产sobel垂直梯度模板，这里用的是自定义的模板

hy=hx’; %生产sobel水平梯度模板

gradx=filter2(hx,c,‘same’);

gradx=abs(gradx); %计算图像的sobel垂直梯度

grady=filter2(hy,c,‘same’);

grady=abs(grady); %计算图像的sobel水平梯度

grad=gradx+grady; %得到图像的sobel梯度subplot(2,5,4);imshow(grad,[ ]);title('图像的sobel梯度');

%canny算子--------------------------------------------------------

[m,n]=size(G);

src_gaussion=my_gaussian(G,3,9); %高斯滤波

w=fspecial(‘sobel’); %这里的sobel算子是与上面不同的算法

src_sobel_x=imfilter(src_gaussion,w,‘replicate’); %求横边缘

w=w’;

src_sobel_y=imfilter(src_gaussion,w,‘replicate’); %求竖边缘

src_sobel_x=double(src_sobel_x);

src_sobel_y=double(src_sobel_y);

src_sobel=sqrt(src_sobel_x.2+src_sobel_y.2); %注意平方和在开方。

% figure;

subplot(2,5,5);

imshow(uint8(src_sobel_x))

title(‘这是sobel进行梯度处理后的图像’);

new_angel_img=zeros(m,n);

for i=1:m

for j=1:n

Mx=src_sobel_x(i,j);

My=src_sobel_y(i,j);

if My~=0 o=atan(Mx/My);%边缘的法线弧度 elseif My==0 && Mx>0 o=pi/2; else o=-pi/2; end if( ( o>=( -22.5/180*pi ) && o<( 22.5/180*pi ) ) || ( o>=( 157.5/180*pi ) && o<=pi ) || ( o<=( -157.5/180*pi ) && o>=-pi ) ) d1=0; elseif ( ( o>= ( 22.5/180*pi ) && o<( 67.5/180*pi ) ) || ( o>=( -157.5/180*pi ) && o<(-112.5/180*pi ) ) ) d1=-45; elseif ( ( o>=(67.5/180*pi) && o<(112.5/180*pi) ) || ( o>=(-112.5/180*pi) && o<-67.5/180*pi ) )d1=90; else d1=45; end new_angel_img(i,j)=d1; end

end

canny_dst=src_sobel;

%%%%%%%%%%%%%%%%%%非最大值抑制%%%%%%%%%%%%%%%%%%%%%%%5

for i=3:m-2

for j=3:n-2

if( new_angel_img(i,j)0 )

A=[canny_dst(i-1,j),canny_dst(i+1,j)];

if( abs(canny_dst(i,j))< max( abs(A) ) )

canny_dst(i,j)=0;

end

elseif( new_angel_img(i,j)-45 )

A=[canny_dst(i-1,j-1),canny_dst(i+1,j+1)];

if( abs(canny_dst(i,j))< max( abs(A) ))

canny_dst(i,j)=0;

end

elseif( new_angel_img(i,j)==90 )

A=[canny_dst(i,j-1),canny_dst(i,j+1)];

if( abs(canny_dst(i,j))< max( abs(A) ) )

canny_dst(i,j)=0;

end

elseif( new_angel_img(i,j)==45 )

A=[canny_dst(i-1,j+1),canny_dst(i+1,j-1)];

if( abs(canny_dst(i,j))< max( abs(A) ) )

canny_dst(i,j)=0;

end

end

end

end

% figure(2);

subplot(2,5,6);

imshow(canny_dst);

title(‘这是非最大值抑制11’);

gh=zeros(m,n);

gl=zeros(m,n);

gl1=zeros(m,n);

for i=1:m

for j=1:n

if( canny_dst(i,j)>=100 )

gh(i,j)=1;

else

gh(i,j)=0;

end

if( canny_dst(i,j)>=30)

gl(i,j)=1;

else

gl(i,j)=0;

end

end

end

gl1=gl;

gl1=logical(gl1);

gh=logical(gh);

% figure(3);

subplot(2,5,7);

imshow(gl1);

title(‘这是非最大值抑制22’);

% figure(4);

subplot(2,5,8);

imshow(gh);

title(‘这是递归连接后的33’);

%%%%%%%%%%%%%%%%%%%%%%%%递归连接%%%%%%%%%%%%%%%%%%%

up=100; %上阈值

low=33; %下阈值

set(0,‘RecursionLimit’,10000); %设置最大递归深度

for i=1:m

for j=1:n

if canny_dst(i,j)>up &&canny_dst(i,j)~=255 %判断上阈值

canny_dst(i,j)=255;

canny_dst=connect(canny_dst,i,j,low);

end

end

end

% figure(5);

subplot(2,5,9);

imshow(canny_dst==255);

title(‘这是递归连接后的’);

BW2 = edge(src_gaussion,‘canny’);

% figure(6);

subplot(2,5,10);

imshow(BW2);

title(‘这是递归连接后的2’);

%所需函数————————————————————————————————————

function d=my_gaussian(src,n,k)

n1=(n+1)/2; %%%%%%%%%%%%%计算高斯滤波模版中心

[m,l]=size(src);

b=zeros(n,n); %%%%%%%%%%%%%空模版

I=zeros(m,l);

d=zeros(m,l);

img=zeros(m,l);

%%%%%%%%%%%%%%%%%计算高斯滤波模版的各个系数

for i=1:n

for j=1:n

b(i,j)=(exp(-((i-n1)^2 +(j-n1)^2)/(2k)))/(2pi*k);

end

end

b = b/sum(b( : )); %%%%%%%%%%%%%%%%%%%%%%对系数进行归一化 （这个地方必须归一化，不然整体灰度会底） 删除括号里面的空格

I=double(src);

img=conv2(I,b,‘same’); %%%%%%%%%%%%图像与高斯滤波模版卷积，输出与输入图像大小相同的图像

d=uint8(img); %%%%%%%%%%转为无符号整型，用来显示

end

function nedge=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

function d=avefilt(x,n)

a(1:n,1:n)=1; %a即n×n模板,元素全是1

p=size(x); %输入图像是p×q的,且p>n,q>n

x1=double(x);

x2=x1;

%A(a:b,c:d)表示A矩阵的第a到b行,第c到d列的所有元素

for i=1:p(1)-n+1

for j=1:p(2)-n+1

c=x1(i:i+(n-1),j:j+(n-1)).a; %取出x1中从(i,j)开始的n行n列元素与模板相乘

s=sum(sum©); %求c矩阵(即模板)中各元素之和

x2(i+(n-1)/2,j+(n-1)/2)=s/(nn); %将模板各元素的均值赋给模板中心位置的元素

end

end

% %未被赋值的元素取原值

d=uint8(x2);

end