OpenLayers 高德地图偏移问题

OpenLayers 高德地图偏移问题

1. 偏移产生的原因

地球上同一个地理位置的经纬度，在不同的坐标系中，会有少许偏移，国内目前常见的坐标系主要分为三种：

1. 地球坐标系——WGS84：常见于 GPS 设备，Google 地图等国际标准的坐标体系。
2. 火星坐标系——GCJ-02：中国国内使用的被强制加密后的坐标体系，高德坐标就属于该种坐标体系。
3. 百度坐标系——BD-09：百度地图所使用的坐标体系，是在火星坐标系的基础上又进行了一次加密处理。

2. 解决思路

web地图坐标系常用的一般是：WGS 84(EPSG:4326)、WGS 84 / Pseudo-Mercator(EPSG:3857)、China Geodetic Coordinate System 2000(EPSG:4490) 等。
既然高德地图偏移产生的原因是由于坐标系不一致，那么将坐标系统一一下即可解决。

3. 代码

OpenLayers地图坐标以EPSG:3857为例：

1. 注册坐标系，定义坐标转换方法，gcj02Mecator.js

// 导入proj控件,使用其方法注入gcj02坐标系
import {
  Projection,
  addProjection,
  addCoordinateTransforms,
  getTransform,
} from 'ol/proj'

const PI = Math.PI
const AXIS = 6378245.0
const OFFSET = 0.00669342162296594323 // (a^2 - b^2) / a^2

function delta(wgLon, wgLat){
  let dLat = transformLat(wgLon - 105.0, wgLat - 35.0)
  let dLon = transformLon(wgLon - 105.0, wgLat - 35.0)
  const radLat = (wgLat / 180.0) * PI
  let magic = Math.sin(radLat)
  magic = 1 - OFFSET * magic * magic
  const sqrtMagic = Math.sqrt(magic)
  dLat = (dLat * 180.0) / (((AXIS * (1 - OFFSET)) / (magic * sqrtMagic)) * PI)
  dLon = (dLon * 180.0) / ((AXIS / sqrtMagic) * Math.cos(radLat) * PI)
  return [dLon, dLat]
}

function outOfChina(lon, lat){
  if (lon < 72.004 || lon > 137.8347) {
    return true
  }
  if (lat < 0.8293 || lat > 55.8271) {
    return true
  }
  return false
}

function transformLat(x, y){
  let ret = -100.0 + 2.0 * x + 3.0 * y + 0.2 * y * y + 0.1 * x * y + 0.2 * Math.sqrt(Math.abs(x))
  ret += ((20.0 * Math.sin(6.0 * x * PI) + 20.0 * Math.sin(2.0 * x * PI)) * 2.0) / 3.0
  ret += ((20.0 * Math.sin(y * PI) + 40.0 * Math.sin((y / 3.0) * PI)) * 2.0) / 3.0
  ret += ((160.0 * Math.sin((y / 12.0) * PI) + 320 * Math.sin((y * PI) / 30.0)) * 2.0) / 3.0
  return ret
}

function transformLon(x, y){
  let ret = 300.0 + x + 2.0 * y + 0.1 * x * x + 0.1 * x * y + 0.1 * Math.sqrt(Math.abs(x))
  ret += ((20.0 * Math.sin(6.0 * x * PI) + 20.0 * Math.sin(2.0 * x * PI)) * 2.0) / 3.0
  ret += ((20.0 * Math.sin(x * PI) + 40.0 * Math.sin((x / 3.0) * PI)) * 2.0) / 3.0
  ret += ((150.0 * Math.sin((x / 12.0) * PI) + 300.0 * Math.sin((x / 30.0) * PI)) * 2.0) / 3.0
  return ret
}

const forEachPoint = function(func){
  return function(input, opt_output, opt_dimension){
    const len = input.length
    const dimension = opt_dimension || 2
    let output
    if (opt_output) {
      output = opt_output
    } else if (dimension !== 2) {
      output = input.slice()
    } else {
      output = new Array(len)
    }
    for (let offset = 0; offset < len; offset += dimension) {
      func(input, output, offset)
    }
    return output
  }
}

const gcj02 = {
  toWGS84: forEachPoint(function(input, output, offset){
    let lng = input[offset]
    let lat = input[offset + 1]
    if (!outOfChina(lng, lat)) {
      const deltaD = delta(lng, lat)
      lng = lng - deltaD[0]
      lat = lat - deltaD[1]
    }
    output[offset] = lng
    output[offset + 1] = lat
  }),
  fromWGS84: forEachPoint(function(input, output, offset){
    let lng = input[offset]
    let lat = input[offset + 1]
    if (!outOfChina(lng, lat)) {
      const deltaD = delta(lng, lat)
      lng = lng + deltaD[0]
      lat = lat + deltaD[1]
    }
    output[offset] = lng
    output[offset + 1] = lat
  }),
}

const projzh = {
  ll2smerc: getTransform('EPSG:4326', 'EPSG:3857'),
  smerc2ll: getTransform('EPSG:3857', 'EPSG:4326'),
}

// wgs84 -> gcj02
projzh.ll2gmerc = function(input, opt_output, opt_dimension){
  const output = gcj02.fromWGS84(input, opt_output, opt_dimension)
  return projzh.ll2smerc(output, output, opt_dimension)
}
// gcj02 -> wgs84
projzh.gmerc2ll = function(input, opt_output, opt_dimension){
  const output = projzh.smerc2ll(input, input, opt_dimension)
  return gcj02.toWGS84(output, opt_output, opt_dimension)
}
// 3857 -> gcj02
projzh.smerc2gmerc = function(input, opt_output, opt_dimension){
  let output = projzh.smerc2ll(input, input, opt_dimension)
  output = gcj02.fromWGS84(output, output, opt_dimension)
  return projzh.ll2smerc(output, output, opt_dimension)
}
// gcj02 -> 3857
projzh.gmerc2smerc = function(input, opt_output, opt_dimension){
  let output = projzh.smerc2ll(input, input, opt_dimension)
  output = gcj02.toWGS84(output, output, opt_dimension)
  return projzh.ll2smerc(output, output, opt_dimension)
}

// 定义GCJ02墨卡托投影坐标系
const gcj02Mecator = new Projection({
  code: 'GCJ-02',
  // extent不能缺少，OpenLayers frame 渲染的时候需要
  extent: [-20037508.342789244, -20037508.342789244, 20037508.342789244, 20037508.342789244],
  units: 'm',
})

// 将GCJ02墨卡托投影坐标系注册进OpenLayers
addProjection(gcj02Mecator)

// 覆盖默认的转换方法
addCoordinateTransforms('EPSG:4326', gcj02Mecator, projzh.ll2gmerc, projzh.gmerc2ll)
addCoordinateTransforms('EPSG:3857', gcj02Mecator, projzh.smerc2gmerc, projzh.gmerc2smerc)

export default gcj02Mecator

2. 创建Layer

import { Map, View } from 'ol'
import { Tile as TileLayer } from 'ol/layer'
import { XYZ as XYZSource } from 'ol/source'
import gcj02Mecator from './gcj02Mecator'

const amapLayer = new TileLayer({
  source: new XYZSource({
    projection: gcj02Mecator,
    url: 'http://webst0{1-4}.is.autonavi.com/appmaptile?style=6&x={x}&y={y}&z={z}&key=你的key',
  }),
})

const map = new Map({
  target: 'mapContainerId',
  view: new View({
    center: [0, 0],
    zoom: 1,
  }), // projection 默认是 'EPSG:3857
  layers: [amapLayer],
})