前方后方空间交会实验报告

单像空间后方交会实习报告摘要本报告旨在总结并评估笔者在单像空间后方交会实习中的经验和收获。

首先，报告介绍了单像空间后方交会实习的目的和背景。

接着，报告详细描述了实习期间所进行的实验和操作步骤。

在实习过程中，笔者遇到了一些挑战，但通过团队合作和专业指导取得了成功。

最后，报告总结了实习对于个人职业发展的重要性，并提出了改进实习体验的建议。

1. 引言单像空间后方交会是测量和分析地球或其他星球上的点的空间坐标的方法之一。

该方法通过将来自不同位置的图像投影到一个共同的平面上，并在该平面上对图像进行测量和分析，以确定点的坐标。

本实习旨在将现实生活中的实地测量和图像处理技术相结合，通过实际操作了解和掌握单像空间后方交会的原理和应用。

2. 实习过程本次实习分为三个步骤：图像获取、图像处理和空间坐标计算。

2.1 图像获取首先，为了进行后续的图像处理和分析，我们需要获取一组具有不同视角的图像。

为了实现这个目标，我们选择了一片公共景区进行实地测量。

在测量过程中，我们使用了专业的测量设备和相机，并按照一定的间隔和角度拍摄了一组图像。

这些图像将被用于后续的图像处理和分析。

2.2 图像处理在图像处理阶段，我们使用了专业的图像处理软件对获取到的图像进行处理。

首先，我们使用了相机标定算法对相机内外参数进行校准，以保证后续测量的精度和准确度。

然后，我们对每张图像进行了特征点提取和匹配，以建立图像之间的对应关系。

最后，根据所获得的对应关系，我们重建了图像场景的三维模型，并将其用于后续的空间坐标计算。

2.3 空间坐标计算在空间坐标计算阶段，我们使用了单像空间后方交会的原理，计算了每个图像特征点的空间坐标。

首先，我们将图像场景的三维模型与图像上的特征点进行对应，以确定特征点在三维空间中的位置。

然后，我们利用三角测量原理计算出特征点的三维坐标。

最后，通过对所有图像特征点的计算，我们可以得到目标点的空间坐标。

3. 实习挑战与解决在实习过程中，我们遇到了一些挑战，如图像质量、算法调优和测量误差等。

《摄 影 测 量 学》单像空间后方交会实习报告班 级： XXXX姓 名： X X X学 号：XXXXXXXXXXXXX指导教师： X X X一、实习目的1、掌握空间后方交会的定义和实现算法；2、了解摄影测量平差的基本过程；3、熟练MATLAB 等程序编写。

二、实习原理利用至少三个已知地面控制点的坐标),,(A A A Z Y X A 、),,(B B B Z Y X B 、),,(C C C Z Y X C ，与其影像上对应的三个像点的影像坐标),(a a y x a 、),(b b y x b 、),(c c y x c ，根据共线方程，反求该像片的外方位元素κωϕ、、、、、S S S Z Y X 。

共线条件方程式：将共线式线性化并取一次小值项得：三、解算过程①获取已知数据。

包括影像比例尺1/m，平均摄影距离（航空摄影的航高）H，内方位元素x0、y0、f，控制点的空间坐标X、Y、Z。

②量测控制点的像点坐标并进行必要的影像坐标系统误差改正，得到像点坐标。

③确定未知数的初始值。

单像空间后方交会必须给出待定参数的初始值，在竖直航空摄影且地面控制点大体对称分布的情况下，Xs0和Ys0为均值，Zs0为航高，φ、ω、κ的初值都设为0。

或者κ的初值可在航迹图上找出或根据控制点坐标通过坐标正反变换求出。

④计算旋转矩阵R。

利用角元素近似值计算方向余弦值，组成R阵。

⑤逐点计算像点坐标的近似值。

利用未知数的近似值按共线条件式计算控制点像点坐标的近似值（x），(y)。

⑥逐点计算误差方程式的系数和常数项，组成误差方程式。

⑦计算法方程的系数矩阵ATA与常数项ATL，组成法方程式。

⑧解求外方位元素。

根据法方程，解求外方位元素改正数，并与相应的近似值求和，得到外方位元素新的近似值。

⑨检查计算是否收敛。

将所求得的外方位元素的改正数与规定的限差比较，通常对φ，ω，κ的改正数△φ，△ω，△κ给予限差，通常为0.000001弧度，当3个改正数均小于0.000001弧度时，迭代结束。

摄影测量学实验报告实验一、单像空间后方交会学院：建测学院班级：测绘082姓名：肖澎学号： 15一．实验目的1.深入了解单像空间后方交会的计算过程；2.加强空间后方交会基本公式和误差方程式，法线方程式的记忆；3.通过上机调试程序加强动手能力的培养。

二．实验原理以单幅影像为基础，从该影像所覆盖地面范围内若干控制点和相应点的像坐标量测值出发，根据共线条件方程，求解该影像在航空摄影时刻的相片外方位元素。

三．实验内容1.程序图框图2.实验数据（1）已知航摄仪内方位元素f＝153.24mm，Xo＝Yo＝0。

限差0.1秒（2）已知4对点的影像坐标和地面坐标：3.实验程序using System;using System.Collections.Generic;using System.Linq;using System.Text;namespace ConsoleApplication3{class Program{static void Main(){//输入比例尺，主距，参与平参点的个数Console.WriteLine("请输入比例尺分母m：\r");string m1 = Console.ReadLine();double m = (double)Convert.ToSingle(m1);Console.WriteLine("请输入主距f：\r");string f1 = Console.ReadLine();double f = (double)Convert.ToSingle(f1);Console.WriteLine("请输入参与平差控制点的个数n：\r");string n1 = Console.ReadLine();int n = (int)Convert.ToSingle(n1);//像点坐标的输入代码double[] arr1 = new double[2 * n];//1.像点x坐标的输入for (int i = 0; i < n; i++){Console.WriteLine("请输入已进行系统误差改正的像点坐标的x{0}值：\r", i+1);string u = Console.ReadLine();for (int j = 0; j < n; j += 2){arr1[j] = (double)Convert.ToSingle(u);}}//2.像点y坐标的输入for (int i = 0; i < n; i++){Console.WriteLine("请输入已进行系统误差改正的像点坐标的y{0}值：\r", i+1);string v = Console.ReadLine();for (int j = 1; j < n; j += 2){arr1[j] = (double)Convert.ToSingle(v);}}//控制点的坐标输入代码double[,] arr2 = new double[n, 3];//1.控制点X坐标的输入for (int j = 0; j < n; j++){Console.WriteLine("请输入控制点在地面摄影测量坐标系的坐标的X{0}值：\r", j+1);string u = Console.ReadLine();arr2[j , 0] = (double)Convert.ToSingle(u);}//2.控制点Y坐标的输入for (int k = 0; k < n; k++){Console.WriteLine("请输入控制点在地面摄影测量坐标系的坐标的Y{0}值：\r", k+1);string v = Console.ReadLine();arr2[k , 1] = (double)Convert.ToSingle(v);}//3.控制点Z坐标的输入for (int p =0; p < n; p++){Console.WriteLine("请输入控制点在地面摄影测量坐标系的坐标的Z{0}值：\r", p+1);string w = Console.ReadLine();arr2[p , 2] = (double)Convert.ToSingle(w);}//确定外方位元素的初始值//1.确定Xs的初始值：double Xs0 = 0;double sumx = 0;for (int j = 0; j < n; j++){double h = arr2[j, 0];sumx += h;}Xs0 = sumx / n;//2.确定Ys的初始值：double Ys0 = 0;double sumy = 0;for (int j = 0; j < n; j++){double h = arr2[j, 1];sumy += h;}Ys0 = sumy / n;//3.确定Zs的初始值：double Zs0 = 0;double sumz = 0;for (int j = 0; j <= n - 1; j++){double h = arr2[j, 2];sumz += h;}Zs0 = sumz / n;doubleΦ0 = 0;doubleΨ0 = 0;double K0 = 0;Console.WriteLine("Xs0,Ys0,Zs0,Φ0,Ψ0,K0的值分别是：{0},{1},{2},{3},{4},{5}", Xs0, Ys0, Zs0, 0, 0, 0);//用三个角元素的初始值按（3-4-5）计算各方向余弦值，组成旋转矩阵,此时的旋转矩阵为单位矩阵I：double[,] arr3 = new double[3, 3];for (int i = 0; i < 3; i++)arr3[i, i] = 1;}double a1 = arr3[0, 0]; double a2 = arr3[0, 1]; double a3 = arr3[0, 2];double b1 = arr3[1, 0]; double b2 = arr3[1, 1]; double b3 = arr3[1, 2];double c1 = arr3[2, 0]; double c2 = arr3[2, 1]; double c3 = arr3[2, 2];/*利用线元素的初始值和控制点的地面坐标，代入共线方程（3-5-2），* 逐点计算像点坐标的近似值*///1.定义存放像点近似值的数组double[] arr4 = new double[2 * n];//----------近似值矩阵//2.逐点像点坐标计算近似值//a.计算像点的x坐标近似值（x）for (int i = 0; i < 2 * n; i += 2){for (int j = 0; j < n; j++){arr4[i] = -f * (a1 * (arr2[j, 0] - Xs0) + b1 * (arr2[j, 1] - Ys0) + c1 * (arr2[j, 2] - Zs0)) / (a3 * (arr2[j, 0] - Xs0) + b3 * (arr2[j, 1] - Ys0) + c3 * (arr2[j, 2] - Zs0)); }}//b.计算像点的y坐标近似值（y）for (int i = 1; i < 2 * n; i += 2){for (int j = 0; j < n; j++){arr4[i] = -f * (a2 * (arr2[j, 0] - Xs0) + b2 * (arr2[j, 1] - Ys0) + c2 * (arr2[j, 2] - Zs0)) / (a3 * (arr2[j, 0] - Xs0) + b3 * (arr2[j, 1] - Ys0) + c3 * (arr2[j, 2] - Zs0)); }}//逐点计算误差方程式的系数和常数项，组成误差方程：double[,] arr5 = new double[2 * n, 6]; //------------系数矩阵（A）//1.计算dXs的系数for (int i = 0; i < 2 * n; i += 2){arr5[i, 0] = -1 / m; //-f/H == -1/m}//2.计算dYs的系数for (int i = 1; i < 2 * n; i += 2){arr5[i, 1] = -1 / m; //-f/H == -1/m}//3.a.计算误差方程式Vx中dZs的系数for (int i = 0; i < 2 * n; i += 2)arr5[i, 2] = -arr1[i] / m * f;}//3.b.计算误差方程式Vy中dZs的系数for (int i = 1; i < 2 * n; i += 2){arr5[i, 2] = -arr1[i] / m * f;}//4.a.计算误差方程式Vx中dΦ的系数for (int i = 0; i < 2 * n; i += 2){arr5[i, 3] = -f * (1 + arr1[i] * arr1[i] / f * f);}//4.a.计算误差方程式Vy中dΦ的系数for (int i = 1; i < 2 * n; i += 2){arr5[i, 3] = -arr1[i - 1] * arr1[i] / f;}//5.a.计算误差方程式Vx中dΨ的系数for (int i = 0; i < 2 * n; i += 2){arr5[i, 4] = -arr1[i] * arr1[i + 1] / f;}//5.b.计算误差方程式Vy中dΨ的系数for (int i = 1; i < 2 * n; i += 2){arr5[i, 4] = -f * (1 + arr1[i] * arr1[i] / f * f);}//6.a.计算误差方程式Vx中dk的系数for (int i = 0; i < 2 * n; i += 2){arr5[i, 5] = arr1[i + 1];}//6.b.计算误差方程式Vy中dk的系数for (int i = 1; i < 2 * n; i += 2){arr5[i, 5] = -arr1[i - 1];}//定义外方位元素组成的数组double[] arr6 = new double[6];//--------------------外方位元素改正数矩阵（X）//定义常数项元素组成的数组double[] arr7 = new double[2 * n];//-----------------常数矩阵（L）//计算lx的值for (int i = 0; i < 2 * n; i += 2)arr7[i] = arr1[i] - arr4[i]; //将近似值矩阵的元素代入}//计算ly的值for (int i = 1; i <= 2 * (n - 1); i += 2){arr7[i] = arr1[i] - arr4[i]; //将近似值矩阵的元素代入}/* 对于所有像点的坐标观测值，一般认为是等精度量测，所以权阵P为单位阵.所以X=(ATA)-1ATL *///1.计算ATdouble[,] arr5T = new double[6, 2 * n];for (int i = 0; i < 6; i++){for (int j = 0; j < 2 * n; j++){arr5T[i, j] = arr5[j, i];}}//A的转置与A的乘积,存放在arr5AA中double[,] arr5AA = new double[6, 6];for (int i = 0; i < 6; i++){for (int j = 0; j < 6; j++){arr5AA[i, j] = 0;for (int l = 0; l < 2 * n; l++){arr5AA[i, j] += arr5T[i, l] * arr5[l, j];}}}nijuzhen(arr5AA);//arr5AA经过求逆后变成原矩阵的逆矩阵//arr5AA * arr5T存在arr5AARATdouble[,] arr5AARAT = new double[6, 2 * n];for (int i = 0; i < 6; i++){for (int j = 0; j < 2 * n; j++){arr5AARAT[i, j] = 0;for (int p = 0; p < 6; p++){arr5AARAT[i, j] += arr5AA[i, p] * arr5T[p, j];}}}//计算arr5AARAT x L，存在arrX中double[] arrX = new double[6];for (int i = 0; i < 6; i++){for (int j = 0; j < 1; j++){arrX[i] = 0;for (int vv = 0; vv < 6; vv++){arrX[i] += arr5AARAT[i, vv] * arr7[vv];}}}//计算外方位元素值double Xs, Ys, Zs, Φ, Ψ, K;Xs = Xs0 + arrX[0];Ys = Ys0 + arrX[1];Zs = Zs0 + arrX[2];Φ = Φ0 + arrX[3];Ψ = Ψ0 + arrX[4];K = K0 + arrX[5];for (int i = 0; i <= 2; i++){Xs += arrX[0];Ys += arrX[1];Zs += arrX[2];Φ += arrX[3];Ψ += arrX[4];K += arrX[5];}Console.WriteLine("Xs,Ys,Zs,Φ,Ψ,K的值分别是：{0},{1},{2},{3},{4},{5}", Xs0, Ys0, Zs0, Φ, Ψ, K);Console.Read();}//求arr5AA的逆矩public static double[,] nijuzhen(double[,] a) {double[,] B = new double[6, 6];int i, j, k;int row = 0;int col = 0;double max, temp;int[] p = new int[6];for (i = 0; i < 6; i++){p[i] = i;B[i, i] = 1;}for (k = 0; k < 6; k++){//找主元max = 0; row = col = i;for (i = k; i < 6; i++){for (j = k; j < 6; j++){temp = Math.Abs(a[i, j]);if (max < temp){max = temp;row = i;col = j;}}}//交换行列，将主元调整到k行k列上if (row != k){for (j = 0; j < 6; j++){temp = a[row, j];a[row, j] = a[k, j];a[k, j] = temp;temp = B[row, j];B[row, j] = B[k, j];B[k, j] = temp;i = p[row]; p[row] = p[k]; p[k] = i; }if (col != k){for (i = 0; i < 6; i++){temp = a[i, col];a[i, col] = a[i, k];a[i, k] = temp;}}//处理for (j = k + 1; j < 6; j++){a[k, j] /= a[k, k];}for (j = 0; j < 6; j++){B[k, j] /= a[k, k];a[k, k] = 1;}for (j = k + 1; j < 6; j++){for (i = 0; j < k; i++){a[i, j] -= a[i, k] * a[k, j];}for (i = k + 1; i < 6; i++){a[i, j] -= a[i, k] * a[k, j];}}for (j = 0; j < 6; j++){for (i = 0; i < k; i++){B[i, j] -= a[i, k] * B[k, j];}for (i = k + 1; i < 6; i++){B[i, j] -= a[i, k] * B[k, j];}for (i = 0; i < 6; i++) {a[i, k] = 0;a[k, k] = 1;}}//恢复行列次序for (j = 0; j < 6; j++){for (i = 0; i < 6; i++) {a[p[i], j] = B[i, j]; }}for (i = 0; i < 6; i++){for (j = 0; j < 6; j++) {a[i, j] = a[i, j];}}return a;}4．实验结果四．实验总结此次实验让我深入了解单像空间后方交会的计算过程，加强了对空间后方交会基本公式和误差方程式，法线方程式的记忆。

摄影测量学单像空间后方交会编程实习报告实习背景在本次实习中，我们学习了摄影测量学单像空间后方交会的编程实现。

这是一种通过计算影像中各点的坐标来确定被摄物的三维坐标的方法，应用广泛于测绘、地理信息、建筑等领域。

本次实习采用 MATLAB 软件进行编程，目的是将理论知识应用到实际操作中，让我们更深入理解摄影测量学单像空间后方交会的原理和应用。

实习内容理论部分首先，我们在工作室进行了理论部分的学习。

老师讲解了单像空间后方交会的原理，以及如何通过影像坐标、相机外方位元素、像点坐标和像平面坐标等参数来计算被摄物的三维坐标。

在理论部分的学习过程中，我们通过公式的推导和实例分析，更加深入地理解了单像空间后方交会的原理。

实践部分实践部分是本次实习的重头戏。

我们利用 MATLAB 软件进行了单像空间后方交会的编程实现，具体步骤如下：1.输入相机外方位元素通过读取文本文件，将相机外方位元素（相机在拍摄时的姿态、位置等参数）输入到 MATLAB 中。

2.输入影像坐标通过读取文本文件，将影像中的像点坐标输入到 MATLAB 中。

3.计算像平面坐标利用相机内定标参数，将像点坐标转化为像平面坐标。

4.计算被摄物三维坐标根据单像空间后方交会的原理，利用相机外方位元素、像平面坐标和像点坐标等参数，计算被摄物的三维坐标。

5.输出结果将计算结果输出到文本文件中，以便后续的数据处理和分析。

在实际操作中，我们首先编写了 MATLAB 脚本文件，根据上述步骤逐步实现了单像空间后方交会的计算过程。

然后，我们利用自己拍摄的实际照片进行实验，将相机外方位元素和像点坐标输入到程序中，最终得到了被摄物的三维坐标结果。

实习收获通过本次实习，我从理论到实践，更深入地理解了摄影测量学单像空间后方交会的原理和应用，同时也掌握了 MATLAB 的编程技能。

在实践中，我遇到了许多问题，包括数据的输入输出、代码的调试和结果的分析等等。

通过和同学的讨论和老师的指导，我不仅解决了这些问题，还对摄影测量学的应用有了更深入的认识。

空间后方—前方交会的原理
以空间后方—前方交会的原理为题，我来为大家描述一下。

空间后方—前方交会是一种用于确定目标位置的方法，常用于航空、导航、测绘等领域。

它利用人眼的立体视觉和视差效应，通过观察目标在不同视角下的位置变化，来推断目标的实际位置。

这种方法可以较精确地确定目标的距离和方位，尤其适用于远距离观测。

在进行空间后方—前方交会时，我们首先需要选择两个观测点，它们之间的距离应足够远，以便产生明显的视差效应。

然后，我们分别在这两个观测点上观察目标，并记录下目标在两个观测点的位置。

接下来，我们需要测量观测点之间的距离，并确定观测点与目标之间的夹角。

这些数据将用于计算目标的实际位置。

通过对两个观测点的位置和距离进行几何分析，我们可以得到目标相对于观测点的位移向量。

然后，我们再将这个位移向量与观测点之间的夹角结合起来，就可以计算出目标相对于观测点的实际位置。

空间后方—前方交会的原理基于视差效应，即当我们观察远处的目标时，由于两只眼睛的视角不同，目标在两只眼睛中的位置也会有所不同。

通过比较这两个位置的差异，我们就可以推断出目标的实际位置。

总的来说，空间后方—前方交会是一种利用视差效应来确定目标位
置的方法。

它可以在远距离观测中提供较为准确的测量结果，具有广泛的应用前景。

空间后交-前交程序设计(实验报告)姓名：班级：学号：时间：空间后交-前交程序设计一、实验目的用 C 、VB或MATLAB语言编写空间后方交会-空间前方交会程序⑴提交实习报告：程序框图、程序源代码、计算结果、体会⑵计算结果：像点坐标、地面坐标、单位权中误差、外方位元素及其精度二、实验数据f=150.000mm，x0=0，y0=0三、实验思路1.利用空间后方交会求左右像片的外方位元素(1).获取m(于像片中选取两点，于地面摄影测量坐标系中选取同点，分别计算距离，距离比值即为m)，x,y,f,X,Y,Z(2).确定未知数初始值Xs,Ys,Zs,q,w,k(3).计算旋转矩阵R(4).逐点计算像点坐标的近似值（x），（y）(5).组成误差方程式(6).组成法方程式(7).解求外方位元素(8).检查是否收敛，即将求得的外方位元素的改正数与规定限差比较，小于限差即终止；否则用新的近似值重复步骤（3）-（7）2.利用求出的外方位元素进行空间前交，求出待定点地面坐标(1).用各自像片的角元素计算出左、右像片的方向余弦值，组成旋转矩阵R1,R2(2).根据左、右像片的外方位元素，计算摄影基线分量Bx，By，Bz(3).计算像点的像空间辅助坐标（X1,Y1,Z1）和（X2,Y2,Z2）(4).计算点投影系数N1和N2(5).计算未知点的地面摄影测量坐标四、实验过程⑴程序框图函数AandL%求间接平差时需要的系数%%%已知%a=像点坐标x,b=像点坐标y，f内方位元素主距%φ=q,ψ=w,κ=k%像空间坐标系X,Y,Z%地面摄影测量坐标系Xs,Ys,Zsfunction [A1,L1,A2,L2]=AandL(a,b,f,q,w,k,X,Y,Z,Xs,Ys,Zs) %%%%%%%%%%%选择矩阵元素a1=cos(q)*cos(k)-sin(q)*sin(w)*sin(k);a2=-cos(q)*sin(k)-sin(q)*sin(w)*cos(k);a3=-sin(q)*cos(w);b1=cos(w)*sin(k);b2=cos(w)*cos(k);b3=-sin(w);c1=sin(q)*cos(k)+cos(q)*sin(w)*sin(k);c2=-sin(q)*sin(k)+cos(q)*sin(w)*cos(k);c3=cos(q)*cos(w);%%%%%%%共线方程的分子分母X_=a1*(X-Xs)+b1*(Y-Ys)+c1*(Z-Zs);Y_=a2*(X-Xs)+b2*(Y-Ys)+c2*(Z-Zs);Z_=a3*(X-Xs)+b3*(Y-Ys)+c3*(Z-Zs);%%%%%%%近似值x=-f*X_/Z_;y=-f*Y_/Z_;%%%%%%%A组成L组成a11=1/Z_*(a1*f+a3*x);a12=1/Z_*(b1*f+b3*x);a13=1/Z_*(c1*f+c3*x);a21=1/Z_*(a2*f+a3*y);a22=1/Z_*(b2*f+b3*y);a23=1/Z_*(c2*f+c3*y);a14=y*sin(w)-(x/f*(x*cos(k)-y*sin(k))+f*cos(k))*cos(w);a15=-f*sin(k)-x/f*(x*sin(k)+y*cos(k)); a16=y;a24=-x*sin(w)-(y/f*(x*cos(k)-y*sin(k))-f*sin(k))*cos(w);a25=-f*cos(k)-y/f*(x*sin(k)+y*cos(k));a26=-x;lx=a-x;ly=b-y;%%%%%%%%%组成一个矩阵，并返回A1=[a11,a12,a13,a14,a15,a16];A2=[a21,a22,a23,a24,a25,a26];L1=lx;L2=ly;函数deg2dms%%%%%%%%角度转度分秒function y=deg2dms(x)a=floor(x);b=floor((x-a)*60);c=(x-a-b/60)*3600;y=a+(b/100)+(c/10000);函数dms2deg%%%%%度分秒转度function y=dms2deg(x)a=floor(x);b=floor((x-a)*100);c=(x-a-b/100)*10000;y=a+b/60+c/3600;函数ok%%%%%%%%%%%%%%目的是为了保证各取的值的有效值%%xy为n*1,a为1*nfunction result=ok(xy,a)format short gi=size(xy,1);for n=1:io=xy(n)-floor(xy(n,1));o=round(o*(10^a(n)))/(10^a(n));xy(n,1)=floor(xy(n,1))+o;endformat long gresult=xy;函数rad2dmsxy%%%%求度分秒表现形式的三个外方位元素，三个角度function xydms=rad2dmsxy(xy)[a,b,c,d,e,f]=testvar(xy);d=deg2dms(rad2deg(d));e=deg2dms(rad2deg(e));f=deg2dms(rad2deg(f));xydms=[a,b,c,d,e,f]';函数spacehoujiao%%%%%%%空间后交%%% f%%输入p（2*n,1）%%像点坐标x,y,X,Y,Z,均为（n,1）function [xy,m,R]=spacehoujiao(p,x,y,f,X,Y,Z)format long;%%%%%权的矢量化,这是等精度时的，如果非，将函数参数改为P P=diag(p);%%求nj=size(X,2);%%初始化Xs=0;Ys=0;Zs=0;for n=1:jXs=Xs+X(n);Ys=Ys+Y(n);Zs=Zs+Z(n);endSx=sqrt((x(2)-x(1))^2+(y(2)-y(1))^2);%%%%两像点之间距离Sd=sqrt((X(2)-X(1))^2+(Y(2)-Y(1))^2);%%%%两地面控制点之间距离m=Sd/Sx; %%%%图像比例系数Xs=Xs/j;Ys=Ys/j;Zs=m*f+Zs/j;m0=0;q=0;w=0;k=0;i=0;a=rand(2*j,6);l=rand(2*j,1);%%%%for n=1:j[a(2*n-1,:),l(2*n-1,1),a(2*n,:),l(2*n,1)]=AandL(x(n),y(n),f,q,w,k,X(n),Y(n),Z(n ),Xs,Ys,Zs);enddet=inv(a'*P*a)*transpose(a)*P*l;%%%%%%%%%循环体while 1%%%%%%%%%%%%%%%%[dXs,dYs,dZs,dq,dw,dk]=testvar(det);detXs=abs(dXs);detYs=abs(dYs);detZs=abs(dZs);detq=abs(dq);detw=abs(dw);detk=abs(dk);%%%%%%%%%if((detXs<0.01)&&(detYs<0.01)&&(detZs<0.01)&&(detq<pi/648000)&&(detw<pi/648000)&& (detk<pi/648000))break;elseV=(a*det-l);Q=inv(a'*P*a);m0=m0+sqrt((V'*P*V)/(2*j-6));%%m0需要每次的改正数算出来相加%%%Xs=Xs+dXs;Ys=Ys+dYs;Zs=Zs+dZs;q=q+dq;w=w+dw;k=k+dk;%%%for n=1:j[a(2*n-1,:),l(2*n-1,1),a(2*n,:),l(2*n,1)]=AandL(x(n),y(n),f,q,w,k,X(n),Y(n),Z(n ),Xs,Ys,Zs);enddet=inv(a'*P*a)*transpose(a)*P*l;i=i+1;%%%%end%%%end[dXs,dYs,dZs,dq,dw,dk]=testvar(det);detXs=abs(dXs);detYs=abs(dYs);detZs=abs(dZs);detq=abs(dq);detw=abs(dw);detk=abs(dk);V=(a*det-l);Q=inv(a'*P*a);m0=m0+sqrt((V'*P*V)/(2*n-6));%%%Xs=Xs+dXs;Ys=Ys+dYs;Zs=Zs+dZs;q=q+dq;w=w+dw;k=k+dk;%%%%%%%%%%%%%可以输出迭代次数的i%%%%%%%%%%%%Xs,Ys,Zs,q,w,k,i,dXs,dYs,dZs,dq,dw,dk,detXs,detYs,detZs %%%%%%%%%%%精度mo=m0*sqrt(Q);m=[mo(1,1),mo(2,2),mo(3,3),mo(4,4),mo(5,5),mo(6,6)]';[mXs,mYs,mZs,mq,mw,mk]=testvar(m);%%%%%%%%%输出xy=[Xs,Ys,Zs,q,w,k]';%%输出（6,1）的外方位元素m=[m0,mXs,mYs,mZs,mq,mw,mk]';%%单位误差，各元素中误差R=xyR(xy);%%旋转矩阵函数spaceqianjiao%%空间前交%输入f%输入x1,y1,x2,y2,R1,R2,xy1,xy2 （n,1）%输出X,Y,Z （n,1）function [X,Y,Z]=spaceqianjiao(x1,y1,x2,y2,f,R1,R2,xy1,xy2) i=size(x1,2);[Xs1,Ys1,Zs1,q1,w1,k1]=testvar(xy1);[Xs2,Ys2,Zs2,q2,w2,k2]=testvar(xy2);for n=1:i[X1(n),Y1(n),Z1(n)]=testvar(R1*[x1(n),y1(n),-f]');[X2(n),Y2(n),Z2(n)]=testvar(R2*[x2(n),y2(n),-f]');Bx=Xs2-Xs1;By=Ys2-Ys1;Bz=Zs2-Zs1;N1=(Bx*Z2(n)-Bz*X2(n))/(X1(n)*Z2(n)-X2(n)*Z1(n));N2=(Bx*Z1(n)-Bz*X1(n))/(X1(n)*Z2(n)-X2(n)*Z1(n));X(n)=Xs1+N1*X1(n);Z(n)=Zs1+N1*Z1(n);Y(n)=0.5*((Ys1+N1*Y1(n))+(Ys2+N2*Y2(n)));end函数testvar%分割矩阵。

单片空间后方交会实习报告范文程序设计实习报告班级：学号：姓名：实习任务：用C或VC++语言实现单片后方交汇的计算。

实习目的：1、深入理解单片空间后方交会的原理，2、在有多余观测情况下，用最小二乘平差方法编程实现解求影像外方位元素的过程。

3．通过上机调试程序加强动手能力的培养实习原理以单幅影象为基础，从该影象所覆盖地面范围内若干控制点的已知地面坐标和相应点的像坐标量测值出发，根据共线条件方程，求解该影象在航空摄影时刻的像片外方位元素某，Y，Z,ф,ω,κ.共线条件方程如下：某-某0=-f某[a1(某-某)+b1(Y-Y)+c1(Z-Z)]/[a3(某-某)+b3(Y-Y)+c3(Z-Z)]y-y0=-f某[a2(某-某)+b2(Y-Y)+c2(Z-Z)]/[a3(某-某)+b3(Y-Y)+c3(Z-Z)]某,y为像点的像平面坐标；某0,y0,f为影像的外方位元素；某,Y,Z为摄站点的物方空间坐标；某,Y,Z为物方点的物方空间坐标；.实习主要仪器：电脑（创天中文VC++）程序框图：输入原始数据归算像点坐标某,,y计算和确定初值某0,Y0,Z0,φ0,ω0,κ0组成旋转矩阵R计算（某），（y）和l某，ly逐点组成误差方程并法化迭代次数小于限差否？所有点完否？解法方程，求未知数改正数计算改正后外方位元素未知数改正数尺为1:50000；4个地面控制点的地面坐标及其对应像点的像片坐标：结果输出：已知条件:像点坐标某,y:-53.482.21-14.78-76.6310.4664.43153.24-86.15地点坐标某a,Ya,Za:37631.131324.5728.6939101249352386.540426.530319.8757.31-68.9936589.425273.3内方位元素:某0=y0=0f=153.24mm计算结果:旋转矩阵：0.9977090.06753340.00398399-0.06752540.997715-0.00211178-0.00411750.001837920.99999像点坐标位:(单位：mm)-86.15-68.99-53.4182.21-14.78-76.6310.4764.43单位权中误差：7.2602632e-006外方位元素：某=39795.435精度为:1.1254813Y=27476.479精度为:1.2437625Z=7572.6929精度为:0.48380521q=-0.0039840098精度为:0.00018182003w=0.0021117837精度为:0.00015959235k=-0.067576934精度为:7.2440432e-005迭代次数：4Preanykeytocontinue实习总结：通过这次实习，在最初的程序中遇到了很多地困难，但最终还是克服了它。

空间后方交会—空间前方交会程序编程实验一．实验目的要求掌握运用空间后方交会-空间前方交会求解地面点的空间位置.学会运用空间后方交会的原理，根据所给控制点的地面摄影测量坐标系坐标以及相应的像平面坐标系中的坐标，利用计算机编程语言实现空间后方交会的过程，完成所给像对中两张像片各自的外方位元素的求解。

然后根据空间后方交会所得的两张像片的内外方位元素，利用同名像点在左右像片上的坐标，求解其对应的地面点在摄影测量坐标系中的坐标，并完成精度评定过程，利用计算机编程语言实现此过程.二．仪器用具计算机、编程软件(MATLAB）三．实验数据实验数据包含四个地面控制点（GCP）的地面摄影测量坐标及在左右像片中的像平面坐标。

此四对坐标运用最小二乘法求解左右像片的外方位元素，即完成了空间后方的过程.另外还给出了5对地面点在左右像片中的像平面坐标和左右像片的内方位元素。

实验数据如下：内方位元素：f=152。

000mm，x0=0，y0=0 四．实验框图此过程完成空间后方交会求解像片的外方位元素，其中改正数小于限差（0。

00003，相当于0。

1'的角度值)为止。

在这个过程中采用迭代的方法，是外方位元素逐渐收敛于理论值,每次迭代所得的改正数都应加到上一次的初始值之中。

在空间后方交会中运用的数学模型为共线方程确定Xs，Ys，Zs的初始值时，对于左片可取地面左边两个GCP的坐标的平均值作为左片Xs 和Ys的初始值，取右边两个GCP的坐标平均值作为右片Xs 和Ys的初始值。

Zs可取地面所有GCP的Z坐标的平均值再加上航高.空间前方交会的数学模型为：五．实验源代码function Main_KJQHFJH（）global R g1 g2 m G a c b1 b2;m=10000；a=5;c=4；feval（＠shuru)；％调用shuru(）shurujcp（)函数完成像点及feval（＠shurujcp)；%CCP有关数据的输入XYZ=feval（＠MQZqianfangjh）; %调用MQZqianfangjh（）函数完成空间前方、%%%％%% 单位权中误差％%%％%后方交会计算解得外方位元素global V1 V2；%由于以上三个函数定义在外部文件中故需VV=[］; ％用feval（）完成调用过程for i=1:2＊cVV(i)=V1（i)；VV(2*i+1）=V2(i)；endm0=sqrt（VV＊（VV’)/（2＊c-6））；disp('单位权中误差m0为正负：’）；disp(m0）; ％计算单位权中误差并将其输出显示输入GCP像点坐标及地面摄影测量坐标系坐标的函数和输入所求点像点坐标函数：function shurujcp（)global c m;m=input（’摄影比例尺:')；%输入GCP像点坐标数据函数并分别将其c=input（'GCP的总数=')；％存入到不同的矩阵之中disp（'GCP左片像框标坐标：')；global g1;g1=zeros（c，2）;i=1;while i<=cm=input('x='）；n=input（'y='）;g1(i，1)=m;g1（i，2)=n；i=i+1；enddisp('GCP右片像框标坐标：’);global g2；g2=zeros(c,2);i=1;while i〈=cm=input（'x=’）;n=input('y=’）;g2（i，1)=m；g2（i，2）=n;i=i+1;end％%%%％%%%%%%%%％％%%％％%％％％％%％%%％％%%％％%％％%%％% function shuru（)global a；a=input（'计算总像对点数='); %完成想计算所需的像平面坐标global b1；％坐标输入,存入不同的矩阵中b1=zeros(a,2)；disp('左片像点坐标：'）i=1；while i〈=am=input（'x=’）；n=input(’y=’);b1(i，1)=m；b1（i,2)=n;i=i+1;end%%global b2;b2=zeros（a,2）；disp（’右片像点坐标：')i=1;while i〈=am=input('x=’）;n=input（'y=’）；b2(i,1）=m；b2(i,2)=n；i=i+1;end％%global c;c=input(’GCP的总数='）；disp('GCP摄影测量系坐标：’)global G;G=zeros(3,c)；i=1；while i〈=cm=input（’X='）;n=input(’Y='）;v=input（’Z=')；G(i，1）=m;G（i，2)=n;G（i,3)=v；i=i+1;end %%％％％％%%％%%%%%%％％％％％％％％%%％％%％%%%％％%%%空间前方交会和后方交会函数:function XYZ=MQZqianfangjh（）global R1 R2 a f b1 b2 Ra Rb；global X1 X2;R1=Ra；R2=Rb;R1=zeros（3,3）;R2=zeros（3,3)；global g1 g2 G V1 V2 V WF c QXX QXX1 QXX2；xs0=（G（1，1)+G(3,1))/2；ys0=（G（1，2）+G(3，2)）/2;［Xs1，Ys1，Zs1,q1,w1，k1 R]=houfangjh(g1，xs0,ys0）；％对左片调用后方交会函数R1=R；V1=V；WF1=WF;QXX1=QXX;save '左片外方位元素为。

任务已知f=153.24mm，m=10000，限差0.1’各点坐标点号像点坐标地面坐标x（mm）y（mm）X（m）Y（m）Z（m）1 -86.15 -68.99 36589.41 25273.32 2195.172 -53.40 82.21 37631.08 31324.51 728.693 -14.78 -76.63 39100.97 24934.98 2386.504 10.46 64.43 40426.54 30319.81 757.31求近似垂直摄影情况下后方交会解设计任务1、确定未知数的初始值：Φ0 =ω0 =К0 = 0 , 内方位元素，，f=153.24mm。

； = 38437m ；= 27963.16m2、计算旋转矩阵R利用角元素的近似值计算方向余弦值，组成R阵根据《摄影测量学》P32中的公式（3-9），初步计算R阵R[0][0]=cos(Φ)*cos(K)-sin(Φ)*sin(W)*sin(K);R[0][1]=-cos(Φ)*sin(K)-sin(Φ)*sin(W)*cos(K);R[0][2]=-sin(Φ)*cos(W);R[1][0]=cos(W)*sin(K);R[1][1]=cos(W)*cos(K);R[1][2]=-sin(W);R[2][0]=sin(Φ)*cos(K)+cos(Φ)*sin(W)*sin(K);R[2][1]=-sin(Φ)*sin(K)+cos(Φ)*sin(W)*cos(K);R[2][2]=cos(Φ)*cos(W);得初始R阵3、逐点计算近似值（x），（y）：带入《摄影测量学》P61的公式（5-1）；得4、组成误差方程式：按（5-8）；（5-9b）、（5-4）式逐点计算误差方程式的系数和常数项根据Lx=x-(x)；Ly=y-(y)得解得A阵为根据公式解除X阵，检验Φ1 、ω1 、К1 并不在限差0.1’之内;利用解得的Φ1 、ω1 、К1，带入2中循环；知道解出合格结果！5、结果，通过C++程序解得结果6、源程序代码#include <iostream>#include <fstream>#include <cmath>#include <string>#include <iomanip>using namespace std;const int n=6;void inverse (double c[n][n]);template<typename T1,typename T2>void transpose (T1*mat1,T2*mat2,int a,int b);template<typename T1,typename T2>void multi(T1*mat1,T2 * mat2,T2 * result,int a,int b,int c);template<typename T>void input (T*mat,int a,int b);template<typename T>void output(T*mat,char*s,int a,int b);int main(){ofstream outFile;cout.precision(5);double x0=0.0, y0=0.0; double fk=0.15324; //内方位元素double m=10000; //估算比例尺double B[4][5]={0.0},R[3][3],XG[6][1],AT[6][8],ATA[6][6],ATL[6][1];input (B,4,5); //从文件中读取控制点的影像坐标和地面坐标，存入数组Bdouble Xs=0.0, Ys=0.0, Zs=0.0,Q=0.0,W=0.0,K=0.0;double X,Y,Z,L[8][1],A[8][6];//确定未知数的出始值for(int i=0;i<4;i++){Xs=Xs+B[i][2];Ys=Ys+B[i][3];Zs=Zs+B[i][4];}Xs=Xs/4; Ys=Ys/4; Zs=Zs/4+m*fk;//cout<<Xs<<" "<<Ys<<" "<<Zs<<endl;//int f=0;do//迭代计算{f++;//组成旋转矩阵R[0][0]=cos(Q)*cos(K)-sin(Q)*sin(W)*sin(K);R[0][1]=-cos(Q)*sin(K)-sin(Q)*sin(W)*cos(K);R[0][2]=-sin(Q)*cos(W);R[1][0]=cos(W)*sin(K);R[1][1]=cos(W)*cos(K);R[1][2]=-sin(W);R[2][0]=sin(Q)*cos(K)+cos(Q)*sin(W)*sin(K);R[2][1]=-sin(Q)*sin(K)+cos(Q)*sin(W)*cos(K);R[2][2]=cos(Q)*cos(W);for(int t=0;t<3;t++){for(int j=0;j<3;j++)cout<<" "<<R[t][j];cout<<endl;}cout<<endl;//计算系数阵和常数项for(int i=0,k=0,j=0;i<=3;i++,k++,j++){X=R[0][0]*(B[i][2]-Xs)+R[1][0]*(B[i][3]-Ys)+R[2][0]*(B[i][4]-Zs);Y=R[0][1]*(B[i][2]-Xs)+R[1][1]*(B[i][3]-Ys)+R[2][1]*(B[i][4]-Zs);Z=R[0][2]*(B[i][2]-Xs)+R[1][2]*(B[i][3]-Ys)+R[2][2]*(B[i][4]-Zs);//L[j][0]=B[i][0]-(x0-fk*X/Z); //l=x-(x)//插入cout<<L[j][0]<<" ";//L[j+1][0]=B[i][1]-(y0-fk*Y/Z); //l=y-(y)//插入cout<<L[j+1][0]<<endl;//j++;A[k][0]=(R[0][0]*fk+R[0][2]*(B[i][0]-x0))/Z;A[k][1]=(R[1][0]*fk+R[1][2]*(B[i][0]-x0))/Z;A[k][2]=(R[2][0]*fk+R[2][2]*(B[i][0]-x0))/Z;A[k][3]=(B[i][1]-y0)*sin(W)-((B[i][0]-x0)*((B[i][0]-x0)*cos(K)-(B[i][1]-y0) *sin(K))/fk+fk*cos(K))*cos(W);A[k][4]=-fk*sin(K)-(B[i][0]-x0)*((B[i][0]-x0)*sin(K)+(B[i][1]-y0)*cos(K))/f k;A[k][5]=B[i][1]-y0;A[k+1][0]=(R[0][1]*fk+R[0][2]*(B[i][1]-y0))/Z;A[k+1][1]=(R[1][1]*fk+R[1][2]*(B[i][1]-y0))/Z;A[k+1][2]=(R[2][1]*fk+R[2][2]*(B[i][1]-y0))/Z;A[k+1][3]=-(B[i][0]-x0)*sin(W)-((B[i][1]-y0)*((B[i][0]-x0)*cos(K)-(B[i][1]-y0)*sin(K))/fk-fk*sin(K))*cos(W);A[k+1][4]=-fk*cos(K)-(B[i][1]-y0)*((B[i][0]-x0)*sin(K)+(B[i][1]-y0)*cos(K)) /fk;A[k+1][5]=-(B[i][0]-x0);k++;}transpose(A,AT,6,8);multi(AT,A,ATA,6,8,6);inverse(ATA);multi(AT,L,ATL,6,8,1);multi(ATA,ATL,XG,6,6,1);Xs=Xs+XG[0][0]; Ys=Ys+XG[1][0]; Zs=Zs+XG[2][0];Q=Q+XG[3][0]; W=W+XG[4][0]; K=K+XG[5][0];//cout<<Q<<"**"<<W<<"**"<<K<<endl;//}while(XG[3][0]>=6.0/206265.0||XG[4][0]>=6.0/206265.0||XG[5][0]>=6.0/206265 .0);cout<<"迭代次数为："<<f<<endl;//精度评定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][0];transpose (V,VT,1,8);multi(VT,V,VTV,1,8,1);m0=VTV[0][0]/2;for(i=0;i<6;i++)for(int j=0;j<6;j++)D[i][j]=m0*ATA[i][j];//屏幕输出误差方程系数阵、常数项、改正数output(A,"误差方程系数阵A为：",8,6);output(L,"常数项L为：",8,1);output(XG,"改正数为：",6,1);out("aim.txt",ios::app); //打开并添加aim.txt文件out(10);//以文件的形式输出像片外方位元素、旋转矩阵、方差阵outFile<<"一、像片的外方位元素为："<<endl<<endl;outFile<<setw(10)<<"Xs="<<Xs<<setw(10)<<"Ys="<<Ys<<setw(10)<<"Zs="<<Zs<<endl;outFile<<setw(20)<<"航向倾角为："<<Q<<setw(10)<<"旁向倾角为："<<W<<setw(10)<<"像片旋角为："<<K<<endl;outFile<<endl<<"二、旋转矩阵R为:"<<endl<<endl;for( i=0;i<3;i++){for(int j=0;j<3;j++)outFile<<setw(25)<<R[i][j]<<setw(25);outFile<<endl;}outFile<<endl;outFile<<setw(0)<<"三、精度评定结果为："<<endl;out(5);for(i=0;i<6;i++){for(int j=0;j<6;j++)outFile<<setw(14)<<D[i][j]<<setw(14);outFile<<endl;}out();return 0;}template<typename T1,typename T2>void transpose(T1*mat1,T2*mat2,int a,int b){ int i,j;for(i=0;i<b;i++)for(j=0;j<a;j++)mat2[j][i]=mat1[i][j];return;}template<typename T1,typename T2>void multi(T1*mat1,T2 * mat2,T2 * result,int a,int b,int c){ int i,j,k;for(i=0;i<a;i++){for(j=0;j<c;j++){result[i][j]=0;for(k=0;k<b;k++)result[i][j]+=mat1[i][k]*mat2[k][j];}}return;}template <typename T>void input (T*mat,int a,int b){ ifstream inFile;in("控制点坐标.txt");while(!in()){for (int i=0;i<a;i++)for(int j=0;j<b;j++)inFile>>mat[i][j];}in();return;}template<typename T>void output(T*mat,char*s,int a,int b){ cout<<setw(15)<<s<<endl;for(int i=0;i<a;i++){for(int j=0;j<b;j++)cout<<setw(13)<<mat[i][j];cout<<endl;}return;}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;elseq[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;j<12;j++){ q[i][j]*=p;q[i][j]-=q[k][j];}}for(h=k=n-1;k>0;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++)//将对角线上数据化为1{ 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];}7、感想这次摄影测量作业感觉非常难！一开始审题就出现了很大问题，读不懂，不知道要求什么东西。

摄影测量学单像空间后方交会编程实习报告本次实习中，我使用编程语言进行了单像空间后方交会的实现，并取得了一定的成果。

首先，我了解了单像空间后方交会的基本原理。

根据光线在透镜上的成像规律，可以推导出物体点在像平面上的坐标与图像点在像平面上的坐标之间的关系式。

通过已知的摄像机内外方位元素和图像点坐标，可以反求得物体点的坐标。

在程序编写过程中，我采用了Python编程语言。

首先，我定义了一个类，用于存储摄像机的内外方位元素和图像点坐标。

然后，我编写了一个函数，用于计算物体点的坐标。

该函数根据已知的内外方位元素和图像点坐标，使用逆向投影的方式反求物体点的坐标。

最后，我编写了一个主函数，通过读取输入文件中的数据，调用计算函数，并将结果保存到输出文件中。

在实现的过程中，我遇到了一些问题。

首先，由于摄像机的内外方位元素需要提前获取，因此我通过测量方法获得了实际的内外方位元素。

然而，测量的过程中存在一定的误差，因此在计算物体点坐标时可能存在一定的误差。

其次，图像坐标与物体点坐标之间的关系式中存在一些参数，如焦距、主点坐标等，这些参数也需要提前获取。

在程序中，我将这些参数作为输入参数，通过外部文件进行输入。

在实习的过程中，我充分运用了自己所学的摄影测量学知识，并将其与编程技能相结合。

在实现过程中，我遇到了一些难题，但通过查阅资料和与老师的讨论，最终得以解决。

通过编程实习，我深入理解了单像空间后方交会的原理，并通过实际操作提高了自己的计算能力。

总的来说，本次实习使我对摄影测量学有了更深入的认识，也提升了我的计算和编程能力。

通过此次实习，我对摄影测量学的兴趣更加浓厚，也更加期待在今后的学习和研究中能够进一步深入探索。

一、实习目的本次实习的主要目的是通过单张影像的空间后方交会算法，掌握摄影测量学中的基本原理和方法，培养实际操作能力和理论联系实际的能力。

通过实习，要求学生了解空间后方交会的概念、原理和计算方法，熟练运用相关软件进行空间后方交会的实际操作，并对结果进行分析和评定。

二、实习内容1. 了解空间后方交会的基本概念：空间后方交会是一种摄影测量方法，通过分析单张影像上若干个控制点的地面坐标和像坐标，利用共线条件方程，求解影像的外方位元素。

2. 学习空间后方交会的原理：影像上的控制点与相机镜头中心形成共线方程，通过最小二乘法求解影像的外方位元素。

外方位元素包括六个参数：三个表示相机在三维空间中的位置，三个表示相机在三维空间中的姿态。

3. 掌握空间后方交会的计算方法：利用共线条件方程，通过最小二乘法求解外方位元素。

计算过程中需要注意控制点的选择、权重的分配和迭代求解的精度。

4. 熟悉空间后方交会的实际操作：使用相关软件（如MATLAB、ERDAS IMAGINE等）进行空间后方交会的实际操作，包括数据的输入、参数的设置和结果的输出。

5. 结果分析和评定：对空间后方交会的结果进行分析和评定，包括计算精度、误差分析和结果的可信度。

三、实习过程1. 理论学习和讨论：在学习空间后方交会的基本概念和原理的基础上，进行小组讨论，理解并掌握空间后方交会的计算方法。

2. 软件操作练习：在老师的指导下，使用相关软件进行空间后方交会的实际操作，熟悉软件的操作界面和功能。

3. 数据处理和计算：根据实习提供的数据，进行空间后方交会的计算，包括控制点的选择、权重的分配和迭代求解的过程。

4. 结果分析和评定：对空间后方交会的结果进行分析和评定，计算精度、误差分析和结果的可信度。

四、实习总结通过本次实习，我对空间后方交会的基本概念、原理和计算方法有了更深入的了解，掌握了相关软件的操作方法，并能够进行空间后方交会的实际操作。

在实习过程中，我学会了如何选择控制点、分配权重和判断迭代求解的精度。

摄影测量学实验报告实验名称：空间后方交会班级：测绘工程12-1一．实验目的掌握运用空间后方交会求解外方元素的过程。

学会运用空间后方交会的原理，根据所给控制点的地面摄影测量坐标系坐标以及相应的像平面坐标系中的坐标，利用已知公式和计算机编程语言实现空间后方交会的过程，完成外方位元素的求解。

二．实验工具个人电脑，visual basic 6.0三．实验源代码Private Sub Command1_Click()Dim c(), kx(), ky(), l(), aa(), a11(), a12(), a13(), a14(), a15(), a16(), a21(), a22(), a23(), a24(), a25()Dim a26(), v(), mm(5)nn = Val(Text2.Text)nnt = nn - 1ReDim c(nnt, 5), kx(nnt), ky(nnt), l(2 * nnt + 1, 0), aa(2 * nnt + 1, 5), a11(nnt), a12(nnt), a13(nnt), a14(nnt), a15(nnt), a16(nnt), a21(nnt), a22(nnt), a23(nnt), a24(nnt), a25(nnt)ReDim a26(nnt)f = Val(Text1.Text)fei = 0#w = 0#k = 0#If Label2.Caption = "" ThenMsgBox "尚未读取文本", , "提示信息"Elsea = Split(Label2.Caption, vbCrLf)For i = 0 To nnthh = Split(a(i), " ")For j = 0 To 5c(i, j) = Val(hh(j))Next jNext im1 = (c(0, 1) - c(1, 1)) ^ 2 + (c(0, 2) - c(1, 2)) ^ 2 m2 = (c(0, 3) - c(1, 3)) ^ 2 + (c(0, 4) - c(1, 4)) ^ 2 m = Sqr(m2 / m1)zs = m * fxs = 0ys = 0For i = 0 To nntxs = xs + c(i, 3)ys = xs + c(i, 4)Next ixs = xs / nnys = ys / nnDoa1 = Cos(fei) * Cos(k) - Sin(fei) * Sin(w) * Sin(k)a2 = -Cos(fei) * Sin(k) - Sin(fei) * Sin(w) * Cos(k)a3 = -Sin(fei) * Cos(w)b1 = Cos(w) * Sin(k)b2 = Cos(w) * Cos(k)b3 = -Sin(w)c1 = Sin(fei) * Cos(k) + Cos(fei) * Sin(w) * Sin(k)c2 = -Sin(fei) * Sin(k) + Cos(fei) * Sin(w) * Cos(k)c3 = Cos(fei) * Cos(w)For i = 0 To nntkx(i) = -f * (a1 * (c(i, 3) - xs) + b1 * (c(i, 4) - ys) + c1 * (c(i, 5) - zs)) / (a3 * (c(i, 3) - xs) + b3 * (c(i, 4) - ys) + c3 * (c(i, 5) - zs))ky(i) = -f * (a2 * (c(i, 3) - xs) + b2 * (c(i, 4) - ys) + c2 * (c(i, 5) - zs)) / (a3 * (c(i, 3) - xs) + b3 * (c(i, 4) - ys) + c3 * (c(i, 5) - zs))Next iFor i = 0 To nntl(i * 2, 0) = c(i, 1) - kx(i)l(i * 2 + 1, 0) = c(i, 2) - ky(i)Next iFor i = 0 To nntzb = (a3 * (c(i, 3) - xs) + b3 * (c(i, 4) - ys) + c3 * (c(i, 5) - zs))a11(i) = (a1 * f + a3 * c(i, 1)) / zba12(i) = (b1 * f + b3 * c(i, 1)) / zba13(i) = (c1 * f + c3 * c(i, 1)) / zba21(i) = (a2 * f + a3 * c(i, 2)) / zba22(i) = (b2 * f + b3 * c(i, 2)) / zba23(i) = (c2 * f + c3 * c(i, 2)) / zba14(i) = c(i, 2) * Sin(w) - (c(i, 1) * (c(i, 1) * Cos(k) - c(i, 2) * Sin(k)) / f + f * Cos(k)) * Cos(w)a15(i) = -f * Sin(k) - c(i, 1) / f * (c(i, 1) * Sin(k) + c(i, 2) * Cos(k))a16(i) = c(i, 2)a24(i) = -c(i, 1) * Sin(w) - (c(i, 2) * (c(i, 1) * Cos(k) - c(i, 2) * Sin(k)) / f - f * Sin(k)) * Cos(w)a25(i) = -f * Cos(k) - c(i, 2) / f * (c(i, 1) * Sin(k) + c(i, 2) * Cos(k))a26(i) = -c(i, 1)Next iFor i = 0 To nntaa(2 * i, 0) = a11(i)aa(2 * i, 1) = a12(i)aa(2 * i, 2) = a13(i)aa(2 * i, 3) = a14(i)aa(2 * i, 4) = a15(i)aa(2 * i, 5) = a16(i)aa(2 * i + 1, 0) = a21(i)aa(2 * i + 1, 1) = a22(i)aa(2 * i + 1, 2) = a23(i)aa(2 * i + 1, 3) = a24(i)aa(2 * i + 1, 4) = a25(i)aa(2 * i + 1, 5) = a26(i)Next iReDim at(5, 2 * nnt + 1), ata(5, 5), nata(5, 5), nataat(5, 2 * nnt + 1), detx(5, 0)at = zhuanzhi(aa)ata = cheng(at, aa)nata = qiuni(ata)nataat = cheng(nata, at)detx = cheng(nataat, l)xs = xs + detx(0, 0)ys = ys + detx(1, 0)zs = zs + detx(2, 0)fei = fei + detx(3, 0)w = w + detx(4, 0)k = k + detx(5, 0)Loop Until Abs(detx(0, 0)) < 0.000001 And Abs(detx(1, 0)) < 0.000001 And Abs(detx(2, 0)) < 0.000001 And Abs(detx(3, 0)) <0.000001 And Abs(detx(4, 0)) < 0.000001 And Abs(detx(5, 0)) < 0.000001ReDim v(2 * nnt + 1, 0), ax(2 * nnt + 1, 0)ax = cheng(aa, detx)v = jian(ax, l)vv = 0For i = 0 To 2 * nnt + 1vv = v(i, 0) ^ 2 + vvNext imo = Sqr(vv / (2 * nn - 6))For i = 0 To 5mm(i) = mo * Sqr(nata(i, i))Next itxtshow.Text = txtshow.Text & "外方元素为中误差为" & vbCrLf & "Xs " & xs & " " & mm(0) & vbCrLf & "Ys " & ys & " " & mm(1) & vbCrLf & "Zs " & zs & " " & mm(2) & vbCrLf & "Φ" & fei & " " & mm(3) & vbCrLf & "W " & w & " " & mm(4) & vbCrLf & "K " & k & " " & mm(5) & vbCrLf & vbCrLftxtshow.Text = txtshow.Text & " 报告时间：" & Year(Date) & Format(Month(Date), "00") & Format(Day(Date), "00") & vbCrLftxtshow.Text = txtshow.Text & " 审批人签字："End IfEnd SubPrivate Sub Command2_Click()CommonDialog1.Filter = "文档|*.txt"CommonDialog1.ShowOpenLabel1.Caption = CommonDialog1.FileNameEnd SubPrivate Sub Command3_Click()Label2.Caption = ""Dim mylineIf Label1.Caption <> "" ThenOpen Label1.Caption For Input As #1Do While Not EOF(1)Line Input #1, mylineLabel2.Caption = Label2.Caption + myline + vbCrLfLoopClose #1End SubPrivate Sub Command4_Click()CommonDialog1.DialogTitle = "保存文件" CommonDialog1.Filter = "文档|*.doc" CommonDialog1.Action = 2If CommonDialog1.FileName <> "" ThenOpen CommonDialog1.FileName For Output As #1Print #1, txtshow.TextClose #1End IfEnd SubPrivate Sub Command5_Click()txtshow.Text = ""End SubPrivate Sub Command6_Click()If MsgBox("确认退出程序？", 68, "提示信息") = 6 Then EndEnd IfEnd SubPublic Function cheng(a(), b())ReDim s(UBound(a, 1), UBound(b, 2))If UBound(a, 2) <> UBound(b, 1) Then MsgBox "不符合乘法要求"EndElseFor i = 0 To UBound(s, 1)For j = 0 To UBound(s, 2)s(i, j) = 0For Q = 0 To UBound(a, 2)s(i, j) = Val(s(i, j)) + Val(a(i, Q)) * Val(b(Q, j)) Next QNext jNext icheng = sEnd IfEnd FunctionPublic Function zhuanzhi(a())Dim sReDim s(UBound(a, 2), UBound(a, 1))For i = 0 To UBound(a, 1)For j = 0 To UBound(a, 2)s(j, i) = a(i, j)Next jNext izhuanzhi = sEnd FunctionPublic Function qiuni(m())Dim s()If UBound(m, 1) <> UBound(m, 2) Then MsgBox ("求逆错误")EndElsen = UBound(m, 1)Dim e(), m1, m2, m3ReDim e(n, 2 * n + 1)ReDim s(n, 2 * n + 1)For i = 0 To nFor j = 0 To ne(i, j) = m(i, j)Next jNext iFor i = 0 To nFor j = n + 1 To 2 * n + 1e(i, j) = 0Next jNext iFor i = 0 To ne(i, n + i + 1) = 1 Next iFor t = 0 To nIf e(t, t) = 0 ThenFor i = t + 1 To nFor j = 0 To 2 * n + 1 s(i, j) = e(i, j)e(i, j) = e(t, j)e(t, j) = s(t, j)Next jIf e(t, t) <> 0 Then GoTo daima1End IfNext iIf e(t, t) = 0 Then MsgBox ("求逆错误") EndGoTo lastlineEnd Ifdaima1:m1 = e(t, t)For j = 0 To 2 * n + 1e(t, j) = e(t, j) / m1Next jFor i = t + 1 To nm2 = e(i, t)For j = 0 To 2 * n + 1e(i, j) = e(i, j) - m2 * e(t, j) Next jNext iNext tFor t = 0 To n - 1For i = t + 1 To nm3 = e(t, i)For j = 0 To 2 * n + 1e(t, j) = e(t, j) - m3 * e(i, j) Next jNext iNext tReDim c(n, n)For i = 0 To nFor j = 0 To nc(i, j) = Round(e(i, j + n + 1), 2)Next jNext iqiuni = cEnd Iflastline:End FunctionPublic Function jian(a(), b())Dim s, tReDim s(UBound(a, 1), UBound(a, 2))If UBound(a, 1) <> UBound(b, 1) Or UBound(a, 2) <> UBound(b, 2) ThenMsgBox "不符合减法要求"EndElseFor i = 0 To UBound(a, 1)For j = 0 To UBound(a, 2)s(i, j) = Val(a(i, j)) - Val(b(i, j))Next jNext ijian = sEnd IfEnd FunctionPrivate Sub Form_Load()Form1.Image1.Top = 0Form1.Image1.Left = 0Form1.Image1.Width = Form1.Width Form1.Image1.Height = Form1.Height End SubPrivate Sub Form_Resize()Form1.Image1.Top = 0Form1.Image1.Left = 0Form1.Image1.Width = Form1.Width Form1.Image1.Height = Form1.Height End Sub四．实验框图获取已知数据量测控制点的像点坐标确定未知数初始值五．实验数据1 -0.08615 -0.06899 36589.41 25273.32 2195.172 -0.05340 0.08221 37631.08 31324.51 728.693 -0.01478 -0.07663 39100.97 24934.98 2386.54 0.01046 0.06443 40426.54 30319.81 757.31 六．实验截图七．实验心得此次实验让我更加了解空间后方交会的计算过程，加强了对空间后方交会基本公式和误差方程式，法线方程式的记忆。

立体像对前方交会实验报告一、实验目的在掌握前方交会原理的基础上，自己编写前方交会程序，在计算机上调试，输出计算结果并对计算结果进行检验。

通过上机调试程序加强动手能力的培养，通过对实验过程的掌握以及对实验结果的分析，增强学生综合运用所学知识解决实际问题的能力。

二、实验仪器计算机，VB6.0三、实验数据1.模拟像片一对：左片号23，右片号24；2.航摄机主距：f=150mm；3.左片23号片外方位元素：φ=−0°25′00″ω=−1°00′00″k=−0°10′00″Xs=103007.006117 Ys=139998.994849 Zs=4801.9989994 （m）右片24号片外方为元素：φ=1°39′59″ω=−0°10′00″k=0°40′00″Xs=106002.023762 Ys=140005.002780 Zs=4797.009648 （m）待求像点坐标如下表：四、实验内容利用所给立体像对两张像片的内、外方位元素，编写空间前方交会程序，根据所给像对中若干同名像点在左右像片上的坐标，解求其对应的地面点的物方坐标，实现空间前方交会的过程。

五、实验成果程序流程图：程序设计界面：程序运行界面：运行结果：（注：表上显示地面点坐标依次是：7,9,4,6,5,8）附：excel进行角度转换：六、程序如下：Dim f#, x1#, y1#, x2#, y2#, i%, j%, u1#, u2#, v1#, v2#, w1#, w2#, fai1#, kab1#, omg1#, fai2#, kab2#, omg2# Dim a12#, a13#, b11#, b12#, b13#, c11#, c12#, c13#, a21#, a22#, a23#, b21#, b22#, b23#, c21#, c22#, c23# Dim n1#, n2#, bu#, bv#, bw#Dim xs1#, xs2#, ys1#, ys2#, zs1#, zs2#Private Sub Command1_Click()fai1 = Val(Text1.Text): kab1 = Val(Text3.Text): omg1 = Val(Text2.Text)fai2 = Val(Text11.Text): kab2 = Val(Text10.Text): omg2 = Val(Text7.Text)f = Val(Text13.Text)xs1 = Val(Text4.Text): xs2 = Val(Text9.Text): ys1 = Val(Text5.Text): ys2 = Val(Text8.Text): zs1 = Val(Text6.Text): zs2 = Val(Text12.Text)End SubPrivate Sub Command3_Click()x1 = Val(InputBox("输入23片坐标x的值"))y1 = Val(InputBox("输入23片坐标y的值"))x2 = Val(InputBox("输入24片坐标x的值"))y2 = Val(InputBox("输入24片坐标y的值"))End SubPrivate Sub Command2_Click()Text1.Visible = FalseText2.Visible = FalseText3.Visible = FalseText4.Visible = FalseText5.Visible = FalseText6.Visible = FalseText7.Visible = FalseText8.Visible = FalseText9.Visible = FalseText10.Visible = FalseText11.Visible = FalseText12.Visible = FalseText13.Visible = FalseLabel1.Visible = FalseLabel2.Visible = FalseLabel3.Visible = FalseLabel4.Visible = FalseLabel5.Visible = FalseLabel6.Visible = FalseLabel7.Visible = FalseLabel8.Visible = FalseLabel9.Visible = FalseLabel10.Visible = FalseLabel11.Visible = FalseLabel12.Visible = FalseLabel13.Visible = FalseLabel14.Visible = Falsea11 = Cos(fai1) * Cos(kab1) - Sin(fai1) * Sin(omg1) * Sin(kab1)a12 = -Cos(fai1) * Sin(kab1) - Sin(fai1) * Sin(omg1) * Cos(kab1)a13 = -Sin(fai1) * Cos(omg1)b11 = Cos(omg1) * Sin(kab1)b12 = Cos(omg1) * Cos(kab1)b13 = -Sin(omg1)c11 = Sin(fai1) * Cos(kab1) + Cos(fai1) * Sin(omg1) * Sin(kab1)c12 = -Sin(fai1) * Sin(kab1) + Cos(fai1) * Sin(omg1) * Cos(kab1)c13 = Cos(fai1) * Cos(omg1)a21 = Cos(fai2) * Cos(kab2) - Sin(fai2) * Sin(omg2) * Sin(kab2)a22 = -Cos(fai2) * Sin(kab2) - Sin(fai2) * Sin(omg2) * Cos(kab2)a23 = -Sin(fai2) * Cos(omg2)b21 = Cos(omg2) * Sin(kab2)b22 = Cos(omg2) * Cos(kab2)b23 = -Sin(omg2)c21 = Sin(fai2) * Cos(kab2) + Cos(fai2) * Sin(omg2) * Sin(kab2)c22 = -Sin(fai2) * Sin(kab2) + Cos(fai2) * Sin(omg2) * Cos(kab2)c23 = Cos(fai2) * Cos(omg2)u1 = a11 * x1 + a12 * y1 - a13 * fv1 = b11 * x1 + b12 * y1 - b13 * fw1 = c11 * x1 + c12 * y1 - c13 * fu2 = a21 * x2 + a22 * y2 - a23 * fv2 = b21 * x2 + b22 * y2 - b23 * fw2 = c21 * x2 + c22 * y2 - c23 * fbu = xs2 - xs1bv = ys2 - ys1bw = zs2 - zs1n1 = (bu * w2 - bw * u2) / (u1 * w2 - u2 * w1)n2 = (bu * w1 - bw * u1) / (u1 * w2 - u2 * w1)Dim xx#, yy#, yy1#, yy2#, zz#xx = Fix((xs1 + u1 * n1) * 1000 + 0.5) / 1000yy1 = Val(Text5.Text) + v1 * n1yy2 = Val(Text8.Text) + v2 * n2yy = Fix(((yy1 + yy2) / 2) * 1000 + 0.5) / 1000zz = Fix((Val(Text6.Text) + w1 * n1) * 1000 + 0.5) / 1000Form1.Print "地面坐标为:"Form1.Print "xx="; xx, "yy="; yy; "zz="; zzEnd Sub七、实验心得通过本次前方交会实验，熟悉掌握了前方交会原理，并利用计算机程序对其进行了解决，收益颇多！。

中南大学本科生课程设计（实践）任务书、设计报告（摄影测量与遥感概论）题目：空间后方交会-前交院系：地球科学与信息物理学院班级：测绘1201班********学号：***********名：***二零一四年十一月一、实验目的通过对数字影像空间后交前交的程序设计实验，要求我们进一步理解和掌握影像外方位元素的有关理论、原理和方法。

利用计算机程序设计语言编写摄影测量空间交会软件进行快速确定影响的外方位元素及其精度，然后通过求得的外方位元素求解位置点的地面摄影测量坐标，达到通过摄影测量量测地面地理坐标的目的。

二、实验要求1.用C、VB、C++或MA TLAB语言编写空间后方交会-空间前方交会程序2.提交实习报告：程序框图、程序源代码、计算结果、体会3.计算结果：像点坐标、地面坐标、单位权中误差、外方位元素及其精度4.完成时间：11月11日前完成三、实验数据四、实验思路➢利用后方交会得出两张相片各自的外方位元素1）获取已知数据影响比例尺m,，内方位元素x0 、y0 、f ，控制点的地面摄影测量坐标Xtp, Ytp, Ztp2）量测控制点左片和右片的像点坐标 x，y3）确定未知数初值 Xs0, Ys0, Zs0, ω，φ，κ（线元素可用控制点均值代替，角元素可用0初始化），即：∑=Xtp X 41s0，∑=Ytp Y 41s0，f Z *m s =ω=0，φ=0，κ=0 4）计算旋转矩阵R5）利用共线方程逐点计算像点坐标的近似值 6）组成误差方程式并法化 7）解求外方位元素改正数8）检查迭代是否收敛（改正值是否小于某一特定常数） ➢ 利用解出的外方位元素进行前方交会1）获取已知数据x0 , y0 , f, XS1, YS1, ZS1,φ1,ω1,κ1 , XS2, YS2, ZS2，φ2，ω2，κ22）量测像点坐标 x1,y1 ,x2,y23）由外方位线元素计算基线分量BX, BY, BZ4）由外方位角元素计算像空间辅助坐标 X1, Y1, Z1 , X2, Y2, Z2 5）计算点投影系数 N1 , N2 6）计算地面坐标 XA, YA, ZA五、实验过程➢ 程序流程图此过程完成空间后方交会求解像片的外方位元素，其中改正数小于限差（0.00003，相当于0.1的角度值）为止。