admin管理员组

文章数量:1398975

I have two bodies in space composed by 3 points.

A = [-50.0098 -181.0869  -37.5881];
B = [  -49.8410 -180.8097  -41.1383];
C = [-49.5887 -180.3952  -46.4462];

A = [-50.199    -180.93 -37.249];
B = [-50.256    -180.823    -33.686];
C = [-50.341    -180.657    -28.359];

I want to plot a plane through one of the 3D points and calculate the distance between the other point and this plane.

Additionally, I want to plot it.

I used the following function:

%% ====== FUNCTION ======
function [normal, d, X, Y, Z] = plot_line_deviazione(p1, p2, p3, x_min, x_max, y_min, y_max)
% This function plots a line from three points. 
% I/P arguments: 
%   p1, p2, p3 eg, p1 = [x y z]
% 
%
% O/P is: 
% normal: it contains a,b,c coeff , normal = [a b c]
% d : coeff
normal = cross(p1 - p2, p1 - p3);
d = p1(1)*normal(1) + p1(2)*normal(2) + p1(3)*normal(3);
d = -d;
x = x_min:x_max; y = y_min:y_max;
[X,Y] = meshgrid(x,y);
Z = (-d - (normal(1)*X) - (normal(2)*Y))/normal(3);
mesh(X,Y,Z)

In the following way:

   figure
plot3(A(:,1),A(:,2),A(:,3),'r.','Markersize',25)
hold on
plot3(B(:,1),B(:,2),B(:,3),'r.','Markersize',25)
plot3(C(:,1),C(:,2),C(:,3),'r.','Markersize',25)
[normal0, d0, X0, Y0, Z0] = plot_line_deviazione(A, B, C, -100, 100, -100, 100);  % new plane

dot(X-d0,normal0)

The "dot" would be the distance that I was looking for. Is this correct? How can I plo it?

Thanks!

I have two bodies in space composed by 3 points.

A = [-50.0098 -181.0869  -37.5881];
B = [  -49.8410 -180.8097  -41.1383];
C = [-49.5887 -180.3952  -46.4462];

A = [-50.199    -180.93 -37.249];
B = [-50.256    -180.823    -33.686];
C = [-50.341    -180.657    -28.359];

I want to plot a plane through one of the 3D points and calculate the distance between the other point and this plane.

Additionally, I want to plot it.

I used the following function:

%% ====== FUNCTION ======
function [normal, d, X, Y, Z] = plot_line_deviazione(p1, p2, p3, x_min, x_max, y_min, y_max)
% This function plots a line from three points. 
% I/P arguments: 
%   p1, p2, p3 eg, p1 = [x y z]
% 
%
% O/P is: 
% normal: it contains a,b,c coeff , normal = [a b c]
% d : coeff
normal = cross(p1 - p2, p1 - p3);
d = p1(1)*normal(1) + p1(2)*normal(2) + p1(3)*normal(3);
d = -d;
x = x_min:x_max; y = y_min:y_max;
[X,Y] = meshgrid(x,y);
Z = (-d - (normal(1)*X) - (normal(2)*Y))/normal(3);
mesh(X,Y,Z)

In the following way:

   figure
plot3(A(:,1),A(:,2),A(:,3),'r.','Markersize',25)
hold on
plot3(B(:,1),B(:,2),B(:,3),'r.','Markersize',25)
plot3(C(:,1),C(:,2),C(:,3),'r.','Markersize',25)
[normal0, d0, X0, Y0, Z0] = plot_line_deviazione(A, B, C, -100, 100, -100, 100);  % new plane

dot(X-d0,normal0)

The "dot" would be the distance that I was looking for. Is this correct? How can I plo it?

Thanks!

Share Improve this question edited Mar 25 at 14:39 topzera2 asked Mar 25 at 12:34 topzera2topzera2 93 bronze badges 2
  • 1 Help me understand your question better. So, you have 2 3-point clusters(the two clusters you listed are identical, which I assume is an error). Are you asking for the distance between the planes created by each cluster? If so, note that this is only applicable when the planes are parallel, because the distance between planes is non-constant. Or, are you trying to find the distance that points in the other cluster are from the plane? If so, do you want all of the points or just a specific one? I ask because you use X in your final, dot() call, but that is not defined anywhere. – Jonathan F. Commented Mar 25 at 13:06
  • Sorry, I corrected the points. Each A,B,C point is located in the same object. I want to use one of the groups of three points to create a plane. and then from point A of the other group, calculate the smallest distance to the plane and plot it. – topzera2 Commented Mar 25 at 14:59
Add a comment  | 

2 Answers 2

Reset to default 1

The distance between a point and a plane will be the length of a line normal to the plane which intersects the point. Given that, this question has three distinct parts, as I see it.

  1. Create a plane based on 3 clustered points
  2. Calculate lines normal to that plane which also go through specified points
  3. Get the length of those lines between the points and plane intersections
  4. Plotting

A great resource for this problem can be found here. Aspects of this answer are modified snippets from that article, and I recommend giving it a read.

The function you used is well suited for plotting plane meshes. The normal calculation, as you included, is also very useful. However, since we need to calculate intersection points on the plane, having a function version of the plane is helpful. To do this, I have used scalar symbolic vars and functions.

function [normal, planefunction, solvedplanefunction, P] = generate_plane(p1, p2, p3)
    normal = cross(p1 - p2, p1 - p3); % get the normal
    
    syms x y z
    P = [x, y, z];
    planefunction = dot(normal, P - p1); % define a function version of the plane
    solvedplanefunction = solve(planefunction, z); % solve this function for z (used later)
end

The next step is to calculate the lines normal to the plane which intersect given points. Since we have the equation for the plane, and a version of it solved for z, this is pretty straightforward. First, we need to create a line that starts at our point, and has the same normal as our plane. Then we just need to solve for when that line's equation is equal to our plane's equation.

function [closest_point] = generate_closest_planar_point(normal, planefunction, P, point)
    syms t
    line = point + t*normal; % create line equation
    
    tempfunction = subs(planefunction, P, line);
    t0 = solve(tempfunction); % get value of t where the equations are equal
    closest_point = subs(line,t,t0); % plug in t0 to line equation to get point
end

With that, most of our problem is solved. Essentially, we need to calculate the norm of each of our calculated line segments, and plot the planes and line segments. So, setting the problem up with arbitrary points that comprise objects m and n:

% defining the objects
ObjectM_A = [1 1 1];
ObjectM_B = [2 0 2];
ObjectM_C = [1 3 3];

ObjectN_A = [3 3 3];
ObjectN_B = [2 2 9];
ObjectN_C = [1 4 2];

% creating the plane for m
figure()
[normalM, planeM, solvedplaneM, PM] = generate_plane(ObjectM_A, ObjectM_B, ObjectM_C);
fsurf(solvedplaneM, [-5 5 -5 5]), hold on % plot the mesh by its function

% creating the plane for n
[normalN, planeN, solvedplaneN, PN] = generate_plane(ObjectN_A, ObjectN_B, ObjectN_C);
fsurf(solvedplaneN, [-5 5 -5 5])

% for each pt in m, find the closest point on n
% print the distance
for pt = {ObjectM_A, ObjectM_B, ObjectM_C}
    cur = pt{1};
    cur_closest = double(generate_closest_planar_point(normalN, planeN, PN, cur));
    plot3([cur_closest(1) cur(1)], [cur_closest(2) cur(2)], [cur_closest(3) cur(3)], 'r','LineWidth', 2); hold on
    
    fprintf(strcat("Object M, point: [", num2str(cur), "] -> [", num2str(cur_closest), "](length = ", string(norm(cur-cur_closest)),")\n"))
end

% for each pt in n, find the closest point on m
% print the distance
for pt = {ObjectN_A, ObjectN_B, ObjectN_C}
    cur = pt{1};
    cur_closest = double(generate_closest_planar_point(normalM, planeM, PM, cur));
    plot3([cur_closest(1) cur(1)], [cur_closest(2) cur(2)], [cur_closest(3) cur(3)], 'b','LineWidth', 2); hold on
    
    fprintf(strcat("Object N, point: [", num2str(cur), "] -> [", num2str(cur_closest), "](length = ", string(norm(cur-cur_closest)),")\n"))
end

% plotting points
plot3(ObjectM_A(1),ObjectM_A(2),ObjectM_A(3),'r.','Markersize',25)
plot3(ObjectM_B(1),ObjectM_B(2),ObjectM_B(3),'r.','Markersize',25)
plot3(ObjectM_C(1),ObjectM_C(2),ObjectM_C(3),'r.','Markersize',25)
plot3(ObjectN_A(1),ObjectN_A(2),ObjectN_A(3),'b.','Markersize',25)
plot3(ObjectN_B(1),ObjectN_B(2),ObjectN_B(3),'b.','Markersize',25)
plot3(ObjectN_C(1),ObjectN_C(2),ObjectN_C(3),'b.','Markersize',25)

zlim([-5 5]);

With the following 2 support functions

1.-

function [a b c d]=plane_points_to_coeffs(A,B,C)
% input : 3 plane points
% output : [a b c d] plane coefficients
% no check whether 3 points may be parallel
x1=A(1);y1=A(2);z1=A(3);
x2=B(1);y2=B(2);z2=B(3);
x3=C(1);y3=B(2);z3=C(3);
a=(y2-y1)*(z3-z1)-(z2-z1)*(y3-y1);
b=(z2-z1)*(x3-x1)-(x2-x1)*(z3-z1);
c=(x2-x1)*(y3-y1)-(y2-y1)*(x3-x1);
d=-(a*x1+b*y1+c*z1);

2.-

function L=dist_point_plane(P,a,b,c,d)
% input point P = [x y z]
% input plane coefficients [a b c d]
% L : output distance point to plane
L=abs(a*P(1)+b*P(2)+c*P(3)-d)/(a^2+b^2+c^2)^.5;

and script

% 1st object
A1 = [-50.0098 -181.0869 -37.5881];
B1 = [-49.8410 -180.8097 -41.1383];
C1 = [-49.5887 -180.3952 -46.4462];
% 2nd object
A2 = [-50.199 -180.93 -37.249];
B2 = [-50.256 -180.823 -33.686];
C2 = [-50.341 -180.657 -28.359];
% 1st plane
[a1 b1 c1 d1]=plane_points_to_coeffs(A1,B1,C1);
% 2nd plane
[a2 b2 c2 d2]=plane_points_to_coeffs(A2,B2,C2);

% distances 1st plane to points on 2nd object
L_plane1_A2=dist_point_plane(A2,a1,b1,c1,d1)
L_plane1_B2=dist_point_plane(B2,a1,b1,c1,d1)
L_plane1_C2=dist_point_plane(C2,a1,b1,c1,d1)

% distance2 2nd plane to points on 1st object
L_plane2_A1=dist_point_plane(A1,a2,b2,c2,d2)
L_plane2_B1=dist_point_plane(B1,a2,b2,c2,d2)
L_plane2_C1=dist_point_plane(C1,a2,b2,c2,d2)

I obtain

L_plane1_A2 =
     1.035855768604593e+02
L_plane1_B2 =
     1.034733625214198e+02
L_plane1_C2 =
     1.033053734283317e+02
L_plane2_A1 =
     1.018878096607852e+02
L_plane2_B1 =
     1.017752913940908e+02
L_plane2_C1 =
     1.016071380599225e+02

I have corrected the error kindly pointed by Jonathan in plane_points_to_coeff

Now, I got the distance(point,plane) formula from

https://en.wikipedia./wiki/Distance_from_a_point_to_a_plane

In this source d has a negative sign.

I'd like to agree with Jonathan again upon another error here, that -d should be +d

because it makes more sense to calculate with ax+by+cy+d

but I am going to keep both in this my answer, for the question originator to find out whether

the formula in Wikipedia is the source of the error, or there's a 20dB distance error and d should remain -d .

With +d the results agree with Jonathan's

L_plane1_A2 =
   0.172913873214268
L_plane1_B2 =
   0.060699534174707
L_plane1_C2 =
   0.107289558913367
L_plane2_A1 =
   0.183549819138508
L_plane2_B1 =
   0.296068085832844
L_plane2_C1 =
   0.464221420001122

本文标签: plotMATLAB How to calculate the distance between a plane and a pointStack Overflow