From: "Robert Phillips" <>
Newsgroups: comp.soft-sys.matlab
Subject: Re: Shortest distance from point to ellipsoid surface
Date: Sat, 24 Sep 2011 20:02:14 +0000 (UTC)
Organization: Embry-Riddle Aeronautical University
Lines: 113
Message-ID: <j5ld06$blu$>
References: <ivd95k$6h2$> <ivep55$3b$> <ivfb8h$rmv$> <ivfd2l$43v$> <j5e0fb$7es$> <j5fed5$f6h$>
Reply-To: "Robert Phillips" <>
Content-Type: text/plain; charset=UTF-8; format=flowed
Content-Transfer-Encoding: 8bit
X-Trace: 1316894534 11966 (24 Sep 2011 20:02:14 GMT)
NNTP-Posting-Date: Sat, 24 Sep 2011 20:02:14 +0000 (UTC)
X-Newsreader: MATLAB Central Newsreader 2598083
Xref: comp.soft-sys.matlab:744126

Thank you for the help; I ended up defining the coefficients as calculated numeric values and then using a for-loop to take ROOTS for each point (I'm calculating minimum distances from points in a Nx3 matrix to the ellipsoid surface). 

I've checked my coefficients and roots, and they seem to be correct. However, I am in fact having trouble with the transformation procedure...

N is the number of points that are used for calculating minimum distances to the ellipsoid surface. N = size(astrd_int_2R,1).

a_sqr, b_sqr, & c_sqr are a^2, b^2, & c^2, respectively, and a>b>c.

rm is the "ellipsoid rotation matrix", where its 3 rows are normalized and orthogonal semiaxis vectors, in order a; b; c. So essentially it's [ax,ay,az; bx,by,bz; cx,cy,cz].

ctrd is the "ellipsoid centroid", [Ex,Ey,Ez].

I have the coefficients C1 through C7 (which match exactly to Roger's I believe), all of which are Nx1 size, which are used in my loop. They are ordered such that C1 corresponds to t^6, C2 to t^5, and so on until C7=1 for t^0.

    % For all interior points inside 2*max(R(E-1)):
    for j=1:N
        % Calculate t_roots.
        t_roots = roots([C1(j,1),C2(j,1),C3(j,1),C4(j,1),...
        % Take real t_roots.
        t_roots = t_roots(imag(t_roots)==0);
        % Calculate ellipsoid surface xyz coordinates corresponding to
        % t_roots.
        x_e_rel = (a_sqr*astrd_int_2R(j,1))./(a_sqr-t_roots);
        y_e_rel = (b_sqr*astrd_int_2R(j,2))./(b_sqr-t_roots);
        z_e_rel = (c_sqr*astrd_int_2R(j,3))./(c_sqr-t_roots);
        % Assign e_rel to xyz coordinates corresponding to point with
        % minimum distance.
        [~,idx_mn] = min(sum([x_e_rel - astrd_int_2R(j,1),...
                                          y_e_rel - astrd_int_2R(j,2),...
                                          z_e_rel - astrd_int_2R(j,3)].^2,2));
        e_rel(j,:) = [x_e_rel(idx_mn),...
    end % End all interior points inside 2*max(R(E-1)).

And I perform a transformation, *which I'm not sure is appropriate / correct*:
e_new = [sum(bsxfun(@times,e_rel,(rm(:,1)')),2) + ctrd(1,1)...
               sum(bsxfun(@times,e_rel,(rm(:,2)')),2) + ctrd(1,2)...
               sum(bsxfun(@times,e_rel,(rm(:,3)')),2) + ctrd(1,3)];

which is (or at least should be) equivalent to:

x_e_new = e_rel(:,1)*rm(1,1) +...
                 e_rel(:,2)*rm(2,1) +...
                 e_rel(:,3)*rm(3,1) + ctrd(1,1);
y_e_new = e_rel(:,1)*rm(1,2) +...
                 e_rel(:,2)*rm(2,2) +...
                 e_rel(:,3)*rm(3,2) + ctrd(1,2);
z_e_new = e_rel(:,1)*rm(1,3) +...
                 e_rel(:,2)*rm(2,3) +...
                 e_rel(:,3)*rm(3,3) + ctrd(1,3);
e_new = [x_e_new, y_e_new, z_e_new];

When I use scatter3 to visualize the e_rel points and e_new points (along with other objects), they look (very) wrong. The e_new points do not seem to lie on the surface of the ellipsoid (whose surface points are properly visualized in the same figure). Also the e_rel and e_new points are very small in magnitude, too.

I figured that, since the ellipsoid surface points are (at least interpreted to be) oriented along the coordinate system's axes and centered at the coordinate system's origin, then it is appropriate to transform the coordinates e_new according to the ellipsoid rotation matrix and centroid.

Is there something that I'm doing wrongly here?

For reference, the real t_roots for the first for-loop iteration are
and the t_roots selected for minimum distance was -444.64

"Steven_Lord" <> wrote in message <j5fed5$f6h$>...
> "Robert Phillips" <> wrote in message 
> news:j5e0fb$7es$
> > I've finally had the time to work more on this. I've gone through 
> > step-by-step and obtained the polynomial equation and its seven 
> > coefficients (although it seems my coefficients are the negatives of the 
> > ones provided by Roger Stafford). I wasn't familiar with how Lagrange 
> > multipliers work so it took a bit to absorb it.
> >
> > Before defining the polynomial's seven coefficients and taking roots etc., 
> > I'm simply "interpreting" each ellipsoid as being translated and rotated. 
> > In my head it wouldn't make sense to literally transform it because 
> > there's technically "nothing" to transform. Am I correct in this thought 
> > process?
> >
> >
> > Also, I'm having difficulty with the roots command (I haven't worked much 
> > with the symbolic toolbox). I receive an error:
> > "??? Undefined function or method 'isfinite' for input
> > arguments of type 'sym'.
> >
> > Error in ==> roots at 27
> > if ~all(isfinite(c))"
> > Naturally, I have c1,c2,c3,c4,c5,c6,c7 defined as 'sym'. Am I missing a 
> > step here? Am I supposed to define the coefficients in some other way?
> ROOTS is intended to be used with _numeric_ vectors of polynomial 
> coefficients; it will not work with _symbolic_ vectors of polynomial 
> coefficients. Since you have a symbolic vector of polynomial coefficients, 
> use POLY2SYM and pass the resulting polynomial into SOLVE.
> *snip*
> -- 
> Steve Lord
> To contact Technical Support use the Contact Us link on