Canny边缘检测源码

#include <iostream>
#include <opencv2/core.hpp>
#include <opencv2/highgui.hpp>
#include <opencv2/imgproc.hpp>
#include<vector>
#include<Eigen/Core>
#include <corecrt_math_defines.h>
using namespace std;
using namespace cv;

class canny
{
public:
	Mat src;
	// 高斯滤波算子
	double** guassArray(int size, double sigma);
	// 高斯滤波
	Mat guassFilter(Mat& src,int size,double sigma);
	// 像素判定
	uchar pixelJudge(double pixel);
	// sobel梯度计算:梯度方向与边缘方向垂直
	Mat sobel(Mat& src);
	// 非极大值抑制（NMS）
	/*
	非极大值抑制是一种将宽边缘带细化，只保留边缘峰值的处理。主要的手法是根据上一步获取
	的梯度方向，将图像像素点沿梯度方向或逆方向的邻域像素点进行比较，如果为最大值则保留,
	否则抑制，即设置像素点为0。
	*/
	Mat NMS(Mat& src);
	/*双阈值连接
	设置两个阈值，一个高阈值和一个低阈值。图像像素点如果高于高阈值则表示是强边缘点，是
	真实边缘点，高于低阈值但小于高阈值表示是弱边缘点，其它不是边缘，抑制。之后需要对弱
	边缘点进行处理，若其邻域存在强边缘点则表示是真实边缘点，否则不是。这样设置两个阈值
	可以滤除图像中噪声，改善图像质量。弱边缘的处理原则是因为真实边缘的弱边缘点都存在强边缘点。
	*/
	Mat doubleThresh(Mat& src, uchar high, uchar low);
	//测试
	Mat cannyTest(Mat& src, int size, double sigma, uchar high, uchar low);
};

Mat canny::cannyTest(Mat& src, int size, double sigma, uchar high, uchar low)
{
	Mat dst = guassFilter(src, size, sigma);
	Mat edge = sobel(dst);
	edge = NMS(edge);
	vector<Mat> grade;
	split(edge, grade);
	imshow("最大值抑制", grade[0]);
	edge = doubleThresh(edge, 30, 100);
	imshow("双阈值连接", edge);
	return edge;
}

uchar canny::pixelJudge(double pixel)
{
	if (pixel < 0)
	{
		return 0;
	}
	else if (pixel > 255)
	{
		return 255;
	}
	else
	{
		return static_cast<uchar>(pixel);
	}
}

double** canny::guassArray(int size, double sigma)
{
	// [1] 初始化权值数组
	double** array = new double* [size];
	for (int i = 0; i < size; i++) {
		array[i] = new double[size];
	}
	// [2] 高斯分布计算
	int center_i, center_j;
	center_i = center_j = size / 2;
	double pi = M_PI;
	double sum = 0;
	for (int i = 0; i < size; i++)
	{
		for (int j = 0; j < size; j++)
		{
			array[i][j] =
				//后面进行归一化，这部分可以不用
				//0.5f *pi*(sigma*sigma) * 
				exp(-(1.0f) * (((i - center_i) * (i - center_i) + (j - center_j) * (j - center_j)) /
					(2.0f * sigma * sigma)));
			sum += array[i][j];
			//cout << array[i][j] << endl;
		}
	}
	// [2-2] 归一化求权值
	for (int i = 0; i < size; i++)
	{
		for (int j = 0; j < size; j++)
		{
			array[i][j] /= sum;
			//cout << array[i][j] << "  ";
		}
	}

	return array;
};

Mat canny::guassFilter(Mat& src, int size, double sigma)
{
	Mat dst = src.clone();
	double** array = guassArray(size, sigma);
	for (int i = 0; i < src.rows; i++)
	{
		for (int j = 0; j < src.cols; j++)
		{
			double sum = 0;
			for (int m = -size / 2; m <= size / 2; m++)
			{
				for (int n = -size / 2; n <= size / 2; n++)
				{
					// 边缘补零
					if (i + m < 0 || i + m >= src.rows || j + n < 0 || j + n >= src.cols)
					{
						sum += 0;
					}
					else
					{
						sum += src.at<uchar>(i + m, j + n)*array[m+size/2][n+size/2];
					}
				}
			}
			dst.at<uchar>(i, j) = pixelJudge(sum);
		}
	}
	return dst;
}

Mat canny::sobel(Mat& src)
{
	Mat gradeValue = src.clone();
	Mat gradeDirect = src.clone();
	// x方向边缘的y模板
	double model1[9] = { -1, -2, -1, 0, 0, 0, 1, 2, 1 };
	// y方向边缘的x模板
	double model2[9] = { -1, 0, 1, -2, 0, 2, -1, 0, 1 };

	for (int i = 0; i < src.rows; i++)
	{
		for (int j = 0; j < src.cols; j++)
		{
			int p = 0;
			double sum_x = 0;
			double sum_y = 0;
			for (int m = -1; m <= 1; m++)
			{
				for (int n = -1; n <= 1; n++)
				{
					// 边缘补零
					if (i + m < 0 || i + m >= src.rows || j + n < 0 || j + n >= src.cols)
					{
						p++;
					}
					else
					{
						sum_y += src.at<uchar>(i + m, j + n) * model1[p];
						sum_x += src.at<uchar>(i + m, j + n) * model2[p];
						p++;
					}
				}
			}
			gradeValue.at<uchar>(i, j) = pixelJudge(pow(pow(sum_x, 2) + pow(sum_y, 2), 1.0 / 2));
			// 0,7垂直方向， 1,2斜对角方向\, 3,4水平方向, 5,6斜对角方向/
			// 梯度方向与边缘垂直，因此model1检测出横向边缘，梯度方向为竖向
			// 同理model2检测出竖向边缘，梯度方向为横向，因此梯度方向为y/x
			//arctanx的范围是-pi/2到pi/2
			gradeDirect.at<uchar>(i, j) = static_cast<uchar>(atan(sum_y / (sum_x + 0.001)) * 8 / CV_PI + 4);
			//cout << (int)gradeDirect.at<uchar>(i, j) << endl;
		}
	}

	Mat temp[] = { gradeValue, gradeDirect };
	Mat grad;
	// 将梯度幅值与梯度方向合并,通道数为2，类似RGB通道
	// 但是不能显示，因为只能显示通道数1、3、4的图像
	// Mat dst1 = Mat::zeros(src.rows, src.cols * 2, CV_8UC2);
	// imshow(dst1);
	merge(temp, 2, grad);

	//对多通道像素进行访问
	//int num = 0;
	//for (int i = 0; i <=0; i++)
	//{
	//	for (int j = 0; j <= 1; j++)
	//	{
	//		cout << (int)grad.at<Vec2b>(i, j)[0] << endl;
	//		cout << (int)grad.at<Vec2b>(i, j)[1] << endl;
	//	}
	//}

	//显示梯度图
	//imshow("edge", gradeValue);
	return grad;
}

Mat canny::NMS(Mat& src)
{
	// 梯度通道和方向通道
	vector<Mat> grad;
	split(src, grad);
	Mat gradValue = grad[0].clone();
	//imshow("edge", grade[0]);
	for (int i = 0; i < gradValue.rows; i++)
		for (int j = 0; j < gradValue.cols; j++)
		{
			//将宽边缘带细化：对于梯度方向是垂直方向，边缘是水平方向，想要细化边缘，就要比较上下gradValue
			//梯度方向是垂直方向，比较上下gradValue
			if (grad[1].at<uchar>(i, j) == 0 || grad[1].at<uchar>(i, j) == 7)
			{
				if (i < gradValue.rows - 1)
					if (gradValue.at<uchar>(i, j) <= gradValue.at<uchar>(i + 1, j))
					{
						grad[0].at<uchar>(i, j) = 0;
						continue;
					}
				if (i > 0)
					if (gradValue.at<uchar>(i, j) <= gradValue.at<uchar>(i - 1, j))
					{
						grad[0].at<uchar>(i, j) = 0;
						continue;
					}
			}
			//梯度方向是“\”方向，比较对角gradValue
			if (grad[1].at<uchar>(i, j) == 1 || grad[1].at<uchar>(i, j) == 2)
			{
				if (i < gradValue.rows - 1 && j < gradValue.cols - 1)
					if (gradValue.at<uchar>(i, j) <= gradValue.at<uchar>(i + 1, j + 1))
					{
						grad[0].at<uchar>(i, j) = 0;
						continue;
					}
				if (i > 0 && j > 0)
					if (gradValue.at<uchar>(i, j) <= gradValue.at<uchar>(i - 1, j - 1))
					{
						grad[0].at<uchar>(i, j) = 0;
						continue;
					}
			}
			//梯度方向是水平方向，比较左右gradValue
			if (grad[1].at<uchar>(i, j) == 3 || grad[1].at<uchar>(i, j) == 4)
			{
				if (j < gradValue.cols - 1)
					if (gradValue.at<uchar>(i, j) <= gradValue.at<uchar>(i, j + 1))
					{
						grad[0].at<uchar>(i, j) = 0;
						continue;
					}
				if (j > 0)
					if (gradValue.at<uchar>(i, j) <= gradValue.at<uchar>(i, j - 1))
					{
						grad[0].at<uchar>(i, j) = 0;
						continue;
					}
			}
			//梯度方向是“/”方向，比较对角gradValue
			if (grad[1].at<uchar>(i, j) == 5 || grad[1].at<uchar>(i, j) == 6)
			{
				if (i < gradValue.rows - 1 && j > 0)
					if (gradValue.at<uchar>(i, j) <= gradValue.at<uchar>(i + 1, j - 1))
					{
						grad[0].at<uchar>(i, j) = 0;
						continue;
					}
				if (i > 0 && j < gradValue.cols - 1)
					if (gradValue.at<uchar>(i, j) <= gradValue.at<uchar>(i - 1, j + 1))
					{
						grad[0].at<uchar>(i, j) = 0;
						continue;
					}
			}
		}
	//imshow("非极大值处理", grad[0]);
	merge(grad, src);
	return src;
}

Mat canny::doubleThresh(Mat& src, uchar high, uchar low)
{
	vector<Mat> grad;
	// 梯度通道和方向通道
	split(src, grad);
	for (int i = 0; i < grad[0].rows; i++)
		for (int j = 0; j < grad[0].cols; j++)
		{
			if (grad[0].at<uchar>(i, j) > high)
			{
				grad[0].at<uchar>(i, j) = 255;
				// 标记强边缘点的位置
				grad[1].at<uchar>(i, j) = 2;
			}
			else if (grad[0].at<uchar>(i, j) > low)
			{
				grad[0].at<uchar>(i, j) = 0;
				// 标记弱边缘点的位置
				grad[1].at<uchar>(i, j) = 1;
			}
			else
			{
				grad[0].at<uchar>(i, j) = 0;
				grad[1].at<uchar>(i, j) = 0;
			}
		}
	// 真实的边缘会在弱边缘点的邻域内存在强边缘点
	for (int i = 0; i < grad[0].rows; i++)
		for (int j = 0; j < grad[0].cols; j++)
		{
			if (grad[1].at<uchar>(i, j) == 1)
			{
				for (int n = -1; n <= 1; n++)
					for (int m = -1; m <= 1; m++)
					{
						if (i + n >= 0 && j + m >= 0 && i + n < src.rows && j + m < src.cols && grad[1].at<uchar>(i + n, j + m) == 2)
							grad[0].at<uchar>(i, j) = 255;
					}
			}
		}
	return grad[0];
}

int main() 
{
	Mat imgLeft = imread("D:\\VC\\c++\\opencv源码\\opencv源码\\Stereo_Sample\\ImageL1.jpg", IMREAD_GRAYSCALE);
	canny c;
	c.src = imgLeft;
	c.cannyTest(c.src, 5, 1, 40, 100);
	cv::waitKey(0);
	return 0;
}