开发电子海图肯定要涉及到墨卡托计算,分享墨卡托转平面坐标的算法
struct GeoPoint
{
double lon;
double lat;
};
double m_R;
//计算海图的基准维度,一般取30°
void getR(double lat)
{
lat = lat*PI / 180;
double N = 6378137 / sqrt(1 - 0.08189*0.08189*sin(lat)*sin(lat));
m_R = N*cos(lat);
}
//将墨卡托海图经纬度坐标转换为平面坐标,单位 :米
void getMokato(GeoPoint *point)
{
point->lon = point->lon / m_icoord*PI / 180;
point->lat = point->lat / m_icoord*PI / 180;
double Q = log(tan(0.785398163375 + point->lat / 2)) - 0.040945*log((1 + 0.08189*sin(point->lat)) / (1 - 0.08189*sin(point->lat)));
point->lon = point->lon*m_R ;
point->lat = m_R*Q ;
}
double getLat(double pointx, double x)
{
double d_e = 0.08189;
double G = d_e / 2 * log((1 + d_e *sin(x)) / (1 - d_e * sin(x))) + pointx / m_R;
return 2 * atan(exp(G)) - PI / 2;
}
//将屏幕坐标转换为墨卡托海图坐标,进行多次迭代计算,
void getLatlon(GeoPoint *point)
{
point->lat = point->lat;
double d_x = 0.0;
d_x = getLat(point->lat, d_x);
d_x = getLat(point->lat, d_x);
d_x = getLat(point->lat, d_x);
d_x = getLat(point->lat, d_x);
point->lat= d_x * 180 / PI;
point->lon = point->lon / m_R * 180 / PI;
}
联系作者及参考文章:http://www.sailxy.com