计算坐标之间的距离(以千米为单位)
Calculating distance in kilometers between coordinates
我正在尝试使用半正弦公式计算两个地理坐标之间的公里距离。
代码:
Dim dbl_dLat As Double
Dim dbl_dLon As Double
Dim dbl_a As Double
dbl_P = WorksheetFunction.Pi / 180
dbl_dLat = dbl_P * (dbl_Latitude2 - dbl_Latitude1)
dbl_dLon = dbl_P * (dbl_Longitude2 - dbl_Longitude1)
dbl_a = Sin(dbl_dLat / 2) * Sin(dbl_dLat / 2) + Cos(dbl_Latitude1 * dbl_P) * Cos(dbl_Latitude2 * dbl_P) * Sin(dbl_dLon / 2) * Sin(dbl_dLon / 2)
dbl_Distance_KM = 6371 * 2 * WorksheetFunction.Atan2(Sqr(dbl_a), Sqr(1 - dbl_a))
我正在使用这些坐标进行测试:
dbl_Longitude1 = 55.629178
dbl_Longitude2 = 29.846686
dbl_Latitude1 = 37.659466
dbl_Latitude2 = 30.24441
而代码returns20015.09,显然是错误的。根据 Yandex 地图应该是 642 公里。
我哪里错了?经纬度格式有误吗?
据我所知,问题在于 atan2() 的参数顺序因语言而异。以下对我有用*:
Option Explicit
Public Sub Distance()
Dim dbl_Longitude1 As Double, dbl_Longitude2 As Double, dbl_Latitude1 As Double, dbl_Latitude2 As Double
dbl_Longitude1 = 55.629178
dbl_Longitude2 = 29.846686
dbl_Latitude1 = 37.659466
dbl_Latitude2 = 30.24441
Dim dbl_dLat As Double
Dim dbl_dLon As Double
Dim dbl_a As Double
Dim dbl_P As Double
dbl_P = WorksheetFunction.Pi / 180
dbl_dLat = dbl_P * (dbl_Latitude2 - dbl_Latitude1) 'to radians
dbl_dLon = dbl_P * (dbl_Longitude2 - dbl_Longitude1) 'to radians
dbl_a = Sin(dbl_dLat / 2) * Sin(dbl_dLat / 2) + _
Cos(dbl_Latitude1 * dbl_P) * Cos(dbl_Latitude2 * dbl_P) * Sin(dbl_dLon / 2) * Sin(dbl_dLon / 2)
Dim c As Double
Dim dbl_Distance_KM As Double
c = 2 * WorksheetFunction.Atan2(Sqr(1 - dbl_a), Sqr(dbl_a)) ' *** swapped arguments to Atan2
dbl_Distance_KM = 6371 * c
Debug.Print dbl_Distance_KM
End Sub
*输出:2507.26205401321
,虽然 gcmap.com 说答案是 2512 公里。这可能是一个精度问题——我认为它已经足够接近可以算作工作了。 (编辑 gcmap 也可能使用局部地球半径而不是平均半径;我不确定。)
说明
我找到了 this description 大圆距离的半正弦公式,这就是你正在实现的。该页面上的 JavaScript 实现给出了 c
:
的计算
var c = 2 * Math.atan2(Math.sqrt(a), Math.sqrt(1-a));
在JavaScript中,atan2() takes parameters y
, x
. However, in Excel VBA, WorksheetFunction.Atan2
取参数x
,y
.您的原始代码将 Sqr(dbl_a)
作为第一个参数传递,就像在 JavaScript 中一样。但是,Sqr(dbl_a)
需要是ExcelVBA中的second参数。
对命名的评论
基于@JohnColeman 的观点,有很多方法可以命名变量。在这种情况下,我建议使用单位前缀而不是类型前缀:例如 deg_Latitude1
、RadPerDeg = Pi/180
和 rad_dLat = RadPerDeg * (deg_Latitude2 - deg_Latitude1)
。我个人认为这有助于避免 unit-conversion mishaps.
我的VBA代码,returns答案以英尺为单位;然而 'd' 是以公里为单位的答案。
Imports System.Math
Module Haversine
Public Function GlobalAddressDistance(sLat1 As String, sLon1 As String, sLat2 As String, sLon2 As String) As String
Const R As Integer = 6371
Const cMetersToFeet As Single = 3.2808399
Const cKiloMetersToMeters As Integer = 1000
Dim a As Double = 0, c As Double = 0, d As Double = 0
'Convert strings to numberic double values
Dim dLat1 As Double = Val(sLat1)
Dim dLat2 As Double = Val(sLat2)
Dim dLatDiff As Double = DegreesToRadians(CDbl(sLat2) - CDbl(sLat1))
Dim dLonDiff As Double = DegreesToRadians(CDbl(sLon2) - CDbl(sLon1))
a = Pow(Sin(dLatDiff / 2), 2) + Cos(DegreesToRadians(dLat1)) * Cos(DegreesToRadians(dLat2)) * Pow(Sin(dLonDiff / 2), 2)
c = 2 * Atan2(Sqrt(a), Sqrt(1 - a))
d = R * c
'Convert kilometers to feet
Return Format((d * cKiloMetersToMeters * cMetersToFeet), "0.##").ToString
End Function
Private Function DegreesToRadians(ByVal dDegrees As Double) As Double
Return (dDegrees * PI) / 180
End Function
模块结束
我正在尝试使用半正弦公式计算两个地理坐标之间的公里距离。
代码:
Dim dbl_dLat As Double
Dim dbl_dLon As Double
Dim dbl_a As Double
dbl_P = WorksheetFunction.Pi / 180
dbl_dLat = dbl_P * (dbl_Latitude2 - dbl_Latitude1)
dbl_dLon = dbl_P * (dbl_Longitude2 - dbl_Longitude1)
dbl_a = Sin(dbl_dLat / 2) * Sin(dbl_dLat / 2) + Cos(dbl_Latitude1 * dbl_P) * Cos(dbl_Latitude2 * dbl_P) * Sin(dbl_dLon / 2) * Sin(dbl_dLon / 2)
dbl_Distance_KM = 6371 * 2 * WorksheetFunction.Atan2(Sqr(dbl_a), Sqr(1 - dbl_a))
我正在使用这些坐标进行测试:
dbl_Longitude1 = 55.629178
dbl_Longitude2 = 29.846686
dbl_Latitude1 = 37.659466
dbl_Latitude2 = 30.24441
而代码returns20015.09,显然是错误的。根据 Yandex 地图应该是 642 公里。
我哪里错了?经纬度格式有误吗?
据我所知,问题在于 atan2() 的参数顺序因语言而异。以下对我有用*:
Option Explicit
Public Sub Distance()
Dim dbl_Longitude1 As Double, dbl_Longitude2 As Double, dbl_Latitude1 As Double, dbl_Latitude2 As Double
dbl_Longitude1 = 55.629178
dbl_Longitude2 = 29.846686
dbl_Latitude1 = 37.659466
dbl_Latitude2 = 30.24441
Dim dbl_dLat As Double
Dim dbl_dLon As Double
Dim dbl_a As Double
Dim dbl_P As Double
dbl_P = WorksheetFunction.Pi / 180
dbl_dLat = dbl_P * (dbl_Latitude2 - dbl_Latitude1) 'to radians
dbl_dLon = dbl_P * (dbl_Longitude2 - dbl_Longitude1) 'to radians
dbl_a = Sin(dbl_dLat / 2) * Sin(dbl_dLat / 2) + _
Cos(dbl_Latitude1 * dbl_P) * Cos(dbl_Latitude2 * dbl_P) * Sin(dbl_dLon / 2) * Sin(dbl_dLon / 2)
Dim c As Double
Dim dbl_Distance_KM As Double
c = 2 * WorksheetFunction.Atan2(Sqr(1 - dbl_a), Sqr(dbl_a)) ' *** swapped arguments to Atan2
dbl_Distance_KM = 6371 * c
Debug.Print dbl_Distance_KM
End Sub
*输出:2507.26205401321
,虽然 gcmap.com 说答案是 2512 公里。这可能是一个精度问题——我认为它已经足够接近可以算作工作了。 (编辑 gcmap 也可能使用局部地球半径而不是平均半径;我不确定。)
说明
我找到了 this description 大圆距离的半正弦公式,这就是你正在实现的。该页面上的 JavaScript 实现给出了 c
:
var c = 2 * Math.atan2(Math.sqrt(a), Math.sqrt(1-a));
在JavaScript中,atan2() takes parameters y
, x
. However, in Excel VBA, WorksheetFunction.Atan2
取参数x
,y
.您的原始代码将 Sqr(dbl_a)
作为第一个参数传递,就像在 JavaScript 中一样。但是,Sqr(dbl_a)
需要是ExcelVBA中的second参数。
对命名的评论
基于@JohnColeman 的观点,有很多方法可以命名变量。在这种情况下,我建议使用单位前缀而不是类型前缀:例如 deg_Latitude1
、RadPerDeg = Pi/180
和 rad_dLat = RadPerDeg * (deg_Latitude2 - deg_Latitude1)
。我个人认为这有助于避免 unit-conversion mishaps.
我的VBA代码,returns答案以英尺为单位;然而 'd' 是以公里为单位的答案。
Imports System.Math
Module Haversine
Public Function GlobalAddressDistance(sLat1 As String, sLon1 As String, sLat2 As String, sLon2 As String) As String
Const R As Integer = 6371
Const cMetersToFeet As Single = 3.2808399
Const cKiloMetersToMeters As Integer = 1000
Dim a As Double = 0, c As Double = 0, d As Double = 0
'Convert strings to numberic double values
Dim dLat1 As Double = Val(sLat1)
Dim dLat2 As Double = Val(sLat2)
Dim dLatDiff As Double = DegreesToRadians(CDbl(sLat2) - CDbl(sLat1))
Dim dLonDiff As Double = DegreesToRadians(CDbl(sLon2) - CDbl(sLon1))
a = Pow(Sin(dLatDiff / 2), 2) + Cos(DegreesToRadians(dLat1)) * Cos(DegreesToRadians(dLat2)) * Pow(Sin(dLonDiff / 2), 2)
c = 2 * Atan2(Sqrt(a), Sqrt(1 - a))
d = R * c
'Convert kilometers to feet
Return Format((d * cKiloMetersToMeters * cMetersToFeet), "0.##").ToString
End Function
Private Function DegreesToRadians(ByVal dDegrees As Double) As Double
Return (dDegrees * PI) / 180
End Function
模块结束