单向后方交会实验报告C.docx
《单向后方交会实验报告C.docx》由会员分享，可在线阅读，更多相关《单向后方交会实验报告C.docx（10页珍藏版）》请在冰豆网上搜索。
单向后方交会实验报告C
单向后方交会实验报告C++
实验报告
班级：
测绘一班
学号:
日期：
2016.5.5
目录
一、计算原理.........................................3
二、算法流程.........................................4
三、源程序...........................................5
四、计算结果........................................13
五、结果分析........................................13
六、心得体会........................................13
一、计算原理
已知条件
摄影机主距f=153.24mm，x0=0.01mm，y0=0.02mm,像片比例尺为1:
40000，有四对点的像点坐标与相应的地面坐标如下表。
以单像空间后方交会方法，求解该像片的外方位元素。
二、算法流程
（1）获取已知数据。
从航摄资料中差取平均航高与摄影机主距；获取控制点的地面测量坐标并转换为地面摄影坐标。
（2）量测控制点的像点坐标并作系统误差改正。
（3）确定未知数的初始值。
在竖直摄影且地面控制点大体对称分布的情况下，按如下方法确定初始值，即
X0SX?
，Yn0SY?
，Zn0
S?
mf?
1?
Zn
?
0?
?
0?
?
0?
0
（4）用三个角元素的初始值按下式，计算各个方向余弦值，组成旋转矩阵R
a1?
cos?
cos?
?
sin?
sin?
sin?
a2?
?
cos?
sin?
?
sin?
sin?
cos?
a3?
?
sin?
cos?
b1?
cos?
sin?
b2?
cos?
cos?
b3?
?
sin?
c1?
sin?
cos?
?
cos?
sin?
sin?
c2?
?
sin?
sin?
?
cos?
sin?
cos?
c3?
cos?
cos?
（5）逐点计算像点坐标的近似值。
利用未知数的近似值和控制点的地面坐标；带入共线方程式，逐点近似像点坐标的近似值（x）、（y）。
（6）逐点计算误差方程式的系数和常数项，组成误差方程式。
（7）计算法方程的系数矩阵ATA和常数项ALL，组成法方程式。
（8）解法方程，求得外方位元素的改正数dXS、dYS、dZS、d?
、d?
、d?
。
（9）用前次迭代取得的近似值，加本次迭代的改正数，计算外方位元素的新值。
KK?
1KKK?
1K
XS?
XS?
dXS,YSK?
YSK?
1?
dYSK,ZS?
ZS?
dZS
?
K?
?
K?
1?
d?
K,?
K?
?
K?
1?
d?
K,?
K?
?
K?
1?
d?
K
（10）将求得的外方位元素改正数与规定的限差比较，若小于限差，则迭代结束。
负责用新的近似值重复（4）-（9），直到满足要求为止。
用共线方程进行空间后方交会的程序框如图所示。
否
二．源程序
//2014113214徐福辉
#include<iostream>
#include<fstream>
#include<iomanip>
usingnamespacestd;
constintn=6;
voidinverse（doublec[n][n]）;
template<typenameT1,typenameT2>voidtranspose（T1*mat1,T2*mat2,inta,intb）;template<typenameT1,typenameT2>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{cout<<"y"<<i+1<<"="<<setw（10）<<x[i][j]<<endl;}}cout<<"已知地面四个点的坐标为：
\n";for（i=0;i<4;i++）for（j=0;j<3;j++）{if（j==0）{cout<<"X"<<i+1<<"="<<X[i][j]<<"";}elseif（j==1）{cout<<"Y"<<i+1<<"="<<X[i][j]<<"";}else{cout<<"Z"<<i+1;cout<<"=";cout<<X[i][j]<<endl;}}cout<<endl;for（j=0;j<3;j++）
for（i=0;i<4;i++）
sumXYZ[j]+=X[i][j];
for（i=0;i<2;i++）
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
升级会员 