同上一篇，还是为了体现FPGA的强大功能，实现实时的边缘检测能力！这一部分简单的可以用Sobel实现，如果想做的好，可以用高斯+Canny来实现，总是，只要你想做，FPGA没有什么做不到的，只要你静得下心来！！！！以下仅仅介绍Soebl算法的Matlab算法实现，希望大家各自努力！

1. Sobel算子的边缘检测实现

1.1. 边缘检测概念

所谓边缘是指其周围像素灰度急剧变化的那些象素的集合，它是图像最基本的特征。边缘存在于目标、背景和区域之间，所以，它是图像分割所依赖的最重要的依据。由于边缘是位置的标志，对灰度的变化不敏感，,因此，边缘也是图像匹配的重要的特征。

边缘检测和区域划分是图像分割的两种不同的方法，二者具有相互补充的特点。在边缘检测中，是提取图像中不连续部分的特征，根据闭合的边缘确定区域。而在区 域划分中，是把图像分割成特征相同的区域，区域之间的边界就是边缘。由于边缘检测方法不需要将图像逐个像素地分割，因此更适合大图像的分割。
边缘大致可以分为两种，一种是阶跃状边缘，边缘两边像素的灰度值明显不同；另一种为屋顶状边缘，边缘处于灰度值由小到大再到小的变化转折点处。
边缘检测的主要工具是边缘检测模板。

边缘检测的有很多，典型的有索贝尔算子，普里维特算子，罗伯茨交叉边缘检测等边缘检测技术，在Matlab中有现成的IPT函数，提供边缘检测，如下，Sobel边缘检测：

IMG1 = imread('D:\Matlab_Project\BMP\Lenna.jpg');    % 读取RGB文件，Lanna PlayBoy

subplot(1,3,1)

imshow(IMG1);

title('原图像');

IMG1 = rgb2gray(IMG1);

[m,n] = size(IMG1); %用Sobel微分算子进行边缘检测

IMG2 = edge(IMG1,'sobel');

subplot(1,3,2);

imshow(IMG2);

title('Sobel边缘检测得到的图像');

但效果不佳，灵活性也不高。此处我主要介绍Sobel算子的使用：

1.2. Sobel算法实现

1.2.1. Sobel算法分析

索贝尔算子（Sobel operator）主要用作边缘检测，在技术上，它是一离散性差分算子，用来运算图像亮度函数的灰度之近似值。在图像的任何一点使用此算子，将会产生对应的灰度矢量或是其法矢量。

Sobel卷积因子为：

该算子包含两组3x3的矩阵，分别为横向及纵向，将之与图像作平面卷积，即可分别得出横向及纵向的亮度差分近似值。如果以A代表原始图像，Gx及Gy分别代表经横向及纵向边缘检测的图像灰度值，其公式如下：

图像的每一个像素的横向及纵向灰度值通过以下公式结合，来计算该点灰度的大小：

通常，为了提高效率 使用不开平方的近似值,但这样做会损失精度

如果梯度G大于某一阀值 则认为该点(x,y)为边缘点。

        if(temp3 > THRESHOLD)

            IMG_Sobel(i,j) = 0;     %Black

        else

            IMG_Sobel(i,j) = 255;   %White

        end

然后可用以下公式计算梯度方向(当然只要检测边缘，则不用计算方向)：

1.2.2. Sobel算子Matlab算法的实现

Sobel算子根据像素点上下、左右邻点灰度加权差，在边缘处达到极值这一现象检测边缘。对噪声具有平滑作用，提供较为精确的边缘方向信息，边缘定位精度不够高。当对精度要求不是很高时，是一种较为常用的边缘检测方法。
Soble算子操作的是灰度图像，因此必须首先对图像进行灰度化。这一部分在第一章中详细介绍过，不再做具体分析，一下算法直接应用灰度图像。

（1）Sobel算子卷积

% -----------------------------------------------------------------------

%         Gx                  Gy                  Pixel

% [   -1  0   +1  ]   [   +1  +2   +1 ]     [   P1  P2   P3 ]

% [   -2  0   +2  ]   [   0   0    0  ]     [   P4  P5   P6 ]

% [   -1  0   +1  ]   [   -1  -2   -1 ]     [   P7  P8   P9 ]

Sobel_X = [-1, 0, 1, -2, 0, 2, -1, 0, 1];   % Mask x

Sobel_Y = [1, 2, 1, 0, 0, 0, -1, -2, -1];   % Mask y

（2）点对点卷积运算实现梯度的计算

IMG_Gray = double(IMG2);    %将图片转换为双精度类型

IMG_Sobel = true(h,w);     %新建一个二值矩阵

THRESHOLD =90;

for i = 2 : h-1     %舍弃了边缘信息

    for j = 2 : w-1

        temp1 = Sobel_X(1) * IMG_Gray(i-1,j-1)  + Sobel_X(2) * IMG_Gray(i-1,j)  + Sobel_X(3) * IMG_Gray(i-1,j+1) +...

                Sobel_X(4) * IMG_Gray(i,j-1)    + Sobel_X(5) * IMG_Gray(i,j)    + Sobel_X(6) * IMG_Gray(i,j+1) +...

                Sobel_X(7) * IMG_Gray(i+1,j-1)  + Sobel_X(8) * IMG_Gray(i+1,j)  + Sobel_X(9) * IMG_Gray(i+1,j+1);

        temp2 = Sobel_Y(1) * IMG_Gray(i-1,j-1)  + Sobel_Y(2) * IMG_Gray(i-1,j)  + Sobel_Y(3) * IMG_Gray(i-1,j+1) +...

                Sobel_Y(4) * IMG_Gray(i,j-1)    + Sobel_Y(5) * IMG_Gray(i,j)    + Sobel_Y(6) * IMG_Gray(i,j+1) +...

                Sobel_Y(7) * IMG_Gray(i+1,j-1)  + Sobel_Y(8) * IMG_Gray(i+1,j)  + Sobel_Y(9) * IMG_Gray(i+1,j+1);

        temp3 = sqrt(temp1^2 + temp2^2);

        %temp3 = abs(temp1) + abs(temp2);   %just for speed

        if(temp3 > THRESHOLD)

            IMG_Sobel(i,j) = 0; %Black

        else

            IMG_Sobel(i,j) = 1; %White

        end

    end

end

（3）THRESHOLD为手动设定，此处为90。根据图像的质量来调节，针对于边缘检测而言，不需要进行梯度方向计算，因此直接进行梯度计算，与阀值对比，既可以实现边缘检测。

1.2.3. 效果图

THRESHOLD = 90下的图像

最后晒个图吧，呵呵FPGA移植的，绝对不是抄袭 “无双OO”，他那个写的实在是太恶心了（无意冒犯前辈，不好意思）。。。。。。

实际上简化后的算法如下：

//---------------------------------------
//Caculate vertical Grade with |abs|
reg	[9:0]	Gy_temp1;	//postive result
reg	[9:0]	Gy_temp2;	//negetive result
reg	[9:0]	Gy_data;	//Vertical grade data
always@(posedge clk or negedge rst_n)
begin
	if(!rst_n)
		begin
		Gy_temp1 <= 0;
		Gy_temp2 <= 0;
		Gy_data <= 0;
		end
	else
		begin
		Gy_temp1 <= image_p11 + (image_p12 << 1) + image_p13;	//postive result
		Gy_temp2 <= image_p31 + (image_p32 << 1) + image_p33;	//negetive result
		Gy_data <= (Gy_temp1 >= Gy_temp2) ? Gy_temp1 - Gy_temp2 : Gy_temp2 - Gy_temp1;
		end
end

//---------------------------------------
//Caculate the square of distance = (Gx^2 + Gy^2)
reg	[20:0]	Gxy_square;
always@(posedge clk or negedge rst_n)
begin
	if(!rst_n)
		Gxy_square <= 0;
	else
		Gxy_square <= Gx_data * Gx_data + Gy_data * Gy_data;
end

//---------------------------------------
//Caculate the distance of P5 = (Gx^2 + Gy^2)^0.5
wire	[10:0]	Dim;
SQRT	u_SQRT
(
	.radical	(Gxy_square),
	.q			(Dim),
	.remainder	()
);