MATLAB圆拟合

查过很多材料，圆拟合的代码主要有两种，但是其实原理相同只是从两种方式表达而已。开始时只查到代码，并不知其原理也没查到，偶然一个机会下载一个资源得知其原理。

以下主要从三个方面总结MATLAB圆拟合

代码（含圆拟合、画圆）

原理

注意(这个最重要，见证了我沙雕的智商~)

1.代码

第一种：利用行列式，短小精悍

function [xc,yc,R,a] = circfit(x,y)
%圆拟合函数
%CIRCFIT Fits a circle in x,y plane
% [XC, YC, R, A] = CIRCFIT(X,Y)
% Result is center point (yc,xc) and radius R.A is an
% optional output describing the circle’s equation:
%
% x^2+y^2+a(1)*x+a(2)*y+a(3)=0
% by Bucher izhak 25/oct/1991
n=length(x); xx=x.*x; yy=y.*y; xy=x.*y;
A=[sum(x) sum(y) n;sum(xy) sum(yy)...
sum(y);sum(xx) sum(xy) sum(x)];
B=[-sum(xx+yy) ; -sum(xx.*y+yy.*y) ; -sum(xx.*x+xy.*y)];
a=A\B;            %x = A\B 用来求解线性方程 A*x = B.  A 和 B 的行数一致.
xc = -.5*a(1);
yc = -.5*a(2);
R = sqrt((a(1)^2+a(2)^2)/4-a(3));
 theta=0:0.1:2*pi;  
    Circle1=xc+R*cos(theta);  
    Circle2=yc+R*sin(theta);   
    plot(Circle1,Circle2,'g','linewidth',1);  
    axis equal  
    end  

第二种：循环计算

function [R,A,B]=circ(x,y)
%circfit函数的易看版本，原理相同，只是处理数据用循环
x1 = 0;
x2 = 0;
x3 = 0;
y1 = 0;
y2 = 0;
y3 = 0;
x1y1 = 0;
x1y2 = 0;
x2y1 = 0;
N=length(x);
for i = 1 : N
x1 = x1 + x(i);
x2 = x2 + x(i)*x(i);
x3 = x3 + x(i)*x(i)*x(i);
y1 = y1 + y(i);
y2 = y2 + y(i)*y(i);
y3 = y3 + y(i)*y(i)*y(i);
x1y1 = x1y1 + x(i)*y(i);
x1y2 = x1y2 + x(i)*y(i)*y(i);
x2y1 = x2y1 + x(i)*x(i)*y(i);
end
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;
A = a/(-2); %x 坐标
B = b/(-2); %y 坐标
R = sqrt(a * a + b * b - 4 * c)/2;
end

如何画圆：

[R,A,B]=circ(X_new,Y_new);
%画结果图
theta=0:pi/20:2*pi;  
    Circle1=A+R*cos(theta);  
    Circle2=B+R*sin(theta);
    
subplot(122);imshow(Image);hold on
 plot(Circle1,Circle2,'m','linewidth',1.5);hold off

2.原理

最小二乘法(least squares analysis)是一种数学优化技术，它通过

最小化误差的平方和找到一组数据的最佳函数匹配。最小二乘法是用最

简的方法求得一些绝对不可知的真值，而令误差平方之和为最小。 最

小二乘法通常用于曲线拟合 (least squares fitting) 。

这里有拟合圆曲线 的公式推导过程和vc实现。

 

 

 

注意

MATLAB的直角坐标系起点在左下角，而图像的坐标系起点在左上角。如果为了表示图像在直角坐标系的位置，那么要把图像倒过来。代码： seg_smooth_1=flipud(seg_smooth);两种坐标系的图像不能在同一坐标系混用，比方说将图像颠倒以便正着在直角坐标系显示苹果轮廓，那么由颠倒图像得出的拟合圆在数字图像坐标系显示会出现错位，如图3，刚刚好圆和苹果是对称的，惊不惊喜，意不意外！之前减二值图代码里就出现这个问题，我不知道怎么回事，归因于玄学，看吧还是学艺不精~三个月不到就想明白了，所以下次代码出错要好好理解。