VC++ MFC高斯平均引数大地主题正反算
高斯投影正反算c代码

高斯投影正反算程序设计一.程序设计流程本程序(de)设计思路如下:(1),程序采用VS08版本作为开发平台,并采用C语言作为开发语言,设计为WindowsForm窗体程序形式.(2),程序主要(de)算法来自于教材.但是本程序为了更加实用,添加了更多(de)解算基准,包括:WGS-84,国际椭球1975,克氏椭球,和2000国家大地坐标系.(3),程序为了更方便(de)读取数据和输出数据,故需要自己定义了固定(de)数据输入格式和数据输出格式或形式,请老师注意查看.二.代码using System;using Systusing SystemponentModel;using System.Data;using System.Drawing;using System.Text;namespace Gauss{public partial class Form1 : Form{//大地坐标//Geodetic Coordinatepublic struct CRDGEODETIC{public double dLongitude;public double dLatitude;public double dHeight;}//笛卡尔坐标//Cartesian Coordinatepublic struct CRDCARTESIAN{public double x;public double y;public double z;}public Form1(){InitializeComponent();}private void button1_Click(object sender, EventArgs e) {double ee = 0;double a = 0;string tt;try{}catch{MessageBox.Show("Gauss Inverse: Choose datum error");return;}if (ttpareTo("克氏椭球")==0){a = 6378245.00;}if (ttpareTo("WGS-84") == 0){a = 6378.00;}if (ttpareTo("1975国际椭球") == 0){a = 6378140.00;ee =}if (ttpareTo("2000国家大地坐标系") == 0){a = 6378.0;}const double pai = 3.1415926;double b = Math.Sqrt(a a (1 - ee ee));double c = a a / b;double epp = Math.Sqrt((a a - b b) / b / b);CRDGEODETIC pcrdGeo;CRDCARTESIAN pcrdCar;double midlong;//求纬度string[] temp;double[] tempradius = new double[3];for (int i = 0; i < 3; i++){tempradius[i] = Convert.ToDouble(temp[i]);}pcrdGeo.dLatitude = tempradius[0] / 180.0 pai + tempradius[1] / 180.0 / 60.0 pai + tempradius[2] / 180 / 60.0 / 60 pai;//求经度for (int i = 0; i < 3; i++){tempradius[i] = Convert.ToDouble(temp[i]);}pcrdGeo.dLongitude = tempradius[0] / 180.0 pai + tempradius[1] / 180.0 / 60.0 pai + tempradius[2] / 180 / 60.0 / 60 pai;int deglon = Convert.ToInt32(pcrdGeo.dLongitude180 / pai);//求中央经度int num; //带号midlong = 0; //默认值,需要制定分带try{}catch{MessageBox.Show("Choose 3/6 error");return;}if (ttpareTo("3度带") == 0){num = Convert.ToInt32(deglon / 6 + 1);midlong = (6 num - 3) / 180.0 pai;}if (ttpareTo("6度带") == 0){num = Convert.ToInt32((deglon + 1.5) / 3);midlong = num 3 pai / 180;}double lp=pcrdGeo.dLongitude - midlong;double N = c / Math.Sqrt(1 + epp eppMath.Cos(pcrdGeo.dLatitude) Math.Cos(pcrdGeo.dLatitude));double M = c / Math.Pow(1 + epp eppMath.Cos(pcrdGeo.dLatitude) Math.Cos(pcrdGeo.dLatitude), 1.5);double ita = epp Math.Cos(pcrdGeo.dLatitude);double t = Math.Tan(pcrdGeo.dLatitude);double Nscnb = N Math.Sin(pcrdGeo.dLatitude) Math.Cos(pcrdGeo.dLatitude);double Ncosb = N Math.Cos(pcrdGeo.dLatitude);double cosb = Math.Cos(pcrdGeo.dLatitude);double X;double m0, m2, m4, m6, m8;double a0, a2, a4, a6, a8;m0 = a (1 - ee ee);m2 = 3.0 / 2.0 m0 ee ee;m4 = 5.0 / 4.0 ee ee m2;m6 = 7.0 / 6.0 ee ee m4;m8 = 9.0 / 8.0 ee ee m6;a0 = m0 + m2 / 2.0 + 3.0 / 8.0 m4 + 5.0 / 16.0 m6 + 35.0 / 128.0 m8;a2 = m2 / 2 + m4 / 2 + 15.0 / 32.0 m6 + 7.0 / 16.0 m8;a4 = m4 / 8.0 + 3.0 / 16.0 m6 + 7.0 / 32.0 m8;a6 = m6 / 32.0 + m8 / 16.0;a8 = m8 / 128.0;double B = pcrdGeo.dLatitude;double sb = Math.Sin(B);double cb = Math.Cos(B);double s2b = sb cb 2;double s4b = s2b (1 - 2 sb sb) 2;double s6b = s2b Math.Sqrt(1 - s4b s4b) + s4b Math.Sqrt(1 - s2b s2b);X = a0 B - a2 / 2.0 s2b + a4 s4b / 4.0 - a6 / 6.0 s6b;pcrdCar.x = Nscnb lp lp / 2.0 + Nscnb cosb cosb Math.Pow(lp, 4) (5 - t t + 9 ita ita + 4 Math.Pow(ita, 4)) / 24.0 + Nscnb Math.Pow(cosb, 4) Math.Pow(lp, 6) (61 - 58 t t + Math.Pow(t, 4)) / 720.0 + X;pcrdCar.y = Ncosb Math.Pow(lp, 1) + Ncosb cosb cosb (1 - t t + ita ita) / 6.0 Math.Pow(lp, 3) + NcosbMath.Pow(lp, 5) Math.Pow(cosb, 4) (5 - 18 t t+ Math.Pow(t, 4) + 14 ita ita - 58 ita ita t t) / 120.0 ;if (pcrdCar.y < 0)pcrdCar.y += 500000;richTextBox1.Text = "Results:\nX:\t" +Convert.ToString(pcrdCar.x) +"\nY:\t"+ Convert.ToString(pcrdCar.y);}private void button2_Click(object sender, EventArgs e) {double ee = 0;double a = 0;string tt;int num; //带号string ytext; //利用y值求带号和中央经线try{}catch{MessageBox.Show("Gauss Inverse: Choose datumerror");return;}if (ttpareTo("克氏椭球") == 0){a = 6378245.00;}if (ttpareTo("WGS-84") == 0){a = 6378.00;}if (ttpareTo("1975国际椭球") == 0){a = 6378140.00;}if (ttpareTo("2000国家大地坐标系") == 0){a = 6378.0;}double b = Math.Sqrt(a a (1 - ee ee));double c = a a / b;double epp = Math.Sqrt((a a - b b) / b / b); CRDGEODETIC pcrdGeo;CRDCARTESIAN pcrdCar;double midlong = 0;//求X,Y和带号pcrdCar.x = Convert.ToDouble(textBox4.Text); ytext = textBox5.Text;string temp = ytext.Substring(0, 2);num = Convert.ToInt32(temp);ytext = ytext.Remove(0, 2);pcrdCar.y = Convert.ToDouble(ytext) - 500000; try{}catch{MessageBox.Show("Choose 3/6 error");return;}if (ttpareTo("3度带") == 0){midlong = num 3 pai / 180;}if (ttpareTo("6度带") == 0){midlong = (6 num - 3) pai / 180;}b = Math.Sqrt(a a (1 - ee ee));c = a a / b;epp = Math.Sqrt(a a - b b) / b;double m0, m2, m4, m6, m8;double a0, a2, a4, a6, a8;m0 = a (1 - ee ee);m2 = 3.0 / 2.0 m0 ee ee;m4 = 5.0 / 4.0 ee ee m2;m6 = 7.0 / 6.0 ee ee m4;m8 = 9.0 / 8.0 ee ee m6;a0 = m0 + m2 / 2.0 + 3.0 / 8.0 m4 + 5.0 / 16.0 m6 + 35.0 / 128.0 m8;a2 = m2 / 2 + m4 / 2 + 15.0 / 32.0 m6 + 7.0 / 16.0 m8;a4 = m4 / 8.0 + 3.0 / 16.0 m6 + 7.0 / 32.0 m8;a6 = m6 / 32.0 + m8 / 16.0;a8 = m8 / 128.0;double Bf, B;Bf = pcrdCar.x / a0;B = 0.0;while (Math.Abs(Bf - B) > 1E-10){B = Bf;double sb = Math.Sin(B);double cb = Math.Cos(B);double s2b = sb cb 2;double s4b = s2b (1 - 2 sb sb) 2;double s6b = s2b Math.Sqrt(1 - s4b s4b) + s4b Math.Sqrt(1 - s2b s2b);Bf = (pcrdCar.x - (-a2 / 2.0 s2b + a4 / 4.0 s4b - a6 / 6.0 s6b)) / a0;}double itaf, tf, Vf, Nf;itaf = epp Math.Cos(Bf);tf = Math.Tan(Bf);Vf = Math.Sqrt(1 + epp epp Math.Cos(Bf)Math.Cos(Bf));Nf = c / Vf;double ynf = pcrdCar.y / Nf;pcrdGeo.dLatitude = Bf - 1.0 / 2.0 Vf Vf tf (ynf ynf - 1.0 / 12.0 Math.Pow(ynf, 4) (5 + 3 tf tf + itaf itaf - 9 Math.Pow(itaf tf, 2)) +1.0 / 360.0 (61 + 90 tf tf + 45 Math.Pow(tf, 4)) Math.Pow(ynf, 6));pcrdGeo.dLongitude = (ynf / Math.Cos(Bf) - (1 + 2 tf tf + itaf itaf) Math.Pow(ynf, 3) / 6.0 / Math.Cos(Bf) +(5 + 28 tf tf + 24 Math.Pow(tf, 4) + 6 itaf itaf + 8 Math.Pow(itaf tf, 2)) Math.Pow(ynf, 5) / 120.0 /Math.Cos(Bf));pcrdGeo.dLongitude = pcrdGeo.dLongitude + midlong;//pcrdGeo.dLatitude = pcrdGeo.dLatitude;richTextBox2.Text = "Results:\nLatitude: " + Convert.ToString(pcrdGeo.dLatitude) + "\nLongtitude: " +Convert.ToString(pcrdGeo.dLongitude);}private void label13_Click(object sender, EventArgs e) {}}}三.程序运行结果分析通过选取书上(de)具体实例进行测试,本程序(de)精度大体满足要求,一般正算(de)精度在0.01米和0.001米之间,反算(de)精度在0.0001秒左右.以下是程序运行(de)截图.。
高斯投影正算与反算的理论方法与实

高斯投影正算与反算的理论方法与实现代码高斯投影是正形投影的一种,同一坐标系中的高斯投影换带计算公式是根据正形投影原理推导出的两个高斯坐标系间的显函数式。
在同一大地坐标系中(例如1954北京坐标系或1980西安坐标系),如果两个高斯坐标系只是主子午线的经度不同,那么显函数式前的系数可以根据坐标系使用的椭球元素和主子午线经度唯一确定。
但如果两个高斯坐标系除了主子午线的经度不同以外,还存在其他线性系,则将线性变换公式代入换带计算的显函数式中,仍然可以得到严密的坐标变换公式。
此时显函数式前的系数等价于使用两个坐标系主子午线的经度和线性变换参数联合求解得到的,可以唯一确定。
//6度带宽54北京坐标系//高斯投影由大地坐标(Unit:Metres)反算经纬度(Unit:DD)void GaussProjInvCal(double X, double Y, double *longitude, double *latitude){int ProjNo; int ZoneWide; ////带宽double longitude1,latitude1, longitude0,latitude0, X0,Y0, xval,yval;double e1,e2,f,a, ee, NN, T,C, M, D,R,u,fai, iPI;iPI = 0.0174532925199433; ////3.1415926535898/180.0;a = 6378245.0; f = 1.0/298.3; //54年北京坐标系参数////a=6378140.0; f=1/298.257; //80年西安坐标系参数ZoneWide = 6; ////6度带宽ProjNo = (int)(X/1000000L) ; //查找带号longitude0 = (ProjNo-1) * ZoneWide + ZoneWide / 2;longitude0 = longitude0 * iPI ; //中央经线X0 = ProjNo*1000000L+500000L;Y0 = 0;xval = X-X0; yval = Y-Y0; //带内大地坐标e2 = 2*f-f*f;e1 = (1.0-sqrt(1-e2))/(1.0+sqrt(1-e2));ee = e2/(1-e2);M = yval;u = M/(a*(1-e2/4-3*e2*e2/64-5*e2*e2*e2/256));fai = u+(3*e1/2-27*e1*e1*e1/32)*sin(2*u)+(21*e1*e1/16-55*e1*e1*e1*e1/32)*sin( 4*u)+(151*e1*e1*e1/96)*sin(6*u)+(1097*e1*e1*e1*e1/512)*sin(8*u);C = ee*cos(fai)*cos(fai);T = tan(fai)*tan(fai);NN = a/sqrt(1.0-e2*sin(fai)*sin(fai));R = a*(1-e2)/sqrt((1-e2*sin(fai)*sin(fai))*(1-e2*sin(fai)*sin(fai))*(1-e2*sin(fai)*sin(fai)));D = xval/NN;//计算经度(Longitude) 纬度(Latitude)longitude1 =longitude0+(D-(1+2*T+C)*D*D*D/6+(5-2*C+28*T-3*C*C+8*ee+24*T*T)*D*D*D*D*D/120)/cos(fai);latitude1 = fai -(NN*tan(fai)/R)*(D*D/2-(5+3*T+10*C-4*C*C-9*ee)*D*D*D*D/24+(61+90*T+298*C+45*T*T-256*ee-3*C*C)*D*D*D*D*D*D/720);//转换为度DD*longitude = longitude1 / iPI;*latitude = latitude1 / iPI;}//高斯投影由经纬度(Unit:DD)正算平面坐标(含带号,Unit:Metres)void GaussProjCal(double longitude, double latitude, double *X, double *Y) {int ProjNo=0; int ZoneWide; ////带宽double longitude1,latitude1, longitude0,latitude0, X0,Y0, xval,yval; double a,f, e2,ee, NN, T,C,A, M, iPI;iPI = 0.0174532925199433; ////3.1415926535898/180.0;ZoneWide = 6; ////6度带宽a=6378245.0; f=1.0/298.3; //54年北京坐标系参数////a=6378140.0; f=1/298.257; //80年西安坐标系参数ProjNo = (int)(longitude / ZoneWide) ;longitude0 = ProjNo * ZoneWide + ZoneWide / 2;longitude0 = longitude0 * iPI ;latitude0=0;longitude1 = longitude * iPI ; //经度转换为弧度latitude1 = latitude * iPI ; //纬度转换为弧度e2=2*f-f*f;ee=e2*(1.0-e2);NN=a/sqrt(1.0-e2*sin(latitude1)*sin(latitude1));T=tan(latitude1)*tan(latitude1);C=ee*cos(latitude1)*cos(latitude1);A=(longitude1-longitude0)*cos(latitude1);M=a*((1-e2/4-3*e2*e2/64-5*e2*e2*e2/256)*latitude1-(3*e2/8+3*e2*e2/32+45*e2 *e2*e2/1024)*sin(2*latitude1)+(15*e2*e2/256+45*e2*e2*e2/1024)*sin(4*latitude1)-(3 5*e2*e2*e2/3072)*sin(6*latitude1));xval = NN*(A+(1-T+C)*A*A*A/6+(5-18*T+T*T+72*C-58*ee)*A*A*A*A*A/120); yval = M+NN*tan(latitude1)*(A*A/2+(5-T+9*C+4*C*C)*A*A*A*A/24+(61-58*T+T*T+600*C-330*ee)*A*A*A*A*A*A/720);X0 = 1000000L*(ProjNo+1)+500000L;Y0 = 0;xval = xval+X0; yval = yval+Y0;*X = xval;*Y = yval;}NN卯酉圈曲率半径,测量学里面用N表示M为子午线弧长,测量学里用大X表示fai为底点纬度,由子午弧长反算公式得到,测量学里用Bf表示R为底点所对的曲率半径,测量学里用Nf表示。
高斯投影正反算原理

高斯投影正反算原理高斯投影是一种常用于地图制图的投影方式,也被广泛应用于其他领域的空间数据处理。
高斯投影正反算是对于已知的地球坐标系上的位置(经纬度),通过计算得到该点的平面坐标(东、北坐标),或者对于已知的平面坐标(东、北坐标),通过计算得到该点的地球坐标系上的位置(经纬度)的过程。
本文将详细介绍高斯投影正反算的原理。
一、高斯投影简介高斯投影是一种圆锥投影,其投影面在地球表面的某个经线上,也就是说,投影面是以该经线为轴的圆锥面。
经过对圆锥体的调整后,使其切于地球椭球面,在该经线上进行投影,同时保持沿该经线方向的比例尺一致,从而达到地图上各点在包括该经线的垂直面上映射的目的。
这种投影方式在某一特定区域内得到高精度的结果,因此广泛应用于地图制图。
二、高斯投影数学模型对于高斯投影正反算,需要先建立高斯投影坐标系与地球坐标系的转换模型。
1.高斯投影坐标系的建立高斯投影坐标系的建立需要确定圆锥面的基本参数,首先需要确定其所处的中央子午线,再确定该子午线上的经度为零点,并利用该经线上某一点的经度和该点的高度来确定该点所在的圆锥体。
圆锥体的底面包括所有与地球椭球面相切的圆面,通过对这些圆面进行调整,使得圆锥体转动后能够在中央子午线上进行投影。
在此基础上,可建立高斯投影坐标系,其中投影面为圆锥面,且中央子午线与投影面的交点称为该投影坐标系的中心,投影面的上端点和下端点分别对应正北方向和正南方向。
2.地球坐标系的建立地球坐标系是以地球椭球体为基础建立的,其坐标系原点确定为地球椭球体上的一个特定点。
在已知该点经纬度和高度的前提下,可确定以该点为中心的地球椭球体,并可根据它与地球坐标系之间的转换关系得到平面坐标系。
3.高斯投影坐标系与地球坐标系之间的转换关系由于高斯投影坐标系与地球坐标系存在不同的坐标体系和基准面,因此需要通过数学关系式来建立它们之间的转换关系。
(1)高斯投影坐标系转地球坐标系:已知高斯投影坐标系中任意一点的东北坐标(N,E),以及所属的中央子午线经度λ0、椭球参数a和e,则可通过以下公式求出该点的地球坐标系经纬度(φ,λ)和高度H:A0为以地球椭球体中心为原点,高斯投影坐标系中心投影坐标为(0,0)的点到椭球面的距离。
高斯投影正反算公式_新

高斯投影坐大地坐标系由大地基准面和地图投影确定,由地图投影到特定椭圆柱面后在南北两极剪开展开而成,是对地球表面的逼近,各国或地区有各自的大地基准面,我国目前主要采用的基准面为:1.WGS84基准面,为GPS基准面,17届国际大地测量协会上推荐,椭圆柱长半轴a=6378137m,短半轴b=6356752.3142451m;2.西安80坐标系,1975年国际大地测量协会上推荐,椭圆柱长半轴a=6378140m,短半轴b=6356755.2881575m;3.北京54坐标系,参照前苏联克拉索夫斯基椭球体建立,椭圆柱长半轴6度0-6度,角。
值为正,Y增加500公里;反算则是由高斯平面坐标(X,Y)求解大地坐标(L,B)。
二、计算模/***************************************本文直接依据空间立体三角函数关系得出结果。
*****/(一)正算由图表1,由方程式(1),令,可得在图表2中,,则由椭圆方程,令可知:(三、程序代doubleL=(m_L-6.0*L0//换算成弧度doublexita=atan(b*b*tan(B)/a/a/cos(L));doubledxita=0.000001;doublexi=dxita;x=0.0;doublec=a*a/b/b;while(xi<xita){x+=dxita/sqrt(c*sin(xi)*sin(xi)+cos(xi)*cos(xi));xi+=dxita;坐标*&B,do ubledoubledxi=0.000001;doublexi=dxi;doubleX=0.0;doublec=a*a/b/b;while(X<x/a){X+=dxi/sqrt(c*sin(xi)*sin(xi)+cos(xi)*cos(xi));xi+=dxi;}doubler=a/sqrt(c*sin(xi)*sin(xi)+cos(xi)*cos(xi));。
高斯正反算及空间直角坐标与大地地理坐标转换

高斯正反算及空间直角坐标与大地地理坐标转换一、实验目的与要求1.对以上理论内容的验证与应用。
2.通过学习掌握测绘软件开发过程与方法,初步具备测绘软件开发基本技能。
3.熟练掌握Visual C++编程环境的使用,了解其特点与程序开发过程,掌软件调试、测试的技术方法。
4.分析测绘程序设计技术课程中相关软件的结构和模块功能,掌握结构化程序设计方法和技术,掌握测绘数据处理问题的基本特点。
5.开发相关程序功能模块,独立完成相关问题概念结构分析、程序结构设计、模块设计、代码编写、调试、测试等工作。
二、实验安排1.实验时数12学时。
2.每实验小组可以由3~4人组成,或独立完成。
若由几个人完成程序设计,应进行合理的分工。
三、实验步骤和要点1.熟悉程序设计任务书的基本内容,调查了解软件需求状况,进行需求分析;2.进行总体设计。
根据所调查收集的资料和任务书的要求,对系统的硬件资源进行初步设计,提出硬件配置计划;进行软件总体设计,设计出软件程序功能的模块;3.根据总体设计的结果,进行详细设计,进行数据存储格式设计、算法等,写出逻辑代码;4.编写程序代码,调试运行;5.程序试运行。
最后同学们可根据自己的选题,写出软件开发设计书一份,打印程序代码和运行结果。
四实验原理高斯正反算:高斯正反算包括两部分内容:高斯正算和高斯反算。
简单的说就是大地地理坐标系坐标(B,L)与其对应的高斯平面直角坐标系坐标(x,y)之间的转换。
若已知大地地理坐标系坐标(B,L)解求对应的高斯平面直角坐标系坐标(x,y)称为高斯正算;反之,则为高斯反算。
空间直角坐标与大地地理坐标转换:地球表面可用一个椭球面表示。
设空间直角坐标系为OXYZ,当椭球的中心与空间直角坐标系原点重合,空间坐标系Z 轴与地球旋转重合(北极方向为正),X 轴正向经度为零时,就可以确定空间直角坐标系与大地地理坐标系的数学关系。
⎪⎩⎪⎨⎧+-=+=+=B H e N Z LB H N Y L B H N X sin ])1([sin cos )(cos cos )(2 式中 N 为卯酉圈曲率半径,B e a N 22sin 1-=; e 为椭球偏心率,222a b a e -=(a ,b 为椭球长半轴和短半轴)。
高斯正反算程序

高斯正反算程序
高斯正反算程序是一种用于计算地理坐标和直角坐标之间转换的程序。
它基于高斯投影原理,将地球表面上的点投影到平面上,实现地理坐标和直角坐标之间的转换。
高斯正算程序是将地理坐标(经度和纬度)转换为直角坐标(X 和Y)的程序。
其计算步骤包括:
1. 将经度和纬度转换为弧度;
2. 计算经度余弦、纬度正弦和纬度余弦;
3. 计算X和Y坐标。
高斯反算程序是将直角坐标(X和Y)转换为地理坐标(经度和纬度)的程序。
其计算步骤包括:
1. 计算X和Y坐标的平方和;
2. 计算平方差;
3. 计算经度和纬度的弧度值;
4. 将弧度值转换为度数。
以上是高斯正反算程序的基本步骤,具体的计算公式和细节可能因不同的投影方法和地区而有所不同。
在实际应用中,需要根据具体情况选择合适的投影方法和参数进行计算。
高斯投影正反算编程

高斯投影正反算编程一.高斯投影正反算基本公式(1)高斯正算基本公式(2)高斯反算基本公式以上主要通过大地测量学基础课程得到.这不进行详细的推导.只是列出基本公式指导编程的进行。
二.编程的基本方法和流程图(1)编程的基本方法高斯投影正反算基本上运用了所有的编程基本语句.本文中是利用C++语言进行基本的设计。
高斯正算中对椭球参数和带宽的选择主要运用了选择语句。
而高斯反算中除了选择语句的应用.在利用迭代算法求底点纬度还应用了循环语句。
编程中还应特别注意相关的度分秒和弧度之间的相互转换.这是极其重要的。
(2)相关流程图1)正算2)反算三.编程的相关代码(1)正算# include "stdio.h"# include "stdlib.h"# include "math.h"# include "assert.h"#define pi (4*atan(1.0))int i;struct jin{double B;double L;double L0;};struct jin g[100];main(int argc, double *argv[]) {FILE *r=fopen("a.txt","r"); assert(r!=NULL);FILE *w=fopen("b.txt","w"); assert(r!=NULL);int i=0;while(fscanf(r,"%lf %lf %lf",&g[i].B,&g[i].L,&g[i].L 0)!=EOF){double a,b;int zuobiao;printf("\n请输入坐标系:北京54=1.西安80=2.WGS84=3:");scanf("%d",&zuobiao);getchar();if(zuobiao==1){a=6378245;b=6356863.0187730473;}if(zuobiao==2){a=6378140;b=6356755.2881575287;}if(zuobiao==3){a=6378137;b=6356752.3142;} //选择坐标系//double f=(a-b)/a;double e,e2;e=sqrt(2*f-f*f);e2=sqrt((a/b)*(a/b)-1);//求椭球的第一.第二曲率//double m0,m2,m4,m6,m8;double a0,a2,a4,a6,a8;m0=a*(1-e*e);m2=3*e*e*m0/2;m4=5*e*e*m2/4;m6=7*e*e*m4/6;m8=9*e*e*m6/8;a0=m0+m2/2+3*m4/8+5*m6/16+35*m8/128;a2=m2/2+m4/2+15*m6/32+7*m8/16;a4=m4/8+3*m6/16+7*m8/32;a6=m6/32+m8/16;a8=m8/128;double Bmiao,Lmiao, L0miao;Bmiao=(int)(g[i].B)*3600.0+(int)((g[i].B-(int)(g[i].B) )*100.0)*60.0+(g[i].B*100-(int)(g[i].B*100))*100.0;Lmiao=(int)(g[i].L)*3600.0+(int)((g[i].L-(int)(g[i].L) )*100.0)*60.0+(g[i].L*100-(int)(g[i].L*100))*100.0;L0miao=(int)(g[i].L0)*3600.0+(int)((g[i].L0-(int)(g[i] .L0))*100.0)*60.0+(g[i].L0*100-(int)(g[i].L0*100))*100 .0;double db;db=pi/180.0/3600.0;double B1,L1,l;B1=Bmiao*db;L1= Lmiao*db;l=L1-L0miao*db;//角度转化为弧度//double T=tan(B1)*tan(B1);double n=e2*e2*cos(B1)*cos(B1);double A=l*cos(B1);double X,x,y;X=a0*(B1)-a2*sin(2*B1)/2+a4*sin(4*B1)/4-a6*sin(6*B1)/6 +a8*sin(8*B1)/8;//求弧长//double N=a/sqrt(1-e*e*sin(B1)*sin(B1));int Zonewide;int Zonenumber;printf("\n请输入带宽:3度带或6度带Zonewide=");scanf("%d",&Zonewide);getchar();if(Zonewide==3){Zonenumber=(int)((g[i].L-Zonewide/2)/Zonewide+1);}else if(Zonewide==6){Zonenumber=(int)g[i].L/Zonewide+1;}else{printf("错误");exit(0);}//选择带宽//doubleFE=Zonenumber*1000000+500000;//改写为国家通用坐标//y=FE+N*A+A*A*A*N*(1-T*T+n*n)/6+A*A*A*A*A*N*(5-18*T*T+T *T*T*T+14*n*n-58*n*n*T*T)/120;x=X+tan(B1)*N*A*A/2+tan(B1)*N*A*A*A*A*(5-T*T+9*n*n+4*n *n*n*n)/24+tan(B1)*N*A*A*A*A*A*A*(61-58*T*T+T*T*T*T)/7 20;printf("\n所选坐标系的转换结果:x=%lf y=%lf\n",x,y);fprintf(w,"%lf %lf\n",x,y);//输出结果到文本文件//}fclose(r);fclose(w);system("pause");return 0;}(2)反算# include "stdio.h"# include "stdlib.h"# include "math.h"# include "assert.h"#define pi (4*atan(1.0))double X,Y,B1,B2,B3,F,t;double m0,m2,m4,m6,m8;double a0,a2,a4,a6,a8,a1,b1;double BB,LL,Bf;double e,e1;int d,m,s,i,zuobiao;double sort(double,double);struct jin{double x;double y;double L0;};struct jin g[100];//x,y,L0为输入量:x,y坐标和中央子午线经度//main(int argc, double *argv[]){FILE *r=fopen("c.txt","r");assert(r!=NULL);FILE *w=fopen("d.txt","w");assert(r!=NULL);int i=0;while(fscanf(r,"%lf %lf %lf",&g[i].x,&g[i].y,&g[i].L 0)!=EOF)//文件为空.无法打开//{double a1=6378245.0000000000;//克拉索夫斯基椭球参数//double b1=6356863.0187730473;double a75=6378140.0000000000;//1975国际椭球参数//double b75=6356755.2881575287;double a84=6378137.0000000000;//WGS-84系椭球参数//double b84=6356752.3142000000;double M,N;//mouyou圈曲率半径.子午圈曲率半径// double t,n;double A,B,C;double BB,LL,Bf,LL0,BB0;double a,b;printf("\n选择参考椭球:1=克拉索夫斯基椭球.2=1975国际椭球.3=WGS-84系椭球:");scanf("%d",&zuobiao);getchar();if(zuobiao==1){a=a1;b=b1;}if(zuobiao==2){a=a75;b=b75;}if(zuobiao==3){a=a84;b=b84;}//选择参考椭球.求解第一偏心率e,第二偏心率e1// Bf=sort(a,b);//调用求解底点纬度的函数//double q=sqrt(1-e*e*sin(Bf)*sin(Bf));double G=cos(Bf);M=a*(1-e*e)/(q*q*q);N=a/q;double H,I;A=g[i].y/N;H=A*A*A;I=A*A*A*A*A;t=tan(Bf);n=e1*cos(Bf);B=t*t;C=n*n;BB0=Bf-g[i].y*t*A/(2*M)+g[i].y*t*H/(24*M)*(5+3*B+C-9*B *C)-g[i].y*t*I/(720*M)*(61+90*B+45*B*B);LL0=g[i].L0*pi/180.0+A/G-H/(6*G)*(1.0+2*B+C)+I/(120*G)*(5.0+28*B+24*B*B+6*C+8*B*C);//利用公式求解经纬度// int Bdu,Bfen,Ldu,Lfen;double Bmiao,Lmiao;Ldu=int(LL0/pi*180);Lfen=int((LL0/pi*180)*60-Ldu*60);Lmiao=LL0/pi*180*3600-Ldu*3600-Lfen*60;Bdu=int(BB0/pi*180);Bfen=int((BB0/pi*180)*60-Bdu*60);Bmiao=BB0/pi*180*3600-Bdu*3600-Bfen*60;//将弧度转化为角度//printf("\n所选坐标系的转换结果:%d度%d分%lf秒 %d 度%d分%lf秒 \n",Bdu,Bfen,Bmiao,Ldu,Lfen,Lmiao);fprintf(w,"%d°%d’%lf”%d°%d’%lf”\n",Bdu,Bfen,Bmiao,Ldu,Lfen,Lmiao);//将结果输出到文本文件//}fclose(r);fclose(w);system("pause");return 0;}double sort(double a,double b){double e,e1;e=sqrt(1-(b/a)*(b/a));e1=sqrt((a/b)*(a/b)-1);double m0,m2,m4,m6,m8;double a0,a2,a4,a6,a8;m0=a*(1-e*e);m2=3*e*e*m0/2;m4=5*e*e*m2/4;m6=7*e*e*m4/6;m8=9*e*e*m6/8;a0=m0+m2/2+3*m4/8+5*m6/16+35*m8/128;a2=m2/2+m4/2+15*m6/32+7*m8/16;a4=m4/8+3*m6/16+7*m8/32;a6=m6/32+m8/16;a8=m8/128;B1=g[i].x/a0;do{F=-a2*sin(2*B1)/2+a4*sin(4*B1)/4-a6*sin(6*B1)/6+a8*s in(8*B1)/8;B2=(g[i].x-F)/a0;B3=B1;B1=B2;} while(fabs(B3-B2)>10e-10);//利用迭代算法求解底点纬度//return B2;}。
高斯投影坐标反算公式

d B 0
x dx dl l y dy dl l
x dx l tan dy y l
x N sin B cos 3 N N sin B cos Bl (5 t 2 9 2 4 4 )l 3 l 6 N sin B cos 5 B (61 58t 2 t 4 )l 5 120 y cos 2 B cos 4 B 2 2 2 N cos B{1 (1 t )l (5 18t 2 t 4 )l 4 } l 2 24
121212角子午线收敛角方大地方位角坐标方位对于克拉索夫斯基椭球1根据高斯坐标确定带号计算中央子午线经度六度带24114923计算中央子午线经度2迭代法求取大地纬度迭代开始时设以后每次迭代按下式计算为止3计算反算公式中的各符号的值2428cos120掌握子午线收敛角的定义及作用
上节回顾
高斯投影坐标正算要求: 1、掌握公式中各符号的意义,公式不要求记住; 2、每个同学必须会计算,一般计算两遍。
1 sin B cos 4 Bl 5 (2 t 2 ) 15
1 3
说明:⑴ 为 l 的奇函数,而且 l 越大, 越大; ⑵ 有正负,当描写点在中央子午线以东时, 为正, 以西时, 为负; ⑶ 当 l 不变时,则 随纬度增加而增大。 (2)由平面坐标f B Bf y2 2M f N f cos B f 1 1 2 2 l y (1 2t f n f ) y 3 3 N f cos B f 6 N f cos B f
算例:带第38带 3 x 5586514 .369 m y 4374 .724 m X 5586514 .369 m Y 38504374 .724 m
- 1、下载文档前请自行甄别文档内容的完整性,平台不提供额外的编辑、内容补充、找答案等附加服务。
- 2、"仅部分预览"的文档,不可在线预览部分如存在完整性等问题,可反馈申请退款(可完整预览的文档不适用该条件!)。
- 3、如文档侵犯您的权益,请联系客服反馈,我们会尽快为您处理(人工客服工作时间:9:00-18:30)。
地球科学与环境工程学院实验报告书课程名:学号:姓名:指导老师:日期:目录一、目的与要求 (1)二、实验容 (1)三、计算公式整理 (1)四、程序代码 (4)五、计算结果 (15)六、实验体会 (16)一、目的与要求参考椭球面是测量计算的基准面。
坐标是椭球面上的基本坐标系,根据测量的观测成果(如距离与方向),从原点出发,逐点计算在椭球面上的坐标;或根据两点的坐标,计算它们之间的线长度和方位角,这类计算称为问题解算(或称为主题解算)。
问题解算的用途是多方面的,随着现代空间技术和航空航天、航海等领域的发展,问题解算(尤其是反算)有着更为重要的作用,因此需要熟练掌握其计算。
二、实验容在《测量学基础》教材中,介绍了高斯平均引数法与白塞尔方法的计算过程、步骤。
鉴于此,需要熟练掌握高斯平均引数法与白塞尔方法解主题问题的基本方法与原理。
采用所熟悉的计算机语言编程计算。
计算时采用克拉索夫椭球参数,至少完成其中一种方反算,按照数据序号选取不同的已知数据,在计算结果中注明所选取的数据序号,选取其它数据作为无效数据处理。
三、计算公式整理3.1、高斯平均引数正算计算公式(S< 200 km)3.2、高斯平均引数正算计算公式(S< 200 km)四、程序代码4.1、角度转换类的头文件:#pragma onceconst double Pi=3.141592653589793;class AngleTrans{public:AngleTrans(void);~AngleTrans(void);double D,F,M,DFM,Rad,Ten;double trans1(double DFM), //度分秒形式的角度转换为弧度形式trans2(double Rad), //弧度形式的角度转换为度分秒形式trans3(double D); //十进制度转化为弧度};4.2、角度转换类的源文件:#include"StdAfx.h"#include"AngleTrans.h"#include<cmath>AngleTrans::AngleTrans(void){}AngleTrans::~AngleTrans(void){}//度分秒转换为弧度double AngleTrans::trans1(double DFM){ D=floor(DFM);F=floor((DFM-D)*100);M=(DFM-D-F/100)*10000;Ten=D+F/60+M/3600;Rad=Ten/180*Pi;return Rad;}//弧度转换为度分秒double AngleTrans::trans2(double Rad){ Ten=Rad/Pi*180;D=floor(Ten);F=(Ten-D)*60;M=(F-floor(F))*60;F=floor(F);DFM=D+F/100+M/10000;return DFM;}//十进制度转化为弧度double AngleTrans::trans3(double D){ Rad=D/180*Pi;return Rad;}4.3、正反算类的头文件:#pragma onceclass ZhengFanSuan{public:ZhengFanSuan(void);~ZhengFanSuan(void);double zB1,zL1,zA12,zS,fB1,fL1,fB2,fL2;double ZhengSuanB(double zB1,double zL1,double zA12,double zS),ZhengSuanL(double zB1,double zL1,double zA12,double zS),ZhengSuanA(double zB1,double zL1,double zA12,double zS);double FanSuanA12(double fB1,double fL1,double fB2,double fL2),FanSuanS(double fB1,double fL1,double fB2,double fL2),FanSuanA21(double fB1,double fL1,double fB2,double fL2); };4.3、正反算类的源文件:#include"StdAfx.h"#include"ZhengFanSuan.h"#include"AngleTrans.h"#include<cmath>ZhengFanSuan::ZhengFanSuan(void)}ZhengFanSuan::~ZhengFanSuan(void){}AngleTrans _AngleTrans;const double e1= 0.0066934216622966,e2=0.006738525414683,a=6378245.0000,b= 6356863.01877,temp=pow(10.0, -10);//精度要求double Calc_M(double z) //计算Mm{ double x=a*(1-pow(e1,2))/sqrt(pow(1-pow(e1,2)*pow(sin(z),2),3));return x;}double Calc_N(double z) //计算Nm{ double x=a/sqrt(1-pow(e1,2)*pow(sin(z),2));return x;}double Calc_t(double z) //计算tm{ double x=tan(z);return x;}double Calc_yita(double z) //计算yitam{ double x=pow(e2,2)*pow(cos(z),2);return x;}//正算纬度double ZhengFanSuan::ZhengSuanB(double zB1,double zL1,double zA12,double zS) { double M[10000],N[10000],t[10000],B[10000],Bm[10000],L[10000],Lm[10000],A[10000],Am[10000],yita[10000];double_zB1=_AngleTrans.trans1(zB1),_zL1=_AngleTrans.trans1(zL1),_zA12=_AngleTrans.trans1(zA12),final;t[1]=tan(_zB1);yita[1]=pow(e2,2)*pow(cos(_zB1),2);N[1]=a/sqrt(1-pow(e1,2)*pow(sin(_zB1),2));M[1]=a*(1-pow(e1,2))/sqrt(pow(1-pow(e1,2)*pow(sin(_zB1),2),3));B[0]=_AngleTrans.trans3(206265/M[1]*zS*cos(_zA12)/3600),L[0]=_AngleTrans.trans3(206265*zS*sin(_zA12)/(N[1]*cos(_zB1))/3600),A[0]=L[0]*sin(_zB1);Bm[1]=_zB1+1/2*B[0];Lm[1]=_zL1+1/2*L[0];Am[1]=_zA12+1/2*A[0];int i=1;B[1]=_AngleTrans.trans3((206265/M[1]*zS*cos(Am[1])*(1+pow(zS,2)/(24*pow(N[1],2))*(pow(sin(Am[1]),2))*(2+3*pow(t[1],2)+2*pow(yita[1],2))+3*pow(cos(Am[1]),2)*pow(yita[1],2)*(pow(t[1],2)-1-pow(yita[1],2)-4*pow(yita[1],2)*pow(t[1],2))))/3600);L[i]=_AngleTrans.trans3((206265/(N[1]*cos(Bm[1]))*zS*sin(Am[i])*(1+pow(zS,2)/(24*pow(N[i],2))*(pow(t[i] ,2)*pow(sin(Am[i]),2)-pow(cos(Am[i]),2)*(1+pow(yita[i],2)-9*pow(yita[i],2)*pow(t[i],2)))))/3600);A[i]=_AngleTrans.trans3((206265/N[i]*zS*sin(Am[i])*t[i]*(1+pow(zS,2)/(24*pow(N[i],2))*(pow(cos(Am[i]),2 )*(2+7*pow(yita[i],2)+9*pow(yita[i],2)*pow(t[i],2)+5*pow(yita[i],4))+pow(sin(Am[i]),2)*(2+pow(t[i],2)+2*pow(yita[i],2)))))/3600);while(B[i]-B[i-1]>=temp){B[i+1]=_AngleTrans.trans3((206265/M[i]*zS*cos(Am[i])*(1+pow(zS,2)/(24*pow(N[i],2))*(pow(sin(Am[i]),2) )*(2+3*pow(t[i],2)+2*pow(yita[i],2))+3*pow(cos(Am[i]),2)*pow(yita[i],2)*(pow(t[i],2)-1-pow(yita[i],2)-4*pow(yita[i],2)*pow(t[i],2))))/3600);L[i+1]=_AngleTrans.trans3((206265/(N[i]*cos(Bm[i]))*zS*sin(Am[i])*(1+pow(zS,2)/(24*pow(N[i],2))*(pow(t[ i],2)*pow(sin(Am[i]),2)-pow(cos(Am[i]),2)*(1+pow(yita[i],2)-9*pow(yita[i],2)*pow(t[i],2)))))/3600);A[i+1]=_AngleTrans.trans3((206265/N[i]*zS*sin(Am[i])*t[i]*(1+pow(zS,2)/(24*pow(N[i],2))*(pow(cos(Am[i] ),2)*(2+7*pow(yita[i],2)+9*pow(yita[i],2)*pow(t[i],2)+5*pow(yita[i],4))+pow(sin(Am[i]),2)*(2+pow(t[i],2)+2*pow(yita[i],2)))))/3600);Bm[i+1]=_zB1+1/2*B[i];Lm[i+1]=_zL1+1/2*L[i];Am[i+1]=_zA12+1/2*A[i];M[i+1]=Calc_M(B[i]);N[i+1]=Calc_N(N[i]);yita[i+1]=Calc_yita(yita[i]);t[i+1]=Calc_t(t[i]);final=B[i];i++;}double Final=final+_zB1;return _AngleTrans.trans2(Final);}//正算经度double ZhengFanSuan::ZhengSuanL(double zB1,double zL1,double zA12,double zS){ double M[10000],N[10000],t[10000],B[10000],Bm[10000],L[10000],Lm[10000],A[10000],Am[10000],yita[10000];double_zB1=_AngleTrans.trans1(zB1),_zL1=_AngleTrans.trans1(zL1),_zA12=_AngleTrans.trans1(zA12),final;t[1]=tan(_zB1);yita[1]=pow(e2,2)*pow(cos(_zB1),2);N[1]=a/sqrt(1-pow(e1,2)*pow(sin(_zB1),2));M[1]=a*(1-pow(e1,2))/sqrt(pow(1-pow(e1,2)*pow(sin(_zB1),2),3));B[0]=_AngleTrans.trans3(206265/M[1]*zS*cos(_zA12)/3600),L[0]=_AngleTrans.trans3(206265*zS*sin(_zA12)/(N[1]*cos(_zB1))/3600),A[0]=L[0]*sin(_zB1);Bm[1]=_zB1+1/2*B[0];Lm[1]=_zL1+1/2*L[0];Am[1]=_zA12+1/2*A[0];int i=1;B[1]=_AngleTrans.trans3((206265/M[1]*zS*cos(Am[1])*(1+pow(zS,2)/(24*pow(N[1],2))*(pow(sin(Am[1]),2))*(2+3*pow(t[1],2)+2*pow(yita[1],2))+3*pow(cos(Am[1]),2)*pow(yita[1],2)*(pow(t[1],2)-1-pow(yita[1],2)-4*pow(yita[1],2)*pow(t[1],2))))/3600);L[i]=_AngleTrans.trans3((206265/(N[1]*cos(Bm[1]))*zS*sin(Am[i])*(1+pow(zS,2)/(24*pow(N[i],2))*(pow(t[i],2)*pow(sin(Am[i]),2)-pow(cos(Am[i]),2)*(1+pow(yita[i],2)-9*pow(yita[i],2)*pow(t[i],2)))))/3600);A[i]=_AngleTrans.trans3((206265/N[i]*zS*sin(Am[i])*t[i]*(1+pow(zS,2)/(24*pow(N[i],2))*(pow(cos(Am[i]),2 )*(2+7*pow(yita[i],2)+9*pow(yita[i],2)*pow(t[i],2)+5*pow(yita[i],4))+pow(sin(Am[i]),2)*(2+pow(t[i],2)+2*pow(yita[i],2)))))/3600);if(L[i]-L[i-1]>=temp){ while(L[i]-L[i-1]>=temp){B[i+1]=_AngleTrans.trans3((206265/M[i]*zS*cos(Am[i])*(1+pow(zS,2)/(24*pow(N[i],2))*(pow(sin(Am[i]),2) )*(2+3*pow(t[i],2)+2*pow(yita[i],2))+3*pow(cos(Am[i]),2)*pow(yita[i],2)*(pow(t[i],2)-1-pow(yita[i],2)-4*pow(yita[i],2)*pow(t[i],2))))/3600);L[i+1]=_AngleTrans.trans3((206265/(N[i]*cos(Bm[i]))*zS*sin(Am[i])*(1+pow(zS,2)/(24*pow(N[i],2))*(pow(t[ i],2)*pow(sin(Am[i]),2)-pow(cos(Am[i]),2)*(1+pow(yita[i],2)-9*pow(yita[i],2)*pow(t[i],2)))))/3600);A[i+1]=_AngleTrans.trans3((206265/N[i]*zS*sin(Am[i])*t[i]*(1+pow(zS,2)/(24*pow(N[i],2))*(pow(cos(Am[i] ),2)*(2+7*pow(yita[i],2)+9*pow(yita[i],2)*pow(t[i],2)+5*pow(yita[i],4))+pow(sin(Am[i]),2)*(2+pow(t[i],2)+2*pow(yita[i],2)))))/3600);Bm[i+1]=_zB1+1/2*B[i];Lm[i+1]=_zL1+1/2*L[i];Am[i+1]=_zA12+1/2*A[i];M[i+1]=Calc_M(B[i]);N[i+1]=Calc_N(N[i]);yita[i+1]=Calc_yita(yita[i]);t[i+1]=Calc_t(t[i]);final=L[i];i++;}}elsefinal=L[i];double Final=final+_zL1;return _AngleTrans.trans2(Final);}//正算方位角double ZhengFanSuan::ZhengSuanA(double zB1,double zL1,double zA12,double zS){ double M[100],N[100],t[100],B[100],Bm[100],L[100],Lm[100],A[100],Am[100],yita[10000];double_zB1=_AngleTrans.trans1(zB1),_zL1=_AngleTrans.trans1(zL1),_zA12=_AngleTrans.trans1(zA12),final;t[1]=tan(_zB1);yita[1]=pow(e2,2)*pow(cos(_zB1),2);N[1]=a/sqrt(1-pow(e1,2)*pow(sin(_zB1),2));M[1]=a*(1-pow(e1,2))/sqrt(pow(1-pow(e1,2)*pow(sin(_zB1),2),3));B[0]=_AngleTrans.trans3(206265/M[1]*zS*cos(_zA12)/3600),L[0]=_AngleTrans.trans3(206265*zS*sin(_zA12)/(N[1]*cos(_zB1))/3600),A[0]=L[0]*sin(_zB1);Bm[1]=_zB1+1/2*B[0];Lm[1]=_zL1+1/2*L[0];Am[1]=_zA12+1/2*A[0];int i=1;B[1]=_AngleTrans.trans3((206265/M[1]*zS*cos(Am[1])*(1+pow(zS,2)/(24*pow(N[1],2))*(pow(sin(Am[1]),2))*(2+3*pow(t[1],2)+2*pow(yita[1],2))+3*pow(cos(Am[1]),2)*pow(yita[1],2)*(pow(t[1],2)-1-pow(yita[1],2)-4*pow(yita[1],2)*pow(t[1],2))))/3600);L[i]=_AngleTrans.trans3((206265/(N[1]*cos(Bm[1]))*zS*sin(Am[i])*(1+pow(zS,2)/(24*pow(N[i],2))*(pow(t[i] ,2)*pow(sin(Am[i]),2)-pow(cos(Am[i]),2)*(1+pow(yita[i],2)-9*pow(yita[i],2)*pow(t[i],2)))))/3600);A[i]=_AngleTrans.trans3((206265/N[i]*zS*sin(Am[i])*t[i]*(1+pow(zS,2)/(24*pow(N[i],2))*(pow(cos(Am[i]),2 )*(2+7*pow(yita[i],2)+9*pow(yita[i],2)*pow(t[i],2)+5*pow(yita[i],4))+pow(sin(Am[i]),2)*(2+pow(t[i],2)+2*pow(yita[i],2)))))/3600);while(A[i]-A[i-1]>=temp){B[i+1]=_AngleTrans.trans3((206265/M[i]*zS*cos(Am[i])*(1+pow(zS,2)/(24*pow(N[i],2))*(pow(sin(Am[i]),2) )*(2+3*pow(t[i],2)+2*pow(yita[i],2))+3*pow(cos(Am[i]),2)*pow(yita[i],2)*(pow(t[i],2)-1-pow(yita[i],2)-4*pow(yita[i],2)*pow(t[i],2))))/3600);L[i+1]=_AngleTrans.trans3((206265/(N[i]*cos(Bm[i]))*zS*sin(Am[i])*(1+pow(zS,2)/(24*pow(N[i],2))*(pow(t[ i],2)*pow(sin(Am[i]),2)-pow(cos(Am[i]),2)*(1+pow(yita[i],2)-9*pow(yita[i],2)*pow(t[i],2)))))/3600);A[i+1]=_AngleTrans.trans3((206265/N[i]*zS*sin(Am[i])*t[i]*(1+pow(zS,2)/(24*pow(N[i],2))*(pow(cos(Am[i] ),2)*(2+7*pow(yita[i],2)+9*pow(yita[i],2)*pow(t[i],2)+5*pow(yita[i],4))+pow(sin(Am[i]),2)*(2+pow(t[i],2)+2*pow(yita[i],2)))))/3600);Bm[i+1]=_zB1+1/2*B[i];Lm[i+1]=_zL1+1/2*L[i];Am[i+1]=_zA12+1/2*A[i];M[i+1]=Calc_M(B[i]);N[i+1]=Calc_N(N[i]);yita[i+1]=Calc_yita(yita[i]);t[i+1]=Calc_t(t[i]);final=A[i];i++;}double Final;if(_zA12>Pi){ Final=final+_zA12-Pi;}else{ Final=final+_zA12+Pi;}return _AngleTrans.trans2(Final);}//反算Sdouble ZhengFanSuan::FanSuanS(double fB1,double fL1,double fB2,double fL2){ double A12,A21,S;double _fB1=_AngleTrans.trans1(fB1),_fL1=_AngleTrans.trans1(fL1),_fB2=_AngleTrans.trans1(fB2),_fL2=_AngleTrans.trans1(fL2);double Bm=(_fB1+_fB2)/2,_B=_fB2-_fB1,_L=_fL2-_fL1;double Nm=Calc_N(Bm),tm=Calc_t(Bm),yitam=Calc_yita(Bm);double Vm=sqrt(1+pow(e2,2)*pow(cos(Bm),2)),r01=Nm*cos(Bm),r03=(-Nm*pow(cos(Bm),3)*pow(tm,2))/24,r21=(Nm*cos(Bm)*(1+pow(yitam,2)-9*pow(yitam,2)*pow(tm,2)*pow(yitam,4))/pow(Vm,4))/24,s10=Nm/pow(Vm,2),s12=(Nm/pow(Vm,2)*pow(cos(Bm),2)*(2+3*pow(tm,2)+2*pow(yitam,2)))/24,s30=(Nm/pow(Vm,6)*(pow(yitam,2)-pow(yitam,2)*pow(tm,2)))/8,t01=tm*cos(Bm),t21=cos(Bm)*tm/(24*pow(Vm,4))*(2+7*pow(yitam,2)+9*pow(yitam,2)*pow(tm,2)+5*pow(yitam,4)),t03=(tm*pow(cos(Bm),3)*(2+pow(tm,2)+2*pow(yitam,2)))/24;double U=r01*_L+r21*pow(_B,2)*_L+r03*pow(_L,3),V=s10*_B+s12*_B*pow(_L,2)+s30*pow(_B,3),_A=t01*_L+t21*pow(_B,2)*_L+t03*pow(_L,3);double Am=atan(U/V);S=V/cos(Am);return S;}//反算A12double ZhengFanSuan::FanSuanA12(double fB1,double fL1,double fB2,double fL2){ double T,A12,A21,S;double _fB1=_AngleTrans.trans1(fB1),_fL1=_AngleTrans.trans1(fL1),_fB2=_AngleTrans.trans1(fB2),_fL2=_AngleTrans.trans1(fL2);double Bm=(_fB1+_fB2)/2,_B=_fB2-_fB1,_L=_fL2-_fL1;double Nm=Calc_N(Bm),tm=Calc_t(Bm),yitam=Calc_yita(Bm);double Vm=sqrt(1+pow(e2,2)*pow(cos(Bm),2)),r01=Nm*cos(Bm),r03=(-Nm*pow(cos(Bm),3)*pow(tm,2))/24,r21=(Nm*cos(Bm)*(1+pow(yitam,2)-9*pow(yitam,2)*pow(tm,2)*pow(yitam,4))/pow(Vm,4))/24,s10=Nm/pow(Vm,2),s12=(Nm/pow(Vm,2)*pow(cos(Bm),2)*(2+3*pow(tm,2)+2*pow(yitam,2)))/24,s30=(Nm/pow(Vm,6)*(pow(yitam,2)-pow(yitam,2)*pow(tm,2)))/8,t01=tm*cos(Bm),t21=cos(Bm)*tm/(24*pow(Vm,4))*(2+7*pow(yitam,2)+9*pow(yitam,2)*pow(tm,2)+5*pow(yitam,4)),t03=(tm*pow(cos(Bm),3)*(2+pow(tm,2)+2*pow(yitam,2)))/24;double U=r01*_L+r21*pow(_B,2)*_L+r03*pow(_L,3),V=s10*_B+s12*_B*pow(_L,2)+s30*pow(_B,3),_A=t01*_L+t21*pow(_B,2)*_L+t03*pow(_L,3);double Am=atan(U/V);A12=Am-0.5*_A;return _AngleTrans.trans2(A12);}//反算A21double ZhengFanSuan::FanSuanA21(double fB1,double fL1,double fB2,double fL2){ double T,A12,A21,S;double _fB1=_AngleTrans.trans1(fB1),_fL1=_AngleTrans.trans1(fL1),_fB2=_AngleTrans.trans1(fB2),_fL2=_AngleTrans.trans1(fL2);double Bm=(_fB1+_fB2)/2,_B=_fB2-_fB1,_L=_fL2-_fL1;double Nm=Calc_N(Bm),tm=Calc_t(Bm),yitam=Calc_yita(Bm);double Vm=sqrt(1+pow(e2,2)*pow(cos(Bm),2)),r01=Nm*cos(Bm),r03=(-Nm*pow(cos(Bm),3)*pow(tm,2))/24,r21=(Nm*cos(Bm)*(1+pow(yitam,2)-9*pow(yitam,2)*pow(tm,2)*pow(yitam,4))/pow(Vm,4))/24,s10=Nm/pow(Vm,2),s12=(Nm/pow(Vm,2)*pow(cos(Bm),2)*(2+3*pow(tm,2)+2*pow(yitam,2)))/24,s30=(Nm/pow(Vm,6)*(pow(yitam,2)-pow(yitam,2)*pow(tm,2)))/8,t01=tm*cos(Bm),t21=cos(Bm)*tm/(24*pow(Vm,4))*(2+7*pow(yitam,2)+9*pow(yitam,2)*pow(tm,2)+5*pow(yitam,4)),t03=(tm*pow(cos(Bm),3)*(2+pow(tm,2)+2*pow(yitam,2)))/24;double U=r01*_L+r21*pow(_B,2)*_L+r03*pow(_L,3),V=s10*_B+s12*_B*pow(_L,2)+s30*pow(_B,3),_A=t01*_L+t21*pow(_B,2)*_L+t03*pow(_L,3);double Am=atan(U/V);A12=Am-0.5*_A;if(A12>Pi){ A21=Am+0.5*_A-Pi;}else{ A21=Am+0.5*_A+Pi;}return _AngleTrans.trans2(A21);}4.4、正算的计算按钮代码:void C主题高斯引数正反算Dlg::OnBnClickedButton2(){ UpdateData(true);ZhengFanSuan _ZhengFanSuan;double B1=_wtof(zB1),L1=_wtof(zL1),A12=_wtof(zA12),S=_wtof(zS),B2=_ZhengFanSuan.ZhengSuanB(B1,L1,A12,S),L2=_ZhengFanSuan.ZhengSuanL(B1,L1,A12,S),A21=_ZhengFanSuan.ZhengSuanA(B1,L1,A12,S);zB2.Format(_T("%.7f"),B2);zL2.Format(_T("%.7f"),L2),zA21.Format(_T("%.7f"),A21);UpdateData(false);}4.5、反算的计算按钮代码:void C主题高斯引数正反算Dlg::OnBnClickedButton1(){ UpdateData(true);ZhengFanSuan _ZhengFanSuan;double B1=_wtof(fB1),L1=_wtof(fL1),B2=_wtof(fB2),L2=_wtof(fL2),S=_ZhengFanSuan.FanSuanS(B1,L1,B2,L2),A12=_ZhengFanSuan.FanSuanA12(B1,L1,B2,L2),A21=_ZhengFanSuan.FanSuanA21(B1,L1,B2,L2);fS.Format(_T("%.7f"),S);fA12.Format(_T("%.7f"),A12);fA21.Format(_T("%.7f"),A21);UpdateData(false);}4.6、清零按钮代码:void C主题高斯引数正反算Dlg::OnBnClickedButton3(){ UpdateData(true);fB1.Format(_T("%.0f"),0);fL1.Format(_T("%.0f"),0);fB2.Format(_T("%.0f"),0);fL2.Format(_T("%.0f"),0);fA12.Format(_T("%.0f"),0);fA21.Format(_T("%.0f"),0);fS.Format(_T("%.0f"),0);UpdateData(false);}void C主题高斯引数正反算Dlg::OnBnClickedButton4(){ UpdateData(true);zB1.Format(_T("%.0f"),0);zL1.Format(_T("%.0f"),0);zB2.Format(_T("%.0f"),0);zL2.Format(_T("%.0f"),0);zA12.Format(_T("%.0f"),0);zA21.Format(_T("%.0f"),0);zS.Format(_T("%.0f"),0);UpdateData(false);}五、计算结果数据组号:18六、实习体会实验一开始本来考虑用控制台做,后来发现用控制台做出来可视化效果不好,而且程序中要用很多cin和cout,写代码十分麻烦,所以最后选择了用MFC做。