【图像处理】（2.1）霍夫直线变换

1、算法思想

边缘检测比如canny算子可以识别出图像的边缘，但是实际中由于噪声和光照不均匀等因素，很多情况下获得的边缘点是不连续的，必须通过边缘连接将他们转换为有意义的边缘。Hough变化是一个重要的检测间断点边界形状的方法，它通过将图像坐标空间变化到参数空间来实现直线和曲线的拟合。

霍夫变换于1962年由Paul Hough 首次提出，后于1972年由Richard Duda和Peter Hart推广使用，经典霍夫变换用来检测图像中的直线，后来霍夫变换扩展到任意形状物体的识别，多为圆和椭圆。

Hough变换是图像处理中从图像中识别几何形状的基本方法之一。Hough直线检测的基本原理在于利用点与线的对偶性，在我们的直线检测任务中，即图像空间中的直线与参数空间中的点是一一对应的，参数空间中的直线与图像空间中的点也是一一对应的。这意味着我们可以得出两个非常有用的结论：

1）图像空间中的每条直线在参数空间中都对应着单独一个点来表示；

2）图像空间中的直线上任何一部分线段在参数空间对应的是同一个点。

因此Hough直线检测算法就是把在图像空间中的直线检测问题转换到参数空间中对点的检测问题，通过在参数空间里寻找峰值来完成直线检测任务，也即把检测整体特性转化为检测局部特性。

2、算法原理

1）图像空间和参数空间

霍夫变换的数学理解是“换位思考”，比如一条直线y=k*x+q有两个参数，在给定坐标系下，这条直线就可以用k和q进行完整的表述。如果我们把x和y看作参数，把a和b看作变量的话，那么图像空间下的坐标点（x1，y1）对应着参数空间里的一条直线q=-x1*k+y1, 图像空间直线上的点（x1，y1）就是参数空间的斜率和截距，其中k，q为参数空间的自变量。

2）参数空间转换过程

下面用不同空间下的点和线的变换过程示例说明。

一条直线可由两个点A=(X1,Y1)和B=(X2,Y2)确定(笛卡尔坐标)。

 另一方面，y=kx+q也可以写成关于（k,q）的函数表达式（霍夫空间）：

对应的变换可以通过图形直观表示：

变换后的空间成为霍夫空间。即：笛卡尔坐标系中一条直线，对应霍夫空间的一个点。

反过来同样成立（霍夫空间的一条直线，对应笛卡尔坐标系的一个点）：

再来看看A、B两个点，对应霍夫空间的情形：

 再看一下三个点共线的情况：

可以看出如果笛卡尔坐标系的点共线，这些点在霍夫空间对应的直线交于一点：这也是必然，共线只有一种取值可能。

如果不止一条直线呢？再看看多个点的情况（有两条直线）：

其实（3，2）与（4，1）也可以组成直线，只不过它有两个点确定，而图中A、B两点是由三条直线汇成，这也是霍夫变换的后处理的基本方式：选择由尽可能多直线汇成的点。

  到这里问题似乎解决了，已经完成了霍夫变换的求解，但是如果像下图这种情况呢？

k=∞是不方便表示的，而且q怎么取值呢，这样不是办法。这里就需要给k和q换一种表示方式，让xy对应的笛卡尔坐标系与新的参数空间进行转换。

对于这里我有一个补充，这里有人对这种或者表现形式是否使用了极坐标产生了不同的理解（这里的和是等价的符号，由于公式是由其他博主那搬运，本人没有进行更改）：

（1）在“会飞的吴克”课程3.hough transform(霍夫变换直线检测)_哔哩哔哩_bilibili中，该视频作者认为这并没有引入极坐标，而是仅仅引入了新的变量和新的表示方法，逻辑如下图所示：

“会飞的吴克”课程中的介绍更容易理解，因为他是从编写代码的角度出发，来剖析原理。

 

 后续的点与线的变换如上图所示。

（2）考虑到上述问题，有人将笛卡尔坐标系换为：极坐标表示。（参考文件里大佬博客里面的图错了，下图是正确的极坐标表示方法，并且给出了辅助线几何解释）

在极坐标系下，其实是一样的：极坐标的点→霍夫空间的直线，只不过霍夫空间不再是[k,q]的参数，而是[ρ, θ]的参数，给出对比图：

这里我针对上面这2种不同理解说一说我的看法，我认为没有引入极坐标点。

如果认为引入了极坐标，上述的第（2）种方式极坐标的点→霍夫空间的直线就有问题，极坐标的点是有和theta表示，怎么对应上的笛卡尔坐标系中的xy呢，而且参数空间变换产生了变化，与实际意义产生了脱离。 只有第（1）种方式中描述的，通过引入一个角和一个边来替换k和q，生成新的公式，这里唯一的缺憾就是先引入的中的sin(theta)有可能为0，也就是theta取90度时，该公式的条件就是theta不能为90度。通过公式变换和移项得到，该公式的theta可以取0，这也符合我们的预期目标，就是上述提到的斜率无穷大，但这变相的扩宽了条件，这种操作方式是不是合理。

从上面可以看到，参数空间的每个点(ρ,θ)都对应了图像空间的一条直线，或者说图像空间的一个点在参数空间中就对应为一条曲线。这样就把在图像空间中检测直线的问题转化为在参数空间中找通过点(ρ,θ)的最多正弦曲线数的问题。霍夫空间中，曲线的交点次数越多，所代表的参数越确定，画出的图形越饱满，其实本质就是在笛卡尔坐标系中的多个点都经过参数空间中点(ρ,θ)对应在笛卡尔空间的直线。

霍夫直线检测就是把图像空间中的直线变换到参数空间中的点，通过统计特性来解决检测问题。具体来说，如果一幅图像中的像素构成一条直线，那么这些像素坐标值（x, y）在参数空间对应的曲线一定相交于一个点，所以我们只需要将图像中的所有像素点（坐标值）变换成参数空间的曲线，并在参数空间检测曲线交点就可以确定直线了。

3）那么怎么确定交点呢？

从数学角度出发，2个未知数2个方程可以求出交点。但是针对图像需要做一点小变动，由于theta表示与直线垂直的线与x轴的夹角，那么他的取值范围就是-pi到pi，但由于是直线，只需要选取0到pi就可以表示所有直线。r表示原点到直线的距离，他的取值范围应该是0到sqrt(图像宽的平方+图像高的平方)。

针对这个有限范围，我们可以让theta在0到pi取样，假设步长为1度，那么可以去180个点，分别为[1度，2度...，180度]，针对每一个点都会有一个r值，这时的r值可能不能像根据数学公式计算出的结果那样没有偏差，所以只要两者差值小于一定的范围就认为该点是两条曲线的交点。

基于以上思想，我们需要进一步方便计算，如下图所示，假设对theta取得间隔为45，那么可以的取值为[45,90,135,180]，假设r最长为9，间隔步长为1，这里对于r的间隔步长表示上面提及的差值范围（这里的theta和r的间隔步长可以自己设置，我这里只是为了方便举例）。对于图像中的每一个点都带入公式中，此时的x和y就变成了常量，再对theta分别取[45,90,135,180]这几个值分别得到对应r值，只要r值落在对应的区间范围内，在下图矩阵对应位置加1。如此重复过程，那么当矩阵中值越大时，表明越多原图像中的点落在由theta和r构成的直线上。

上面的理论解释我觉得没什么问题，我自己也会有其他的一个更直观的解释方法。

如上图所示，共有4个点，每个点根据上述举例的[45,90,135,180]4个theta值得到对应的r值，有了r和theta那么就能确定一条直线，如下所示每个点都有4条直线，分别就是对应的theta值和r值。观察A和B他两在theta等于90时，A和B在横向方向上在同一条直线上，假设此时的r在2-3这个范围内，那么如下图所示A和B的theta和r值都在红框所在区域，那么他的值就从0变成了2。如此类推即可。

 

4）非极大值抑制

当我们取得theta角度间隔比较小时，那么就会出现如上图所示情况，最终确定的直线在一个较小的角度内很密集，这些直线我们称之为候选直线。而我们需要的是一条直线，而不是多条直线，，也就是从这几条直线中选出一条作为最终的直线。

那就是在theta相近的直线内，找一个投票数最多的一条之间就行，方法就是在候选直线中先选一条作为处理直线，在其附近画一个窗口，认定窗口内的就是与他重叠度较高的直线，如下图的绿框，如果在绿框内处理直线为最高得票数，否则该处理直线被删除。接着从候选直线中选取下一条直线，重复上面的过程。知道所有候选直线都被处理完毕即可得到最终的曲线集合。

5）总结：使用霍夫变换检测直线具体步骤

3、代码测试

实现代码（3.hough transform(霍夫变换直线检测)_哔哩哔哩_bilibili）

#coding:utf-8

import numpy as np
import matplotlib.pyplot as plt
import os
import cv2

def lines_detector_hough(edge,ThetaDim = None,DistStep = None,threshold = None,halfThetaWindowSize = 2,halfDistWindowSize = None):
    '''
    :param edge: 经过边缘检测得到的二值图
    :param ThetaDim: hough空间中theta轴的刻度数量(将[0,pi)均分为多少份),反应theta轴的粒度,越大粒度越细
    :param DistStep: hough空间中dist轴的划分粒度,即dist轴的最小单位长度
    :param threshold: 投票表决认定存在直线的起始阈值
    :return: 返回检测出的所有直线的参数(theta,dist)
    @author: bilibili-会飞的吴克
    '''
    imgsize = edge.shape
    if ThetaDim == None:
        ThetaDim = 90
    if DistStep == None:
        DistStep = 1
    MaxDist = np.sqrt(imgsize[0]**2 + imgsize[1]**2)
    DistDim = int(np.ceil(MaxDist/DistStep))

    if halfDistWindowSize == None:
        halfDistWindowSize = int(DistDim/50)
    accumulator = np.zeros((ThetaDim,DistDim)) # theta的范围是[0,pi). 在这里将[0,pi)进行了线性映射.类似的,也对Dist轴进行了线性映射

    sinTheta = [np.sin(t*np.pi/ThetaDim) for t in range(ThetaDim)]
    cosTheta = [np.cos(t*np.pi/ThetaDim) for t in range(ThetaDim)]

    for i in range(imgsize[0]):
        for j in range(imgsize[1]):
            if not edge[i,j] == 0:
                for k in range(ThetaDim):
                    accumulator[k][int(round((i*cosTheta[k]+j*sinTheta[k])*DistDim/MaxDist))] += 1

    M = accumulator.max()

    if threshold == None:
        threshold = int(M*2.3875/10)
    result = np.array(np.where(accumulator > threshold)) # 阈值化
    temp = [[],[]]
    for i in range(result.shape[1]):
        eight_neiborhood = accumulator[max(0, result[0,i] - halfThetaWindowSize + 1):min(result[0,i] + halfThetaWindowSize, accumulator.shape[0]),
                           max(0, result[1,i] - halfDistWindowSize + 1):min(result[1,i] + halfDistWindowSize, accumulator.shape[1])]
        if (accumulator[result[0,i],result[1,i]] >= eight_neiborhood).all():
            temp[0].append(result[0,i])
            temp[1].append(result[1,i])

    result = np.array(temp)    # 非极大值抑制

    result = result.astype(np.float64)
    result[0] = result[0]*np.pi/ThetaDim
    result[1] = result[1]*MaxDist/DistDim

    return result

def drawLines(lines,edge,color = (255,0,0),err = 3):
    if len(edge.shape) == 2:
        result = np.dstack((edge,edge,edge))
    else:
        result = edge
    Cos = np.cos(lines[0])
    Sin = np.sin(lines[0])

    for i in range(edge.shape[0]):
        for j in range(edge.shape[1]):
            e = np.abs(lines[1] - i*Cos - j*Sin)
            if (e < err).any():
                result[i,j] = color

    return result


if __name__=='__main__':
    pic_path = './HoughImg/'
    pics = os.listdir(pic_path)

    for i in pics:
        if i[-5:] == '.jpeg' or i[-4:] == '.jpg':
            img = plt.imread(pic_path+i)
            blurred = cv2.GaussianBlur(img, (3, 3), 0)
            plt.imshow(blurred,cmap='gray')
            plt.axis('off')
            plt.show()

            if not len(blurred.shape) == 2:
                gray = cv2.cvtColor(blurred, cv2.COLOR_RGB2GRAY)
            else:
                gray = blurred
            edge = cv2.Canny(gray, 50, 150)   #  二值图 (0 或 255) 得到 canny边缘检测的结果


            lines = lines_detector_hough(edge)
            final_img = drawLines(lines,blurred)
            

            plt.imshow(final_img,cmap='gray')
            plt.axis('off')
            plt.show()

opencv代码

函数原型：

CV_EXPORTS_W void HoughLines( InputArray image, OutputArray lines,
                              double rho, double theta, int threshold,
                              double srn = 0, double stn = 0,
                              double min_theta = 0, double max_theta = CV_PI );

 测试代码如下：

int main() {
    Mat src_img;
    Mat dst_img, cdst, cdst2;

    src_img = imread("D:\\WORK\\5.OpenCV\\LeanOpenCV\\pic_src\\pic18.bmp");
    imshow("原图", src_img);
    
    Canny(src_img, dst_img, 50, 200, 3);
    imshow("Canny图", dst_img);


    cdst = src_img.clone();
    cdst2 = dst_img.clone();
    vector<Vec2f> lines;
    HoughLines(dst_img, lines, 1, CV_PI / 180, 100, 0, 0);

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

    imshow("detected lines", cdst);
    imshow("detected lines2", cdst2);

    waitKey(0);
}

测试效果图如下，canny边缘检测有不连续的边缘，霍夫变换直线检测可以连接不连续的直线边缘。

4、参考文献

1、《数字图像处理与机器视觉》

第二版。 张铮、徐超、任淑霞、韩海玲等编著。

2、霍夫变换直线检测（Line Detection）原理及示例

霍夫变换直线检测（Line Detection）原理及示例_leonardohaig的博客-CSDN博客_霍夫直线检测

3、霍夫变换

霍夫变换 - 疯狂奔跑 - 博客园

4、霍夫变换-----特征提取

霍夫变换-----特征提取_程序猿的视界-CSDN博客_霍夫变换曲线

5、霍夫线变换

霍夫线变换 — OpenCV 2.3.2 documentation

6、霍夫圆变换

霍夫圆变换 — OpenCV 2.3.2 documentation

7、经典霍夫变换（Hough Transform）

经典霍夫变换（Hough Transform）_YuYunTan的专栏-CSDN博客_hough变换论文

8、霍夫变换（Hough Transform）的原理以及代码（Matlab&C）实现

霍夫变换（Hough Transform）的原理以及代码（Matlab&C）实现_ljwcdtj的博客-CSDN博客_hough transform

9、 hough变换算法 - 啊哈彭 - 博客园 

10、3.hough transform(霍夫变换直线检测)_哔哩哔哩_bilibili