水准网间接平差程序设计C++.docx
《水准网间接平差程序设计C++.docx》由会员分享，可在线阅读，更多相关《水准网间接平差程序设计C++.docx（16页珍藏版）》请在冰点文库上搜索。
水准网间接平差程序设计C++
////////////////////////////////////////////////////
//visualC++6.0编译通过//
////////////////////////////////////////////////////
///////////////////////////////////////////////////
//参考资料//
//部分网络资料//
//宋力杰《测量平差程序设计》//
//姚连壁《基于matlab的控制网平差程序设计》//
///////////////////////////////////////////////////
#include
#include
#include
#include
#include
usingnamespacestd;
//////////////////////////////////////////////////////////////////////////class
classSZWPC
{
private:
intgcz_zs;//高差总数
intszd_zs;//总点数
intyz_szd_zs;//已知点数
doublem_pvv;//[pvv]
int*qsd_dh;//高差起点号
int*zd_dh;//高差终点号
char**dm;//点名地址数组
double*gcz;//观测值数组
double*szd_gc;//高程值数组
double*P;//观测值的权
double*ATPA,*ATPL;//法方程系数矩阵与自由项
double*dX;//高程改正数、平差值
double*V;//残差
doublem_mu;//单位权中误差
public:
SZWPC（）;
~SZWPC（）;
intij（inti,intj）;//对称矩阵下标计算函数
boolinverse（doublea[],intn）;//对称正定矩阵求逆（仅存下三角元素）（参考他人）
voidinputdata（char*datafile）;//输入原始数据函数
intdm_dh（char*name）;//点名转点号
voidca_H0（）;//近似高程计算函数
voidca_ATPA（）;//法方程组成函数
voidca_dX（）;//高程平差值计算函数
voidprintresult（char*resultfile）;//精度估计与平差值输出函数
doubleca_V（）;//残差计算函数
voidzxecpc（char*resultfile）;//最小二乘平差函数
};
//////////////////////////////////////////////////////////////////////
//构造函数
SZWPC:
:
SZWPC（）
{
gcz_zs=0;
szd_zs=0;
yz_szd_zs=0;
}
//////////////////////////////////////////////////////////////////////
//析构函数
SZWPC:
:
~SZWPC（）
{
if（gcz_zs>0）
{
delete[]qsd_dh;
delete[]zd_dh;
delete[]gcz;
delete[]P;
delete[]V;
}
if（szd_zs>0）
{
delete[]szd_gc;
delete[]ATPA;
delete[]ATPL;
delete[]dX;
for（inti=0;iif（dm[i]!
=NULL）delete[]（dm[i]）;
delete[]dm;
}
}
//////////////////////////////////////////////////////////////////////////
//对称矩阵下标计算函数
intSZWPC:
:
ij（inti,intj）
{
return（i>=j）?
i*（i+1）/2+j:
j*（j+1）/2+i;
}
//////////////////////////////////////////////////////////////////////////
//对称正定矩阵求逆（仅存下三角元素）（参考他人）
boolSZWPC:
:
inverse（doublea[],intn）
{
double*a0=newdouble[n];
for（intk=0;k{
doublea00=a[0];
if（a00+1.0==1.0）
{
delete[]a0;
returnfalse;
}
for（inti=1;i{
doubleai0=a[i*（i+1）/2];
if（i<=n-k-1）a0[i]=-ai0/a00;
elsea0[i]=ai0/a00;
for（intj=1;j<=i;j++）
{
a[（i-1）*i/2+j-1]=a[i*（i+1）/2+j]+ai0*a0[j];
}
}
for（i=1;i{
a[（n-1）*n/2+i-1]=a0[i];
}
a[n*（n+1）/2-1]=1.0/a00;
}
delete[]a0;
returntrue;
}
///////////////////////////////////////////////////////////////////////原始数据输入函数
voidSZWPC:
:
inputdata（char*datafile）
{
ifstreaminfile（datafile,ios:
:
in）;
if（!
infile）
{
cerr<<"Openerror!
"<}
infile>>gcz_zs>>szd_zs>>yz_szd_zs;
intunPnumber=szd_zs-yz_szd_zs;
szd_gc=newdouble[szd_zs];
dX=newdouble[szd_zs];
ATPA=newdouble[szd_zs*（szd_zs+1）/2];
ATPL=newdouble[szd_zs];
qsd_dh=newint[gcz_zs];
zd_dh=newint[gcz_zs];
gcz=newdouble[gcz_zs];
V=newdouble[gcz_zs];
P=newdouble[gcz_zs];
dm=newchar*[szd_zs];
for（inti=0;i{
dm[i]=NULL;//dm_dh函数根据dm[i]是否为NULL确定dm[i]是否为点名地址
}
charbuffer[128];//临时数组，保存从文件中读到的点名
for（i=0;i<=yz_szd_zs-1;i++）//读取已知高程数据
{
infile>>buffer;
intc=dm_dh（buffer）;
infile>>szd_gc[i];
}
for（i=0;i{
infile>>buffer;//读取高程起点名
qsd_dh[i]=dm_dh（buffer）;
infile>>buffer;//读取高程终点
zd_dh[i]=dm_dh（buffer）;
infile>>gcz[i]>>P[i];//读取高差值与路线长度
P[i]=1.0/P[i];//线路长转化为观测值的权
}
infile.close（）;
}
////////////////////////////////////////////////////////////////////点名转点号，返回点名对应的点号
intSZWPC:
:
dm_dh（char*name）
{
for（inti=0;i{
if（dm[i]!
=NULL）
{
if（strcmp（name,dm[i]）==0）returni;//将待查点名与已经存入点名数组的点名比较，若存在返回点号
}
else
{
intlen=strlen（name）;//判断点名长度
dm[i]=newchar[len+1];//为点名申请存储空间
strcpy（dm[i],name）;//待查点名是一个新的点名，将新点名的地址放到dm数组中
returni;//返回点号
}
}
return-1;//dm数组已经存满，且没有待查点名
}
////////////////////////////////////////////////////////////////////////////高程近似值计算
voidSZWPC:
:
ca_H0（）
{
for（inti=yz_szd_zs;ifor（intj=1;;j++）
{
intk=0;//计算出近似高程的点数
for（i=0;i{
intk1=qsd_dh[i];//高差起点号
intk2=zd_dh[i];//高差终点号
if（szd_gc[k1]>-10000.0&&szd_gc[k2]<-10000.0）//k1点高程或高程近似值已知，k2点高程或高程近似值未知
{
szd_gc[k2]=szd_gc[k1]+gcz[i];//计算近似高程
k++;
}
else
if（szd_gc[k1]<-10000.0&&szd_gc[k2]>-10000.0）//k2点高程或高程近似值已知，k1点高程或高程近似值未知
{
szd_gc[k1]=szd_gc[k2]-gcz[i];//计算近似高程
k++;
}
}
if（k==（szd_zs-yz_szd_zs））break;//所有的近似高程计算完成，退出
}
}
//////////////////////////////////////////////////////////////////////////
//组成法方程
voidSZWPC:
:
ca_ATPA（）
{
//intt=szd_zs;
for（inti=0;ifor（i=0;ifor（intk=0;k{
inti=qsd_dh[k];//获取点号
intj=zd_dh[k];//获取点号
doublePk=P[k];//获取权值
doublelk=gcz[k]-（szd_gc[j]-szd_gc[i]）;//获得第k个自由项
ATPL[i]-=Pk*lk;//获得法方程自由项
ATPL[j]+=Pk*lk;
ATPA[ij（i,i）]+=Pk;//获得法方程系数矩阵
ATPA[ij（j,j）]+=Pk;
ATPA[ij（i,j）]-=Pk;
}
}
//////////////////////////////////////////////////////////////////////////
//高程平差值计算
voidSZWPC:
:
ca_dX（）
{
for（inti=0;iif（!
inverse（ATPA,szd_zs））//矩阵求逆
{
cerr<<"法方程系数矩阵降秩!
"<exit（0）;//退出程序
}
for（i=0;i{
doublexi=0.0;
for（intj=0;j{
xi+=ATPA[ij（i,j）]*ATPL[j];
}
dX[i]=xi;
szd_gc[i]+=xi;//计算高程平差值
}
}
//////////////////////////////////////////////////////////////////////////
//残差计算
doubleSZWPC:
:
ca_V（）
{
doublepvv=0.0;
for（inti=0;i<=gcz_zs-1;i++）
{
intk1=qsd_dh[i];
intk2=zd_dh[i];
V[i]=szd_gc[k2]-szd_gc[k1]-gcz[i];
pvv+=V[i]*V[i]*P[i];
}
return（pvv）;
}
//////////////////////////////////////////////////////////////////////////
//原始数据和平差值输出
voidSZWPC:
:
printresult（char*resultfile）
{
doublepvv=ca_V（）;//残差计算
ofstreamoutfile（resultfile,ios:
:
out）;//以输出方式打开文件，若文件不存在，创建文件
//输出原始观测数据
outfile<"<"<outfile<<""<<"已知点数：
"<outfile<for（inti=0;i<=yz_szd_zs-1;i++）
{
outfile<<""<outfile<:
fixed）;
outfile<}
outfile<outfile<<"起始点名"<<""<<"终点点名"<<""<<"高差观测值（m）"<<""<<"两点间距离（km）"<for（i=0;i<=gcz_zs-1;i++）
{
outfile<<""<outfile<:
fixed）;
outfile<outfile<:
fixed）;
outfile<}
m_mu=sqrt（pvv/（gcz_zs-（szd_zs-yz_szd_zs）））;//计算单位权中误差
outfile<outfile<outfile<outfile<<"点名近似高程改正数高程平差值中误差"<for（i=0;i{
outfile<（2）<doubledx=dX[i];
doubleqii=ATPA[ij（i,i）];
outfile<:
fixed）;
outfile<outfile<:
fixed）;
outfile<outfile<:
fixed）;
outfile<outfile<:
fixed）;
outfile<}
//输出观测值平差值及其精度
outfile<outfile<<"起点终点观测高差ｖ"<<"高差平差值观测权中误差"<for（i=0;i<=gcz_zs-1;i++）
{
intk1=qsd_dh[i];
intk2=zd_dh[i];
doubleqii=ATPA[ij（k1,k1）];
doubleqjj=ATPA[ij（k2,k2）];
doubleqij=ATPA[ij（k1,k2）];
doubleml=sqrt（qii+qjj-2.0*qij）*m_mu;
outfile.width
（2）;
outfile<outfile.width（7）;
outfile<outfile<:
fixed）;
outfile<outfile<:
fixed）;
outfile<outfile<:
fixed）;
outfile<outfile<:
fixed）;
outfile<outfile<:
fixed）;
outfile<}
outfile.close（）;
}
//////////////////////////////////////////////////////////////////////////
//水准网最小二乘平差
voidSZWPC:
:
zxecpc（char*resultfile）
{
ca_H0（）;//近似高程计算
ca_ATPA（）;//组成法方程
ca_dX（）;//高程平差值计算
}
///////////////////////////////////////////////////////////////////////
intmain（）
{
char*datafile="算例\\Data.txt";//原始数据文件存储地址指针
char*resultfile="算例\\Result.txt";//平差结果输出地址指针
cout<cout<<"原数据文件位置：
"<cout<<"平差结果文件位置：
"<SZWPCnew_net;//定义新的对象
new_net.inputdata（datafile）;//输入原始数据
new_net.zxecpc（resultfile）;//最小二乘平差计算
new_net.printresult（resultfile）;//输出平差结果
return0;
}
升级会员 