近景单张像片空间后方交会
【武汉大学-摄影测量学-单张相片解析】3.5.5单片空间后方交会
cos
0
s
in
0 0 1 0 0 0
1 0 0
X
YZ
R 1 R
R1
X
Y
Z
Xs
Ys Zs
0 R1 0
1
0 0 0
1X X s
0 0
Y Z
Ys Zs
武汉大学
摄影测量基础
偏导数－2-2
X
YZ
R 1
0 0
1
0 0 0
1 X
0 0
R
YZ
c1
X s
X)
f Z 2 (a1Z a3 X )
1
X
Z
(a1 f
f
Z
a3 )
1 Z
a1 f
a3 (x
x0 )
武汉大学
摄影测量基础
偏导数－1
x X s
1 Z
a1 f
a3(x x0 )
x Ys
1 Z
b1 f
b3(x x0 )
x Z s
1 Z
c1 f
c3(x x0 )
y X s
1 Z
c2 c3
0 0 0
a1 a1
a2 a3
bc11
a2 b2 c2
a3 b3 c3
X Y Z
0
aa23cc11
a1c2 a1c3
a1c2 a2c1 0
a3c2 a2c3
a1c3 a2c3
0
a3c1 a3c2
X
YZ
0 bb32
b3 0
b1
b2 b1
0
X
武汉大学
摄影测量基础
误差方程的建立
☺ 已知值 x0 , y0 , f , m, X, Y, Z ☺ 观测值 x, y
单像空间后方交会名词解释
单像空间后方交会名词解释
单像空间后方交会是摄影测量学中的一个重要概念，它是指利用单个影像进行地物测量和定位的方法。
在单像空间后方交会中，通过对单张影像进行分析，可以确定地面上物体的位置和形状。
这个过程涉及到对影像中的特征点进行识别和匹配，然后利用相机内外参数以及影像上的像点坐标来计算地物的三维坐标。
单像空间后方交会的过程包括以下几个步骤，首先是对影像进行预处理，包括去畸变、影像配准等操作；然后是特征点的提取和匹配，这一步是通过计算机视觉算法来实现的，可以利用角点、边缘等特征来进行匹配；接下来是相机内外参数的标定，这一步是为了将像素坐标转换为实际世界坐标而进行的；最后是利用已知的相机参数和像点坐标来计算地物的三维坐标。
单像空间后方交会在航空摄影、遥感影像解译和地图制图等领域有着广泛的应用。
它可以通过对单张影像的处理，实现对地物的测量和定位，为地理信息系统和地图制图提供了重要的数据基础。
同时，随着计算机视觉和图像处理技术的不断发展，单像空间后方交会的精度和效率也在不断提高，为各种应用领域提供了更加可靠和精确的地物信息。
第五讲 单片空间后方交会
x12 − f (1 + 2 ) f xy − 1 1 f
2 x2 − f (1 + 2 ) f
−
x1 y1 f
y12 − f (1 + 2 ) f − x2 y2 f
x y − 2 2 f
2 x3 − f (1 + 2 ) f
2 y2 − f (1 + 2 ) f
−
x3 y3 f
xy − 3 3 f
Y B
A
C X
利用航摄像片上三个以上像点坐标和对应像 点坐标和对应地面点坐标，计算像片外方位元 素的工作，称为单张像片的空间后方交会。 进行空间后方交会运算，常用的一个基本公 式是前面提到的共线方程。式中的未知数，是 六个外方位元素。由于一个已知点可列出两个 方程式，如有三个不在一条直线上的已知点， 就可列出六个独立的方程式，解求六个外方位 元素。由于共线条件方程的严密关系式是非线 性函数，不便于计算机迭代计算。为此，要由 严密公式推导出一次项近似公式，即变为线性 函数。
（5） 用所取未知数的初始值和控制点的地面坐标，代入共线方程式，逐 ） 用所取未知数的初始值和控制点的地面坐标，代入共线方程式， 点计算像点坐标的近似值 ( x), ( y ) 并计算 lx , l y a ( X − X S ) + b1 (Y − YS ) + c1 ( Z − Z S ) x=−f 1 a3 ( X − X S ) + b3 (Y − YS ) + c3 ( Z − Z S ) a ( X − X S ) + b2 (Y − YS ) + c2 ( Z − Z S ) y=−f 2 a3 ( X − X S ) + b3 (Y − YS ) + c3 ( Z − Z S ) （6） 组成误差方程式。 ） 组成误差方程式。 7） 计算法方程式的系数矩阵与常数项，组成法方程式。 （7） 计算法方程式的系数矩阵与常数项，组成法方程式。 （8） 解算法方程，迭代求得未知数的改正数。 ） 解算法方程，迭代求得未知数的改正数。
单张相片的空间后方交会在两部固定相机相对定位中的应用
计算。直到满足 要求为止． 三、 计算结 果
用上述 方法编写 vc++程序， ．求得了算 例中已知 的一个固定 相机的六 个外 方位元 素为
以=-302．9879mm乓=一L3426mm磊=3．9376mm尹=91．4045。
[ 关键词】后方交会相对定位坐标转换 中图分类号：0 29 文献标识码：^ 文章编 号： 16 71- - 759 7( 20 08) 1220 107 一01
【x，Lzr =剐工，乃z】1
0
0
R=撇 =雕 习 1 1刭 R=月 肭 =I 1 0 I‘ 0
’l
cos r
A如 屯I
【s i I l 矿c∞尹J 0 si nw cos oJ J 【 0 0 J 【q c2 c，J
降sitnx"=0。雕I=100
弘- ，畿XA蔫Xs篇盛嬲慨- ， 。码( 一 )+岛(匕一b)+c3(乙一z；)
‘J ’l ’
，，：～，．垒! 当 二墨! ±垒! 匕二垦2±生! 幺=圣2 (3．2，
7
。码(爿_一五) +岛(匕一乓)确定 将( 3) 式线性化得到误差方程为
∞=0．2752’J r =- 2．25 29。
( 下转第155页)
囫
；i - ●
历 史 教学 中 合 理运 用 “讨论法”
刘敬顺 ( 山东省昌邑市北孟二中 山东昌邑261318) 中圈分类号：G42 文献标识码z ^ 文章 编号 ：1671 - - 759 7( 200 8) 1220 155一 01
按上述过程分别求解得到两个相机的外方位元素即12 个元素。完成了在 物坐标 系下，两部 固定相机相 对位置的确 定．
单像空间后方交会实习报告
单像空间后方交会实习报告一、实习目的单像空间后方交会是摄影测量中确定像片外方位元素的重要方法。
通过本次实习，旨在深入理解单像空间后方交会的基本原理和计算过程，熟练掌握相关软件的操作，提高对摄影测量数据处理的实践能力，并培养解决实际问题的思维和方法。
二、实习原理单像空间后方交会的目的是利用像片上的像点坐标以及相应的地面控制点坐标，通过数学模型求解像片的外方位元素（三个线元素 Xs、Ys、Zs 和三个角元素φ、ω、κ）。
其基本原理基于共线条件方程，即摄影中心、像点和相应的地面点位于同一条直线上。
共线条件方程可以表示为：＼＼begin{align}x x_0&＝ f\frac{a_1(X X_s) ＋ b_1(Y Y_s) ＋ c_1(Z Z_s)｝｛a_3(X X_s) ＋ b_3(Y Y_s) ＋ c_3(Z Z_s)｝＼＼y y_0&＝ f\frac{a_2(X X_s) ＋ b_2(Y Y_s) ＋ c_2(Z Z_s)｝｛a_3(X X_s) ＋ b_3(Y Y_s) ＋ c_3(Z Z_s)｝＼end{align}＼其中，＼（（x,y)＼）为像点坐标，＼（（x_0,y_0)＼）为主点坐标，＼（f\）为摄影机焦距，＼（（X,Y,Z)＼）为地面点的物方空间坐标，＼（（X_s,Y_s,Z_s)＼）为摄影中心的物方空间坐标，＼（（a_1,b_1,c_1)，（a_2,b_2,c_2)，（a_3,b_3,c_3)＼）为由角元素φ、ω、κ 构成的旋转矩阵的元素。
三、实习数据本次实习使用了一组航空像片，像片比例尺为 1:5000，焦距为152mm，像主点坐标为\（（x_0,y_0)＝（5000mm,5000mm)＼）。
同时，提供了 6 个均匀分布在像片范围内的地面控制点的物方空间坐标和像点坐标。
四、实习步骤1、数据准备整理地面控制点的物方空间坐标和像点坐标，确保数据的准确性。
输入像片的基本参数，如像主点坐标、焦距等。
单像空间后方交会
单像空间后方交会测绘学院 成晓倩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 的函数。
第2章后方交会
则
f sin H f a22 cos H y y0 a23 H ( x x0 )( y y0 ) ( y y0 ) 2 a24 cos ( f )sin f f a21 ( y y0 ) 2 ( x x0 )( y y0 ) a25 ( f ) cos sin f f
sin sin 0 0 cos cos
cos 0 sin
X Y 1 X X s 0 0 1 X X s Z R R 1 R 1 Y Y R 0 0 0 Y Y s s Z Z s 1 0 0 Z Z s
《摄影测量学》
第2 章
单片空间后方交会
主要内容
一、定义
二、误差方程和法方程 三、计算过程
一、定义 z
y s(Xs, Ys, Zs) Z
b a
x
c
根据影像覆盖范 围内一定数量的 分布合理的地面 控制点（已知其 像点和地面点的 坐标），利用共 线条件方程求解 像片外方位元素 C X
Y
B
A
注：如果我们有每张像片的六个外方位元素，就能恢复航摄像片 与被摄地面之间的几何关 R 1 Y Y s Z Z Z s
0 a 2 c1 a1c2 a c a c 3 1 1 3 0 b3 b 2 b3 0 b1
a1c2 a 2 c1 0 a 3 c 2 a 2 c3 b2 X b1 Y 0 Z
x f X X Z ( ) Z Z f b2 Z b3Y Z Y b2 f b3 f f Z X (b1Y b2 X ) Z X Y X (b1 b2 ) Z Z Z
单张像片的空间后方交会原理
单像空间后方交会原理非常简单，稍稍懂一点线性代数和微分学知识就很容易理解。
但它的意义很大，所有的基于视觉技术的定位（测量）都是从它扩展来的，没有它就没有后面的复杂算法的演进，以面向更复杂的场景以及需求。
这块内容一般来说就是在已知地面上若干点的地面坐标以后，反求该相应摄影光束的外参数，当用以作摄影机的标定时，还可借以同时求出摄影的内参数。
具体的说就是一个多元线性模型，根据不同的需求可以基于它进一步构建各种平差模型。
摄影测量学空间后方交会实验报告
摄影测量学实验报告实验一、单像空间后方交会学院：建测学院班级：测绘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．实验结果四．实验总结此次实验让我深入了解单像空间后方交会的计算过程，加强了对空间后方交会基本公式和误差方程式，法线方程式的记忆。
《摄影测量学》第10讲-空间后方交会
0 0 Fx ( X S ,YS0 , Z S ,ϕ 0 , ω0 ,κ 0 ) → Fx0
0 (XS − XS ) +
(YS − YS0 ) +
0 (Z S − Z S ) +
(ω − ω0 ) +
∂Fx0 ∂κ
0 0 (κ − κ 0 ) + Fx ( X S , YS0 , Z S ,ϕ 0 , ω0 ,κ 0 )
内 容 安 排
• 单像空间后方交会概述 • 共线方程的线性化（难点） 共线方程的线性化（难点） • 利用共线条件方程解算像片的外方位元 点） （ 点）
[一]概述
1、什么叫单像空间后方交会 什么叫单像空间后方交会 利用地面控制点及其在片像上的像点， 利用地面控制点及其在片像上的像点，确定一 张像片外方位元素的方法。 张像片外方位元素的方法。
2
(
)
求：a = ?
取初值
任取a0＝0： da = 6 由于：da = a − a0，a = a0 + da = 6 da = −36 / 13 = −2.8 取a0＝6： 由于：da = a − a0，a = a0 + da = 3.2 da = −1 取a0＝3.2 由于：da = a − a0，a = a0 + da = 2.2
S S
) + b2 ( Y − Y S ) + c 2 ( Z − Z S ) ) + b3 ( Y − Y S ) + c 3 ( Z − Z S )
a1 ( X − X S ) + b1 (Y − YS ) + c1 ( Z − Z S ) Fx = x + f =0 a3 ( X − X S ) + b3 (Y − YS ) + c3 ( Z − Z S ) Fy = y + f a 2 ( X − X S ) + b2 (Y − YS ) + c 2 ( Z − Z S ) =0 a3 ( X − X S ) + b3 (Y − YS ) + c3 ( Z − Z S )
- 1、下载文档前请自行甄别文档内容的完整性，平台不提供额外的编辑、内容补充、找答案等附加服务。
- 2、"仅部分预览"的文档,不可在线预览部分如存在完整性等问题,可反馈申请退款(可完整预览的文档不适用该条件!)。
- 3、如文档侵犯您的权益，请联系客服反馈,我们会尽快为您处理(人工客服工作时间：9:00-18:30)。
实验一近景单张像片空间后方交会一、实验目的通过单张像片空间后方交会算法编程，掌握空间后方交会的基本原理和基本步骤，为近景摄影测量解析处理方法打下理论基础。
二、实验内容利用C++编程平台，通过给定的单片像点的像点坐标、内方位元素和控制点物方空间坐标，计算出像片的外方位元素。
三、实验步骤：1、编程流程图：2、编程代码：#include "iostream"#include"stdio.h"#include "stdlib.h"#include<math.h>#define N 10using namespace std;void mult(double *m1,double *m2,double *result,int i_1,int j_12,int j_2)//矩阵相乘{int i,j,k;for(i=0;i<i_1;i++)for(j=0;j<j_2;j++){result[i*j_2+j]=0.0;for(k=0;k<j_12;k++)result[i*j_2+j]+=m1[i*j_12+k]*m2[j+k*j_2];}return;}/* ######################矩阵求逆##################*/int invers_matrix(double *m1,int n){int *is,*js;int i,j,k,l,u,v;double temp,max_v;is=(int *)malloc(n*sizeof(int));js=(int *)malloc(n*sizeof(int));if(is==NULL||js==NULL){printf("out of memory!\n");return(0);}for(k=0;k<n;k++){max_v=0.0;for(i=k;i<n;i++)for(j=k;j<n;j++){temp=fabs(m1[i*n+j]);if(temp>max_v){max_v=temp; is[k]=i; js[k]=j;}}if(max_v==0.0){free(is); free(js);printf("invers is not availble!\n");return(0);}if(is[k]!=k)for(j=0;j<n;j++){u=k*n+j; v=is[k]*n+j;temp=m1[u]; m1[u]=m1[v]; m1[v]=temp;}if(js[k]!=k)for(i=0;i<n;i++){u=i*n+k; v=i*n+js[k];temp=m1[u]; m1[u]=m1[v]; m1[v]=temp;}l=k*n+k;m1[l]=1.0/m1[l];for(j=0;j<n;j++)if(j!=k){u=k*n+j;m1[u]*=m1[l];}for(i=0;i<n;i++)if(i!=k)for(j=0;j<n;j++)if(j!=k){u=i*n+j;m1[u]-=m1[i*n+k]*m1[k*n+j];}for(i=0;i<n;i++)if(i!=k){u=i*n+k;m1[u]*=-m1[l];}}for(k=n-1;k>=0;k--){if(js[k]!=k)for(j=0;j<n;j++){u=k*n+j; v=js[k]*n+j;temp=m1[u]; m1[u]=m1[v]; m1[v]=temp;}if(is[k]!=k)for(i=0;i<n;i++){u=i*n+k; v=i*n+is[k];temp=m1[u]; m1[u]=m1[v]; m1[v]=temp;}}free(is); free(js);return(1);}/*------------------------------矩阵转置----------------------------------*/void transpose(double *m1,double *m2,int m,int n){int i,j;for(i=0;i<m;i++)for(j=0;j<n;j++)m2[j*m+i]=m1[i*n+j];return;}void main(){double Xs,Ys,Zs,q,w,k;//外方位元素//double a[3],b[3],c[3];//旋转矩阵//double x0,y0,f;//内方位元素//double x[N],y[N];double X[N],Y[N],Z[N];double x1[N],y1[N];double m;//中误差//double L[2*N];double XX[6];double A[2*N][6];double X0[N],Y0[N],Z0[N],At[6][2*N],result1[6][6],result2[6][1];//外方位元素初始置//int i,n=0;//迭代次数//double sum=0,m0;/*for(i=0;i<N;i++){printf("请输入第%d个点的地面坐标：",i+1); //----------输入点地面坐标-scanf_s("%lf%lf%lf",&X[i],&Y[i],&Z[i]);} */X[0]=1617.712 ;Y[0]=1883.870 ;Z[0]=-4917.308;X[1]=1615.572 ;Y[1]= 1687.514 ; Z[1]= -4913.750;X[2]=1614.640 ;Y[2]= 1315.066 ;Z[2]=-4916.724;X[3]=1615.784 ;Y[3]= 1001.358 ;Z[3]= -4916.691;X[4]=1617.197 ;Y[4]= 690.040 ;Z[4]= -4917.070;X[5]=1615.933 ;Y[5]= 387.492;Z[5]= -4917.234;X[6]=1618.946 ; Y[6]= 145.491;Z[6]= -4917.807;X[7]=1616.753 ;Y[7]= -186.640 ;Z[7]= -4918.531;X[8]=1617.714 ;Y[8]= -482.264;Z[8]= -4918.561;X[9]=1617.859 ;Y[9]= -680.767;Z[9]= -4918.799;/*for(i=0;i<N;i++){printf("请输入第%d个点的像片坐标：",i+1); //---------输入点像片坐标scanf_s("%lf%lf",&x[i],&y[i]);}cout<<endl;*/x[0]=1105.312313;y[0]=3367.979242;x[1]=1111.483438;y[1]=3203.998009;x[2]=1128.936009;y[2]=2886.146374;x[3]=1145.729428;y[3]=2615.385301;x[4]=1163.579521;y[4]=2343.603584;x[5]=1179.483899;y[5]=2076.672738;x[6]=1196.664025;y[6]=1861.971885;x[7]=1216.048311;y[7]=1565.046087;x[8]=1236.196347;y[8]=1299.275413;x[9]=1249.998721;y[9]=1120.52425;/*-----------------设定外方位元素初始值--------------*/x0=48.6401;y0=11.4765;f=4548.5949;Xs=3391.658;Ys=-135.645;Zs=97.250;q=-0.0102;w=0.0620;k=-0.0887;XX[3]=1;/*------------------迭代计算--------------------------*/while((XX[3]>0.00001 || XX[4]>0.00001 || XX[5]>0.00001)&&n<100) {/*----------------旋转矩阵R-----------------------*/a[0]=cos(q)*cos(k)-sin(q)*sin(w)*sin(k);a[1]=-cos(q)*sin(k)-sin(q)*sin(w)*cos(k);a[2]=-sin(q)*cos(w);b[0]=cos(w)*sin(k);b[1]=cos(w)*cos(k);b[2]=-sin(w);c[0]=sin(q)*cos(k)+cos(q)*sin(w)*sin(k);c[1]=-sin(q)*sin(k)+cos(q)*sin(w)*cos(k);c[2]=cos(q)*cos(w);/*-----------------像点坐标计算值------------------*/for(i=0;i<N;i++){X0[i]=a[0]*(X[i]-Xs)+b[0]*(Y[i]-Ys)+c[0]*(Z[i]-Zs);Y0[i]=a[1]*(X[i]-Xs)+b[1]*(Y[i]-Ys)+c[1]*(Z[i]-Zs);Z0[i]=a[2]*(X[i]-Xs)+b[2]*(Y[i]-Ys)+c[2]*(Z[i]-Zs);x1[i]=x0-f*X0[i]/Z0[i];y1[i]=y0-f*Y0[i]/Z0[i];}/*-------------误差方程中各偏导数的值--------------*/for(i=0;i<N;i++){A[2*i][0]=((a[0]*f+a[2]*(x[i]-x0)))/Z0[i];A[2*i][1]=((b[0]*f+b[2]*(x[i]-x0)))/Z0[i];A[2*i][2]=((c[0]*f+c[2]*(x[i]-x0)))/Z0[i];A[2*i][3]=(y[i]-y0)*sin(w)-((x[i]-x0)*((x[i]-x0)*cos(k)-y[i]*sin( k))/f+f*cos(k))*cos(w);A[2*i][4]=-f*sin(k)-(x[i]-x0)*((x[i]-x0)*sin(k)+(y[i]-y0)*cos(k))/f;A[2*i][5]=y[i]-y0;L[2*i]=x[i]-x1[i];A[1+2*i][0]=((a[1]*f+a[2]*(y[i]-y0)))/Z0[i];A[1+2*i][1]=((b[1]*f+b[2]*(y[i]-y0)))/Z0[i];A[1+2*i][2]=((c[1]*f+c[2]*(y[i]-y0)))/Z0[i];A[1+2*i][3]=-(x[i]-x0)*sin(w)-((y[i]-y0)*((x[i]-x0)*cos(k)-(y[i]-y0)*sin(k))/f-f*sin(k))*cos(w);A[1+2*i][4]=-f*cos(k)-(y[i]-y0)*((x[i]-x0)*sin(k)+(y[i]-y0)*cos(k ))/f;A[1+2*i][5]=-x[i]+x0;L[1+2*i]=y[i]-y1[i];}/*-------------------解法方程--------------------*/transpose(&A[0][0],&At[0][0],2*N,6);mult(&At[0][0],&A[0][0],&result1[0][0],6,2*N,6);invers_matrix(&result1[0][0],6);mult(&At[0][0],L,&result2[0][0],6,2*N,1);mult(&result1[0][0],&result2[0][0],&XX[0],6,6,1);Xs+=XX[0];Ys+=XX[1];Zs+=XX[2];q+=XX[3];w+=XX[4];k+=XX[5];n++;}/*----------------旋转矩阵R-----------------------*/a[0]=cos(q)*cos(k)-sin(q)*sin(w)*sin(k);a[1]=-cos(q)*sin(k)-sin(q)*sin(w)*cos(k);a[2]=-sin(q)*cos(w);b[0]=cos(w)*sin(k);b[1]=cos(w)*cos(k);b[2]=-sin(w);c[0]=sin(q)*cos(k)+cos(q)*sin(w)*sin(k);c[1]=-sin(q)*sin(k)+cos(q)*sin(w)*cos(k);c[2]=cos(q)*cos(w);cout<<"迭代次数为:"<<n<<endl;printf("\n像片外方位元素的解\n");cout<<"航向顷角q:"<<q<<endl;cout<<"旁向倾角w:"<<w<<endl;cout<<"像片旋角k:"<<k<<endl;printf("Xs=%lf\t\tYs=%lf\t\tZs=%lf\n",Xs,Ys,Zs);cout<<endl;printf("旋转矩阵R:\n");for(i=0;i<3;i++)printf("%lf\t",a[i]);printf("\n");for(i=0;i<3;i++)printf("%lf\t",b[i]);printf("\n");for(i=0;i<3;i++)printf("%lf\t",c[i]);printf("\n");/*-------------------计算单位权中误差---------------*/ for(i=0;i<2*N;i++)sum+=L[i]*L[i];m0=sqrt(sum/(2*N-6));cout<<"单位权中误差m0="<<m0<<endl;cout<<"测量精度:"<<endl;cout<<"Xs的精度为:"<<m0*sqrt(result1[0][0])<<endl;cout<<"Ys的精度为:"<<m0*sqrt(result1[1][1])<<endl;cout<<"Zs的精度为:"<<m0*sqrt(result1[2][2])<<endl;cout<<"q的精度为:"<<m0*sqrt(result1[3][3])<<endl;cout<<"w的精度为:"<<m0*sqrt(result1[4][4])<<endl;cout<<"k的精度为:"<<m0*sqrt(result1[5][5])<<endl;system("pause");}四、实验结果分析五、实验心得通过这次实验，使我掌握了在C++及Microsoft visual studio 环境下单张像片空间后方交会算法编程及实现，掌握空间后方交会的基本原理和基本步骤，并且锻炼了我编程的能力，看到自己编出的程序成功实现后，感到很有成就感，这也是我有了更大的自信，也为以后的工作打下了更好的基础，感谢老师的这次实习。