Fitting an ellipsoid to 3D data points

I have a large set of 3D data points to which I want to fit to an ellipsoid.

My maths is pretty poor, so I'm having trouble implementing the least squares method without any math libraries.

Does anyone know of or have a piece of code that can fit an ellipsoid to data which I can plug straight into my project? In C would be best, but it should be no problem for me to convert from C++, Java, C#, python etc.

EDIT: Just being able to find the cent开发者_JAVA技巧re would be a huge help too. Note that the points aren't evenly spaced so taking the mean won't result in the centre.

here you go:

This paper describes fitting an ellipsoid to multiple dimensions AS WELL AS finding the center for the ellipois. Hope this helps,

http://www.physics.smu.edu/~scalise/SMUpreprints/SMU-HEP-10-14.pdf

(btw, I'm assuming this answer is a bit late, but I figured I would add this solution for anyone who stumbles across your question in search for the same thing :)

If you want the minimum-volume enclosing ellipsoid, check out this SO answer for a bounding ellipsoid.

If you want the best fitting ellipse in a least-squares sense, check out this MATLAB code for error ellipsoids where you find the covariance matrix of your mean-shifted 3D points and use that to construct the ellipsoid.

Least Squares data fitting is probably a good methodology give the nature of the data you describe. The GNU Scientific Library contains linear and non-linear least squares data fitting routines. In your case, you may be able to transform your data into a linear space and use linear least-squares, but that would depend on your actual use case. Otherwise, you'll need to use non-linear methods.

I could not find a good Java based algorithm for fitting an ellipsoid, so I ended up writing it myself. There were some good algorithms for an ellipse with 2D points, but not for an ellipsoid with 3D points. I experimented with a few different MATLAB scripts and eventually settled on Yury Petrov's Ellipsoid Fit. It fits an ellipsoid to the polynomial Ax^2 + By^2 + Cz^2 + 2Dxy + 2Exz + 2Fyz + 2Gx + 2Hy + 2Iz = 1. It doesn't use any constraints to force an ellipsoid, so you have to have a fairly large number of points to prevent a random quardic from being fit instead of the ellipsoid. Other than that, it works really well. I wrote a small Java library using Apache Commons Math that implements Yury Petrov's script in Java. The GIT repository can be found at https://github.com/BokiSoft/EllipsoidFit.

We developed a set of Matlab and Java codes to fit ellipsoids here: https://github.com/pierre-weiss

You can also check our open-source Icy plugin. The following tutorial can be helpful: https://www.youtube.com/endscreen?video_referrer=watch&v=nXnPOG_YCxw

Note: most of the existing codes fit a generic quadric and do not impose an ellipsoidal shape. To get more robustness, you need to go to convex programming rather than just linear algebra. This is what is done in the indicated sources.

Cheers, Pierre

Here is unstrict solution with fast and simple random search approach*. Best side - no heavy linear algebra library required**. Seems it worked fine for mesh collision detection.

Is assumes that ellipsoid center matches cloud center and then uses some sort of mirrored average to search for main axis.

Full working code is slightly bigger and placed on git, idea of main axis search is here:

np.random.shuffle(pts)

pts_len = len(pts)
pt_average =  np.sum(pts, axis = 0) / pts_len

vec_major = pt_average * 0
minor_max, major_max = 0, 0

# may be improved with overlapped pass, 
for pt_cur in pts:
    vec_cur = pt_cur - pt_average
    proj_len, rej_len = proj_length(vec_cur, vec_major)

    if proj_len < 0:
        vec_cur = -vec_cur
    vec_major += (vec_cur - vec_major) / pts_len

    major_max = max(major_max, abs(proj_len))
    minor_max = max(minor_max, rej_len)

It can be improved/optimized even more at some points. Examples what it will produce:

And full experiment code with plots

*_{i.e. adjusting code lines randomly until they work}
**_{was actually reason to figure out this solution}

I have an idea. Approximately solution, not the best but will keep points inside. In XY plane find the radius R1 that will obtain all points. Same do for the XZ plane (R2) and YZ plane (R3). Then use the maximums on each axes. A=max(R1,R2), B=max(R1,R3) and C=max(R2,R3). But, first of all find the average (center) of all points and align it to origin.

I have just gone through the same process. Here is a python module which is based on work by Nima Moshtagh. Referenced in many places but also in this question about a Bounding ellipse

This module also handles plotting of the final ellipsoid. Enjoy!

https://github.com/minillinim/ellipsoid/blob/master/ellipsoid.py

I ported Yury Petrov's least-squares Matlab fitter to Java some time ago, it only needs JAMA: https://github.com/mdoube/BoneJ/blob/master/src/org/doube/geometry/FitEllipsoid.java