案例 1:霍夫变换
在计算机视觉问题解决中,我们经常在边缘信息中搜索圆圈。这个问题的特点是有很多数据点,可能来自许多不同的圆圈,存在大量噪音。
解决此问题的常用方法是 霍夫变换 https://en.wikipedia.org/wiki/Circle_Hough_Transform。基本思想是对可以通过每个点 (x, y) 的圆的证据求和。
我们创建了一个名为 Hough [a, b, r] 的整数数组,它参数化了所有可能通过您的点 (x,y) 的圆。这相当于在以(x,y)为中心的r=1平面上画一个半径为1的圆;以 (x,y) 等为中心的 r=2 平面中半径为 2 的圆。
每次通过 [a, b, r] 中的一个点绘制一个圆时,我们将相应的值加 1。有些点积累了很多的证据。这些点对应于感兴趣的圆圈。
来自 cis.rit.edu 的图片说明了其中一个 r 平面中发生的情况。
对每个点 (x,y) 执行此操作将生成指向 [a,b,r] 中与您寻找的圆圈相对应的每个点的证据。因此,只需扫描此数组即可找到具有最大证据的点。那是你的圈子。
霍夫变换示例
知道圆的半径可以将这个问题从 O(n^3) 问题减少到 O(n^2) 问题,因为只需要构建和扫描一个平面。我在对数空间中绘制半径以给出(不太准确)O(n^2 log n)算法也有很好的结果。
案例 2:圆形拟合
如果已知点位于单个圆的边界附近,和/或如果点不是很多和/或我们非常确定噪声非常小,则霍夫变换是一个糟糕的解决方案因为它是计算密集型的,内存占用很大,并且由于累加器数组的光栅特性可能不是很准确。
在这种情况下,我们可能希望通过类似于使用线性回归的线拟合技术来拟合圆。圆拟合技术的讨论可以在https://pdfs.semanticscholar.org/faac/44067f04abf10af7dd583fca0c35c5937f95.pdf找到。
这个算法的一个(相当简单的)实现如下所示。
#include <stdio.h>
#include <stdlib.h>
#include <math.h>
typedef struct {
float x;
float y;
} point;
/*
* Function using a modified least squares approach to circle fitting.
*
* Reference :
*
* Umbach, D. and Jones, K. N. "A few methods for fitting circles to data",
* IEEE Trans Instrumentation and Measurement
* vol XX(Y) 2000
*
* https://pdfs.semanticscholar.org/faac/44067f04abf10af7dd583fca0c35c5937f95.pdf
*
* NOTES
*
* The code below has not been checked for numerical stability or conditioning.
*/
int circle_fit_MLS (point P[], int n, double *x_pos, double *y_pos, double *radius)
{
int i;
double sum_x=0.0, sum_y=0.0, sum_xx=0.0, sum_xy=0.0, sum_yy=0.0, sum_xyy=0.0, sum_yxx=0.0, sum_xxx=0.0, sum_yyy=0.0;
double A, B, C, D, E;
double x2, y2, xy, F, xdif, ydif;
for (i=0; i<n; i++)
{
sum_x += P[i].x;
sum_y += P[i].y;
x2 = P[i].x * P[i].x;
y2 = P[i].y * P[i].y;
xy = P[i].x * P[i].y;
sum_xx += x2;
sum_yy += y2;
sum_xy += xy;
sum_xyy += xy*P[i].y;
sum_yxx += P[i].y*x2;
sum_xxx += x2*P[i].x;
sum_yyy += y2*P[i].y;
}
A = n * sum_xx - sum_x * sum_x;
B = n * sum_xy - sum_x * sum_y;
C = n * sum_yy - sum_y * sum_y;
D = 0.5 * ( n * (sum_xyy + sum_xxx) - sum_x * sum_yy - sum_x * sum_xx);
E = 0.5 * ( n * (sum_yxx + sum_yyy) - sum_y * sum_xx - sum_y * sum_yy);
F = A*C - B*B;
*x_pos = (D*C - B*E) / F;
*y_pos = (A*E - B*D) / F;
*radius = 0;
for (i=0; i<n; i++)
{
xdif = P[i].x - *x_pos;
ydif = P[i].y - *y_pos;
*radius += sqrt(xdif * xdif + ydif * ydif);
}
*radius /= n;
return 0;
}
下面的主程序可用于测试代码。请发回任何结果/观察/建议以改进 cmets。
int main()
{
point *P;
int n, i;
double xpos, ypos, radius;
printf ("Please enter the number of points \n> ");
scanf ("%d", &n);
P = malloc (n * sizeof(point));
for (i=0; i<n; i++)
{
printf ("x%d = ", i);
scanf ("%f", &P[i].x);
printf ("y%d = ", i);
scanf ("%f", &P[i].y);
}
circle_fit_MLS (P, n, &xpos, &ypos, &radius);
printf (" a = %f\n", xpos);
printf (" b = %f\n", ypos);
printf (" r = %f\n", radius);
}