作业4--空间后方交会
空间后方交会的解算
空间后方交会的解算一． 空间后方交会的目的摄影测量主要利用摄影的方法获取地面的信息，主要是是点位信息，属性信息，因此要对此进行空间定位和建模，并首先确定模型的参数，这就是空间后方交会的目的，用以求出模型外方位元素。
二． 空间后方交会的原理空间后方交会的原理是共线方程。
共线方程是依据相似三角形原理给出的，其形式如下111333222333()()()()()()()()()()()()A S A S A S A S A S A S AS A S A S A S A S A S a X X b Y Y c Z Z x f a X X a Y Y a Z Z a X X b Y Y c Z Z y f a X X a Y Y a Z Z -+-+-=--+-+--+-+-=--+-+-上式成为中心投影的构线方程，我们可以根据几个已知点，来计算方程的参数，一般需要六个方程，或者要三个点，为提高精度，可存在多余观测，然后利用最小二乘求其最小二乘解。
将公式利用泰勒公式线性化，取至一次项，得到其系数矩阵A ；引入改正数（残差）V ，则可将其写成矩阵形式：V AX L =-其中111333222333[,]()()()()()()()()()()()()()()Tx y A S A S A S x A S A S A S A S A S A S y A S A S A S L l l a X X b Y Y c Z Z l x x x fa X X a Y Y a Z Z a X Xb Y Yc Z Z l y y y fa X X a Y Y a Z Z =-+-+-=-=+-+-+--+-+-=-=+-+-+- 则1()T T X A A A L -=X 为外方位元素的近似改正数，由于采用泰勒展开取至一次项，为减少误差，要将的出的值作为近似值进行迭代，知道小于规定的误差三． 空间后方交会解算过程1. 已知条件近似垂直摄影00253.24mmx y 0f ===2. 解算程序流程图MATLAB 程序format long;s1=xlsread('data.xls');%读取数据a1=s1(1:4,1:2);%影像坐标b1=s1(1:4,3:5);%地面摄影测量坐标a2=s1.*10^-3;%影像坐标单位转化j1=a2(1,:)-a2(2,:);j2=j1(1,1)^2+j1(1,2)^2;lengh_a1=sqrt(j2); %相片某一长度j1=b1(1,:)-b1(1,:);j2=j1(1,1)^2+j1(1,2)^2;lengh_b1=sqrt(j2); %地面对应的长度m=lengh_b1/lengh_a1;%求出比例尺n0=0;p0=0;q0=0;x0=mean(b1(:,1));y0=mean(b1(:,2));f=153.24*10^-3;z0=m*f;x001={x0,x0,x0,x0};X0=cell2mat(x001)';y001={y0,y0,y0,y0};Y0=cell2mat(y001)';z001={z0,z0,z0,z0};Z0=cell2mat(z001)';%初始化外方位元素的值aa1=cos(n0)*cos(q0)-sin(n0)*sin(p0)*sin(q0);aa2=-sin(q0)*cos(n0)-sin(n0)*sin(p0)*cos(q0);aa3=-sin(n0)*cos(p0);bb1=sin(q0)*cos(p0);bb2=cos(q0)*cos(p0);bb3=-sin(p0);cc1=sin(n0)*cos(q0)+sin(p0)*cos(n0)*sin(q0);cc2=-sin(n0)*sin(q0)+sin(p0)*cos(q0)*cos(n0);cc3=cos(n0)*cos(p0);%计算改正数XX1=aa1.*(b1(:,1)-X0)+bb1.*(b1(:,2)-Y0)+cc1.*(b1(:,3)-Z0); XX2=aa2.*(b1(:,1)-X0)+bb2.*(b1(:,2)-Y0)+cc2.*(b1(:,3)-Z0); XX3=aa3.*(b1(:,1)-X0)+bb3.*(b1(:,2)-Y0)+cc3.*(b1(:,3)-Z0); lx=a1(:,1)+f.*(XX1./XX3);ly=a1(:,2)+f.*(XX2./XX3);l={lx',ly'};L=cell2mat(l)';%方程系数A=[-3.969*10^-5 0 2.231*10^-5 -0.2 -0.04 -0.06899;0 -3.969*10^-5 1.787*10^-5 -0.04 -0.18 0.08615;-2.88*10^-5 0 1*10^-5 -0.17 0.03 0.08211;0 -2.88*10^-5 -1.54*10^-5 0.03 -0.2 0.0534;-4.14*10^-5 0 4*10^-6 -0.15 -7.4*10^-3 -0.07663;0 -4.14*10^-5 2.07*10^-5 -7.4*10^-3 -0.19 0.01478;-2.89*10^-5 0 -1.98*10^-6 -0.15 -4.4*10^-3 0.06443;0 -2.89*10^-5 -1.22*10^-5 -4.4*10^-3 -0.18 0.01046];%L=[-1.28 3.78 -3.02 -1.45 -4.25 4.98 -4.72 -0.385]'.*10^-2; %第一次迭代X=inv(A'*A)*A'*L;3.结果X=1492.41127406195-554.4015671761941425.68660973544-0.0383847815608609 0.00911624039769785 -0.105416434087641S=1492.41127406195-554.401567176194 1425.68660973544 38436.9616152184 27963.1641162404-0.105416434087641。
单像空间后方交会名词解释
单像空间后方交会名词解释
单像空间后方交会是摄影测量学中的一个重要概念，它是指利用单个影像进行地物测量和定位的方法。
在单像空间后方交会中，通过对单张影像进行分析，可以确定地面上物体的位置和形状。
这个过程涉及到对影像中的特征点进行识别和匹配，然后利用相机内外参数以及影像上的像点坐标来计算地物的三维坐标。
单像空间后方交会的过程包括以下几个步骤，首先是对影像进行预处理，包括去畸变、影像配准等操作；然后是特征点的提取和匹配，这一步是通过计算机视觉算法来实现的，可以利用角点、边缘等特征来进行匹配；接下来是相机内外参数的标定，这一步是为了将像素坐标转换为实际世界坐标而进行的；最后是利用已知的相机参数和像点坐标来计算地物的三维坐标。
单像空间后方交会在航空摄影、遥感影像解译和地图制图等领域有着广泛的应用。
它可以通过对单张影像的处理，实现对地物的测量和定位，为地理信息系统和地图制图提供了重要的数据基础。
同时，随着计算机视觉和图像处理技术的不断发展，单像空间后方交会的精度和效率也在不断提高，为各种应用领域提供了更加可靠和精确的地物信息。
摄影测量空间后方交会
摄影测量空间后方交会以单张影像空间后方交会方法，求解该像的外方位元素一、实验数据与理论基础：1、实验数据：航摄仪内方位元素f＝153.24mm，x0＝y0＝0，以及4对点的影像坐标和相应的地面坐标：影像坐标地面坐标x(mm)y(mm)X(m)Y(m)Z(m)1-86.15-68.9936589.4125273.322195.172-53.4082.2137631.0831324.51728.693-14.78-76.6339100.9724934.982386.50410.4664.4340426.5430319.81757.312、理论基础(1) 空间后方交会是以单幅影像为基础，从该影像所覆盖地面范围内若干控制点的已知地面坐标和相应点的像坐标量测值出发，根据共线条件方程，解求该影像在航空摄影时刻的外方位元素Xs，Ys，Zs，φ，ω，κ。
(2) 每一对像方和物方点可列出2个方程，若有3个已知地面坐标的控制点，可列出6个方程，求取外方位元素改正数△Xs，△Ys，△Zs，△φ，△ω，△κ。
二、数学模型和算法公式1、数学模型：后方交会利用的理论模型为共线方程。
共线方程的表达公式为：)()()()()()(333111S A S A S A S A S A S A Z Z c Y Y b X X a Z Z c Y Y b X X a fx -+-+--+-+--=)()()()()()(333222S A S A S A S A S A S A Z Z c Y Y b X X a Z Z c Y Y b X X a fy -+-+--+-+--=其中参数分别为：κωϕκϕsin sin sin cos cos 1-=aκωϕκϕsin sin sin sin cos 2--=a ωϕcos sin 3-=aκωsin cos 1=b κωcos cos 2=b ωsin 3-=bκωϕκϕsin sin cos cos sin 1+=c κωϕκϕcos sin cos sin sin 2+-=c ωϕcos cos 3=c旋转矩阵R 为⎥⎥⎥⎦⎤⎢⎢⎢⎣⎡=321321321c c c b b b a a a R2、 由于外方位元素共有6个未知数，根据上述公式可知，至少需要3个不在一条直线上的已知地面点坐标就可以求出像片的外方位元素。
后方交会法计算推导公式
后方交会法计算推导公式后方交会法是一种用于计算物体在空间中的坐标和距离的方法。
它基于两个观测者在不同位置观测同一个物体的现象。
假设有两个观测者A和B，在空间中观测同一个物体P。
观测者A 和B的位置分别为A(xA, yA, zA)和B(xB, yB, zB)。
物体P在观测者A和B的朝向上的投影分别为a和b，它们的长度分别为dA和dB。
根据几何关系，可以推导出以下公式：dA = sqrt((xA - xP)^2 + (yA - yP)^2 + (zA - zP)^2)dB = sqrt((xB - xP)^2 + (yB - yP)^2 + (zB - zP)^2)其中，(xP, yP, zP)是物体P的坐标。
如果已知dA、dB和相关观测者位置的坐标，可以使用这些公式来计算物体P的坐标(xP, yP, zP)。
同时，如果已知物体P在两个观测者朝向上的投影长度a和b，也可以利用这些公式计算物体P到观测者A和B的距离。
需要注意的是，后方交会法在实际应用中可能会受到观测误差的影响，因此在计算时需要考虑这些误差，并采取合适的数据处理和精度控制方法。
拓展：后方交会法是测量和定位的重要方法之一，广泛应用于地理测量、摄影测量、建筑工程等领域。
它可以通过精确的测量和计算，确定物体在三维空间中的准确位置和形状，对于工程设计、地理信息系统等具有重要的实际应用价值。
除了后方交会法，还有其他一些方法可以用于测量和定位物体的坐标和距离，比如三角测量法、三角高程测量法、全站仪测量法等。
每种方法都有其适用的场景和局限性，根据具体的测量需求和条件选择合适的方法是非常重要的。
此外，随着科技的进步和发展，新的测量和定位技术不断涌现，为实现更精确和高效的测量和定位提供了更多的选择。
摄影测量学空间后方交会实验报告
摄影测量学实验报告实验一、单像空间后方交会学院：建测学院班级：测绘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．实验结果四．实验总结此次实验让我深入了解单像空间后方交会的计算过程，加强了对空间后方交会基本公式和误差方程式，法线方程式的记忆。
单像空间后方交会实验报告(c++版)
单像空间后方交会实验报告(c++版)单像空间后方交会姓名：学号：时间：目录一、作业任务..................................................... - 4 -二、计算原理..................................................... - 4 -三、算法流程..................................................... - 8 -四、源程序....................................................... - 9 -五、计算结果..................................................... - 9 -六、结果分析..................................................... - 9 -七、心得与体会................................................... - 9 -八、附页......................................................... - 9 -1.c++程序.................................................. - 10 -2.C++程序截图.............................................. - 17 -3.matlb程序 ............................................... - 17 -一、 作业任务已知条件：摄影机主距f=153.24mm ，x0=0，y0=0, 像片比例尺为1:40000，有四对点的像点坐标与相应的地面坐标如下表。
摄影测量学空间后方交会实验报告13页
摄影测量学空间后方交会实验报告13页报告摘要：本实验以三张已知高度物体的相片为样本，运用摄影测量学的空间后方交会算法，通过MATLAB编程实现对物体三维坐标的计算，并加以检验与探讨，得出结论相对合理。
同时，本实验也通过对MATLAB编程的应用，掌握了空间后方交会算法的理论及实践方法。
1. 实验目的（1）学习和应用摄影测量学中的空间后方交会算法，掌握其计算方法和程序实现过程。
（2）了解数字像相机的相关特性，并掌握其使用方法。
（3）通过对样本数据的处理，熟悉和掌握MATLAB编程的应用技巧。
2. 实验器材数字相机一部，尺子一把，样本图像三张，MATLAB软件。
3. 实验原理在精密测量领域中，采用摄影测量学的空间后方交会算法，可实现三维坐标的测量和重建。
这种方法是将已知物体的照片，通过对像点的提取及校准，得出像点坐标系下的物体的三维坐标系下的坐标。
这种算法的基本思路是：利用像点坐标系下的物体三维坐标系下的坐标关系，构建一个误差最小的方程组，通过矩阵的求解，得到物体三维坐标系下的坐标。
数字相机是一种基于CCD或CMOS成像器材料的成像设备，根据数字信号的处理能力，合成电子图像。
数字相机的性能主要包括分辨率、感光度、曝光控制、焦距、光圈等参数。
使用数字相机拍摄时，应根据拍摄对象的光线条件、距离、尺寸、景深等因素，进行调节。
4. 实验过程（1）利用数字相机拍摄三张已知高度物体的照片，并在样本上面贴标记，用尺子测算实际高度。
（2）利用图像处理软件MATLAB，对照片进行像点识别和校准，得到像点坐标系下的坐标。
（3）根据相片中已知物体的测高值及像点坐标系下的坐标值，通过MATLAB编写空间后方交会的程序算法，得出物体的三维坐标系下的坐标。
（4）对得出的坐标值进行检验及探讨，分析误差来源及部分工具库的使用方法。
5. 实验结果与分析本实验所得出的三维坐标值，原本应是在一个确定点之间展开的点集。
知道参数计算不全或精度不够是常有的事情（尤其在没有精密测量器具的条件下），这种情况我们应该考虑从给出的角度和图像来看和计算得到位置。
空间后方交会原理
空间后方交会原理.txtゅ你不用一上线看见莪在线，就急着隐身，放心。
莪不会去缠你。
说好的不离不弃现在反而自己却做不到╮单幅航空影像空间后方交会程序设计来源: 摘要：航空影像单像空间后方交会是利用影像所覆盖地面范围内若干控制点的已知地面坐标和相应点的像点坐标，根据共线方程，在最小二乘法的原则下，计算航空影像的外方位元素。
本文采用 编译平台和C#语言编写了单像空间后方交会的参数计算程序。
利用相关文献中的实验数据与结果，证明本程序运行良好，参数解算精确，可以满足工程项目要求。
关键词：摄影测量与遥感；外方位元素；共线方程；最小二乘原则；空间后方交会；0 引言确定影像或摄影光束在摄影瞬间的空间位置和姿态的参数，称为影像的外方位元素。
一幅影像的外方位元素包括6 个参数，其中3 个参数线参数，用于描述摄影中心S 相对于物方空间坐标系的位置（X，Y，Z）；另外3 个是角元素，用于描述影像面在摄影瞬间的空中姿态。
角元素有三种不同的表达形式：（1）以Y 轴为主轴的φ-ω-κ系统（主轴是在旋转过程中，空间方向不变的一个固定轴）：以Y 为主轴旋转φ角，然后绕X 轴旋转ω角，最后围绕Z 轴旋转κ角。
（2）以X 轴为主轴的φ?-ω?-κ?系统：以X 为主轴旋转ω?角，然后绕轴旋转φ?角，最后围绕Z 轴旋转κ?角。
（3）以Z 轴为主轴的A-α-κ系统：以Z 为主轴旋转A 角，然后绕Y 轴旋转α角，最后围绕Z 轴旋转κ角。
本文的角元素系统选择使用第一种表达方式，即φ-ω-κ[1]。
如果知道了每幅影像的6 个外方位元素，就能确定被摄物体与航摄影像的关系。
因此，如何获取影像的外方位元素，一直是摄影测量工作者所探讨的问题。
目前，从技术方面来说，解决这个问题主要有两种手段：一种是利用雷达、全球定位系统（GPS），惯性导航系统（INS）以及星相摄影机来获取影像的外方位元素；另一种是利用影像覆盖范围内一定数量的控制点的空间坐标与摄像坐标，根据共线条件方程，反求该影像的外方位元素，这种方法称为单幅影像的空间后方交会。
摄影测量学单像空间后方交会程序设计作业
using System;using System.Collections.Generic;using System.Linq;using System.Text;namespace 单像空间后方交会{class Program{static void Main(string[] args){int x0, y0, i, j; double f, m;Console.Write("请输入像片比例尺：");m = double.Parse(Console.ReadLine());Console.Write("请输入像片的内方位元素x0:");//均以毫米为单位x0 = int.Parse(Console.ReadLine());Console.Write("请输入像片的内方位元素y0:");y0 = int.Parse(Console.ReadLine());Console.Write("请输入摄影机主距f:");f = double.Parse(Console.ReadLine());Console.WriteLine();//输入坐标数据double[,] zuobiao = new double[4, 5];for (i = 0; i < 4; i++){for (j = 0; j < 5; j++){if (j < 3){Console.Write("请输入第{0}个点的第{1}个地面坐标：", i + 1, j + 1);zuobiao[i, j] =double.Parse(Console.ReadLine());}else{Console.Write("请输入第{0}个点的第{1}个像点坐标：", i + 1, j - 2);zuobiao[i, j] =double.Parse(Console.ReadLine());}} Console.WriteLine();}//归算像点坐标for (i = 0; i < 4; i++){for (j = 3; j < 5; j++){if (j == 3)zuobiao[i, j] = zuobiao[i, j] - x0;elsezuobiao[i, j] = zuobiao[i, j] - y0;}}//计算和确定初值double zs0 = m * f, xs0 = 0, ys0 = 0;for (i = 0; i < 4; i++){xs0 = xs0 + zuobiao[i, 0];ys0 = ys0 + zuobiao[i, 1];}xs0 = xs0 / 4;ys0 = ys0 / 4;//逐点计算误差方程系数double[,] xishu = new double[8, 6];for (i = 0; i < 8; i += 2){double x, y;x = zuobiao[i / 2, 3]; y = zuobiao[i / 2, 4];xishu[i, 0] = xishu[i + 1, 1] = -1 / m; xishu[i, 1] = xishu[i + 1, 0] = 0; xishu[i, 2] = -x / (m * f); xishu[i, 3] = -f * (1 + x * x / (f * f));xishu[i, 4] = xishu[i + 1, 3] = -x * y / f; xishu[i, 5] = y; xishu[i + 1, 2] = -y / (m * f); xishu[i + 1, 4] = -f * (1 + y * y / (f * f)); xishu[i + 1, 5] = -x;}//计算逆阵double[,] dMatrix =matrixChe(matrixTrans(xishu), xishu);double[,] dReturn = ReverseMatrix(dMatrix, 6);Console.WriteLine("逆矩阵为：");if (dReturn != null){matrixOut(dReturn);}//求解过程double phi0 = 0, omega0 = 0, kappa0 = 0; int q = 0;double[,] r = new double[3, 3];double[,] jinsi = new double[4, 2];double[] chazhi = new double[8];double[] jieguo = new double[6];double[,] zhong = matrixChe(dReturn,matrixTrans(xishu));do{ //计算旋转矩阵rr[0, 0] = Math.Cos(phi0) * Math.Cos(kappa0) - Math.Sin(phi0) * Math.Sin(omega0) * Math.Sin(kappa0);r[0, 1] = -Math.Cos(phi0) * Math.Sin(kappa0) - Math.Sin(phi0) * Math.Sin(omega0) * Math.Cos(kappa0);r[0, 2] = -Math.Sin(phi0) * Math.Cos(omega0);r[1, 0] = Math.Cos(omega0) * Math.Sin(kappa0);r[1, 1] = Math.Cos(omega0) * Math.Cos(kappa0);r[1, 2] = -Math.Sin(omega0);r[2, 0] = Math.Sin(phi0) * Math.Cos(kappa0) + Math.Cos(phi0) * Math.Sin(omega0) * Math.Sin(kappa0);r[2, 1] = -Math.Sin(phi0) * Math.Sin(kappa0) + Math.Cos(phi0) * Math.Sin(omega0) * Math.Cos(kappa0);r[2, 2] = Math.Cos(phi0) * Math.Cos(omega0);//计算x,y的近似值for (i = 0; i < 4; i++){jinsi[i, 0] = -f * (r[0, 0] * (zuobiao[i, 0] - xs0) + r[1, 0] * (zuobiao[i, 1] - ys0) + r[2, 0] * (zuobiao[i, 2] - zs0)) / (r[0, 2] * (zuobiao[i, 0] - xs0) + r[1, 2] * (zuobiao[i, 1] - ys0) + r[2, 2] * (zuobiao[i, 2] - zs0));jinsi[i, 1] = -f * (r[0, 1] * (zuobiao[i, 0] - xs0) + r[1, 1] * (zuobiao[i, 1] - ys0) + r[2, 1] * (zuobiao[i, 2] - zs0)) / (r[0, 2] * (zuobiao[i, 0] - xs0) + r[1, 2] * (zuobiao[i, 1] - ys0) + r[2, 2] * (zuobiao[i, 2] - zs0));}for (i = 0; i < 8; i += 2){chazhi[i] = zuobiao[i / 2, 3] - jinsi[i / 2, 0];chazhi[i + 1] = zuobiao[i / 2, 4] - jinsi[i / 2, 1];}for (i = 0; i < zhong.GetLength(0); i++){double k = 0;for (j = 0; j < zhong.GetLength(1); j++){k = k + zhong[i, j] * chazhi[j];}jieguo[i] = k;}//求新的近似值xs0 += jieguo[0]; ys0 += jieguo[1]; zs0 += jieguo[2];phi0 += jieguo[3]; omega0 += jieguo[4]; kappa0 += jieguo[5];q++;if (q > 1000)break;} while ((Math.Abs(jieguo[0]) > 0.020 ||Math.Abs(jieguo[1]) > 0.020) || Math.Abs(jieguo[2]) > 0.020);Console.WriteLine("共进行了{0}次运算", q);Console.WriteLine("旋转矩阵为");matrixOut(r);for (i = 0; i < jieguo.GetLength(0); i++){Console.Write("第{0}个外方位元素为：{1}", i + 1, jieguo[i]);}}//矩阵转置public static double[,] matrixTrans(double[,] X){double[,] A = X;double[,] C = new double[A.GetLength(1),A.GetLength(0)];for (int i = 0; i < A.GetLength(1); i++)for (int j = 0; j < A.GetLength(0); j++){C[i, j] = A[j, i];}return C;}//矩阵输出public static void matrixOut(double[,] X){double[,] C = X;for (int i = 0; i < C.GetLength(0); i++){for (int j = 0; j < C.GetLength(1); j++){Console.Write(" {0}", C[i, j]);}Console.Write("\n");}}//二维矩阵相乘public static double[,] matrixChe(double[,] X, double[,] Y){int i, j, n; double m;double[,] C = X; double[,] D = Y;double[,] E = new double[C.GetLength(0),C.GetLength(0)];for (i = 0; i < C.GetLength(0); i++){for (n = 0; n < C.GetLength(0); n++){m = 0;for (j = 0; j < C.GetLength(1); j++){m = m + C[i, j] * D[j, n];}E[i, n] = m;}}return E;}//计算行列式的值public static double MatrixValue(double[,] MatrixList, int Level){double[,] dMatrix = new double[Level, Level];for (int i = 0; i < Level; i++)for (int j = 0; j < Level; j++)dMatrix[i, j] = MatrixList[i, j];double c, x;int k = 1;for (int i = 0, j = 0; i < Level && j < Level; i++, j++){if (dMatrix[i, j] == 0){int m = i;for (; dMatrix[m, j] == 0; m++) ;if (m == Level)return 0;else{for (int n = j; n < Level; n++){c = dMatrix[i, n];dMatrix[i, n] = dMatrix[m, n];dMatrix[m, n] = c;}k *= (-1);}}for (int s = Level - 1; s > i; s--){x = dMatrix[s, j];for (int t = j; t < Level; t++)dMatrix[s, t] -= dMatrix[i, t] * (x / dMatrix[i, j]);}}double sn = 1;for (int i = 0; i < Level; i++){if (dMatrix[i, i] != 0)sn *= dMatrix[i, i];elsereturn 0;}return k * sn;}//计算逆阵public static double[,] ReverseMatrix(double[,] dMatrix, int Level){double dMatrixValue = MatrixValue(dMatrix, Level);if (dMatrixValue == 0) return null;double[,] dReverseMatrix = new double[Level, 2 * Level];double x, c;for (int i = 0; i < Level; i++){for (int j = 0; j < 2 * Level; j++){if (j < Level)dReverseMatrix[i, j] = dMatrix[i, j];elsedReverseMatrix[i, j] = 0;}dReverseMatrix[i, Level + i] = 1;}for (int i = 0, j = 0; i < Level && j < Level; i++, j++){if (dReverseMatrix[i, j] == 0){int m = i;for (; dMatrix[m, j] == 0; m++) ;if (m == Level)return null;else{for (int n = j; n < 2 * Level; n++)dReverseMatrix[i, n] += dReverseMatrix[m, n];}}x = dReverseMatrix[i, j];if (x != 1){for (int n = j; n < 2 * Level; n++)if (dReverseMatrix[i, n] != 0)dReverseMatrix[i, n] /= x;}for (int s = Level - 1; s > i; s--){x = dReverseMatrix[s, j];for (int t = j; t < 2 * Level; t++)dReverseMatrix[s, t] -= (dReverseMatrix[i, t] * x);}}for (int i = Level - 2; i >= 0; i--){for (int j = i + 1; j < Level; j++)if (dReverseMatrix[i, j] != 0){c = dReverseMatrix[i, j];for (int n = j; n < 2 * Level; n++)dReverseMatrix[i, n] -= (c * dReverseMatrix[j, n]);}}double[,] dReturn = new double[Level, Level];for (int i = 0; i < Level; i++)for (int j = 0; j < Level; j++)dReturn[i, j] = dReverseMatrix[i, j + Level];return dReturn;}}}。
摄影测量实验报告(空间后方交会—前方交会)
空间后方交会—空间前方交会程序编程实验一．实验目的要求掌握运用空间后方交会-空间前方交会求解地面点的空间位置.学会运用空间后方交会的原理，根据所给控制点的地面摄影测量坐标系坐标以及相应的像平面坐标系中的坐标，利用计算机编程语言实现空间后方交会的过程，完成所给像对中两张像片各自的外方位元素的求解。
然后根据空间后方交会所得的两张像片的内外方位元素，利用同名像点在左右像片上的坐标，求解其对应的地面点在摄影测量坐标系中的坐标，并完成精度评定过程，利用计算机编程语言实现此过程.二．仪器用具计算机、编程软件(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 '左片外方位元素为。
单像空间后方交会
单像空间后方交会测绘学院 成晓倩1 概述1.1 定义利用一定数量的地面控制点和对应像点坐标求解单张像片外方位元素的方法称为空间后方交会。
1.2 所需控制点个数与分布共线条件方程的一般形式为：⎪⎪⎩⎪⎪⎨⎧-+-+--+-+--=--+-+--+-+--=-)()()()()()()()()()()()(33322203331110S S S S S S S S S S S S Z Z c Y Y b X X a Z Z c Y Y b X X a f y y Z Z c Y Y b X X a Z Z c Y Y b X X a f x x （1）式中包含有六个外方位元素，即κωϕ、、、、、S S S Z Y X ，只有确定了这六个外方位元素的值，才能利用共线条件方程真正确定一张像片的任一像点与对应地面点的坐标关系。
个数：对任一控制点，我们已知其地面坐标)(i i i Z Y X 、、和对应像点坐标)(i i y x 、，代入共线条件方程可以列出两个方程式，因此，只少需要3个控制点才能解算出六个外方位元素。
在实际应用中，为了避免粗差，应有多余检查点，因此，一般需要4～6个控制点。
分布：为了最有效地控制整张像片，控制点应均匀分布于像片边缘，如下图所示。
由于共线条件方程是非线性的，直接答解十分困难，所以首先将共线方程改化为线性形式，然后再答解最为简单的线性方程组。
2 空间后方交会的基本思路分布合理 分布合理 分布不合理2.1 共线条件方程线性化的基本思路在共线条件方程中，令)()()()()()()()()(333222111S S S S S S S S S Z Z c Y Y b X X a Z Z Z c Y Y b X X a Y Z Z c Y Y b X X a X -+-+-=-+-+-=-+-+-= （2） 则共线方程变为⎪⎪⎩⎪⎪⎨⎧-=--=-ZY fy y Z Xf x x 00 （3） 对上式两侧同乘Z ，并移至方程同侧，则有⎩⎨⎧=-+=-+0)(0)(00Z y y Y f Z x x X f （4） 令⎩⎨⎧-+=-+=Zy y Y f Fy Zx x X f Fx )()(00 （5） 由于上式是共线方程的变形，因此，Fy Fx 、是κωϕ、、、、、S S S Z Y X 的函数。
(空间后方交会的计算过程)空间后方交会
将上述偏导数代入，可以求得其余的系数如下
x ( x cos k y sin k ) f cos k ] cos f x a15 f sin k ( x sin k y cos k ) f a16 y a14 y sin [ x ( x cos k y sin k f sin k ) f sin k ] cos f y a25 f cos k ( x sin k y cos k f a26 x a24 x sin [
计算中，通常将地面控制点的坐标认为是真值，而把相应的像点 Vy 列 坐标认为是观测值，加入相应的改正数 Vx ,Vy ，得 x Vx , y ， 出如下的每个点的误差方程式为：
x x x x x x V dX dY dZ d d dk ( x) x S S S x X Y Z k S S S V y dX y dY y dZ y d y d y dk ( y ) y y S S S X Y Z k S S S
当竖直投影时，角元素都是小角（小于3度），此时可近似认为 k 0, Z A Z S H ，各个系数的表达式可以得到简化。
空间后方交会计算中的误差方程和法方程 由于有六个未知数，所以至少需要知道三个 已知的地面控制点，为了能够平差，通常在 像片的四个角选取四个或更多的地面控制点。
1 4 YS 0 Ytpi 4 i 1
4) 计算旋转矩阵R：利用角元素的近似值计算 方向元素，组成旋转矩阵R。 5）逐点计算像点坐标的近似值：利用未知数 的近似值按照共线方程计算控制点像点坐 标的近似值(x),(y)； 6) 组成误差方程式 7) 组成法方程式 8）解求外方位元素 9）检查计算是否收敛：将求得外方位元素的 改正数与规定的限差比较，小于限差则计 算终止，否则迭代计算。
空间后方交会的计算例题
一、选择题1. 在空间后方交会中，如果两条视线交会于点A，且已知视线与水平面的夹角为30度和45度，计算交点A的高度。
- A. 50米- B. 75米- C. 100米- D. 125米2. 若在空间后方交会中，已知两点的坐标为A(10, 20, 30)和B(15, 25, 35)，计算这两点之间的距离。
- A. 10米- B. 15米- C. 20米- D. 25米3. 在空间后方交会中，假设三点的视线交会分别为0.5°、1.0°和1.5°，计算其交会点的相对误差。
- A. 0.1米- B. 0.2米- C. 0.3米- D. 0.4米4. 若一个目标点在空间中，观测点1的坐标为(100, 200, 300)，观测点2的坐标为(150, 250, 350)，计算两观测点到目标点的视线夹角为30度和60度。
目标点的高度为多少？- A. 250米- B. 300米- C. 350米- D. 400米5. 在空间后方交会中，已知一个点A的坐标为(5, 5, 5)，另一个点B的坐标为(10, 10, 10)，计算两点之间的直线距离。
- A. 7.07米- B. 10米- C. 12.25米- D. 15米6. 如果在空间后方交会中，已知三个点A、B、C的坐标分别为(0, 0, 0)、(4, 0, 0)和(0, 3, 0)，计算点C到点AB的距离。
- A. 3米- B. 4米- C. 5米- D. 6米7. 在空间后方交会中，若两条视线分别与X轴和Y轴成30度和45度角，计算这两条视线交会点的Z坐标。
- A. 50米- B. 60米- C. 70米- D. 80米8. 一个空间点P的坐标为(20, 30, 40)，在三条视线交会中，视线与X、Y、Z轴的夹角分别为30度、45度和60度。
计算点P的相对坐标误差。
- A. 1.5米- B. 2.0米- C. 2.5米- D. 3.0米9. 在空间后方交会中，如果已知视线A和视线B的夹角为90度，且两条视线交会点的坐标为(0, 0, 0)和(100, 100, 0)，计算第二个交会点的坐标。
空间后方交会基本原理
F y X F y S 0 d S X F Y y S 0 d S Y F Z y S 0 d S Z F y 0 d F y 0 d F y 0 d F y 0
因为
Fx xf
a1(XXS)b1(YYS)c1(ZZS) 0 a3(XXS)b3(YYS)c3(ZZS)
④在摄影过程中直接获取。
h
4
内
• 单像空间后方交会概述
容
• 共线方程的线性化（难点） • 利用共线条件方程解算像片的外方位元
素（重点）
安
排
[一]概述
1、什么叫单像空间后方交会 利用地面控制点及其在片像上的像点，确定一
张像片外方位元素的方法。
S (XS 、 YS 、 ZS)
c
Z
Y
b
a
C
A X 地面控B制点(Ground Control Point, GCP)
第十讲 空间后方交会
已学过的主要内容
绪论 摄影与航空摄影 第一部分 单张航摄像片的解析 第二部分 立体像对的基本知识 第三部分 作业基本理论
问题的引出
XS，YS，ZS
YXYXss((ZZZZs)s)ba1cx11xxb2acy22yyb3ca3f3ff
若 ai，bi，ci
f x，y ，Z
已知
X，Y
h
6
[一]概述
1、什么叫单像空间后方交会 利用地面控制点及其在像片上的像点，确定一
张像片外方位元素的方法。 2、单像空间后方交会的基本方法
角锥体法
S
cb a
ZT
C B
YT
A
D
XT
h
7
角锥体法
S
c
b
空间后方交会实验报告
空间后方交会实验报告1. 引言空间后方交会是一种常用的测量手段，用于确定目标在空间中的坐标位置。
本实验旨在通过空间后方交会实验，了解空间测量的基本原理和方法，并锻炼我们的观察、计算和推理能力。
2. 实验目的- 学习和掌握空间后方交会实验的基本原理和步骤；- 熟悉使用测量仪器和处理数据的方法；- 锻炼团队合作和问题解决能力。
3. 实验原理空间后方交会利用多个测量设备，通过测量目标物体在不同位置上的观测数据，计算出目标物体的空间坐标。
主要基于以下原理：- 多点定位原理：通过多次测量目标物体在不同位置上的观测数据，可以确定目标物体的坐标位置；- 观测误差校正原理：由于观测仪器的误差存在，需要对观测数据进行误差校正，以提高测量精度。
4. 实验步骤4.1 实验准备- 确定实验区域，搭建固定的测量基准点；- 使用全站仪对测量基准点进行校准，确保测量精度；- 准备至少三个可移动的晶体棱镜，用于固定在目标物体上；- 准备笔记本电脑和相应的测量软件，用于数据处理和计算。
4.2 实验操作1. 将晶体棱镜固定在目标物体上，并确定测量的起始位置；2. 使用全站仪测量目标物体在不同位置上的观测数据，包括水平角、垂直角和斜距；3. 将观测数据导入计算软件，并进行误差校正；4. 根据误差校正后的观测数据，计算目标物体的空间坐标；5. 重复步骤2至4，获得更多的观测数据，以提高计算精度；6. 对多次观测结果进行平均计算，得到最终的目标物体空间坐标；7. 将计算结果进行分析和评估，确定测量精度和可靠度。
5. 实验结果与讨论经过实验操作和数据处理，获得了目标物体的空间坐标。
通过对测量结果的分析和比较，可以得出以下结论：- 目标物体的空间坐标可以通过空间后方交会实验进行测量和计算；- 观测数据的误差校正对测量结果的精度和可靠度具有重要影响；- 重复观测和平均计算可以提高测量结果的准确性。
6. 实验总结空间后方交会实验是一种常用的测量手段，可以用于测量目标物体在空间中的坐标位置。
空间后方交会
空间后方交会
科技名词定义
中文名称：空间后方交会
英文名称：space resection
定义：利用航摄像片上三个以上不在一条直线上的已知点按构像方程计算该像片外方位元
素的方法。
空间后方交会，是指利用航摄像片上三个以上不在一条直线上的控制点按构像方程计算该像片外方位元素的方法。
内容
仅在待定点上设站，向三个已知控制点观测两个水平夹角a、b，从而计算待定点的坐标，称为后方交会。
交会测量是加密控制点常用的方法，它可以在数个已知控制点上设站，分别向待定点观测方向或距离，也可以在待定点上设站向数个已知控制点观测方向或距离，而后计算待定点的坐标。
常用的交会测量方法有前方交会、后方交会、侧边交会和自由测站法。
后方交会法首先出现于测绘地形图工作中，测量上称为“三点题”，是用图解法作为加密图根点之用。
后来随着解析法、公式法的出现，在工程建设控制测量中也经常被采用。
比如隧道工程控制网往往由于隧道开工前测设完成，而洞口土石方施工完毕后，需补设洞口投点，以便控制隧道轴线，测设投点就要用到后方交会法；深水桥墩放样测量中的墩心定位也可以应用此法，还可用来测定施工控制导线的始终点等。
应用范围之广说明了此法的实用性很强。
其代表图形如下图所示。
图中三角形ABC是控制网中的一个三角形，P点即为后方交会点（需确定坐标的待定点），只要置棱镜于P点，用全圆测回法测定a，b，r 三个角值，即可应用解析法公式算出待定点P的坐标。
此法内外业工作量小，只要P点的点位精度符合施工放样要求或作为洞口投点的精度要求，就可以成为广大测绘科技工作者所乐意选用的方法之一。
(空间后方交会的计算过程)空间后方交会 PPT
大家好
11
计算中，通常将地面控制点的坐标认为是真值，而把相应的像点 坐标认为是观测值，加入相应的改正数 V x ,V y ，得 xVx,y，Vy 列
出如下的每个点的误差方程式为：
V x X xSdS X Y x SdSY Z xSdSZ xd xd k xd k(x)x V y X ySdS X Y y SdSY Z ySdSZ yd yd k yd k(y)y
式中，x,y为像点坐标的观测值，(x),(y)为用控制点的物方坐标及 外方位元素的近似值大家代好入中心投影方程求得的像点坐标近12似值。
舍弃二次项，使之线性化得：
x y ((xy)) X X xySSd dX S X S Y Y xySSd dSS Y Y Z Z xySSd dS Z S Z x yd d xyd d k xk yd dkk
式中，(x),(y)为函数的近似值。
k
f
sin
k)
f
sin
k ] cos
a25
f
cos k
y f
(x sin k
y cos k
a26 x
当竖直投影时，角元素都是小角（小于3度），此时可近似认为
k 0 ,Z A Z S H ，各个系数的表达式可以得到简化。
大家好
10
空间后方交会计算中的误差方程和法方程
由于有六个未知数，所以至少需要知道三个 已知的地面控制点，为了能够平差，通常在 像片的四个角选取四个或更多的地面控制点。
第三章 单张航摄像片解析
大家好
1
单张像片的空间后方交会
如果我们有每张像片的六个外方位元素，就 能恢复航摄像片与被摄地面之间的几何关系， 重建地面的立体模型。
空间后方交会原理
式2-6，2-7其中：X 和Y 是控制点的地面坐标；和y'是控制点在摄影测量坐标系下的坐标值。
由此可以解算出参数X Y a,b,T ,T ，据此确定参数的初始值：
κ式2-8，注意：κ 值的确定采用的是全圆反正切，这和初等数学中的反正切图像不一样。该函数是计算机系统中的一个数学函数，将最后所得值归算到-180°到+180°之间。
文章来自中国硕士论文网，更多工程硕士论文就在中国硕士论文网，如需转载请保留一个链接：
参考文献
张剑清等著．摄影测量学［M］.武汉：武汉大学出版社
李德仁等著．摄影测量与遥感概论［M］.北京：测绘出版社
张凤举，张华海等编著．控制测量学［M］.北京：煤炭工业出版社出版社
计算方向余弦矩阵方向余弦矩阵的形式如下：
式方向余弦的各个元素的计算方法求解影像的外方位元素的参数根据控制点的个数和坐标值，确定法方程，以其中一点A 为例，解释如下：各个系数的值即是2-2 中的各个改正数的系数值。
如果有四个已知的控制点，那么矩阵B就是8× 6的矩阵。使用四个以上的控制点，可以解出Δ的值（Δ = (BT B)?1 (BTε )）。然后重复步骤和1.3 进行一次迭代运算，即可完成单像空间后方交会的工作。
- 1、下载文档前请自行甄别文档内容的完整性，平台不提供额外的编辑、内容补充、找答案等附加服务。
- 2、"仅部分预览"的文档,不可在线预览部分如存在完整性等问题,可反馈申请退款(可完整预览的文档不适用该条件!)。
- 3、如文档侵犯您的权益，请联系客服反馈,我们会尽快为您处理(人工客服工作时间：9:00-18:30)。
作业报告空间后方交会专业：测绘工程班级：2008级(1)班姓名：陈闻亚指导教师：陈强2010 年 4 月16 日1 作业任务------------------------------------------------------------------------------------ 32 作业思想--------------------------------------------------------------------------------------- 3 3 作业条件及数据-------------------------------------------------------------------- 34 作业过程--------------------------------------------------------------------------- 35 源程序----------------------------------------------------------------------------- 46 计算结果---------------------------------------------------------------------------177心得体会与建议----------------------------------------------------------------------------- 171 作业任务计算近似垂直摄影情况下后方交会解。
即利用摄影测量空间后方交会的方法，获取相片的6个外方位元素。
限差为0.1。
2作业思想利用摄影测量空间后方交会的方法求解。
该方法的基本思想是利用至少三个一直地面控制点的坐标A（X A，Y A，Z A）、B（X B，Y B，Z B）C（X C，Y C，Z C），与其影像上对应的三个像点的影像坐标a（x a，y a）、b（x b，y b）、c（x c，y c），根据共线方程，反求该相片的外方位元素X S、Y S、Z S、φ、ω、κ。
3作业条件及数据已知摄影机主距f=153.24mm，四对点的像点坐标与相应的地面坐标列入下表：4作业过程4．1 获取已知数据相片比例尺1/m=1：10000，内方位元素f=153.24mm，x0，y0；获取控制点的地面测量坐标X t、Y t、Z t。
4．2 量测控制点的像点坐标：本次作业中为已知。
见表1。
4．3 确定未知数的初始值：在近似垂直摄影情况下，胶原素的初始值为0，即φ0 = ω0 = κ0=0；线元素中，Z S0=H=mf=1532.4m，X S0、Y S0的取值可用四个控制点坐标的平均值，即：X S0= 144ii1X=∑=38437.00Y S0= 144ii1Y=∑=89106.624．4 计算旋转矩阵R：利用胶原素的近似值计算方向余弦值，组成R阵。
4．5 逐点计算像点坐标的近似值：利用未知数的近似值按共线方程式计算控制点像点坐标的近似值（x）（y）。
4．6 组成误差方程：逐点计算误差方程式的系数和常数项。
4．7 组成法方程式：计算法方程的系数矩阵A T A与常数项A T L。
4．8 求解外方位元素：根据法方程，由式X=（A t A）-1 A T L解求外方位元素改正数，并与相应的近似值求和，得到外方位元素新的近似值。
4．9 求解外方位元素：将求得的外方位元素的改正数与规定的限差（0.1）比较，小于限差则计算终止，否则用新的近似值重复第4.4至4.8步骤的计算，知道满足要求为止。
5 源程序#include <stdio.h>#include <stdlib.h>#include <math.h>const double PRECISION=1e-5;typedef double DOUBLE[5];int InputData(int &Num, DOUBLE *&Data,double &m,double &f);int Resection(const int &Num,const DOUBLE *&Data,const double &m,const double &f);int InverseMatrix(double *matrix,const int &row);int main(int argc, char* argv[]){DOUBLE *Data=NULL;int Num;double f(0),m(0);if(InputData(Num,Data,m,f)){if (Data!=NULL){delete []Data;}return 1;}if(Resection(Num,Data,m,f)){if (Data!=NULL){delete []Data;}return 1;}if (Data!=NULL){delete []Data;}printf("解算完毕...\n");do{printf("计算结果保存于\"结果.txt\"文件中\n""请选择操作（输入P打开结果数据，R打开原始数据，其它退出程序）：");fflush(stdin); //刷新输入流char order=getchar();if ('P'==order || 'p'==order){system("结果.txt");}else if ('R'==order || 'r'==order){system("data.txt");}elsebreak;system("cls");}while(1);system("PAUSE");return 0;}/********************************************** *函数名：InputData*函数介绍：从文件（data.txt）中读取数据,*文件格式如下：*点数m(未知写作0)* 内方位元素(f x0 y0)*编号x y X Y Z*实例：4 0153.24 0 01 -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*参数：(in/out)Num(点数)，*(in/out)Data(存放数据),m,f,x0,y0*返回值：int ，0成功，1文件打开失败，2控制点个*数不足，3文件格式错误**********************************************/int InputData(int &Num, DOUBLE *&Data,double &m,double &f) {double x0,y0;FILE *fp_input;if (!(fp_input=fopen("data.txt","r"))){return 1;}fscanf(fp_input,"%d%lf",&Num,&m);if (Num<4){return 2;}fscanf(fp_input,"%lf%lf%lf",&f,&x0,&y0); f/=1000;if (m<0 || f<0){return 3;}Data=new DOUBLE[Num];double *temp= new double[Num-1]; double scale=0;int i;for (i=0;i<Num;i++){//读取数据，忽略编号if(fscanf(fp_input,"%*d%lf%lf%lf%lf%lf",&Data[i][0],&Data[i][1],&Data[i][2],&Data[i][3],&Data[i][4])!=5){return 3;}//单位换算成mData[i][0]/=1000.0;Data[i][1]/=1000.0;}//如果m未知则归算其值if (0==m){for (i=0;i<Num-1;i++){temp[i]=(Data[i][2]-Data[i+1][2])/(Data[i][0]-Data[i+1][0])+ (Data[i][3]-Data[i+1][3])/(Data[i][1]-Data[i+1][1]);scale+=temp[i]/2.0;}m=scale/(Num-1);}fclose(fp_input);delete []temp;return 0;}/***********************************************函数名：MatrixMul*函数介绍：求两个矩阵的积,*参数：Jz1(第一个矩阵),row(第一个矩阵行数),*Jz2(第二个矩阵),row(第二个矩阵列数),com(第一个*矩阵列数),(out)JgJz(存放结果矩阵)*返回值：void**********************************************/void MatrixMul(double *Jz1,const int &row,double *Jz2, const int &line,const int &com,double *JgJz) {for (int i=0;i<row;i++){for (int j=0;j<line;j++){double temp=0;for (int k=0;k<com;k++){temp+=*(Jz1+i*com+k)*(*(Jz2+k*line+j));}*(JgJz+i*line+j)=temp;}}}/***********************************************函数名：OutPut*函数介绍：向结果.txt文件输出数据*参数：Q协因数阵，m精度，m0单位权中误差，6个外*方位元素，旋转矩阵*返回值：int，0成功，1失败**********************************************/int OutPut(const double *&Q,const double *&m,const double &m0, const double &Xs,const double &Ys,const double &Zs,const double &Phi,const double &Omega,const double &Kappa,const double *R){FILE *fp_out;if (!(fp_out=fopen("结果.txt","w"))){return 1;}FILE *fp_input;if (!(fp_input=fopen("data.txt","r"))){return 1;}fprintf(fp_out,"**************************************""**************************************""**************************************""*********************************\n");fprintf(fp_out,"\n空间后方交会程序(C\\C++)\n测绘一班\n""学号：20080729\n姓名：陈闻亚\n\n");fprintf(fp_out,"**************************************""**************************************""**************************************""*********************************\n");fprintf(fp_out,"已知数据：\n\n已知点数：");int num;double temp,x,y;fscanf(fp_input,"%d%lf",&num,&temp);fprintf(fp_out,"%d\n",num);fprintf(fp_out,"摄影比例尺（0表示其值位置）：");fprintf(fp_out,"%10.0lf\n",temp);fprintf(fp_out,"内方位元素(f x0 y0)：");fscanf(fp_input,"%lf%lf%lf",&temp,&x,&y);fprintf(fp_out,"%10lf\t%10lf\t%10lf\n",temp,x,y);for (int i=0;i<num;i++){double temp[5];fscanf(fp_input,"%*d%lf%lf%lf%lf%lf",&temp[0],&temp[1],&temp[2],&temp[3],&temp[4]);fprintf(fp_out,"%3d\t%10lf\t%10lf\t%10lf\t%10lf\t%10lf\n", i+1,temp[0],temp[1],temp[2],temp[3],temp[4]);}fclose(fp_input);fprintf(fp_out,"**************************************""**************************************""**************************************""*********************************\n");fprintf(fp_out,"计算结果如下：\n\n外方位元素：\n"); fprintf(fp_out,"\tXs=%10lf\n",Xs);fprintf(fp_out,"\tYs=%10lf\n",Ys);fprintf(fp_out,"\tZs=%10lf\n",Zs);fprintf(fp_out,"\tPhi=%10lf\n",Phi);fprintf(fp_out,"\tOmega=%10lf\n",Omega);fprintf(fp_out,"\tKappa=%10lf\n\n",Kappa);fprintf(fp_out,"旋转矩阵：\n");for (i=0;i<3;i++){fprintf(fp_out,"\t");for (int j=0;j<3;j++){fprintf(fp_out,"%10lf\t",*(R+i*3+j));}fprintf(fp_out,"\n");}fprintf(fp_out,"\n单位权中误差：%10lf\n\n",m0); fprintf(fp_out,"协因数阵：\n");for (i=0;i<6;i++){fprintf(fp_out,"\t");for (int j=0;j<6;j++){fprintf(fp_out,"%20lf\t",*(Q+i*6+j));}fprintf(fp_out,"\n");}fprintf(fp_out,"\n外方位元素精度：");for (i=0;i<6;i++){fprintf(fp_out,"%10lf\t",m[i]);}fprintf(fp_out,"\n");fprintf(fp_out,"**************************************""**************************************""**************************************""*********************************\n");fclose(fp_out);return 0;}/***********************************************函数名：Resection*函数介绍：计算*参数：Num(点数)，Data(数据)，m,,f(焦距),x0,y0*返回值：int，0成功，其它失败**********************************************/int Resection(const int &Num,const DOUBLE *&Data,const double &m, const double &f){double Xs=0,Ys=0,Zs=0;int i,j;//设置初始值for (i=0;i<Num;i++){Xs+=Data[i][2];Ys+=Data[i][3];}Xs/=Num;Ys/=Num;Zs=m*f;double Phi(0),Omega(0),Kappa(0);double R[3][3]={0.0};double *L=new double[2*Num];typedef double Double6[6]; Double6 *A=new Double6[2*Num]; double *AT=new double[2*Num*6]; double *ATA=new double[6*6]; double *ATL=new double[6]; double *Xg=new double[6];//迭代计算do{//旋转矩阵R[0][0]=cos(Phi)*cos(Kappa)-sin(Phi)*sin(Omega)*sin(Kappa);R[0][1]=-cos(Phi)*sin(Kappa)-sin(Phi)*sin(Omega)*cos(Kappa);R[0][2]=-sin(Phi)*cos(Omega);R[1][0]=cos(Omega)*sin(Kappa);R[1][1]=cos(Omega)*cos(Kappa);R[1][2]=-sin(Omega);R[2][0]=sin(Phi)*cos(Kappa)+cos(Phi)*sin(Omega)*sin(Kappa);R[2][1]=-sin(Phi)*sin(Kappa)+cos(Phi)*sin(Omega)*cos(Kappa);R[2][2]=cos(Phi)*cos(Omega);for (i=0;i<Num;i++){double X=R[0][0]*(Data[i][2]-Xs)+R[1][0]*(Data[i][3]-Ys)+R[2][0]*(Data[i][4]-Zs);double Y=R[0][1]*(Data[i][2]-Xs)+R[1][1]*(Data[i][3]-Ys)+R[2][1]*(Data[i][4]-Zs);double Z=R[0][2]*(Data[i][2]-Xs)+R[1][2]*(Data[i][3]-Ys)+R[2][2]*(Data[i][4]-Zs);double xxx,yyy;xxx=-f*X/Z;yyy=-f*Y/Z;//常数项L[2*i]=Data[i][0]-(-f*X/Z);L[2*i+1]=Data[i][1]-(-f*Y/Z);A[2*i][0]=(R[0][0]*f+R[0][2]*(xxx))/Z;A[2*i][1]=(R[1][0]*f+R[1][2]*(xxx))/Z;A[2*i][2]=(R[2][0]*f+R[2][2]*(xxx))/Z;A[2*i][3]=(yyy)*sin(Omega)-(((xxx)/f)* ((xxx)*cos(Kappa)-(yyy)*sin(Kappa))+f*cos(Kappa))*cos(Omega);A[2*i][4]=-f*sin(Kappa)-((xxx)/f)*((xxx)* sin(Kappa)+(yyy)*cos(Kappa));A[2*i][5]=(yyy);A[2*i+1][0]=(R[0][1]*f+R[0][2]*(yyy))/Z; A[2*i+1][1]=(R[1][1]*f+R[1][2]*(yyy))/Z; A[2*i+1][2]=(R[2][1]*f+R[2][2]*(yyy))/Z; A[2*i+1][3]=-(xxx)*sin(Omega)-(((yyy)/f)*((xxx)*cos(Kappa)-(yyy)*sin(Kappa))-f*sin(Kappa))*cos(Omega);A[2*i+1][4]=-f*cos(Kappa)-((yyy)/f)*((xxx)* sin(Kappa)+(yyy)*cos(Kappa));A[2*i+1][5]=-(xxx);}//求矩阵A的转置矩阵ATfor (i=0;i<2*Num;i++){for (j=0;j<6;j++){*(AT+j*2*Num+i)=A[i][j];}}//求ATAMatrixMul(AT,6,&A[0][0],6,2*Num,ATA);if(InverseMatrix(ATA,6))return 1;MatrixMul(AT,6,L,1,2*Num,ATL);MatrixMul(ATA,6,ATL,1,6,Xg);Xs+=Xg[0];Ys+=Xg[1];Zs+=Xg[2];Phi+=Xg[3];Omega+=Xg[4];Kappa+=Xg[5];} while(fabs(Xg[0])>=PRECISION ||fabs(Xg[1])>=PRECISION || fabs(Xg[2])>=PRECISION ||fabs(Xg[3])>=PRECISION ||fabs(Xg[4])>=PRECISION || (Xg[5])>=PRECISION);//注：协因数阵，旋转矩阵等计算本应该使用最后外方位元素值，//由于变换很小忽略double *Q=ATA;double *V=new double[2*Num];MatrixMul(&A[0][0],2*Num,Xg,1,6,V);double VTV=0;for(i=0;i<2*Num;i++){V[i]-=L[i];VTV+=V[i]*V[i];}double m0=sqrt(VTV/(2*Num-6));double *mm=new double[6];for (i=0;i<6;i++){mm[i]=sqrt(*(Q+i*6+i))*m0;}OutPut(Q,mm,m0,Xs,Ys,Zs,Phi,Omega,Kappa,&R[0][0]);delete []L;delete []A;delete []AT;delete []ATA;delete []ATL;delete []Xg;delete []mm;delete []V;return 0;}void swap(double &a,double &b){double temp=a;a=b;b=temp;}/***********************************************函数名：InverseMatrix*函数介绍：求矩阵的逆（高斯-约当法）*输入参数：(in/out)matrix（矩阵首地址），*（in）row（矩阵阶数）*输出参数：matrix（原矩阵的逆矩阵）*返回值：int ，0成功，1失败*调用函数：swap(double&,double&)**********************************************/int InverseMatrix(double *matrix,const int &row) {double *m=new double[row*row];double *ptemp,*pt=m;int i,j;ptemp=matrix;for (i=0;i<row;i++){for (j=0;j<row;j++){*pt=*ptemp;ptemp++;pt++;}}int k;int *is=new int[row],*js=new int[row];for (k=0;k<row;k++){double max=0;//全选主元//寻找最大元素for (i=k;i<row;i++){for (j=k;j<row;j++){if (fabs(*(m+i*row+j))>max) {max=*(m+i*row+j);is[k]=i;js[k]=j;}}}if (0 == max){return 1;}//行交换if (is[k]!=k){for (i=0;i<row;i++){swap(*(m+k*row+i),*(m+is[k]*row+i)); }}//列交换if (js[k]!=k){for (i=0;i<row;i++){swap(*(m+i*row+k),*(m+i*row+js[k])); }}*(m+k*row+k)=1/(*(m+k*row+k));for (j=0;j<row;j++){if (j!=k){*(m+k*row+j)*=*((m+k*row+k));}}for (i=0;i<row;i++){if (i!=k){for (j=0;j<row;j++){if(j!=k){*(m+i*row+j)-=*(m+i*row+k)**(m+k*row+j);}}}}for (i=0;i<row;i++){if(i!=k){*(m+i*row+k)*=-(*(m+k*row+k));}}int r;//恢复行列for (r=row-1;r>=0;r--){if (js[r]!=r){for (j=0;j<row;j++){swap(*(m+r*row+j),*(m+js[r]*row+j));}}if (is[r]!=r){for (i=0;i<row;i++){swap(*(m+i*row+r),*(m+i*row+is[r]));}}}ptemp=matrix;pt=m;for (i=0;i<row;i++) {for (j=0;j<row;j++){*ptemp=*pt;ptemp++;pt++;}}delete []is;delete []js;delete []m;return 0;}6 计算结果外方位元素：Xs=39795.452297Ys=27476.462210Zs=7572.685927φ= - 0.003987ω= 0.002114κ= - 0.067578旋转矩阵：0.997709 0.067534 0.003987-0.067526 0.997715 -0.002114-0.004121 0.001840 0.999990单位权中误差：0.0000077 心得体会与建议有挑战才会有激情，成功了便更有满足感。