vba - excel vba中大型数据集中纬度/经度之间的最近距离

Question

初学者在这里...我正在研究这个井距项目，该项目着眼于纬度/经度并确定下一个最近的井。我想我可能正在创建一个无限循环，或者程序只是永远运行(它循环通过 15,000 行)。我的主要斗争是试图确保将每个位置与数据集中的每个位置进行比较。从那里我取第二低的距离(因为与自身相比，最低距离为零)。

Sub WellSpacing()
Dim r As Integer, c As Integer, L As Integer, lastrow As Integer
Dim lat1 As Double, lat2 As Double, long1 As Double, long2 As Double
Dim distance As Double, d1 As Double, d2 As Double, d3 As Double
Dim PI As Double

PI = Application.WorksheetFunction.PI()
L = 2
r = 3
c = 10
lastrow = Sheets("Test").Cells(Rows.Count, "J").End(xlUp).Row

For L = 2 To lastrow
    For r = 2 To lastrow
        lat1 = Sheets("Test").Cells(L, c)
        long1 = Sheets("Test").Cells(L, c + 1)
        lat2 = Sheets("Test").Cells(r, c)
        long2 = Sheets("Test").Cells(r, c + 1)
        d1 = Sin((Abs((lat2 - lat1)) * PI / 180 / 2)) ^ 2 + Cos(lat1 * PI / 180) * Cos(lat2 * PI / 180) * Sin(Abs(long2 - long1) * PI / 180 / 2) ^ 2
        d2 = 2 * Application.WorksheetFunction.Atan2(Sqr(1 - d1), Sqr(d1))
        d3 = 6371 * d2 * 3280.84
        Sheets("Working").Cells(r - 1, c - 9) = d3
    Next r

    Sheet2.Activate
    Range("A:A").Sort Key1:=Range("A1"), Order1:=xlAscending
    distance = Sheet2.Range("A2")
    Sheets("Test").Cells(L, c + 2) = distance
    Sheet2.Range("A:A").Clear
    Sheet1.Activate

Next L
End Sub

Answer 1

我最近一直在使用地理位置数学(又名 coordinate geometry)，并编写了一个子程序来完成您正在寻找的几乎相同的事情。

您的代码可能没有创建无限循环，但计算数千个坐标之间的距离可以是 处理器密集型即使是对代码的微小更改也会对处理时间产生巨大影响。

计算最近的坐标对:蛮力法

有许多算法可用于确定最近点，但最容易编码(因此可能最适合一次性使用)被称为 。蛮力法 .

For p1 = 1 to numPoints
    For p2 = p1 + 1 to numPoints
        ...calculate {distance}
        ...if {distance} < minDistance then minDist = {distance}
    Next p2
Next p1

使用此方法，将计算 之间的距离。 x * ( n - 1 ) / 2 点对。

例如， 的列表5分需要 10 个比较 :

1. Point 1 ↔ Point 2
2. Point 1 ↔ Point 3
3. Point 1 ↔ Point 4
4. Point 1 ↔ Point 5
5. Point 2 ↔ Point 3
6. Point 2 ↔ Point 4
7. Point 2 ↔ Point 5
8. Point 3 ↔ Point 4
9. Point 3 ↔ Point 5
10. Point 4 ↔ Point 5

由于额外的点会以指数方式增加执行时间，因此这种方法会产生一些冗长的处理时间，尤其是在速度较慢的机器上或点数过多的情况下。

我用于 的方法计算点之间的距离对于 比较点列表之间的距离远非 [代码重] 最有效的替代方案，但它们可以满足我的“一次性”需求。

根据我的目的，我将在 Excel 和 Access 之间切换(几乎相同的代码)，但 Access 更快，因此您可能希望将列表移动到表格中并这样做。

我比较的点列表之一是 252 项 , 这需要 31,628 个人比较 使用这种“简单代码”方法。在 Excel ，该过程在 中完成1.12 秒 ，即 访问 只需 0.16 秒 .

在我们开始处理更长的点列表之前，这似乎没有太大区别:我的另一个列表(接近你的大小)大约有 12000点 , 这需要 71,994,000 次计算 使用蛮力方法。在 访问 ，该过程在 中完成8.6 分钟 ，所以我估计需要 大约一个小时在 Excel .

当然，所有这些时间都基于我的操作系统、处理能力、Office 版本等。VBA 不适合这种级别的计算，你可以做的一切减少代码长度都会产生很大的不同，包括注释掉状态栏更新、即时窗口输出、关闭屏幕更新等。

这段代码有点凌乱且没有注释，因为我出于自己的目的将它放在一起，但它对我有用。如果您对它的工作原理有任何疑问，请告诉我。所有计算均以公制为单位，但可以轻松转换。

Sub findShortestDist_Excel()

    Const colLatitude = "C" ' Col.C = Lat, Col.D = Lon
    Dim pointList As Range, pointCount As Long, c As Range, _
        arrCoords(), x As Long, y As Long
    Dim thisDist As Double, minDist As Double, minDist_txt As String
    Dim cntCurr As Long, cntTotal As Long, timerStart As Single

    timerStart = Timer
    Set pointList = Sheets("Stops").UsedRange.Columns(colLatitude)
    pointCount = WorksheetFunction.Count(pointList.Columns(1))

    'build array of numbers found in Column C/D
    ReDim arrCoords(1 To 3, 1 To pointCount)
    For Each c In pointList.Columns(1).Cells
        If IsNumeric(c.Value) And Not IsEmpty(c.Value) Then
            x = x + 1
            arrCoords(1, x) = c.Value
            arrCoords(2, x) = c.Offset(0, 1).Value
        End If
    Next c

    minDist = -1
    cntTotal = pointCount * (pointCount + 1) / 2

    'loop through array
    For x = 1 To pointCount
        For y = x + 1 To pointCount
            If (arrCoords(1, x) & arrCoords(2, x)) <> (arrCoords(1, y) & arrCoords(2, y)) Then
                cntCurr = cntCurr + 1
                thisDist = Distance(arrCoords(1, x), arrCoords(2, x), _
                    arrCoords(1, y), arrCoords(2, y))
                'check if this distance is the smallest yet
                If ((thisDist < minDist) Or (minDist = -1)) And thisDist > 0 Then
                    minDist = thisDist
                    'minDist_txt = arrCoords(1, x) & "," & arrCoords(2, x) & " -> " & arrCoords(1, y) & "," & arrCoords(2, y)
                End If
                'Application.StatusBar = "Calculating Distances: " & Format(cntCurr / cntTotal, "0.0%")
            End If
        Next y
        'DoEvents
    Next x

    Debug.Print "Minimum distance: " & minDist_txt & " = " & minDist & " meters"
    Debug.Print "(" & Round(Timer - timerStart, 2) & "sec)"
    Application.StatusBar = "Finished.  Minimum distance: " & minDist_txt & " = " & minDist & "m"

End Sub

请注意 上述过程取决于以下 (Access 与 Excel 的版本略有不同):

Excel:计算点之间的距离

Public Function Distance(ByVal lat1 As Double, ByVal lon1 As Double, _
    ByVal lat2 As Double, ByVal lon2 As Double) As Double
'returns Meters distance in Excel (straight-line)
    Dim theta As Double: theta = lon1 - lon2
    Dim Dist As Double: Dist = Math.Sin(deg2rad(lat1)) * Math.Sin(deg2rad(lat2)) + Math.Cos(deg2rad(lat1)) * Math.Cos(deg2rad(lat2)) * Math.Cos(deg2rad(theta))
    Dist = rad2deg(WorksheetFunction.Acos(Dist))
    Distance = Dist * 60 * 1.1515 * 1.609344 * 1000
End Function

Function deg2rad(ByVal deg As Double) As Double
    deg2rad = (deg * WorksheetFunction.PI / 180#)
End Function

Function rad2deg(ByVal rad As Double) As Double
    rad2deg = rad / WorksheetFunction.PI * 180#
End Function

...以及 Microsoft Access 的替代代码:

访问:最短距离

Sub findShortestDist_Access()

    Const tableName = "Stops"
    Dim pointCount As Long, arrCoords(), x As Long, y As Long
    Dim thisDist As Double, minDist As Double
    Dim cntCurr As Long, cntTotal As Long, timerStart As Single
    Dim rs As Recordset

    timerStart = Timer

    Set rs = CurrentDb.OpenRecordset("SELECT * FROM " & tableName)
    With rs
        .MoveLast
        .MoveFirst
        pointCount = .RecordCount

        'build array of numbers found in Column C/D
        ReDim arrCoords(1 To 2, 1 To pointCount)
        Do While Not .EOF
            x = x + 1
            arrCoords(1, x) = !stop_lat
            arrCoords(2, x) = !stop_lon
            .MoveNext
        Loop
        .Close
    End With

    minDist = -1
    cntTotal = pointCount * (pointCount + 1) / 2
    SysCmd acSysCmdInitMeter, "Calculating Distances:", cntTotal

    'loop through array
    For x = 1 To pointCount
        For y = x + 1 To pointCount
                cntCurr = cntCurr + 1
                thisDist = Distance(arrCoords(1, x), arrCoords(2, x), _
                    arrCoords(1, y), arrCoords(2, y))
                'check if this distance is the smallest yet
                If ((thisDist < minDist) Or (minDist = -1)) And thisDist > 0 Then
                    minDist = thisDist
                End If
                SysCmd acSysCmdUpdateMeter, cntCurr
        Next y
        DoEvents
    Next x
    SysCmd acSysCmdRemoveMeter
    Debug.Print "Minimum distance: " & minDist_txt & " = " & minDist & " meters"
    Debug.Print "(" & Round(Timer - timerStart, 2) & "sec)"

End Sub

请注意 上述过程取决于以下 ...(Access 可以更快地处理大量计算，但我们必须自己构建一些内置在 Excel 中的函数)

访问:计算点之间的距离

Const pi As Double = 3.14159265358979

Public Function Distance(ByVal lat1 As Double, ByVal lon1 As Double, _
    ByVal lat2 As Double, ByVal lon2 As Double) As Double
'returns Meters distance in Access (straight-line)
    Dim theta As Double: theta = lon1 - lon2
    Dim dist As Double
    dist = Math.Sin(deg2rad(lat1)) * Math.Sin(deg2rad(lat2)) + Math.Cos(deg2rad(lat1)) _
        * Math.Cos(deg2rad(lat2)) * Math.Cos(deg2rad(theta))
    dist = rad2deg(aCos(dist))
    Distance = dist * 60 * 1.1515 * 1.609344 * 1000
End Function

Function deg2rad(ByVal deg As Double) As Double
    deg2rad = (deg * pi / 180#)
End Function

Function rad2deg(ByVal rad As Double) As Double
    rad2deg = rad / pi * 180#
End Function

Function aTan2(x As Double, y As Double) As Double
    aTan2 = Atn(y / x)
End Function

Function aCos(x As Double) As Double
    On Error GoTo aErr
    If x = 0 Or Abs(x) = 1 Then
        aCos = 0
    Else
        aCos = Atn(-x / Sqr(-x * x + 1)) + 2 * Atn(1)
    End If
    Exit Function
aErr:
    aCos = 0
End Function

平面案例

另一种计算更接近点的方法称为 平面案例 .我还没有看到任何现成的代码示例，而且我不需要它足够糟糕来编写它，但它的要点是:



阅读更多关于 的信息Closest pair of points problem 在维基百科上。