wgs84 gc02 bd09 坐标系转换excel宏，方便大家获取

坐标系转换excel vba宏，根据python代码改编

原python代码地址：https://github.com/wandergis/coordTransform_py

由于vba精度问题，转换精度没有python 高，但基本符合要求

百度->腾讯/高德 BD09toGCJ02lang BD09toGCJ02Lat
腾讯/高德->百度 GCJ02toBD09Lng GCJ02toBD09lat
google->百度 WGS84toGCJ02Lng WGS84toGCJ02Lat

Private Const PI As Double = 3.14159265358979
Private Const x_pi = 3.14159265358979 * 3000# / 180#
Private Const X_EARTH_RADIUS = 6378245
Private Const EE = 6.69342162296594E-03

Private Type Point
lng As Double
lat As Double
End Type
Private Function transform(lng As Double, lat As Double) As Point
dlat = transformLat(lng – 105#, lat – 35#)
dlng = transformlng(lng – 105#, lat – 35#)
radlat = lat / 180# * PI
magic = Math.Sin(radlat)
magic = 1 – EE * magic * magic
sqrtmagic = Sqr(magic)
dlat = (dlat * 180#) / ((X_EARTH_RADIUS * (1 – EE)) / (magic * sqrtmagic) * PI)
dlng = (dlng * 180#) / (X_EARTH_RADIUS / sqrtmagic * Math.Cos(radlat) * PI)
mgLat = lat + dlat
mglng = lng + dlng
Dim ret As Point
ret.lng = mglng
ret.lat = mgLat
transform = ret
End Function
Private Function transformlng(x As Double, y As Double) As Double
ret = 300# + x + 2# * y + 0.1 * x * x + 0.1 * x * y + 0.1 * Sqr(Math.Abs(x))
ret = ret + (20# * Math.Sin(6# * x * PI) + 20# * Math.Sin(2# * x * PI)) * 2# / 3#
ret = ret + (20# * Math.Sin(x * PI) + 40# * Math.Sin(x / 3# * PI)) * 2# / 3#
ret = ret + (150# * Math.Sin(x / 12# * PI) + 300# * Math.Sin(x / 30# * PI)) * 2# / 3#
transformlng = ret
End Function
Private Function transformLat(x As Double, y As Double) As Double
ret = -100# + 2# * x + 3# * y + 0.2 * y * y + 0.1 * x * y + 0.2 * Sqr(Math.Abs(x))
ret = ret + (20# * Math.Sin(6# * x * PI) + 20# * Math.Sin(2# * x * PI)) * 2# / 3#
ret = ret + (20# * Math.Sin(y * PI) + 40# * Math.Sin(y / 3# * PI)) * 2# / 3#
ret = ret + (160# * Math.Sin(y / 12# * PI) + 320 * Math.Sin(y * PI / 30#)) * 2# / 3#
transformLat = ret
End Function
Private Function BD09toGCJ02(lng As Double, lat As Double) As Point
Dim p As Point
x = lng – 0.0065
y = lat – 0.006
Z = Sqr(x * x + y * y) – 0.00002 * Math.Sin(y * x_pi)
theta = WorksheetFunction.Atan2(x, y) – 0.000003 * Math.Cos(x * x_pi)
p.lng = Z * Math.Cos(theta)
p.lat = Z * Math.Sin(theta)
BD09toGCJ02 = p
End Function
Private Function GCJ02toBD09(lng As Double, lat As Double) As Point
Dim p As Point
Z = Sqr(lng * lng + lat * lat) + 0.00002 * Math.Sin(lat * x_pi)
theta = WorksheetFunction.Atan2(lng, lat) + 0.000003 * Math.Cos(lng * x_pi)
p.lng = Z * Math.Cos(theta) + 0.0065
p.lat = Z * Math.Sin(theta) + 0.006
GCJ02toBD09 = p
End Function

Private Function WGS84toGCJ02(lng As Double, lat As Double) As Point
Dim p As Point
dlat = transformLat(lng – 105#, lat – 35#)
dlng = transformlng(lng – 105#, lat – 35#)
radlat = lat / 180# * PI
magic = Math.Sin(radlat)
magic = 1 – EE * magic * magic
sqrtmagic = Sqr(magic)
dlat = (dlat * 180#) / ((X_EARTH_RADIUS * (1 – EE)) / (magic * sqrtmagic) * PI)
dlng = (dlng * 180#) / (X_EARTH_RADIUS / sqrtmagic * Math.Cos(radlat) * PI)
p.lat = lat + dlat
p.lng = lng + dlng
WGS84toGCJ02 = p
End Function

Public Function BD09toGCJ02Lng(lng As Double, lat As Double) As Double
Dim p As Point
p = BD09toGCJ02(lng, lat)
BD09toGCJ02Lng = p.lng
End Function
Public Function BD09toGCJ02Lat(lng As Double, lat As Double) As Double
Dim p As Point
p = BD09toGCJ02(lng, lat)
BD09toGCJ02Lat = p.lat
End Function

Public Function GCJ02toBD09Lat(lng As Double, lat As Double) As Double
Dim p As Point
p = GCJ02toBD09(lng, lat)
GCJ02toBD09Lat = p.lat
End Function

Public Function GCJ02toBD09Lng(lng As Double, lat As Double) As Double
Dim p As Point
p = GCJ02toBD09(lng, lat)
GCJ02toBD09Lng = p.lng
End Function

Public Function WGS84toGCJ02Lng(lng As Double, lat As Double) As Double
Dim p As Point
p = WGS84toGCJ02(lng, lat)
WGS84toGCJ02Lng = p.lng
End Function

Public Function WGS84toGCJ02Lat(lng As Double, lat As Double) As Double
Dim p As Point
p = WGS84toGCJ02(lng, lat)
WGS84toGCJ02Lat = p.lat
End Function

Public Function GCJ02toWGS84Lng(lng As Double, lat As Double) As Double
Dim p As Point
p = transform(lng, lat)
mglng = lng * 2 – p.lng
mgLat = lat * 2 – p.lat
GCJ02toWGS84Lng = mglng
End Function
Public Function GCJ02toWGS84Lat(lng As Double, lat As Double) As Double
Dim p As Point
p = transform(lng, lat)
mglng = lng * 2 – p.lng
mgLat = lat * 2 – p.lat
GCJ02toWGS84Lat = mgLat
End Function

Public Function WGS84toBD09Lat(lng As Double, lat As Double) As Double
Dim p As Point
p = WGS84toGCJ02(lng, lat)
WGS84toBD09Lat = GCJ02toBD09Lat(p.lng, p.lat)
End Function

Public Function WGS84toBD09Lng(lng As Double, lat As Double) As Double
Dim p As Point
p = WGS84toGCJ02(lng, lat)
WGS84toBD09Lng = GCJ02toBD09Lng(p.lng, p.lat)
End Function

Public Function BD09toWGS84Lat(lng As Double, lat As Double) As Double
Dim p As Point
p = BD09toGCJ02(lng, lat)
BD09toWGS84Lat = GCJ02toWGS84Lat(p.lng, p.lat)
End Function

Public Function BD09toWGS84Lng(lng As Double, lat As Double) As Double
Dim p As Point
p = BD09toGCJ02(lng, lat)
BD09toWGS84Lng = GCJ02toWGS84Lng(p.lng, p.lat)
End Function