单像空间后方交会程序《Java版》
本科实验报告
课程名称：单像空间后方交会编程实现姓名：李恒松
学院：地学学院
专业：测绘工程一班
学号：20133294
指导老师：陈强
2015年5 月20 日
目录
一、作业任务 ............................................................................................................... - 2 -
二、计算原理 ............................................................................................................... - 2 -
三、算法流程 ............................................................................................................... - 6 -
四、源程序 ................................................................................................................... - 7 -
五、计算结果 ............................................................................................................... - 7 -
六、结果分析 ............................................................................................................... - 7 -
七、心得与体会 ........................................................................................................... - 7 -
八、附页 ....................................................................................................................... - 7 -
Java程序 ................................................................................................................................................................ - 7 -
一、
作业任务
已知条件：
摄影机主距f=153.24mm ，x0=0，y0=0, 像片比例尺为1:40000，有四对点的像点坐标与相应的地面坐标如下表。
以单像空间后方交会方法，求解该像片的外方位元素。
二、
计算原理
1. 获取已知数据。从航摄资料中查取平均航高与摄影机主距；获取控制点的地面测量坐标并转换为地面摄测坐标。
2. 测量控制点的像点坐标并作系统误差改正。
3. 确定未知数的初始值。在竖直摄影且地面控制点大体对称分布的情况下，按如下方法确定初始值，即
1
,,S
S
S
X Y X
Y
Z
mf Z n
n
n
===+
∑∑∑ 0000?ωκ===
式中：m 为摄影比例尺分母；n 为控制点个数
4. 用三个角元素的初始值按下式计算各方向余弦值，组成旋转矩阵R
12
31
231
2
3a a a R b b b c c c ????=??????
其中：
12
3123
123cos cos sin sin sin cos sin sin sin cos sin cos cos sin cos cos sin sin cos cos sin sin sin sin cos sin cos cos cos a a a b b b c c c ?κ?ωκ?κ?ωκ?ωωκωκ
ω?κ?ωκ?κ?ωκ?ω
=-??=--??=-?
=??
=??=-??=+?
=-+??=? 5. 逐点计算像点坐标的近似值。利用未知数的近似值和控制点的地面坐标，带入以下共线方程式，
111333222333()()()
()()()()()()()()()
A S A S A S A S A S A S A S A S A S A S A S A S a X X b Y Y c Z Z x f a X X b Y Y c Z Z a X X b Y Y c Z Z y f a X X b Y Y c Z Z -+-+-?
=-?-+-+-??
-+-+-?=-?-+-+-?
逐点计算像点坐标的近似值()x 、()y ，将题目所给数据代入上式可得像点坐标的近似数据。
6.
逐点计算误差方程式的系数和常数项，组成误差方程式。
??
?
??
?
????????????????????=868584
83
82
81
7675747372716665646362615655545352514645444342413635343332312625242322
211615
14131211
6
*8A a a a a a a a a a a a a a a a a a a a a a a a a a a a a a a a a a a a a a a a a a a a a a a a a
其中
????
????
??
?
???????
?????
?
????
??
?
-≈+-≈-≈-≈-
≈≈≈-
≈+-≈-
≈≈-≈x
a f
y f a f xy a H
y a H f a a y a f xy a f x f a H x a a H
f a 2622
252423
222116152
2
14131211
)1(0)1(0
其他项类似可得。
由常数项计算公式：
111333222333()()()()()()()()()()()()(1,2,3,4)i S i S i S i
xi i S i S i S i yi i S i S i S i
i S i S i S a X X b Y Y c Z Z x f l a X X b Y Y c Z Z L l a X X b Y Y c Z Z y f a X X b Y Y c Z Z i -+-+-?
?+?
?-+-+-???
?==???-+-+-?
??+??-+-+-??
=
将题目所给数据代入上式，可得常数项。
7. 计算法方程的是系数矩阵T A A 和常数项T A L ，组成法方程式。
L AX V -=
8. 解法方程，求得外方位元素的改正数,,,,,S S S dX dY dZ d d d ?ωκ。
L A A A X T T 1)(-=
9. 用前次迭代取得的近似值，加本次迭代的改正数，计算外方位元素的新值。
111111,,,,K K K K K K K K K
S S S S S S S S S K
K K K K K K K K S S
S S S
S S S
S
X X dX Y Y dY Z Z dZ d d d ??
?ωω
ωκκ
κ
------=+=+=+=+=+=+
式中：K 代表迭代次数。
10. 将求得的外方位元素改正数与规定的限差比较，若小于限差，则迭代结束。否则用新的近似值重复4~9满足要求为止。 11. 由误差方程的计算式：
1112131415160
212223242526x s s s y s s s v a X a Y a Z a a a x x
v a X a Y a Z a a a y y
?ωκ?ωκ?=?+?+?+?+?+?+-??=?+?+?+?+?+?+-?? 以及单位权中误差的计算式
0σ=
（式中：r 表示多余观测数） 以及平差中6个参数的协因数阵
最终得到参数i X 的中误差为
i
x σσ=
1
()T xx Q A A -=
源程序代码请见附页
三、计算结果
在经过三次迭代之后得到的最终成果如表1所示
表1，最终计算结果
39795.45225 75952727476.462385
240564
7572.685987
778016
-0.00398692355
9975288
0.00211388806
645731
-0.0675779705
3971112
四、结果分析
由计算结果可知在拍摄照片瞬间，摄影中心在地面摄影测量坐标系中的坐标为（39795.452258，27476.462385，7572.685988）（单位：米），航向倾角ψ为-0.003987弧度，旁向倾角ω为0.002114弧度，像片旋角κ为0.067578弧度。
表2，精度评定结果
参数中误差
Xs 1.15327223063957
Ys 1.301244230747318
Zs 0.5083803266176986
ψ 1.8602199035719645E-4
ω 1.6814981476239593E-4
κ7.501994950982504E-5
五、心得与体会
通过这次大作业，我的心得体会以下几点：
第一：编程确实是一件非诚艰难的事情，即使整个单像空间后方交会实验的原理非常简单，清晰明了，但是要用编程把它表示出来还是一件非常艰难的事情。我遇到的最大的难点是在矩阵求逆的步骤上，为了解决这个问题。我用了两天的时间，还请教了数学学院的同学，程序里求逆矩阵的方法是来自那位同学的。这次编程也让我对于线性代数中的求逆矩阵的问题理解的更加深刻。
第二：在发现问题，解决问题的过程中，我也认识到了自己基础知识的不足，对于Java语言的不熟悉让我在很多常见的bug面前无能为力。以后要加强编程语言的学习。
第三：很多东西看似简单，其实细节处才体现出能力。把小的细节处理好，才能有大的作品。
六、附页
1.Java程序
以下第一部分为主程序
import java.util.*;
import java.math.*;
public class li {
public static void main(String[] args){
heng hengFind=new heng();
double
x[][]={{-0.08615,-0.06899},{-0.05340,0.08221},{-0.01478,-0.07663},{0.01046,0.06443}};
double
X[][]={{36589.41,25273.32,2195.17},{37631.08,31324.51,728.69},{39100.97,24934.98,2386.5 0},{40426.54,30319.81,757.31}};
int i,j,num=0;
double X0[]=new double[6];
double f=0.15324;
double a=1/40000.0;
double R[][]=new double[3][3];
double evaluate_x[]=new double[8];
double A[][]=new double[8][6];
double A_T[][]=new double[6][8];
double L[]=new double[8];
double pi=3.1415926535897932;
double AT_A[][]=new double[6][6];
double inverseAT_A[][]=new double[6][6];//inverseAT_A=hengFind.inverse(AT_A,6);
double result_v[]=new double[8];//用来存放用于计算单位权中误差的v值
double standard[]=new double[6];
double resultF=0;//用来存放v的转置矩阵
double result1[][]=new double[6][8];
double result2[]=new double[6];
double sumXYZ[]=new double[3];
System.out.println("已知像点坐标为： ");
for(i=0;i<4;i++){
for(j=0;j<2;j++)
{
if(j==0)
{
System.out.println("x"+(i+1)+"= "+"\t"+x[i][j]);
}
else
{
System.out.println("y"+(i+1)+"= "+"\t"+x[i][j]);
}
}
}
System.out.println("已知地面四个点的坐标为：");
for(i=0;i<4;i++){
for(j=0;j<3;j++)
{
if(j==0)
{
System.out.println("X"+(i+1)+"= "+"\t"+X[i][j]);
}
else
if(j==1)
{
System.out.println("Y"+(i+1)+"= "+"\t"+X[i][j]);
}
else
{
System.out.println("Z"+(i+1)+"= "+"\t"+X[i][j]);
}
}
}
for(j=0;j<3;j++){
for(i=0;i<4;i++){
sumXYZ[j]+=X[i][j];}
}
for(i=0;i<2;i++){
X0[i]=sumXYZ[i]/4;
}
X0[i]=1/a*f+sumXYZ[2]/4.0;
do{
R[0][0]=Math.cos(X0[3])*Math.cos(X0[5])-Math.sin(X0[3])*Math.sin(X0[4])*Math.sin(X0 [5]);
R[0][1]=-Math.cos(X0[3])*Math.sin(X0[5])-Math.sin(X0[3])*Math.sin(X0[4])*Math.cos(X 0[5]);
R[0][2]=-Math.sin(X0[3])*Math.cos(X0[4]);
R[1][0]=Math.cos(X0[4])*Math.sin(X0[5]);
R[1][1]=Math.cos(X0[4])*Math.cos(X0[5]);
R[1][2]=-Math.sin(X0[4]);
R[2][0]=Math.sin(X0[3])*Math.cos(X0[5])+Math.cos(X0[3])*Math.sin(X0[4])*Math.sin(X0 [5]);
R[2][1]=-Math.sin(X0[3])*Math.sin(X0[5])+Math.cos(X0[3])*Math.sin(X0[4])*Math.cos(X 0[5]);
R[2][2]=Math.cos(X0[3])*Math.cos(X0[4]);
evaluate_x[0]=-f*(R[0][0]*(X[0][0]-X0[0])+R[1][0]*(X[0][1]-X0[1])+R[2][0]*(X[0][2]-X0[2]))/(R[0][2]*(X[0][0]-X0[0])+R[1][2]*(X[0][1]-X0[1])+R[2][2]*(X[0][2]-X0[2]));
evaluate_x[1]=-f*(R[0][1]*(X[0][0]-X0[0])+R[1][1]*(X[0][1]-X0[1])+R[2][1]*(X[0][2]-X0[2]))/(R[0][2]*(X[0][0]-X0[0])+R[1][2]*(X[0][1]-X0[1])+R[2][2]*(X[0][2]-X0[2]));
evaluate_x[2]=-f*(R[0][0]*(X[1][0]-X0[0])+R[1][0]*(X[1][1]-X0[1])+R[2][0]*(X[1][2]-
X0[2]))/(R[0][2]*(X[1][0]-X0[0])+R[1][2]*(X[1][1]-X0[1])+R[2][2]*(X[1][2]-X0[2]));
evaluate_x[3]=-f*(R[0][1]*(X[1][0]-X0[0])+R[1][1]*(X[1][1]-X0[1])+R[2][1]*(X[1][2]-X0[2]))/(R[0][2]*(X[1][0]-X0[0])+R[1][2]*(X[1][1]-X0[1])+R[2][2]*(X[1][2]-X0[2]));
evaluate_x[4]=-f*(R[0][0]*(X[2][0]-X0[0])+R[1][0]*(X[2][1]-X0[1])+R[2][0]*(X[2][2]-X0[2]))/(R[0][2]*(X[2][0]-X0[0])+R[1][2]*(X[2][1]-X0[1])+R[2][2]*(X[2][2]-X0[2]));
evaluate_x[5]=-f*(R[0][1]*(X[2][0]-X0[0])+R[1][1]*(X[2][1]-X0[1])+R[2][1]*(X[2][2]-X0[2]))/(R[0][2]*(X[2][0]-X0[0])+R[1][2]*(X[2][1]-X0[1])+R[2][2]*(X[2][2]-X0[2]));
evaluate_x[6]=-f*(R[0][0]*(X[3][0]-X0[0])+R[1][0]*(X[3][1]-X0[1])+R[2][0]*(X[3][2]-X0[2]))/(R[0][2]*(X[3][0]-X0[0])+R[1][2]*(X[3][1]-X0[1])+R[2][2]*(X[3][2]-X0[2]));
evaluate_x[7]=-f*(R[0][1]*(X[3][0]-X0[0])+R[1][1]*(X[3][1]-X0[1])+R[2][1]*(X[3][2]-X0[2]))/(R[0][2]*(X[3][0]-X0[0])+R[1][2]*(X[3][1]-X0[1])+R[2][2]*(X[3][2]-X0[2]));
for(i=0;i<4;i++)
{
/*a11*/A[2*i][0]=(R[0][0]*f+R[0][2]*evaluate_x[2*i])/(R[0][2]*(X[i][0]-X0[0])+R[1][ 2]*(X[i][1]-X0[1])+R[2][2]*(X[i][2]-X0[2]));
/*a12*/A[2*i][1]=(R[1][0]*f+R[1][2]*evaluate_x[2*i])/(R[0][2]*(X[i][0]-X0[0])+R[1][ 2]*(X[i][1]-X0[1])+R[2][2]*(X[i][2]-X0[2]));
/*a13*/A[2*i][2]=(R[2][0]*f+R[2][2]*evaluate_x[2*i])/(R[0][2]*(X[i][0]-X0[0])+R[1][ 2]*(X[i][1]-X0[1])+R[2][2]*(X[i][2]-X0[2]));
/*a21*/A[2*i+1][0]=(R[0][1]*f+R[0][2]*evaluate_x[2*i+1])/(R[0][2]*(X[i][0]-X0[0])+R [1][2]*(X[i][1]-X0[1])+R[2][2]*(X[i][2]-X0[2]));
/*a22*/A[2*i+1][1]=(R[1][1]*f+R[1][2]*evaluate_x[2*i+1])/(R[0][2]*(X[i][0]-X0[0])+R [1][2]*(X[i][1]-X0[1])+R[2][2]*(X[i][2]-X0[2]));
/*a23*/A[2*i+1][2]=(R[2][1]*f+R[2][2]*evaluate_x[2*i+1])/(R[0][2]*(X[i][0]-X0[0])+R [1][2]*(X[i][1]-X0[1])+R[2][2]*(X[i][2]-X0[2]));
/*a14*/A[2*i][3]=evaluate_x[2*i+1]*Math.sin(X0[4])-(evaluate_x[2*i]/f*(evaluate_x[2 *i]*Math.cos(X0[5])-evaluate_x[2*i+1]*Math.sin(X0[5]))+f*Math.cos(X0[5]))*Math.cos(X0[4 ]);
/*a15*/A[2*i][4]=-f*Math.sin(X0[5])-evaluate_x[2*i]/f*(evaluate_x[2*i]*Math.sin(X0[ 5])+evaluate_x[2*i+1]*Math.cos(X0[5]));
/*a16*/A[2*i][5]=evaluate_x[2*i+1];
/*a24*/A[2*i+1][3]=-1*evaluate_x[2*i]*Math.sin(X0[4])-(evaluate_x[2*i+1]/f*(evaluat e_x[2*i]*Math.cos(X0[5])-evaluate_x[2*i+1]*Math.sin(X0[5]))-f*Math.sin(X0[5]))*Math.cos
(X0[4]);
/*a25*/A[2*i+1][4]=-1*f*Math.cos(X0[5])-evaluate_x[2*i+1]/f*(evaluate_x[2*i]*Math.s in(X0[5])+evaluate_x[2*i+1]*Math.cos(X0[5]));
/*a26*/A[2*i+1][5]=-evaluate_x[2*i];
}
//进行常数项的初始化
for(i=0;i<4;i++)
{
L[2*i]=x[i][0]-evaluate_x[2*i];
L[2*i+1]=x[i][1]-evaluate_x[2*i+1];
}
//A的转置矩阵
for(i=0;i<8;i++)
for(j=0;j<6;j++)
{
A_T[j][i]=A[i][j];
}
//实现A与AT相乘
int k=0;
for(i=0;i<6;i++)
for(j=0;j<6;j++)
AT_A[i][j]=0;
for(i=0;i<6;i++)
for(k=0;k<6;k++)
for(j=0;j<8;j++)
AT_A[i][k]+=A_T[i][j]*A[j][k];
//得到AT*A的逆矩阵存放在inverseAT_A中,此处采用面向对象的方法
inverseAT_A=hengFind.inverse(AT_A,6);
//实现矩阵AT_A[6][6]与AT[6][8]的相乘,结果存放在result1[6][8]中
for(i=0;i<6;i++)
for(j=0;j<8;j++)
result1[i][j]=0;
for(i=0;i<6;i++)
for(k=0;k<8;k++)
for(j=0;j<6;j++)
result1[i][k]+=inverseAT_A[i][j]*A_T[j][k];
//实现result1[6][8]与l[8]的相乘，得到结果放在result2[6]中
for(i=0;i<6;i++)
result2[i]=0;
for(i=0;i<6;i++)
for(j=0;j<8;j++)
result2[i]+=result1[i][j]*L[j];
num++;
for(i=0;i<6;i++)
{
X0[i]=X0[i]+result2[i];
}
System.out.println("\n"+"进行第"+num+"次迭代得到Xs,Ys,Zs, ψ，ω，κ改正数分别为：");
for(i=0;i<6;i++)
{
System.out.print(result2[i]+"\t");
}
}while(Math.abs(result2[3]*206265.0)>6||Math.abs(result2[4]*206265.0)>6||Math.abs(r esult2[5]*206265.0)>6);
System.out.print("\n"+"满足限差条件结束循环，最终结果为："+"\n");
System.out.println("Xs"+"\t"+"Ys"+"\t"+"Zs"+"\t"+"ρ"+"\t"+"ω"+"\t"+"k");
for(i=0;i<6;i++)
{
System.out.print(X0[i]+"\t");
}
// 之前代码所求为外方位元素的六个值
// 以下代码所求为单位权中误差
for(i=0;i<8;i++){
result_v[i]=0;
}
for(i=0;i<8;i++){//这儿可能会出错
for(j=0;j<6;j++){
result_v[i]+=A[i][j]*result2[j];
}
}
for(int q=0;q<6;q++){
result_v[q]-=L[q];
}
for(int m=0;m<8;m++){
resultF+=(result_v[m]*result_v[m]);
}
resultF=Math.sqrt(resultF/2);
for(int p=0;p<6;p++){
standard[p]=(Math.sqrt(inverseAT_A[p][p]))*resultF;
}
System.out.println("\n"+"各所求量的中误差为");
System.out.print("Xs"+"\t"+"Ys"+"\t"+"Zs"+"\t"+"ρ"+"\t"+"ω"+"\t"+"k"+"\n");
System.out.print(standard[0]+"\t"+standard[1]+"\t"+standard[2]+"\t"+standard[3]+"\t "+standard[4]+"\t"+standard[5]+"\n");
}
}
以下第二部分为计算逆矩阵的面向对象程序
public class heng{
public double[][] inverse(double[][] a,int n){
double b[][]=new double[n][2*n];
double c[][]=new double[n][n];
double det1,yinzhi,bb;
int i,j,k,u;
for(i=0;i for(j=0;j<2*n;j++){ b[i][j]=0; } } for(i=0;i for(j=0;j b[i][j]=a[i][j]; } } for(j=0;j b[j][n+j]=1; } for(i=0;i if(b[i][i]==0){ for(j=i;j if(b[j][i]!=0){ this.temp(b[i],b[j],n); } } } for(k=i+1;k yinzhi=(-1)*b[k][i]/b[i][i]; for(u=0;u<2*n;u++){ b[k][u]=b[k][u]+b[i][u]*yinzhi; } } } det1=this.fun(a,n); if(det1==0){ System.out.println("此矩阵不可逆"); } if(det1!=0){ for(i=0;i bb=b[i][i]; for(j=0;j<2*n;j++){ b[i][j]=b[i][j]/bb; } } for(i=n-1;i>0;i--){ for(k=0;k bb=b[k][i]; for(u=0;u<2*n;u++){ b[k][u]=b[k][u]-bb*b[i][u]; } } } } for(i=0;i for(j=0;j c[i][j]=b[i][j+n]; } } return c; } public void temp(double aa[],double bb[],int n ){ int i; double temp1; for(i=0;i temp1=aa[i];aa[i]=bb[i];bb[i]=temp1; } } public double fun(double array[][],int n){ int ii,jj,k,u; double det1=1,yin; for(ii=0;ii if(array[ii][ii]==0){ for(jj=ii;jj if(array[jj][ii]!=0){ temp(array[ii],array[jj],n); } } } for(k=ii+1;k yin=(-1)*array[k][ii]/array[ii][ii]; for(u=0;u array[k][u]=array[k][u]+array[ii][u]*yin; } } for(ii=0;ii det1=det1*array[ii][ii]; } } return(det1); } } 实验程序结果如下图： 空间后方交会编程实习报告 一实习目的 用程序设计语言（Visual C++或者C语言）编写一个完整的单片空间后方交会程序，通过对提供的试验数据进行计算，输出像片的外方位元素并评定精度。本实验的目的在于让学生深入理解单片空间后方交会的原理，体会在有多余观测情况下，用最小二乘平差方法编程实现解求影像外方位元素的过程。通过上机调试程序加强动手能力的培养，通过对实验结果的分析，增强学生综合运用所学知识解决实际问题的能力。 二实习内容 利用一定数量的地面控制点，根据共线条件方程求解像片外方位元素。 三实习数据 已知航摄仪的内方位元素：f k =153.24mm，x =y =0.0mm，摄影比例尺为1:50000； 4个地面控制点的地面坐标及其对应像点的像片坐标： 四实习原理 如果我们知道每幅影像的6个外方位元素，就能确定被摄物体与航摄影像的关系。因此，如何获取影像的外方位元素，一直是摄影测量工作者所探讨的问题。可采取的方法有：利用雷达、全球定位系统（GPS）、惯性导航系统（INS）以及星相摄影机来获取影像的外方位元素；也可以利用影像覆盖范围内一定数量的控制点的空间坐标与摄影坐标，根据共线条件方程，反求该影像的外方位元素，这种方法称为单幅影像的空间后方交会。 单像空间后方交会的基本思想是：以单幅影像为基础，从该影像所覆盖地面范围内若干控制点的已知地面坐标和相应点的像坐标量测值出发，根据共线条件方程，解求该影像在航空摄影时刻的外方位元素Xs，Ys，Zs，t，w，k。 五实习流程 （1）获取已知数据。从摄影资料中查取影像比例尺1/m，平均摄影距离（航空摄影的航高、内方位元素x0，y0，f；获取控制点的空间坐标Xt，Yt，Zt。 （2）量测控制点的像点坐标并进行必要的影像坐标系统误差改正，得到像点坐标。 （3）确定未知数的初始值。单像空间后方交会必须给出待定参数的初始值，在竖直航空摄影且地面控制点大体对称分布的情况下，可按如下方法确定初始值： 空间后方交会的解算 一． 空间后方交会的目的 摄影测量主要利用摄影的方法获取地面的信息，主要是是点位信息，属性信息，因此要对此进行空间定位和建模，并首先确定模型的参数，这就是空间后方交会的目的，用以求出模型外方位元素。 二． 空间后方交会的原理 空间后方交会的原理是共线方程。 共线方程是依据相似三角形原理给出的，其形式如下 111333222333()()() ()()() ()()()()()()A S A S A S A S A S A S A S 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[,]()()()()()()()()()()()()()()T x 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 f a X X a Y Y a Z Z a X X b Y Y c Z Z l y y y f a X X a Y Y a Z Z =-+-+-=-=+-+-+--+-+-=-=+-+-+- 则1()T T X A A A L -= X 为外方位元素的近似改正数， 由于采用泰勒展开取至一次项，为减少误差，要将的出的值作为近似值进行迭代，知道小于规定的误差 三． 空间后方交会解算过程 1. 已知条件 近似垂直摄影 实验一近景单张像片空间后方交会 一、实验目的 通过单张像片空间后方交会算法编程，掌握空间后方交会的基本原理和基本步骤，为近景摄影测量解析处理方法打下理论基础。 二、实验内容 利用C++编程平台，通过给定的单片像点的像点坐标、内方位元素和控制点物方空间坐标，计算出像片的外方位元素。 三、实验步骤： 1、编程流程图： 2、编程代码： #include "iostream" #include"stdio.h" #include "stdlib.h" #include 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 单像空间后方交会和双像解析空间后方-前 方交会的算法程序实现 遥感科学与技术 摘要：如果已知每张像片的6个外方位元素，就能确定被摄物体与航摄像片的关系。因此，利用单像空间后方交会的方法，可以迅速的算出每张像片的6个外方位元素。而前方交会的计算，可以算出像片上点对应于地面点的三维坐标。基于这两点，利用计算机强大的运算能力，可以代替人脑快速的完成复杂的计算过程。 关键词：后方交会，前方交会，外方位元素，C++编程 0.引言： 单张像片空间后方交会是摄影测量基本问题之一,是由若干控制点及其相应像点坐标求解摄站参数(X S,Y S,ZS,ψ、ω、κ)。单像空间后方交会主要有三种方法:基于共线条件方程的平差解法、角锥法、基于直接线性变换的解法。而本文将介绍第一种方法，基于共线条件方程反求象片的外方位元素。 而空间前方交会先以单张像片为单位进行空间后方交会，分别求出两张像片的外方位元素，再根据待定点的一对像点坐标，用空间前方交会的方法求解待定点的地面坐标。可以说，这种求解地面点的坐标的方法是以单张像片空间后方交会为基础的，因此，单张像片空间后方交会成为解决这两个问题以及算法程序实现的关键。 1.单像空间后方交会的算法程序实现： （1）空间后方交会的基本原理：对于遥感影像,如何获取像片的外方位元素,一直是摄影测量工作者探讨的问题,其方法有:利用雷达(Radar)、全球定位系统(GPS)、惯性导航系统(I N S)以及星像摄影机来获取像片的外方位元素;也可以利用一定数量的地面控制点,根据共线方程,反求像片的外方位元素,这种方法称为单像空间后方交会(如图1所示)。 图中,地面坐标X i、Yi、Zi和对应的像点坐标x i、yi是已知的,外方位元素XS、Y S、ZS(摄站点坐标),ψ、ω、κ(像片姿态角)是待求的。 （2）空间后方交会数学模型：空间后方交会的数学模型是共线方程, 即中心投影的构像方程： 式中X、Y、Z是地面某点在地面摄影测量坐标系中的坐标,x,y是该地面点在像片上的构像点的像片坐标,对 于空间后方交会而言它们是已知的,还有主距f是已知的。而9个方向余弦a 1,a 2,a3;b1,b 2,b 3;c 1,c2,c 3是未知的，具体表达式可以取 摄影测量原理 单张影像后方交会实习 目录 一实习目的 (3) 二实习原理 (3) 1. 间接平差 (3) 2. 共线方程 (3) 3. 单向空间后方交会 (4) 三计算流程 (4) 1. 求解步骤 (4) 2.计算机框图 (4) 四程序实现 (5) 五结果分析 (6) 1.外方位元素 (6) 2.误差 (6) 3.旋转矩阵R (7) 六实习体会 (7) 1. 平台的选择 (7) 2.问题的解决 (7) 3.心得体会 (8) 七代码展示 (8) 一实习目的 为了增强同学们对后方交会公式的理解，培养同学们对迭代循环编程的熟悉感，本次摄影测量课间实习内容定为用C语言或其他程序编写单片空间后方交会程序，最终输出像点坐标、地面坐标、单位权中误差、外方位元素及其精度。 已知四对点的影像坐标和地面坐标如下。内方位元素fk=153.24mm，x0=y0=0。 本次实习，我使用了matlab2014进行后方交会程序实现。结果与参考答案一致，精度良好。 二实习原理 题干中有四个控制点在地面摄影测量坐标系中的坐标和对应的像点坐标，由此可列出8个误差方程，存在2个多余观测（n=2）。故可利用间接平差的最小二乘法则求解。 由于共线方程是非线性函数模型，为了方便计算，需要将其“线性化”。但如果仅取泰勒级数展开式的一次项，未知数的近似值改正是不精确的。因此必须采用迭代趋近法计算，直到外方位元素的改正值小于限差。 1.间接平差 间接平差为平差计算最常用的方法。在确定多个未知量的最或然值时，选择它们之间不存在任何条件关系的独立量作为未知量组成用未知量表达测量的函数关系、列出误差方程式，按最小二乘法原理求得未知量的最或然值的平差方法。 在一个间接平差问题中，当所选的独立参数X个数与必要观测值t个数相等时，可将每个观测值表达成这t个参数的函数，组成观测方程。 函数模型为：L = BX + d。 2.共线方程 共线方程是中心投影构像的数学基础，也是各种摄影测量处理方法的重要理论基础。 式中： x，y 为像点的像平面坐标； x0，y0，f 为影像的内方位元素； XS，YS，ZS 为摄站点的物方空间坐标； XA，YA，ZA 为物方点的物方空间坐标； ai，bi，ci （i = 1,2,3）为影像的3 个外方位角元素组成的9 个方向余弦。 一． 实验目的： 掌握摄影测量空间后方交会的原理，利用计算机编程语言实现空间后方交会外方位元素的解算。 二． 仪器用具及已知数据文件： 计算机windows xp 系统，编程软件（VISUAL C++6.0），地面控制点在摄影测量坐标系中的坐标及其像点坐标文件shuju.txt 。 三． 实验内容： 单张影像的空间后方交会：利用已知地面控制点数据及相应像点坐标根据共线方程反求影像的外方位元素。 数学模型：共线条件方程式： )(3)(3)(3)(1)(1)(1Zs Z c Ys Y b Xs X a Zs Z c Ys Y b Xs X a f x -+-+--+-+--= )(3)(3)(3)(2)(2)(2Zs Z c Ys Y b Xs X a Zs Z c Ys Y b Xs X a f y -+-+--+-+--= 求解过程： （1）获取已知数据。从航摄资料中查取平均航高与摄影机主距；获取控制点的地面测量坐标并转换为地面摄影测量坐标。 （2）量测控制点的像点坐标并做系统改正。 （3）确定未知数的初始值。在竖直摄影且地面控制点大致分布均匀的情况下，按如下方法确定初始值，即： n X X S ∑=0,n Y Y S ∑=0,n Z mf Z S ∑=0 φ =ω=κ=0 式中；m 为摄影比例尺分母；n 为控制点个数。 （4）用三个角元素的初始值，计算个方向余弦，组成旋转矩阵R 。 （5）逐点计算像点坐标的近似值。利用未知数的近似值和控制点的地面 坐标代入共线方程式，逐点计算像点坐标的近似值（x ）、（y ）。 （6）逐点计算误差方程式的系数和常数项，组成误差方程式。 （7）计算法方程的系数矩阵A A T 和常数项l A T ，组成法方程式。 （8）解法方程，求得外方位元素的改正数dXs ,S dY ，s dZ ，d φ,d ω,d κ。 （9）用前次迭代取得的近似值，加本次迭代的改正数，计算外方位元素 的新值。 作业报告 空间后方交会 专业：测绘工程 班级：2008级(1)班姓名：陈闻亚 指导教师：陈强 2010 年 4 月16 日 1 作业任务------------------------------------------------------------------------------------ 3 2 作业思想 --------------------------------------------------------------------------------------- 3 3 作业条件及数据 -------------------------------------------------------------------- 3 4 作业过程--------------------------------------------------------------------------- 3 5 源程序----------------------------------------------------------------------------- 4 6 计算结果--------------------------------------------------------------------------- 17 7心得体会与建议----------------------------------------------------------------------------- 17 1 作业任务 计算近似垂直摄影情况下后方交会解。即利用摄影测量空间后方交会的方法，获取相片的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。 单像空间后方交会 姓名： 学号： 时间： 目录 一、作业任务 ............................................................................................................... - 3 - 二、计算原理 ............................................................................................................... - 3 - 三、算法流程 ............................................................................................................... - 7 - 四、源程序 ................................................................................................................... - 8 - 五、计算结果 ............................................................................................................... - 8 - 六、结果分析 ............................................................................................................... - 8 - 七、心得与体会 ........................................................................................................... - 8 - 八、附页 ....................................................................................................................... - 8 - 1.c++程序 ........................................................................................................... - 8 - 2.C++程序截图.................................................................................................. - 15 - 3.matlb程序..................................................................................................... - 16 - 空间后方交会程序 ————————————————————————————————作者：————————————————————————————————日期： ? 一． 实验目的: 掌握摄影测量空间后方交会的原理,利用计算机编程语言实现空间 后方交会外方位元素的解算。 二． 仪器用具及已知数据文件： 计算机wind ｏws xp 系统,编程软件(VI ＳＵA Ｌ C ＋+6．0)，地面控 制点在摄影测量坐标系中的坐标及其像点坐标文件shu ｊu.txt 。 三． 实验内容: 单张影像的空间后方交会:利用已知地面控制点数据及相应像点坐标根据 共线方程反求影像的外方位元素。 数学模型：共线条件方程式: ) (3)(3)(3) (1)(1)(1Zs Z c Ys Y b Xs X a Zs Z c Ys Y b Xs X a f x -+-+--+-+--= ) (3)(3)(3)(2)(2)(2Zs Z c Ys Y b Xs X a Zs Z c Ys Y b Xs X a f y -+-+--+-+--= 求解过程: （1)获取已知数据。从航摄资料中查取平均航高与摄影机主距;获取 控制点的地面测量坐标并转换为地面摄影测量坐标。 （2)量测控制点的像点坐标并做系统改正。 (３)确定未知数的初始值。在竖直摄影且地面控制点大致分布均匀 的情况下,按如下方法确定初始值,即: n X X S ∑=0,n Y Y S ∑=0,n Z mf Z S ∑=0 φ =ω=κ=0 式中;ｍ为摄影比例尺分母;ｎ为控制点个数。 （４)用三个角元素的初始值,计算个方向余弦，组成旋转矩阵R 。 (5)逐点计算像点坐标的近似值。利用未知数的近似值和控制点的地面坐标代入共 线方程式，逐点计算像点坐标的近似值（x ）、（y ）。 （6)逐点计算误差方程式的系数和常数项,组成误差方程式。 (7）计算法方程的系数矩阵A A T 和常数项l A T ,组成法方程式。 （8）解法方程，求得外方位元素的改正数dXs ,S dY ，s dZ ,d φ,ｄω，d κ。 （９）用前次迭代取得的近似值，加本次迭代的改正数，计算外方位元素的新值。 空间后方交会的直接解 空间后方交会，即由物方已知若干个控制点以及相应的像点坐标，解求摄站的坐标与影像的方位，这是一个摄影测量的基本问题。通常采用最小二乘解算，由于原始的观测值方程是非线性的，因此，一般空间后方交会必须已知方位元素的初值，且解算过程是个迭代解算过程。但是，在实时摄影测量的某些情况下，影像相对于物方坐标系的方位是任意的，且没有任何初值可供参考。这时常规的空间后方交会最小二乘算法就无法处理，而必须建立新的空间后方交会的直接解法。 直接解法的基本思想是将它分成两步：先求出三个已知点i P 到摄站S 的距离i S ；然后求出摄站S 的坐标和影像方位。 物方一已知点()i i i i ,Z ,Y X P 在影像上的成像()i i i ,y x p ，根据影像已知的内方位元素()0 ,y f,x 可求得从摄站()S S S S ,Z ,Y X 到已知点i P 的观测方向i ,βαi 。 () ??? ????-+-= -=2 020 tan tan x x f y y βf x x αi i i i i (1) 距离方程组可以写成如下形式： ?? ??? =+++=+++=+++020202312 1133123232 3322322122 2211221b x x x a x b x x x a x b x x x a x (2) 其中()j ;i ,,i,j S ,b a ij ij ij ij ≠===321cos ?。因此，解算摄站S 到三个 控制点的距离问题，被归结为解算一个三元二次联立方程组的问题。这个方程组的解算方法选用迭代法。 迭代计算公式可写成： 主方法： private void Cal_Click(object sender, EventArgs e) { string[] lines = RichText.Text.Split('\n'); long m = lines.Length; m = m - 1;//真实数据行数 double[] Coor_x = new double[m];//已知点x坐标 double[] Coor_y = new double[m];//已知点x坐标 double[] Coor_X = new double[m];//已知点X坐标 double[] Coor_Y = new double[m];//已知点Y坐标 double[] Coor_Z = new double[m];//已知点Z坐标 ///赋值 for (int i = 0; i < m; i++) { string[] FJstring = Regex.Split(lines[i+1], ","); Coor_x[i] = 0.001*(Convert.ToDouble(FJstring[0])); Coor_y[i] = 0.001 *( Convert.ToDouble(FJstring[1])); Coor_X[i] = Convert.ToDouble(FJstring[2]); Coor_Y[i] = Convert.ToDouble(FJstring[3]); Coor_Z[i] = Convert.ToDouble(FJstring[4]); } if (textBox_m.Text == "") { MessageBox.Show("请输入参数！"); } if (textBox_m.Text != "") { double M = double.Parse(textBox_m.Text);//比例尺 double f = 0.001 * (double.Parse(textBox_f.Text));//焦距 double x0 = 0.001 * double.Parse(textBox_x0.Text);//内方位元素x0 double y0 = 0.001 * double.Parse(textBox_y0.Text);//内方位元素y0 double X0 = 0, Y0 = 0, Z0 = 0;//外方位坐标元素初始值 double min = (double.Parse(textBox_k.Text));//焦距 double angle1 = 0, angle2 = 0, angle3 = 0;//外方位角元素初始值 for (int i = 0; i < m; i++) {//累加 X0 = Coor_X[i] + X0; Y0 = Coor_Y[i] + Y0; 单像空间后方交会程序报告 指导老师：刘老师 班级：测绘101 姓名：尚锋 学号： 19号 1、应用程序的主入口部分的代码： using System; using System.Collections.Generic; using System.Linq; using System.Windows.Forms; namespace单像空间后方交会 { static class Program { /// public double Z; } private double f; //主距 private double u; //u为外方位元素，下面5个相同 private double w; private double k; private double Xs; private double Ys; private double Zs; private image_point[] p = new image_point[4]; //四个控制点 private double[] R = new double[9]; //旋转矩阵 private double[] a = new double[8]; //像点坐标近似值 private double[,] A = new double[8, 6]; //误差方程式系数 private double[] L = new double[8]; //误差方程式常数项 private int count = 0; //统计代次数 public Tongyong(double g, double[] q) //构造函数，初始化各变量,单位m { f = g; for (int i = 0; i < 4; i++) { int j = i * 5; p[i].x = q[j]; p[i].y = q[j + 1]; p[i].X = q[j + 2]; p[i].Y = q[j + 3]; p[i].Z = q[j + 4]; } double ave = 0, sum = 0; //求比例尺分母 for (int i = 0; i < 3; i++) { for (int j = i + 1; j < 4; j++) { sum += Math.Sqrt(Math.Pow(p[i].X - p[j].X, 2) + Math.Pow(p[i].Y - p[j].Y, 2)) / Math.Sqrt(Math.Pow(p[i].x - p[j].x, 2) + Math.Pow(p[i].y - p[j].y, 2)); } } ave = sum / 6; u = 0; //给定外方位元素的初始值,角度均设置为0 w = 0; k = 0; Xs = (p[0].X + p[1].X + p[2].X + p[3].X) / 4; //Xs为四个控制点X的平均值，Ys类似 输入文件形式如下： C++源程序如下： #include 任务已知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.3 2 2195.17 2 -53.40 82.21 37631.08 31324.51 728.69 3 -14.78 -76.63 39100.97 24934.98 2386.50 4 10.46 64.43 40426.54 30319.81 757.31 求近似垂直摄影情况下后方交会解 设计任务 1、确定未知数的初始值： Φ0 =ω0 =К0 = 0 , 内方位元素，，f=153.24mm。 ； = 38437m ； = 27963.16m 2、计算旋转矩阵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阵为 #include // for(i=k-1;i>=0;i--) // { // if(q[i][h]==0) // continue; // p=q[k][h]/q[i][h]; // // p=q[i][h]/q[k][h]; // for(j=11;j>0;j--) // { // q[i][j]*=p; // q[i][j]-=q[k][j]; // } // } // for(i=0;i<6;i++)//将对角线上数据化为1 // { // p=1.0/q[i][i]; // for(j=0;j<12;j++) // q[i][j]*=p; // } // for(i=0;i<6;i++) //提取逆矩阵 // for(j=0;j 《摄影测量基础》答卷（A） 一、填空题（20分，每空1分） 1、摄影测量中常用的坐标系有像平面直角坐标系、像空间直角坐标系、像空间辅助坐标系、地面摄影测量坐标系、地面测量坐标系。 2、解求单张像片的外方位元素最少需要 3 个平高地面控制点。 3、GPS辅助空中三角测量的作用是大量减少甚至完全免除地面控制点，缩短成图周期，提高生产效率，降低生产成本。 4、两个空间直角坐标系间的坐标变换最少需要 2 个平高和 1 个高程地面控制点。 5、摄影测量加密按平差范围可分为单模型、单航带和区域网三种方法。 6、摄影测量的发展经历了模拟摄影测量、解析摄影测量和数字摄影测量三个 阶段。 7、恢复立体像对左右像片的相互位置关系依据的是共面条件方程。 二、名词解释（20分，每个4分） 1、内部可靠性：一定假设下，平差系统所能发现的模型误差的最小值。 2、绝对定向元素：确定模型在地面空间坐标系中的绝对位置和姿态的参数。 3、像主点：相机主光轴与像平面的交点。 4、带状法方程系数矩阵的带宽：带状法方程系数矩阵的主对角线元素沿某行（列）到最远非零 元素间所包含未知数的个数。 5、自检校光束法区域网平差：选用若干附加参数构成系统误差模型，在光束法区域网平差中同时解求这些附加参数，从而在平差过程中自行检定和消除系统误差影响的区 域网平差。 《摄影测量基础》答卷（A） 一、填空题（20分，每空1分） 1、表示航摄像片的外方位角元素可以采用以Y轴为主轴的κ ω ?? ? 、以X轴为主轴的''' κ?ω ??和以Z轴为主轴的kaA ?? 三种转角系统。 2、航摄像片是所覆盖地物的中心投影。 3、摄影测量加密按数学模型可分为航带法、独立模型法和光束法三种方法。 4、从航摄像片上量测的像点坐标可能带有摄影材料变形、摄影机物镜畸变、大气折光误差和地球曲率误差四种系统误差。 5、要将地物点在摄影测量坐标系中的模型坐标转换到地面摄影测量坐标系，最少需要2 个平高和 1 个高程地面控制点。 6、带状法方程系数矩阵的带宽是指法方程系数矩阵中主对角线元素起沿某一行到最远处的非零元素间所包含的未知数个数。 7、人眼观察两幅影像能产生立体视觉的基本条件是在不同摄站获取的具有一定重叠的两幅影像、观察时每只眼睛只能看一张像片、两幅影像的摄影比例尺尽量一致和两幅影像上相同地物的连线与眼基线尽量平行。 二、名词解释（20分，每个4分） 1、内部可靠性：一定的假设条件下，平差系统所能发现的模型误差的最小值。 2、解析相对定向：根据摄影时同名光线位于一个核面的条件，利用共面条件方程解算立体像对中两张像片的相互关系参数，使同名光线对对相交。 3、GPS辅助空中三角测量：利用载波相位差分GPS动态定位技术获取摄影时刻摄影中心的三维坐标，将其作为带权观测值引入摄影测量区域网平差中，整体确定物方点坐标和像片外方位元素并对其质量进行评定的理论和方法。 4、主合点：地面上一组平行于摄影方向线的光束在像片上的构像。 5、单片空间后方交会：在单张像片上，利用一定数量的地面控制点及其对应的像点坐标，根据共线条件方程求解像片的6个外方位元素。 %读取包含点号的像点坐标、控制点坐标分别到矩阵IP、CP% IP=load('f:'); CP=load('f:'); %删除像点坐标、控制点坐标矩阵中的点号列% IP(:,1)=[]; CP(:,1)=[]; %像片大小% H=4008; W=5344; IP1=IP(:,1)-W/2; ￥ IP2=H/2-IP(:,2); IP=[IP1,IP2]; %焦距% fx=; fy=; %内方位元素% x0=; y0=; %畸变参数% k1=; { k2=; p1=; p2=; %外方位元素初值% Xs=3000; Ys=-100; Zs=100; Phi=0; Omega=0; Kappa=0; ! m=size(IP,1); AIP=zeros(m,2); L=zeros(2*m,1); A=zeros(2*m,6); X=ones(6,1); while X(4,1)>=||X(5,1)>=||X(6,1)>= %旋转矩阵的系数% a1=cos(Phi)*cos(Kappa)-sin(Phi)*sin(Omega)*sin(Kappa); a2=-cos(Phi)*sin(Kappa)-sin(Phi)*sin(Omega)*cos(Kappa); a3=-sin(Phi)*cos(Omega); 》 b1=cos(Omega)*sin(Kappa); b2=cos(Omega)*cos(Kappa); b3=-sin(Omega); c1=sin(Phi)*cos(Kappa)+cos(Phi)*sin(Omega)*sin(Kappa); c2=-sin(Phi)*sin(Kappa)+cos(Phi)*sin(Omega)*cos(Kappa); c3=cos(Phi)*cos(Omega); for i=1:m %求像点坐标常数项% r=sqrt((IP(i,1)-x0)^2+(IP(i,2)-y0)^2); Dx=(IP(i,1)-x0)*(k1*(r^2)+k2*(r^4))+p1*(r^2+2*((IP(i,1)-x0)^2))+2*p2*(IP(i,1)-x0)*(IP(i,2)-y0); 、 Dy=(IP(i,2)-y0)*(k1*(r^2)+k2*(r^4))+p2*(r^2+2*((IP(i,2)-y0)^2))+2*p1*(IP(i,1)-x0)*(IP(i,2)-y0); X_=a1*(CP(i,1)-Xs)+b1*(CP(i,2)-Ys)+c1*(CP(i,3)-Zs); Y_=a2*(CP(i,1)-Xs)+b2*(CP(i,2)-Ys)+c2*(CP(i,3)-Zs); Z_=a3*(CP(i,1)-Xs)+b3*(CP(i,2)-Ys)+c3*(CP(i,3)-Zs); %求像点坐标的近似值% AIP(i,1)=-fx*X_/Z_+Dx+x0; AIP(i,2)=-fy*Y_/Z_+Dy+y0; %求误差方程式系数% a11=(a1*fx+a3*IP(i,1))/Z_; a12=(b1*fx+b3*IP(i,1))/Z_; 【 a13=(c1*fx+c3*IP(i,1))/Z_; a14=sin(Omega)*IP(i,2)-(IP(i,1)*(IP(i,1)*cos(Kappa)-IP(i,2)*sin(Kappa))/fx+fx*cos(Kappa))*cos(Om ega); a15=-fx*sin(Kappa)-IP(i,1)*(IP(i,1)*sin(Kappa)+IP(i,2)*cos(Kappa))/fx; a16=IP(i,2); a21=(a2*fy+a3*IP(i,2))/Z_; a22=(b2*fy+b3*IP(i,2))/Z_; a23=(c2*fy+c3*IP(i,2))/Z_; a24=-sin(Omega)*IP(i,1)-(IP(i,2)*(IP(i,1)*cos(Kappa)-IP(i,2)*sin(Kappa))/fy-fy*sin(Kappa))*cos(Om ega); a25=-fy*cos(Kappa)-IP(i,2)*(IP(i,1)*sin(Kappa)+IP(i,2)*cos(Kappa))/fy; a26=-IP(i,1); 。 %求常数项% L(2*i-1,1)=IP(i,1)-AIP(i,1); L(2*i,1)=IP(i,2)-AIP(i,2); %求误差方程式系数阵% A(2*i-1,1)=a11;空间后方交会编程实习报告
空间后方交会的解算
近景单张像片空间后方交会
单像空间后方交会和双像解析空间后方-前方交会的算法程序实现
摄影测量学后方交会matlab实习报告
空间后方交会程序
作业4--空间后方交会
单像空间后方交会实验报告(c++版)
空间后方交会程序
空间后方交会的直接解
单片空间后方交会C#源代码
单像空间后方交会程序报告
摄影测量单片空间后方交会程序《C++》
空间后方交会报告
C语言空间后方交会源代码
摄影测量
单向空间后方交会matlab编程