Canny边缘检测算法的实现

图像边缘信息主要集中在高频段，一般说图像锐化或检测边缘，实质就是高频滤波。咱们知道微分运算是求信号的变化率，具备增强高频份量的做用。在空域运算中来讲，对图像的锐化就是计算微分。因为数字图像的离散信号，微分运算就变成计算差分或梯度。图像处理中有多种边缘检测（梯度）算子，经常使用的包括普通一阶差分，Robert算子（交叉差分），Sobel算子等等，是基于寻找梯度强度。拉普拉斯算子（二阶差分）是基于过零点检测。经过计算梯度，设置阀值，获得边缘图像。

Canny边缘检测算子是一种多级检测算法。1986年由John F. Canny提出，同时提出了边缘检测的三大准则：

1. 低错误率的边缘检测：检测算法应该精确地找到图像中的尽量多的边缘，尽量的减小漏检和误检。
2. 最优定位：检测的边缘点应该精确地定位于边缘的中心。
3. 图像中的任意边缘应该只被标记一次，同时图像噪声不该产生伪边缘。

Canny算法出现之后一直是做为一种标准的边缘检测算法，此后也出现了各类基于Canny算法的改进算法。时至今日，Canny算法及其各类变种依旧是一种优秀的边缘检测算法。并且除非前提条件很适合，你很难找到一种边缘检测算子能显著地比Canny算子作的更好。

关于各类差分算子，还有Canny算子的简单介绍，这里就不罗嗦了，网上均可以找获得。直接进入Canny算法的实现。Canny算法分为几步。

1. 高斯模糊。

这一步很简单，相似于LoG算子（Laplacian of Gaussian）做高斯模糊同样，主要做用就是去除噪声。由于噪声也集中于高频信号，很容易被识别为伪边缘。应用高斯模糊去除噪声，下降伪边缘的识别。可是因为图像边缘信息也是高频信号，高斯模糊的半径选择很重要，过大的半径很容易让一些弱边缘检测不到。

  

             Lena原图                                       Lena高斯模糊，半径2

2. 计算梯度幅值和方向。

图像的边缘能够指向不一样方向，所以经典Canny算法用了四个梯度算子来分别计算水平，垂直和对角线方向的梯度。可是一般都不用四个梯度算子来分别计算四个方向。经常使用的边缘差分算子（如Rober，Prewitt，Sobel）计算水平和垂直方向的差分Gx和Gy。这样就能够以下计算梯度模和方向：

梯度角度θ范围从弧度-π到π，而后把它近似到四个方向，分别表明水平，垂直和两个对角线方向（0°,45°,90°,135°）。能够以±iπ/8（i=1,3,5,7）分割，落在每一个区域的梯度角给一个特定值，表明四个方向之一。

这里我选择Sobel算子计算梯度。Sobel算法很简单，处处均可以找到，就不列出代码来了。相对于其余边缘算子，Sobel算子得出来的边缘粗大明亮。

 

 

 

下图是对上面半径2的高斯模糊图像L通道（HSL）应用Sobel算子的梯度模图，没有施加任何阀值。

           Sobel算子，无阀值

3. 非最大值抑制。

非最大值抑制是一种边缘细化方法。一般得出来的梯度边缘不止一个像素宽，而是多个像素宽。就像咱们所说Sobel算子得出来的边缘粗大而明亮，从上面Lena图的Sobel结果能够看得出来。所以这样的梯度图仍是很“模糊”。而准则3要求，边缘只有一个精确的点宽度。非最大值抑制能帮助保留局部最大梯度而抑制全部其余梯度值。这意味着只保留了梯度变化中最锐利的位置。算法以下：

1. 比较当前点的梯度强度和正负梯度方向点的梯度强度。
2. 若是当前点的梯度强度和同方向的其余点的梯度强度相比较是最大，保留其值。不然抑制，即设为0。好比当前点的方向指向正上方90°方向，那它须要和垂直方向，它的正上方和正下方的像素比较。

注意，方向的正负是不起做用的，好比东南方向和西北方向是同样的，都认为是对角线的一个方向。前面咱们把梯度方向近似到水平，垂直和两个对角线四个方向，因此每一个像素根据自身方向在这四个方向之一进行比较，决定是否保留。这一部分的代码也很简单，列出以下。pModule，pDirection分别记录了上一步梯度模值和梯度方向。

pmoddrow = pModule + Width + 1; pdirdrow = pDirection + Width + 1; pstrongdrow = pStrong + Width + 1; for (i = 1; i < Hend - 1; i++) { 　　pstrongd = pstrongdrow; 　　pmodd = pmoddrow; 　　pdird = pdirdrow; 　　for (j = 1; j < Wend - 1; j++) 　　　　{ 　 switch (*pdird) { case 0:        // x direction
            case 4: if (*pmodd > *(pmodd - 1) && *pmodd > *(pmodd + 1)) *pstrongd = 255; break; case 1:        // northeast-southwest direction. Notice the data order on y direction of bmp data
            case 5: if (*pmodd > *(pmodd + Width + 1) && *pmodd > *(pmodd - Width - 1)) *pstrongd = 255; break; case 2:        // y direction
            case 6: if (*pmodd > *(pmodd - Width) && *pmodd > *(pmodd + Width)) *pstrongd = 255; break; case 3:        // northwest-southeast direction. Notice the data order on y direction of bmp data
            case 7: if (*pmodd > *(pmodd + Width - 1) && *pmodd > *(pmodd - Width + 1)) *pstrongd = 255; break; default: ASSERT(0); break; } pstrongd++; pmodd++; pdird++; } pstrongdrow += Width; pmoddrow += Width; pdirdrow += Width; }

View Code

下图是非最大值抑制的结果。可见边缘宽度已经大大减少。可是这个图像中由于没有应用任何阀值，还含有大量小梯度模值的点，也就是图中很暗的地方。下面，阀值要上场了。

             非最大值抑制结果

4. 双阀值。

通常的边缘检测算法用一个阀值来滤除噪声或颜色变化引发的小的梯度值，而保留大的梯度值。Canny算法应用双阀值，即一个高阀值和一个低阀值来区分边缘像素。若是边缘像素点梯度值大于高阀值，则被认为是强边缘点。若是边缘梯度值小于高阀值，大于低阀值，则标记为弱边缘点。小于低阀值的点则被抑制掉。这一步算法很简单。

5. 滞后边界跟踪。

至此，强边缘点能够认为是真的边缘。弱边缘点则多是真的边缘，也多是噪声或颜色变化引发的。为获得精确的结果，后者引发的弱边缘点应该去掉。一般认为真实边缘引发的弱边缘点和强边缘点是连通的，而又噪声引发的弱边缘点则不会。所谓的滞后边界跟踪算法检查一个弱边缘点的8连通领域像素，只要有强边缘点存在，那么这个弱边缘点被认为是真是边缘保留下来。

这个算法搜索全部连通的弱边缘，若是一条连通的弱边缘的任何一个点和强边缘点连通，则保留这条弱边缘，不然抑制这条弱边缘。搜索时能够用广度优先或者深度优先算法，我在这里实现了应该是最容易的深度优先算法。一次连通一条边缘的深度优先算法以下：

1. 准备一个栈s，一个队列q，设联通指示变量connected为假。从图像的第一个点开始，进入2。
2. 若是这个点是弱边界点而且没有被标记，把它标记，并把它做为第一个元素放入栈s中，同时把它放入记录连通曲线的队列q，进入3。若是这个点不是弱边界或者已经被标记过，到图像的下一个点，重复2。
3. 从栈s中取出一个元素，查找它的8像素领域。若是一个领域像素是弱边界而且没有被标记过，把这个领域像素标记，并加入栈s中，同时加入队列q。同时查找领域对应的强边界图，若是有一个像素是强边界，表示这条弱边界曲线和强边界联通，设置connected为真。重复3直到栈中没有元素了。若是connected为假，则依次从队列q中取出每一个元素，清空标记。若是connected为真，保留标记。
4. 清空队列q，设置connected为假，移动到图像的下一个点，到2。

// 5. Edge tracking by hysteresis
    stack<CPoint> s; queue<CPoint> q; BOOL connected = FALSE; long row_idx = Width; for (i = 1; i < Height - 1; i++, row_idx += Width) { for (j = 1; j < Width - 1; j++) { pweakd = pWeak + row_idx + j; if (*pweakd == 255) { s.push(CPoint(j, i)); q.push(CPoint(j, i)); *pweakd = 1;        // Label it

                while (!s.empty()) { CPoint p = s.top(); s.pop(); // Search weak edge 8-point neighborhood
                    pweakd = pWeak + p.y*Width + p.x; if (*(pweakd - Width - 1) == 255) { CPoint np = CPoint(p.x - 1, p.y - 1); s.push(np); q.push(np); *(pweakd - Width - 1) = 1;        // Label it
 } if (*(pweakd - Width) == 255) { CPoint np = CPoint(p.x, p.y - 1); s.push(np); q.push(np); *(pweakd - Width) = 1;        // Label it
 } if (*(pweakd - Width + 1) == 255) { CPoint np = CPoint(p.x + 1, p.y - 1); s.push(np); q.push(np); *(pweakd - Width + 1) = 1;        // Label it
 } if (*(pweakd - 1) == 255) { CPoint np = CPoint(p.x - 1, p.y); s.push(np); q.push(np); *(pweakd - 1) = 1;        // Label it
 } if (*(pweakd + 1) == 255) { CPoint np = CPoint(p.x + 1, p.y); s.push(np); q.push(np); *(pweakd + 1) = 1;        // Label it
 } if (*(pweakd + Width - 1) == 255) { CPoint np = CPoint(p.x - 1, p.y + 1); s.push(np); q.push(np); *(pweakd + Width - 1) = 1;        // Label it
 } if (*(pweakd + Width) == 255) { CPoint np = CPoint(p.x, p.y + 1); s.push(np); q.push(np); *(pweakd + Width) = 1;        // Label it
 } if (*(pweakd + Width + 1) == 255) { CPoint np = CPoint(p.x + 1, p.y + 1); s.push(np); q.push(np); *(pweakd + Width + 1) = 1;        // Label it
 } // Search strong edge 8-point neighborhood
                    if (connected == FALSE) { pstrongd = pStrong + p.y*Width + p.x; for (int m = -1; m <= 1; m++) { for (int n = -1; n <= 1; n++) { if (*(pstrongd + m*Width + n) == 255) { connected = TRUE; goto next; } } } } next: continue; } // No more element in the stack
                if (connected == FALSE) { // The weak edge is not connected to any strong edge. Suppress it.
                    while (!q.empty()) { CPoint p = q.front(); q.pop(); pWeak[p.y*Width + p.x] = 0; } } else { // Clean the queue
                    while (!q.empty()) q.pop(); connected = FALSE; } } } } // Add the connected weak edges (labeled) into strong edge image. // All strong edge pixels are labeled 255, otherwise 0.
    for (i = 0; i < len; i++) { if (pWeak[i] == 1) pStrong[i] = 255; }

View Code

下面是对Lena图计算Canny边缘检测的梯度模图和二值化图，高斯半径2，高阀值100，低阀值50。

  

             Canny检测梯度模图                        Canny检测梯度二值图

做为对比，下面是用一阶差分和Sobel算子对原图计算的结果，阀值100。因为一阶差分的梯度值相对较小，我对一阶差分的梯度值放大了必定倍数，使得它和Sobel的梯度值保持一样的水平。

  

               一阶差分梯度模图                        一阶差分梯度二值图

  

                 Sobel梯度模图                          Sobel梯度二值图

很明显，Canny边缘检测的效果是很显著的。相比普通的梯度算法大大抑制了噪声引发的伪边缘，并且是边缘细化，易于后续处理。对于对比度较低的图像，经过调节参数，Canny算法也能有很好的效果。

     

                       原图                            Canny梯度模，高斯半径2，低阀值30，高阀值100       Canny梯度二值化图，高斯半径2，低阀值30，高阀值100

     

                      原图                         Canny梯度模，高斯半径1，低阀值40，高阀值80  Canny梯度二值化图，高斯半径1，低阀值40，高阀值80

参考

Canny edge detector

canny 算法