本文主要讲述地理坐标系统的原理以及怎么利用Python进行地理坐标系统转换，主要内容包含以下几块：

1. 什么是地理坐标系统？
2. 常用的地理坐标系统有哪些？常用地图产品分别是什么地理坐标系？
3. 怎么样利用Python实现地理坐标系的转换？
4. 介绍一个能够不写代码就进行地理坐标系转换的软件。

1. 什么是地理坐标系统

用一张图来解释什么是地理坐标系统。

上图中红色实线圈为地球自然表面，由于地球自然表面坑坑洼洼、凹凸不平，所以在进行精准测量和定位的时候是不能直接用数据公式来拟合地球表面的。

这时，用一个可以近似表示地球表面的规则的椭圆（如上图中的蓝色虚线框所示）来进行地球表面的定位和测量，这个规则的三维球面就是地球椭球体。但是，地球表面高低起伏不平，造成地球椭球体与地球自然表面在有的地方能够十分贴合，在有的地方则不是很贴合，这会造成测量及定位地误差。

所以，需要在地球椭球体以及地球自然表面之间加一层中间层，使其能够尽量贴合真实的地球表面，那这中间层就是大地基准面（如图中绿色实线圈所示），其可以理解为在特定区域内与地球自然表面极为吻合的椭球体。地球上各个地方的地理位置及环境都是独一无二的，为了满足当地测量及定位的精度要求，全世界的地理学家们建立了无数个适合当地的大地基准面。

地理坐标系统就是由大地基准面衍生而来的，其是使用三维球面来定义地球表面位置，以实现通过经纬度对地球表面点位引用的坐标系。一个地理坐标系包括角度测量单位、本初子午线和基准面三部分。一个大地基准面可以对应多个地理坐标系统，而一个地理坐标系统只能对应一个大地基准面。（这边为了方便理解，没有加入大地水准面的概念）

同一个坐标点在不同地理坐标系的地图上，会落在不同的区域；同一个地点获取不同地理坐标系下的坐标数据，值不相同。

因此，在地图制图以及空间分析之前，先要了解坐标点以及地图的地理坐标系统。

2. 常用的地理坐标系统有哪些？常用地图产品分别是什么地理坐标系？

常用的地理坐标系统有：

1. 地球坐标系（WGS-84）

WGS-84是国际通用坐标系，也叫地球坐标系，GPS系统就是采用的WGS-84。WGS-84对于具体地方的位置描述可能不如当地坐标系来的准确，但是它对全球范围内的位置估计更准确。在中国范围内，一般不直接使用WGS-84，而是使用由国家测绘局在WGS-84基础上加密的火星坐标系（GCJ-02），或者使用企业在GCJ-02上二次加密的坐标系，例如百度坐标系（BD-09）、搜狗坐标系等。

2. 北京54坐标系（BJ-54）

BJ-54是建国初期提出的地理坐标系，因此在早期有比较广泛的运用，有一定比例的数据使用的是BJ-54。从现代的眼光看，它并不能十分准确地表达我国国境内的空间位置。

3. 西安80坐标系(XIAN-80)

XIAN-80由于后期意识到BJ-54的不足，我国1978年4月在西安召开全国天文大地网平差会议，确定重新定位，建立的我国新地理坐标系，它在中国经济建设、国防建设和科学研究中发挥了巨大作用。

4. 2000国家大地坐标系(CGCS-2000)

CGCS-2000我国当前最新的国家大地坐标系。2018年，我国国土资源系统全面采用CGCS-2000，并要求各类国土资源数据向CGCS-2000进行转换。

5. 地方独立坐标系

许多城市、矿区基于实用、方便与科学的目的，建立了地方坐标系。

国内地图产品的地理坐标系集中于CGCS-2000、GCJ-02以及一些对GCJ-02加密的商业坐标系，国外地图产品的地理坐标系则基本是WGS-84，下表罗列了常用地图产品的地理坐标系，供大家参考。

3. 利用Python实现地理坐标系统的转换

本节主要介绍怎么利用Python实现常用地理坐标系统之间的转换，并对转换结果进行精度检验。

话不多说，先上Python代码。

# -*- coding: utf-8 -*-
import math
# 设置常量
pi = 3.141592653589793234  # π
r_pi = pi * 3000.0 / 180.0
la = 6378245.0  # 长半轴
ob = 0.00669342162296594323  # 扁率

# gcj02 -> bd09
# lng为gcj02的经度
# lat为gcj02的纬度
# 返回值为转换后坐标的列表形式，[经度, 纬度]
def gcj02_bd09(lon_gcj02, lat_gcj02):
    b = math.sqrt(lon_gcj02 * lon_gcj02 + lat_gcj02 * lat_gcj02) + 0.00002 * math.sin(lat_gcj02 * r_pi)
    o = math.atan2(lat_gcj02, lon_gcj02) + 0.000003 * math.cos(lon_gcj02 * r_pi)
    lon_bd09 = b * math.cos(o) + 0.0065
    lat_bd09 = b * math.sin(o) + 0.006
    return [lon_bd09, lat_bd09]

# bd09 -> gcj02
# lng为bd09的经度
# lat为bd09的纬度
# 返回值为转换后坐标的列表形式，[经度, 纬度]
def bd09_gcj02(lon_bd09, lat_bd09):
    m = lon_bd09 - 0.0065
    n = lat_bd09 - 0.006
    c = math.sqrt(m * m + n * n) - 0.00002 * math.sin(n * r_pi)
    o = math.atan2(n, m) - 0.000003 * math.cos(m * r_pi)
    lon_gcj02 = c * math.cos(o)
    lat_gcj02 = c * math.sin(o)
    return [lon_gcj02, lat_gcj02]

# wgs84 -> gcj02
# lng为wgs84的经度
# lat为wgs84的纬度
# 返回值为转换后坐标的列表形式，[经度, 纬度]
def wgs84_gcj02(lon_wgs84, lat_wgs84):
    if judge_China(lon_wgs84, lat_wgs84):  # 判断是否在国内
        return [lon_wgs84, lat_wgs84]
    tlat = transformlat(lon_wgs84 - 105.0, lat_wgs84 - 35.0)
    tlng = transformlng(lon_wgs84 - 105.0, lat_wgs84 - 35.0)
    rlat = lat_wgs84 / 180.0 * pi
    m = math.sin(rlat)
    m = 1 - ob * m * m
    sm = math.sqrt(m)
    tlat = (tlat * 180.0) / ((la * (1 - ob)) / (m * sm) * pi)
    tlng = (tlng * 180.0) / (la / sm * math.cos(rlat) * pi)
    lat_gcj02 = lat_wgs84 + tlat
    lon_gcj02 = lon_wgs84 + tlng
    return [lon_gcj02, lat_gcj02]

# gcj02 -> wgs84
# lng为gcj02的经度
# lat为gcj02的纬度
# 返回值为转换后坐标的列表形式，[经度, 纬度]
def gcj02_wgs84(lon_gcj02, lat_gcj02):
    if judge_China(lon_gcj02, lat_gcj02):
        return [lon_gcj02, lat_gcj02]
    tlat = transformlat(lon_gcj02 - 105.0, lat_gcj02 - 35.0)
    tlng = transformlng(lon_gcj02 - 105.0, lat_gcj02 - 35.0)
    rlat = lat_gcj02 / 180.0 * pi
    m = math.sin(rlat)
    m = 1 - ob * m * m
    sm = math.sqrt(m)
    tlat = (tlat * 180.0) / ((la * (1 - ob)) / (m * sm) * pi)
    tlng = (tlng * 180.0) / (la / sm * math.cos(rlat) * pi)
    lat_wgs84 = 2 * lat_gcj02 - (lat_gcj02 + tlat)
    lon_wgs84 = 2 * lon_gcj02 - (lon_gcj02 + tlng)
    return [lon_wgs84, lat_wgs84]

# wgs84 -> bd09
# lng为wgs84的经度
# lat为wgs84的纬度
# 返回值为转换后坐标的列表形式，[经度, 纬度]
def wgs84_bd09(lon_wgs84, lat_wgs84):
    # 先把wgs84坐标系的坐标转换为火星坐标系
    tmpList_gcj02 = wgs84_gcj02(lon_wgs84, lat_wgs84)
    # 然后把火星坐标系的坐标转换为百度坐标系
    return gcj02_bd09(tmpList_gcj02[0], tmpList_gcj02[1])

# bd09 -> wgs84
# lng为bd09的经度
# lat为bd09的纬度
# 返回值为转换后坐标的列表形式，[经度, 纬度]
def bd09_wgs84(lon_bd09, lat_bd09):
    # 先把百度坐标系的经纬度转换为火星坐标系
    tmpList_gcj02 = bd09_gcj02(lon_bd09, lat_bd09)
    # 然后把火星坐标系的坐标转换为百度坐标系
    return gcj02_wgs84(tmpList_gcj02[0], tmpList_gcj02[1])

# 经纬度计算功能类
def transformlat(lon, lat):
    r = -100.0 + 2.0 * lon + 3.0 * lat + 0.2 * lat * lat + 0.1 * lon * lat + 0.2 * math.sqrt(math.fabs(lon))
    r += (20.0 * math.sin(6.0 * lon * pi) + 20.0 * math.sin(2.0 * lon * pi)) * 2.0 / 3.0
    r += (20.0 * math.sin(lat * pi) + 40.0 * math.sin(lat / 3.0 * pi)) * 2.0 / 3.0
    r += (160.0 * math.sin(lat / 12.0 * pi) + 320 * math.sin(lat * pi / 30.0)) * 2.0 / 3.0
    return r

def transformlng(lon, lat):
    r = 300.0 + lon + 2.0 * lat + 0.1 * lon * lon + 0.1 * lon * lat + 0.1 * math.sqrt(math.fabs(lon))
    r += (20.0 * math.sin(6.0 * lon * pi) + 20.0 * math.sin(2.0 * lon * pi)) * 2.0 / 3.0
    r += (20.0 * math.sin(lon * pi) + 40.0 * math.sin(lon / 3.0 * pi)) * 2.0 / 3.0
    r += (150.0 * math.sin(lon / 12.0 * pi) + 300.0 * math.sin(lon / 30.0 * pi)) * 2.0 / 3.0
    return r

# 简单判断坐标点是否在中国
# 不在的话返回True
# 在的话返回False
def judge_China(lon, lat):
    if lon < 70 or lon > 140:
        return True
    if lat < 0 or lat > 55:
        return True
    return False

# 坐标系转换的main函数
# 0->gcj02
# 1->wgs84
# 2->bd09
def main(lon, lat, fromCoord, toCoord):
    fromCoord = int(fromCoord)
    toCoord = int(toCoord)

    if fromCoord - toCoord != 0:
        # 新建变量
        # 进行坐标转换
        if fromCoord == 0 and toCoord == 1:
            temp = gcj02_wgs84(lon, lat)
        elif fromCoord == 0 and toCoord == 2:
            temp = gcj02_bd09(lon, lat)
        elif fromCoord == 1 and toCoord == 0:
            temp = wgs84_gcj02(lon, lat)
        elif fromCoord == 1 and toCoord == 2:
            temp = wgs84_bd09(lon, lat)
        elif fromCoord == 2 and toCoord == 0:
            temp = bd09_gcj02(lon, lat)
        elif fromCoord == 2 and toCoord == 1:
            temp = bd09_wgs84(lon, lat)
        return temp
    else:
        return [lon, lat]

if __name__ == '__main__':
    # 原坐标
    lon = 116.50987588293626
    lat = 39.82122830702036
    print("原坐标：", lon, ",", lat)
    # 将WGS-84坐标转换为GCJ-02
    result1 = wgs84_gcj02(lon, lat)
    print("WGS-84 -> GCJ-02：", result1)
    # 将GCJ-02坐标转换为WGS-84
    result2 = gcj02_wgs84(lon, lat)
    print("GCJ-02 -> WGS-84：", result2)
    # 将WGS-84坐标转换为BD-09
    result3 = wgs84_bd09(lon, lat)
    print("WGS-84 -> BD-09：", result3)
    # 将BD-09坐标转换为WGS-84
    result4 = bd09_wgs84(lon, lat)
    print("BD-09 -> WGS-84：", result4)
    # 将BD-09坐标转换为GCJ-02
    result5 = bd09_gcj02(lon, lat)
    print("BD-09 -> GCJ-02：", result5)
    # 将GCJ-02坐标转换为BD-09
    result6 = gcj02_bd09(lon, lat)
    print("GCJ-02 -> BD-09：", result6)

接下来，对上述地理坐标系转换代码的精确度进行验证，具体的验证步骤为：

1. 在A坐标系的地图产品下拾取5个坐标点
2. 利用代码对5个坐标点进行A坐标系到B坐标系的转换
3. 在B坐标系地图产品下同样拾取这5个位置的坐标点
4. 计算转换后坐标点与真值坐标点的偏移距离，并计算平均偏移距离

其中，WGS-84坐标系的地图产品是OpenStreetMap，GCJ-02坐标系的地图产品是高德地图，BD-09坐标系的地图产品是百度地图。

以下为选择的5个坐标点

以下为精确度计算过程表

通过上表可以得到：

• WGS-84 => GCJ-02的误差均值是8.9米
• WGS-84 => BD-09的误差均值是7.6米
• GCJ-02 => WGS-84的误差均值是8.2米
• GCJ-02 => BD-09的误差均值是3.9米
• BD-09 =>WGS-84的误差均值是7.2米
• BD-09 => GCJ-02的误差均值是3.9米

三个坐标系之间转换公式的误差均值都在10米以下，在能够使用的合理范围之内。由于不同地图产品之间存在一定的差异以及在拾取坐标的时候存在一定的误差，真实的误差均值可能比看到的误差还要再小。

介绍一个零代码实现地理坐标系转换的工具

很多小伙伴在工作和学习中会遇见地理坐标系转换的问题，不知道转换的原理也不知道应该用什么工具来转换。因此，我利用年终技术总结的时间把自己在工作中使用到的坐标系转换算法进行落地，开发了一个能够进行帮助小伙伴们进行地理坐标系转换的工具。

扫描下面二维码，关注地学大数据公众号，输入地理坐标转换，即可获取零代码实现地理坐标转换的工具

如果是手机用户可以保存或截图二维码 并在微信扫一扫功能中识别二维码。

该工具的目的是帮助使用者进行地理坐标系的转换，目前支持WGS-84、GCJ-02以及BD-09三个坐标系之间的相互转换。软件目前支持”复制粘贴“的方式进行数据输入，可以通过”复制粘贴“以及”输出CSV“两种形式进行输出。

本文对于地理坐标系统原理以及转换的介绍就到这里，如果有小伙伴对于其中某些知识点有疑问或者有问题探讨，请加关注私聊。

如果觉得我这篇文章对你学习有帮助，请点赞、收藏并加关注，谢谢~