专栏首页数据科学学习手札用Python实现WGS84、火星坐标系、百度坐标系、web墨卡托四种坐标相互转换

用Python实现WGS84、火星坐标系、百度坐标系、web墨卡托四种坐标相互转换

一、简介

  主流被使用的地理坐标系并不统一,常用的有WGS84、GCJ02(火星坐标系)、BD09(百度坐标系)以及百度地图中保存矢量信息的web墨卡托,本文利用Python编写相关类以实现4种坐标系统之间的互相转换。

二、代码及说明

import math

class LngLatTransfer():

    def __init__(self):
        self.x_pi = 3.14159265358979324 * 3000.0 / 180.0
        self.pi = math.pi  # π
        self.a = 6378245.0  # 长半轴
        self.es = 0.00669342162296594323  # 偏心率平方
        pass

    def GCJ02_to_BD09(self, gcj_lng, gcj_lat):
        """
        实现GCJ02向BD09坐标系的转换
        :param lng: GCJ02坐标系下的经度
        :param lat: GCJ02坐标系下的纬度
        :return: 转换后的BD09下经纬度
        """
        z = math.sqrt(gcj_lng * gcj_lng + gcj_lat * gcj_lat) + 0.00002 * math.sin(gcj_lat * self.x_pi)
        theta = math.atan2(gcj_lat, gcj_lng) + 0.000003 * math.cos(gcj_lng * self.x_pi)
        bd_lng = z * math.cos(theta) + 0.0065
        bd_lat = z * math.sin(theta) + 0.006
        return bd_lng, bd_lat


    def BD09_to_GCJ02(self, bd_lng, bd_lat):
        '''
        实现BD09坐标系向GCJ02坐标系的转换
        :param bd_lng: BD09坐标系下的经度
        :param bd_lat: BD09坐标系下的纬度
        :return: 转换后的GCJ02下经纬度
        '''
        x = bd_lng - 0.0065
        y = bd_lat - 0.006
        z = math.sqrt(x * x + y * y) - 0.00002 * math.sin(y * self.x_pi)
        theta = math.atan2(y, x) - 0.000003 * math.cos(x * self.x_pi)
        gcj_lng = z * math.cos(theta)
        gcj_lat = z * math.sin(theta)
        return gcj_lng, gcj_lat


    def WGS84_to_GCJ02(self, lng, lat):
        '''
        实现WGS84坐标系向GCJ02坐标系的转换
        :param lng: WGS84坐标系下的经度
        :param lat: WGS84坐标系下的纬度
        :return: 转换后的GCJ02下经纬度
        '''
        dlat = self._transformlat(lng - 105.0, lat - 35.0)
        dlng = self._transformlng(lng - 105.0, lat - 35.0)
        radlat = lat / 180.0 * self.pi
        magic = math.sin(radlat)
        magic = 1 - self.es * magic * magic
        sqrtmagic = math.sqrt(magic)
        dlat = (dlat * 180.0) / ((self.a * (1 - self.es)) / (magic * sqrtmagic) * self.pi)
        dlng = (dlng * 180.0) / (self.a / sqrtmagic * math.cos(radlat) * self.pi)
        gcj_lng = lat + dlat
        gcj_lat = lng + dlng
        return gcj_lng, gcj_lat


    def GCJ02_to_WGS84(self, gcj_lng, gcj_lat):
        '''
        实现GCJ02坐标系向WGS84坐标系的转换
        :param gcj_lng: GCJ02坐标系下的经度
        :param gcj_lat: GCJ02坐标系下的纬度
        :return: 转换后的WGS84下经纬度
        '''
        dlat = self._transformlat(gcj_lng - 105.0, gcj_lat - 35.0)
        dlng = self._transformlng(gcj_lng - 105.0, gcj_lat - 35.0)
        radlat = gcj_lat / 180.0 * self.pi
        magic = math.sin(radlat)
        magic = 1 - self.es * magic * magic
        sqrtmagic = math.sqrt(magic)
        dlat = (dlat * 180.0) / ((self.a * (1 - self.es)) / (magic * sqrtmagic) * self.pi)
        dlng = (dlng * 180.0) / (self.a / sqrtmagic * math.cos(radlat) * self.pi)
        mglat = gcj_lat + dlat
        mglng = gcj_lng + dlng
        lng = gcj_lng * 2 - mglng
        lat = gcj_lat * 2 - mglat
        return lng, lat


    def BD09_to_WGS84(self, bd_lng, bd_lat):
        '''
        实现BD09坐标系向WGS84坐标系的转换
        :param bd_lng: BD09坐标系下的经度
        :param bd_lat: BD09坐标系下的纬度
        :return: 转换后的WGS84下经纬度
        '''
        lng, lat = self.BD09_to_GCJ02(bd_lng, bd_lat)
        return self.GCJ02_to_WGS84(lng, lat)


    def WGS84_to_BD09(self, lng, lat):
        '''
        实现WGS84坐标系向BD09坐标系的转换
        :param lng: WGS84坐标系下的经度
        :param lat: WGS84坐标系下的纬度
        :return: 转换后的BD09下经纬度
        '''
        lng, lat = self.WGS84_to_GCJ02(lng, lat)
        return self.GCJ02_to_BD09(lng, lat)


    def _transformlat(self, lng, lat):
        ret = -100.0 + 2.0 * lng + 3.0 * lat + 0.2 * lat * lat + \
              0.1 * lng * lat + 0.2 * math.sqrt(math.fabs(lng))
        ret += (20.0 * math.sin(6.0 * lng * self.pi) + 20.0 *
                math.sin(2.0 * lng * self.pi)) * 2.0 / 3.0
        ret += (20.0 * math.sin(lat * self.pi) + 40.0 *
                math.sin(lat / 3.0 * self.pi)) * 2.0 / 3.0
        ret += (160.0 * math.sin(lat / 12.0 * self.pi) + 320 *
                math.sin(lat * self.pi / 30.0)) * 2.0 / 3.0
        return ret


    def _transformlng(self, lng, lat):
        ret = 300.0 + lng + 2.0 * lat + 0.1 * lng * lng + \
              0.1 * lng * lat + 0.1 * math.sqrt(math.fabs(lng))
        ret += (20.0 * math.sin(6.0 * lng * self.pi) + 20.0 *
                math.sin(2.0 * lng * self.pi)) * 2.0 / 3.0
        ret += (20.0 * math.sin(lng * self.pi) + 40.0 *
                math.sin(lng / 3.0 * self.pi)) * 2.0 / 3.0
        ret += (150.0 * math.sin(lng / 12.0 * self.pi) + 300.0 *
                math.sin(lng / 30.0 * self.pi)) * 2.0 / 3.0
        return ret

    def WGS84_to_WebMercator(self, lng, lat):
        '''
        实现WGS84向web墨卡托的转换
        :param lng: WGS84经度
        :param lat: WGS84纬度
        :return: 转换后的web墨卡托坐标
        '''
        x = lng * 20037508.342789 / 180
        y = math.log(math.tan((90 + lat) * self.pi / 360)) / (self.pi / 180)
        y = y * 20037508.34789 / 180
        return x, y

    def WebMercator_to_WGS84(self, x, y):
        '''
        实现web墨卡托向WGS84的转换
        :param x: web墨卡托x坐标
        :param y: web墨卡托y坐标
        :return: 转换后的WGS84经纬度
        '''
        lng = x / 20037508.34 * 180
        lat = y / 20037508.34 * 180
        lat = 180 / self.pi * (2 * math.atan(math.exp(lat * self.pi / 180)) - self.pi / 2)
        return lng, lat

  整个模块的使用方式可用下面的导图概括,其中每个函数都只需要传入经纬度坐标信息:

  以上就是本文的全部内容,如有笔误之处望指出!

本文参与腾讯云自媒体分享计划,欢迎正在阅读的你也加入,一起分享。

我来说两句

0 条评论
登录 后参与评论

推荐阅读

  • 「网安夜校」开课啦!多门网络安全课程开启限时优惠报名

    众志成城,共抗疫情。腾讯安全联合腾讯云大学、腾讯课堂启动「网安夜校」,为大家提供限时优惠的网络安全课程。欢迎网络安全从业者和信息安全专业学生报名参加学习,快速充电提升自我。

    腾讯安全
    安全培训腾讯云大学
  • Flink源码走读(一):Flink工程目录

    导语 | Flink已经成为未来流计算趋势,目前在很多大厂已经有了大规模的使用。最近在学习Flink源码,就想把自己学习的过程分享出来,希望能帮助到志同道合的朋友。开始阅读源码,说明读者已经对flink的基本概念有一些了解,这里就不再重复介绍Flink了。本文作为学习过程的第一章,首先对Flink的工程目录做一个解读,了解了工程下各个模块的作用,才能在遇到问题时准确定位到代码,进一步学习。

    2011aad
    大数据解决方案
  • Flink源码走读(二):Flink+Kafka实现端到端Exactly Once语义

    Flink通过Checkpoint机制实现了消息对状态影响的Exactly Once语义,即每条消息只会影响Flink内部状态有且只有一次。但无法保证输出到Sink中的数据不重复。以图一所示为例,Flink APP收到Source中的A消息,将其转化为B消息输出到Sink,APP在处理完A1后做了一次Checkpoint,假设APP在处理到A4时发生错误重启,APP将会重新从A2开始消费并处理数据,就会导致B2和B3重复输出到Sink中两次。

    2011aad
    大数据解决方案Kafka
  • kubernetes系列教程(十九)使用metric-server让HPA弹性伸缩愉快运行

    kubernetes监控指标大体可以分为两类:核心监控指标和自定义指标,核心监控指标是kubernetes内置稳定可靠监控指标,早期由heapster完成,现由metric-server实现;自定义指标用于实现核心指标的扩展,能够提供更丰富的指标支持,如应用状态指标,自定义指标需要通过Aggregator和k8s api集成,当前主流通过promethues实现。

    HappyLau谈云计算
    Kubernetes容器微服务微服务架构腾讯微服务平台 TFS
  • 三分钟入坑指北 🔜 Docsify + Serverless Framework 快速创建个人博客系统

    之前由于学摄影的关系,为了提高自己的审美,顺便锻炼下自己的英文能力,翻译了不少国外艺术类的 文章。最近一直想搭一个个人博客来存放这些内容,又懒得折腾建站,遂一直搁置。

    Aceyclee
    ServerlessHTML网站GitGitHub
  • NVM作为主存上对数据库管理系统的影响

    implications of non-volatile memory as primary storage for database management systems

    yzsDBA
    存储缓存数据库数据结构SQL
  • DevOps平台架构演进

    附最新架构图https://www.processon.com/view/5cbd897de4b0bab90962c435

    我思故我在
    DevOps 解决方案微服务架构架构设计
  • 【腾讯云AI小程序大赛】中山大学作品《小耳朵天使》

    ----------------------------------------------------------------------------------

    陈华山
    小程序 · 云开发小程序语音识别文字识别对话机器人
  • Kona JDK 在腾讯大数据领域内的实践与发展

    经常听人谈到 OpenJDK,那它到底是什么呢?相信大家都听说过 Java SE、ME、EE等规范, 通常意义上对 Open JDK 的定义指:Java SE规范的一个免费和开源参考实现。

    腾小云
    JDKJavaJVM大数据Oracle
  • 公告丨腾讯安全产品更名通知

    为了更好地为政企客户的安全保驾护航,腾讯安全即日起更新旗下身份安全、网络安全、终端安全、应用安全、数据安全、业务安全、安全管理、安全服务等八类安全产品的命名,致力于打造全栈安全产品“货架”,让客户选购安全产品/服务更加便捷,更快地找到合适的安全产品,从而对自身的安全建设“对症下药”。

    腾讯安全
    DDoS 防护应用安全 MS验证码(业务安全)应用安全(移动安全)漏洞扫描服务

扫码关注云+社区

领取腾讯云代金券