图像处理--霍夫变换

  霍夫变换是一个非常重要的检测间断点边界形状的方法，他通过将图像坐标空间变换到参数空间，来实现直线和曲线的拟合。

一、直线检测：

1.直角坐标参数空间：

一条直线可以用数学表达式y = mx + c 或者 = x cos θ + y sinθ 表示。ρ是从原点到直线的垂直距离，θ 是直线的垂线与横轴顺时针方向的夹角（如果你使用的坐标系不同方向也可能不同，我是按OpenCV 使用的坐标系描述的）。如下图所示：
所以如果一条线在原点下方经过，θ的值就应该大于0，角度小于180。但是如果从原点上方经过的话，角度不是大于180，也是小于180，但ρ 的值小于0。垂直的线角度为0 度，水平线的角度为90 度。

霍夫变换是如何工作？每一条直线都可以用(ρ; θ) 表示。所以首先创建一个2D 数组（累加器），初始化累加器，所有的值都为0。行表示ρ，列表示θ。这个数组的大小决定了最后结果的准确性。如果你希望角度精确到1 度，你就需要180 列。对于ρ，最大值为图片对角线的距离。所以如果精确度要达到一个像素的级别，行数就应该与图像对角线的距离相等。

想象一下我们有一个大小为100x100 的直线位于图像的中央。取直线上的第一个点，我们知道此处的（x，y）值。把x 和y 带入上边的方程组，然后遍历θ 的取值：0，1，2，3，…，180。分别求出与其对应的 的值，这样我们就得到一系列(ρ; θ) 的数值对，如果这个数值对在累加器中也存在相应的位置，就在这个位置上加1。所以现在累加器中的（50，90）=1。（一个点可能存在与多条直线中，所以对于直线上的每一个点可能是累加器中的多个值同时加1）。

使用opencv实现上述过程：

#include<opencv2\opencv.hpp>
#include<stdio.h>
#include<iostream>

using namespace cv;
using namespace std;

int main(int arg, char* argv[])
{
	//在如图片并显示；
	Mat img = imread("H:\\MATLAB代码\\霍夫变换的理解\\线条.png");
	//namedWindow("原图");
	//imshow("原图", img);

	destroyAllWindows();
	//将图像灰度化并显示
	Mat grayImage; //创建无初始化矩阵
	cvtColor(img, grayImage, CV_RGB2GRAY);
	//namedWindow("灰度图");
	//imshow("灰度图", grayImage);
	//边缘检测
	Mat edge;
	Canny(grayImage, edge, 3, 9, 3);

	int i, j;

	//行列
	int row = grayImage.rows;
	int col = grayImage.cols;

	//极径最大值为对角线+宽
	//cvCeil 取整，返回不小于参数的整数。
	int max_r = col + cvCeil(sqrt(double(row*row + col*col)));

	//累加器 三角函数
	int *line_cnt[180];  //指针数组，  是一个普通的数组，数组中的每一个元素都是指针。
	double sin_[180], cos_[180], rad_ = CV_PI / 180;
	for (i = 0; i < 180; i++)
	{
		//初始化累加器为0
		line_cnt[i] = new int[max_r]();   //列为180列，列表示theta，行表示r,r最大值为max_r;
		//初始化三角函数
		sin_[i] = sin(i*rad_);
		cos_[i] = cos(i*rad_);
	}

	//极径，极角
	int r = 0;
	int theta = 0;

	//遍历图像，判断并进行累加。
	uchar *p;
	for (i = 0; i < row; i++)
	{
		//用指针遍历图像中的元素。
		//每一行图像的指针；
		p = edge.ptr<uchar>(i);
		for (j = 0; j < col; j++)
		{
			if (p[j] != 0)
			{
				//不等于0，证明是有用的边缘像素点
				//由该点的参数方程 更改累加器
				for (theta = 0; theta < 180; theta++)
				{
					//极坐标 直线方程
					r = cvRound(j*cos_[theta] + i*sin_[theta]);
					//偏移量， 因为累计器（霍夫矩阵的索引没有负值）
					r = r + col;
					//累加器完成累加
					line_cnt[theta][r]++;
				}
			}
		}
	}

	//存放取出最长的n条线
	int n = 6;
	int *line_n[3];            //三行分别用来存储极角、极径、共点数目
	line_n[0] = new int[n]();  //极角
	line_n[1] = new int[n]();  //极径
	line_n[2] = new int[n]();  //共点数目

	int tt = 0, rr = 0, cnt = 0;
	//寻找累加器中最大的值 
	for (theta = 0; theta < 180; theta++)
	{
		for (r = 0; r < max_r; r++)
		{
			//最少共点 < 这条直线的共点 则替换 并尝试进行冒泡
			if (line_n[2][n-1] < line_cnt[theta][r])
			{
				line_n[0][n-1] = theta;
				//累计的时候偏移过，将偏移修改回来
				line_n[1][n-1] = r - col;
				line_n[2][n-1] = line_cnt[theta][r];
				//冒泡排序
				for (i = n - 1; i > 0; i--)
				{
					//如果大于 则交换
					if (line_n[2][i] > line_n[2][i - 1])
					{
						tt = line_n[0][i];
						rr = line_n[1][i];
						cnt = line_n[2][i];
						line_n[0][i] = line_n[0][i - 1];
						line_n[1][i] = line_n[1][i - 1];
						line_n[2][i] = line_n[2][i - 1];
						line_n[0][i - 1] = tt;
						line_n[1][i - 1] = rr;
						line_n[2][i - 1] = cnt;
					}
					else
						break;
				}
			}
		}
	}

	//画出线段
	for (i = 0; i < n; i++)
	{
		Point pt1, pt2;
		double a = cos_[line_n[0][i]], b = sin_[line_n[0][i]];
		double x0 = a*line_n[1][i], y0 = b*line_n[1][i];
		//用点画出检测到的直线。
		//下面四个式子为纯直线几何推导。
		pt1.x = cvRound(x0 + max_r*(-b));
		pt1.y = cvRound(y0 + max_r*(a));
		pt2.x = cvRound(x0 - max_r*(-b));
		pt2.y = cvRound(y0 - max_r*(a));

		//绿线
		line(img, pt1, pt2, Scalar(0, 255, 0), 1, CV_AA);
	}
	cvNamedWindow("直线检测图");
	imshow("直线检测图",img);
	waitKey(0);
	//----------------------释放内存-----------------
	for (i = 0; i < 180; i++)
	{
		delete[]line_cnt[i];
	}
	delete[]line_n[0];
	delete[]line_n[1];
	delete[]line_n[2];
	return 0;
}