最小二乘法拟合圆公式推导及vc实现[r]

         最小二乘法(least squares analysis)是一种 数学 优化 技术,它通过 最小化 误差 的平方和找到一组数据的最佳 函数 匹配。 最小二乘法是用最简的方法求得一些绝对不可知的真值,而令误差平方之和为最小。 最小二乘法通常用于 曲线拟合 (least squares fitting) 。这里有 拟合圆曲线 的公式推导过程 和 vc实现。








VC实现的代码:
void CViewActionImageTool::LeastSquaresFitting()
{
    
if (m_nNum<3)
    
{
        
return;
    }


    
int i=0;

    
double X1=0;
    
double Y1=0;
    
double X2=0;
    
double Y2=0;
    
double X3=0;
    
double Y3=0;
    
double X1Y1=0;
    
double X1Y2=0;
    
double X2Y1=0;

    
for (i=0;i<m_nNum;i++)
    
{
        X1 
= X1 + m_points[i].x;
        Y1 
= Y1 + m_points[i].y;
        X2 
= X2 + m_points[i].x*m_points[i].x;
        Y2 
= Y2 + m_points[i].y*m_points[i].y;
        X3 
= X3 + m_points[i].x*m_points[i].x*m_points[i].x;
        Y3 
= Y3 + m_points[i].y*m_points[i].y*m_points[i].y;
        X1Y1 
= X1Y1 + m_points[i].x*m_points[i].y;
        X1Y2 
= X1Y2 + m_points[i].x*m_points[i].y*m_points[i].y;
        X2Y1 
= X2Y1 + m_points[i].x*m_points[i].x*m_points[i].y;
    }


    
double C,D,E,G,H,N;
    
double a,b,c;
    N 
= m_nNum;
    C 
= N*X2 - X1*X1;
    D 
= N*X1Y1 - X1*Y1;
    E 
= N*X3 + N*X1Y2 - (X2+Y2)*X1;
    G 
= N*Y2 - Y1*Y1;
    H 
= N*X2Y1 + N*Y3 - (X2+Y2)*Y1;
    a 
= (H*D-E*G)/(C*G-D*D);
    b 
= (H*C-E*D)/(D*D-G*C);
    c 
= -(a*X1 + b*Y1 + X2 + Y2)/N;

    
double A,B,R;
    A 
= a/(-2);
    B 
= b/(-2);
    R 
= sqrt(a*a+b*b-4*c)/2;

    m_fCenterX 
= A;
    m_fCenterY 
= B;
    m_fRadius 
= R;

    
return;
}


工程下载
编译运行后随便打开一个图片,当然最好是全白的图片,然后就点吧,大于三个点后就会开始拟合。红线画的圆为拟合的圆,深蓝的点为鼠标点击设置的样本点。单击鼠标右键清空样本集。

路漫漫其修远兮 吾将上下而求索

posted @ 2006-10-09 21:37 .Live 阅读(4730) 评论(25)  编辑 收藏 网摘 所属分类: 图形图像

  回复  引用  查看    
#1楼 2006-10-09 22:53 | GouGou      
VC6.0的项目~
  回复  引用    
#2楼 2007-04-05 09:53 | 路不平 [未注册用户]
公式错误!

最小二乘法的意义是距离差值的平方(di-R)^2,不是距离平方差的平方(di^2-R^2)^2

两者是不同的!事实上

(di^2-R^2)=(di-R)(di+R)

因此,以上结果将使圆半径“趋小”!

  回复  引用  查看    
#3楼 [楼主]2007-04-05 13:37 | .Live      
@路不平
其实是一样的。

最小化 误差 的平方和,将距离平方差认为是 误差,可以表征和真实值得误差大小,使之最小也可以得到一个相对较优的拟合值。同时简化公式的推导。

如果使用距离的平方差,则平方和中会出现开根号,公式就很难推导了。
  回复  引用    
#4楼 2007-04-06 09:18 | 路不平 [未注册用户]
回:dotLive

我编制过关于圆最小二乘的程序。起先我以为两者是一样的,就是:评价(di-R)^2与评价(di^2-R^2)^2是一样的。
但正如你所说,后者是把R^2作为真值进行考虑的。当然如果散点恰能构成一个圆,两者是一样的。因为:
(a^2-b^2)=(a+b)(a-b)=0
但我实际编程后(采用最速下降法)发现,后者的半径偏小(散点越不接近于圆越明显)!为什么呢?因为,最小二乘的几何原理考察的是“距离”的误差(di-R)。那么变成平方差 (di^2-R^2)后,事实上是给(di-R)乘上了一个(di+R约=2R)的权重,最终实际的结果,当然是半径偏小。
事实上,你可以通过编程发现两者的差异!
当然上述推导依然是有价值的,因为这个结果作为迭代法的初值是再合适不过了!
实际的推导是不“可解”的,否则,《数学手册》就会出现相应的方法,就象线性回归和抛物线回归一样!



  回复  引用  查看    
#5楼 [楼主]2007-04-06 10:30 | .Live      
@路不平
你的分析是对的,对我很有启发。我没有考虑到偏差的程度是如何变化的。你的分析应该说验证了这个推导在整体圆检测算法中的可行性。
事实上这个拟合公式确实迭代检测圆的一个步骤,我采取的方法是每次迭代删除方差大于某个阈值的点,认为它们是非圆上的点。
要使用这个迭代拟合的方法,还需要前面粗检测来限定范围,以保证进行拟合的点在圆周的一个特定范围内;最后删除点的迭代次数,也就是迭代停止条件也需要考虑和分析。最后我们希望是可以得到全部都是圆上的点,或者都是离真实边缘很接近的点,那么拟合的圆稍微有所偏差也是可以接受的。
还有一点是,这个检测的算法是假设在很高的分辨率下进行的,因此需要采用一些牺牲少量精度,大量降低计算量的方法,检测精度已经由高分辨的采集装置保证了一部分。
  回复  引用    
#6楼 2007-04-20 16:24 | Merlin [未注册用户]
不错。我在工程上验证了这个算法是可行的。
我给你了我邮件地址,
  回复  引用    
#7楼 2007-05-13 17:20 | 谁来帮帮我 [未注册用户]
1 熟悉圆度测量的基本知识和工作原理 2。编制计算机软件,实现以下基本功能 用最下二乘法来评定圆度误差。
设计圆度测量的偏心校正算法,并在计算机上予以仿真,验证。
根据实际测试数据计算某一个实际工件的偏心值和圆度误差
3设计程序流程图和编制软件
``两位大哥 都那么厉害 `谁能做这个?
我的毕业设计 我给钱
  回复  引用    
#8楼 2007-05-13 21:51 | 谁来帮帮我 [未注册用户]
我的QQ:55847142
  回复  引用    
#9楼 2007-10-11 11:29 | vb [未注册用户]
楼主,那您最后的解决方案是什么呢?
@.Live

  回复  引用  查看    
#10楼 [楼主]2007-10-15 08:58 | .Live      
@vb
可以在07年11月份的仪器仪表学报上参考《多分辨率快速圆检测在抑菌圈测量系统中的应用》,或者检索学位论文《基于机器视觉的抗生素残留量检测技术研究》
  回复  引用    
#11楼 2007-11-10 16:26 | ShiZD [未注册用户]
请问楼主和“路不平”,楼主列出的算法存在拟和出来的半径R可能偏小的问题,那么拟和得到的圆心坐标A和B是否准确呢?如果也存在偏差,偏差的趋势是怎样的?
  回复  引用  查看    
#12楼 [楼主]2007-11-13 09:19 | .Live      
@ShiZD
圆心的偏差怎么衡量呢,半径有大小问题,圆心偏上下左右都没有衡量的方向的问题,而且最小二乘法就是求得一个最小化误差的方法,如果样本点确定的圆是唯一的,那么得到的结果也肯定是精确的。当样本确定的圆是不确定的,那么要分析算法结果就需要确定一个比较的标准真值。

这个算法的结果半径“偏小”是指和使用 距离差值的平方 的方法得到结果比较而言。因为使用 距离差值的平方 的方法显然会科学一点,就认为它更接近真值。

就实践上来说,可以使用人工确定的真值来统计结果的准确率,在相对固定的应用环境下达到一定标准就可以认为是可行的了。
  回复  引用    
#13楼 2008-02-23 13:16 | 球王 [未注册用户]
非常不错 现在正需要这方面的源码参考学习 谢谢
  回复  引用    
#14楼 2008-02-23 14:18 | 球王 [未注册用户]
请教一下:由于不是数学专业,对最小二乘法不是很懂,请问拟合椭圆和拟合圆有多大不同呢?最好能有源码来学习一下 谢谢
  回复  引用    
#15楼 2008-03-05 13:24 | dotlive [未注册用户]
@球王
@球王
参数个数不同,你要解的方程组就不同,跟阶数也有关,椭圆要5个参数,更复杂了点,不好解。能解出来的话,源码没什么不同的。建议查找相关的论文看看。
  回复  引用    
#16楼 2008-03-20 19:32 | 和规范的话 [未注册用户]
谢谢了,两位达淫!
  回复  引用    
#17楼 2008-05-27 09:06 | Pterosaur [未注册用户]
只想确认的一点是,这个拟和算法对于一个确定的目标圆,是不是会无穷趋近于真实值?
谢谢楼主了!
  回复  引用    
#18楼 2008-05-28 21:39 | Linchang [未注册用户]
其实我觉得楼主的做法是可以的,你可以把平方的差认为是面积之差,所以拟合出来的使面积差的平方和最小,应该还是说得通的


  回复  引用  查看    
#19楼 [楼主]2008-05-31 00:43 | .Live      
@Pterosaur

对于在同一个圆上的点,只要不在一条直上的3个点就可以拟合出圆心和半径了

  回复  引用    
#20楼 2008-06-01 11:58 | Pterosaur [未注册用户]
楼主理解错我的意思了,我是想知道这个估计的无偏性和一致性...
按照仿真的数据来看,半径拟和的无偏性不能保证,一致性不清楚.
还没有具体的证明
  回复  引用  查看    
#21楼 [楼主]2008-06-15 12:26 | .Live      
@Pterosaur

这个真的很难回答,因为 “确定的目标圆”带有主观因素的,不同的应用领域对目标的要求是不同的,在程序中或许可以控制一些参数来使结果偏向应用的要求。但是从公式上去证明就有点困难了。我也没有很多关注公式,更多关心应用吧。
  回复  引用    
#22楼 2008-06-16 09:13 | Pterosaur [未注册用户]
换个说法,我已经知道了这种拟和算法的半径是有偏估计,偏差与拟和数据所带有的噪声的方差有关(噪声的均值为0),但是我现在不清楚加入的噪声具体有多大,在这种情况下,这个半径该如何估计?
一个很重要的问题是噪声不能忽略.
  回复  引用  查看    
#23楼 [楼主]2008-06-21 22:16 | .Live      
@Pterosaur
明白了,"噪声不能忽略"是关键问题,由此带来偏差。但在我的所碰到的实际应用中,前期做了很多工作来消除噪声。而且更重要的是,最后的进行统计的并不是点,而是线,而且是在很低的分辨率下经过HOUGH算法过滤后的线,然后再在高分辨率进行一个数据的精确统计。只是因为图像数据是离散的,所以表现出来的都是点坐标而已。因此在我的应用实例中,“噪声被忽略” :-|。
噪声确实会产生偏差,如果能分析噪声源的规律,或许可以得到修正的公式。直观的可以看到,如果所有噪声点均匀分布,那么对圆的影响还跟圆的空间位置有关的。因此,还是要跟实际应用结合吧。
  回复  引用    
#24楼 2008-06-23 22:35 | Pterosaur [未注册用户]
感谢楼主的热情回复
的确,消噪可以带来很好的效果,而且我在应用的前期也进行了很多的消噪过程.但是精度仍然达不到工程的需求.最大的问题就是出在半径的估计上,因为这个半径的估计误差在后期的运算过程中会被放大大约100倍(并不是估计出半径就结束的).
很头疼的问题是,即使引入的噪声是高斯白噪声,我也不清楚对半径的估计影响会有多大.
另外,如果拟和这个圆的所有点被集中在一个很小的范围之内(大于45度),拟和结果会很差,这个问题该如何解决?
  回复  引用    
#25楼 2008-08-30 14:57 | 小金鱼 [未注册用户]
--引用--------------------------------------------------
路不平: 回:dotLive
<br>
<br> 我编制过关于圆最小二乘的程序。起先我以为两者是一样的,就是:评价(di-R)^2与评价(di^2-R^2)^2是一样的。
<br> 但正如你所说,后者是把R^2作为真值进行考虑的。当然如果散点恰能构成一个圆,两者是一样的。因为:
<br> (a^2-b^2)=(a+b)(a-b)=0
<br> 但我实际编程后(采用最速下降法)发现,后者的半径偏小(散点越不接近于圆越明显)!为什么呢?因为,最小二乘的几何原理考察的是“距离”的误差(di-R)。那么变成平方差 (di^2-R^2)后,事实上是给(di-R)乘上了一个(di+R约=2R)的权重,最终实际的结果,当然是半径偏小。
<br> 事实上,你可以通过编程发现两者的差异!
<br> 当然上述推导依然是有价值的,因为这个结果作为迭代法的初值是再合适不过了!
<br> 实际的推导是不“可解”的,否则,《数学手册》就会出现相应的方法,就象线性回归和抛物线回归一样!
<br>
<br>
<br>
--------------------------------------------------------
请问“路不平”所说的:“但我实际编程后(采用最速下降法)发现,”中最速下降法指的是什么?

标题  
姓名  
主页
Email (博主才能看到) 
验证码 *  看不清,换一张 [登录][注册]
内容(请不要发表任何与政治相关的内容)  
  登录  使用高级评论  新用户注册  返回页首  恢复上次提交      
该文被作者在 2006-10-12 17:30 编辑过
"五向定位"职业成长路线公开课(上海、南京、大连)
Google站内搜索


相关链接: