C#图像处理程序实现--Canny实现 笔记整理（仅作相关知识点了解）

于Canny算法的笔记
以下内容参考本篇博文，看了许多的介绍，这篇个人感觉诗写得比较详细的，本篇参考这个并结合C#的实现过程进行整理
l
Canny边缘检测算法可以分为以下5个步骤：
1)使用高斯滤波器，以平滑图像，滤除噪声。
2)计算图像中每个像素点的梯度强度和方向。
3)应用非极大值（Non-Maximum Suppression）抑制，以消除边缘检测带来的杂散响应。
4)应用双阈值（Double-Threshold）检测来确定真实的和潜在的边缘。
5)通过抑制孤立的弱边缘最终完成边缘检测。
下面详细介绍每一步的实现思路。
3.1高斯平滑滤波
为了尽可能减少噪声对边缘检测结果的影响，所以必须滤除噪声以防止由噪声引起的错误检测。为了平滑图像，使用高斯滤波器与图像进行卷积，该步骤将平滑图像，以减少边缘检测器上明显的噪声影响。大小为(2k+1)x(2k+1)的高斯滤波器核的生成方程式由下式给出：

若图像中一个3x3的窗口为A，要滤波的像素点为e，则经过高斯滤波之后，像素点e的亮度值为：

其中" * "为卷积符号，sum表示矩阵中所有元素相加求和。关于卷积的实现过程可以参考这个
重要的是需要理解，高斯卷积核大小的选择将影响Canny检测器的性能。尺寸越大，检测器对噪声的敏感度越低，但是边缘检测的定位误差也将略有增加。一般5x5是一个比较不错的trade off。
3.2 计算梯度强度和方向
图像中的边缘可以指向各个方向，因此Canny算法使用四个算子来检测图像中的水平、垂直和对角边缘。边缘检测的算子（如Roberts，Prewitt，Sobel等）返回水平Gx和垂直Gy方向的一阶导数值，由此便可以确定像素点的梯度G和方向theta 。

其中G为梯度强度， theta表示梯度方向，arctan为反正切函数。下面以Sobel算子为例讲述如何计算梯度强度和方向。
x和y方向的Sobel算子分别为：

其中Sx表示x方向的Sobel算子，用于检测y方向的边缘； Sy表示y方向的Sobel算子，用于检测x方向的边缘（边缘方向和梯度方向垂直）。在直角坐标系中，Sobel算子的方向如下图所示。

若图像中一个3x3的窗口为A，要计算梯度的像素点为e，则和Sobel算子进行卷积之后，像素点e在x和y方向的梯度值分别为：

其中*为卷积符号，sum表示矩阵中所有元素相加求和。根据公式（3-2）便可以计算出像素点e的梯度和方向。参考我这里之前写的实现Sobel的过程
这里单独拎出来看一下，上面的链接讲了实现梯度的算法，这里需要单独提取出X和Y方向的梯度，其实就是用卷积的方式，先将Sobel的X方向卷积得出Gx，然后在得出Gy。很简单的，实现过程记得都把每个方法写成一个单独额函数来调用。会极大地方便代码管理。

在这里插入代码片
​static byte[,] XgradientOfSobel(Bitmap bmpToLaplace)
        {
            int[,] filterMatrix = new int[3, 3] { { -1, 0, 1 }, { -2, 0, 2 }, { -1, 0, 1 } };//in X direction
            if ((filterMatrix.GetLength(0) & 1) == 0 || (filterMatrix.GetLength(1) & 1) == 0 || filterMatrix.GetLength(0) != filterMatrix.GetLength(1))
                return null;
            int height = bmpToLaplace.Height;
            int width = bmpToLaplace.Width;
            int filterOffset = (filterMatrix.GetLength(1) - 1) / 2;
            BitmapData bmpData = bmpToLaplace.LockBits(new Rectangle(0, 0, bmpToLaplace.Width, bmpToLaplace.Height), ImageLockMode.ReadWrite, bmpToLaplace.PixelFormat);
            int channels = Image.GetPixelFormatSize(bmpToLaplace.PixelFormat) / 8;
            int stride = bmpData.Stride;

            int offset = stride - width * channels;
            byte[,] greyPixelArray = new byte[height, width];

            GCHandle handle = GCHandle.Alloc(greyPixelArray, GCHandleType.Pinned);
            IntPtr destination = handle.AddrOfPinnedObject();
            unsafe
            {

                byte* src = (byte*)bmpData.Scan0.ToPointer() + filterOffset * stride;
                byte* dst = (byte*)destination + filterOffset * width;
                for (int i = filterOffset; i < height - filterOffset; i++)
                {
                    src += filterOffset * channels;
                    dst += filterOffset;
                    for (int j = filterOffset; j < width - filterOffset; j++, src += channels, dst++)
                    {
                        int b = 0, g = 0, r = 0, bc = 1, gc = 1, rc = 1;
                        if (src[0] == 0)
                            bc = -8;
                        if (src[1] == 0)
                            gc = -8;
                        if (src[2] == 0)
                            rc = -8;

                        for (int fi = -filterOffset; fi <= filterOffset; fi++)
                            for (int fj = -filterOffset; fj <= filterOffset; fj++)
                            {
                                int cpos = fi * stride + fj * channels;
                                int fv = filterMatrix[fi + filterOffset, fj + filterOffset];
                                b += src[cpos] * fv;
                                g += src[cpos + 1] * fv;
                                r += src[cpos + 2] * fv;
                            }
                        int blue = b / bc, green = g / gc, red = r / rc;
                        blue = blue > 255 ? 255 : ((blue < 0) ? 0 : blue);
                        green = green > 255 ? 255 : ((green < 0) ? 0 : green);
                        red = red > 255 ? 255 : ((red < 0) ? 0 : red);
                        *dst = (byte)(blue * 0.11 + green * 0.59 + red * 0.3);
                    }
                    src += filterOffset * channels + offset;
                    dst += filterOffset;
                }
            }
            bmpToLaplace.UnlockBits(bmpData);
            handle.Free();
            return greyPixelArray;
        }​

3.3 非极大值抑制
讲到这里，先补习一下之前的一些小知识：

导数和偏导数
导数与偏导数本质是一致的，都是当自变量的变化量趋于0时，函数值的变化量与自变量变化量比值的极限。直观地说，偏导数也就是函数在某一点上沿坐标轴正方向的的变化率。
区别在于：
导数，指的是一元函数中，函数y=f(x)在某一点处沿x轴正方向的变化率；
偏导数，指的是多元函数中，函数y=f(x1,x2,…,xn)在某一点处沿某一坐标轴（x1,x2,…,xn）正方向的变化率。
导数与方向导数：
在前面导数和偏导数的定义中，均是沿坐标轴正方向讨论函数的变化率。那么当我们讨论函数沿任意方向的变化率时，也就引出了方向导数的定义，即：某一点在某一趋近方向上的导数值。
通俗的解释是：我们不仅要知道函数在坐标轴正方向上的变化率（即偏导数），而且还要设法求得函数在其他特定方向上的变化率。而方向导数就是函数在其他特定方向上的变化率。
导数与梯度
梯度的提出只为回答一个问题：
　函数在变量空间的某一点处，沿着哪一个方向有最大的变化率？
　梯度定义如下：
　函数在某一点的梯度是这样一个向量，它的方向与取得最大方向导数的方向一致，而它的模为方向导数的最大值。
　这里注意三点：
　1）梯度是一个向量，即有方向有大小；
　2）梯度的方向是最大方向导数的方向；
　3）梯度的值是最大方向导数的值。
导数与向量
向量的定义是有方向（direction）有大小（magnitude）的量。
　从前面的定义可以这样看出，偏导数和方向导数表达的是函数在某一点沿某一方向的变化率，也是具有方向和大小的。因此从这个角度来理解，我们也可以把偏导数和方向导数看作是一个向量，向量的方向就是变化率的方向，向量的模，就是变化率的大小。
　那么沿着这样一种思路，就可以如下理解梯度：
　梯度即函数在某一点最大的方向导数，函数沿梯度方向函数有最大的变化率

极大值抑制是一种边缘稀疏技术，非极大值抑制的作用在于瘦边。对图像进行梯度计算后，仅仅基于梯度值提取的边缘仍然很模糊。对于标准3，对边缘有且应当只有一个准确的响应。而非极大值抑制则可以帮助将局部最大值之外的所有梯度值抑制为0，对梯度图像中每个像素进行非极大值抑制的算法是：

1. 将当前像素的梯度强度与沿正负梯度方向上的两个像素进行比较。
2. 如果当前像素的梯度强度与另外两个像素相比最大，则该像素点保留为边缘点，否则该像素点将被抑制。
   通常为了更加精确的计算，在跨越梯度方向的两个相邻像素之间使用线性插值来得到要比较的像素梯度，现举例如下：
   
   如图3-2所示，将梯度分为8个方向，分别为E、NE、N、NW、W、SW、S、SE，其中0代表00_{45o,1代表450}90o，2代表-900~-45o， 3代表-450~0o。像素点P的梯度方向为theta，则像素点P1和P2的梯度线性插值为：
   
   因此非极大值抑制的伪代码描写如下：
   
   需要注意的是，如何标志方向并不重要，重要的是梯度方向的计算要和梯度算子的选取保持一致。
   这里补习一下正切函数，如下图（来自百度百科）：
   
   以上做些补充：来自于这里
   实际数字图像中的像素点是离散的二维矩阵，所以处在真正中心位置C处的梯度方向两侧的点是不一定存在的，或者说是一个亚像素(sub pixel)点，而这个不存在的点， 以及这个点的梯度值就必须通过对其两侧的点进行插值来得到。
   对于上面的代码，如果|gy|>|gx|，这说明该点的梯度方向更靠近Y轴方向，所以g2和g4则在C的上下，我们可以用下面来说明这两种情况（方向相同和方向不同）：
   
   上图中，C表示中心位置点，斜的直线表示梯度方向（非极大值抑制是在梯度方向上的极大值），左边的一副表示gy与gx的方向相同，而右边的这幅这表示gy与gx的方向相反（注意原点在左上角），而权重则为weight = |gx|/|gy|，因此则根据此种情况，其插值表示则为：

在这里插入代码片
​dTemp1 = weight*g1 + (1-weight)*g2;
dTemp2 = weight*g3 + (1-weight)*g4;	​

同理，我们可以得到|gx|>|gy|的情况，此时说明该点的梯度方向更靠近X轴方向，g2和g4则在水平方向，我们可以用下图来说明该种情况：

上图中，C表示中心位置点，斜的直线表示梯度方向（非极大值抑制是在梯度方向上的极大值），左边的一副表示gy与gx的方向相同，而右边的这幅这表示gy与gx的方向相反（注意原点在左上角），而权重则为weight = |gy|/|gx|，因此则根据此种情况，其插值表示则为：

在这里插入代码片
​dTemp1 = weight*g1 + (1-weight)*g2;
dTemp2 = weight*g3 + (1-weight)*g4;	​

通过上面的分析，我们可以了解Canny算子中的非极大值抑制之前的准备工作，也即进行必要的插值。插值的原因再啰嗦下：①由于在Canny算子中采用的简化方法来进行边缘方向的确定，自然图像中边缘梯度方向的不一定沿着该四个方向，因此为了找出在一个像素点上最能吻合其所在梯度方向的两侧的像素值，就必须进行必要的插值；② 也由于实际数字图像中的像素点是离散的二维矩阵，处在真正中心位置C处的梯度方向两侧的点是不一定存在的，或者说是一个亚像素(sub pixel)点，而这个不存在的点， 以及这个点的梯度值就必须通过对其两侧的点进行插值来得到。

一定要注意插值法在不同角度下，权重的计算是不一样)

3.4 双阈值检测
在施加非极大值抑制之后，剩余的像素可以更准确地表示图像中的实际边缘。然而，仍然存在由于噪声和颜色变化引起的一些边缘像素。为了解决这些杂散响应，必须用弱梯度值过滤边缘像素，并保留具有高梯度值的边缘像素，可以通过选择高低阈值来实现。如果边缘像素的梯度值高于高阈值，则将其标记为强边缘像素；如果边缘像素的梯度值小于高阈值并且大于低阈值，则将其标记为弱边缘像素；如果边缘像素的梯度值小于低阈值，则会被抑制。阈值的选择取决于给定输入图像的内容。
双阈值检测的伪代码描写如下：

3.5 抑制孤立低阈值点
到目前为止，被划分为强边缘的像素点已经被确定为边缘，因为它们是从图像中的真实边缘中提取出来的。然而，对于弱边缘像素，将会有一些争论，因为这些像素可以从真实边缘提取也可以是因噪声或颜色变化引起的。为了获得准确的结果，应该抑制由后者引起的弱边缘。通常，由真实边缘引起的弱边缘像素将连接到强边缘像素，而噪声响应未连接。为了跟踪边缘连接，通过查看弱边缘像素及其8个邻域像素，只要其中一个为强边缘像素，则该弱边缘点就可以保留为真实的边缘。
抑制孤立边缘点的伪代码描述如下：

4 总结
通过以上5个步骤即可完成基于Canny算法的边缘提取