后方交会MATLAB程序实习报告 本文关键词：交会，实习报告，程序，MATLAB
后方交会MATLAB程序实习报告 本文简介：《摄影测量学》单像空间后方交会实习报告班级：XXXX00000姓名：XXX000学号：XXXXXXXXXXXXX指导教师：XXX1、实习目的1、掌握空间后方交会的定义和实现算法；2、了解摄影测量平差的基本过程；3、熟练MATLAB等程序编写。2、实习原理利用至少三个已知地面控制点的坐标、、，与其影像
后方交会MATLAB程序实习报告 本文内容：
《摄
影
测
量
学》
单像空间后方交会
实习报告
班
级：
XXXX00000
姓
名：
X
X
X000
学
号：XXXXXXXXXXXXX
指导教师：
X
X
X
1、
实习目的
1、掌握空间后方交会的定义和实现算法；
2、了解摄影测量平差的基本过程；
3、熟练MATLAB等程序编写。
2、
实习原理
利用至少三个已知地面控制点的坐标、、，与其影像上对应的三个像点的影像坐标、、，根据共线方程，反求该像片的外方位元素。
共线条件方程式：
将共线式线性化并取一次小值项得：
3、
解算过程
①获取已知数据。包括影像比例尺1/m，平均摄影距离（航空摄影的航高）H，内方位元素x0、y0、f，控制点的空间坐标X、Y、Z。
②量测控制点的像点坐标并进行必要的影像坐标系统误差改正，得到像点坐标。
③确定未知数的初始值。单像空间后方交会必须给出待定参数的初始值，在竖直航空摄影且地面控制点大体对称分布的情况下，Xs0和Ys0为均值，Zs0为航高，φ、ω、κ的初值都设为0。或者κ的初值可在航迹图上找出或根据控制点坐标通过坐标正反变换求出。
④计算旋转矩阵R。利用角元素近似值计算方向余弦值，组成R阵。
⑤逐点计算像点坐标的近似值。利用未知数的近似值按共线条件式计算控制点像点坐标的近似值（x），(y)。
⑥逐点计算误差方程式的系数和常数项，组成误差方程式。
⑦计算法方程的系数矩阵ATA与常数项ATL，组成法方程式。
⑧解求外方位元素。根据法方程，解求外方位元素改正数，并与相应的近似值求和，得到外方位元素新的近似值。
⑨检查计算是否收敛。将所求得的外方位元素的改正数与规定的限差比较，通常对φ，ω，κ的改正数△φ，△ω，△κ给予限差，通常为0.000001弧度，当3个改正数均小于0.000001弧度时，迭代结束。否则用新的近似值重复④～⑧步骤的计算，直到满足要求为止。
四、程序框图
五、主要代码
[fn,pn,fi]=uigetfile(.txt,选择文件
);
I=load([pn
fn]);
x=I(:,1)
;y=I(:,2)
;X=I(:,3)
;Y=I(:,4)
;Z=I(:,5)
;
s=0;S=0;
for
i=1:3;
j=i+1;
sij=sqrt((x(i)-x(j))^2+(y(i)-y(j))^2);
Sij=sqrt((X(i)-X(j))^2+(Y(i)-Y(j))^2);
s=s+sij;
S=S+Sij;
end
m=S*1000/s;
f=153.24;
Xs0=(X(1)+X(2)+X(3)+X(4))/4;
Ys0=(Y(1)+Y(2)+Y(3)+Y(4))/4;
Zs0=m*f;
fai0=0;omig0=0;ka0=0;
for
v=1:+inf;
a1=cos(fai0)*cos(ka0)-sin(fai0)*sin(omig0)*sin(ka0);
a2=-cos(fai0)*sin(ka0)-sin(fai0)*sin(omig0)*cos(ka0);
a3=-sin(fai0)*cos(omig0);
b1=cos(omig0)*sin(ka0);
b2=cos(omig0)*cos(ka0);
b3=-sin(omig0);
c1=sin(fai0)*cos(ka0)+cos(fai0)*sin(omig0)*sin(ka0);
c2=-sin(fai0)*sin(ka0)+cos(fai0)*sin(omig0)*cos(ka0);
c3=cos(fai0)*cos(omig0);
R=[a1
a2
a3;b1
b2
b3;c1
c2
c3];
l=[];A=[];
for
h=1:4;
O=a1*(X(h)-Xs0)+b1*(Y(h)-Ys0)+c1*(Z(h)-Zs0);
P=a2*(X(h)-Xs0)+b2*(Y(h)-Ys0)+c2*(Z(h)-Zs0);
Q=a3*(X(h)-Xs0)+b3*(Y(h)-Ys0)+c3*(Z(h)-Zs0);
x1=-f*O/Q;
y1=-f*P/Q;
a11=(a1*f+a3*x(h))/Q;
a12=(b1*f+b3*x(h))/Q;
a13=(c1*f+c3*x(h))/Q;
a14=y(h)*sin(omig0)-(x(h)/f*(x(h)*cos(ka0)-y(h)*sin(ka0))+f*cos(ka0))*cos(omig0);
a15=-f*sin(ka0)-x(h)/f*(x(h)*sin(ka0)+y(h)*cos(ka0));
a16=y(h);
a21=(a2*f+a3*y(h))/Q;
a22=(b2*f+b3*y(h))/Q;
a23=(c2*f+c3*y(h))/Q;
a24=-x(h)*sin(omig0)-(y(h)/f*(x(h)*cos(ka0)-y(h)*sin(ka0))-f*sin(ka0))*cos(omig0);
a25=-f*cos(ka0)-y(h)/f*(x(h)*sin(ka0)+y(h)*cos(ka0));
a26=-x(h);
lx=x(h)-x1;
ly=y(h)-y1;
lh=[lx
ly]
;
Ah=[a11,a12,a13,a14,a15,a16;a21,a22,a23,a24,a25,a26];
A=[A;Ah];
l=[l;lh];
end
XX=inv(AA)*Al;
Xs0=Xs0+XX(1);
Ys0=Ys0+XX(2);
Zs0=Zs0+XX(3);
fai0=fai0+XX(4);
omig0=omig0+XX(5);
ka0=ka0+XX(6);
R=[a1,a2,a3;b1,b2,b3;c1,c2,c3]
Xs=Xs0
Ys=Ys0
Zs=Zs0
fai=fai0
omig=omig0
ka=ka0
If
abs(XX(4))<0.0000291b1,b2,b3;c1,c2,c3];
V=A*XX-l;
Qii=inv(AA);
m0=sqrt(VV/2);
mi=m0*sqrt(Qii);
m=diag(mi);
m=[m(1),m(2),m(3),m(4)/pi*180*3600,m(5)/pi*180*3600,m(6)/pi*180*3600];
mXs=m(1);mYs=m(2);mZs=m(3);mfai=m(4);momig=m(5);mka=m(6);
fp=fopen(
后方交会计算结果.txt,wt
);
fprintf(fp,迭代次数:%d/n,v);
fprintf(fp,/n旋转矩阵R：/n
);
[m,n]=size(R);
for
i=1:1:m
for
j=1:1:n
if
j==n
fprintf(fp,%g/n,R(i,j));
else
fprintf(fp,%g/t,R(i,j));
end
end
end
fprintf(fp,/n外方位元素解：/n
);
fprintf(fp,Xs=%g/t,Xs);fprintf(fp,Ys=%g/t,Ys);fprintf(fp,Zs=%g/n,Zs);
fprintf(fp,φ=%g/t,fai);fprintf(fp,ω=%g/t,omig);fprintf(fp,κ=%g/n,ka);
fprintf(fp,/n单位权中误差的值：mo=%g/n,m0);
fprintf(fp,/n外方位元素中误差为：/n
);
fprintf(fp,mXs=%g米/t,mXs);fprintf(fp,mYs=%g米/t,mYs);fprintf(fp,mZs=%g米/n,mZs);
fprintf(fp,mφ=%g秒/t,mfai);fprintf(fp,mω=%g秒/t,momig);fprintf(fp,mκ=%g秒/n,mka);
6、
实习结果
1、起算数据
2、运行结果
3、几次迭代的数据