理解霍夫变换及圆检测

近日在做小球动力特性识别过程中，用到了opencv中的霍夫变换圆检测来识别图像中的小球，从而确定小球的位置。但是在调节参数的过程中不能明确各参数的真正含义，无法达到最优的识别效果。
所以想深入理解霍夫圆检测，以及了解opencv对于该算法做了哪些改进，从而更好地使用之。

谁是霍夫?

这里的霍夫不是那个CPU之父特德霍夫，这里的霍夫是Paul hough。一位喜欢发专利的数学家。

什么是霍夫变换

霍夫变换于1962年由Paul Hough 首次提出，后于1972年由Richard Duda和Peter Hart推广使用，经典霍夫变换用来检测图像中的直线，后来霍夫变换扩展到任意形状物体的识别，多为圆和椭圆。
霍夫变换 ( Hough Transform, HT)是模式识别 领域中对二值图像进行直线检测的有效方法 .参考 图 1,在标准参数化方式下 ,图像空间中的直线 l表 达为：
d= x cosθ + y sinθ ,　d≥0, 0≤θ<π
右边的图意思是在一个点上，不同的θ 对应不同的ρ，θ是穿过该点的直线的斜率，ρ是穿过该点的直线与原点的距离。这样穿过每个点作无数条直线，参数ρ和θ就能绘出一个正弦图像。
而如果1,2,3在一条直线上就可以找到一个θ使得过1，2，3点的直线的ρ很接近。也就是有图上相交的过程。而通常使用的方法是累加法，求得累加次数最多的（θ，ρ）从而确定这条直线。

如何用于圆检测

在检测圆时，附上沈子恒的解释辅助说明：（他写的十分清晰）
3.1如何表示一个圆？
与使用（r,theta）来表示一条直线相似，使用（a,b,r）来确定一个圆心为（a,b）半径为r的圆。
某个圆过点（x1,y1），则有：(x1-a1)^2 + (y1-b1)^2 = r1^2 。
那么过点（x1,y1）的所有圆可以表示为（a1(i),b1(i),r1(i)），其中r1∈（0，无穷），每一个 i 值都对应一个不同的圆，（a1(i),b1(i),r1(i)）表示了无穷多个过点（x1,y1）的圆。
3.2 如何确定多个点在同一个圆上？
如3.2中说明，过点（x1,y1）的所有圆可以表示为（a1(i),b1(i),r1(i)），过点（x2,y2）的所有圆可以表示为（a2(i),b2(i),r2(i)），过点（x3,y3）的所有圆可以表示为（a3(i),b3(i),r3(i)），如果这三个点在同一个圆上，那么存在一个值（a0,b0,r0），使得 a0 = a1(k)=a2(k)=a3(k) 且b0 = b1(k) = b2(k) = b3(k) 且r0=r1(k)=r2(k)=r3(k)，即这三个点同时在圆（a0,b0,r0)上。
从下图可以形象的看出：

首先，分析过点（x1,y1）的所有圆（a1(i),b1(i),r1(i)），当确定r1(i)时 ，（a1(i),b1(i)）的轨迹是一个以（x1,y1,r1(i)）为中心半径为r1(i)的圆。那么，所有圆（a1(i),b1(i),r1(i)）的组成了一个以（x1,y1,0）为顶点，锥角为90度的圆锥面。三个圆锥面的交点A 既是同时过这三个点的圆。累加次数最多的A点就是一个圆心。

然而这种方法计算量巨大。速度很慢，opencv中使用的是较简便的算法。

opencv是怎么实现圆检测的

OpenCV实现的是一个比标准霍夫圆变换更为灵活的检测方法——霍夫梯度法，该方法运算量相对于标准霍夫圆变换大大减少。其检测原理是依据圆心一定是在圆上的每个点的模向量上，这些圆上点模向量的交点就是圆心，霍夫梯度法的第一步就是找到这些圆心，这样三维的累加平面就又转化为二维累加平面。第二步是根据所有候选中心的边缘非0像素对其的支持程度来确定半径。
霍夫梯度法的原理是这样的。
【1】首先对图像应用边缘检测，比如用canny边缘检测。

【2】然后，对边缘图像中的每一个非零点，考虑其局部梯度，即用Sobel（）函数计算x和y方向的Sobel一阶导数得到梯度。

【3】利用得到的梯度，由斜率指定的直线上的每一个点都在累加器中被累加，这里的斜率是从一个指定的最小值到指定的最大值的距离。

【4】同时，标记边缘图像中每一个非0像素的位置。

【5】然后从二维累加器中这些点中选择候选的中心，这些中心都大于给定阈值并且大于其所有近邻。这些候选的中心按照累加值降序排列，以便于最支持像素的中心首先出现。

【6】接下来对每一个中心，考虑所有的非0像素。

【7】这些像素按照其与中心的距离排序。从到最大半径的最小距离算起，选择非0像素最支持的一条半径。
【8】如果一个中心收到边缘图像非0像素最充分的支持，并且到前期被选择的中心有足够的距离，那么它就会被保留下来。
霍夫梯度法的缺点：
<1>在霍夫梯度法中，我们使用Sobel导数来计算局部梯度，那么随之而来的假设是，其可以视作等同于一条局部切线，并这个不是一个数值稳定的做法。在大多数情况下，这样做会得到正确的结果，但或许会在输出中产生一些噪声。

<2>在边缘图像中的整个非0像素集被看做每个中心的候选部分。因此，如果把累加器的阈值设置偏低，算法将要消耗比较长的时间。第三，因为每一个中心只选择一个圆，如果有同心圆，就只能选择其中的一个。

<3>因为中心是按照其关联的累加器值的升序排列的，并且如果新的中心过于接近之前已经接受的中心的话，就不会被保留下来。且当有许多同心圆或者是近似的同心圆时，霍夫梯度法的倾向是保留最大的一个圆。可以说这是一种比较极端的做法，因为在这里默认Sobel导数会产生噪声，若是对于无穷分辨率的平滑图像而言的话，这才是必须的。

调用API

HoughCircles(image, method, dp, minDist, circles=None, param1=None, param2=None, minRadius=None, maxRadius=None)

其中参数含义如下：

image-源图像，灰度图需要先进行灰度化

method-检测方法选择参数，目前只支持 CV_HOUGH_GRADIENT这种方法，即霍夫梯度法

dp-累加器图像分辨率的反比。要理解这个，必须举例子，比如若dp=1,那么累加器的分辨率与原图像一样，如果dp=2，那么累加器的分辨率只有原图像的一半。可见，这个值越大，累加器的分辨率反而越低，所以这是个反比值。

minDist-检测到圆的圆心之间的最小距离，显然，这个值越小，伪圆可能越多，而这个值越大，则有越多的圆被漏检。

param1-函数HoughCircles是用到了Canny作边缘检测的，而Canny算法中的滞后阈值法要求设定高低两个阈值，这个参数就是设定这个高阈值的，而低阈值为这个高阈值的二分之一。阈值太高细微边缘不易被检测到。

param2 -圆心阈值参数，我前面已经说了（或者大家可以参考利用霍夫变换做直线检测的原理及OpenCV代码实现），圆心是通过投票得出的，那么多少票才算是圆心呢？这个值就是确定这个问题的。

代码案例

frame_lwpCV = cv2.resize(frame_lwpCV, None, fx=0.5, fy=0.5, interpolation=cv2.INTER_CUBIC)
cimage = cv2.cvtColor(frame_lwpCV, cv2.COLOR_RGB2GRAY)
circles = cv2.HoughCircles(cimage, cv2.HOUGH_GRADIENT, 1, 20, param1=50, param2=30, minRadius=16, maxRadius=24)
for i in circles[0, :]:
	if len(circles) > 1:
		print('不止一个圆')
	else:
		cv2.circle(frame_lwpCV, (i[0], i[1]), i[2], (0, 0, 255), 2)  # 画圆
		cv2.circle(frame_lwpCV, (i[0], i[1]), 2, (0, 0, 255), 2)  # 画圆心
		circle_number += 1
#circle.shpe为（1，1，3）（未知，圆序号，列（圆心X坐标，圆心y坐标，半径））

参考文章：

基于直径累积的霍夫变换检测圆算法_屈汉章
快速霍夫变换算法_孙丰荣
霍夫变换算法在圆心视觉定位中的应用研究_焦圣喜
Hough变换检测原理
Python+OpenCV图像处理（十五）—— 圆检测
OpenCV霍夫变换：霍夫线变换，霍夫圆变换