前言
程序适用于二维离散点在平面内的圆拟合,其中离散点为nx2的数列,n>=3。以下步骤1和步骤2选择任一即可。调用函数由最小二乘法推导得来。
代码
1.读取外部文本
代码如下(示例):
clear;clc;close all;
[FileNmae,PathName]=uigetfile({'*.*'});
str=[PathName,FileNmae];
D=importdata(str);
D=D.data;
2.创建离散点
代码如下(示例):
x=[1-sqrt(2)/2 1 1+sqrt(2)/2];
y=[1-sqrt(2)/2 0 1-sqrt(2)/2];
D=[x' y'];
3.主程序
代码如下(示例):
[X0,Y0,R]=Circlefit(D);
theta=0:0.1:2*pi;
Circle1=X0+R*cos(theta);
Circle2=Y0+R*sin(theta);
plot(Circle1,Circle2,'g','linewidth',1);
xlabel('x轴');ylabel('y轴');
axis equal
grid on
4.调用函数
代码如下(示例):
function [X0,Y0,R]=Circlefit(D0)
[row,~]=size(D0);
count_x=0;count_y=0;
count_xx=0;count_yy=0;
count_xy=0;
count_xxx=0;count_xyy=0;
count_xxy=0;count_yyy=0;
count_x_y=0;
for i=1:row
count_xx=count_xx+D0(i,1)^2;
count_yy=count_yy+D0(i,2)^2;
count_x=count_x+D0(i,1);
count_y=count_y+D0(i,2);
count_xy=count_xy+D0(i,1)*D0(i,2);
count_xxx=count_xxx+D0(i,1)^3;
count_yyy=count_yyy+D0(i,2)^3;
count_xyy=count_xyy+D0(i,1)*D0(i,2)^2;
count_xxy=count_xxy+D0(i,1)^2*D0(i,2);
count_x_y=count_x_y+D0(i,1)^2+D0(i,2)^2;
end
C=row*count_xx-count_x*count_x;
D=row*count_xy-count_x*count_y;
E=row*count_xxx+row*count_xyy-count_x_y*count_x;
G=row*count_yy-count_y*count_y;
H=row*count_xxy+row*count_yyy-count_x_y*count_y;
a=(D*H-E*G)/(C*G-D^2);
b=(D*E-C*H)/(C*G-D^2);
c=-(count_x_y+a*count_x+b*count_y)/row;
X0=-a/2;Y0=-b/2;
R=sqrt(a^2+b^2-4*c)/2;
end