摄影测量作业3.docx

1、摄影测量作业3地球科学与环境工程学院摄影测量实验报告书 课程名：摄影测量测量学 学 号： * * 日 期： 2016年4月17日一实验目的 掌握摄影测量空间后方交会的原理，利用计算机编程语言实现空间后方交会外方位元素的解算。二实验数据及实验设备 计算机windows7系统，编程软件VS2015，地面控制点在摄影测量坐标系中的坐标及其像点坐标文件data.txt。已知条件摄影机主距f=153.24mm，x0=0.01mm，y0=-0.02mm, 像片比例尺为1:40000，像片上四个像点的框标坐标及其对应的地面坐标如下表。 点号像点坐标地面坐标x(mm)y(mm)X(m)Y(m)Z(m)1-86

2、.15-68.9936589.4125273.322195.172-53.4082.2137631.0831324.51728.693-14.78-76.6339100.9724934.982386.50410.4664.4340426.5430319.81757.31以单像空间后方交会方法，求解该像片的外方位元素。三实验所用到的数学公式及程序计算步骤。单张影像的空间后方交会：利用已知地面控制点数据及相应像点坐标 根据共线方程反求影像的外方位元素。 数学模型：共线条件方程式：求解过程：（1）获取已知数据。从航摄资料中查取平均航高与摄影机主距；获取控制点的地面测量坐标并转换为地面摄影测量坐标。（

3、2）量测控制点的像点坐标并做系统改正。（3）确定未知数的初始值。在竖直摄影且地面控制点大致分布均匀的情况下，按如下方法确定初始值，即： q=w=k=0 式中:m为摄影比例尺分母；n为控制点个数。 （4）用三个角元素的初始值，计算个方向余弦，组成旋转矩阵R。（5）逐点计算像点坐标的近似值。利用未知数的近似值和控制点的地面坐标代入共线方程式，逐点计算像点坐标的近似值（x）、（y）。（6）逐点计算误差方程式的系数和常数项，组成误差方程式。（7）计算法方程的系数矩阵A和常数项L，组成法方程式。（8）解法方程，求得外方位元素的改正数dXs，dYs，dZs，dq，dw，dk。（9）用前次迭代取得的近似值，

4、加本次迭代的改正数，计算外方位元素的新值。（10）将求得的外方位元素改正数与规定的限差比较，若小于限差则迭代结束。否则用新的近似值重复（4）（9），直到满足要求为止。四程序流程图。五程序的主要源代码如下所示： #include stdafx.h#include KongJianHouFangJiaoHui.h#include KongJianHouFangJiaoHuiDlg.h#include afxdialogex.h#include CeHuiLei.h#include Matrix.h#ifdef _DEBUG#define new DEBUG_NEW#endifCMatrix xyX

5、YZ;int ButtonClickeOk = 1;void CKongJianHouFangJiaoHuiDlg:OnBnClickedOk()/导入测量数据按钮代码 CFileDialog dlgOpenFile(TRUE, _T(txt), NULL, OFN_FILEMUSTEXIST, _T(文本文件)|*.txt|(所有文件）|*.*)|); if (dlgOpenFile.DoModal() = IDCANCEL) return;/如果选择取消按钮，则退出 CString strPathName = dlgOpenFile.GetPathName();/获取选择的文件的完整路径

6、CString strFileTitle = dlgOpenFile.GetFileTitle();/获取文件名 CString strExtName = dlgOpenFile.GetFileExt();/获取文件扩展名 CStdioFile sf; /创建文件对象 /以读的形式打开文件，如果打开失败 if (!sf.Open(strPathName, CFile:modeRead) MessageBox(_T(读取文件出错！); return; CString strLine;/存放每一行文本 CStringArray strSplit; CeHuiLei SplitString; int

7、 hang = 0; while (sf.ReadString(strLine) if (hang = 1) hang+; else int n = SplitString.SplitStringArray(strLine, , strSplit); int index1 = m_list1.GetItemCount(); CString sno1 = _T(); sno1.Format(_T(%d), index1 + 1); int nR = m_list1.InsertItem(index1, sno1); for (int i = 0; i 5; i+) m_list1.SetItem

8、Text(index1, i + 1, strSpliti); sf.Close(); int Line = m_list1.GetItemCount(); CMatrix listXYZ(Line, 5); for (int i = 0; i Line; i+) for (int j = 0; j 5;j+) listXYZ(i,j)=_wtof(m_list1.GetItemText(i, j + 1); xyXYZ = listXYZ; ButtonClickeOk = 0;/表示导入数据按钮已被点击并且已导入数据CMatrix CKongJianHouFangJiaoHuiDlg:Ge

9、tX(CMatrix A, CMatrix L)/计算外方位元素改正数矩阵 CMatrix X,_A,_AA,N_AA; _A = A;/A的转置 _AA = _A*A; N_AA = _AA.Inv();/_AA的逆矩阵 X = N_AA*_A*L; return X;CMatrix CKongJianHouFangJiaoHuiDlg:GetA(CMatrix xyXYZ, double f, CMatrix XX)/计算系数矩阵A int iRow = xyXYZ.Row(); CMatrix A(2 * iRow, 6); double XS = XX(0, 0); double YS

10、 = XX(0, 1); double ZS = XX(0, 2); double w = XX(0, 4); double k = XX(0, 5); double q = XX(0, 3); double a1 = cos(q)*cos(k) - sin(q)*sin(w)*sin(k); double a2 = -cos(q)*sin(k) - sin(q)*sin(w)*cos(k); double a3 = -sin(q)*cos(w); double b1 = cos(w)*sin(k); double b2 = cos(w)*cos(k); double b3 = -sin(w)

11、; double c1 = sin(q)*cos(k) + cos(q)*sin(w)*sin(k); double c2 = -sin(q)*sin(k) + cos(q)*sin(w)*cos(k); double c3 = cos(q)*cos(w); f = f / 1000.0; for (int i = 0; i iRow; i+) double ZA = xyXYZ(i, 4); double XA = xyXYZ(i, 2); double YA = xyXYZ(i, 3); double x = (xyXYZ(i, 0)-_wtof(strx0)/1000.0; double

12、 y = (xyXYZ(i, 1)-_wtof(stry0)/1000.0; double _Z = a3*(XA - XS) + b3*(YA - YS) + c3*(ZA - ZS); A(2*i, 0) = (a1*f + a3*x) / _Z; A(2*i, 1) = (b1*f + b3*x) / _Z; A(2*i, 2) = (c1*f + c3*x) / _Z; A(2*i, 3) = y*sin(w) - (x*(x*cos(k) - y*sin(k) / f + f*cos(k)*cos(w); A(2*i, 4) = -f*sin(k) - x*(x*sin(k) + y

13、*cos(k) / f; A(2*i, 5) = y; A(2*i+1, 0) = (a2*f + a3*y) / _Z; A(2 * i + 1, 1) = (b2*f + b3*y) / _Z; A(2 * i + 1, 2) = (c2*f + c3*y) / _Z; A(2 * i + 1, 3) = -x*sin(w) - (y*(x*cos(k) - y*sin(k) / f - f*sin(k)*cos(w); A(2 * i + 1, 4) = -f*cos(k) - y/ f*(x*sin(k) + y*cos(k); A(2 * i + 1, 5) = -x; return

14、 A;CMatrix CKongJianHouFangJiaoHuiDlg:GetL(CMatrix xyXYZ, double f, CMatrix XX)/计算L矩阵 int iRow = xyXYZ.Row(); CMatrix L(2 * iRow, 1); double XS = XX(0, 0); double YS = XX(0, 1); double ZS = XX(0, 2); double w = XX(0, 4); double k = XX(0, 5); double q = XX(0, 3); double a1 = cos(q)*cos(k) - sin(q)*si

15、n(w)*sin(k); double a2 = -cos(q)*sin(k) - sin(q)*sin(w)*cos(k); double a3 = -sin(q)*cos(w); double b1 = cos(w)*sin(k); double b2 = cos(w)*cos(k); double b3 = -sin(w); double c1 = sin(q)*cos(k) + cos(q)*sin(w)*sin(k); double c2 = -sin(q)*sin(k) + cos(q)*sin(w)*cos(k); double c3 = cos(q)*cos(w); f = f

16、 / 1000.0; for (int i = 0; i iRow; i+) double ZA = xyXYZ(i, 4); double XA = xyXYZ(i, 2); double YA = xyXYZ(i, 3); double x = (xyXYZ(i, 0) - _wtof(strx0) / 1000.0; double y = (xyXYZ(i, 1) - _wtof(stry0) / 1000.0; double x0 = -f*(a1*(XA - XS) + b1*(YA - YS) + c1*(ZA - ZS) / (a3*(XA - XS) + b3*(YA - YS

17、) + c3*(ZA - ZS); double y0 = -f*(a2*(XA - XS) + b2*(YA - YS) + c2*(ZA - ZS) / (a3*(XA - XS) + b3*(YA - YS) + c3*(ZA - ZS); L(2*i, 0) =x - x0; L(2*i+1, 0) =y - y0; return L;void CKongJianHouFangJiaoHuiDlg:OnBnClickedButton1()/计算按钮 UpdateData(true); if (ButtonClickeOk = 1) MessageBox(_T(请先点击“导入已知坐标”按

18、钮导入数据！), _T(提示); return; else if (strf = _T() | strm = _T() | strx0 = _T() | stry0 = _T() MessageBox(_T(请输入数据！), _T(提示); return; else int c = strm.GetLength(); double m = _wtof(strm.Right(c - 2); int iRow = xyXYZ.Row(); double sumX = 0; double sumY = 0; double sumZ = 0; for (int i = 0; i iRow; i+) s

19、umX += xyXYZ(i, 2); sumY += xyXYZ(i, 3); sumZ += xyXYZ(i, 4); double XS = sumX / iRow; double YS = sumY / iRow; double ZS = sumZ / iRow + m*_wtof(strf) / 1000; double w = 0; double k = 0; double q = 0; CMatrix XX(1, 6); XX(0, 0) = XS; XX(0, 1) = YS; XX(0, 2) = ZS; XX(0, 3) = q; XX(0, 4) = w; XX(0, 5

20、) = k; CMatrix dX(6, 1); /迭代计算： do CMatrix A = GetA(xyXYZ, _wtof(strf), XX); CMatrix L = GetL(xyXYZ, _wtof(strf), XX); dX = GetX(A, L); for (int n = 0; n = 2.908882087e-5 | dX(4, 0) = 2.908882087e-5 | dX(5, 0) = 2.908882087e-5);/各角元素迭代计算至其改正值小于6秒,6s=2.908882087e-5 弧度。 CString strarray6;/输出外方位元素计算结果

21、strarray0 = _T(XS); strarray1 = _T(YS); strarray2 = _T(ZS); strarray3 = _T(q); strarray4 = _T(w); strarray5 = _T(k); CString strOutPut = _T(外方位元素为：rn); CString str = _T(); for (int i = 0; i 6; i+) str.Format(_T(%s=%.4f ), strarrayi, XX(0, i); strOutPut += str; if (i = 2) strOutPut += _T(rn); MessageBox(strOutPut, _T(提示), 1); 六程序运行的结果为：（1）导入数据和输入数据（2）点击计算按钮显示计算结果为：七实验总结： 通过这次的实验我学到了很多的东西，通过编程加深了对摄影测量空间后方交会相关知识的理解。在老师的严格要求下，我翻阅了很多的C+语言的书籍，看了很多的算法，这个程序最难的地方是矩阵的计算及迭代计算，以及如何将文本文档中的数据赋给程序中的矩阵，这些在以前的学习中都没有学到，感谢老师的严格要求。