摄影测量后方交会法求解外方位元素.doc

摄影测量后方交会求外方位元素09地信一班 肖明梅 解题思路：定义要用到的变量并初始化，定义一个函数用于求解旋转矩阵R，系数矩阵A,近似坐标矩阵JSZB,常数矩阵L;矩阵的转置，逆，矩阵相乘，相减，求外方位元素的近似值初值以及结果输出都定义为相应的函数。最后定义一个用于循环求解的函数（程序中xhqiujie()）,在该函数中调用之前定义的函数，求出外方位元素近似值初值，改正数Dv6,1,误差V8,1，用dowhile语句进行循环，使精度达到20m,限差低于20m，然后调用结果输出函数用于输出达到要求的结果。在主函数中创建对象的实例，引用该实例的方法即xhqiujie()函数，就可以求出外方位元素。代码：using System;using System.Collections.Generic;using System.Linq;using System.Text;namespace 摄影测量后方交会求外方位元素 class Program double , , , Xs, Ys, Zs, m, f, Sx = 0, Sy = 0; double n = 206265 / 3600;/*角元素的单位从m到度的转换系数*/ double, zuobiao = -86.15, -53.40, -14.78, 10.46 , -68.99, 82.21,-76.63,64.43,36589.41,37631.08,39100.97,40426.54, 25273.32,31324.51,24934.98,30319.81 , 2195.17,728.69,2386.50,757.31 ; double a = new double3;/*存储a1,a2.a3*/ double b = new double3;/*存储b1,b2,b3*/ double c = new double3;/*存储c1,c2,c3*/ double XX = new double4; double YY = new double4; double ZZ = new double4; double, L = new double8, 1; double, JSZB = new double2, 4;/*像点近似坐标*/ double, A = new double8, 6;/*系数阵*/ double, AT = new double6, 8; double, AB = new double6, 6; double, AN = new double6, 6; double, AL = new double6, 1; double, Dv = new double6, 1; double, AX = new double8, 1;/*系数阵与改正数矩阵的乘积*/ double, V = new double8, 1;/*中误差矩阵*/ private void Qchuzhi() for (int i = 0; i 4; i+) Sx += zuobiao2, i; for (int i = 0; i 4; i+) Sy += zuobiao3, i; for (int i = 0; i 4; i+) zuobiao0, i = zuobiao0, i / 1000; zuobiao1, i = zuobiao1, i / 1000; /*求外方位元素初始值*/ = = = 0; m = 50000; f = 0.15324; Xs = Sx / 4; Ys = Sy / 4; Zs = m * f; /*求旋转矩阵R,像点坐标近似值，常数项矩阵L,系数矩阵A*/ private void R() /旋转矩阵R a0 = Math.Cos() * Math.Cos() - Math.Sin() * Math.Sin() * Math.Sin(); a1 = -Math.Cos() * Math.Sin() - Math.Sin() * Math.Sin() * Math.Cos(); a2 = -Math.Sin() * Math.Cos(); b0 = Math.Cos() * Math.Sin(); b1 = Math.Cos() * Math.Cos(); b2 = -Math.Sin(); c0 = Math.Sin() * Math.Cos() + Math.Cos() * Math.Sin() * Math.Sin(); c1 = -Math.Sin() * Math.Sin() + Math.Cos() * Math.Sin() * Math.Cos(); c2 = Math.Cos() * Math.Cos(); for (int i = 0; i 4; i+) XXi = a0 * (zuobiao2, i - Xs) + b0 * (zuobiao3, i - Ys) + c0 * (zuobiao4, i - Zs); YYi = a1 * (zuobiao2, i - Xs) + b1 * (zuobiao3, i - Ys) + c1 * (zuobiao4, i - Zs); ZZi = a2 * (zuobiao2, i - Xs) + b2 * (zuobiao3, i - Ys) + c2 * (zuobiao4, i - Zs); /由共线条件方程式求得xo,yo的近似值 JSZB0, i = -f * (XXi) / (ZZi); JSZB1, i = -f * (YYi) / (ZZi); /常数项矩阵 Li * 2, 0 = zuobiao0, i - JSZB0, i; Li * 2 + 1, 0 = zuobiao1, i - JSZB1, i; for (int i = 0; i 4; i+) A2 * i, 0 = (a0 * f + a2 * zuobiao0, i) / ZZi; A2 * i, 1 = (b0 * f + b2 * zuobiao0, i) / ZZi; A2 * i, 2 = (c0 * f + c2 * zuobiao0, i) / ZZi; A2 * i, 3 = zuobiao1, i * Math.Sin() - (zuobiao0, i * (zuobiao0, i * Math.Cos() - zuobiao1, i * Math.Sin() / f + f * Math.Cos() * Math.Cos(); A2 * i, 4 = -f * Math.Sin() - zuobiao0, i * (zuobiao0, i * Math.Sin() + zuobiao1, i * Math.Cos() / f; A2 * i, 5 = zuobiao1, i; A2 * i + 1, 0 = (a1 * f + a2 * zuobiao1, i) / ZZi; A2 * i + 1, 1 = (b1 * f + b2 * zuobiao1, i) / ZZi; A2 * i + 1, 2 = (c1 * f + c2 * zuobiao1, i) / ZZi; A2 * i + 1, 3 = -zuobiao0, i * Math.Sin() - (zuobiao1, i * (zuobiao0, i * Math.Cos() - zuobiao1, i * Math.Sin() / f - f * Math.Sin() * Math.Cos(); A2 * i + 1, 4 = -f * Math.Cos() - zuobiao1, i * (zuobiao0, i * Math.Sin() + zuobiao1, i * Math.Cos() / f; A2 * i + 1, 5 = -zuobiao0, i; /*定义一个函数用来求矩阵的转置*/ private void T(double, A, double, B, int a, int b) for (int i = 0; i b; i+) for (int j = 0; j a; j+) Bi, j = Aj, i; /*定义一个函数用来求矩阵的乘积*/ private void Multiply(double, A, double, B, double, C, int a, int b, int c) for (int i = 0; i a; i+) for (int j = 0; j c; j+) Ci, j = 0; for (int m = 0; m b; m+) Ci, j = Ai, m * Bm, j + Ci, j; private void substract(double, A, double, B, double, C, int a, int b) for(int i=0;ia;i+) for (int j = 0; j b; j+) Ci, j = Ai, j - Bi, j; /*定义一个函数用来求矩阵的逆*/ /矩阵的求逆函数 public double, ReverseMatrix(double, dMatrix, int Level) double dMatrixValue = MatrixValue(dMatrix, Level); if (dMatrixValue = 0) return null; double, dReverseMatrix = new doubleLevel, 2 * Level; double x, c; for (int i = 0; i Level; i+) for (int j = 0; j 2 * Level; j+) if (j Level) dReverseMatrixi, j = dMatrixi, j; else dReverseMatrixi, j = 0; dReverseMatrixi, Level + i = 1; for (int i = 0, j = 0; i Level & j Level; i+, j+) if (dReverseMatrixi, j = 0) int m = i; for (; dMatrixm, j = 0; m+) ; if (m = Level - 1) return null; else / Add i-row with m-row for (int n = j; n 2 * Level; n+) dReverseMatrixi, n += dReverseMatrixm, n; / Format the i-row with 1 start x = dReverseMatrixi, j; if (x != 1) for (int n = j; n i; s-) x = dReverseMatrixs, j; for (int t = j; t = 0; i-) for (int j = i + 1; j Level; j+) if (dReverseMatrixi, j != 0) c = dReverseMatrixi, j; for (int n = j; n 2 * Level; n+) dReverseMatrixi, n -= (c * dReverseMatrixj, n); double, dReturn = new doubleLevel, Level; for (int i = 0; i Level; i+) for (int j = 0; j Level; j+) dReturni, j = dReverseMatrixi, j + Level; return dReturn; /求得矩阵行列式的值 public double MatrixValue(double, MatrixList, int Level) double, dMatrix = new doubleLevel, Level; for (int i = 0; i Level; i+) for (int j = 0; j Level; j+) dMatrixi, j = MatrixListi, j; double c, x; int k = 1; for (int i = 0, j = 0; i Level & j Level; i+, j+) if (dMatrixi, j = 0) int m = i; for (; dMatrixm - 2, j = 0; m+) ; if (m = Level) return 0; else / Row change between i-row and m-row for (int n = j; n i; s-) x = dMatrixs, j; for (int t = j; t Level; t+) dMatrixs, t -= dMatrixi, t * (x / dMatrixi, j); double sn = 1; for (int i = 0; i Level; i+) if (dMatrixi, i != 0) sn *= dMatrixi, i; else return 0; return k * sn; private void output() Console.WriteLine(外方位元素Xs：t0m, Xs); Console.WriteLine(外方位元素Ys：t0m, Ys);