单向后方交会实验报告CDOCWord文件下载.docx
《单向后方交会实验报告CDOCWord文件下载.docx》由会员分享，可在线阅读，更多相关《单向后方交会实验报告CDOCWord文件下载.docx（19页珍藏版）》请在冰点文库上搜索。
36589.41
25273.32
2195.17
2
-53.40
82.21
37631.08
31324.51
728.69
3
-14.78
-76.63
39100.97
24934.98
2386.50
4
10.46
64.43
40426.54
30319.81
757.31
以单像空间后方交会方法，求解该像片的外方位元素。
二、算法流程
（1）获取已知数据。
从航摄资料中差取平均航高与摄影机主距；
获取控制点的地面测量坐标并转换为地面摄影坐标。
（2）量测控制点的像点坐标并作系统误差改正。
（3）确定未知数的初始值。
在竖直摄影且地面控制点大体对称分布的情况下，按如下方法确定初始值，即
，
（4）用三个角元素的初始值按下式，计算各个方向余弦值，组成旋转矩阵R
（5）逐点计算像点坐标的近似值。
利用未知数的近似值和控制点的地面坐标；
带入共线方程式，逐点近似像点坐标的近似值（x）、（y）。
（6）逐点计算误差方程式的系数和常数项，组成误差方程式。
（7）计算法方程的系数矩阵
和常数项
，组成法方程式。
（8）解法方程，求得外方位元素的改正数
。
（9）用前次迭代取得的近似值，加本次迭代的改正数，计算外方位元素的新值。
（10）将求得的外方位元素改正数与规定的限差比较，若小于限差，则迭代结束。
负责用新的近似值重复（4）-（9），直到满足要求为止。
用共线方程进行空间后方交会的程序框如图所示。
输入原始数据
像点坐标计算，系统误差正
确定外方位因素初始值
组成旋转矩阵R
逐点组成误差方程式并法化
所有像点完否
是
解法方程，求外方位元素改正数
计算改正后的外方位元素
外方位元素改正数是否小于限差
输出计算成果，计算并结束
结束并显示错误信息
迭代次数小于n
是
否
否
2．源程序
//2014113214徐福辉
#include<
iostream>
fstream>
iomanip>
usingnamespacestd;
constintn=6;
voidinverse（doublec[n][n]）;
template<
typenameT1,typenameT2>
voidtranspose（T1*mat1,T2*mat2,inta,intb）;
voidmulti（T1*mat1,T2*mat2,T2*result,inta,intb,intc）;
intmain（）
{
doublex[4][2]={-0.08616,-0.06897,-0.05341,0.08223,-0.01479,-0.07661,0.01045,0.06445};
doubleX[4][3]={36589.41,25273.32,2195.17,37631.08,31324.51,728.69,39100.97,24934.98,2386.50,40426.54,30319.81,757.31};
inti,j,m=1;
//n为迭代次数
doubleX0[6]={0};
//设定未知数（XS,YS,ZS,ψ，ω，κ）初始值
doublef=0.15324;
//摄影机主距f=153.24mm
doublea=1/40000.0;
//像片比例尺为1:
40000
doubleR[3][3]={0};
//初始化旋转矩阵R
doubledayue_x[8]={0};
//用于存放像点估计值
doubleA[8][6]={0};
//系数阵
doubleAT[6][8]={0};
//A的转置矩阵
doubleL[8]={0};
//存放常数项
constdoublepi=3.1415926535897932;
doubleAsum[6][6]={0};
doublejieguo2[6]={0};
doublejieguo1[6][8]={0};
doublesumXYZ[3]={0};
cout.precision（5）;
cout<
<
"
已知像点坐标为：
\n"
;
for（i=0;
i<
4;
i++）
for（j=0;
j<
2;
j++）
{
cout<
fixed;
if（j==0）
{
cout<
x"
<
i+1<
="
setw（10）<
x[i][j]<
"
}
else
y"
endl;
}
已知地面四个点的坐标为：
3;
X"
X[i][j]<
if（j==1）
Y"
Z"
i+1;
sumXYZ[j]+=X[i][j];
X0[i]=sumXYZ[i]/4;
//X0,Y0初始化
X0[i]=1/a*f+sumXYZ[2]/4.0;
//Z0初始化
do{
R[0][0]=cos（X0[3]）*cos（X0[5]）-sin（X0[3]）*sin（X0[4]）*sin（X0[5]）;
R[0][1]=-cos（X0[3]）*sin（X0[5]）-sin（X0[3]）*sin（X0[4]）*cos（X0[5]）;
R[0][2]=-sin（X0[3]）*cos（X0[4]）;
R[1][0]=cos（X0[4]）*sin（X0[5]）;
R[1][1]=cos（X0[4]）*cos（X0[5]）;
R[1][2]=-sin（X0[4]）;
R[2][0]=sin（X0[3]）*cos（X0[5]）+cos（X0[3]）*sin（X0[4]）*sin（X0[5]）;
R[2][1]=-sin（X0[3]）*sin（X0[5]）+cos（X0[3]）*sin（X0[4]）*cos（X0[5]）;
R[2][2]=cos（X0[3]）*cos（X0[4]）;
//第一个像点的估计值，其坐标位于X[0][0],X[0][1],X[0][2]
dayue_x[0]=-f*（R[0][0]*（X[0][0]-X0[0]）+R[1][0]*（X[0][1]-X0[1]）+R[2][0]*（X[0][2]-X0[2]））/（R[0][2]*（X[0][0]-X0[0]）+R[1][2]*（X[0][1]-X0[1]）+R[2][2]*（X[0][2]-X0[2]））;
dayue_x[1]=-f*（R[0][1]*（X[0][0]-X0[0]）+R[1][1]*（X[0][1]-X0[1]）+R[2][1]*（X[0][2]-X0[2]））/（R[0][2]*（X[0][0]-X0[0]）+R[1][2]*（X[0][1]-X0[1]）+R[2][2]*（X[0][2]-X0[2]））;
//第二个像点的估计值,其坐标位于X[1][0],X[1][1],X[1][2]
dayue_x[2]=-f*（R[0][0]*（X[1][0]-X0[0]）+R[1][0]*（X[1][1]-X0[1]）+R[2][0]*（X[1][2]-X0[2]））/（R[0][2]*（X[1][0]-X0[0]）+R[1][2]*（X[1][1]-X0[1]）+R[2][2]*（X[1][2]-X0[2]））;
dayue_x[3]=-f*（R[0][1]*（X[1][0]-X0[0]）+R[1][1]*（X[1][1]-X0[1]）+R[2][1]*（X[1][2]-X0[2]））/（R[0][2]*（X[1][0]-X0[0]）+R[1][2]*（X[1][1]-X0[1]）+R[2][2]*（X[1][2]-X0[2]））;
//第三个像点的估计值,其坐标位于X[2][0],X[2][1],X[2][2]
dayue_x[4]=-f*（R[0][0]*（X[2][0]-X0[0]）+R[1][0]*（X[2][1]-X0[1]）+R[2][0]*（X[2][2]-X0[2]））/（R[0][2]*（X[2][0]-X0[0]）+R[1][2]*（X[2][1]-X0[1]）+R[2][2]*（X[2][2]-X0[2]））;
dayue_x[5]=-f*（R[0][1]*（X[2][0]-X0[0]）+R[1][1]*（X[2][1]-X0[1]）+R[2][1]*（X[2][2]-X0[2]））/（R[0][2]*（X[2][0]-X0[0]）+R[1][2]*（X[2][1]-X0[1]）+R[2][2]*（X[2][2]-X0[2]））;
//第四个像点的估计值,其坐标位于X[3][0],X[3][1],X[3][2]
dayue_x[6]=-f*（R[0][0]*（X[3][0]-X0[0]）+R[1][0]*（X[3][1]-X0[1]）+R[2][0]*（X[3][2]-X0[2]））/（R[0][2]*（X[3][0]-X0[0]）+R[1][2]*（X[3][1]-X0[1]）+R[2][2]*（X[3][2]-X0[2]））;
dayue_x[7]=-f*（R[0][1]*（X[3][0]-X0[0]）+R[1][1]*（X[3][1]-X0[1]）+R[2][1]*（X[3][2]-X0[2]））/（R[0][2]*（X[3][0]-X0[0]）+R[1][2]*（X[3][1]-X0[1]）+R[2][2]*（X[3][2]-X0[2]））;
for（i=0;
//第i个像点估计值放在dayue_x[2*（i-1）]
A[2*i][0]=（R[0][0]*f+R[0][2]*dayue_x[2*i]）/（R[0][2]*（X[i][0]-X0[0]）+R[1][2]*（X[i][1]-X0[1]）+R[2][2]*（X[i][2]-X0[2]））;
A[2*i][1]=（R[1][0]*f+R[1][2]*dayue_x[2*i]）/（R[0][2]*（X[i][0]-X0[0]）+R[1][2]*（X[i][1]-X0[1]）+R[2][2]*（X[i][2]-X0[2]））;
A[2*i][2]=（R[2][0]*f+R[2][2]*dayue_x[2*i]）/（R[0][2]*（X[i][0]-X0[0]）+R[1][2]*（X[i][1]-X0[1]）+R[2][2]*（X[i][2]-X0[2]））;
A[2*i+1][0]=（R[0][1]*f+R[0][2]*dayue_x[2*i+1]）/（R[0][2]*（X[i][0]-X0[0]）+R[1][2]*（X[i][1]-X0[1]）+R[2][2]*（X[i][2]-X0[2]））;
A[2*i+1][1]=（R[1][1]*f+R[1][2]*dayue_x[2*i+1]）/（R[0][2]*（X[i][0]-X0[0]）+R[1][2]*（X[i][1]-X0[1]）+R[2][2]*（X[i][2]-X0[2]））;
A[2*i+1][2]=（R[2][1]*f+R[2][2]*dayue_x[2*i+1]）/（R[0][2]*（X[i][0]-X0[0]）+R[1][2]*（X[i][1]-X0[1]）+R[2][2]*（X[i][2]-X0[2]））;
A[2*i][3]=dayue_x[2*i+1]*sin（X0[4]）-（dayue_x[2*i]/f*（dayue_x[2*i]*cos（X0[5]）-dayue_x[2*i+1]*sin（X0[5]））+f*cos（X0[5]））*cos（X0[4]）;
A[2*i][4]=-f*sin（X0[5]）-dayue_x[2*i]/f*（dayue_x[2*i]*sin（X0[5]）+dayue_x[2*i+1]*cos（X0[5]））;
A[2*i][5]=dayue_x[2*i+1];
A[2*i+1][3]=-1*dayue_x[2*i]*sin（X0[4]）-（dayue_x[2*i+1]/f*（dayue_x[2*i]*cos（X0[5]）-dayue_x[2*i+1]*sin（X0[5]））-f*sin（X0[5]））*cos（X0[4]）;
A[2*i+1][4]=-1*f*cos（X0[5]）-dayue_x[2*i+1]/f*（dayue_x[2*i]*sin（X0[5]）+dayue_x[2*i+1]*cos（X0[5]））;
A[2*i+1][5]=-dayue_x[2*i];
//初始化常数项
L[2*i]=x[i][0]-dayue_x[2*i];
L[2*i+1]=x[i][1]-dayue_x[2*i+1];
//A的转置矩阵
8;
for（j=0;
6;
AT[j][i]=A[i][j];
//实现A与AT相乘
intk=0;
Asum[i][j]=0;
for（k=0;
k<
k++）
Asum[i][k]+=AT[i][j]*A[j][k];
//得到AT*A的逆矩阵存放在inverseAsum[6][6]中
inverse（Asum）;
//实现矩阵Asum[6][6]与AT[6][8]的相乘,结果存放在result1[6][8]中
jieguo1[i][j]=0;
jieguo1[i][k]+=Asum[i][j]*AT[j][k];
//实现result1[6][8]与l[8]的相乘，得到结果放在result2[6]中；
jieguo2[i]=0;
jieguo2[i]+=jieguo1[i][j]*L[j];
X0[i]=X0[i]+jieguo2[i];
ofstreamf7（"
d:
\\A.txt"
）;
f7<
std:
:
fixed;
进行第"
m<
次迭代带得到Xs,Ys,Zs,ψ，ω，κ改正数结果为：
setw（12）<
jieguo2[i];
f7<
endl<
f7.close（）;
getchar（）;
m+=1;
}while（abs（jieguo2[3]*206265.0）>
6||abs（jieguo2[4]*206265.0）>
6||abs（jieguo2[5]*206265.0）>
6）;
\n满足条件的结果为\n"
Xs"
Ys"
Zs"
ρ"
ω"
κ"
ofstreamf7（"
f7<
cout.precision（4）;
X0[i];
setw（16）<
f7.close（）;
doubleXG[6][1];
XG[i][0]=jieguo2[i];
doubleAXG[8][1],V[8][1],VT[1][8],VTV[1][1],m0,D[6][6];
multi（A,XG,AXG,8,6,1）;
i++）//计算改正数
V[i][0]=AXG[i][0]-L[i];
transpose（V,VT,1,8）;
multi（VT,V,VTV,1,8,1）;
m0=VTV[0][0]/2;
ofstreamf6（"
\\what.txt"
cout<
评定完精度的结果是"
endl;
for（intj=0;
D[i][j]=m0*Asum[i][j];
D[i][j];
f6<
f6<
所得中误差为"
sqrt（D[i][i]）<
f6.close（）;
getchar（）;
return0;
}
voidinverse（doublec[n][n]）
inti,j,h,k;
doublep;
doubleq[n][12];
n;
q[i][j]=c[i][j];
for（j=n;
12;
if（i+6==j）
q[i][j]=1;
q[i][j]=0;
for（h=k=0;
n-1;
k++,h++）
for（i=k+1;
if（q[i][h]==0）
continue;
p=q[k][h]/q[i][h];
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--）
p=q[k][h]/q[i
升级会员 