通用水准网形秩亏自由网和拟稳平差程序设计..doc

MATLAB设计任何网形的秩亏自由网平差和拟稳平差程序设计 长安大学王省超 2015.5 程序介绍：程序适合于任何网形的水准网平差，原始数据输入到连个excel表格 程序界面： 原始数据录入表格： （1）DH表 （2）GXLB表 程序代码： function xsz=XSZ(num1,num2)%函数功能提取误差系数A [m1,n1]=size(num1); [m2,n2]=size(num2); n=0; for i=1:m1 %用来判断参数个数 if num1(i,2)==1 n=n+0; else n=n+1; end end xsz=zeros(m2,n); %建立系数阵，全为零 for i=1:m2 %提取系数阵 q=num2(i,1); z=num2(i,2); xsz(i,z)=1; xsz(i,q)=-1; end end %---------------常数项L-------------------------------------------------------------- function l=L(num1,num2) [m1,n1]=size(num1); [m2,n2]=size(num2); l=zeros(m2,1); for i=1:m2 %计算l q=num2(i,1); z=num2(i,2); l(i,1)=num2(i,4)-num1(z,4)+num1(q,4); end l=l*1000; %把l从米换算为毫米 end %-------------求平差权阵P------------------------------------------------------------ function p=P1(num1,num2) [m1,n1]=size(num1); [m2,n2]=size(num2); p=zeros(m2,m2); for i=1:m2 p(i,i)=1/num2(i,5); end end %------------秩亏自由网平差--------------------------------------------------------------------- function[v,Qxx]=PTPC clear all clc global H1; global K1; [num1]=xlsread(DH); K=num1(:,4); K1=K; [m1,n1]=size(num1); %用来判断num1得行列 [num2]=xlsread(GXLB); [m2,n2]=size(num2); A=XSZ(num1,num2) %误差方程系数 r=rank(A); %矩阵的秩 d=1; %秩亏数 l=L(num1,num2) %误差方程常数项单位mm P=P1(num1,num2); %权 N=A*P*A %法方程系数 W=A*P*l NN=N*N N0=NN(1:r,1:r); NN_=blkdiag(inv(N0),zeros(d,d));%广义逆求逆 Nm=N*NN_ x=Nm*W %-----------高程平差值------ H=x/1000+num1(:,4)%单位统一为M H1=H; %----精度评定------- v=A*x-l %单位mm Qxx=N*NN_*N*NN_*N msgbox(秩亏自由网普通平差完成) %---------------拟稳平差--------------------------------------------------------------------------------------- function[v,Qxx]= NWPC clear all global H2; global K1; [num1]=xlsread(DH); K=num1(:,4); K1=K;; [m1,n1]=size(num1); %用来判断num1得行列 [num2]=xlsread(GXLB); [m2,n2]=size(num2); A=XSZ(num1,num2) %误差方程系数 r=rank(A);