摄影测量作业3-空间后方交会计算
空间后方交会与前方交会求解地面点坐标的计算方法

双像解析计算的空间后交-前交方法当我们通过航空摄影,获得地面的一个立体像对时,采用双像解析计算的空间后交-前交方法计算地面点的空间点位。
这种方法首先由单片后方交会求出左、右像片的外方位元素,再用空间前方交会公式求出待定点坐标,其具体的作业步骤如下:● 像片野外控制测量一个立体像对采用空间后方交会-前方交会法计算点的地面坐标时,像对内必须具有一定数量的地面控制点坐标。
一般情况下,在一个像对的重叠范围四个角上,找出四个明显地物点,在野外判识出地面的实际位置,并准确地在像片上刺出各点的位置,要求在像片的背面绘出各点与周围地物关系的点位略图,加注记说明。
然后用普通测量计算方法,求出四个控制点的地面坐标X,Y,Z 。
● 用立体坐标量测仪测像点的坐标像片在仪器上归心定向后,测出四个控制点的像片坐标()11y ,x 与()22y ,x ,然后测出所需要解求的地面点坐标()11y ,x 和()22y ,x 。
● 空间后方交会法计算像片外方位元素利用控制点分别计算每个像片的六个外方位元素,包括:S1X ,S1Y ,S1Z ,1ϕ,1ω,1κ和S2X ,S2Y ,S2Z ,2ϕ,2ω,2κ。
● 空间前方交会计算所求点的地面坐标1. 用各自像片的角元素,计算出左、右像片的旋转矩阵1R 与2R 。
2. 根据左、右像片的外方位线元素计算摄影基线分量X B ,Y B ,Z B :⎪⎭⎪⎬⎫-=-=-=S1S2Z S1S2Y S1S2X Z Z B Y Y B X X B3. 逐点计算像点的像空间辅助坐标:⎥⎥⎥⎦⎤⎢⎢⎢⎣⎡-=⎥⎥⎥⎦⎤⎢⎢⎢⎣⎡f y x R Z Y X 111111,⎥⎥⎥⎦⎤⎢⎢⎢⎣⎡-=⎥⎥⎥⎦⎤⎢⎢⎢⎣⎡f y x R Z Y X 1111114. 计算点投影系数:⎪⎪⎭⎪⎪⎬⎫--=--=12211121221221Z X Z X X B Z B N Z X Z X X B Z B N Z X Z X5. 计算所求点的地面摄影测量坐标:⎥⎥⎥⎦⎤⎢⎢⎢⎣⎡+⎥⎥⎥⎦⎤⎢⎢⎢⎣⎡+⎥⎥⎥⎦⎤⎢⎢⎢⎣⎡=⎥⎥⎥⎦⎤⎢⎢⎢⎣⎡+⎥⎥⎥⎦⎤⎢⎢⎢⎣⎡=⎥⎥⎥⎦⎤⎢⎢⎢⎣⎡+⎥⎥⎥⎦⎤⎢⎢⎢⎣⎡=⎥⎥⎥⎦⎤⎢⎢⎢⎣⎡222222111222222222111111111Z N Y N X N B B B Z Y X Z N Y N X N Z Y X Z N Y N X N Z Y X Z Y X Z Y X S S S S S S S S S A A A 6. 重复3-5步骤完成所有点地面坐标的计算。
空间后方交会的解算

空间后方交会的解算一. 空间后方交会的目的摄影测量主要利用摄影的方法获取地面的信息,主要是是点位信息,属性信息,因此要对此进行空间定位和建模,并首先确定模型的参数,这就是空间后方交会的目的,用以求出模型外方位元素。
二. 空间后方交会的原理空间后方交会的原理是共线方程。
共线方程是依据相似三角形原理给出的,其形式如下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个不在一条直线上的已知地面点坐标就可以求出像片的外方位元素。
摄影测量-空间前交、后交【精选文档】

空间后交—前交程序设计(实验报告)姓名:班级:学号:时间:空间后交-前交程序设计一、实验目的用 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;%%%%%权的矢量化,这是等精度时的,如果非,将函数参数改为PP=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。
摄影测量学空间后方交会实验报告

摄影测量学实验报告实验一、单像空间后方交会学院:建测学院班级:测绘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.实验结果四.实验总结此次实验让我深入了解单像空间后方交会的计算过程,加强了对空间后方交会基本公式和误差方程式,法线方程式的记忆。
(完整版)摄影测量作业部分答案

1.什么是摄影测量学,摄影测量发展的三个阶段摄影测量学是对研究的物体进行摄影,测量和解译所获得的影像,获取被摄物体的几何信息和物理信息的一门科学。
三个阶段:模拟、解析、数字摄影测量。
2.根据成图的需要,规定了摄影比例尺后,如何选择航空摄影机与航摄高度?答:采用特宽角航摄机,航高值就小,采用常角或窄角航摄机,航高值就大。
(航高的大小将决定飞机实际升限和最小安全高度的限制,另外,测图的高程测定精度与航高有关(高程中误差与航高成正比))大比例尺单像测图,应选用常角或窄角航摄机,小比例尺立体测图应选用特宽角航摄机。
3.什么是航摄像片的内外方位元素,各有何作用?答:内方位元素包括三个参数,即摄影中心S到像片的垂距(主距)f及像主点o在像框标坐标系中的坐标00,y x,用其来恢复摄影光束。
确定摄影光束在摄影瞬间的空间位置和姿态的参数,称为外方位元素,一张的外方位元素包括六个参数,其中有三个是直线元素,用于描述摄影中心的空问坐标值;另外三个是角元素,用于表达像片面的空间姿态。
4.摄影测量中常用的坐标系有哪些,各有何作用?答:摄影测量中常用的坐标系有两大类。
一类是用于描述像点的位置,称为像方空间坐标系;另——类是用于描述地面点的位置.称为物方空间坐标系。
(1)像方空间坐标系①像平面坐标系像平面坐标系用以表示像点在像平面上的位置,通常采用右手坐标系,xy轴的选择按需要而定.在解析和数字摄影测量中,常根据框标来确定像平面坐标系,称为像框标坐标系。
②像空间坐标系为了便于进行空间坐标的变换,需要建立起描述像点在像空间位置的坐标系,即像空间坐标系。
以摄影中心S为坐标原点,xy,轴与像平面坐标系的xy轴平行,z轴与主光轴重合,形成像空间右手直角坐标系s-xyz③像空间辅助坐标系像点的像空间坐标可直接以像平面坐标求得,但这种坐标的待点是每张像片的像空间坐标系不统一,这给计算带来困难。
为此,需要建立一种相对统一的坐标系.称为像空间辅助坐标系,用s-XYZ表示。
摄影测量学空间后方交会实验报告13页

摄影测量学空间后方交会实验报告13页报告摘要:本实验以三张已知高度物体的相片为样本,运用摄影测量学的空间后方交会算法,通过MATLAB编程实现对物体三维坐标的计算,并加以检验与探讨,得出结论相对合理。
同时,本实验也通过对MATLAB编程的应用,掌握了空间后方交会算法的理论及实践方法。
1. 实验目的(1)学习和应用摄影测量学中的空间后方交会算法,掌握其计算方法和程序实现过程。
(2)了解数字像相机的相关特性,并掌握其使用方法。
(3)通过对样本数据的处理,熟悉和掌握MATLAB编程的应用技巧。
2. 实验器材数字相机一部,尺子一把,样本图像三张,MATLAB软件。
3. 实验原理在精密测量领域中,采用摄影测量学的空间后方交会算法,可实现三维坐标的测量和重建。
这种方法是将已知物体的照片,通过对像点的提取及校准,得出像点坐标系下的物体的三维坐标系下的坐标。
这种算法的基本思路是:利用像点坐标系下的物体三维坐标系下的坐标关系,构建一个误差最小的方程组,通过矩阵的求解,得到物体三维坐标系下的坐标。
数字相机是一种基于CCD或CMOS成像器材料的成像设备,根据数字信号的处理能力,合成电子图像。
数字相机的性能主要包括分辨率、感光度、曝光控制、焦距、光圈等参数。
使用数字相机拍摄时,应根据拍摄对象的光线条件、距离、尺寸、景深等因素,进行调节。
4. 实验过程(1)利用数字相机拍摄三张已知高度物体的照片,并在样本上面贴标记,用尺子测算实际高度。
(2)利用图像处理软件MATLAB,对照片进行像点识别和校准,得到像点坐标系下的坐标。
(3)根据相片中已知物体的测高值及像点坐标系下的坐标值,通过MATLAB编写空间后方交会的程序算法,得出物体的三维坐标系下的坐标。
(4)对得出的坐标值进行检验及探讨,分析误差来源及部分工具库的使用方法。
5. 实验结果与分析本实验所得出的三维坐标值,原本应是在一个确定点之间展开的点集。
知道参数计算不全或精度不够是常有的事情(尤其在没有精密测量器具的条件下),这种情况我们应该考虑从给出的角度和图像来看和计算得到位置。
《摄影测量学》第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)。
方法确定初始值,即:
q=w=k=0
式中:m 为摄影比例尺分母;n 为控制点个数。 (4)用三个角元素的初始值,计算个方向余弦,组成旋转矩阵 R。 (5)逐点计算像点坐标的近似值。利用未知数的近似值和控制点的地面 坐标代入共线方程式,逐点计算像点坐标的近似值(x)、(y)。 (6)逐点计算误差方程式的系数和常数项,组成误差方程式。 (7)计算法方程的系数矩阵 A 和常数项 L,组成法方程式。 (8)解法方程,求得外方位元素的改正数 dXs,dYs,dZs,dq,dw,dk。 (9)用前次迭代取得的近似值,加本次迭代的改正数,计算外方位元素 的新值。 (10)将求得的外方位元素改正数与规定的限差比较,若小于限差则迭代结束。否则用 新的近似值重复(4)~(9),直到满足要求为止。
4
10.46 64.43 40426.54 30319.81 757.31
以单像空间后方交会方法,求解该像片的外方位元素。
作业要求: 1、用一门计算机语言(如 C,C++,C#,VB)编写单像空间后方交会程序,各角
元素迭代计算至其改正值小于 6 秒。 2、提交正式的课程作业报告。 3、作业报告包括:封面、目录、正文等,其中正文部分包括:作业任务、计算
2
一.实验目的
掌握摄影测量空间后方交会的原理,利用计算机编程语言实现空间后方交 会外方位元素的解算。
二.实验数据及实验设备
计算机 windows10 系统,编程软件 VS2015,地面控制点在摄影测量坐标系 中的坐标及其像点坐标文件 data.txt。
已知条件
摄影机主距 f=153.24mm,x0=0.01mm,y0=-0.02mm, 像片比例尺为 1:40000,像片上四个像点的框标坐标及其对应的地面坐标如下表。
CString strPathName = dlgOpenFile.GetPathName();//获取选择的文件的完整路径
CString strFileTitle = dlgOpenFile.GetFileTitle();//获取文件名 CString strExtName = dlgOpenFile.GetFileExt();//获取文件扩展名 CStdioFile sf; //创建文件对象
点号
像点坐标
地面坐标
x(mm) y(mm) X(m)
Y(m)
Z(m)
1
-86.15 -68.99 36589.41 25273.32 2195.17
2
-53.40 82.21 37631.08 31324.51 728.69
3
-14.78 -76.63 39100.97 24934.98 2386.50
double _Z = a3*(XA - XS) + b3*(YA - YS) + c3*(ZA - ZS);
A(2*i, 0) = (, 1) = (b1*f + b3*x) / _Z; A(2*i, 2) = (c1*f + c3*x) / _Z;
5
CFileDialog dlgOpenFile(TRUE, _T("txt"), NULL, OFN_FILEMUSTEXIST, _T("(文本文件)|*.txt|(所有文件)|*.*)||"));
if (dlgOpenFile.DoModal() == IDCANCEL) return;//如果选择取消按钮,则退出
for (int j = 0; j < 5;j++) listXYZ(i,j)=_wtof(m_list1.GetItemText(i, j + 1)); xyXYZ = listXYZ; ButtonClickeOk = 0;//表示导入数据按钮已被点击并且已导入数据 }
CMatrix CKongJianHouFangJiaoHuiDlg::GetX(CMatrix A, CMatrix L)//计算外方位元素改正数矩 阵 {
double w = XX(0, 4); double k = XX(0, 5); double q = XX(0, 3);; double a1 = cos(q)*cos(k) - sin(q)*sin(w)*sin(k); double a2 = -cos(q)*sin(k) - sin(q)*sin(w)*cos(k); double a3 = -sin(q)*cos(w); double b1 = cos(w)*sin(k); double b2 = cos(w)*cos(k); double b3 = -sin(w); double c1 = sin(q)*cos(k) + cos(q)*sin(w)*sin(k);
m_list1.SetItemText(index1, i + 1, strSplit[i]);
} }
} sf.Close();
6
int Line = m_list1.GetItemCount(); CMatrix listXYZ(Line, 5); for (int i = 0; i < Line; i++)
CMatrix X,_A,_AA,N_AA; _A = ~A;//A 的转置 _AA = _A*A; N_AA = _AA.Inv();//_AA 的逆矩阵 X = N_AA*_A*L; return X; }
CMatrix CKongJianHouFangJiaoHuiDlg::GetA(CMatrix xyXYZ, double f, CMatrix XX)//计算系数矩 阵A {
//以读的形式打开文件,如果打开失败 if (!sf.Open(strPathName, CFile::modeRead)) {
MessageBox(_T("读取文件出错!")); return; } CString strLine;//存放每一行文本 CStringArray strSplit; CeHuiLei SplitString; int hang = 0; while (sf.ReadString(strLine)) { if (hang <= 1) {
原理、算法流程、源程序、计算结果、结果分析、心得体会等。
三.实验所用到的数学公式及程序计算步骤。
单张影像的空间后方交会:利用已知地面控制点数据及相应像点坐标 根据共线方程反 求影像的外方位元素。 数学模型:共线条件方程式:
3
求解过程: (1)获取已知数据。从航摄资料中查取平均航高与摄影机主距;获取控制点的地面测
四.程序流程图。
4
五.程序的主要源代码如下所示:
#include "stdafx.h" #include "KongJianHouFangJiaoHui.h" #include "KongJianHouFangJiaoHuiDlg.h" #include "afxdialogex.h" #include "CeHuiLei.h" #include "Matrix.h" #ifdef _DEBUG #define new DEBUG_NEW #endif CMatrix xyXYZ; int ButtonClickeOk = 1; void CKongJianHouFangJiaoHuiDlg::OnBnClickedOk()//导入测量数据按钮代码 {
hang++; } else {
int n = SplitString.SplitStringArray(strLine, ',', strSplit); int index1 = m_list1.GetItemCount(); CString sno1 = _T(""); sno1.Format(_T("%d"), index1 + 1); int nR = m_list1.InsertItem(index1, sno1); for (int i = 0; i < 5; i++) {
A(2*i, 3) = y*sin(w) - (x*(x*cos(k) - y*sin(k)) / f + f*cos(k))*cos(w); A(2*i, 4) = -f*sin(k) - x*(x*sin(k) + y*cos(k)) / f; A(2*i, 5) = y; A(2*i+1, 0) = (a2*f + a3*y) / _Z; A(2 * i + 1, 1) = (b2*f + b3*y) / _Z; A(2 * i + 1, 2) = (c2*f + c3*y) / _Z; A(2 * i + 1, 3) = -x*sin(w) - (y*(x*cos(k) - y*sin(k)) / f - f*sin(k))*cos(w); A(2 * i + 1, 4) = -f*cos(k) - y/ f*(x*sin(k) + y*cos(k)); A(2 * i + 1, 5) = -x; } return A; }
地球科学与环境工程学院 摄影测量实验报告书
课 程 名: 《摄影测量学》 学 号: 姓 名: 指导老师: 日 期: 2017 年 4 月 20 日
1
目录
一、目的与要求 ................................................ 3 二、实验数据与实验准备 ........................................ 3 三、程序流程图 ................................................ 4 四、实验代码 .................................................. 5 五、实验结果截图 ............................................. 10 六、实验总结 ................................................. 11