opencv HoughLine 理解

参考:

1. opencv HoughLine Transform Tutorial
2. https://guiqing.blog.csdn.net/article/details/8058336
3. https://www.cnblogs.com/lancer2015/p/6852488.html

opencv HoughLine结果分析

opencv中，HoughLine的结果，是找到的所有直线的集合，其中每条直线的表示是(ρ, θ)的形式，过图像左上角(0,0)点作找到的直线的垂线，ρ即是(0,0)点到垂足的距离，θ是指垂线与图像水平方向的夹角，如下图所示：

从图中可以看出，用一对这样的(ρ,θ)就能确定一条唯一的直线。ρ的计算值始终是非负值，但opencv进行了处理，当直线与Y轴的交点位于(0,0)点下方时，θ的范围是(0,π)，ρ，θ值不变，当直线与Y轴的交点位于(0,0)点上方时，θ的范围为(π,2π)，此时，将θ值变为θ-π，变换到(0,π)范围，而ρ值取-ρ。

这样变换的依据是直线的另一种表示方式，即用直线上的一个点坐标加上直线的角度来表示，此处选择的点为上述的垂足，其坐标可表示为

$$(\rho \ast cos(\theta ),\rho \ast sin(\theta ))=(-\rho \ast cos(\theta -\pi ),-\rho \ast sin(\theta -\pi ))$$

故将θ变为θ-π，ρ变为-ρ，计算得到的垂足是一样的，由于θ为法线角度，直线与法线垂直，θ变为θ-π即法线旋转180度，直线的角度不变，故变换前后表示的是同一条直线。

而且用垂线与水平轴的夹角和垂足的坐标来表示直线的话，在绘制该直线时比较方便，直接从垂足向直线两头分别移动一段距离得到两个点，即可用来进行绘制，计算方式如下：

$$cx=\rho \ast cos(\theta )$$
$$cy=\rho \ast sin(\theta )$$
$$x1=cx+1000\ast sin(\theta )$$
$$y1=cy-1000\ast cos(\theta )$$
$$x2=cx-1000\ast sin(\theta )$$
$$y2=cy+1000\ast cos(\theta )$$

计算原理如下图所示：

opencv HoughLine 原理分析

对于Hough变换，很多文章都会提到图像空间和参数空间的对应关系，参考文章[3]中有如下转换(稍加了修改)：

$$y=mx+b(图像空间)$$
$$b=-xm+y(参数空间)$$

在图像空间中，点的坐标是(x,y)，而x,y在参数空间中是方程的系数，一对系数就确定参数空间中的一条直线，因此说图像空间中的一个点，对应参数空间中的一条直线，反过来，在参数空间中，点的坐标是(m,b)，而m,b在图像空间中是系数，一对系数确定图像空间中的一条直线，因此说参数空间的一个点对应图像空间的一条直线。

但是这对于算法的意义何在？其实HoughLine的本质，是针对图像中所有可能出现的直线的情况，对于每条直线涵盖的图像上的点进行加和，然后看每条线加和的值是否达到阈值，关键就在于对每一对(m,b)，确定图像上的哪些点在此对(m,b)确定的直线上，这实现起来比较麻烦，所以转换了另一种思路，对于图像上的每个点(x,y)，用参数空间方程计算每个可能的m值对应的b值，由此得到一对(m,b)，因为参数空间方程和图像空间方程只是一个简单的形式变换，这对(x,y)和这对(m,b)是由参数空间方程计算得来的，也就必然满足图像空间方程，也就是说，该图像空间的点(x,y)必然在参数空间的点(m,b)确定的图像空间直线上，在同一条直线上的点(x,y)，计算得到的(m,b)也必然相同，这样，遍历整个图像，将每个点(x,y)的亮度值分别累加到每条经过该点的直线(m,b)上，就相当于把每条直线涵盖的图像上的点进行了累加，简化了实现，现在关键的问题在于，所有可能的直线是连续的，而算法需要进行离散化，就需要对参数空间中的参数(m,b)进行细分，而这种细分是要做到均匀取值的。

在图像空间方程中，系数b是截距，也就是x为0时，y的值，其取值范围是无穷的，对截距进行离散化无法实现，系数m是斜率，其取值范围也是无穷的，对其进行等离散化也无法实现，故难以直接用该形式的参数空间来进行HoughLine计算，通过上面对opencv HoughLine结果的分析，我们知道可以用一对(ρ,θ)来表示一条直线，且ρ的范围为[0,sqrt(w*w + h*h)]，θ的范围为[0,π)，我们将参数空间方程看作以下形式：

$$\rho =f((x,y),\theta )$$

这时，问题转换到对图像上的每个点(x,y)，θ取[0,π)范围内的每个离散值，计算对应的ρ的离散值，然后将该像素的亮度累加到(ρ,θ)表示的直线中，最后每条直线的累加值和阈值比较得到符合条件的直线。

现在的问题就是如何由((x,y), θ)计算ρ？考虑上图中的点(x,y)，有如下关系：

$$x=r\ast cos(\alpha )$$
$$y=r\ast sin(\alpha )$$
$$\rho =r\ast cos(\theta -\alpha )=r\ast cos(\theta )\ast cos(\alpha )+r\ast sin(\theta )\ast sin(\alpha )=x\ast cos(\theta )+y\ast sin(\theta )$$

这就是最终的参数空间方程。

到此为止，所有的谜题都皆晓了，最根本的原理，就是对所有可能的情况进行累加计算，再比较每个情况的累加值与阈值，玄妙的空间变换，其根本，只是为了方便算法的实现。

HoughLine python实现

# -*- coding: utf-8 -*-
import cv2 as cv
import time
import numpy as np
from matplotlib import pyplot as plt


def HoughLines(img, rho, theta, threshold, srn=1, stn=1, min_theta=0, max_theta=np.pi):
    """HoughLine 实现
    # Return
        找到的直线列表 每条直线表示为[[rho, theta]]
    # Arguments
        img: 输入的二值化图像
        rho: 极径分辨率
        theta: 极角分辨率
        threshold: 阈值
        srn: 多尺度计算时，精细极径分辨率=rho/srn   此处直接实现为精细版
        stn: 多尺度计算时，精细极角分辨率=theta/stn 此处直接实现为精细版
        min_theta: 最小极角
        max_theta: 最大极角
    """
    rho_acc = rho / srn
    theta_acc = theta / stn
    # 小于180度的 rho为正
    # 大于180度的 rho取负 角度-180
    # rho 以 rho_range//2 为0点
    rho_range = int((img.shape[0] + img.shape[1]) / rho_acc) * 2
    # 直接计算从[0,pi)范围, 再最后取结果时再过滤
    pi_range = int(np.pi / theta_acc)

    hough = np.zeros([pi_range, rho_range], dtype='int32')

    for y in range(img.shape[0]):
        for x in range(img.shape[1]):
            if img[y, x] < 128:
                continue
            for theta in range(pi_range*2):
                ct = np.cos(theta * theta_acc)
                st = np.sin(theta * theta_acc)
                # 参数空间方程
                rho = int((y * st + x * ct) / rho_acc)
                if theta < pi_range:
                    the = theta
                    rho = rho_range // 2 + rho
                    hough[the, rho] += 1
                else:
                    the = theta - pi_range
                    rho = rho_range // 2 - rho
                    hough[the, rho] += 1
    # plt.imshow(np.array(hough))
    # plt.show()
    res = []
    for theta in range(pi_range * 2):
        if theta * theta_acc < min_theta or theta * theta_acc > max_theta:
            continue
        if theta >= pi_range:
            theta -= pi_range
        for rho in range(rho_range):
            if hough[theta, rho] > threshold:
                res.append([[(rho - rho_range // 2) * rho_acc, theta * theta_acc]])

    return np.array(res)


def test(img, method=0):
    gray = cv.cvtColor(img, cv.COLOR_BGR2GRAY)
    if method == 0:
        lines = HoughLines(gray, 1, np.pi / 90, 250)
    else:
        lines = cv.HoughLines(gray, 1, np.pi / 90, 250)
    print('lines:', np.array(lines).shape)
    finds = []
    for line in lines:
        rho, theta = line[0]
        ct = np.cos(theta)
        st = np.sin(theta)
        x0 = ct * rho
        y0 = st * rho
        # x0, y0 为垂足坐标 下面的两个点是垂足延直线向两头各移动1000个单位
        x1 = int(x0 - 1000 * st)
        y1 = int(y0 + 1000 * ct)
        x2 = int(x0 + 1000 * st)
        y2 = int(y0 - 1000 * ct)

        # 过滤重合的线
        for find in finds:
            a1 = np.array([rho, theta * 180 / np.pi])
            a2 = np.array([x1, y1, x2, y2])
            b1 = np.array(find[:2])
            b2 = np.array(find[2:])
            d1 = np.abs(a1 - b1)
            d2 = np.abs(a2 - b2)
            if np.all(d1 < 20) or np.all(d2 < 20):
                # print('line', line[0], [x1,y1,x2,y2], 'find', find, 'overlap')
                break
        else:
            finds.append([rho, theta * 180 / np.pi, x1, y1, x2, y2])
            # print('line', line[0], [x1,y1,x2,y2])
            cv.line(img, (x1, y1), (x2, y2), (0,0,255), 2)
    print('finds', len(finds))
    return img


def main():
    img = cv.imread('d:/0.png')

    time1 = time.time()
    t1 = test(img, method=0)
    time2 = time.time()
    print('self time:', time2 - time1)

    time1 = time.time()
    t2 = test(img, method=1)
    time2 = time.time()
    print('opencv time:', time2 - time1)

    t = np.hstack([t1, t2])
    cv.imshow('self <-----> opencv', t)
    cv.waitKey()


if __name__ == '__main__':
    main()

对比手动实现版本和opencv版本，结果基本一致，不过python实现版本的效率相比opencv差了很多，如下图所示，左边是手动实现版本，右边是opencv版本：

以上就是本人对opencv HoughLine的理解，由于水平有限，可能有理解错误的地方，欢迎大家交流和讨论。