Hough变换原理及应用

霍夫变换原理介绍可参见冈萨勒斯的《数字图像处理》第二版476页。

在xy平面上的一点(x0,y0)要转换为极坐标下的(θ,ρ)：

x0 = ρcosθ

y0 = ρsinθ

那么

x0cosθ = ρcosθcosθ

y0sinθ = ρsinθsinθ

可得

x0cosθ + y0sinθ = ρ

那么xy坐标平面到极坐标θρ平面的转换公式为：

xcosθ + ysinθ = ρ

可以看出xy平面上的一点在极坐标θρ平面上的投影为一条正弦曲线，那么n个点对应n条正弦曲线，这n条曲线如果有m条曲线线交于一点(θ0,ρ0)，那么这m条曲线对应在xy平面下的m个点就在同一直线上，而θ0就是这条直线的法线与x轴的夹角，ρ0为原点到直线的距离。

如果理解了霍夫变换的原理，接下来就可以看一下变换源码了。

下面是OpenCV自带的标准霍夫变换函数：

static void
icvHoughLinesStandard( const CvMat* img, float rho, float theta,
                       int threshold, CvSeq *lines, int linesMax )
{
    cv::AutoBuffer<int> _accum, _sort_buf;
    cv::AutoBuffer<float> _tabSin, _tabCos;

    const uchar* image;
    int step, width, height;
    int numangle, numrho;
    int total = 0;
    int i, j;
    float irho = 1 / rho;
    double scale;

    CV_Assert( CV_IS_MAT(img) && CV_MAT_TYPE(img->type) == CV_8UC1 );

    image = img->data.ptr;
    step = img->step;
    width = img->cols;
    height = img->rows;

    numangle = cvRound(CV_PI / theta);
    numrho = cvRound(((width + height) * 2 + 1) / rho);

    _accum.allocate((numangle+2) * (numrho+2));
    _sort_buf.allocate(numangle * numrho);
    _tabSin.allocate(numangle);
    _tabCos.allocate(numangle);
    int *accum = _accum, *sort_buf = _sort_buf;
    float *tabSin = _tabSin, *tabCos = _tabCos;

    memset( accum, 0, sizeof(accum[0]) * (numangle+2) * (numrho+2) );

    float ang = 0;
    for(int n = 0; n < numangle; ang += theta, n++ )
    {
        tabSin[n] = (float)(sin((double)ang) * irho);
        tabCos[n] = (float)(cos((double)ang) * irho);
    }

    // stage 1. fill accumulator
    for( i = 0; i < height; i++ )
        for( j = 0; j < width; j++ )
        {
            if( image[i * step + j] != 0 )
                for(int n = 0; n < numangle; n++ )
                {
                    int r = cvRound( j * tabCos[n] + i * tabSin[n] );
                    r += (numrho - 1) / 2;
                    accum[(n+1) * (numrho+2) + r+1]++;
                }
        }

    // stage 2. find local maximums
    for(int r = 0; r < numrho; r++ )
        for(int n = 0; n < numangle; n++ )
        {
            int base = (n+1) * (numrho+2) + r+1;
            if( accum[base] > threshold &&
                accum[base] > accum[base - 1] && accum[base] >= accum[base + 1] &&
                accum[base] > accum[base - numrho - 2] && accum[base] >= accum[base + numrho + 2] )
                sort_buf[total++] = base;
        }

    // stage 3. sort the detected lines by accumulator value
    icvHoughSortDescent32s( sort_buf, total, accum );

    // stage 4. store the first min(total,linesMax) lines to the output buffer
    linesMax = MIN(linesMax, total);
    scale = 1./(numrho+2);
    for( i = 0; i < linesMax; i++ )
    {
        CvLinePolar line;
        int idx = sort_buf[i];
        int n = cvFloor(idx*scale) - 1;
        int r = idx - (n+1)*(numrho+2) - 1;
        line.rho = (r - (numrho - 1)*0.5f) * rho;
        line.angle = n * theta;
        cvSeqPush( lines, &line );
    }
}

步骤为：

1、将θ离散值的sin、cos的计算结果存在数组中，便于下一步使用。

2、将原图中灰度值不为0的点变换到θρ空间的正弦曲线，而θρ空间相当于一个累加器。

3、将累加器的值进行排序，使用的是最大邻域值法，可以看出邻域很小，邻域小计算速度就快。

4、输出累加器中值最大的几个值。

参见上面的源码写出自己的代码与之比较。

原图：

源码：

<pre name="code" class="cpp">#include "stdafx.h"
#include <cv.h>
#include <opencv2\highgui\highgui.hpp>

using namespace cv;

int _tmain(int argc, _TCHAR* argv[])
{
	Mat src = imread("E:\\VC++Projects\\road.jpg");
	Mat gray = Mat::zeros(src.size(),CV_8UC1);
	cvtColor(src, gray, CV_BGR2GRAY);

	// 求得x和y方向的一阶微分
	Mat sobelx;
	Mat sobely;
	Sobel(gray, sobelx, CV_32F, 1, 0, 3);
	Sobel(gray, sobely, CV_32F, 0, 1, 3);
	gray.release();

	// 求得梯度和方向
	Mat norm;
	Mat dir;
	cartToPolar(sobelx, sobely, norm, dir);
	sobelx.release();
	sobely.release();

	// 转换为8位单通道图像
	double normMax;
	minMaxLoc(norm, NULL, &normMax);
	Mat grad;
	norm.convertTo(grad, CV_8UC1, 255.0/normMax, 0);

	// OpenCV自带大津阈值函数处理
	Mat cvOstu;
	threshold(grad, cvOstu, 0, 255, CV_THRESH_BINARY|CV_THRESH_OTSU);
	norm.release();
	dir.release();
	grad.release();

	//vector<Vec2f> lines;
	//HoughLines( cvOstu, lines, 1, CV_PI/180, 100 );

	//for( size_t i = 0; i < 10/*lines.size()*/; i++ )
	//{
	//	float rho = lines[i][0];
	//	float theta = lines[i][1];
	//	double a = cos(theta), b = sin(theta);
	//	double x0 = a*rho, y0 = b*rho;
	//	Point pt1(cvRound(x0 + 1000*(-b)),
	//		cvRound(y0 + 1000*(a)));
	//	Point pt2(cvRound(x0 - 1000*(-b)),
	//		cvRound(y0 - 1000*(a)));
	//	line( src, pt1, pt2, Scalar(0,0,255), 1, 8 );
	//}

	vector<Vec4i> lines2;
	HoughLinesP( cvOstu, lines2, 1, CV_PI/180, 80, 30, 10 );
	for( size_t i = 0; i < 10/*lines2.size()*/; i++ )
	{
		line( src, Point(lines2[i][0], lines2[i][1]),
			Point(lines2[i][2], lines2[i][3]), Scalar(0,255,0), 1, 8 );
	}

	namedWindow("src");
	imshow("src", src);
	namedWindow("binary");
	imshow("binary", cvOstu);
	waitKey(0);
	return 0;
}

使用标准霍夫变换结果：

使用概率霍夫变换的结果：

看到上面的结果，会发现虽然找出了直线，但是可能结果并不是想要的，比如在几条直线过于接近(θ，ρ的差值很小)的只显示一条。那么必须扩大最大值的搜素邻域，当然扩大搜索邻域会增加搜索时间。根据自己的要求，采用合适的方法！