绝对定向作业报告Word文档格式.docx
《绝对定向作业报告Word文档格式.docx》由会员分享，可在线阅读，更多相关《绝对定向作业报告Word文档格式.docx（14页珍藏版）》请在冰豆网上搜索。
4作业过程
4．1模型点坐标计算
根据上次相对定向作业，正确求解出相对定向元素后，利用空间前方交会共识计算出模型点的坐标。
任一模型点坐标：
4．2求模型点的摄影测量坐标：
为了后续计算，把上步求得坐标平移到摄影测量坐标系中，同时放大模型比例尺，使之接近实地大小。
4．3根据所给绝对定向元素重新计算旋转矩阵R，以及上述坐标共同代入绝对定向基本关系式，求解各点的地面摄影测量坐标（X,Y,Z）。
4．4将所求得得各点的地面摄影测量坐标和模型点的摄影测量坐标作为已知条件，绝对定向元素视为未知，进行绝对定向。
4．5确定相对定向元素的初始值；
ψ0=ω0=κ0=0，Xs0=Ys0=Zs0=0，λ0=1
4．6根据确定的初始值（或新的近似值），计算出误差方程式的常数项。
4．7逐点组成误差方程并法化，逐点法化。
4．8解求法方程，得七个绝对定向元素的改正数。
4．8计算绝对定向元素的新值。
4．9判断绝对定向元素的改正数是否小于限值0.00003rad，如满足条件，则结束相对定向计算。
否则重复4.5~4.9。
4．10根据求得的绝对定向元素，将所有模型点的摄测坐标转换为地面摄测坐标。
5源程序
#include<
iostream.h>
fstream.h>
stdlib.h>
iomanip.h>
math.h>
constintN=7;
inti,j,k;
doubleb,x[6][3],y[6][3],z[6][3]={0};
//x[6][3]和y[6][3]分别为左右片像点像空间坐标，z[6][3]为右片各点像空间辅助坐标
//求转置矩阵
template<
typenameT1,typenameT2>
voidTranspose（T1*mat1,T2*mat2,inta,intb）
{
for（i=0;
i<
a;
i++）
for（j=0;
j<
b;
j++）
mat2[j][i]=mat1[i][j];
return;
}
//求矩阵的乘积
voidArray_mul（T1*mat1,T2*mat2,T2*result,inta,intb,intc）
{inti,j,k;
{for（j=0;
c;
{result[i][j]=0;
for（k=0;
k<
k++）
result[i][j]+=mat1[i][k]*mat2[k][j];
}
}
return;
//求逆矩阵
voidswap（double*a,double*b）{doublec;
c=*a;
*a=*b;
*b=c;
};
Inverse（doubleA[N][N],intn）
inti,j,k;
doubled;
intJS[N],IS[N];
for（k=0;
n;
{
d=0;
for（i=k;
for（j=k;
j++）{
if（fabs（A[i][j]）>
d）{
d=fabs（A[i][j]）;
IS[k]=i;
JS[k]=j;
};
if（d+1.0==1.0）return0;
if（IS[k]!
=k）for（j=0;
j++）swap（&
A[k][j],&
A[IS[k]][j]）;
if（JS[k]!
=k）for（i=0;
i++）swap（&
A[i][k],&
A[i][JS[k]]）;
A[k][k]=1/A[k][k];
for（j=0;
j++）if（j!
=k）A[k][j]=A[k][j]*A[k][k];
for（i=0;
i++）if（i!
=k）A[i][j]=A[i][j]-A[i][k]*A[k][j];
=k）A[i][k]=-A[i][k]*A[k][k];
for（k=n-1;
k>
=0;
k--）
j++）if（JS[k]!
=k）swap（&
A[JS[k]][j]）;
i++）if（IS[k]!
A[i][IS[k]]）;
return1;
}
//原始数据导入
voidInput（）
{doublem;
j=0;
i=0;
ifstreamf1（"
左片各点像空间坐标.txt"
）;
if（!
f1）
cerr<
<
"
左片各点像空间坐标.txtfilenotopen!
endl;
exit
（1）;
while（f1>
>
m）
x[i][j++]=m;
if（j>
=3）
{i++;
f1.close（）;
cout<
左片各点像空间坐标坐标为:
6;
3;
{
if（j%3==0）
cout<
setw（8）<
x[i][j];
i=0;
ifstreamf2（"
右片各点像空间坐标.txt"
f2）
右片像各点空间坐标.txtfilenotopen!
while（f2>
y[i][j++]=m;
f2.close（）;
右片各点像空间坐标坐标为:
y[i][j];
b=0.000527;
voidmain（）
doublea[5]={0},d[5]={0},R[3][3],N[2][6],A[6][7]={0},AT[7][6]={0},l[6][1]={0},ATA[7][7]={0},ATl[7][1]={0},DG[7][1]={0};
intt=0;
Input（）;
do{
t++;
5;
a[i]=a[i]+d[i];
//计算旋转矩阵
R[0][0]=cos（a[0]）*cos（a[2]）-sin（a[0]）*sin（a[1]）*sin（a[2]）;
R[0][1]=-cos（a[0]）*sin（a[2]）-sin（a[0]）*sin（a[1]）*cos（a[2]）;
R[0][2]=-sin（a[0]）*cos（a[1]）;
R[1][0]=cos（a[1]）*sin（a[2]）;
R[1][1]=cos（a[1]）*cos（a[2]）;
R[1][2]=-sin（a[1]）;
R[2][0]=sin（a[0]）*cos（a[2]）+cos（a[0]）*sin（a[1]）*sin（a[2]）;
R[2][1]=-sin（a[0]）*sin（a[2]）+cos（a[0]）*sin（a[1]）*cos（a[2]）;
R[2][2]=cos（a[0]）*cos（a[1]）;
//计算右片各点空间辅助坐标
for（j=0;
z[i][j]=R[j][0]*y[i][0]+R[j][1]*y[i][1]+R[j][2]*y[i][2];
N[0][i]=（b*z[i][2]-b*a[4]*z[i][0]）/（x[i][0]*z[i][2]-z[i][0]*x[i][2]）;
N[1][i]=（b*x[i][2]-b*a[4]*x[i][0]）/（x[i][0]*z[i][2]-z[i][0]*x[i][2]）;
l[i][0]=N[0][i]*x[i][1]-N[1][i]*z[i][1]-b*a[3];
i++）
A[i][0]=-z[i][0]*z[i][1]/z[i][2]*N[1][i];
A[i][1]=-（z[i][2]+z[i][1]*z[i][1]/z[i][2]）*N[1][i];
A[i][2]=z[i][0]*N[1][i];
A[i][3]=b;
A[i][4]=-z[i][1]*b/z[i][2];
Transpose（A,AT,6,5）;
Array_mul（AT,A,ATA,5,6,5）;
Inverse（ATA,5）;
Array_mul（AT,l,ATl,5,6,1）;
Array_mul（ATA,ATl,DG,5,5,1）;
d[i]=DG[i][0];
}while（（fabs（d[0]）>
0.00003）||（fabs（d[1]）>
0.00003）||（fabs（d[2]）>
0.00003）||（fabs（d[3]）>
0.00003）||（fabs（d[4]）>
0.00003））;
相对定向元素计算迭代次数为:
t<
未知数的解为:
ψ="
a[0]<
"
ω="
a[1]<
κ="
a[2]<
;
μ="
a[3]<
γ="
a[4]<
."
ofstreamf3（"
计算结果.txt"
f3）
{cerr<
计算结果.txtfilenotopen!
exit
（1）;
f3<
迭代次数为:
相对定向元素的值为:
bv="
b*a[3]<
bw="
b*a[4]<
//计算出个像点的模型坐标mzb[6][3]、摄影测量坐标mszb[6][3]（m=37000），地面摄影测量坐标ds[6][3]
//采用如下绝对定向元素：
缩放系数k＝.00156,ψ＝.0527,ω＝.1426,κ=0.2478,X0=6385.067,Y0=1954.325,Z0=724.215，计算出各点的地面摄影测量坐标；
doublemzb[6][3]={0},mszb[6][3]={0},ds[6][3]={0},m=37000.0,f=0.024,k=1.00156,X0=6385.067,Y0=1954.325,Z0=724.215;
a[0]=0.0527;
a[1]=0.1426;
a[2]=0.2478;
{
mzb[i][0]=N[0][i]*x[i][0];
mzb[i][1]=0.5*（N[0][i]*x[i][1]+N[1][i]*y[i][1]+b*a[3]）;
mzb[i][2]=N[0][i]*x[i][2];
mszb[i][0]=mzb[i][0]*m;
mszb[i][1]=mzb[i][1]*m;
mszb[i][2]=m*f+mzb[i][2]*m;
6个像点的模型坐标分别为：
mzb[i][0]<
"
mzb[i][1]<
mzb[i][2]<
6个像点的摄影测量坐标分别为：
mszb[i][0]<
mszb[i][1]<
mszb[i][2]<
//计算各点的地面摄影测量坐标
ds[i][0]=k*（R[0][0]*mszb[i][0]+R[0][1]*mszb[i][1]+R[0][2]*mszb[i][2]）+X0;
ds[i][1]=k*（R[1][0]*mszb[i][0]+R[1][1]*mszb[i][1]+R[1][2]*mszb[i][2]）+Y0;
ds[i][2]=k*（R[2][0]*mszb[i][0]+R[2][1]*mszb[i][1]+R[2][2]*mszb[i][2]）+Z0;
各点的地面摄影测量坐标为:
ds[i][0]<
ds[i][1]<
ds[i][2]<
计算得各点的地面摄影测量坐标为:
XYZ"
f3.close（）;
//根据地面摄影坐标和模型点的摄影坐标，实现绝对定向，计算出绝对定向元素和各点的地面摄影测量坐标
doublejX[7][1],dX[7],jA[18][7],jAT[7][18],jATA[7][7],jL[18][1],jATL[7][1],jG[7][1],jV[18][1],jAX[18][1];
//重心化坐标
doubleXp=0,Yp=0,Zp=0,Xmp=0,Ymp=0,Zmp=0;
Xp=Xp+mszb[i][0];
Yp=Yp+mszb[i][1];
Zp=Zp+mszb[i][2];
Xmp=Xmp+ds[i][0];
Ymp=Ymp+ds[i][1];
Zmp=Zmp+ds[i][2];
Xp=Xp/6;
Yp=Yp/6;
Zp=Zp/6;
Xmp=Xmp/6;
Ymp=Ymp/6;
Zmp=Zmp/6;
mszb[i][0]=mszb[i][0]-Xp;
mszb[i][1]=mszb[i][1]-Yp;
mszb[i][2]=mszb[i][2]-Zp;
ds[i][0]=ds[i][0]-Xmp;
ds[i][1]=ds[i][1]-Ymp;
ds[i][2]=ds[i][2]-Zmp;
t=0;
do{
7;
jX[i][0]=jX[i][0]+dX[i];
//计算旋转矩阵
R[0][0]=cos（jX[4][0]）*cos（jX[6][0]）-sin（jX[4][0]）*sin（jX[5][0]）*sin（jX[6][0]）;
R[0][1]=-cos（jX[4][0]）*sin（jX[6][0]）-sin（jX[4][0]）*sin（jX[5][0]）*cos（jX[6][0]）;
R[0][2]=-sin（jX[4][0]）*cos（jX[5][0]）;
R[1][0]=cos（jX[5][0]）*sin（jX[6][0]）;
R[1][1]=cos（jX[5][0]）*cos（jX[6][0]）;
R[1][2]=-sin（jX[5][0]）;
R[2][0]=sin（jX[4][0]）*cos（jX[6][0]）+cos（jX[4][0]）*sin（jX[5][0]）*sin（jX[6][0]）;
R[2][1]=-sin（jX[4][0]）*sin（jX[6][0]）+cos（jX[4][0]）*sin（jX[5][0]）*cos（jX[6][0]）;
R[2][2]=cos（jX[4][0]）*cos（jX[5][0]）;
//计算常数项L
jL[i*3][0]=ds[i][0]-jX[3][0]*（R[0][0]*mszb[i][0]+R[0][1]*mszb[i][1]+R[0][2]*mszb[i][2]）-jX[0][0];
jL[i*3+1][0]=ds[i][1]-jX[3][0]*（R[1][0]*mszb[i][0]+R[1][1]*mszb[i][1]+R[1][2]*mszb[i][2]）-jX[1][0];
jL[i*3+2][0]=ds[i][2]-jX[3][0]*（R[2][0]*mszb[i][0]+R[2][1]*mszb[i][1]+R[2][2]*mszb[i][2]）-jX[2][0];
//构造矩阵A
jA[i*3][0]=1;
jA[i*3][1]=0;
jA[i*3][2]=0;
jA[i*3][3]=mszb[i][0];
jA[i*3][4]=-mszb[i][2];
jA[i*3][5]=0;
jA[i*3][6]=-mszb[i][1];
jA[i*3+1][0]=0;
jA[i*3+1][1]=1;
jA[i*3+1][2]=0;
jA[i*3+1][3]=mszb[i][1];
jA[i*3+1][4]=0;
jA[i*3+1][5]=-mszb[i][2];
jA[i*3+1][6]=mszb[i][1];
jA[i*3+2][0]=0;
jA[i*3+2][1]=0;
jA[i*3+2][2]=1;
jA[i*3+2][3]=mszb[i][2];
jA[i*3+2][4]=mszb[i][0];
jA[i*3+2][5]=mszb[i][1];
jA[
升级会员 