Discover MakerZone

MATLAB and Simulink resources for Arduino, LEGO, and Raspberry Pi

Learn more

Discover what MATLAB® can do for your career.

Opportunities for recent engineering grads.

Apply Today

Thread Subject:
Coplanar Points

Subject: Coplanar Points

From: Ammar T. Al-Sayegh

Date: 9 Jun, 2001 11:33:47

Message: 1 of 11

Hi Everyone,

Is there a function that can determine whether a set of points
in 3D space is coplanar?


-ammar

Subject: Coplanar Points

From: roger.elliexyzzy@mindspring.com.invalid (Roger Stafford)

Date: 9 Jun, 2001 14:02:57

Message: 2 of 11

In article <vQnU6.169369$f85.25050450@typhoon.nyroc.rr.com>, "Ammar T.
Al-Sayegh" <ammar@twcny.rr.com.spam.not> wrote:

> Hi Everyone,
>
> Is there a function that can determine whether a set of points
> in 3D space is coplanar?
>
> -ammar

------

  Ammar, probably the most insightful method of solving your problem is to
find the mean square orthogonal distance from the set of points to the
best-fitting plane in the least squares sense. If that distance is zero,
then the points are necessarily coplanar, and otherwise not.

  Let x, y , and z be n x 1 column vectors of the three coordinates of the
point set. Subtract from each, their respective mean values to get V, and
form from it the positive definite matrix A,

V = [x-mean(x),y-mean(y),z-mean(z)];
A = (1/n)*V'*V;

Then from

[U,D] = eig(A);

select the smallest eigenvalue in the diagonal matrix D. This is the mean
square orthogonal distance of the points from the best fitting plane and
the corresponding eigenvector of U gives the coefficients in the equation
of that plane, along with the fact that it must contain the mean point
(mean(x),mean(y),mean(z)).

  Note: you can use matlab's svd instead of eig and the smallest
eigenvalue will always be located in the lower right corner of D.

Roger
--
(Remove "xyzzy" and ".invalid" to send me email.)
Roger Stafford

Subject: Coplanar Points

From: Brett Shoelson

Date: 9 Jun, 2001 22:12:13

Message: 3 of 11

Hi Ammar,
This might not be as elegant as Roger's solution, but it's a good approach
nonetheless. It uses the fact that the triple scalar product of 3 vectors
vanishes iff the vectors are coplanar.

function copl=iscoplanar(x,y,z)

%Let vectors v1 and v2 determine the test plane:
v1=[x(1),y(1),z(1)];
v2=[x(2),y(2),z(2)];
crossprod=cross(v1,v2);

%The triple product A.(BxC) vanishes iff vectors A,B, and C are coplanar
for i=3:length(x)
    v3=[x(i),y(i),z(i)];
    if dot(crossprod,v3)~=0
       copl = 0;
       return;
    end
end
copl = 1;

Regards,
Brett


"Ammar T. Al-Sayegh" <ammar@twcny.rr.com.spam.not> wrote in message
news:vQnU6.169369$f85.25050450@typhoon.nyroc.rr.com...
> Hi Everyone,
>
> Is there a function that can determine whether a set of points
> in 3D space is coplanar?
>
>
> -ammar
>
>

Subject: Coplanar Points

From: roger.elliexyzzy@mindspring.com.invalid (Roger Stafford)

Date: 9 Jun, 2001 16:15:57

Message: 4 of 11

In article <1bxU6.22538$e34.7521897@typhoon.southeast.rr.com>, "Brett
Shoelson" <b@cox.rr.com> wrote:

> Hi Ammar,
> This might not be as elegant as Roger's solution, but it's a good approach
> nonetheless. It uses the fact that the triple scalar product of 3 vectors
> vanishes iff the vectors are coplanar.
>
> function copl=iscoplanar(x,y,z)
>
> %Let vectors v1 and v2 determine the test plane:
> v1=[x(1),y(1),z(1)];
> v2=[x(2),y(2),z(2)];
> crossprod=cross(v1,v2);
>
> %The triple product A.(BxC) vanishes iff vectors A,B, and C are coplanar
> for i=3:length(x)
> v3=[x(i),y(i),z(i)];
> if dot(crossprod,v3)~=0
> copl = 0;
> return;
> end
> end
> copl = 1;
>
> Regards,
> Brett

  Brett, with reference to Ammar's problem, your triple scalar product
approach here is valid in principle, but it is necessary to base your
vectors on one of the points in the set, not the origin. Otherwise you
would find iscoplanar declaring, for example, most sets of just three
points as being non-coplanar.

Roger
--
(Remove "xyzzy" and ".invalid" to send me email.)
Roger Stafford

Subject: Coplanar Points

From: Ammar T. Al-Sayegh

Date: 9 Jun, 2001 23:53:54

Message: 5 of 11

Roger and Brett,

Thank you for your elaborate follow-ups.

The method that I am using now is similar to Brett's
method, where I determine the vector orthogonal to the
lines connecting three of the points using the cross
product, then checking whether the dot product of this
orthogonal vector with the lines connecting the other
points is zero.

The method outlined by Roger looks elegant, but
unfortunately I cannot fully understand it. Perhaps
someone can rewrite it in the form of ISCOPLANAR
function, as Brett illustrated his method?


-ammar

Subject: Coplanar Points

From: Brett Shoelson

Date: 9 Jun, 2001 23:45:16

Message: 6 of 11

You're absolutely correct. New, improved code below. Thanks.

function copl=iscoplanar(x,y,z)
if length(x)~=length(y) | length(y)~=length(z)
    error('Invalid number of vertices.');
end
if length(x)<3
    copl=1;
    return;
end

% Let vectors v1 and v2 determine the test plane:
v1=[x(1)-x(2),y(1)-y(2),z(1)-z(2)];
v2=[x(3)-x(2),y(3)-y(2),z(3)-z(2)];
crossprod=cross(v1,v2);

% The triple product A.(BxC) vanishes iff vectors A,B, and C are coplanar
for i=3:length(x)
   v3=[x(i)-x(2),y(i)-y(2),z(i)-z(2)];
   if dot(crossprod,v3)~=0
      copl = 0;
      return
   end
end
copl = 1;


"Roger Stafford" <roger.elliexyzzy@mindspring.com.invalid> wrote in message
news:roger.elliexyzzy-0906011615570001@pool0108.cvx15-bradley.dialup.earthli
nk.net...
> In article <1bxU6.22538$e34.7521897@typhoon.southeast.rr.com>, "Brett
> Shoelson" <b@cox.rr.com> wrote:
>
> > Hi Ammar,
> > This might not be as elegant as Roger's solution, but it's a good
approach
> > nonetheless. It uses the fact that the triple scalar product of 3
vectors
> > vanishes iff the vectors are coplanar.
> >
> > function copl=iscoplanar(x,y,z)
> >
> > %Let vectors v1 and v2 determine the test plane:
> > v1=[x(1),y(1),z(1)];
> > v2=[x(2),y(2),z(2)];
> > crossprod=cross(v1,v2);
> >
> > %The triple product A.(BxC) vanishes iff vectors A,B, and C are coplanar
> > for i=3:length(x)
> > v3=[x(i),y(i),z(i)];
> > if dot(crossprod,v3)~=0
> > copl = 0;
> > return;
> > end
> > end
> > copl = 1;
> >
> > Regards,
> > Brett
>
> Brett, with reference to Ammar's problem, your triple scalar product
> approach here is valid in principle, but it is necessary to base your
> vectors on one of the points in the set, not the origin. Otherwise you
> would find iscoplanar declaring, for example, most sets of just three
> points as being non-coplanar.
>
> Roger
> --
> (Remove "xyzzy" and ".invalid" to send me email.)
> Roger Stafford

Subject: Coplanar Points

From: roger.elliexyzzy@mindspring.com.invalid (Roger Stafford)

Date: 9 Jun, 2001 19:25:59

Message: 7 of 11

In article <mGyU6.2844$3y3.512815@typhoon.nyroc.rr.com>, "Ammar T.
Al-Sayegh" <ammar@twcny.rr.com.spam.not> wrote:

> The method outlined by Roger looks elegant, but
> unfortunately I cannot fully understand it. Perhaps
> someone can rewrite it in the form of ISCOPLANAR
> function, as Brett illustrated his method?
> -ammar

------

  Ammar, here is a function for the minimum orthogonal distance which I
described earlier as per your request. I will let you convert it into a
boolean function, but you should probably provide some tolerance for round
off errors rather than requiring exact equalities, particularly in view of
the use of matlab's 'svd' function here. It is in setting such a
tolerance level that this algorithm may have a value, since returned
quantity s is in the same units (squared) as the coordinates of the given
point set.

  You will note that the function also returns with the best fitting plane
defined, if you happened to be interested in that.

------
function [s,p] = best_fit_plane(x,y,z);

% x,y,z are n x 1 column vectors of the three coordinates
% of a set of n points in three dimensions. s returns with
% the minimum mean square orthogonal distance to a least
% squares best-fit plane. The four coefficients of that
% plane's equation, Ax + By + Cz + D = 0, are returned in
% row vector p = [A,B,C,D]. A,B,C are normalized: A^2+B^2+C^2=1.

[n,mx] = size(x); [ny,my] = size(y); [nz,mz] = size(z);
if (mx~=my)|(mx~=mz)|(ny~=n)|(nz~=n)
 error('The arguments must be column vectors of the same length.')
end
m = [mean(x),mean(y),mean(z)];
w = [x-m(1),y-m(2),z-m(3)]; % Use "mean" point as base
a = (1/n)*w'*w; % 'a' is a positive definite matrix
[u,d,v] = svd(a); % 'eig' & 'svd' get same eigenvalues for this matrix
s = d(3,3); % Select the smallest eigenvalue
p = u(:,3)'; % Get the corresponding (normalized) eigenvector (A,B,C)
p(4) = -p*m'; % Choose D so plane contains "mean" point
------

Roger
--
(Remove "xyzzy" and ".invalid" to send me email.)
Roger Stafford

Subject: Coplanar Points

From: Ammar T. Al-Sayegh

Date: 10 Jun, 2001 03:31:45

Message: 8 of 11

Roger and Brett, many thanks for the functions.
I'll try them out...


-ammar

Subject: Coplanar Points

From: roger.elliexyzzy@mindspring.com.invalid (Roger Stafford)

Date: 9 Jun, 2001 23:32:22

Message: 9 of 11

In article <gyyU6.23358$e34.7623961@typhoon.southeast.rr.com>, "Brett
Shoelson" <b@cox.rr.com> wrote:

> function copl=iscoplanar(x,y,z)
> if length(x)~=length(y) | length(y)~=length(z)
> error('Invalid number of vertices.');
> end
> if length(x)<3
> copl=1;
> return;
> end
>
> % Let vectors v1 and v2 determine the test plane:
> v1=[x(1)-x(2),y(1)-y(2),z(1)-z(2)];
> v2=[x(3)-x(2),y(3)-y(2),z(3)-z(2)];
> crossprod=cross(v1,v2);
>
> % The triple product A.(BxC) vanishes iff vectors A,B, and C are coplanar
> for i=3:length(x)
> v3=[x(i)-x(2),y(i)-y(2),z(i)-z(2)];
> if dot(crossprod,v3)~=0
> copl = 0;
> return
> end
> end
> copl = 1;

------

  Brett, on considering the matter more carefully, I believe there is
another aspect of your method that could cause trouble in certain cases.
If the user were unlucky and constructed the initial two vectors v1 and v2
from three points which happened to be colinear, the resulting cross
product vector would be of zero length and no further points could
possibly fail the test of being coplanar with the first three. If the
first three points were nearly colinear and if you put some tolerance
allowance in your criterion to allow for round-off errors, then too wide a
variety of points might pass. I think your initial three points should be
selected with some considerable degree of care to make this work properly.

Roger
--
(Remove "xyzzy" and ".invalid" to send me email.)
Roger Stafford

Subject: Coplanar Points

From: Brett Shoelson

Date: 10 Jun, 2001 15:19:40

Message: 10 of 11

Roger, once again, you are correct. That should teach me to try to do
something quickly without thinking it through. So, though my solution gets
less and less elegant, I'm reposting a modified version. (I think I have
sufficiently accounted for coincident and collinear points in this version.)
Thanks for your vigilance. Ammar, sorry this took a few iterations, but I
think I've covered the bases this time around.
Brett

function copl=iscoplanar(x,y,z)
% Tests input points for coplanarity.
% COPL = ISCOPLANAR(X,Y,Z) takes input arguments x,y, and z as column
vectors.
% COPL = ISCOPLANAR(x) takes an n x 3 input argument
% in the form [x1 y1 z1;x2 y2 z2;...;xn yn zn]
if nargin == 0
   error('Requires at least one input argument.');
elseif nargin == 1
   if size(x,2) == 3
      % Matrix of all x,y,z is input
      allpoints = x;
      %Convert input to column vectors
      y = x(:,2);
      z = x(:,3);
      x = x(:,1);
   else
      error('Invalid input.')
   end
else
   % Compile a matrix of all x,y,z
   allpoints = [x y z];
end

if length(x)<=3
   disp('Three or fewer points are necessarily coplanar.');
   copl=1;
   return;
end

% Select first two vectors that span the test plane
% Any two non-coincident points define v1:
[B,I,J] = unique(allpoints,'rows');
if length(I) == 1
   disp('All points are coincident (and hence coplanar).');
   copl = 1;
   return;
else
   v1 = [x(I(2)) - x(I(1)), y(I(2)) - y(I(1)), z(I(2)) - z(I(1))];
end

% If possible, define v2 such that it is not parallel to v1
for i = 1:length(x)
   v2 = [x(i) - x(I(1)), y(i) - y(I(1)),z(i) - z(I(1))];
   crossprod=cross(v1,v2);
   if norm(crossprod) ~= 0
      % v1 and v2 are non-parallel, and thus are sufficient to define a
plane
      break;
   end
end
if crossprod == 0
   disp('All points are collinear (and hence coplanar).');
    copl = 1;
    return;
end

% At this point, v1 and v2 necessarily span a test plane
crossprod=cross(v1,v2);

% %The triple product A.(BxC) vanishes iff vectors A,B, and C are coplanar
for i = 1:length(x)
   v3 = [x(i) - x(I(2)), y(i) - y(I(2)), z(i) - z(I(2))];
   if dot(crossprod,v3) ~= 0
       copl = 0;
       return
   end
end
copl = 1;



> Brett, on considering the matter more carefully, I believe there is
> another aspect of your method that could cause trouble in certain cases.
> If the user were unlucky and constructed the initial two vectors v1 and v2
> from three points which happened to be colinear, the resulting cross
> product vector would be of zero length and no further points could
> possibly fail the test of being coplanar with the first three. If the
> first three points were nearly colinear and if you put some tolerance
> allowance in your criterion to allow for round-off errors, then too wide a
> variety of points might pass. I think your initial three points should be
> selected with some considerable degree of care to make this work properly.
>
> Roger
> --
> (Remove "xyzzy" and ".invalid" to send me email.)
> Roger Stafford

Subject: Coplanar Points

From: Ammar T. Al-Sayegh

Date: 10 Jun, 2001 18:53:08

Message: 11 of 11

"Brett Shoelson" <b@cox.rr.com> wrote in message
news:geMU6.25492$e34.8542512@typhoon.southeast.rr.com...
> Roger, once again, you are correct. That should teach me to try to do
> something quickly without thinking it through. So, though my solution gets
> less and less elegant, I'm reposting a modified version. (I think I have
> sufficiently accounted for coincident and collinear points in this version.)
> Thanks for your vigilance. Ammar, sorry this took a few iterations, but I
> think I've covered the bases this time around.

Don't worry, I'm sure that it would have taken me a lot
more iterations to arrive to the same result. That's why
I came here :)


-ammar

Tags for this Thread

No tags are associated with this thread.

What are tags?

A tag is like a keyword or category label associated with each thread. Tags make it easier for you to find threads of interest.

Anyone can tag a thread. Tags are public and visible to everyone.

Contact us