间接平差程序设计

1、P间接平差程序设计论文闭合导线网目录编程任务3间接平差原理4平差过程9平差结果10结论27课程设计的体会及建议28参考文献29编程任务L1L6L5L4L3L2S5S4S3S2S1T0BP2P4P3P1导线网示意图： 1）起算数据：T0=100°2030 XB = 100.208 m YB = 120.365 m2） 观测角、观测边长及其观测中误差：L1 = 90°2402± 1.0, S1 = 908.018 m ± 3.4 mm L2 = 220°3055± 3.5, S2 = 533.226 m ± 6.2 mmL3 =

2、 286°3556± 3.1, S3 = 623.836 m ± 4.6 mmL4 = 227°0538± 4.0, S4 = 545.382 m ± 6.8 mmL5 = 230°2810± 3.1, S5 = 609.897 m ± 8.8 mmL6 = 295°1918± 2.2。 要求：1） 计算待定三角点Pi（i = 1.4）的坐标平差值（Xi,Yi）及其方差协方差阵DXX；2） 计算待定三角点Pi的点位中误差；3） 计算观测值的平差值及其方差协方差阵DLL；4） 计算各导线

3、边边长的相对中误差；5） 绘出待定点的点位误差椭圆；6） 验算平差结果的正确性。 间接平差原理² 设某一平差问题中,有n个观测值L，已知其协因数阵 , ，必要观测数为t.² 选定t个独立量为参数X，其估量为 ² 观测值L与改正数V之和为: ，称为观测值平差值。² 按具体平差问题,可列出n个平差值方程为:² 令则平方值方程的矩阵形式为: 令：得误差方程为:² 按最小二乘原理: 按求函数自由极值的方法，得转置后得推得间接平差的基础方程: 解此基础方程，得法方程:上式法方程式可简写成: 其纯量形式为: 解之，得 或将求出的 代入误差方程,即

4、可求得观测值的改正数V: 从而求得平差结果为: 误差方程及其系数按角度坐标平差：按两点距离平差单位权中误差误差椭圆三要素计算：点位误差椭圆形象的反映了控制点在不同方向上的位差，可称为点位精度曲线。根据这个图可以找出待定点坐标平差值在各个方向上的位差，进而进行精度评定。F2=0.5*(x2+y2-(x2-y2)2+4*xy2）tan2*=4*xy2/(x2-y2)E2=0.5*(x2+y2+(x2-y2)2+4*xy2)平差过程一， 该平差问题中有11个观测量，8个必要观测值；二， 所以要设8个未知参量，分别为四个待定点的X,Y坐标. 三， 将每一个观测量的平差值分别表达成所选参数的函

5、数，若函数非线性要将其线性化，列出误差方程 四， 由误差方程系数B和自由项组成法方程五， 法方程个数等于参数的个数8； 六， 解算法方程，求出参数x，y，计算参数的平差值 +=   七， 由误差方程计算V，求出观测量平差值VLL+=Ù八， 评定精度九， 计算待定三角点Pi的点位中误差；a) 计算观测值的平差值及其方差协方差阵DLL；b) 计算各导线边边长的相对中误差；c) 绘出待定点的点位误差椭圆；d) 验算平差结果的正确性计算过程各待定点的近似坐标P1P2P3P4X0992.3141326.018748.195253.71

6、1Y0289.611705.509940.651710.592边号方向X0/mY0/m方位角近似值S0/mX0/ S0(cosa)Y0/S0(sina)l=S-S0/mm1BP1892.106169.24610°4432908.0180.982480.1863902P1P2333.704415.89851°1527533.2260.625820.7799703P2P3-577.823235.142157°5123623.836-0.926240.3769304P3P4-494.484-230.059204°570.84545.382-0.90667-0.

7、4218305BP4153.503590.227255°2518.4609.86152730.251700.9678035.27297角度号观测值L近似值L0l=L-L0/方向编号abL190°240290°24020BP110.04234-.223179L2220°3055220°30550P1P220.30171-.242083L3286°3556286°35560P2P330.1246270.306252L4227°0538227°05380P3P44-0.15953790.34291L5230&#

8、176;2810230°2817.45421-7.45421BP450.327327-0.08513L6295°1918295°1913.54794.45421权阵P1000000000000.082000000000000.104000000000000.0625000000000000.104000000000000.207000000000000.087000000000000.026000000000000.047000000000000.022000000000000.013系数阵B-a1-b1000000a1+a2b1+b2-a1-b20000-a2-

9、b2a3+a2b3+b2-a2-b20000-a3-b3a3+a4b3+b4-a4b40000-a4-b4a4-a5b4-b5-a1-b10000a5b5cosa1sina1000000-cosa2-sina2cosa2sina2000000-cosa3-sina3cosa3sina3000000-cosa4-sina4cosa4sina4000000cosa5sina5-0.042340.22318000000 0.344-0.465-0.3020.2420.0000.0000.0000.000 0.0000.0000.4260.064-0.125-0.306-0.3020.242 0.00

10、00.000-0.125-0.306-0.0350.6490.160-0.343 0.0000.0000.0000.0000.160-0.343-0.486860.42804 -0.0420.2230.0000.0000.0000.0000.32733-0.08513 0.9820.1860.0000.0000.0000.0000.0000.000 -0.626-0.7800.6260.7800.0000.0000.0000.000 0.0000.0000.926-0.377-0.9260.3770.0000.000 0.0000.0000.0000.0000.9070.422-0.907-0

11、.422 0.0000.0000.0000.0000.0000.0000.2520.968观测量观测值改正数v=Bx-l观测值的平差值L190°24020.07326/90°2402.07L2220°3055-2.50499220°3052.5L3286°3556-2.80412286°3553.2L4227°05383.323435227°0541.3L5230°28104.51760230°2814.5L6295°19180.30314295°1918.3S1908.01

12、8 m2.00325/mm908.020S2533.226 m6.16387533.23216S3623.836 m-0.28447623.83572S4545.382 m-5.95462545.3761S5609.897 m-8.82192609.88818源程序#include<iostream>#include<cmath>#include"iomanip"#include<fstream>using namespace std;const double PI=3.14159265358979312;const double n=3

13、600*180/PI;const int v=8;#include<windows.h>/*-矩阵相乘-*/ void mult(double *m1,double *m2,double *result,int i_1,int j_12,int j_2) int i,j,k;for(i=0;i<i_1;i+)for(j=0;j<j_2;j+)resulti*j_2+j=0.0;for(k=0;k<j_12;k+)resulti*j_2+j+=m1i*j_12+k*m2j+k*j_2;return;/*-矩阵求逆-*/void inverse(double cvv)

14、 int i,j,h,k; double p; double qv16; for(i=0;i<v;i+) for(j=0;j<v;j+) qij=cij; for(i=0;i<v;i+) for(j=v;j<16;j+) if(i+8=j) qij=1; else qij=0; for(h=k=0;k<v-1;k+,h+) for(i=k+1;i<v;i+) if(qih=0) continue; p=qkh/qih; for(j=0;j<16;j+) qij*=p; qij-=qkj; for(h=k=v-1;k>0;k-,h-) for(i=

15、k-1;i>=0;i-) if(qih=0) continue; p=qkh/qih; for(j=0;j<16;j+) qij*=p; qij-=qkj;for(i=0;i<v;i+) p=1.0/qii; for(j=0;j<16;j+) qij*=p; for(i=0;i<v;i+) for(j=0;j<v;j+) cij=qij+8;void main()double S5,L6,X5,Y5,T6,B118,l111,a5,b5,c5,s5,m11,x81;double BH811,BHP811,BHPB88,BHPl81,P1111;int i,j

16、,k;S0=908.018;S1=533.226;S2=623.836;S3=545.382;S4=609.897;L0=325442;/角度是以秒的单位输入的L1=793855; L2=1031756;L3=817538;L4=829690; L5=1063158;m0=1.0;m1=3.5;m2=3.1;m3=4.0;m4=3.1;m5=2.2;m6=3.4;m7=6.2;m8=4.6;m9=6.8;m10=8.8;X0=100.208;Y0=120.365; T0=361230;for(i=0;i<=5;i+)/把角度单位化为弧度Li=Li/n;T0=T0/n;for(i=1;i&

17、lt;=4;i+)/求出近似坐标及近似方位角Ti=Ti-1+Li-1-PI;Xi=Xi-1+Si-1*cos(Ti);Yi=Yi-1+Si-1*sin(Ti);/*-误差方程系数-*/for(i=0;i<=3;i+)ci=(Xi+1-Xi)/sqrt(pow(Xi+1-Xi,2)+pow(Yi+1-Yi,2);si=(Yi+1-Yi)/sqrt(pow(Xi+1-Xi,2)+pow(Yi+1-Yi,2);ai=0.001*n*si/sqrt(pow(Xi+1-Xi,2)+pow(Yi+1-Yi,2);bi=-0.001*n*ci/sqrt(pow(Xi+1-Xi,2)+pow(Yi+1-

18、Yi,2);c4=(X4-X0)/sqrt(pow(X4-X0,2)+pow(Y4-Y0,2);s4=(Y4-Y0)/sqrt(pow(X4-X0,2)+pow(Y4-Y0,2);a4=0.001*n*s4/sqrt(pow(X4-X0,2)+pow(Y4-Y0,2);b4=-0.001*n*c4/sqrt(pow(X4-X0,2)+pow(Y4-Y0,2);for(i=0;i<=10;i+)for(j=0;j<=7;j+)Bij=0;B00=-a0;B01=-b0;B10=a0+a1;B11=b0+b1;B12=-a1;B13=-b1;B22=a1+a2;B23=b1+b2;B2

19、4=-a2;B25=-b2;B20=-a1;B21=-b1;B32=-a2;B33=-b2;B34=a2+a3;B35=b2+b3;B36=-a3;B37=-b3;B44=-a3;B45=-b3;B46=a3-a4;B47=b3-b4;B50=-a0;B51=-b0;B56=a4;B57=b4;B60=c0;B61=s0;B70=-c1;B71=-s1;B72=c1;B73=s1;B82=-c2;B83=-s2;B84=c2;B85=s2;B94=-c3;B95=-s3;B96=c3;B97=s3;B106=c4;B107=s4;/*-观测值权阵-*/for(i=0;i<=10;i+)f

20、or(j=0;j<=10;j+)Pij=0;for(i=0;i<=10;i+)Pii=1/pow(mi,2);/*-误差方程常数项-*/l00=-PI+L0+T0-atan(Y1-Y0)/(X1-X0);l10=-PI+L1-atan(Y2-Y1)/(X2-X1)+atan(Y1-Y0)/(X1-X0);l20=-2*PI+L2-atan(Y3-Y2)/(X3-X2)+atan(Y2-Y1)/(X2-X1);l30=-PI+L3-atan(Y4-Y3)/(X4-X3)+atan(Y3-Y2)/(X3-X2);l40=L4+atan(Y4-Y3)/(X4-X3)-atan(Y4-Y0

21、)/(X4-X0)-PI;l50=L5+atan(Y4-Y0)/(X4-X0)-atan(Y1-Y0)/(X1-X0)-2*PI;for(i=0;i<=5;i+)li0=n*li0;for(i=6;i<10;i+)li0=Si-6-sqrt(pow(Xi-6-Xi-5,2)+pow(Yi-6-Yi-5,2);l100=S4-sqrt(pow(X4-X0,2)+pow(Y4-Y0,2);/*-解法方程-*/for(i=0;i<11;i+)/系数阵转置for(j=0;j<8;j+)BHji=Bij;mult(&BH00,&P00,&BHP00,8,1

22、1,11);mult(&BHP00,&B00,&BHPB00,8,11,8);inverse(BHPB);mult(&BHP00,&l00,&BHPl00,8,11,1);mult(&BHPB00,&BHPl00,&x00,8,8,1);/*-单位权中误差-*/double e4,F4,tn4, Bx111,v111;double w=0,f=0;mult(&B00,&x00,&Bx00,11,8,1);for(i=0;i<=10;i+)vi0= Bxi0-li0;w=w+Pii*pow(vi

23、0,2);f=w/3;/*-误差椭圆三参数-*/for(i=0;i<=7;i+)for(j=0;j<=7;j+)BHPBij=f*BHPBij;for(i=0;i<=3;i+)ei=sqrt(0.5*(BHPB2*i2*i+BHPB2*i+12*i+1+sqrt(pow(BHPB2*i2*i-BHPB2*i+12*i+1),2)+4*BHPB2*i2*i+1);Fi=sqrt(0.5*(BHPB2*i2*i+BHPB2*i+12*i+1-sqrt(pow(BHPB2*i2*i-BHPB2*i+12*i+1),2)+4*BHPB2*i2*i+1);tni=0.5*atan(BH

24、PB2*i2*i+1/(BHPB2*i2*i-BHPB2*i+12*i+1);for(i=1;i<=4;i+)Xi=Xi+x2*i-20*0.001;Yi=Yi+x2*i-10*0.001;cout<<"待定点坐标X,Y的坐标平差值分别为"<<endl;cout<<fixed;cout.precision(5);for(i=1;i<=4;i+)cout<<Xi<<setw(15)<<Yi<<endl;system("pause");平差结果坐标平差值的方差协方

25、差阵DXX15.159-0.5857.7112.4235.7552.4252.588-0.155-0.58522.505-4.37923.853-7.17210.698-4.3047.3317.711-4.37942.580-7.80735.6626.66410.0980.1302.42323.853-7.80788.235-20.93834.888-1.14113.6465.755-7.17235.662-20.93861.78318.63613.13521.0322.42510.6986.66434.88818.63677.47410.24064.5702.588-4.30410.098-

26、1.14113.13510.24018.52010.750-0.1557.3310.13013.64621.03264.57010.75098.191观测值的平差值及其方差协方差阵DLL1.159-0.9570.030-0.623-1.1470.3560.182-0.197-0.425-0.3020.926-0.9579.241-2.730-1.013-2.666-1.1603.0164.649-2.646-5.670-0.2120.030-2.7306.562-6.1421.622-1.7321.2112.4642.5921.032-1.258-0.623-1.013-6.14212.293-2.2201.116-0.382-2.7203.3605.1494.216-1.147-2.6661.622-2.2208.2792.2850.091-2.365-1.864-0.3594.9180.356-1.160-1.7321.1162.2855.1890.9462.9370.147-1.744-2.9700.1823.0161.211-0.3820.0910.94615.200-2.5200.9392.5731.131-0.1974.6492.464-2.720-2.3652.937-2.52048.6371.3615.3714.727-0.425-2.6462.5923