How to fit Plane (z=ax+by+c) to 3D point cloud data?
57 visualizaciones (últimos 30 días)
Mostrar comentarios más antiguos
Swati Jain
el 22 de Sept. de 2017
Comentada: Stephen
el 29 de Jun. de 2020
Hi, I am trying to do plane fit to 3D point data. Point cloud file is attached. Here is my code I tried using least square method
clc;
clear all;
close all;
%%%%%%%%%%%%%%%%%%%%%%%%%%%reading %%%%%%%%%
arr = step_mask('Step_scan01_ex.xlsx','Sheet1', 'A:C');
%subplot(1,3,1)
x=arr(:,1);
y=arr(:,2);
z=arr(:,3)
plot3(arr(:,1),arr(:,2),arr(:,3),'.'); grid on
xlabel('x(mm)'); ylabel('y(mm)'); zlabel('z(mm)');
title('Masked plot');
%%%%%%%%%%%%%%%%%%%%%%%%%%%%least squares method %%%%%%%%%%%
[xnum,ynum]=size(arr);
%
for i = 1: xnum
xz (i) = x (i) * z (i);
yz (i) = y (i) * z (i);
xy (i) = x (i) * y (i);
end
Sxz = sum (sum (xz));
Syz = sum (sum (yz));
Sz = sum (sum (z));
Sxy = sum (sum (xy));
Sx = ynum*sum (x);
Sy = sum (y);
Sx2 = sum (x.^2);
Sy2 = sum (y.^2);
A = [Sx2 Sxy Sx; Sxy Sy2 Sy; Sx Sy xnum * ynum];
Z = [Sxz; Syz; Sz];
W = A \ Z;
w1 = W (1);
w2 = W (2);
w3 = W (3);
l = - w1/(w1^2 + w2^2 + 1)^0.5;
m = - w2/(w1^2 + w2^2 + 1)^0.5;
n = 1/(w1^2 + w2^2 + 1)^0.5;
p = w3/(w1^2 + w2^2 + 1)^0.5;
%%%%%%%%%%%%%%%%Generating and displaying the obtained plane %%%%%%%%%
z_2=zeros(7340,1);
for i = 1: xnum
z_2(i)=(-l/n)*x(i)+(-m/n)*y(i)+(p/n);
i=i+1;
end
a=-l/n;
b=-m/n;
c=p/n;
averagev = sum(sum(z_2))/(xnum * ynum);
figure;
plot3(x,y,z_2,'.'); grid on
xlabel('x'); ylabel('y'); zlabel('z(mm)');
title('Fitted Plane');
sa = abs(z-z_2);
figure;
plot3(x,y,sa,'.'); grid on
xlabel('x'); ylabel('y'); zlabel('z(mm)');
title('Substracted Plane');
Result of this code is not looking fine to me. Could anyone please suggest what is the best way to fit plane with 3D pint data?
0 comentarios
Respuesta aceptada
Star Strider
el 22 de Sept. de 2017
Try this:
arr = xlsread('Step_scan01_ex.xls');
x=arr(:,1);
y=arr(:,2);
z=arr(:,3);
DM = [x, y, ones(size(z))]; % Design Matrix
B = DM\z; % Estimate Parameters
[X,Y] = meshgrid(linspace(min(x),max(x),50), linspace(min(y),max(y),50));
Z = B(1)*X + B(2)*Y + B(3)*ones(size(X));
figure(1)
plot3(arr(:,1),arr(:,2),arr(:,3),'.')
hold on
meshc(X, Y, Z)
hold off
grid on
xlabel('x(mm)'); ylabel('y(mm)'); zlabel('z(mm)');
title('Masked plot');
grid on
text(-20, 50, 450, sprintf('Z = %.3f\\cdotX %+.3f\\cdotY %+3.0f', B))
5 comentarios
Stephen
el 29 de Jun. de 2020
Is this fitting mechanism using orthogonal distance regression, instead of just minimizing the error in Z direction?
Más respuestas (0)
Ver también
Categorías
Más información sobre Regression en Help Center y File Exchange.
Community Treasure Hunt
Find the treasures in MATLAB Central and discover how the community can help you!
Start Hunting!