近景攝影測(cè)量光束法平差報(bào)告

1、精選優(yōu)質(zhì)文檔-傾情為你奉上近景攝影測(cè)量光束法平差報(bào)告2011 年 6 月 4 日1 作業(yè)目的 -32 外業(yè)控制點(diǎn)的觀測(cè)與解算 -33 近景影像獲取 -44 LPS刺點(diǎn)點(diǎn)位 -45 光束法平差與精度評(píng)定 -56 總結(jié) -111 作業(yè)目的以近景攝影測(cè)量大實(shí)習(xí)為基礎(chǔ)，對(duì)所攝取近景相片解析處理，以外業(yè)控制點(diǎn)的解算成果以及內(nèi)業(yè)LPS平差結(jié)果為依據(jù)，編寫(xiě)光束法平差程序，由22個(gè)控制點(diǎn)的像素坐標(biāo)及5個(gè)“已知控制點(diǎn)”的三維坐標(biāo)求解其余17個(gè)控制點(diǎn)的三維坐標(biāo)，并評(píng)定精度。2 作業(yè)條件及數(shù)據(jù)點(diǎn)號(hào)像素坐標(biāo)列（J）像素坐標(biāo)行（I）XYZ左片：2650.9892114.93497.4532353.7473299.895

2、382792.4912259.531508.8008342.3524298.6832102791.483740.514508.8138342.3548307.0717163928.5592120.49520.2969353.7531300.1146214890.5842130.45527.9857353.5821300.10371648.6242765.5820003660.4521441.4110004728.563816.58500051965.8952557.99600061910.1051210.0700072767.4553044.53100092774.0591493.061000

3、123319.0112665.417000133312.2861986.582000143298.4681284.901000154055.0522705.029000173808.9851539.018000183715.006962.032000193836.444706.426000204883.392691.6510002247541681000234825.4091018.545000右片：2670.9482129.967497.4532353.7473299.895382346.4432264.542508.8008342.3524298.6832102361.448691.079

4、508.8138342.3548307.0717164088.4192115.427520.2969353.7531300.1146205203.4412736.112527.9857353.5821300.10371666.1032764.8820003685.4031472.5740004754.414860.65600051652.4312568.50300061600.0581207.01400072312.0273077.96400092334.4721470.473000123083.1932691.257000133083.0661976.987000143072.4281240

5、.419000154230.5272741.493000173956.4451498.102000183852.033888.02000194075.943613.315000215214.4922119.571000225144.4631507.538000235139.982903.989000外方位元素初始值：(左)Xs1 = 497.9149,Ys1 = 301.2754,Zs1 = 297.2430,Q1=14.9560,W1=4.7765,K1=-0.0308,(右)Xs2 = 509.6501,Ys2 = 301.3448,Zs2 = 297.4727,Q2=2.1777 ,W2

6、=4.5969,K2=0.0791,相片主距：50mm3 平差思想31 共線條件方程本次作業(yè)中，光束法平差基于如下共線條件方程：該共線方程是描述攝影中心S、像點(diǎn)a以及物點(diǎn)A位于一直線上的關(guān)系式。像點(diǎn)a在成像過(guò)程中存在某種系統(tǒng)誤差，其改正數(shù)（x，y）添加在上式左邊，并有：將共線條件方程線性化，其中：32 像點(diǎn)坐標(biāo)誤差方程式由于本作業(yè)采用光束法平差，兩張照片共44個(gè)點(diǎn)（每張22個(gè)），每個(gè)點(diǎn)可列2個(gè)像點(diǎn)坐標(biāo)誤差方程（vx，vy），故總共88個(gè)誤差方程。而每個(gè)誤差方程中：有12個(gè)外方位元素改正（每張相片6個(gè)），6個(gè)內(nèi)方位元素改正（每張相片3個(gè)），51個(gè)加密點(diǎn)坐標(biāo)改正（5個(gè)控制點(diǎn)坐標(biāo)改正為0，故只有有1

7、7個(gè)加密點(diǎn)的三維坐標(biāo)改正：dX，dY,dZ），4個(gè)畸變參數(shù)改正（每張相片2個(gè)，k1，k2）。33 平差由88個(gè)像點(diǎn)坐標(biāo)誤差方程式組成方程組：,其中：V陣為像點(diǎn)坐標(biāo)誤差改正，88行*1列；A陣為系數(shù)陣，88行*73列；X陣偉未知數(shù)陣，73行*1列（未知數(shù)包括12個(gè)外方位元素、6個(gè)內(nèi)方位元素、51個(gè)加密點(diǎn)三維坐標(biāo)、4個(gè)畸變參數(shù)）L陣為常數(shù)項(xiàng)陣，88行*1列由式X=（AtA）-1 ATL解求未知數(shù)即可4 程序#include <iomanip.h>#include <stdlib.h>#include <math.h>#include<fstream.h&g

8、t;#include<iostream.h>const int N=90;/求轉(zhuǎn)置矩陣template<typename T1,typename T2>void Transpose(T1*mat1,T2*mat2,int a,int b)int i,j;for(i=0;i<a;i+)for(j=0;j<b;j+)mat2ji=mat1ij;return;/求矩陣的乘積template<typename T1,typename T2>void Array_mul(T1*mat1,T2 * mat2,T2 * result,int a,int b,i

9、nt c) int i,j,k;for(i=0;i<a;i+)for(j=0;j<c;j+)resultij=0;for(k=0;k<b;k+)resultij+=mat1ik*mat2kj;return;/求逆矩陣inverse(double ANN,int m)int i=0,j=0,k=0;double C2020,B2020;for(i=0;i<2*m;i+) for(j=0;j<2*m;j+) if(i=j) Cij=1.0; else Cij=0.0;for(i=0;i<m;i+) for(j=0;j<m;j+)Bij=Aij;for(i=

10、0;i<m;i+) for(j=m;j<2*m;j+)Bij=Cij-m;cout.precision(5);for(k=0;k<m;k+)for(i=k;i<m;i+)for(j=2*m-1;j>=0;j-)if(Bik!=0)Bij=Bij/Bik; for(i=m-1;i>k;i-)if(Bik=0)continue; for(j=0;j<2*m;j+) Bij=Bij-Bkj; for(k=1;k<m;k+) for(i=0;i<k;i+) for(j=2*m-1;j>=i;j-) Bij=Bij-Bik*Bkj;for(i=

11、0;i<m;i+) for(j=0;j<m;j+)Aji=Bjm+i; return 1; void main()/定義兩張相片共個(gè)控制點(diǎn)和加密點(diǎn)的像素坐標(biāo)（J1,I1,J2,I2）和像平面直角坐標(biāo)（x1,y1,x2,y2），及地面坐標(biāo)（X,Y,Z）double J1N=0.0,I1N=0.0,J2N=0.0,I2N=0.0;double x1N=0.0,x2N=0.0,z1N=0.0,z2N=0.0;double XN=0.0,YN=0.0,ZN=0.0; /導(dǎo)入控制點(diǎn)坐標(biāo)數(shù)據(jù)ifstream infile;infile.open("控制點(diǎn)坐標(biāo)數(shù)據(jù).txt")

12、; if(infile.is_open() while(!infile.eof ()for(int i=0;i<44;i+) infile>>J1i;infile.ignore(1); infile>>I1i;infile.ignore(1); infile>>Xi;infile.ignore(1); infile>>Yi;infile.ignore(1); infile>>Zi;infile.ignore(1); infile.close(); cout<<J11<<" "<&

13、lt;J143<<endl;/像素坐標(biāo)轉(zhuǎn)化為直角坐標(biāo)for (int j = 0; j < 44; j+)x1j = (J1j - 5616/2) * 6. / 1000 ;z1j = (3744/2 - I1j) * 6. / 1000 ;cout<<x11<<" "<<x143<<endl;/左右相片外方位元素初始值及攝影機(jī)主距double Xs1,Ys1,Zs1,Q1,W1,K1,f1,x01,z01,k11,k21;double Xs2,Ys2,Zs2,Q2,W2,K2,f2,x02,z02,k12,

14、k22;Xs1 = 497.9149,Ys1 = 301.2754,Zs1 = 297.2430,Q1=14.9560,W1=4.7765,K1=-0.0308,f1 = 50, x01 = 0, z01 = 0, k11 = 0, k21 = 0;Xs2 = 509.6501,Ys2 = 301.3448,Zs2 = 297.4727,Q2=2.1777 ,W2=4.5969,K2=0.0791,f2 = 50, x02 = 0, z02 = 0, k12 = 0, k22 = 0;double rrN=0.0,dxN=0.0,dzN=0.0;double XX1,YY1,ZZ1,XX2,Y

15、Y2,ZZ2;/定義旋轉(zhuǎn)矩陣，系數(shù)矩陣,常數(shù)項(xiàng)和改正數(shù)double R133=0.0,R233=0.0,A1NN=0.0,l1NN=0.0,dNN=0.0;int t=0;cout<<Xs1<<endl;/組成旋轉(zhuǎn)矩陣 R100=cos(Q1)*cos(K1)-sin(Q1)*sin(W1)*sin(K1); R101=cos(W1)*sin(Q1); R102=-cos(Q1)*sin(K1)-sin(Q1)*sin(W1)*cos(K1); R110=-sin(Q1)*cos(K1)-cos(Q1)*sin(W1)*sin(K1); R111=cos(Q1)*cos

16、(W1); R112=sin(Q1)*sin(K1)-cos(Q1)*sin(W1)*cos(K1); R120=cos(W1)*sin(K1); R121=-sin(W1); R122=cos(W1)*cos(K1); R200=cos(Q2)*cos(K2)-sin(Q2)*sin(W2)*sin(K2); R201=cos(W2)*sin(Q2); R202=-cos(Q2)*sin(K2)-sin(Q2)*sin(W2)*cos(K2); R210=-sin(Q2)*cos(K2)-cos(Q2)*sin(W2)*sin(K2); R211=cos(Q2)*cos(W2); R212=

17、sin(Q2)*sin(K2)-cos(Q2)*sin(W2)*cos(K2); R220=cos(W2)*sin(K2); R221=-sin(W2); R222=cos(W2)*cos(K2); for(int k=0;k<4;k+) XX1=R100*(Xk-Xs1)+R110*(Yk-Ys1)+R120*(Zk-Zs1); YY1=R101*(Xk-Xs1)+R111*(Yk-Ys1)+R121*(Zk-Zs1); ZZ1=R102*(Xk-Xs1)+R112*(Yk-Ys1)+R122*(Zk-Zs1); XX2=R200*(Xk-Xs1)+R210*(Yk-Ys1)+R220

18、*(Zk-Zs1); YY2=R201*(Xk-Xs1)+R211*(Yk-Ys1)+R221*(Zk-Zs1); ZZ2=R202*(Xk-Xs1)+R212*(Yk-Ys1)+R222*(Zk-Zs1); for(int p=0;p<22;p+) rrp=(x1p-x01)*(x1p-x01)+(z1p-z01)*(z1p-z01); dxp=(x1p-x01)*(k11*(x1p-x01)*(x1p-x01)+k21*(x1p-x01)*(x1p-x01)*(x1p-x01)*(x1p-x01); dzp=(z1p-z01)*(k11*(z1p-z01)*(z1p-z01)+k21

19、*(z1p-z01)*(z1p-z01)*(z1p-z01)*(z1p-z01); cout<<rr1<<endl; for(int q=22;q<44;q+) rrq=(x1q-x02)*(x1q-x02)+(z1q-z02)*(z1q-z02); dxq=(x1q-x02)*(k12*(x1q-x02)*(x1q-x02)+k22*(x1q-x02)*(x1q-x02)*(x1q-x02)*(x1q-x02); dzq=(z1q-z02)*(k12*(z1q-z02)*(z1q-z02)+k22*(z1q-z02)*(z1q-z02)*(z1q-z02)*(z

20、1q-z02); cout<<rr22<<endl; /計(jì)算系數(shù)陣(88*73)和常數(shù)項(xiàng) /左片 for(int i=0,H=0;i<22;i+) l1H0=x1i-x01+f1*XX1/YY1-dxi;l1H+10=z1i-z01+f1*ZZ1/YY1-dzi;A1H0=(R101*(x1i-x01)-R100*f1)/YY1;/x/Xs=-x/XA1H1=(R111*(x1i-x01)-R110*f1)/YY1;/x/YsA1H2=(R121*(x1i-x01)-R120*f1)/YY1;/x/ZsA1H3=(z1i-z01)*sin(W1)-(x1i-x01

21、)*(x1i-x01)*cos(K1)-(z1i-z01)*sin(K1)/f1+f1*cos(K1)*cos(W1);/x/QA1H4=-f1*sin(K1)-(x1i-x01)*(x1i-x01)*sin(K1)+(z1i-z01)*cos(K1)/f1;/x/WA1H5=z1i-z01;/x/KA1H6=(x1i-x01)/f1;/x/fA1H7=1;/x/x0A1H8=0;/x/z0A1H9=(x1i-x01)*rri;/x/k1A1H10=(x1i-x01)*rri*rri;/x/k2A1H11=(R201*(x1i-x02)-R200*f2)/YY2;/x/Xs=-x/XA1H12

22、=(R211*(x1i-x02)-R210*f2)/YY2;/x/YsA1H13=(R221*(x1i-x02)-R220*f2)/YY2;/x/ZsA1H14=(z1i-z02)*sin(W2)-(x1i-x02)*(x1i-x02)*cos(K2)-(z1i-z02)*sin(K2)/f2+f2*cos(K2)*cos(W2);/x/QA1H15=-f2*sin(K2)-(x1i-x02)*(x1i-x02)*sin(K2)+(z1i-z02)*cos(K2)/f2;/x/WA1H16=z1i-z02;/x/KA1H17=(x1i-x02)/f2;/x/fA1H18=1;/x/x0A1H1

23、9=0;/x/z0A1H20=(x1i-x02)*rri;/x/k1A1H21=(x1i-x02)*rri*rri;/x/k2A1H+10=(R101*(z1i-z01)-R102*f1)/YY1;/z/Xs=-z/XA1H+11=(R111*(z1i-z01)-R112*f1)/YY1;/z/Ys=-z/YA1H+12=(R121*(z1i-z01)-R122*f1)/YY1;/z/Zs=-z/ZA1H+13=-(x1i-x01)*sin(W1)-(z1i-z01)*(x1i-x01)*cos(K1)-(z1i-z01)*sin(K1)/f1-f1*sin(K1)*cos(W1);/z/QA

24、1H+14=-f1*cos(K1)-(z1i-z01)*(x1i-x01)*sin(K1)+(z1i-z01)*cos(K1)/f1;/z/WA1H+15=-(x1i-x01);/z/KA1H+16=(z1i-z01)/f1;/z/fA1H+17=0;/z/x0A1H+18=1;/z/z0A1H+19=(z1i-z01)*rri;/x/k1A1H+110=(z1i-z01)*rri*rri;/x/k2A1H+111=(R201*(z1i-z02)-R202*f2)/YY2;/z/Xs=-z/XA1H+112=(R211*(z1i-z02)-R212*f2)/YY2;/z/Ys=-z/YA1H+

25、113=(R221*(z1i-z02)-R222*f2)/YY2;/z/Zs=-z/ZA1H+114=-(x1i-x02)*sin(W2)-(z1i-z02)*(x1i-x02)*cos(K2)-(z1i-z02)*sin(K2)/f2-f2*sin(K2)*cos(W2);/z/QA1H+115=-f2*cos(K2)-(z1i-z02)*(x1i-x02)*sin(K2)+(z1i-z02)*cos(K2)/f2;/z/WA1H+116=-(x1i-x02);/z/KA1H+117=(z1i-z02)/f2;/z/fA1H+118=0;/z/x0A1H+119=1;/z/z0A1H+120

26、=(z1i-z02)*rri;/x/k1A1H+121=(z1i-z02)*rri*rri;/x/k2H=H+2; /右片 for(int ii=22,HH=44;ii<44;ii+) l1HH0=x1ii-x01+f1*XX1/YY1-dxii;l1HH+10=z1ii-z01+f1*ZZ1/YY1-dzii;A1HH0=(R101*(x1ii-x01)-R100*f1)/YY1;/x/Xs=-x/XA1HH1=(R111*(x1ii-x01)-R110*f1)/YY1;/x/YsA1HH2=(R121*(x1ii-x01)-R120*f1)/YY1;/x/ZsA1HH3=(z1ii-

27、z01)*sin(W1)-(x1ii-x01)*(x1ii-x01)*cos(K1)-(z1ii-z01)*sin(K1)/f1+f1*cos(K1)*cos(W1);/x/QA1HH4=-f1*sin(K1)-(x1ii-x01)*(x1ii-x01)*sin(K1)+(z1ii-z01)*cos(K1)/f1;/x/WA1HH5=z1ii-z01;/x/KA1HH6=(x1ii-x01)/f1;/x/fA1HH7=1;/x/x0A1HH8=0;/x/z0A1HH9=(x1ii-x01)*rrii;/x/k1A1HH10=(x1ii-x01)*rrii*rrii;/x/k2A1HH11=(R

28、201*(x1ii-x02)-R200*f2)/YY2;/x/Xs=-x/XA1HH12=(R211*(x1ii-x02)-R210*f2)/YY2;/x/YsA1HH13=(R221*(x1ii-x02)-R220*f2)/YY2;/x/ZsA1HH14=(z1ii-z02)*sin(W2)-(x1ii-x02)*(x1ii-x02)*cos(K2)-(z1ii-z02)*sin(K2)/f2+f2*cos(K2)*cos(W2);/x/QA1HH15=-f2*sin(K2)-(x1ii-x02)*(x1ii-x02)*sin(K2)+(z1ii-z02)*cos(K2)/f2;/x/WA1

29、HH16=z1ii-z02;/x/KA1HH17=(x1ii-x02)/f2;/x/fA1HH18=1;/x/x0A1HH19=0;/x/z0A1HH20=(x1ii-x02)*rrii;/x/k1A1HH21=(x1ii-x02)*rrii*rrii;/x/k2A1HH+10=(R101*(z1ii-z01)-R102*f1)/YY1;/z/Xs=-z/XA1HH+11=(R111*(z1ii-z01)-R112*f1)/YY1;/z/Ys=-z/YA1HH+12=(R121*(z1ii-z01)-R122*f1)/YY1;/z/Zs=-z/ZA1HH+13=-(x1ii-x01)*sin(

30、W1)-(z1ii-z01)*(x1ii-x01)*cos(K1)-(z1ii-z01)*sin(K1)/f1-f1*sin(K1)*cos(W1);/z/QA1HH+14=-f1*cos(K1)-(z1ii-z01)*(x1ii-x01)*sin(K1)+(z1ii-z01)*cos(K1)/f1;/z/WA1HH+15=-(x1ii-x01);/z/KA1HH+16=(z1ii-z01)/f1;/z/fA1HH+17=0;/z/x0A1HH+18=1;/z/z0A1HH+19=(z1ii-z01)*rrii;/x/k1A1HH+110=(z1ii-z01)*rrii*rrii;/x/k2A

31、1HH+111=(R201*(z1ii-z02)-R202*f2)/YY2;/z/Xs=-z/XA1HH+112=(R211*(z1ii-z02)-R212*f2)/YY2;/z/Ys=-z/YA1HH+113=(R221*(z1ii-z02)-R222*f2)/YY2;/z/Zs=-z/ZA1HH+114=-(x1ii-x02)*sin(W2)-(z1ii-z02)*(x1ii-x02)*cos(K2)-(z1ii-z02)*sin(K2)/f2-f2*sin(K2)*cos(W2);/z/QA1HH+115=-f2*cos(K2)-(z1ii-z02)*(x1ii-x02)*sin(K2)

32、+(z1ii-z02)*cos(K2)/f2;/z/WA1HH+116=-(x1ii-x02);/z/KA1HH+117=(z1ii-z02)/f2;/z/fA1HH+118=0;/z/x0A1HH+119=1;/z/z0A1HH+120=(z1ii-z02)*rrii;/x/k1A1HH+121=(z1ii-z02)*rrii*rrii;/x/k2HH=HH+2; /計(jì)算左片外方位元素Xs1,Ys1,Zs1,Q,W,Kdouble A1TNN=0,A1TA1NN=0,A1TL1NN=0;Transpose(A1,A1T,N,N);Array_mul(A1T,A1,A1TA1,N,N,N);

33、inverse(A1TA1,N);Array_mul(A1T,l1,A1TL1,N,N,N);Array_mul(A1TA1,A1TL1,d,N,N,N);Xs1=Xs1+d00; Ys1=Ys1+d10; Zs1=Zs1+d20;Q1=Q1+d30; W1=W1+d40; K1=K1+d50;f1=f1+d60; x01=x01+d70; z01=z01+d80; k11=k11+d90; k21=k21+d100;Xs2=Xs2+d110; Ys2=Ys1+d120; Zs2=Zs2+d130;Q2=Q2+d140; W2=W2+d150; K2=K2+d160;f2=f2+d170; x

34、02=x02+d180; z02=z02+d190;k12=k11+d90; k22=k22+d100;cout<<"迭代次數(shù)："<<t<<endl;cout<<"外方位元素的直線元素為"<<endl;cout<<" Xs1 = "<<Xs1<<" , Xs2 = "<<Xs2<<endl;cout<<" Ys1 = "<<Ys1<<&quo

35、t; , Ys2 = "<<Ys2<<endl;cout<<" Zs1 = "<<Zs1<<" , Zs2 = "<<Zs2<<endl;cout<<" Q1 = "<<Q1<<" , Q2 = "<<Q2<<endl;cout<<" W1 = "<<Q1<<" , W2 = "<<W2<<endl; cout<&l