最小二乘法拟合椭圆

转载：http://blog.sina.com.cn/s/blog_471e6c930102x96q.html

为了很好的进行椭圆方程拟合，本文先对椭圆基本知识进行复习，后进行非标准椭圆方程拟合公式推导，最后有matlab代码的实现。

1.       用最小二乘法做椭圆拟合

1.1.  椭圆标准方程

        对椭圆印象最深的就是高中时教过的，一条绳子，两个钉子，一支笔，就可以绘制出一个椭圆。固定两个钉子，让钉子之间的距离小于绳子的长度，然后用绳子的两端分别固定在两个钉子上，放一支笔在绳子的任意位置，拉紧线进行划线，画出来的完成图形就行椭圆。

然后用数学知识来解释，这两个钉子做为两个焦点，绳子的长度为2a（这里用2a，而不是a，仅是为了后续写标准的时候好看而已），钉子之间的距离为2c。为了更好的写方程，建立一个标准方程的坐标系。以两个钉子的连线中点为坐标原点，钉子连线所在的方向为x方向，绘制xoy坐标系。

其中绳子的长度为

两个钉子（焦点）的距离为

两个焦点的坐标分别为（-c,0）和（c,0）。

然后利用绳子长度的定义进行标准方程证明，具体如下：

设椭圆上任意一点M坐标为(x，y)，则根据椭圆定义，写方程如下：

即

继续证明如下：

1.2.  椭圆非标准方程

对于椭圆的非标准方程，可以写成以下形式：

以上，是推导过程。由于是在word中编辑的文档，为了省事，这里直接截图放置在博客中，分享给大家。

以下是matlab封装的最小二乘法拟合椭圆方程的函数
function ellipse_paras = ellipsefit(X,Y,n)
% 椭圆方程：x^2+Axy+By^2+Cx+Dy+E=0
% 采用最小二乘法进行拟合椭圆
% Input: X --- a vector of x measurements
%        Y --- a vector of y measurements
%        n --- the number of measurements
% author:Sigrid,2018-03-19
% M*[A B C D E]' = N
%初始化椭圆方程结果
ellipse_paras.A=0;
ellipse_paras.B=0;
ellipse_paras.C=0;
ellipse_paras.D=0;
ellipse_paras.E=0;
x=X;
y=Y;
xy=x.*y;
x2=x.*x;
y2=y.*y;
x3=x2.*x;
y3=y2.*y;
x2y=x2.*y;
xy2=x.*y2;
x4=x2.*x2;
y4=y2.*y2;
x3y=x3.*y;
xy3=x.*y3;
x2y2=x2.*y2;
% Construct M
M=[ sum(x2y2)    sum(xy3)    sum(x2y)    sum(xy2)    sum(xy);
    sum(xy3)     sum(y4)     sum(xy2)    sum(y3)     sum(y2);
    sum(x2y)     sum(xy2)    sum(x3)     sum(xy)     sum(x) ;
    sum(xy2)     sum(y3)     sum(xy)     sum(y2)     sum(y) ;
    sum(xy)      sum(y2)     sum(x)      sum(y)      n       ];
% Construct N
N = -[sum(x3y) sum(x2y2) sum(x3) sum(x2y) sum(x2)]';
P=M\N;
A=P(1);
B=P(2);
C=P(3);
D=P(4);
E=P(5);