国家大地2000坐标系投影坐标转高德地图经纬度坐标-【高斯投影坐标】

      由于工作需要将大地2000坐标转换后显示到高德地图上。查找了种种资料都没有一篇完美的,走了很多弯路,最终在学习了各种专业知识后,完成了转化,将走的弯路一一记录。
       一开始不懂国家大地2000坐标系,也不懂为什么坐标是这么大数值的。上来就是百度一下大地坐标系2000转wgs84坐标,一堆都是安利proj4的。但是最后实操失败,转化出来的经纬度要么无穷,要么是负的。这里试了很多模型,当时也不大理解各种模型的差异,反正就是浪费了很多时间。类似代码如下,有兴趣的可以自己试一试,网址如下:https://epsg.io/transform 也是大佬安利的,但是我就没有换算成功过。=,= 这里有n多种不同的坐标系。 在这个网站上试错的时候发现了区域的概念和三分度六分度的概念以及经纬度和米制的概念。

        // 定义CGCS2000坐标系
        proj4.defs('EPSG:4490', '+proj=longlat +ellps=GRS80 +no_defs')
        // 定义UTM投影坐标系统
        proj4.defs('UTM', '+proj=utm +zone=50 +datum=WGS84 +units=m +no_defs')
        // 定义高斯-克吕格投影坐标系统
        proj4.defs('EPSG:4522', '+proj=tmerc +lat_0=0 +lon_0=114 +k=1 +x_0=38500000 +y_0=0 +ellps=GRS80 +units=m +no_defs +type=crs')
        // 定义WGS84坐标系
        proj4.defs('EPSG:4326', '+proj=longlat +datum=WGS84 +no_defs')
        const wgs84Coords = proj4('EPSG:4490', 'EPSG:4326', dadilos)
 

     经过查资料得知,我已有的坐标是以m位单位的【投影坐标】(啊西,一开始就搞错了),且y如果有8位的话,前两位表示区域,而且我现有的投影坐标是三分度划分区域的,这个区域的概念之前也是不懂得,这些都需要去看下相关的专业介绍。这里就不说明了,毕竟专业知识有限,只懂了一些皮毛。

      以上的试错完毕以后才感觉自己可能要接近真相了,但是还是没找到合适的转换方式,最后看到了【高斯投影正反算法】(在到高斯投影之前还误入歧途去了UTM投影,后来查到资料说大地2000的投影坐标是用的高斯投影),即将开启正确的大门啊。 原来的博客链接找不到了,写文章的时候找了个类似的链接,他们讲的比较详细:

高斯正反算—投影坐标转大地坐标、大地坐标转投影坐标(附有完整代码及测试结果)-CSDN博客

所以所以,就采用自己计算换算率,但是这里还有一个坑,网上找的案例他们用的参数是克氏椭球(Krasovsky椭球),导致我换算出来的经纬度还是有十几公里误差吧。

再一次研究了一下参数置换的方式,最后的最后终于得到了看着比较正经的经纬度坐标。哦默默,我已经搞了这个玩意儿两天半了!!!!!!贴一下最终的代码!!!!!具体为什么这么算,就自行百度公式了!

getlatlng (utmCoords) {
      let L0 = 0
      let B = 0
      let L = 0
      const x = utmCoords[0]
      let y = utmCoords[1]
      const tmZone = Math.floor(y / 1000000)
      // console.log(tmZone)
      const p = 206264.80625
      L0 = tmZone * 3
      for (let i = 1; utmCoords[1] / i >= 10; i = i * 10) {
        y = utmCoords[1] - Math.floor(utmCoords[1] / i) * i - 500000
        // 中央经线3°×当地带号(适用于1∶1万地形图)。
      }
      const bt = x / 6378137.0 * p
      const BT = x / 6378137.0
      const c3 = Math.cos(BT) * Math.cos(BT)
      const c4 = Math.sin(BT) * Math.cos(BT)
      const Bf = (bt + (50221746 + (293622 + (2350 + 22 * c3) * c3) * c3) * c4 * Math.pow(10, -10) * p) / p
      const c5 = Math.pow(Math.cos(Bf), 2)
      const c6 = Math.sin(Bf) * Math.cos(Bf)
      const Nf = 6399596.652 - (21562.267 - (108.973 - 0.612 * c5) * c5) * c5
      const Z = y / (Nf * Math.cos(Bf))
      const b2 = (0.5 + 0.003369 * c5) * c6
      const b3 = 0.333333 - (0.166667 - 0.001123 * c5) * c5
      const b4 = 0.25 + (0.16161 + 0.00562 * c5) * c5
      const b5 = 0.2 - (0.1667 - 0.0088 * c5) * c5
      const z2 = Math.pow(Z, 2)
      B = (Bf * p - (1 - (b4 - 0.12 * z2) * z2) * z2 * b2 * p) / 3600.0
      L = L0 + ((1 - (b3 - b5 * z2) * z2) * Z * p) / 3600.0
      // console.log('计算高斯投影坐标得大地坐标系:纬度:' + B + '    经度:' + L + '    中央子午线:' + L0)
      // const tmZone = Math.floor((los[0]) / 3) + 1
      return [ L, B ]
    },

最后的最后,这玩意换算出来的是WGS84经纬度坐标,高德地图经纬度坐标还需要再通过高德地图api再换算一下。

AMap.convertFrom(newloslist, 'gps', (status, result) => {
            if (result.info === 'ok') {
              const path = result.locations
              console.log('高德经纬度列表' + result.locations)
              var text = new AMap.Text({
                text: this.info.ProjectName, // 标记显示的文本内容
                anchor: 'left', // 设置文本标记锚点位置
                draggable: true, // 是否可拖拽
                cursor: 'pointer', // 指定鼠标悬停时的鼠标样式。
                angle: 0, // 点标记的旋转角度
                style: {
                  padding: '.75rem 1.25rem',
                  'margin-bottom': '1rem',
                  'border-radius': '.25rem',
                  'background-color': 'white',
                  width: 'auto',
                  'border-width': 0,
                  'box-shadow': '0 2px 6px 0 rgba(114, 124, 245, .5)',
                  'text-align': 'center',
                  'font-size': '18px',
                  color: 'black',
                },
                position: path[0], // 点标记在地图上显示的位置
              })
              // 添加点击事件监听器
              text.on('click', () => {
                this.map.setZoomAndCenter(15, text.getPosition()) // 放大地图并定位到该点
              })
              text.setMap(this.map) // 将文本标记设置到地图上
              geocoder.getLocation(this.info.Location, (status, result) => {
                if (status === 'complete' && result.info === 'OK') {
                  // console.log(result.geocodes[0].location)
                  this.location = [ result.geocodes[0].location.lng, result.geocodes[0].location.lat ]
                  console.log('高德根据地址搜索得坐标' + this.location)
                  const marker = new AMap.Marker({
                    position: this.location,
                    draggable: true,
                  })
                  marker.setTitle(this.info.Location)
                  marker.setMap(this.map)
                }
              })
              const polygon = new AMap.Polygon({
                path: path,
                strokeColor: '#FF33FF', // 描边颜色
                strokeStyle: 'dashed', // 轮廓线样式
                strokeWeight: 6, // 描边宽度
                strokeOpacity: 0.2, // 描边透明度
                fillColor: '#1791fc', // 填充颜色
                fillOpacity: 0.35, // 填充透明度
                zIndex: 100,
              })
              this.map.setCenter(path[0]) // 同时设置地图层级与中心点
              // 将多边形添加到地图
              this.map.add(polygon) 
              this.map.setFitView([ polygon ])
            } 
})

下图是我坐标换算输出的数据对比。

这里贴一个坐标转换验证官网网址:煤炭地质云,这里可以将坐标进行转换验证。通过计算算出来的经纬度坐标与官方换算出来的经纬度基本上相差无几了,精确度在小数点后三位。我已经很满意了!!!!!!!

作者:devil_zhuzhu

物联沃分享整理
物联沃-IOTWORD物联网 » 国家大地2000坐标系投影坐标转高德地图经纬度坐标-【高斯投影坐标】

发表回复