单像空间后方交会实验报告（c++版）

单

空
间

方
交
会





姓名：
学号：
时间：




2013

目录

作业务 3
二 计算原理 3
三 算法流程 7
四 源程序 8
五 计算结果 8
六 结果分析 8
七 心体会 8
八 附页 8
1 c++程序 8
2 C++程序截图 15
3 matlb程序 16
















作业务
已知条件：
摄影机距f15324mmx00y00 片例尺140000四点点坐标相应面坐标表
点号
点坐标
面坐标
x(mm)
y(mm)
X(m)
Y(m)
Z(m)
1
8615
6899
3658941
2527332
219517
2
5340
8221
3763108
3132451
72869
3
1478
7663
3910097
2493498
238650
4
1046
6443
4042654
3031981
75731

单空间方交会方法求解该片外方位元素

二 计算原理

1 获取已知数航摄资料中查取均航高摄影机距获取控制点面测量坐标转换面摄测坐标
2 测量控制点点坐标作系统误差改正
3 确定未知数初始值竖直摄影面控制点体称分布情况方法确定初始值


式中：摄影例尺分母控制点数
4 三角元素初始值式计算方余弦值组成旋转矩阵

矩阵中元素计算公式：

5 逐点计算点坐标似值利未知数似值控制点面坐标带入线方程式

逐点计算点坐标似值

6 逐点计算误差方程式系数常数项组成误差方程式

常数项计算公式：

常数项矩阵计算式：

7 计算法方程系数矩阵常数项组成法方程式
8 解法方程求外方位元素改正数
9 前次迭代取似值加次迭代改正数计算外方位元素新值

式中：代表迭代次数
10 求外方位元素改正数规定限差较限差迭代结束否新似值重复4~9满足求止
11 误差方程计算式：

单位权中误差计算式

（式中：表示余观测数）
差中6参数协数阵


终参数中误差






















三 算法流程


否

否

输入原始数
点坐标计算系统误差改正
确定外方位元素初始值
组成旋转矩阵R
逐点组成误差方程式法化
点完否
计算改正外方位元素
计算中误差输出成果结束
解法方程求外方位元素改正数
外方位元素改正数否限差

结束提示错误信息
迭代次数n


































四 源程序
源程序代码请见附页
五 计算结果
三次迭代终成果表1示
表1终计算结果
Xs（米）
Ys（米）
Zs（米）
ψ（弧度）
ω（弧度）
κ（弧度）
39795452258
27476462385
7572685988
0003987
0002114
0067578

六 结果分析
计算结果知拍摄片瞬间摄影中心面摄影测量坐标系中坐标（39795452258274764623857572685988）（单位：米）航倾角ψ0003987弧度旁倾角ω0002114弧度片旋角κ0067578弧度

表2精度评定结果
　参数
Xs
Ys
Zs
ψ
ω
κ
中误差
11073
12494
04881
00002
00002
00001

七 心体会
通次实验单张片空间方交会计算原理实现程深认识牢记种计算公式推导程基础做熟练应次实验中遇困难解决解决困难程中编程力提高涉知识印象加深
程序身足例迭代程中矩阵相输出等没采编写函数形式逐进行计算程序整体起较臃肿
编程程中确保计算准确性利matlab步编写相功计算单张相片空间方交会迭代程序验算原c++程序正确性
现编写完程序心中豪遗憾想c#编写c#暑期实时隔年没复免学会知识老师想通次实验c#进行复时间原素没达原目标
总感谢老师次锻炼提升机会

八 附页

1 c++程序

#includeiostream
using namespace std
#includefstream
#include
const int n6
void inverse(double c[n][n])
templatevoid transpose (T1*mat1T2*mat2int aint b)
templatevoid multi(T1*mat1T2*mat2T2 * resultint aint bint c)
int main()
{
double x[4][2]{008615006899005340008221001478007663001046006443}
double X[4][3]{365894125273322195173763108313245172869391009724934982386504042654303198175731}
int ijnum0num迭代次数
double X0[6]{0}设定未知数(XSYSZS三角度)初始值
double f015324摄影机距f15324mm
double a1400000片例尺140000
double R[3][3]{0}初始化旋转矩阵R
double approx_x[8]{0}存放点估计值
double A[8][6]{0}定义系数阵
double AT[6][8]{0}定义A转置矩阵
double L[8]{0}定义常数项
const double pi31415926535897932
double Asum[6][6]{0}
double result2[6]{0}
double result1[6][8]{0}
double sumXYZ[3]{0}
coutprecision(5)
cout<<已知点坐标：＼n
for(i0i<4i++)
for(j0j<2j++)
{
cout< if(j0)
{
cout< }
else
{
cout< }
}
cout<<已知面四点坐标：＼n
for(i0i<4i++)
for(j0j<3j++)
{
if(j0)
{
cout< }
else
if(j1)
{
cout< }
else
{
cout< }

}
cout< for(j0j<3j++)
for(i0i<4i++)
sumXYZ[j]+X[i][j]
for(i0i<2i++)
X0[i]sumXYZ[i]4X0Y0初始化
X0[i]1a*f+sumXYZ[2]40Z0进行初始化
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]
approx_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]))
approx_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]))
第二点估计值第2点坐标存放X[1][0]X[1][1]X[1][2]
approx_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]))
approx_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]))
第三点估计值第3点坐标存放X[2][0]X[2][1]X[2][2]
approx_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]))
approx_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]))
第四点估计值第4点坐标存放X[3][0]X[3][1]X[3][2]
approx_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]))
approx_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(i0i<4i++)
{
第i点估计值放approx_x[2*(i1)]approx_x[2*i1]
*a11*A[2*i][0](R[0][0]*f+R[0][2]*approx_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]))
*a12*A[2*i][1](R[1][0]*f+R[1][2]*approx_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]))
*a13*A[2*i][2](R[2][0]*f+R[2][2]*approx_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]))
*a21*A[2*i+1][0](R[0][1]*f+R[0][2]*approx_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]))
*a22*A[2*i+1][1](R[1][1]*f+R[1][2]*approx_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]))
*a23*A[2*i+1][2](R[2][1]*f+R[2][2]*approx_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]))
*a14*A[2*i][3]approx_x[2*i+1]*sin(X0[4])(approx_x[2*i]f*(approx_x[2*i]*cos(X0[5])approx_x[2*i+1]*sin(X0[5]))+f*cos(X0[5]))*cos(X0[4])
*a15*A[2*i][4]f*sin(X0[5])approx_x[2*i]f*(approx_x[2*i]*sin(X0[5])+approx_x[2*i+1]*cos(X0[5]))
*a16*A[2*i][5]approx_x[2*i+1]
*a24*A[2*i+1][3]1*approx_x[2*i]*sin(X0[4])(approx_x[2*i+1]f*(approx_x[2*i]*cos(X0[5])approx_x[2*i+1]*sin(X0[5]))f*sin(X0[5]))*cos(X0[4])
*a25*A[2*i+1][4]1*f*cos(X0[5])approx_x[2*i+1]f*(approx_x[2*i]*sin(X0[5])+approx_x[2*i+1]*cos(X0[5]))
*a26*A[2*i+1][5]approx_x[2*i]
}
进行常数项初始化
for(i0i<4i++)
{
L[2*i]x[i][0]approx_x[2*i]
L[2*i+1]x[i][1]approx_x[2*i+1]
}
A转置矩阵
for(i0i<8i++)
for(j0j<6j++)
{
AT[j][i]A[i][j]
}
实现AAT相
int k0
for(i0i<6i++)
for(j0j<6j++)
Asum[i][j]0
for(i0i<6i++)
for(k0k<6k++)
for(j0j<8j++)
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]中
for(i0i<6i++)
for(j0j<8j++)
result1[i][j]0
for(i0i<6i++)
for(k0k<8k++)
for(j0j<6j++)
result1[i][k]+Asum[i][j]*AT[j][k]
实现result1[6][8]l[8]相结果放result2[6]中
for(i0i<6i++)
result2[i]0
for(i0i<6i++)
for(j0j<8j++)
result2[i]+result1[i][j]*L[j]
num++
for(i0i<6i++)
{
X0[i]X0[i]+result2[i]
}
ofstream f7(d＼＼Atxt)
f7<< stdfixed
cout<<进行第< for(i0i<6i++)
{
cout< f7< }
cout< f7close()
getchar()
}while(abs(result2[3]*2062650)>6||abs(result2[4]*2062650)>6||abs(result2[5]*2062650)>6)
cout<<＼n满足限差条件结束循环终结果：＼n
cout< ofstream f7(d＼＼Atxt)
f7<< stdfixed
coutprecision(4)
for(i0i<6i++)
{
cout< }
f7close()

double XG[6][1]
for(i0i<6i++)
XG[i][0]result2[i]
double AXG[8][1]V[8][1]VT[1][8]VTV[1][1]m0D[6][6]
multi(AXGAXG861)
for( i0i<8i++) 计算改正数
V[i][0]AXG[i][0]L[i]
transpose (VVT18)
multi(VTVVTV181)
m0VTV[0][0]2
cout< ofstream f6(d＼＼whattxt)
f6<< stdfixed
for(i0i<6i++)
{
for(int j0j<6j++)
{
D[i][j]m0*Asum[i][j]
cout< f6< }
cout< f6< }
for(i0i<6i++)
cout< f6close()
屏幕输出误差方程系数阵常数项改正数




getchar()
return 0
}
void inverse(double c[n][n])
{ int ijhk
double p
double q[n][12]
for(i0i for(j0j q[i][j]c[i][j]
for(i0i for(jnj<12j++)
{
if(i+6j)
q[i][j]1
else
q[i][j]0
}
for(hk0k for(ik+1i {
if(q[i][h]0)
continue
pq[k][h]q[i][h]
for(j0j<12j++)
{ q[i][j]*p
q[i][j]q[k][j]
}
}
for(hkn1k>0kh) 消角线数
for(ik1i>0i)
{
if(q[i][h]0)
continue
pq[k][h]q[i][h]
for(j0j<12j++)
{
q[i][j]*p
q[i][j]q[k][j]
}
}
for(i0i {
p10q[i][i]
for(j0j<12j++)
q[i][j]*p
}
for(i0i for(j0j c[i][j]q[i][j+6]
}
templatevoid transpose(T1*mat1T2*mat2int aint b)
{ int ij
for(i0i for(j0j mat2[j][i]mat1[i][j]
return
}
templatevoid multi(T1*mat1T2 * mat2T2 * resultint aint bint c)
{ int ijk
for(i0i {for(j0j {result[i][j]0
for(k0k result[i][j]+mat1[i][k]*mat2[k][j]
}
}
return
}

2 C++程序截图

3 matlb程序

function [XSYSZSphiwkappa]doit
x_originalzeros(81)
Xzeros(43)
x_original[008615006899005340008221001478007663001046006443]
X[365894125273322195173763108313245172869391009724934982386504042654303198175731]
f015324摄影机距f15324mm
a0000025片例尺140000
m40000
XS38437
YS27963155
ZSm*f+15169175
phi0
w0
kappa0
Rzeros(33)
num1
dxzeros(61)
while num1|abs(dx(41)*206265)>6|abs(dx(51)*206265)>6|abs(dx(61)*206265)>6
R(11)cos(phi)*cos(kappa)sin(phi)*sin(w)*sin(kappa)
R(12)cos(phi)*sin(kappa)sin(phi)*sin(w)*cos(kappa)
R(13)sin(phi)*cos(w)
R(21)cos(w)*sin(kappa)
R(22)cos(w)*cos(kappa)
R(23)sin(w)
R(31)sin(phi)*cos(kappa)+cos(phi)*sin(w)*sin(kappa)
R(32)sin(phi)*sin(kappa)+cos(phi)*sin(w)*cos(kappa)
R(33)cos(phi)*cos(w)

x_similarzeros(81)
x_similar(11)f*(R(11)*(X(11)XS)+R(21)*(X(12)YS)+R(31)*(X(13)ZS))(R(13)*(X(11)XS)+R(23)*(X(12)YS)+R(33)*(X(13)ZS))x1
x_similar(21)f*(R(12)*(X(11)XS)+R(22)*(X(12)YS)+R(32)*(X(13)ZS))(R(13)*(X(11)XS)+R(23)*(X(12)YS)+R(33)*(X(13)ZS))y1
x_similar(31)f*(R(11)*(X(21)XS)+R(21)*(X(22)YS)+R(31)*(X(23)ZS))(R(13)*(X(21)XS)+R(23)*(X(22)YS)+R(33)*(X(23)ZS))x2
x_similar(41)f*(R(12)*(X(21)XS)+R(22)*(X(22)YS)+R(32)*(X(23)ZS))(R(13)*(X(21)XS)+R(23)*(X(22)YS)+R(33)*(X(23)ZS))y2
x_similar(51)f*(R(11)*(X(31)XS)+R(21)*(X(32)YS)+R(31)*(X(33)ZS))(R(13)*(X(31)XS)+R(23)*(X(32)YS)+R(33)*(X(33)ZS))x3
x_similar(61)f*(R(12)*(X(31)XS)+R(22)*(X(32)YS)+R(32)*(X(33)ZS))(R(13)*(X(31)XS)+R(23)*(X(32)YS)+R(33)*(X(33)ZS))y3
x_similar(71)f*(R(11)*(X(41)XS)+R(21)*(X(42)YS)+R(31)*(X(43)ZS))(R(13)*(X(41)XS)+R(23)*(X(42)YS)+R(33)*(X(43)ZS))x4
x_similar(81)f*(R(12)*(X(41)XS)+R(22)*(X(42)YS)+R(32)*(X(43)ZS))(R(13)*(X(41)XS)+R(23)*(X(42)YS)+R(33)*(X(43)ZS))y4
azeros(86)
for i14
Z(R(13)*(X(i1)XS)+R(23)*(X(i2)YS)+R(33)*(X(i3)ZS))
a(i*211)1Z*(R(11)*f+R(13)*x_similar(2*i11))
a(i*212)1Z*(R(21)*f+R(23)*x_similar(2*i11))
a(i*213)1Z*(R(31)*f+R(33)*x_similar(2*i11))
a(2*i1)1Z*(R(12)*f+R(13)*x_similar(2*i1))
a(2*i2)1Z*(R(22)*f+R(23)*x_similar(2*i1))
a(2*i3)1Z*(R(32)*f+R(33)*x_similar(2*i1))
a(i*214)x_similar(2*i1)*sin(w)(x_similar(2*i11)f*(x_similar(2*i11)*cos(kappa)x_similar(2*i1)*sin(kappa))+f*cos(kappa))*cos(w)
a(i*215)f*sin(kappa)x_similar(2*i11)f*(x_similar(2*i11)*sin(kappa)+x_similar(2*i1)*cos(kappa))
a(i*216)x_similar(2*i1)
a(2*i4)x_similar(2*i11)*sin(w)(x_similar(2*i1)f*(x_similar(2*i11)*cos(kappa)x_similar(2*i1)*sin(kappa))f*sin(kappa))*cos(w)
a(2*i5)f*cos(kappa)x_similar(2*i1)f*(x_similar(2*i11)*sin(kappa)+x_similar(2*i1)*cos(kappa))
a(2*i6)x_similar(i*211)
end
Lx_originalx_similar
dx((inv(a＇*a))*a＇)*L
num2str(a)
disp(＇改正数＇)
num2str(dx)
size(dx)
XSXS+dx(11)
YSYS+dx(21)
ZSZS+dx(31)
phiphi+dx(41)
ww+dx(51)
kappakappa+dx(61)
numnum+1
num2str(XS)
num2str(YS)
num2str(ZS)
num2str(phi)
num2str(w)
num2str(kappa)
end
文档香网(httpswwwxiangdangnet)户传

《香当网》用户分享的内容，不代表《香当网》观点或立场，请自行判断内容的真实性和可靠性！
该内容是文档的文本内容，更好的格式请下载文档