Something

clear all;
load planet.dat;
x = planet(:,2);
y = planet(:,3);
% plot(x, y);
coeff = fminsearch('planet_fit', [5 4 3 0]);


function E = planet_fit(coeff)
load planet.dat;
x = planet(:,2)
y = planet(:,3)
a = coeff(1)
b = coeff(2)
x0 = x(1)
y0 = y(1)
E = sum(abs(((x - x0).^2 / a.^2) + ((y - y0).^2 / b^.2) - 1).^2);