@[TOC](测绘程序设计基础 实验3(2) CSU)
实验3(2) 函数
（高斯投影正反算程序设计）
(工具：VS2010)
一、 实验目的
• 掌握函数的定义、引用及应用方法
二、实验内容与要求
• 编写高斯投影正算和反算的两个函数，并设计简单界面对函数计算的正确性进行测试。
三、设计与实现：
3.1 设计思路：
3.2 界面设计：
3.3主要代码：
3.3.1文件:<Csupport.cpp>
/***************************************************************************
* 文件名:<Csupport.cpp> *
* *
* 描述:WG84坐标系高斯坐标系 *
* *
* 历史：**日期** **理由** **签名** *
* 2019年3月27日 代码格式修改 *** *
* *
* 外部过程： *
* *
/**************************************************************************/
/***************************************************************************
* 名字：Process_oppDegree(double dDms) *
* *
* 描述：度分秒形式--》十进制的角度值 *
* *
* 历史：**日期** **理由** **签名** *
* 2019年3月8日 创建该函数 *** *
* *
* 参数： 1.double dDms *
* *
* *
* 返回值：double类型数据 dDeg 返回十进制的角度值 *
* *
* 注： *
/**************************************************************************/
double Csupport::Process_oppDegree(double dDms)
{
int iDegree,iMin;
double dSec;
double dDeg;
iDegree=int(dDms);
iMin=int((dDms-iDegree)*100);
dSec=((dDms-iDegree)*100-iMin)*100;
dDeg=iDegree+double(iMin)/60+dSec/3600;
return dDeg;
}
/***************************************************************************
* 名字：double Process_oppD(double dDeg) *
* *
* 描述：十进制的角度值--》弧度式方位角 *
* *
* 历史：**日期** **理由** **签名** *
* 2019年3月27日 创建该函数 *** *
* 参数： 1.double dDeg*PI/180 *
* *
* *
* 返回值：double类型数据 返回弧度值 *
* *
* 注： *
/**************************************************************************/
double Csupport::Process_oppD(double dDeg)
{
return dDeg*Pi/180;
}
/***************************************************************************
* 名字：double Process_Second(double dD) *
* *
* 描述：十进制的角度值--》 秒形式 *
* *
* 历史：**日期** **理由** **签名** *
* 2019年3月8日 创建该函数 *** *
* 2019年3月27日 修改 *** *
* *
* 参数： 1.double dD *
* *
* *
* 返回值：double类型数据 dDegree 返回秒形式 *
* *
* 注： *
/**************************************************************************/
double Csupport::Process_Second(double dD)
{
double dDegree;
int iDegree;
int iMin;
double dSec;
double dTmp;
iDegree=int(dD);
dTmp=(dD-iDegree)*100;
iMin=int(dTmp);
dSec=(dTmp-iMin)*100;
dDegree=iDegree*3600+double(iMin)*60+dSec;
return dDegree/206265;
}
/***************************************************************************
* 名字：GetX(double dB,double e,double a) *
* *
* 描述：求X的初值 *
* *
* 历史：**日期** **理由** **签名** *
* 2019年3月27日 创建该函数 *** *
* 参数： 1.double dB *
* 2.double e *
* 3.double a *
* *
* 返回值：double类型数据 返回X的初值 *
* *
* 注： *
/**************************************************************************/
double Csupport::GetX(double dB,double e,double a)
{
const double A0=(1+e*e*3/4+45/64*pow(e,4)+350/512*pow(e,6)+11025/16384*pow(e,8));
const double A2=-(e*e*3/4+pow(e,4)*60/64+pow(e,6)*525/512+pow(e,8)*17640/16384)/2;
const double A4=(pow(e,4)*15/64+pow(e,6)*210/512+pow(e,8)*8820/16384)/4;
const double A6=-(pow(e,6)*35/512+pow(e,8)*2520/16384)/6;
const double A8=(pow(e,8)*315/16384)/8;
double T;
double T1=1-e*e;
double T2=A0*dB+A2*sin(2*dB)+A4*sin(4*dB)+A6*sin(6*dB)+A8*sin(8*dB);
T=a*T1*T2;
//T=a*(1-e*e)*(A0*dB+A2*sin(2*dB)+A4*sin(4*dB)+A6*sin(6*dB)+A8*sin(8*dB));
return T;
}
/***************************************************************************
* 名字：Zone(double dB,double dL,double t,int &ZoneNum) *
* *
* 描述：依据带号进行计算 *
* *
* 历史：**日期** **理由** **签名** *
* 2019年3月29日 创建该函数 *** *
* 参数： 1.double dB *
* 2.double dL *
* 3.double t *
* 4.int &ZoneNum *
* *
* 返回值：double类型数据 当前带中央经度 并将带号通过引用给与ZoneNum *
* *
* 注： *
/**************************************************************************/
int Csupport::Zone(double dB,double dL,double t,int &ZoneNum)
{
int fZone;
if(t==0)
{
fZone=3;
dL-=1.5;
}
else{
fZone=6;
}
ZoneNum=( int)dL/fZone+1;
if(ZoneNum<=0)
ZoneNum=360/fZone+ZoneNum;
int CML=(fZone==3?(ZoneNum-0.5)*fZone+1.5:CML=6*ZoneNum-3);
return CML;
}
/***************************************************************************
* 名字：void XY(double B,double L,double &dX,double &dY,int dS) *
* *
* 描述：正算主函数，调用了函数Process_oppD 、Process_Second 、GetX三个函数*
* *
* 历史：**日期** **理由** **签名** *
* 2019年3月27日 创建该函数 *** *
* 参数： 1.double B *
* 2.double L *
* 3.double &dX *
* 4.double &dY *
* 5.double dCML *
* *
* 返回值：无返回值，利用引用求解 *
* *
* 注： *
/**************************************************************************/
void Csupport::XY(double B,double L,double &dX,double &dY,int dS)
{
B=Process_oppDegree(B);//转为十进制
//L=Process_oppDegree(L);
int ZoneNum;
int dCML=Zone(B,L,dS,ZoneNum);
//B L 修改
double dB=Process_oppD(B);//转为弧度值
double dL=Process_Second(L-dCML);
//定义某些常量
const double a=6378137;
const double b=6356752.3142451795;
double e1=(sqrt(a*a-b*b)/a);
double e2=(sqrt(a*a-b*b)/b);
double Eta=e2*cos(dB);
double N=(a/sqrt(1-e1*e1*sin(dB)*sin(dB)));
double t=tan(dB);
//计算dX
double t1=N*sin(dB)*cos(dB)*dL*dL/2;
double t2=N*sin(dB)*pow(cos(dB),3)*(5-t*t+9*pow(Eta,2)+4*pow(Eta,4))*pow(dL,4)/24;
double t4=N*sin(dB)*pow(cos(dB),5)*(61-58*t*t+pow(t,4))*pow(dL,6)/720;
double tmp=GetX(dB,e1,a);
dX=tmp+t1+t2+t4;
//计算dY
dY=5e5+N*cos(dB)*dL+N*pow(cos(dB),3)*(1-t*t+Eta*Eta)*pow(dL,3)/6+N*pow(cos(dB),5)*(5-18*t*t+pow(t,4)+14*Eta*Eta-58*pow(Eta*t,2))*pow(dL,5)/120;
dY=(dY>=0?dY+ZoneNum*1e6:dY-ZoneNum*1e6);
}
/***************************************************************************
* 名字：double Process_Degree(double dD) *
* *
* 描述：十进制的角度值--》度分秒形式 *
* *
* 历史：**日期** **理由** **签名** *
* 2019年3月8日 创建该函数 *** *
* *
* 参数： 1.double dD *
* *
* *
* 返回值：double类型数据 dDegree 返回度分秒形式 *
* *
* 注： *
/**************************************************************************/
double Csupport::Process_Degree(double dD)
{
double dDegree;
int iDegree;
int iMin;
double dSec;
double dTmp;
iDegree=int(dD);
dTmp=(dD-iDegree)*60;
iMin=int(dTmp);
dSec=(dTmp-iMin)*60;
dDegree=iDegree+double(iMin)/100+dSec/10000;
return dDegree;
}
/***************************************************************************
* 名字：GetBf(double dX,double e,double a) *
* *
* 描述：求得Bf的值 *
* *
* 历史：**日期** **理由** **签名** *
* 2019年3月27日 创建该函数 *** *
* *
* 参数： 1.double dX *
* 2.double e *
* 3.double a *
* *
* 返回值：double类型数据 Bf 返回度分秒形式 *
* *
* 注： *
/**************************************************************************/
double Csupport::GetBf(double dX,double e,double a)
{
const double A0=(1+e*e*3/4+45/64*pow(e,4)+350/512*pow(e,6)+11025/16384*pow(e,8));
double B0=dX/(a*(1-e*e)*A0);
const double K0=(e*e*3/4+pow(e,4)*45/64+pow(e,6)*350/512+pow(e,8)*11025/16384)/2;
const double K2=-(pow(e,4)*63/64+pow(e,6)*1108/512+pow(e,8)*58239/16384)/3;
const double K4=(pow(e,6)*604/512+pow(e,8)*68484/16384)/3;
const double K6=-(pow(e,8)*26328/16384)/3;
return (B0+sin(2*B0)*(K0+pow(sin(B0),2)*(K2+pow(sin(B0),2)*(K4+K6*pow(sin(B0),2)))));
}
/***************************************************************************
* 名字：Zone2(double dY,double dS) *
* *
* 描述：求得Bf的值 *
* *
* 历史：**日期** **理由** **签名** *
* 2019年3月27日 创建该函数 *** *
* *
* 参数： 1.double dY *
* 2.double dS *
* *
* 返回值：double类型数据 Bf 返回度分秒形式 *
* *
* 注： *
/**************************************************************************/
int Csupport::Zone2(double dY,double dS)
{
int tmp=dY/1e6;
if (dS==0)
{
return (tmp*3);
}
else
return (6*tmp-3);
}
/***************************************************************************
* 名字：BL(double B,double L,double &dX,double &dY,double dCML) *
* *
* 描述：反算主函数，调用了函数GetBf 函数 *
* *
* 历史：**日期** **理由** **签名** *
* 2019年3月27日 创建该函数 *** *
* 参数： 1.double B *
* 2.double L *
* 3.double &dX *
* 4.double &dY *
* 5.double dCML *
* *
* 返回值：无返回值，利用引用求解 *
* *
* 注： *
/**************************************************************************/
void Csupport::BL(double &B,double &L,double dX,double dY,double dS)
{
int L0=Zone2(dY,dS);
int tmp=dY/1e5;
dY=dY-tmp*1e5;
const double a=6378137;
const double b=6356752.3142451795;
const double Alfa=1/298.257223563;
double e1=(sqrt(a*a-b*b)/a);
double e2=(sqrt(a*a-b*b)/b);
double Bf=GetBf(dX,e1,a);
double tf=tan(Bf);
double Etaf=e2*cos(Bf);
double Nf=a/sqrt(1-pow(e1*sin(Bf),2));
double Mf=Nf/(1+pow(e2*cos(Bf),2));
double dB=Bf-tf/(2*Mf*Nf)*dY*dY
+tf/(24*Mf*pow(Nf,3))*(5+3*tf*tf+Etaf*Etaf-9*pow(Etaf*tf,2))*pow(dY,4)
-tf/(720*Mf*pow(Nf,5))*(61+90*tf*tf+45*pow(tf,4))*pow(dY,6);
double l=dY/(Nf*cos(Bf))
-(1+2*tf*tf+Etaf*Etaf)*pow(dY,3)/(6*pow(Nf,3)*cos(Bf))
+(5+28*tf*tf+24*pow(tf,4)+6*Etaf*Etaf+8*pow(Etaf*tf,2))*pow(dY,5)
/(120*pow(Nf,5)*cos(Bf));
B=dB*180/Pi;
L=l*180/Pi+L0;
B=Process_Degree(B);
L=Process_Degree(L);
}
3.3.2文件:< RS_110_zhang_SY3_02Dlg.cpp >
/***************************************************************************
* 文件名:<RS_110_zhang_SY3_02Dlg.cpp> *
* *
* 描述:WG84坐标系高斯坐标系对话框文件 *
* *
* 历史：**日期** **理由** **签名** *
* 2019年3月27日 代码格式修改 *** *
* *
* 外部过程： *
* *
/**************************************************************************/
void CRS_110_zhang_SY3_02Dlg::OnBnClickedButton1Calculationpositive()
{
UpdateData(true);
Csupport k;
k.XY(B,L,X,Y,m_3);
UpdateData(false);
// TODO: 在此添加控件通知处理程序代码
}
void CRS_110_zhang_SY3_02Dlg::OnBnClickedButton2Calculationback()
{
// TODO: 在此添加控件通知处理程序代码
UpdateData(true);
Csupport k;
k.BL(B,L,X,Y,m_3);
UpdateData(false);
}
3.4运行结果
3.5设计技巧：
通过类调用，方便修改补充。直接编写函数，最后调用。
3305
被折叠的 条评论
为什么被折叠?
到【灌水乐园】发言
