Ellipse Fit (Direct method)

Fits an ellipse to a set of points on a plane; returns the coefficients of the ellipse's equation

Description

This is a fast and non-iterative ellipse fit. Usage:

   A = EllipseDirectFit(XY)

Input: XY(n,2) is the array of coordinates of n points
            x(i)=XY(i,1), y(i)=XY(i,2)

 Output: A = [a b c d e f]' is the vector of coefficients
 of the equation of the best fitting ellipse:

   ax^2 + bxy + cy^2 + dx + ey + f = 0,

To convert this vector to the geometric parameters (semi-axes, center, etc.) use standard formulas, e.g., (19) - (24) in Wolfram Mathworld:
       http://mathworld.wolfram.com/Ellipse.html

This ellipse fit was proposed in article
A. W. Fitzgibbon, M. Pilu, R. B. Fisher
"Direct Least Squares Fitting of Ellipses"
IEEE Trans. PAMI, Vol. 21, pages 476-480 (1999)

The authors called it "Direct Ellipse Fit".

My code is based on a numerically stable version
of this fit published by R. Halir and J. Flusser. I only
added data centering to further improve performance.

Note: this fit returns ellipses only! You will get an ellipse
even if points can be better approximated by a hyperbola.

This fit is somewhat biased toward smaller ellipses.

Acknowledgements

The author wishes to acknowledge the following in the creation of this submission:
Ellipse Fit

• MATLAB release: MATLAB 6.0 (R12)

[OH]

??? Input argument "XY" is undefined.

Error in ==> EllipseDirectFit at 17
centroid = mean(XY); % the centroid of the data set
 
??? Input argument "XY" is undefined.

Error in ==> EllipseDirectFit at 17
centroid = mean(XY); % the centroid of the data set

[Nikolai Chernov]

from Author: XY is the Nx2 input array of coordinates of N points,
    x(i)=XY(i,1), y(i)=XY(i,2). To use the program, type
      A = EllipseDirectFit(XY)

[Francisco]

The parameters in your equation (i.e.: a, b, c, d, e and f) do not fit with equation 15 in Wolfram's Mathworld.
This needs further explanation.

[Laurens]

The Author (Fitzgibbon) posted some new code of his implementation at his site:

http://research.microsoft.com/en-us/um/people/awf/ellipse/

might be of interest

[Faten]

I tried to use this code for image and it did not work
Normal jpge image , binary image and gray image all faild

[Ed Shen]

The code works, however you must divde A(2), A(4) and A(5) by 2 to be able to use it with the Mathworld equations. I've verified with points given in http://research.microsoft.com/en-us/um/people/awf/ellipse/ and the results are almost identical

[GAURAV]

I have three points say (a1,b1),(a2,b2),(a3,b3)..
How do I paas this as input?? Please let me know..

Appreciate your help.

Thank You,
Gaurav.

[Shriniket Sarkar]

I tried using an XY array of
575 562
589 201
495 614
799 443
729 540
Using this set it provides exponential values for a,b and c return parameters and we cant use these parameters for a plot.

Is there something i am missing in this ? Can someone help?

[ICESI]

Excellent implementation. Some people want to see the plot of the fit, you just have to add the following code in the end of the script.

hold on

%Convert the A to str
a = num2str(A(1));
b = num2str(A(2));
c = num2str(A(3));
d = num2str(A(4));
e = num2str(A(5));
f = num2str(A(6));

%Equation
eqt= ['(',a, ')*x^2 + (',b,')*x*y + (',c,')*y^2 + (',d,')*x+ (',e,')*y + (',f,')'];
xmin=0.7*min(XY(:,1));
xmax=1.3*max(XY(:,2));
ezplot(eqt,[xmin,xmax])
scatter(XY(:,1),XY(:,2))
hold off