百度墨卡托平面坐标转经纬度坐标--查表法 唯一有效！！！

转自http://www.site-digger.com/html/articles/20150831/98.html。

百度墨卡托平面坐标转经纬度坐标

发布时间：2015-08-31

---

以”钟楼”为例，在百度地图中搜索（结果如上图所示），返回的是12128773.43,4040249.00这种形式的坐标。这是什么坐标呢？这是墨卡托平面坐标（http://developer.baidu.com/map/jsdevelop-6.htm）。

1. 如何测试这种坐标（墨卡托平面坐标）？

作者发现百度的坐标拾取系统（http://api.map.baidu.com/lbsapi/getpoint/index.html）是兼容这种坐标的。简单地说，就是可以在百度的坐标拾取系统里反查这种坐标。如下图所示。

2. 如何将百度地图的墨卡托坐标转为经纬度坐标呢？

百度地图API提供了这个功能，示例页面：

http://www.site-digger.com/tools/mct2latlng.html 转换算法是JS实现的（查看页面源码可见），效果如下图所示。

其值和在百度坐标拾取系统里看到的值（“当前坐标点如下”栏内）一致。

另外，这里有个Python版本的实现(https://gist.github.com/JobsDong/8313339)，亲测准确。

转百度坐标BD09的python代码如下，本人亲测有效：

#!/usr/bin/python
#-*- coding=utf-8 -*-

"""用于和地理信息有关的数据处理
"""

__author__ = ['"wuyadong" <wuyadong@tigerknows.com>']

import math

MCBAND = (12890594.86, 8362377.87, 5591021, 3481989.83, 1678043.12, 0)
MC2LL = ([1.410526172116255e-8, 0.00000898305509648872, -1.9939833816331,
          200.9824383106796, -187.2403703815547, 91.6087516669843, - 23.38765649603339,
          2.57121317296198, -0.03801003308653, 17337981.2],
         [-7.435856389565537e-9, 0.000008983055097726239, -0.78625201886289,
          96.32687599759846, -1.85204757529826, -59.36935905485877, 47.40033549296737,
          -16.50741931063887, 2.28786674699375, 10260144.86],
         [-3.030883460898826e-8, 0.00000898305509983578, 0.30071316287616,
          59.74293618442277, 7.357984074871, -25.38371002664745, 13.45380521110908,
          -3.29883767235584, 0.32710905363475, 6856817.37],
         [-1.981981304930552e-8, 0.000008983055099779535, 0.03278182852591, 40.31678527705744,
          0.65659298677277, -4.44255534477492, 0.85341911805263, 0.12923347998204,
          -0.04625736007561, 4482777.06],
         [3.09191371068437e-9, 0.000008983055096812155, 0.00006995724062, 23.10934304144901,
          -0.00023663490511, -0.6321817810242, -0.00663494467273, 0.03430082397953,
          -0.00466043876332, 2555164.4],
         [2.890871144776878e-9, 0.000008983055095805407, -3.068298e-8, 7.47137025468032,
          -0.00000353937994, -0.02145144861037, -0.00001234426596, 0.00010322952773,
          -0.00000323890364, 826088.5])

X_PI = 3.14159265358979324 * 3000.0 / 180.0


class GISError(Exception):
    """GIS Exception
    """


def convert_MCT_2_BD09(lon, lat):
    """将墨卡托坐标转换成BD09
        Args:
            lon: float, 经度
            lat: float, 维度
        Returns:
            (x, y): tuple, 经过转换的x, y
    """
    ax = None

    # 获取常量ax
    for j in range(len(MCBAND)):
        if lat >= MCBAND[j]:
            ax = MC2LL[j]
            break

    if ax is None:
        raise GISError("error lat:%s" % lat)

    e = ax[0] + ax[1] * abs(lon)
    i = abs(lat) / ax[9]
    aw = ax[2] + ax[3] * i + ax[4] * i * i + ax[5] * i * i * i +\
         ax[6] * i * i * i * i + ax[7] * i * i * i * i * i + ax[8] * i * i * i * i * i * i
    if lon < 0:
        e *= -1
    if lat < 0:
        aw *= -1
    return e, aw


def convert_BD09_2_GCJ03(lon, lat):
    """坐标转换，将BD09转成GCJ03
        Args:
            lon: float, 经度
            lat: float, 维度
        Returns:
            (x, y): tuple, 转换后的结果
    """
    x = lon - 0.0065
    y = lat - 0.006
    z = math.sqrt(x * x + y * y) - 0.00002 * math.sin(y * X_PI)
    theta = math.atan2(y, x) - 0.000003 * math.cos(x * X_PI)
    gg_lon = z * math.cos(theta)
    gg_lat = z * math.sin(theta)
    return gg_lon, gg_lat

if __name__ == '__main__':
    mocator_lng_lat = [13271401.86,2985177.68,13271454.86,2985165.68,13271642.86,2985215.68,13271655.86,2985236.68,13271576.86,2985392.68,13271555.86,2985417.68,13271459.86,2985656.68,13271412.86,2985822.67,13271388.86,2985840.67,13271112.86,2985770.68,13271102.86,2985724.68,13271401.86,2985177.68]
    baidu_lng_lat_lst = []
    for i in range(int(len(mocator_lng_lat)/2)):
        baidu_lng_lat = convert_MCT_2_BD09(mocator_lng_lat[i*2], mocator_lng_lat[i*2+1])
        baidu_lng_lat_lst.append(baidu_lng_lat)
    print(baidu_lng_lat_lst)