File Exchange

image thumbnail

Determine the distance between two ellipses (in 3D)

version 1.2.0.0 (17.6 KB) by Rody Oldenhuis
Determine the minimum distance (and corresponding angles) between two ellipses

2 Downloads

Updated 10 Jun 2018

GitHub view license on GitHub

The problem of finding the geometric (minimum) distance between two arbitrary ellipses is surprisingly difficult. The general problem of finding all stationary points (minimum/maximum/saddle, no less than 12 possible points) has indeed been solved, but that algorithm is horrificly complex and would require thousands if not millions of operations once implemented.
The function distanceEllipseEllipse() is based on a somewhat more practical algorithm which limits itself to minima in the distance function. It is based on repeatedly finding the distance between a point and an ellipse, which can be done analytically (see my other post, distanceEllipsePoints.m). The algorithm by itself is not very robust (it frequently finds a local minimum, which is *not* the true distance).
This function executes the algorithm 4 times, for 4 different initial values, which greatly improves its robustness. Some numerical experimentation (comparing with a brute-force search) has shown that the true minimum distance is returned in more than 95% of the cases. The algorithm is implemented such that MATLAB's JIT-accelerator can accelerate it to the fullest extent, which makes it pretty fast and well-suited to handle large datasets requiring this calculation.
This is an implementation of the algorithm described in
Ik-Sung Kim: "An algorithm for finding the distance between two
% ellipses". Commun. Korean Math. Soc. 21 (2006), No.3, pp.559-567.

A copy-pastable example (also in header of M-file):

% Ellipse1 Ellipse2(=circle)
a = [2.0 1.0];
b = [0.5 1.0];
c = {[0,0,0], [-2,2,0]}; % location of centers
u = {[1,0,0], [1,0,0]}; % both oriented in XY-plane
v = {[0,1,0], [0,1,0]}; % to visualize them more easily

% plot the ellipses
f = 0:0.01:2*pi;
E1 = [a(1)*cos(f) + c{1}(1); b(1)*sin(f) + c{1}(2)];
E2 = [a(2)*cos(f) + c{2}(1); b(2)*sin(f) + c{2}(2)];
figure, hold on
plot(E1(1,:),E1(2,:),'r', E2(1,:),E2(2,:),'b')
axis equal

% run routine
[min_dist, fp_min, fs_min] = ...
distanceEllipseEllipse(a,b,c,u,v)

% plot the minimum distance returned
x = [a(1)*cos(fp_min) + c{1}(1), a(2)*cos(fs_min) + c{2}(1)];
y = [b(1)*sin(fp_min) + c{1}(2), b(2)*sin(fs_min) + c{2}(2)];
line(x,y,'color', 'k')

This should generate the screen shot given here.

If you find this work useful, please consider a small donation:
https://www.paypal.me/RodyO/3.5

Comments and Ratings (4)

@xiaoqing Yes, indeed; I'll update the example to allow for easier tinkering with the input data :)

xiaoqing

okay, so for the example in your .m file to correctly visualize the result, u,v need to be normalized before calculating the coordinates of those two points.

xiaoqing

a = [2.0 1.2];
b = [0.5 1.0];
c = {[0,0,0],[1,3,0]}; % location of centers

u = {[1,1,0], [1,0,0]}; % both oriented in XY-plane
v = {[-1,1,0], [0,1,0]}; % to visualize them more easily

does not seem to calculate a correct minimum distance.

Cihan Ulas

Hi,
I generated two plane with as follows;

par1.a=2;par1.b=2;par1.c=1;par1.d=2;
par2.a=10;par2.b=1;par2.c=2;par2.d=100;

P1=generate_plane(par1,-10:10,-10:10,'r.');
P2=generate_plane(par2,-50:-40,-50:-40,'b.');
---------------------------------------
function P=generate_plane(par,rangex,rangey,color)
[x y]=meshgrid(rangex,rangey);
z=(par.a*x(:)+par.b*y(:)+par.d)/-par.c;
z=3*rand+z;
P=[x(:) y(:) z];
----------------------------------------

Then find the principal axes and orientations as follows,
[U1 E1]=svd(cov(P1));
[U2 E2]=svd(cov(P2));

a=[E1(1,1) E2(1,1)]; %principal axis
b=[E1(2,2) E2(2,2)] ; %minor axis
c={mean(P1), mean(P2)};
u= {U1(:,1)', U2(:,1)'};
v= {U1(:,2)', U2(:,2)'};

However the algorithm can not find the correct distance. Also if two ellipses are intersect it fails. The result should be 0.

For 3D, an explicit example is necessary.

Updates

1.2.0.0

Description update

1.2.0.0

[linked to Github]

1.1.0.0

Updated contact info

MATLAB Release Compatibility
Created with R2009b
Compatible with any release
Platform Compatibility
Windows macOS Linux