单向后方交会实验报告C++

实验报告实验报告 班级 测绘一班班级 测绘一班 学号学号 日期 日期 2016 5 52016 5 5 目录目录 一 计算原理 3 二 算法流程 4 三 源程序 5 四 计算结果 13 五 结果分析 13 六 心得体会 13 一 计算原理 已知条件 摄影机主距 f 153 24mm x0 0 01mm y0 0 02mm 像片比例尺为 1 40000 有四对点的像点坐标与相应的地面坐标如下表 以单像空间后方交会方法 求解该像片的外方位元素 二 算法流程 1 获取已知数据 从航摄资料中差取平均航高与摄影机主距 获取控制点的 地面测量坐标并转换为地面摄影坐标 2 量测控制点的像点坐标并作系统误差改正 3 确定未知数的初始值 在竖直摄影且地面控制点大体对称分布的情况下 按如下方法确定初始值 即 n X XS 0 n Y YS 0 Z n mfZS 1 0 0 000 4 用三个角元素的初始值按下式 计算各个方向余弦值 组成旋转矩阵 R 点号像点坐标地面坐标 x mm y mm X m Y m Z m 1 86 15 68 9936589 4125273 322195 17 2 53 4082 2137631 0831324 51728 69 3 14 78 76 6339100 9724934 982386 50 410 4664 4340426 5430319 81757 31 coscos cossincossinsin sinsincoscossin sin coscos sincos cossin cossinsinsincos sinsinsincoscos 3 2 1 3 2 1 3 2 1 c c c b b b a a a 5 逐点计算像点坐标的近似值 利用未知数的近似值和控制点的地面坐标 带入共线方程式 逐点近似像点坐标的近似值 x y 6 逐点计算误差方程式的系数和常数项 组成误差方程式 7 计算法方程的系数矩阵和常数项 组成法方程式 AATLAL 8 解法方程 求得外方位元素的改正数 dddddd SSS ZYX 9 用前次迭代取得的近似值 加本次迭代的改正数 计算外方位元素的新值 K S K S K S K S K S K S K S K S K S dZZZdYYYdXXX 111 KKKKKKKKK ddd 111 10 将求得的外方位元素改正数与规定的限差比较 若小于限差 则迭代结 束 负责用新的近似值重复 4 9 直到满足要求为止 用共线方程进行空间后方交会的程序框如图所示 输入原始数据 像点坐标计算 系统误差正 确定外方位因素初始值 组成旋转矩阵 R 逐点组成误差方程式并法化 所有像点完否 是 解法方程 求外方位元素改正数 是 否 否否 2 源程序 2014113214 徐福辉 include include include using namespace std const int n 6 void inverse double c n n templatevoid transpose T1 mat1 T2 mat2 int a int b templatevoid multi T1 mat1 T2 mat2 T2 result int a int b int c int main double x 4 2 0 08616 0 06897 0 05341 0 08223 0 01479 0 07661 0 01045 0 06445 double X 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 int i j m 1 n 为迭代次数 double X0 6 0 设定未知数 XS YS ZS 初始值 double f 0 15324 摄影机主距 f 153 24mm double a 1 40000 0 像片比例尺为 1 40000 double R 3 3 0 初始化旋转矩阵 R double dayue x 8 0 用于存放像点估计值 double A 8 6 0 系数阵 double AT 6 8 0 A 的转置矩阵 计算改正后的外方位元素 外方位元素改正数是否小于限差 是 输出计算成果 计算并结束 结束并显示错误信息 迭代次数小于 n double L 8 0 存放常数项 const double pi 3 1415926535897932 double Asum 6 6 0 double jieguo2 6 0 double jieguo1 6 8 0 double sumXYZ 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 else if 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 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 4 i 第 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 初始化常数项 for i 0 i 4 i L 2 i x i 0 dayue x 2 i L 2 i 1 x i 1 dayue x 2 i 1 A 的转置矩阵 for i 0 i 8 i for j 0 j 6 j AT j i A i j 实现 A 与 AT 相乘 int k 0 for i 0 i 6 i for j 0 j 6 j Asum i j 0 for i 0 i 6 i for k 0 k 6 k for j 0 j 8 j 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 i 0 i 6 i for j 0 j 8 j jieguo1 i j 0 for i 0 i 6 i for k 0 k 8 k for j 0 j 6 j jieguo1 i k Asum i j AT j k 实现 result1 6 8 与 l 8 的相乘 得到结果放在 result2 6 中 for i 0 i 6 i jieguo2 i 0 for i 0 i 6 i for j 0 j 8 j jieguo2 i jieguo1 i j L j for i 0 i 6 i X0 i X0 i jieguo2 i ofstream f7 d A txt f7 std fixed cout 进行第 m 次迭代带得到 Xs Ys Zs 改正数结果为 n for i 0 i 6 i cout setw 12 jieguo2 i f7 setw 12 jieguo2 i cout endl 6 abs jieguo2 4 206265 0 6 abs jieguo2 5 206265 0 6 cout n 满足条件的结果为 n cout setw 12 Xs setw 12 Ys setw 12 Zs setw 12 setw 12 setw 12 endl ofstream f7 d A txt f7 std fixed cout precision 4 for i 0 i 6 i cout setw 12 X0 i f7 setw 16 X0 i f7 close double XG 6 1 for i 0 i 6 i XG i 0 jieguo2 i double AXG 8 1 V 8 1 VT 1 8 VTV 1 1 m0 D 6 6 multi A XG AXG 8 6 1 for i 0 i 8 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 cout endl ofstream f6 d what txt cout 评定完精度的结果是 endl for i 0 i 6 i for int j 0 j 6 j D i j m0 Asum i j cout setw 10 D i j f6 setw 10 D i j cout endl f6 endl cout 所得中误差为 endl for i 0 i 6 i cout sqrt D i i endl f6 close getchar return 0 void inverse double c n n int i j h k double p double q n 12 for i 0 i n i for j 0 j n j q i j c i j for i 0 i n i for j n j 12 j if i 6 j q i j 1 else q i j 0 for h k 0 k n 1 k h for i k 1 i n i if q i h 0 continue p q k h q i h for j 0 j0 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 n i p 1 0 q i i for j 0 j 12 j q i j p for i 0 i n i for j 0 j n j c i j q i j 6 t