gusucode.com > matlab程序语言实现的水准网平差程序，使用于测绘人员 > code1/maltab程序语言实现的水准网平差程序，使用于测绘人员/malab水准网平差/adjlevel.m

    function[] = adjlevel(ofile)
%adjlevel: adjustent for leveling network
%INPUT:
%  ofile: output file name
%OUTPUT: None
%
%AUTHOR: Hua WANG
%DATE: 20/05/2009

%structures
pt = struct('id',{},'hgt',{},'flag',{});
edge = struct('sptid',{},'eptid',{},'diffhgt',{},'dist',{});

load pt
load edge

%design matrix
nobs = length(edge);    %observations number
npar = sum([pt.flag]);  %parameters number
npt=length(pt);         %total points number
B=zeros(nobs,npt);
l=zeros(nobs,1);
P=zeros(nobs,1);
for i=1:nobs
  is=find(strcmp([pt.id],edge(i).sptid)==1);
  ie=find(strcmp([pt.id],edge(i).eptid)==1);
  B(i,is)=-pt(is).flag;
  B(i,ie)= pt(ie).flag;
  l(i)=edge(i).diffhgt-(pt(ie).hgt-pt(is).hgt);
end
%remove known points from design matrix
B(:,find([pt.flag]==0))=[];

%weight matrix
P=1./[edge.dist];

%least-squares
[x,stdx,mse] = lscov(B,l,P);
v=B*x-l;

a=find([pt.flag]==1);
for i=1:npar
  pt(a(i)).hgt=x(i);
end

%write output file
fid=fopen(ofile,'w');
if(fid==-1)
  error('Can not open output file!');
end

fprintf(fid,'sigma0: %5.4f \n',sqrt(mse));
fprintf(fid,'\n id  height(m) stdx(m)\n');
sigma=zeros(npt,1);
sigma(a)=stdx;
for i=1:npt
  fprintf(fid,'%s %8.4f %8.4f\n',char(pt(i).id),pt(i).hgt,sigma(i));
end

fprintf(fid,'\nstart end  dh(m)    v(m)   dh_hat(m)\n');
for i=1:nobs
  fprintf(fid,'%4s %4s %8.4f %8.4f %8.4f\n',char(edge(i).sptid),char(edge(i).eptid),edge(i).diffhgt,v(i),edge(i).diffhgt+v(i));
end
  
fclose(fid);
return