Prev: calling value of variable from main program into myfun while using fsolve
Next: Matlab GUI: jump from one callback to another
From: Els on 21 May 2010 08:46 You are correct. The one with the shortest distance to the surface is the one I am interested in. But if that is too much trouble, I can calculate that myself. If I just had the two points, that would make my day. Thanks again.
From: Torsten Hennig on 21 May 2010 04:55 > Where the normal vector is perpendicular to the Line > from P1 to P2 and runs through C and C'. Ok, so you try to find the point C' on the surface of the ellipsoid nearest to C. If c = (cx,cy,cz) and C'=(csx,csy,csz), the problem is min: (cx-csx)^2 + (cy-csy)^2 + (cz-csz)^2 s.t. (csx-xc)^2/rx^2 + (csy-yc)^2/ry^2 + (csz-zc)^2/rz^2 = 1. Introducing a Lagrange parameter lambda, the solution is given by a stationary point of the Lagrange function L(lambda,csx,csy,csz) = (cx-csx)^2 + (cy-csy)^2 + (cz-csz)^2 - lambda* ( (csx-xc)^2/rx^2 + (csy-yc)^2/ry^2 + (csz-zc)^2/rz^2 - 1) So calculate first-order derivatives of L with respect to csx,csy,csz and lambda, set them to zero and solve the resulting system of 4 equations with fsolve, e.g. Best wishes Torsten.
From: Els on 21 May 2010 09:59 Dear Torsten, This sounds great, the solution is near. But I have never calculated a first-order derivative in Matlab before. And the xc,yc and zc, are they the same as cx,cy and cz? Thanks again.
From: Torsten Hennig on 21 May 2010 07:34 > Dear Torsten, > > This sounds great, the solution is near. But I have > never calculated a first-order derivative in Matlab > before. And the xc,yc and zc, are they the same as > cx,cy and cz? > No, (xc,yc,zc) is the centre of the ellipsoid. > Thanks again. L(lambda,csx,csy,csz) = (cx-csx)^2 + (cy-csy)^2 + (cz-csz)^2 - lambda* ( (csx-xc)^2/rx^2 + (csy-yc)^2/ry^2 + (csz-zc)^2/rz^2 - 1) So calculate first-order derivatives of L with respect to csx,csy,csz and lambda, set them to zero and solve the resulting system of 4 equations with fsolve, e.g. The first-order derivatives are given by dL/d(lambda) = -((csx-xc)^2/rx^2 + (csy-yc)^2/ry^2 + (csz-zc)^2/rz^2 - 1) dL/d(csx) = -2*(cx-csx)-2*lambda*(csx-xc)/rx^2 dL/d(csy) = -2*(cy-csy)-2*lambda*(csy-yc)/ry^2 dL/d(csz) = -2*(cz-csz)-2*lambda*(csz-zc)/rz^2 Now solve dL/d(lambda) = 0 dL/d(csx) = 0 dL/d(csy) = 0 dL/d(csz) = 0 for lambda, csx,csy and csz using fsolve. Best wishes Torsten.
From: Bruno Luong on 21 May 2010 13:28
This problem of projection on ellipse in 2D can be solved without optimization toolbox, here is the demo: % Data, generate ellipse centered at origin, x'*Q*x = c L=randn(2); Q = L'*L; % (2 x 2) sdp matrix c = rand; % Given point (to be projected on ellipse), it can be inside or outside P = randn(2,1); % Engine [U S V] = svd(Q); s = diag(S); A = U*diag(sqrt(c./s)); R = [0 1; -1 0]; B = R*Q*A; alpha = -P.'*R*Q*A; rho = sqrt(s(1)/s(2)); beta = c*(rho-1/rho); Poly = beta*[0 -2 0 2 0] + ... alpha(1)*[-1 0 0 0 1] + ... alpha(2)*[0 2 0 2 0]; tprj = roots(Poly); tprj = tprj(imag(tprj)==0); % remove imaginary solution tprj = reshape(tprj, 1, []); % Candidate of the projections xprj = A*[(1-tprj.^2); 2*tprj]; xprj = bsxfun(@rdivide, xprj, (1+tprj.^2)); % Graphical check clf % Generate point on ellipse theta = linspace(0,2*pi); t = tan(theta/2); x = A*[(1-t.^2); 2*t]; x = bsxfun(@rdivide, x, (1+t.^2)); plot(x(1,:),x(2,:)); axis equal; hold on; % Plot P projected to the ellipse plot(P(1),P(2),'or'); for k=1:size(xprj,2) plot([xprj(1,k) P(1)],[xprj(2,k) P(2)],'-r'); end % Bruno |