Hello! How to define on which plane lies a given point with coordinates? I have coordinates of plane and normal

13 visualizzazioni (ultimi 30 giorni)
Hello everyone. I undertsand it might be very simple question, but I could not solve it. I have point for example A(1.5, 1.5, 3.0) and i have 6 planes of cube and I need to determine on which plane lies this point. I found normal for each plane and need to find for which of these planes ( number of plane) rely A.
Thank you so much in advance for your time and consideration . More appreciate help and any advice.
X0 = 1.5;
Y0 = 1.5;
Z0 = 3.0;
planes(:,:,1) = [0 3 3; 0 0 3; 0 3 0; 0 0 0; 0 0 0];
planes(:,:,2) = [0 0 3; 3 0 3; 0 0 0; 3 0 0; 0 0 0];
planes(:,:,3) = [3 0 3; 3 3 3; 3 0 0; 3 3 0; 3 0 0];
planes(:,:,4) = [3 3 3; 0 3 3; 3 3 0; 0 3 0; 0 3 3];
planes(:,:,5) = [0 3 0; 3 3 0; 0 0 0; 3 0 0; 0 0 0];
planes(:,:,6) = [0 3 3; 3 3 3; 0 0 3; 3 0 3; 0 0 3];
location_plane = 6;
for j=1:6 % j is number of plane
j
plane = planes(:,:,j);
p0 = plane(1,:); %p0 is top left point of plane
p1 = plane(2,:); %p1 is top right point of plane
p2 = plane(3,:); %p2 is bottom left point of plane
p3 = plane(4,:); %p3 is bottom right point of plane
V0 = plane(5,:); %point on the plane
% Pi is initial start point on the ray
Pi = [X0 Y0 Z0]; %initial start point
%Ri = [XBar YBar ZBar]; %direction vector with unit length
A = p0-p2; %calculate A and B then
B = p0-p3; %then to calculate Normal of each plane
n=cross(A,B); % Normal for each Plane
n

Risposta accettata

John D'Errico
John D'Errico il 6 Apr 2023
Modificato: John D'Errico il 6 Apr 2023
Does a point lie in a plane? This is quite easy to do.
The equation of a plane is a simple one. In fact, it works in any number of dimensions. I'll write it out using the dot function, to make it clear that it uses a dot product. If P0 is any point in the plane, and N is the normal vector to the plane, then the equation that tells if X lies in the plane is simply:
dot(X - P0,N) == 0
Of course, you NEVER want to test if something is exactly zero in floating point arithmetic, so you might do this instead:
abs(dot(X - P0,N)) < tol
where tol is some sufficiently small number.
First, you need to start working with VECTORS. I see you did that later on, but get used it it.
XYZ = [1.5, 1.5, 3]; % a point to test
planes(:,:,1) = [0 3 3; 0 0 3; 0 3 0; 0 0 0; 0 0 0];
planes(:,:,2) = [0 0 3; 3 0 3; 0 0 0; 3 0 0; 0 0 0];
planes(:,:,3) = [3 0 3; 3 3 3; 3 0 0; 3 3 0; 3 0 0];
planes(:,:,4) = [3 3 3; 0 3 3; 3 3 0; 0 3 0; 0 3 3];
planes(:,:,5) = [0 3 0; 3 3 0; 0 0 0; 3 0 0; 0 0 0];
planes(:,:,6) = [0 3 3; 3 3 3; 0 0 3; 3 0 3; 0 0 3];
% a point in each plane. I'll just use the first point in each plane
% I could also have used squeeze here
PiP = reshape(planes(1,:,:),3,6)';
% normal vectors to each plane
Normals = zeros(6,3);
for i = 1:6;
% subtract off that first point, then apply null
% the result will be a column vector, so transpose it at the end.
Normals(i,:) = null(planes(:,:,i) - PiP(i,:));
end
PiP
PiP = 6×3
0 3 3 0 0 3 3 0 3 3 3 3 0 3 0 0 3 3
Normals
Normals = 6×3
1 0 0 0 -1 0 1 0 0 0 1 0 0 0 1 0 0 1
Now, is the point XYZ in any of those planes? We could get fancy for ALL of the planes, and do the test in one line of code. For example:
dot(XYZ - PiP,Normals,2)
ans = 6×1
1.5000 -1.5000 -1.5000 -1.5000 3.0000 0
And that suggests it is plane number 6 that we want to find. But a tolerance makes sense still. So you should do this:
tol = eps*100; % a reasonably small number
find(abs(dot(XYZ - PiP,Normals,2)) < tol)
ans = 6
We could have done the same thing using a loop of course.
tol = eps*100;
for i = 1:6
if abs(dot(XYZ - PiP(i,:),Normals(i,:))) < tol
loc = i;
disp("plane # "+loc+" was identified")
break
end
end
plane # 6 was identified
As you can see, nothing more than the equation of a plane was ever needed to determine if the point lies in one of those planes. Be careful to make sure the tolerance is not too tight, but 100*eps is fine for this. I suppose, if you wanted to be a little better about it, I might have used a tolerance that would depend on the numbers themselves. So perhaps this:
smartertol = std(PiP(:))*10*eps
smartertol = 3.3415e-15
Note that if your numbers scale in various ways, so the cube gets larger or smaller, then smartertol also scales with them nicely. As well, we could even have been tricky and computed the normal vectors without a loop, using the pagesvd function, but that starts to get a little (unecessarily) deep for this answer.
  3 Commenti
John D'Errico
John D'Errico il 13 Apr 2023
It (loc) tells which plane was identified.
Did you not see the next line?
disp("plane # "+loc+" was identified")
what does that do?

Accedi per commentare.

Più risposte (1)

Davide Masiello
Davide Masiello il 6 Apr 2023
Modificato: Davide Masiello il 6 Apr 2023
Four points can be verified to belong on the same plane if the determinant of the matrix
is equal to zero.
I'd translate that into code in the following way.
x0 = 1.5;
y0 = 1.5;
z0 = 3;
planes(:,:,1) = [0 3 3; 0 0 3; 0 3 0; 0 0 0; 0 0 0];
planes(:,:,2) = [0 0 3; 3 0 3; 0 0 0; 3 0 0; 0 0 0];
planes(:,:,3) = [3 0 3; 3 3 3; 3 0 0; 3 3 0; 3 0 0];
planes(:,:,4) = [3 3 3; 0 3 3; 3 3 0; 0 3 0; 0 3 3];
planes(:,:,5) = [0 3 0; 3 3 0; 0 0 0; 3 0 0; 0 0 0];
planes(:,:,6) = [0 3 3; 3 3 3; 0 0 3; 3 0 3; 0 0 3];
for j = 1:6
x(j,:) = planes(1:3,1,j);
y(j,:) = planes(1:3,2,j);
z(j,:) = planes(1:3,3,j);
end
for j = 1:6
M = [x(j,:),x0;y(j,:),y0;z(j,:),z0;ones(1,4)];
if det(M) == 0
fprintf('The point belongs to face %d.\n',j)
break
end
end
The point belongs to face 6.
Examples with other points
x0 = 1.5;
y0 = 3;
z0 = 1;
for j = 1:6
M = [x(j,:),x0;y(j,:),y0;z(j,:),z0;ones(1,4)];
if det(M) == 0
fprintf('The point belongs to face %d.\n',j)
break
end
end
The point belongs to face 4.
x0 = 1.5;
y0 = 0;
z0 = 2;
for j = 1:6
M = [x(j,:),x0;y(j,:),y0;z(j,:),z0;ones(1,4)];
if det(M) == 0
fprintf('The point belongs to face %d.\n',j)
break
end
end
The point belongs to face 2.
Only thing, for edges 2 answers would be correct, but this code will output only the first of the 2 faces the point belongs to.
  5 Commenti

Accedi per commentare.

Categorie

Scopri di più su Desktop in Help Center e File Exchange

Prodotti


Release

R2021a

Community Treasure Hunt

Find the treasures in MATLAB Central and discover how the community can help you!

Start Hunting!

Translated by