Harris角点 C++

Harris角点检测

  因为国庆节回家的原因，导致很多内容没有更新，国庆假期已过，现就及时更新内容啦！另外国庆节偶然见得几句诗，摘录如下：
————————————————————————————————————————————
  我见青山多妩媚，料青山见我应如是！
————————————————————————————————————————————

正式内容

此处直接贴代码了，原理就不再给予解释了，大家自行搜索即可。

void Test_Harris(Mat &input, Mat &output)
{
	//clone原图像，用来做标记
	output = input.clone();

	//转化为灰度图像进行处理，对于含噪声的图像，需要先进行平滑滤波
	cvtColor(input, input, COLOR_BGR2GRAY);

	//做卷积之前先填充
	Mat imageSource_padded;
	copyMakeBorder(input, imageSource_padded, 1, 1, 1, 1, BORDER_CONSTANT, Scalar::all(0));

	//生成Mat用来存储梯度结果
	Mat imageSobelX = cv::Mat::zeros(imageSource_padded.size(), CV_32FC1);
	Mat imageSobelY = cv::Mat::zeros(imageSource_padded.size(), CV_32FC1);

	//遍历图像，获取Ix,Iy梯度
	for (int i = 1; i < imageSource_padded.rows - 1; i++)
	{
		for (int j = 1; j < imageSource_padded.cols - 1; j++)
		{
			float x1, x2, y1, y2;
			x1 = input.at<uchar>(i - 1, j + 1) + input.at<uchar>(i, j + 1) * 1.414 + input.at<uchar>(i + 1, j + 1);
			x2 = input.at<uchar>(i - 1, j - 1) + input.at<uchar>(i, j - 1) * 1.414 + input.at<uchar>(i + 1, j - 1);
			y1 = input.at<uchar>(i - 1, j - 1) + input.at<uchar>(i - 1, j) * 1.414 + input.at<uchar>(i - 1, j + 1);
			y2 = input.at<uchar>(i + 1, j - 1) + input.at<uchar>(i + 1, j) * 1.414 + input.at<uchar>(i + 1, j + 1);
			imageSobelX.at<float>(i, j) = abs(x1 - x2);
			imageSobelY.at<float>(i, j) = abs(y1 - y2);
		}
	}

	//获取Ix*Ix的值
	Mat SobelAmpXX = cv::Mat::zeros(imageSobelX.size(), CV_32FC1);
	for (int i = 1; i < imageSobelX.rows - 1; i++)
	{
		for (int j = 1; j < imageSobelX.cols - 1; j++)
		{
			SobelAmpXX.at<float>(i, j) = imageSobelX.at<float>(i, j)*imageSobelX.at<float>(i, j);
		}
	}

	//获取Iy*Iy的值
	Mat SobelAmpYY = cv::Mat::zeros(imageSobelY.size(), CV_32FC1);
	for (int i = 1; i < imageSobelY.rows - 1; i++)
	{
		for (int j = 1; j < imageSobelY.cols - 1; j++)
		{
			SobelAmpYY.at<float>(i, j) = imageSobelY.at<float>(i, j)*imageSobelY.at<float>(i, j);
		}
	}

	//获取Ix*Iy的值
	Mat SobelAmpXY = cv::Mat::zeros(imageSobelX.size(), CV_32FC1);
	for (int i = 1; i < imageSobelX.rows - 1; i++)
	{
		for (int j = 1; j < imageSobelX.cols - 1; j++)
		{
			SobelAmpXY.at<float>(i, j) = imageSobelX.at<float>(i, j)*imageSobelY.at<float>(i, j);
		}
	}

	//使用高斯函数进行平滑滤波
	Mat GaussianXX, GaussianYY, GaussianXY;
	GaussianBlur(SobelAmpXX, GaussianXX, Size(3, 3), BORDER_CONSTANT);
	GaussianBlur(SobelAmpYY, GaussianYY, Size(3, 3), BORDER_CONSTANT);
	GaussianBlur(SobelAmpXY, GaussianXY, Size(3, 3), BORDER_CONSTANT);

	//获取Harris响应
	float k = 0.05401;
	Mat Result = cv::Mat::zeros(GaussianXX.size(), CV_32FC1);
	for (int i = 1; i < GaussianXX.rows - 1; i++)
	{
		for (int j = 1; j < GaussianXX.cols - 1; j++)
		{
			float a, b, c;
			a = GaussianXX.at<float>(i, j);
			b = GaussianYY.at<float>(i, j);
			c = GaussianXY.at<float>(i, j);

			//用此式获取Harris响应
			Result.at<float>(i, j) = a * b - c * c - k * (a + b)*(a + b);
		}
	}

	//归一化，便于显示
	normalize(Result, Result, 0, 255, NORM_MINMAX);
	Result.convertTo(Result, CV_8UC1);
	imshow("Test_Harris_Response", Result);

	//int thresh = 60;
	//for (int i = 1; i < Result.rows-1; i++)
	//{
	//	for (int j = 1; j < Result.cols-1; j++)
	//	{
	//		if (Result.at<uchar>(i, j) > thresh)  //阈值判断
	//		{
	//			circle(output, Point(j - 1, i - 1), 5, Scalar(0, 0, 255), 1, 8); //标注角点
	//		}
	//	}
	//}

	//非极大值抑制
	int x[8] = { -1, -1, -1, 0, 0, 1, 1, 1 };
	int y[8] = { -1, 0, 1, -1, 1, -1, 0, 1 };
	for (int i = 1; i < Result.rows - 1; i++)
	{
		for (int j = 1; j < Result.cols - 1; j++)
		{
			int flag = 0;
			for (int k = 0; k < 8; k++)
			{
				if (Result.at<uchar>(i, j) <= Result.at<uchar>(i + x[k], j + y[k]))
				{
					//如果邻域内存在近似像素的值比中心位置的近似像素的值高，则进行非极大值抑制
					flag = 1;
					break;
				}
			}

			if (flag == 0)
			{
				//在原图像中的角点位置做标记
				circle(output, Point(j - 1, i - 1), 5, Scalar(0, 0, 255), 1, 8);

			}
		}
	}

}