using System;
using System.Collections.Generic;
using System.Linq;
using System.Text;
using System.Threading.Tasks;
namespace MRP
{
public class ClassMct
{
static public int __IterativeTimes = ; //反向转换程序中的迭代次数
static public double __IterativeValue = ; //反向转换程序中的迭代初始值
static public double __A = 6378.137; //椭球体长轴,千米
static public double __B = 6356.752314; //椭球体短轴,千米
static public double __B0 = ; //标准纬度,弧度
static public double __L0 = ; //原点经度,弧度
//角度到弧度的转换
static public double DegreeToRad(double degree)
{
return Math.PI * degree / ;
}
//弧度到角度的转换
static public double RadToDegree(double rad)
{
return ( * rad) / Math.PI;
}
//设定__A与__B
static public void SetAB(double a, double b)
{
if (a <= || b <= )
{
return;
}
__A = a;
__B = b;
}
//设定__B0
static public void SetLB0(double pmtL0, double pmtB0)
{
double l0 = DegreeToRad(pmtL0);
if (l0 < -Math.PI || l0 > Math.PI)
{
return;
}
__L0 = l0;
double b0 = DegreeToRad(pmtB0);
if (b0 < -Math.PI / || b0 > Math.PI / )
{
return;
}
__B0 = b0;
}
/*******************************************
经纬度转XY坐标
pmtLB0: 参考点经纬度
pmtLB1: 要转换的经纬度
返回值: 直角坐标,单位:公里
*******************************************/
static public PointXY LBToXY(PointLB pmtLB0, PointLB pmtLB1)
{
SetLB0(pmtLB0.lon, pmtLB0.lat);
double B = DegreeToRad(pmtLB1.lat);
double L = DegreeToRad(pmtLB1.lon);
PointXY xy = new PointXY();
xy.x = ; xy.y = ;
double f/*扁率*/, e/*第一偏心率*/, e_/*第二偏心率*/, NB0/*卯酉圈曲率半径*/, K, dtemp;
double E = Math.Exp();
if (L < -Math.PI || L > Math.PI || B < -Math.PI / || B > Math.PI / )
{
return xy;
}
if (__A <= || __B <= )
{
return xy;
}
f = (__A - __B) / __A;
dtemp = - (__B / __A) * (__B / __A);
if (dtemp < )
{
return xy;
}
e = Math.Sqrt(dtemp);
dtemp = (__A / __B) * (__A / __B) - ;
if (dtemp < )
{
return xy;
}
e_ = Math.Sqrt(dtemp);
NB0 = ((__A * __A) / __B) / Math.Sqrt( + e_ * e_ * Math.Cos(__B0) * Math.Cos(__B0));
K = NB0 * Math.Cos(__B0);
xy.x = K * (L - __L0);
xy.y = K * Math.Log(Math.Tan(Math.PI / + (B) / ) * Math.Pow(( - e * Math.Sin(B)) / ( + e * Math.Sin(B)), e / ));
double y0 = K * Math.Log(Math.Tan(Math.PI / + (__B0) / ) * Math.Pow(( - e * Math.Sin(__B0)) / ( + e * Math.Sin(__B0)), e / ));
xy.y = xy.y - y0;
xy.y = -xy.y;//正常的Y坐标系(向上)转程序的Y坐标系(向下)
return xy;
}
/*******************************************
XY坐标转经纬度
pmtLB0: 参考点经纬度
pmtXY: 要转换的XY坐标,单位:公里
返回值: 经纬度
*******************************************/
static public PointLB XYtoLB(PointLB pmtLB0, PointXY pmtXY)
{
SetLB0(pmtLB0.lon, pmtLB0.lat);
double X = pmtXY.x;
double Y = -pmtXY.y;//程序的Y坐标系(向下)转正常的Y坐标系(向上)
double B = , L = ;
PointLB lb = new PointLB();
lb.lat = ; lb.lon = ;
double f/*扁率*/, e/*第一偏心率*/, e_/*第二偏心率*/, NB0/*卯酉圈曲率半径*/, K, dtemp;
double E = Math.Exp();
if (__A <= || __B <= )
{
return lb;
}
f = (__A - __B) / __A;
dtemp = - (__B / __A) * (__B / __A);
if (dtemp < )
{
return lb;
}
e = Math.Sqrt(dtemp);
dtemp = (__A / __B) * (__A / __B) - ;
if (dtemp < )
{
return lb;
}
e_ = Math.Sqrt(dtemp);
NB0 = ((__A * __A) / __B) / Math.Sqrt( + e_ * e_ * Math.Cos(__B0) * Math.Cos(__B0));
K = NB0 * Math.Cos(__B0);
double y0 = K * Math.Log(Math.Tan(Math.PI / + (__B0) / ) * Math.Pow(( - e * Math.Sin(__B0)) / ( + e * Math.Sin(__B0)), e / ));
Y = Y + y0;
L = X / K + __L0;
B = __IterativeValue;
for (int i = ; i < __IterativeTimes; i++)
{
B = Math.PI / - * Math.Atan(Math.Pow(E, (-Y / K)) * Math.Pow(E, (e / ) * Math.Log(( - e * Math.Sin(B)) / ( + e * Math.Sin(B)))));
}
lb.lon = RadToDegree(L);
lb.lat = RadToDegree(B);
return lb;
}
}
}