Canny的C++实现

好的图像边缘检测应该满足

1. 尽可能的标记处所有的边缘；
2. 标记出的边缘就是实际图像内容的边缘；
3. 图像中的边缘只标记一次。

Canny边缘检测的步骤

1. 图像灰度化；
2. 高斯滤波（不限制）降噪；
3. 使用简单的算子（如Sobel、Prewitt等）检测图像水平g r a d _ x grad_x g r a d _x、垂直梯度g r a d _ y grad_y g r a d _y，并求出梯度g r a d _ x 2 + g r a d _ y 2 \sqrt{grad_x^2 + grad_y^2}g r a d _x 2 +g r a d _y 2 ​和梯度方向（角度）g r a d _ y / g r a d _ x grad_y/grad_x g r a d _y /g r a d _x；
4. 对梯度进行NMS；
5. 使用双阈值法确定强边缘、非边缘、弱边缘，以及弱边缘的二次判定。

C++实现示例

auxiliary.hpp

#pragma once

#include
#include
#include
#include

int GaussianKernel(int kernel_size, std::vector<std::vector<float>> &kernel)
{
    kernel.clear();
    kernel.resize(kernel_size);
    for (auto &it : kernel)
    {
        it.resize(kernel_size);
    }

    std::vector<int> coord_val(kernel_size, -kernel_size / 2);
    for (int i = 1; i < kernel_size; ++i)
    {
        coord_val[i] = coord_val[i - 1] + 1;
    }

    const float kSigma = 0.5f;
    float val1 = 1 / (2 * M_PI * kSigma * kSigma);
    float val2 = -1 / (2 * kSigma * kSigma);
    for (int i = 0; i < kernel_size; ++i)
    {
        for (int j = 0; j < kernel_size; ++j)
        {
            kernel[i][j] = val1 * exp(val2 * (pow(coord_val[i], 2) + pow(coord_val[j], 2)));
        }
    }

    float sum = 0.0f;
    for (int i = 0; i < kernel_size; ++i)
    {
        sum += std::accumulate(kernel[i].begin(), kernel[i].end(), 0.0);
    }

    for (int i = 0; i < kernel_size; ++i)
    {
        for (int j = 0; j < kernel_size; ++j)
        {
            kernel[i][j] /= sum;
        }
    }

    return 0;
}

int GaussianBlur(const std::vector<std::vector<float>> &kernel, const cv::Mat &src, cv::Mat &formula_gaussian)
{
    int kernel_size = kernel.size();

    int pad = kernel_size >> 1;
    cv::Mat formula_pad = cv::Mat(src.rows + 2 * pad, src.cols + 2 * pad, CV_32FC1);
    memset(formula_pad.data, 0, formula_pad.rows * formula_pad.cols * sizeof(float));
    cv::Rect rect(pad, pad, src.cols, src.rows);
    src.copyTo(formula_pad(rect));

    formula_gaussian = cv::Mat(src.rows, src.cols, CV_32FC1);
    memset(formula_gaussian.data, 0, formula_gaussian.rows * formula_gaussian.cols * sizeof(float));
    for (int i = 0; i < src.rows; ++i)
    {
        for (int j = 0; j < src.cols; ++j)
        {
            float tmp = 0.0f;
            for (int m = 0; m < kernel_size; ++m)
            {
                for (int n = 0; n < kernel_size; ++n)
                {
                    tmp += kernel[m][n] * formula_pad.at<float>(i + m, j + n);
                }
            }

            formula_gaussian.at<float>(i, j) = tmp;
        }
    }

    return 0;
}

int DoSobel(const cv::Mat &formula_gaussian, cv::Mat &sobel_xy, cv::Mat &sobel_angle, float &mean_sobel_xy)
{

    const int kSobelKernel = 3;
    int pad = kSobelKernel >> 1;
    std::vector<std::vector<int> > kernel_x{ {-1, 0, 1}, {-2, 0, 2}, {-1, 0, 1} };
    std::vector<std::vector<int> > kernel_y{ {1, 2, 1}, {0, 0, 0}, {-1, -2, -1} };
    cv::Mat gaussian_pad(formula_gaussian.rows + 2 * pad, formula_gaussian.cols + 2 * pad, CV_32FC1);
    memset(gaussian_pad.data, 0, gaussian_pad.rows * gaussian_pad.cols * sizeof(float));
    cv::Rect rect = cv::Rect(pad, pad, formula_gaussian.cols, formula_gaussian.rows);
    formula_gaussian.copyTo(gaussian_pad(rect));

    cv::Mat sobel_x(formula_gaussian.rows, formula_gaussian.cols, CV_32FC1);
    cv::Mat sobel_y(formula_gaussian.rows, formula_gaussian.cols, CV_32FC1);
    sobel_xy = cv::Mat(formula_gaussian.rows, formula_gaussian.cols, CV_32FC1);
    sobel_angle = cv::Mat(formula_gaussian.rows, formula_gaussian.cols, CV_32FC1);
    const float kEps = 0.000001f;
    float sum_sobel_xy = 0.0f;
    mean_sobel_xy = 0.0f;
    for (int i = 0; i < formula_gaussian.rows; ++i)
    {
        for (int j = 0; j < formula_gaussian.cols; ++j)
        {
            float gradient_x = 0.0f;
            float gradient_y = 0.0f;
            float gradient_xy = 0.0f;
            for (int m = 0; m < kSobelKernel; ++m)
            {
                for (int n = 0; n < kSobelKernel; ++n)
                {
                    gradient_x += kernel_x[m][n] * gaussian_pad.at<float>(i + m, j + n);
                    gradient_y += kernel_y[m][n] * gaussian_pad.at<float>(i + m, j + n);
                }
            }
            gradient_xy = sqrt(pow(gradient_x, 2) + pow(gradient_y, 2));
            sobel_x.at<float>(i, j) = gradient_x;
            sobel_y.at<float>(i, j) = gradient_y;
            sobel_xy.at<float>(i, j) = gradient_xy;
            sobel_angle.at<float>(i, j) = gradient_y / (fabs(gradient_x) > kEps ? gradient_x : kEps);

            sum_sobel_xy += gradient_xy;
        }
    }
    mean_sobel_xy = sum_sobel_xy / (sobel_xy.rows * sobel_xy.cols);

    return 0;
}

int DoNMS(const cv::Mat &sobel_xy, const cv::Mat &sobel_angle, cv::Mat &nms_sobel_xy)
{
    nms_sobel_xy = cv::Mat(sobel_xy.rows, sobel_xy.cols, CV_32FC1);
    memset(nms_sobel_xy.data, 0, nms_sobel_xy.rows * nms_sobel_xy.cols * sizeof(float));

    for (int i = 1; i < nms_sobel_xy.rows - 1; ++i)
    {
        for (int j = 1; j < nms_sobel_xy.cols - 1; ++j)
        {

            float en0 = sobel_xy.at<float>(i - 1, j - 1);
            float en1 = sobel_xy.at<float>(i - 1, j);
            float en2 = sobel_xy.at<float>(i - 1, j + 1);
            float en3 = sobel_xy.at<float>(i, j - 1);
            float en4 = sobel_xy.at<float>(i, j);
            float en5 = sobel_xy.at<float>(i, j + 1);
            float en6 = sobel_xy.at<float>(i + 1, j - 1);
            float en7 = sobel_xy.at<float>(i + 1, j);
            float en8 = sobel_xy.at<float>(i + 1, j + 1);

            float grad_inter1 = 0.0f;
            float grad_inter2 = 0.0f;
            float ratio = 0.0f;
            float angle = sobel_angle.at<float>(i, j);
            if (angle >= 0)
            {
                if (angle >= 1)
                {
                    ratio = 1.0f / angle;

                    grad_inter1 = (1 - ratio) * en1 + ratio * en2;
                    grad_inter2 = (1 - ratio) * en7 + ratio * en6;
                }
                else
                {
                    ratio = angle;

                    grad_inter1 = (1 - ratio) * en5 + ratio * en2;
                    grad_inter2 = (1 - ratio) * en3 + ratio * en6;
                }
            }
            else
            {
                if (fabs(angle) >= 1)
                {
                    ratio = 1.0f / fabs(angle);

                    grad_inter1 = (1 - ratio) * en1 + ratio * en0;
                    grad_inter2 = (1 - ratio) * en7 + ratio * en8;
                }
                else
                {
                    ratio = fabs(angle);

                    grad_inter1 = (1 - ratio) * en3 + ratio * en0;
                    grad_inter2 = (1 - ratio) * en5 + ratio * en8;
                }
            }
            if (en4 > grad_inter1 && en4 > grad_inter2)
            {
                nms_sobel_xy.at<float>(i, j) = en4;
            }
        }
    }

    return 0;
}

int DoBinaryThresh(const cv::Mat &nms_sobel_xy, float mean_sobel_xy, cv::Mat &img_canny)
{
    cv::Mat binary_thresh_canny = cv::Mat(nms_sobel_xy.rows, nms_sobel_xy.cols, CV_32FC1);
    memset(binary_thresh_canny.data, 0, binary_thresh_canny.rows * binary_thresh_canny.cols * sizeof(float));

    float low_thresh = mean_sobel_xy * 0.5f;
    float high_thresh = low_thresh * 3.0f;
    for (int i = 1; i < binary_thresh_canny.rows - 1; ++i)
    {
        for (int j = 1; j < binary_thresh_canny.cols - 1; ++j)
        {
            float en0 = nms_sobel_xy.at<float>(i - 1, j - 1);
            float en1 = nms_sobel_xy.at<float>(i - 1, j);
            float en2 = nms_sobel_xy.at<float>(i - 1, j + 1);
            float en3 = nms_sobel_xy.at<float>(i, j - 1);
            float en4 = nms_sobel_xy.at<float>(i, j);
            float en5 = nms_sobel_xy.at<float>(i, j + 1);
            float en6 = nms_sobel_xy.at<float>(i + 1, j - 1);
            float en7 = nms_sobel_xy.at<float>(i + 1, j);
            float en8 = nms_sobel_xy.at<float>(i + 1, j + 1);
            if (en4 >= high_thresh)
            {
                binary_thresh_canny.at<float>(i, j) = 255.0f;
            }
            else if (en4  low_thresh)
            {
                binary_thresh_canny.at<float>(i, j) = 0.0f;
            }
            else
            {
                if (en0 >= high_thresh || en1 >= high_thresh || en2 >= high_thresh
                    || en3 >= high_thresh || en5 >= high_thresh
                    || en6 >= high_thresh || en7 >= high_thresh || en8 >= high_thresh)
                {
                    binary_thresh_canny.at<float>(i, j) = 255.0f;
                }
            }
        }
    }

    binary_thresh_canny.convertTo(img_canny, CV_8UC1);

    return 0;
}

Canny.cpp

#include "auxiliary.hpp"
#include
#include

int main()
{

    cv::Mat img = cv::imread("../data/推导2.jpg", 0);
    cv::Mat imgf;
    img.convertTo(imgf, CV_32FC1, 1 / 255.0);

    int kKernelSize = 5;
    std::vector<std::vector<float>> kernel;
    int ret = GaussianKernel(kKernelSize, kernel);
    if (0 != ret)
    {
        system("pause");
        return ret;
    }

    cv::Mat formula_gaussian;
    ret = GaussianBlur(kernel, imgf, formula_gaussian);
    if (0 != ret)
    {
        system("pause");
        return ret;
    }

    cv::Mat sobel_xy, sobel_angle;
    float mean_sobel_xy = 0.0f;
    ret = DoSobel(formula_gaussian, sobel_xy, sobel_angle, mean_sobel_xy);
    if (0 != ret)
    {
        system("pause");
        return ret;
    }

    cv::Mat nms_sobel_xy;
    ret = DoNMS(sobel_xy, sobel_angle, nms_sobel_xy);
    if (0 != ret)
    {
        system("pause");
        return ret;
    }

    cv::Mat img_canny;
    ret = DoBinaryThresh(nms_sobel_xy, mean_sobel_xy, img_canny);
    if (0 != ret)
    {
        system("pause");
        return ret;
    }

    cv::imshow("img_canny", img_canny);
    cv::waitKey(0);

    return 0;
}

非极大值抑制时双线性插值示意图

图1 非极大值抑制时双线性插值示意图

; 实验图像及结果

以原图 推导2.jpg为示例

图2 推导2.jpg

gaussian.jpg

图3 gaussian.jpg

sobel_x.jpg

图4 sobel_x.jpg

sobel_y.jpg

图5 sobel_y.jpg

sobel_xy.jpg

图6 sobel_xy.jpg

sobel_angle.jpg在NMS时插值计算梯度正反两个方向虚拟像素点梯度时使用。

图7 sobel_angle.jpg

nms_sobel_xy.jpg

图8 nms_sobel_xy.jpg

最终结果 img_canny.jpg

图9 最终结果 img_canny.jpg

Reference

Canny算子中的非极大值抑制（Non-Maximum Suppression）分析

Original: https://blog.csdn.net/liugan528/article/details/122824546
Author: GarryLau
Title: Canny的C++实现