2

我正在尝试使用haversine公式计算两个地理坐标之间的公里距离。

代码：

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

并且代码返回20015.09，这显然是错误的。根据 Yandex 地图，它应该是 642 公里。

我哪里错了？经度和纬度的格式是否错误？

4

2 回答 2

2

据我所知，问题在于 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 使用当地地球半径而不是平均半径；我不确定。）

解释

我发现了大圆距离的半正弦公式的描述，这就是你正在实现的。该页面上的 JavaScript 实现提供了以下计算c：

var c = 2 * Math.atan2(Math.sqrt(a), Math.sqrt(1-a));

在 JavaScript 中，atan2()接受参数y, x. 但是，在 Excel VBA 中，WorksheetFunction.Atan2采用参数x, y. 您的原始代码作为第一个参数传递Sqr(dbl_a)，就像在 JavaScript 中一样。但是，Sqr(dbl_a)需要是 Excel VBA 中的第二个参数。

关于命名的评论

基于@JohnColeman 的观点，有很多方法可以命名变量。在这种情况下，我建议使用单位前缀而不是类型前缀：例如deg_Latitude1，、、RadPerDeg = Pi/180和rad_dLat = RadPerDeg * (deg_Latitude2 - deg_Latitude1)。我个人认为这有助于避免单位转换失误。

于 2017-12-28T13:23:10.610 回答
1

我的 VBA 代码以英尺为单位返回答案；然而，“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

端模块

于 2019-03-25T20:37:37.787 回答