水准网间接平差程序设计(C++)

////////////////////////////////////////////////////

// visual C++6.0 编译通过//

////////////////////////////////////////////////////

///////////////////////////////////////////////////

// 参考资料//

// 部分网络资料//

// 宋力杰《测量平差程序设计》//

//姚连壁《基于matlab的控制网平差程序设计》//

///////////////////////////////////////////////////

#include

#include

#include

#include

#include

using namespace std;

//////////////////////////////////////////////////////////////////////////class

class SZWPC

{

private:

int gcz_zs; //高差总数

int szd_zs; //总点数

int yz_szd_zs; //已知点数

double m_pvv; //[pvv]

int *qsd_dh; //高差起点号

int *zd_dh; //高差终点号

char **dm; //点名地址数组

double *gcz; //观测值数组

double *szd_gc; //高程值数组

double *P; //观测值的权

double *ATPA,*ATPL; //法方程系数矩阵与自由项

double *dX; //高程改正数、平差值

double *V; //残差

double m_mu; //单位权中误差

public:

SZWPC();

~SZWPC();

int ij(int i,int j);//对称矩阵下标计算函数

bool inverse(double a[],int n);//对称正定矩阵求逆(仅存下三角元素)(参考他人）

void inputdata(char *datafile);//输入原始数据函数

int dm_dh(char *name); //点名转点号

void ca_H0(); //近似高程计算函数

void ca_ATPA(); //法方程组成函数

void ca_dX(); //高程平差值计算函数

void printresult(char *resultfile); //精度估计与平差值输出函数double ca_V(); //残差计算函数

void zxecpc(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(int i=0; i

if(dm[i]!=NULL)delete[](dm[i]);

delete []dm;

}

}

////////////////////////////////////////////////////////////////////////// // 对称矩阵下标计算函数

int SZWPC::ij(int i,int j)

{

return (i>=j)? i*(i+1)/2+j :j*(j+1)/2+i;

}

////////////////////////////////////////////////////////////////////////// // 对称正定矩阵求逆(仅存下三角元素)(参考他人）

bool SZWPC::inverse(double a[],int n)

{

double *a0=new double[n];

for(int k=0;k

{

double a00=a[0];

if(a00+1.0==1.0)

{

delete []a0;

return false;

}

for(int i=1;i

{

double ai0 = a[i*(i+1)/2];

if(i<=n-k-1)a0[i]= -ai0/a00;

else a0[i]= ai0/a00;

for(int j=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;

}