受下面两个博客启发,撰写了基于最小二乘法的平面拟合程序,在此感谢两位博主。
https://blog.csdn.net/konglingshneg/article/details/82585868
https://blog.csdn.net/shenziheng1/article/details/51175383
%平面方程的一般表达式为:Ax+By+Cz+D=0
%移项得到:z=(-A/C)*x+(-B/C)*y+(-D/c)
%令:a0=-A/C,a1=-B/C,a2=-D/c
%平面方程表达式变为:z=a0*x+a1*y+a2
%这里利用最小二乘法求出a0,a1,a2,然后再得出平面的法向量和归一化法向量
%% 平面拟合算法
clc,clear,close all
data=importdata('data.txt'); % 读取坐标文件,大家可以根据这个栗子举一反三
x=data(:,1);%定义点的x坐标
y=data(:,2);%定义点的y坐标
z=data(:,3);%定义点的z坐标
x_avr = sum(x)/length(data);% 求取法方程组系数,下同
y_avr = sum(y)/length(data);
z_avr = sum(z)/length(data);
xx_avr = sum(x.*x)/length(data);
yy_avr = sum(y.*y)/length(data);
zz_avr = sum(z.*z)/length(data);
xy_avr = sum(x.*y)/length(data);
xz_avr = sum(x.*z)/length(data);
yz_avr = sum(y.*z)/length(data);
% 法方程组的系数矩阵
A=[xx_avr xy_avr x_avr;
xy_avr yy_avr y_avr;
x_avr y_avr 1];%A矩阵,即系数矩阵
b=[xz_avr;yz_avr;z_avr];%b矩阵
xishu=A^-1*b; % 求取系数
a0=xishu(1);
a1=xishu(2);
a2=xishu(3);
vector=[a0 a1 -1];%平面法向量
nor=norm(vector);%向量的模
normalized_vector=[a0/nor a1/nor -1/nor];%平面法向量归一化
%% 图形绘制
scatter3(x,y,z,'filled')%绘制散点
hold on;
xfit = min(x):0.1:max(x);%x取值
yfit = min(y):0.1:max(y);%y取值
[XFIT,YFIT]= meshgrid (xfit,yfit);%生成X,Y点列
ZFIT=a0*XFIT+a1*YFIT+a2;%计算Z值
mesh(XFIT,YFIT,ZFIT);%绘制平面
生成的图片如下:
data.txt如下:
2.92379616171011 2.24411026297013 13.5799231123306 1.00926844826813 1.70190476178454 8.12425118188989 2.54982092942300 2.02649907973411 12.1791390980483 2.63460644130687 1.80361606750388 11.6800610851254 2.73738941072702 1.15193338338168 9.93057897159909 1.16887169102182 1.47983230710732 7.77724030336559 1.79956529819779 1.24663786967033 8.33904420540658 1.51974080570131 1.36781557656483 8.14292834109712 2.60013696044862 1.47990505132981 10.6399890748866 1.86282765492709 1.83453413816874 10.2292577243604 2.82129518885905 1.09930886065148 9.94051695967255 1.36369405660571 2.80543221983056 12.1436847727031 1.52760583304398 2.88957437944329 12.7239348044178 1.29107796076943 1.98172818493616 9.52734047634735 1.27213711741733 1.97850527680004 9.47979006523477 2.73858441528018 1.67543881964275 11.5034852894886 2.15940917473114 2.80010769283532 13.7191414279683 2.09972040367266 1.73849356224043 10.4149214940666 1.28990959644745 1.22240551058757 7.24703572465763 |