Kirsch边缘检测原理

Kirsch算子由8个3×3窗口模扳组成、每个模板分别代表一个特定的检测方向，其模板算子如图2.3.1所示。

  在进行边缘检测时，把M0-M7所表示的边缘模板(加权矩阵)分别与图像中的一个3×3区域相乘，选取输出值为最大的模板。然后，把这一最大输出值作为中央像素点上的边缘强度，把取得最大值的边缘模板Mk的方向k(k的取值如图2.3.2所示)作为其边缘方向。假设图像中一点P(i，j)及其八邻域的灰度如图3.3.3所示，并设Qk(k＝0，1，…，7)为图像经过kirsch算子第k个模板处理后得到的k方向上的边缘强度，则P(i,j)的边缘强度为s(i,j)=max|qk|(k=0,1，…，7)，而相应的边缘方向D(i，j)={k|qk为最大值}

%读出要处理的图象
clear
clc
close all
bw=imread('e:\11.jpg');

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%对图象进行预处理

%画出原始图象
bw1=rgb2gray(bw);
figure(1)
imshow(bw1)
title('原始图象')

%对图象进行均值滤波处理
bw2=filter2(fspecial('average',3),bw1);
figure(2)
imshow(bw2)
title('均值滤波')

%对图象进行高斯滤波处理
bw3=filter2(fspecial('gaussian'),bw2);
figure(3)
imshow(bw3)
title('高斯滤波')

%利用小波变换对图象进行降噪处理
[thr,sorh,keepapp]=ddencmp('den','wv',bw3);     %获得除噪的缺省参数
bw4=wdencmp('gbl',bw3,'sym4',2,thr,sorh,keepapp);%图象进行降噪处理
figure(4)
imshow(bw4)

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%提取图象边缘
t=1200 ;                        %设定阈值
colormap(gray(256));          %设定调色板
bw5=double(bw4);             %把图变为十进制数
[m,n]=size(bw5);               %得到图像的大小（长和宽）
g=zeros(m,n);              %定义一个大小为S的零矩阵
%利用Sobel算子进行边缘提取
for i=2:m-1
  for j=2:n-1
       d1 =(5*bw5(i-1,j-1)+5*bw5(i-1,j)+5*bw5(i-1,j+1)-3*bw5(i,j-1)-3*bw5(i,j+1)-3*bw5(i+1,j-1)-3*bw5(i+1,j)-3*bw5(i+1,j+1))^2;
       d2 =((-3)*bw5(i-1,j-1)+5*bw5(i-1,j)+5*bw5(i-1,j+1)-3*bw5(i,j-1)+5*bw5(i,j+1)-3*bw5(i+1,j-1)-3*bw5(i+1,j)-3*bw5(i+1,j+1))^2;
       d3 =((-3)*bw5(i-1,j-1)-3*bw5(i-1,j)+5*bw5(i-1,j+1)-3*bw5(i,j-1)+5*bw5(i,j+1)-3*bw5(i+1,j-1)-3*bw5(i+1,j)+5*bw5(i+1,j+1))^2;
       d4 =((-3)*bw5(i-1,j-1)-3*bw5(i-1,j)-3*bw5(i-1,j+1)-3*bw5(i,j-1)+5*bw5(i,j+1)-3*bw5(i+1,j-1)+5*bw5(i+1,j)+5*bw5(i+1,j+1))^2;
       d5 =((-3)*bw5(i-1,j-1)-3*bw5(i-1,j)-3*bw5(i-1,j+1)-3*bw5(i,j-1)-3*bw5(i,j+1)+5*bw5(i+1,j-1)+5*bw5(i+1,j)+5*bw5(i+1,j+1))^2;
       d6 =((-3)*bw5(i-1,j-1)-3*bw5(i-1,j)-3*bw5(i-1,j+1)+5*bw5(i,j-1)-3*bw5(i,j+1)+5*bw5(i+1,j-1)+5*bw5(i+1,j)-3*bw5(i+1,j+1))^2;
       d7 =(5*bw5(i-1,j-1)-3*bw5(i-1,j)-3*bw5(i-1,j+1)+5*bw5(i,j-1)-3*bw5(i,j+1)+5*bw5(i+1,j-1)-3*bw5(i+1,j)-3*bw5(i+1,j+1))^2;
       d8 =(5*bw5(i-1,j-1)+5*bw5(i-1,j)-3*bw5(i-1,j+1)+5*bw5(i,j-1)-3*bw5(i,j+1)-3*bw5(i+1,j-1)-3*bw5(i+1,j)-3*bw5(i+1,j+1))^2; 

       g(i,j)=round(sqrt(d1+d2+d3+d4+d5+d6+d7+d8)); %梯度模取整
    end
end
for i=1:m
  for j=1:n
  if g(i,j)>t
  bw5(i,j)=255;              %将梯度值与阈值比较 ，大于T则把图像的灰度变为255，小于T则把图像的灰度变为0
  else
  bw5(i,j)=0;
    end
  end
end
%显示边缘提取后的图象
figure(5)
imshow(bw5)
title('kirsch边缘检测')

clc
clear all close all
A = imread('lena.jpg');
mask1=[-3,-3,-3;-3,0,5;-3,5,5];  % 建立方向模板
mask2=[-3,-3,5;-3,0,5;-3,-3,5];
mask3=[-3,5,5;-3,0,5;-3,-3,-3];
mask4=[-3,-3,-3;-3,0,-3;5,5,5];
mask5=[5,5,5;-3,0,-3;-3,-3,-3];
mask6=[-3,-3,-3;5,0,-3;5,5,-3];
mask7=[5,-3,-3;5,0,-3;5,-3,-3];
mask8=[5,5,-3;5,0,-3;-3,-3,-3];
B=mat2gray(A);
subplot(121);imshow(B);title('原图');
I = im2double(A);  % 将数据图像转化为双精度
d1 = imfilter(I, mask1);  % 计算8个领域的灰度变化
d2 = imfilter(I, mask2);
d3 = imfilter(I, mask3);
d4 = imfilter(I, mask4);
d5 = imfilter(I, mask5);
d6 = imfilter(I, mask6);
d7 = imfilter(I, mask7);
d8 = imfilter(I, mask8);
dd = max(abs(d1),abs(d2));  % 取差值变化最大的元素组成灰度变化矩阵
dd = max(dd,abs(d3));
dd = max(dd,abs(d4));
dd = max(dd,abs(d5));
dd = max(dd,abs(d6));
dd = max(dd,abs(d7));
dd = max(dd,abs(d8));
%grad = mat2gray(dd);  % 将灰度变化矩阵转化为灰度图像
%level = graythresh(grad);  % 计算灰度阈值
BW = im2bw(grad,0.03);  % 用阈值分割梯度图像
subplot(122); imshow(BW);title('Kirsch算子的处理结果')  % 显示分割后的图像，即边缘图像 title('Kirsch')