?? adjpro0504.cpp
字號:
flag++;
}while(flag<a.obnum);
end.X0=a.Pt[1].X0;end.Y0=a.Pt[1].Y0;
end.X=end.Y=-PI;
// 計算改正后起始坐標方位角與長度
a.L[obi].dist0*=dist(a.Pt[0],a.Pt[1])/dist(a.Pt[0],end);
double ddf=afa(a.Pt[0],a.Pt[1])-afa(a.Pt[0],end);
if(ddf<0) ddf+=2*PI;
a.L[obi].A0+=h_d(ddf);
// 各點近似坐標歸零
for(i=1;i<a.allpnum;i++)
if(i<a.fixpnum){a.Pt[i].X0=a.Pt[i].X;a.Pt[i].Y0=a.Pt[i].Y;}
else a.Pt[i].X0=a.Pt[i].Y0=-PI;
// 近似方位角歸零
for(i=0;i<a.obnum+a.fixafn+a.fixdisn;i++)
if(i!=obi)a.L[i].A0=-PI;
zheng(a.L[obi]); // 坐標正算
// 計算改正后近似方位角
do{flag=0;
statangc(a);
for(int i=0;i<a.obnum;i++)
if(a.L[i].A0==-PI){flag=1;break;}
}while(flag==1);
// 計算改正后各點的近似坐標
do{flag=0;
for(int i=0;i<a.obnum-1;i++)
for(int j=i+1;j<a.obnum;j++)
XY0ang(a.L[i],a.L[j]);
for(i=0;i<a.allpnum;i++)
if(a.Pt[i].X0==-PI)
{ flag=1;break; }
}while(flag==1);
// for(i=0;i<a.allpnum;i++)
// cout<<a.Pt[i].name<<" "<<a.Pt[i].X0<<" "<<a.Pt[i].Y0<<" "
// <<a.Pt[i].X<<" "<<a.Pt[i].Y<<" "<<endl;
return 1;
}
//************************計算控制網未知點的近似坐標**********************************************
int setx0y0(XYnet &a)
{
int n=a.obnum+a.fixdisn+a.fixafn;
// 1.計算近似坐標、近似邊長確定的方位角與邊長
int t(0);
do{
for(int i=0;i<n;i++)
{
if(a.L[i].startp->X0!=-PI && a.L[i].endp->X0!=-PI && a.L[i].A0==-PI)
{
//1.1 近似坐標確定的邊的方位角
a.L[i].A0=h_d(afa(*(a.L[i].startp),*(a.L[i].endp)));
for(int k=0;k<n;k++)
{
if(a.L[i].startp==a.L[k].endp && a.L[i].endp==a.L[k].startp && a.L[k].A0==-PI)
{
if(d_h(a.L[i].A0)-PI>=0)
a.L[k].A0=h_d(d_h(a.L[i].A0)-PI);
else a.L[k].A0=h_d(d_h(a.L[i].A0)+PI);
}
if(a.L[i].startp==a.L[k].startp && a.L[i].endp==a.L[k].endp && a.L[k].A0==-PI)
a.L[k].A0=a.L[i].A0;
}
}
//1.2 邊長觀測值計算的近似邊長
if(a.L[i].style==2 || a.L[i].style==4)
for(int k=0;k<n;k++)
if(a.L[i].startp==a.L[k].endp && a.L[i].endp==a.L[k].startp && a.L[k].dist0==-PI
|| a.L[i].startp==a.L[k].startp && a.L[i].endp==a.L[k].endp && a.L[k].dist0==-PI)
a.L[k].dist0=a.L[i].dist0;
//1.3 由已知近似方位角推算未知近似方位角
if(a.L[i].A0!=-PI)
for(int k=0;k<n;k++)
{
if(a.L[i].startp==a.L[k].endp && a.L[i].endp==a.L[k].startp && a.L[k].A0==-PI)
{
if(d_h(a.L[i].A0)-PI>=0)
a.L[k].A0=h_d(d_h(a.L[i].A0)-PI);
else a.L[k].A0=h_d(d_h(a.L[i].A0)+PI);
}
if(a.L[i].startp==a.L[k].startp && a.L[i].endp==a.L[k].endp && a.L[k].A0==-PI )
a.L[k].A0=a.L[i].A0;
}
}
// 2. 以測站計算觀測值近似方位角
statangc(a);
// 3.近似坐標計算
// 3.1 兩方位角交會法
for(i=0;i<n-1;i++)
for(int j=i+1;j<n;j++)
XY0ang(a.L[i],a.L[j]);
// 3.2 坐標正算法
for(i=0;i<a.obnum+a.fixafn+a.fixdisn;i++)
zheng(a.L[i]);
// 3.3 角度后方交會法
for(i=0;i<a.obnum-2;i++)
for(int j=i+1;j<a.obnum-1;j++)
for(int k=j+1;k<a.obnum;k++)
houj(a.L[i],a.L[j],a.L[k]);
// 3.4邊長后方交會法
for(i=0;i<a.obnum+a.fixafn+a.fixdisn-2;i++)
for(int j=i+1;j<a.obnum+a.fixafn+a.fixdisn-1;j++)
for(int k=j+1;k<a.obnum+a.fixafn+a.fixdisn;k++)
XY0dist(a.L[i],a.L[j],a.L[k]);
// 3.5 無定向導線
Udxdsetx0y0(a);
// 4. 迭代條件的判斷
t=0;
for( i=0;i<a.allpnum;i++)
{
// cout<<a.Pt[i].name<<" "<<a.Pt[i].X0<<" "<<a.Pt[i].Y0<<endl;
if(a.Pt[i].X0==-PI || a.Pt[i].Y0==-PI)
t=1;
}
}while(t==1);
cout<<"---------------------------------------------------------------------------"<<endl;
cout<<" 未知點近似坐標 "<<endl;
cout<<"---------------------------------------------------------------------------"<<endl;
for(int ii=a.fixpnum; ii<a.allpnum;ii++)
cout<<" "<<a.Pt[ii].name<<" "<<setf(a.Pt[ii].X0,3)<<" "<<setf(a.Pt[ii].Y0,3)<<endl;
cout<<"---------------------------------------------------------------------------"<<endl;
return 1;
}
//********************************設置平面網平差的A,P,L*******************************************
int setXYadj(XYnet &a)
{
// 1 確定平差網基本信息
a.aa.m=a.obnum+a.fixafn+a.fixdisn;
a.aa.n=2*a.allpnum+a.statnum;
// 2 確定觀測值的權鎮P
for(int i=0;i<a.aa.m;i++)
{
for(int j=0;j<a.aa.m;j++)
a.aa.P[i][j]=0;
if(a.L[i].style==1) // 設置角度觀測值的權
{
a.aa.P[i][i]=1;
}
if(a.L[i].style==2) // 設置距離觀測值的權
{
double mlml=a.msa+a.msb/1000000*a.L[i].dist;
mlml*=mlml;
a.aa.P[i][i]=a.mangle/rou*a.mangle/rou/mlml;
}
if(a.L[i].style==3 || a.L[i].style==4)// 設置固定距離與方位角觀測值的權
a.aa.P[i][i]=1000000;
}
// matdis(a.aa.P,a.aa.m,a.aa.m);
// 3 確定誤差方程系數陣A
for(i=0;i<a.aa.m;i++)
{
double s=dist(*(a.L[i].startp),*(a.L[i].endp));
// cout<<s<<endl;
double ss=s*s;
double dy=a.L[i].endp->Y0-a.L[i].startp->Y0;
double dx=a.L[i].endp->X0-a.L[i].startp->X0;
double aki=dy/ss; double bki=-1*dx/ss;
for(int j=0;j<a.aa.n;j++)
{a.aa.A[i][j]=0;
if(a.L[i].style==1) // 方向觀測值誤差方程系數
{
a.aa.A[i][2*a.allpnum+a.L[i].sti]=-1; // 測站定向角未知數系數
a.aa.A[i][2*(a.L[i].startp->i)]=aki;
a.aa.A[i][2*(a.L[i].startp->i)+1]=bki;
a.aa.A[i][2*(a.L[i].endp->i)]=-aki;
a.aa.A[i][2*(a.L[i].endp->i)+1]=-bki;
}
if(a.L[i].style==2 || a.L[i].style==4) // 邊長觀測值、固定邊長誤差方程系數
{
a.aa.A[i][2*(a.L[i].startp->i)]=bki*s;
a.aa.A[i][2*(a.L[i].startp->i)+1]=-aki*s;
a.aa.A[i][2*(a.L[i].endp->i)]=-bki*s;
a.aa.A[i][2*(a.L[i].endp->i)+1]=aki*s;
}
if(a.L[i].style==3) // 固定方位角誤差方程系數
{
a.aa.A[i][2*a.L[i].startp->i]=aki;
a.aa.A[i][2*a.L[i].startp->i+1]=bki;
a.aa.A[i][2*a.L[i].endp->i]=-aki;
a.aa.A[i][2*a.L[i].endp->i+1]=-bki;
}
}
}
// matdis(a.aa.A,a.aa.m,a.aa.n);
// 4 確定誤差方程常數項L
for(i=0;i<a.aa.m;i++)
{
a.L[i].A0=h_d(afa(*(a.L[i].startp),*(a.L[i].endp)));
a.L[i].dist0=dist(*(a.L[i].startp),*(a.L[i].endp));
if(a.L[i].style==3)
{
a.aa.l[i][0]=d_h(a.L[i].A-a.L[i].A0);
}
if(a.L[i].style==2 || a.L[i].style==4)
{
a.aa.l[i][0]=a.L[i].dist-a.L[i].dist0;
}
}
//方向觀測值誤差方程常數項
int n=0;
for(int j=0;j<a.statnum;j++)
{
for(int i=n;i<n+a.st[j].aglnum;i++)
{
double scd=d_h(a.L[i].A0)-d_h(a.L[n].A0);
if(scd<0) scd+=2*PI;
a.aa.l[i][0]=d_h(a.L[i].angle)-scd;
}
n+=a.st[j].aglnum+a.st[j].disnum;
}
return 1;
}
//********************************控制網平差計算**************************************************
int doXYadj(XYnet &a)
{
// 6 平差計算
int r0=a.fixafn+a.fixdisn;
if(!doadj(a.aa,2*a.fixpnum,r0))return 0; // 其中第二項為控制點數與測站數之和
cout<<"---------------------------------------------------------------------------"<<endl;
cout<<" 未知數改正數 "<<endl;
cout<<"---------------------------------------------------------------------------"<<endl;
matdis(a.aa.X,a.aa.n);
cout<<"---------------------------------------------------------------------------"<<endl;
cout<<" 驗后單位權中誤差："<<a.aa.m0*206264.8<<" 秒 "<<endl;
cout<<"---------------------------------------------------------------------------"<<endl;
cout<<"---------------------------------------------------------------------------"<<endl;
for(int i=a.fixpnum;i<a.allpnum;i++)
{
a.Pt[i].X=a.Pt[i].X0+a.aa.X[2*i][0];
a.Pt[i].Y=a.Pt[i].Y0+a.aa.X[2*i+1][0];
a.Pt[i].X0+=a.aa.X[2*i][0];
a.Pt[i].Y0+=a.aa.X[2*i+1][0];
cout.precision(11);cout.width(10);
cout<<" "<<a.Pt[i].name<<" "<<setf(a.Pt[i].X,4)<<" "<<setf(a.Pt[i].Y,4)<<endl;
}
cout<<"---------------------------------------------------------------------------"<<endl;
cout<<" 未知點坐標的協方差陣:"<<endl<<endl;
cout<<"---------------------------------------------------------------------------"<<endl;
for(i=2*a.fixpnum;i<2*a.allpnum;i++)
{ cout<<" ";
for(int j=2*a.fixpnum;j<2*a.allpnum;j++)
{
cout<<a.aa.m0*a.aa.m0*a.aa.QXX[i][j]<<" ";
}cout<<endl;
}
cout<<"---------------------------------------------------------------------------"<<endl;
cout<<" 觀測值改正數:"<<endl<<endl;
cout<<"---------------------------------------------------------------------------"<<endl;
matdis(a.aa.V,a.aa.m);
cout<<"---------------------------------------------------------------------------"<<endl;
//7 點位中誤差與誤差橢圓元素計算
for(i=a.fixpnum;i<a.allpnum;i++)
{
// 坐標中誤差計算
a.Pt[i].mX=a.aa.m0*sqrt(a.aa.QXX[2*i][2*i]);
a.Pt[i].mY=a.aa.m0*sqrt(a.aa.QXX[2*i+1][2*i+1]);
// 誤差橢圓參數計算
double K=sqrt((a.aa.QXX[2*i][2*i]-a.aa.QXX[2*i+1][2*i+1])
*(a.aa.QXX[2*i][2*i]-a.aa.QXX[2*i+1][2*i+1])
+4*a.aa.QXX[2*i][2*i+1]*a.aa.QXX[2*i][2*i+1]);
a.Pt[i].mp=sqrt(a.Pt[i].mX*a.Pt[i].mX+a.Pt[i].mY*a.Pt[i].mY);
a.Pt[i].E=sqrt(1/2.0*a.aa.m0*a.aa.m0*(a.aa.QXX[2*i][2*i]+a.aa.QXX[2*i+1][2*i+1]+K));
a.Pt[i].F=sqrt(1/2.0*a.aa.m0*a.aa.m0*(a.aa.QXX[2*i][2*i]+a.aa.QXX[2*i+1][2*i+1]-K));
double F2=atan2(2*a.aa.QXX[2*i][2*i+1],a.aa.QXX[2*i][2*i]-a.aa.QXX[2*i+1][2*i+1]);
if(F2<=0) F2+=2*PI;F2/=2.0;
if(a.aa.QXX[2*i][2*i+1]>0 && F2>PI/2.0 && F2<PI
|| a.aa.QXX[2*i][2*i+1]>0 && F2>PI*3/2.0 && F2<2*PI)
F2+=PI/2.0;
if(a.aa.QXX[2*i][2*i+1]<0 && F2<PI/2.0 && F2>0
|| a.aa.QXX[2*i][2*i+1]<0 && F2>PI && F2<PI*3/2.0)
F2+=PI/2.0;
if(F2>=2*PI) F2-=2*PI;
a.Pt[i].T=h_d(F2);
}
return 1;
}
//*******************************平面網的文件輸入*************************************************
int finXYnet(XYnet &a,char *fname) // 文件輸入平面網函數
{
ifstream in(fname,ios::nocreate); // 建立文件流,并與輸入文件名建立關聯
if(!in) {cout<<fname<<" error: file does not exist! "<<endl; return 0;}
// 文件現實性判斷
// 1 輸入網的基本信息
char name[20];
in>>a.netname;
in>>a.allpnum;
in>>a.fixpnum;
in>>a.statnum;
// in>>a.obnum;
a.obnum=0;
in>>a.fixdisn;
in>>a.fixafn;
in>>a.mangle;
in>>a.msa;
in>>a.msb;
int n(a.fixpnum); // n為已輸入名字的點的個數
// 2輸入控制點信息
for(int i=0;i<a.fixpnum;i++)
{
in>>a.Pt[i].name>>a.Pt[i].X>>a.Pt[i].Y;
a.Pt[i].fixed=1; // 控制點標記
a.Pt[i].X0=a.Pt[i].X;
a.Pt[i].Y0=a.Pt[i].Y;
a.Pt[i].mX=a.Pt[i].mY=0;
a.Pt[i].i=i; // 控制點編號，從0到a.fixpnum-1
}
// 3 輸入未知點相關信息（名字在后面輸入）
for(i=a.fixpnum;i<a.allpnum;i++)
{
a.Pt[i].fixed=0; // 未知點標記
a.Pt[i].X=a.Pt[i].Y=-PI;
a.Pt[i].X0=a.Pt[i].Y0=-PI;
a.Pt[i].mX=a.Pt[i].mY=0;
*(a.Pt[i].name)=0;
a.Pt[i].i=i; // 為未知點編號,從a.fixpnum到a.allpnum-1
}
a.anglenum=0;
a.distnum=0;
// 4 輸入測站及觀測值
int obsnum(0);
for(i=0;i<a.statnum;i++)
{
// 4-1 輸入測站信息
int t=0; // 點名比較標志
in>>name; // 輸入測站名
for(int k=0;k<n;k++)
if(strnicmp(name,a.Pt[k].name,20)==0)
{
a.st[i].stp=&(a.Pt[k]); // 找到同名點，起點指針指向該點
t++; // 找到標志
}
if(t==0) {strcpy(a.Pt[n].name,name);
a.st[i].stp=&(a.Pt[n]); // 找不到同名點，該名輸給新點
n++;}
in>>a.st[i].obnum;
a.st[i].aglnum=0; // 測站角度觀測個數
a.st[i].disnum=0; // 測站距離觀測個數
a.st[i].i=i;
?? 快捷鍵說明
復制代碼
Ctrl + C
搜索代碼
Ctrl + F
全屏模式
F11
切換主題
Ctrl + Shift + D
顯示快捷鍵
?
增大字號
Ctrl + =
減小字號
Ctrl + -