Given a set of points (x, y) this function returns the best fit ellipse.
6 Downloads
Updated 13 Dec 2008
No LicenseUsage:
[semimajor_axis, semiminor_axis, x0, y0, phi] = ellipse_fit(x, y)
Input:
x - a vector of x measurements
y - a vector of y measurements
Output:
semimajor_axis - Magnitude of ellipse longer axis
semiminor_axis - Magnitude of ellipse shorter axis
x0 - x coordinate of ellipse center
y0- y coordinate of ellipse center
phi - Angle of rotation in radians with respect to
the x-axis
Algorithm used:
Given the quadratic form of an ellipse:
a*x^2 + 2*b*x*y + c*y^2 + 2*d*x + 2*f*y + g = 0 (1)
we need to find the best (in the Least Square sense) parameters a,b,c,d,f,g.
To transform this into the usual way in which such estimation problems are presented,
divide both sides of equation (1) by a and then move x^2 to the other side. This gives us:
2*b'*x*y + c'*y^2 + 2*d'*x + 2*f'*y + g' = -x^2 (2)
where the primed parametes are the original ones divided by a. Now the usual estimation technique is used where the problem is presented as:
M * p = b, where M = [2*x*y y^2 2*x 2*y ones(size(x))],
p = [b c d e f g], and b = -x^2. We seek the vector p, given by:
p = pseudoinverse(M) * b.
From here on I used formulas (19) - (24) in Wolfram Mathworld:
http://mathworld.wolfram.com/Ellipse.html
1.1 | Added input and output explanations to description part. |
Inspired by: Circle fit
Inspired: Ellipse Fit (Taubin method), Ellipse Fit (Direct method)
Sumith YD (view profile)
This fitting won’t work well always (N Chernov pointed it out already)..For below data I got erroneous results...
X = {26 26 27 28 29 30 31 32 33 34 35 36 37 38 38 38 39 39 39 40 40 41 42 42 44 44 45 46 47 48 49 50 51 52 56 58 59 60 61 62 63 64 65 65 66 67 68}
Y = {4 5 6 6 6 6 6 6 7 8 8 7 7 5 6 7 4 5 6 4 6 6 5 6 6 7 7 7 6 6 7 7 7 8 8 8 7 7 7 6 7 7 6 7 6 5 4}
Aritra Das (view profile)
Can anyone provide a clear idea about how to implement this practically?
For example if I want to find ellipses in an binary image say bw, how to run this code on the image to get the ellipses?
As I see there is no way to provide the matrix name as a input argument.
And if somebody could explain the input arguments "Input: x,y - a set of points in 2 column vectors. AT LEAST 5 points are needed !"
I mean this statement a bit elaborately it will be very helpful.
duraid deikran (view profile)
i have a central vacancy of point after applying certain function, how can i fit ellipse on the periphery of this vacancy using your function
Marcello (view profile)
Thank you! Great job!
How may I use this function to plot the resulting ellipse?
dingding (view profile)
thank you ! help me so much!
Sreedu (view profile)
Hello,
At least how many points are required to get a reliable result?
Raphael Cautain (view profile)
I tested the function with a true ellipse, successively oriented from 0 to 2pi. All parameters were found
correct, except one :
The orientation was found with a wrong sign each time the axis direction is nearer from Ox than from Oy (that is, 1 time over 2).
Actually, recommendations from Wolfram's pages (results §23) should be - at least - implemented like this :
phi = 0.5 * acot((a-c)/(2*b));
if (a > c)
phi = phi + pi/2;
end
Rather than like the existing line 77 :
% phi = 0.5 * acot((c-a)/(2*b));
Thus the last lines of ellipse_fit are not adapted and must be removed :
% if (a_prime < b_prime)
% phi = pi/2 - phi;
Laura Suad (view profile)
Hello I am using this ellipse fit program but I need to have the: semimajor_axis, semiminor_axis, x0, y0, phi with dispersion (the fitting errors). How can I do to obtain these parameters?
Nikolai Chernov (view profile)
This is an unreliable ellipse fit, as it depends on the choice of the coordinate system (rotating or translating the xy frame would give you a different ellipse). Better fits should be invariant under transformations of coordinates. Besides, this fit may return a hyperbola, instead of ellipse. I have posted two wel-known and well-designed ellipse fits, Taubin Fit and Direct Fit, see files #22683 and #22684.