编写一个程序根据离散点数字地形数据利用移动拟合法计算给定点地面高程.docx
《编写一个程序根据离散点数字地形数据利用移动拟合法计算给定点地面高程.docx》由会员分享，可在线阅读，更多相关《编写一个程序根据离散点数字地形数据利用移动拟合法计算给定点地面高程.docx（15页珍藏版）》请在冰豆网上搜索。
编写一个程序根据离散点数字地形数据利用移动拟合法计算给定点地面高程
编写一个程序，根据离散点数字地形数据利用移动拟合法计算给定点地面高程。
原始数据：
X=102,109,105,103,108,105,115,118,116,113
Y=110,113,115,103,105,108,104,108,113,118
Z=15,18,19,17,21,15,20,15,17,22
数据选择的个数（选择距离最近点）大于等于6
待拟合的点的坐标：
x=110，y=110
程序：
#include
#include
#include
typedefstruct
{
doubleX;
doubleY;
doubleZ;
}Point;
///////////////////////////////////////////////////////////////////读
voidRead（FILE*fp,Pointpoints[10]）
{
inti,j;
doubleDouble[10][3];
if（（fp=fopen（"Raw_Data.txt","rt"））==NULL）
printf（"路径错误,请重新输入!
"）;
for（i=0;i<10;i++）
for（j=0;j<3;j++）
fscanf（fp,"%lf",&Double[i][j]）;
for（i=0;i<10;i++）
{
points[i].X=Double[i][0];
points[i].Y=Double[i][1];
points[i].Z=Double[i][2];
}
fclose（fp）;
}
//////////////////////////////////////////////////////////////////////逆阵
voidMTPMRR（doubleMTPM[][6],doubleMTPMR[][6]）
{
inti,j,h,k,n=6;
doubleP,Q[6][12];
for（i=0;ifor（j=0;jQ[i][j]=MTPM[i][j];
for（i=0;ifor（j=n;j<12;j++）
{
if（i+6==j）
Q[i][j]=1;
else
Q[i][j]=0;
}
for（h=k=0;kfor（i=k+1;i{
if（Q[i][h]==0）
continue;
P=Q[k][h]/Q[i][h];
for（j=0;j<12;j++）
{
Q[i][j]*=P;
Q[i][j]-=Q[k][j];
}
}
for（h=k=n-1;k>0;k--,h--）//消去对角线以上的数据
for（i=k-1;i>=0;i--）
{
if（Q[i][h]==0）
continue;
P=Q[k][h]/Q[i][h];
for（j=0;j<12;j++）
{
Q[i][j]*=P;
Q[i][j]-=Q[k][j];
}
}
for（i=0;i{
P=1.0/Q[i][i];
for（j=0;j<12;j++）
Q[i][j]*=P;
}
for（i=0;ifor（j=0;jMTPMR[i][j]=Q[i][j+6];
}
//////////////////////////////////////////////////////////////////////////////////////写
voidWrite（FILE*fp,Pointpoints[10],doubleW[10][10],doubleXp,doubleYp,doubleZp）
{
inti,j;
if（（fp=fopen（"Out_Data.txt","w"））==NULL）
fprintf（fp,"输出文件不存在!
"）;
fprintf（fp,"参考点坐标是:
\n"）;
for（i=0;i<10;i++）
fprintf（fp,"%lf,%lf,%lf\n",points[i].X,points[i].Y,points[i].Z）;
fprintf（fp,"\n"）;
fprintf（fp,"权矩阵是:
\n"）;
for（i=0;i<10;i++）
{
for（j=0;j<10;j++）
fprintf（fp,"%lf",W[i][j]）;
fprintf（fp,"\n"）;
}
fprintf（fp,"待定点坐标是:
\n"）;
fprintf（fp,"%lf,%lf,%.16f",Xp,Yp,Zp）;
fclose（fp）;
}
voidM_M（doubleA[6][10],doubleB[10][10],doubleC[6][10]）
{
inti,j,k;
for（i=0;i<6;i++）
{
for（j=0;j<10;j++）
{
C[i][j]=0;
for（k=0;k<10;k++）
C[i][j]+=（A[i][k]*B[k][j]）;
}
}
}
voidM_M（doubleA[6][10],doubleB[10][6],doubleC[6][6]）
{
inti,j,k;
for（i=0;i<6;i++）
{
for（j=0;j<6;j++）
{
C[i][j]=0;
for（k=0;k<10;k++）
C[i][j]+=（A[i][k]*B[k][j]）;
}
}
}
voidM_M（doubleA[6][10],doubleB[10][1],doubleC[6][1]）
{
inti,j,k;
for（i=0;i<6;i++）
{
for（j=0;j<1;j++）
{
C[i][j]=0;
for（k=0;k<10;k++）
C[i][j]+=（A[i][k]*B[k][j]）;
}
}
}
////////////////////////////////////////////////////////main（）
voidmain（）
{
FILE*fp;
inti,j,k;
intSelect;
Pointpoints[10];
doubleM[10][6],MT[6][10],W[10][10],Z[10][1],MTP[6][10];
doubleMTPM[6][6],MTPMR[6][6],MTPZ[6][1],XX[6][1],X[10],Y[10];
doubleXp,Yp,Zp;
printf（"\n请选择操作（1--4）:
\n"）;
printf（"=======================================================================\n\n"）;
printf（"1使用移动拟合法求解待定点（110,110）的高程值\n"）;
printf（"2使用移动拟合法求解由您输入的任何待定点（X,Y）的高程值\n"）;
printf（"3使用移动拟合法求解Pending_Data.txt文件中多个待定点的高程值\n"）;
printf（"4退出:
只有选择4时才退出程序!
\n"）;
printf（"=======================================================================\n"）;
loop:
scanf（"%d",&Select）;
switch（Select）
{
case1:
{
printf（"你选择的是:
1使用移动拟合法求解待定点（110,110）的高程值!
\n"）;
Xp=110.0;
Yp=110.0;
Read（fp,points）;
for（i=0;i<10;i++）
{
X[i]=points[i].X-Xp;///////////////////////////////////X,Y存放的是差值
Y[i]=points[i].Y-Xp;
Z[i][0]=points[i].Z;
}
for（i=0;i<10;i++）
{
M[i][0]=X[i]*X[i];
M[i][1]=X[i]*Y[i];
M[i][2]=Y[i]*Y[i];
M[i][3]=X[i];
M[i][4]=Y[i];
M[i][5]=1;
}
for（i=0;i<10;i++）
for（j=0;j<10;j++）
{
if（i==j）
W[i][j]=1/（sqrt（X[i]*X[i]+Y[i]*Y[i]））;
else
W[i][j]=0;
}
for（i=0;i<6;i++）
for（j=0;j<10;j++）
MT[i][j]=M[j][i];
M_M（MT,W,MTP）;
M_M（MTP,M,MTPM）;
MTPMRR（MTPM,MTPMR）;
M_M（MTP,Z,MTPZ）;
for（i=0;i<6;i++）
{
XX[i][0]=0;
for（k=0;k<6;k++）
XX[i][0]+=MTPMR[i][k]*MTPZ[k][0];////MTPMR*MTPZ存在XX中
}
Zp=XX[5][0];
Write（fp,points,W,Xp,Yp,Zp）;
printf（"待定点的高程是:
%.16f\n\n",Zp）;
gotoloop;
}
case2:
{
printf（"你选择的是:
2使用移动拟合法求解由您输入的任何待定点（X,Y）的高程值\n"）;
printf（"请输入所求待定点的Xp坐标:
"）;
scanf（"%lf",&Xp）;
printf（"请输入所求待定点的Yp坐标:
"）;
scanf（"%lf",&Yp）;
Read（fp,points）;
for（i=0;i<10;i++）
{
X[i]=points[i].X-Xp;///////////////////////////////////X,Y存放的是差值
Y[i]=points[i].Y-Xp;
Z[i][0]=points[i].Z;
}
for（i=0;i<10;i++）
{
M[i][0]=X[i]*X[i];
M[i][1]=X[i]*Y[i];
M[i][2]=Y[i]*Y[i];
M[i][3]=X[i];
M[i][4]=Y[i];
M[i][5]=1;
}
for（i=0;i<10;i++）
for（j=0;j<10;j++）
{
if（i==j）
W[i][j]=1/（sqrt（X[i]*X[i]+Y[i]*Y[i]））;
else
W[i][j]=0;
}
for（i=0;i<6;i++）
for（j=0;j<10;j++）
MT[i][j]=M[j][i];
M_M（MT,W,MTP）;
M_M（MTP,M,MTPM）;
MTPMRR（MTPM,MTPMR）;
M_M（MTP,Z,MTPZ）;
for（i=0;i<6;i++）
{
XX[i][0]=0;
for（k=0;k<6;k++）
XX[i][0]+=MTPMR[i][k]*MTPZ[k][0];////MTPMR*MTPZ存在XX中
}
Zp=XX[5][0];
Write（fp,points,W,Xp,Yp,Zp）;
printf（"待定点的高程是:
%.16f\n\n",Zp）;
gotoloop;
}
case3:
{
printf（"你选择的是:
3使用移动拟合法求解Pending_Data.txt文件中多个待定点的高程值\n\n"）;
if（（fp=fopen（"Pending_Data.txt","rt"））==NULL）
{
printf（"待定点数据不存在!
请重新输入!
\n\n"）;
break;
}
intCounts=0,a,b;
doubleP_P[100]={0},PXp[50]={0},PYp[50]={0};
for（a=0;a<100;a++）
fscanf（fp,"%lf",&P_P[a]）;
while（P_P[Counts]!
=0）
{
if（Counts%2）
PYp[Counts/2]=P_P[Counts];
else
PXp[Counts/2]=P_P[Counts];
Counts++;
}
Read（fp,points）;
for（a=0;aif（PXp[a]!
=0&&PYp[a]!
=0）
{
for（i=0;i<10;i++）
{
X[i]=points[i].X-PXp[a];///////////////////////////////////X,Y存放的是差值
Y[i]=points[i].Y-PXp[a];
Z[i][0]=points[i].Z;
}
for（i=0;i<10;i++）
{
M[i][0]=X[i]*X[i];
M[i][1]=X[i]*Y[i];
M[i][2]=Y[i]*Y[i];
M[i][3]=X[i];
M[i][4]=Y[i];
M[i][5]=1;
}
for（i=0;i<10;i++）
for（j=0;j<10;j++）
{
if（i==j）
W[i][j]=1/（sqrt（X[i]*X[i]+Y[i]*Y[i]））;
else
W[i][j]=0;
}
for（i=0;i<6;i++）
for（j=0;j<10;j++）
MT[i][j]=M[j][i];
M_M（MT,W,MTP）;
M_M（MTP,M,MTPM）;
MTPMRR（MTPM,MTPMR）;
M_M（MTP,Z,MTPZ）;
for（i=0;i<6;i++）
{
XX[i][0]=0;
for（k=0;k<6;k++）
XX[i][0]+=MTPMR[i][k]*MTPZ[k][0];////MTPMR*MTPZ存在XX中
}
Zp=XX[5][0];
Write（fp,points,W,PXp[a],PYp[a],Zp）;
printf（"待定点（%lf,%lf）的高程是:
%.16f\n\n",PXp[a],PYp[a],Zp）;
}
gotoloop;
}
case4:
{
printf（"你选择的是4退出!
\n\n\n\n"）;
exit（0）;
}
}
}
输入数据：
建立一个Raw_Data.txt
10211015
10911318
10511519
10310317
10810521
10510815
11510420
11810815
11611317
11311822
输出结果：
参考点坐标是:
102.000000,110.000000,15.000000
109.000000,113.000000,18.000000
105.000000,115.000000,19.000000
103.000000,103.000000,17.000000
108.000000,105.000000,21.000000
105.000000,108.000000,15.000000
115.000000,104.000000,20.000000
118.000000,108.000000,15.000000
116.000000,113.000000,17.000000
113.000000,118.000000,22.000000
权矩阵是:
0.1250000.0000000.0000000.0000000.0000000.0000000.0000000.0000000.0000000.000000
0.0000000.3162280.0000000.0000000.0000000.0000000.0000000.0000000.0000000.000000
0.0000000.0000000.1414210.0000000.0000000.0000000.0000000.0000000.0000000.000000
0.0000000.0000000.0000000.1010150.0000000.0000000.0000000.0000000.0000000.000000
0.0000000.0000000.0000000.0000000.1856950.0000000.0000000.0000000.0000000.000000
0.0000000.0000000.0000000.0000000.0000000.1856950.0000000.0000000.0000000.000000
0.0000000.0000000.0000000.0000000.0000000.0000000.1280370.0000000.0000000.000000
0.0000000.0000000.0000000.0000000.0000000.0000000.0000000.1212680.0000000.000000
0.0000000.0000000.0000000.0000000.0000000.0000000.0000000.0000000.1490710.000000
0.0000000.0000000.0000000.0000000.0000000.0000000.0000000.0000000.0000000.117041
待定点坐标是:
110.000000,110.000000,17.7522235554966700
升级会员 