摄影测量实习报告单片空间后方交会Word文档下载推荐.docx
- 文档编号:21309554
- 上传时间:2023-01-29
- 格式:DOCX
- 页数:10
- 大小:35.43KB
摄影测量实习报告单片空间后方交会Word文档下载推荐.docx
《摄影测量实习报告单片空间后方交会Word文档下载推荐.docx》由会员分享,可在线阅读,更多相关《摄影测量实习报告单片空间后方交会Word文档下载推荐.docx(10页珍藏版)》请在冰豆网上搜索。
a1=cosΦ*cosκ-sinΦ*sinω*sinκ,a2=-cosΦ*sinκ-sinΦ*sinω*cosκ,a3=-sinΦ*cosω
b1=cosω*sinκ,b2=cosω*cosκ,b3=-sinω
c1=sinΦ*cosκ+cosΦ*sinω*sinκ,c2=-sinΦ*sinκ+cosΦ*sinω*cosκ,c3=cosΦ*cosω
4、逐点计算控制点的像点坐标的近似值,共线方程为:
x=(-f*(a1*(X-Xs)+b1*(Y-Ys)+c1*(Z-Zs)))/(a3*(X-Xs)+b3*(Y-Ys)+c3*(Z-Zs))
y=(-f*(a2*(X-Xs)+b2*(Y-Ys)+c2*(Z-Zs)))/(a3*(X-Xs)+b3*(Y-Ys)+c3*(Z-Zs))
(x)=(-f*(a1*(X-Xs)+b1*(Y-Ys)+c1*(Z-Zs)))/(a3*(X-Xs)+b3*(Y-Ys)+c3*(Z-Zs))
(y)=(-f*(a2*(X-Xs)+b2*(Y-Ys)+c2*(Z-Zs)))/(a3*(X-Xs)+b3*(Y-Ys)+c3*(Z-Zs))
5、组成误差方程:
Vx=a11*dXs+a12*dYs+a13*dZs+a14*dΦ+a15*dω+a16*dκ+(x)-x
Vy=a21*dXs+a22*dYs+a23*dZs+a24*dΦ+a25*dω+a26*dκ+(y)-y
系数求解(是共线方程分别外方位元素求导,是共线方程线性化的系数):
变量代换
A=a1*(X-Xs)+b1*(Y-Ys)+c1*(Z-Zs)
B=a2*(X-Xs)+b2*(Y-Ys)+c2*(Z-Zs)
C=a3*(X-Xs)+b3*(Y-Ys)+c3*(Z-Zs)
a11=(a1*f+a3*x)/C,a12=(b1*f+b3*x)/C,a13=(c1*f+c3*f)/C,a14=y*sinω-((x*(x*cosκ-y*sinκ))/f+f*cosκ)cosω,a15=-f*sinκ-(x*(x*sinκ+y*cosκ)/f),a16=y
a21=(a2*f+a3*x)/C,a22=(b2*f+b3*x)/C,a23=(c2*f+c3*x)/C,a24=-xsinω-((x*(x*cosκ-y*sinκ))/f-f*sinκ)cosω,a25=-f*cosκ-(y*(x*sinκ+y*cosκ)/f),a26=-x
误差方程的常系数(是像点坐标观测值与计算的近似值的差值):
lx=x-(x),ly=y-(y)
6、组成法方程,解求外方位元素改正数X=(ATA)-1ATL(A为误差方程的系数矩阵,L为误差方程的常系数矩阵,通过步骤5求得,此处先求ATA再求矩阵的逆矩阵,解得的改正数加上相应的近似值得到外方位元素新的近似值)
7、检查计算是否收敛:
将求得的外方位元素改正数与规定的限差比较,大于限差继续迭代,小于限差则终止。
四、各子函数详细设计的关键技术参数
子函数(输入函数、Input,矩阵求积Matrixmultiply,计算函数Resection,矩阵转置Matrixtranspose,矩阵求逆Matrixinverse,输出函数Output)主要用了Matrixtranspose,矩阵的行变列,列变行,参数为要转置的矩阵,转置后的矩阵,要转置矩阵的行列数,Matrixinverse,求矩阵的代数余子式,参数有要求逆的矩阵和,逆矩阵的行数,Matrixmultiply,一矩阵的行乘以二矩阵的列,参数为一矩阵,二矩阵,所求矩阵,一的行,一的列,二的列
五、像片外方位元素解算结
六、实习体会
从周一开始进行整个实习框架进行构建,编写了实习的思路既步骤;
然后是梳理框架,要用到的子函数,子函数参数,怎么编写,有不懂的东西网上参考资料;
再是自己动手编程,发现编程要求很细腻,出不得一点差错,程序编写出来后,有很多的错误,要求逐一进行改正,修改;
最后才得以运行。
通过本次实习,我深刻理解了单片空间后方交会原理,进一步熟悉理解了C语言编程,认识到了时间的搭配的重要性和参考资料的必要性,当遇到难题是要迎难而上,达到突破,当完成是能感到一丝丝欣慰和成就感。
附录:
源代码
#include"
stdio.h"
math.h"
Matrixmultiply.c"
Matrixtranspose.c"
Matrixinverse.c"
voidmain()
{
inti,j,k,f=0;
doublex0=0.00018,y0=0.00026,fk=0.15324;
//内方位元素
doublem=40000;
//估算比例尺
doubleR[3][3],XG[6][1],AT[6][8],ATA[6][6],ATL[6][1];
doubleXs=0.0,Ys=0.0,Zs=0.0,Q=0.0,W=0.0,K=0.0;
doubleX,Y,Z,L[8][1],A[8][6];
doubleB[4][5]={-0.08615,-0.06899,36589.41,25273.32,2195.17,-0.05340,0.08221,37631.08,31324.51,728.69,-0.01478,-0.07663,39100.97,24934.98,2386.80,0.01046,0.06443,40426.54,30319.81,757.31};
for(i=0;
i<
4;
i++)
{
Xs=Xs+B[i][2];
Ys=Ys+B[i][3];
Zs=Zs+B[i][4];
}
Xs=Xs/4;
Ys=Ys/4;
Zs=m*fk;
//求得外方位线元素的初始值
do//迭代计算
f++;
//迭代次数
//组成旋转矩阵
R[0][0]=cos(Q)*cos(K)-sin(Q)*sin(W)*sin(K);
R[0][1]=-cos(Q)*sin(K)-sin(Q)*sin(W)*cos(K);
R[0][2]=-sin(Q)*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(Q)*cos(K)+cos(Q)*sin(W)*sin(K);
R[2][1]=-sin(Q)*sin(K)+cos(Q)*sin(W)*cos(K);
R[2][2]=cos(Q)*cos(W);
//计算系数阵和常数项
for(i=0,k=0,j=0;
=3;
i++,k++,j++)
{
X=R[0][0]*(B[i][2]-Xs)+R[1][0]*(B[i][3]-Ys)+R[2][0]*(B[i][4]-Zs);
Y=R[0][1]*(B[i][2]-Xs)+R[1][1]*(B[i][3]-Ys)+R[2][1]*(B[i][4]-Zs);
Z=R[0][2]*(B[i][2]-Xs)+R[1][2]*(B[i][3]-Ys)+R[2][2]*(B[i][4]-Zs);
//为了计算简单而进行变量代换
L[j][0]=B[i][0]-(x0-fk*X/Z);
L[j+1][0]=B[i][1]-(y0-fk*Y/Z);
j++;
A[k][0]=(R[0][0]*fk+R[0][2]*(B[i][0]-x0))/Z;
A[k][1]=(R[1][0]*fk+R[1][2]*(B[i][0]-x0))/Z;
A[k][2]=(R[2][0]*fk+R[2][2]*(B[i][0]-x0))/Z;
A[k][3]=(B[i][1]-y0)*sin(W)-((B[i][0]-x0)*((B[i][0]-x0)*cos(K)-(B[i][1]-y0)*sin(K))/fk+fk*cos(K))*cos(W);
A[k][4]=-fk*sin(K)-(B[i][0]-x0)*((B[i][0]-x0)*sin(K)+(B[i][1]-y0)*cos(K))/fk;
A[k][5]=B[i][1]-y0;
A[k+1][0]=(R[0][1]*fk+R[0][2]*(B[i][1]-y0))/Z;
A[k+1][1]=(R[1][1]*fk+R[1][2]*(B[i][1]-y0))/Z;
A[k+1][2]=(R[2][1]*fk+R[2][2]*(B[i][1]-y0))/Z;
A[k+1][3]=-(B[i][0]-x0)*sin(W)-((B[i][1]-y0)*((B[i][0]-x0)*cos(K)-(B[i][1]-y0)*sin(K))/fk-fk*sin(K))*cos(W);
A[k+1][4]=-fk*cos(K)-(B[i][1]-y0)*((B[i][0]-x0)*sin(K)+(B[i][1]-y0)*cos(K))/fk;
A[k+1][5]=-(B[i][0]-x0);
k++;
}
Matrixtranspose(A,AT,8,6);
//此为转置函数的调用,求AT
Matrixmultiply(AT,A,ATA,6,8,6);
//此为矩阵相乘函数的调用,求ATA
Matrixinverse(ATA,6);
//此为求矩阵逆函数的调用,求ATA的逆
Matrixmultiply(AT,L,ATL,6,8,1);
//此为矩阵相乘函数的调用,求ATL
Matrixmultiply(ATA,ATL,XG,6,6,1);
//此为矩阵相乘函数的调用,求外方位元素改正数
Xs=Xs+XG[0][0];
Ys=Ys+XG[1][0];
Zs=Zs+XG[2][0];
Q=Q+XG[3][0];
W=W+XG[4][0];
K=K+XG[5][0];
//初始值加外方位元素改正数进行迭代
}while(XG[3][0]>
=0.00000001||XG[4][0]>
=0.00000001||XG[5][0]>
=0.00000001);
//当限差满足要求时要再一次进行旋转矩阵的求解
R[0][0]=cos(Q)*cos(K)-sin(Q)*sin(W)*sin(K);
printf("
迭代次数:
%d"
f);
//屏幕输出误差方程系数阵、常数项、改正数
\n\n误差方程系数矩阵A为:
\n\n"
);
for(i=0;
6;
for(j=0;
j<
j++)
printf("
%13.5e"
A[i][j]);
printf("
\n"
\n常数项L为:
8;
1;
L[i][j]);
\n改正数XG为:
XG[i][j]);
\n相片的外方位元素为:
Xs=%13.7e,Ys=%13.7e,Zs=%13.7e\n\n"
Xs,Ys,Zs);
Q=%13.7e,W=%13.7e,K=%13.7e\n"
Q,W,K);
\n旋转矩阵R为:
3;
R[i][j]);
//子函数
#include<
stdlib.h>
math.h>
stdio.h>
intMatrixinverse(a,n)
intn;
doublea[];
{int*is,*js,i,j,k,l,u,v;
doubled,p;
is=malloc(n*sizeof(int));
js=malloc(n*sizeof(int));
for(k=0;
k<
=n-1;
k++)
{d=0.0;
for(i=k;
i<
i++)
for(j=k;
j<
j++)
{l=i*n+j;
p=fabs(a[l]);
if(p>
d){d=p;
is[k]=i;
js[k]=j;
}
if(d+1.0==1.0)
{free(is);
free(js);
err**notinv\n"
return(0);
if(is[k]!
=k)
for(j=0;
{u=k*n+j;
v=is[k]*n+j;
p=a[u];
a[u]=a[v];
a[v]=p;
if(js[k]!
for(i=0;
{u=i*n+k;
v=i*n+js[k];
p=a[u];
}
l=k*n+k;
a[l]=1.0/a[l];
for(j=0;
if(j!
{u=k*n+j;
a[u]=a[u]*a[l];
for(i=0;
if(i!
for(j=0;
if(j!
{u=i*n+j;
a[u]=a[u]-a[i*n+k]*a[k*n+j];
}
for(i=0;
if(i!
{u=i*n+k;
a[u]=-a[u]*a[l];
for(k=n-1;
k>
=0;
k--)
{if(js[k]!
{u=k*n+j;
v=js[k]*n+j;
{u=i*n+k;
v=i*n+is[k];
free(is);
return
(1);
voidMatrixmultiply(a,b,c,m,n,k)
intm,n,k;
doublea[],b[],c[];
{
inti,j,l,u;
m;
for(j=0;
k;
u=i*k+j;
c[u]=0.0;
for(l=0;
l<
n;
l++)
c[u]+=a[i*n+l]*b[l*k+j];
return;
voidMatrixtranspose(a,b,m,n)
intm,n;
doublea[],b[];
inti,j,u;
u=j*m+i;
b[u]=0.0;
b[j*m+i]=a[i*n+j];
return;
}
- 配套讲稿:
如PPT文件的首页显示word图标,表示该PPT已包含配套word讲稿。双击word图标可打开word文档。
- 特殊限制:
部分文档作品中含有的国旗、国徽等图片,仅作为作品整体效果示例展示,禁止商用。设计者仅对作品中独创性部分享有著作权。
- 关 键 词:
- 摄影 测量 实习 报告 单片 空间 后方 交会