VB.NET 实现的 Douglas-Peucker 算法，无递归

• 下载源代码 - 1.5 KB

引言

这篇关于 Douglas-Peucker 算法的文章是对 Code Project 上已发表的两篇文章的补充：Polyline Simplification（多线段简化），作者是 Elmar de Koning，以及 A C# Implementation of Douglas Peucker（Douglas Peucker 的 C# 实现），作者是 CraigSelbert。以下代码是对 Craig Selbert 代码的一个大致的转换，该代码简洁明了。

我不太喜欢 Douglas-Peucker 算法中递归的部分，因为它可能会导致程序超出堆栈限制，超出托管代码的控制范围。正如墨菲定律所说：任何可能出错的事情，都会在最糟糕的时刻出错。

因此，除了 Craig 的代码之外，我编写了 VB.NET 代码来实现 Douglas-Peucker 算法，使用迭代而不是递归。我在这里发布代码，不做过多解释。有关该算法的优秀解释，请参阅 Elmar de Koning 的文章。

我在这段代码片段中使用的类 PointD（以及许多其他方法）是一个辅助类（请参阅我的文章：2D Polyline Vertex Smoothing（二维多线段顶点平滑）），代码位于本文末尾。

代码

public 函数 GetSimplifyDouglasPeucker 是迭代的包装器。相对于 Craig 的代码，我用三个方法代替了他的 PerpendicularDistance 方法，分别计算三角形的面积、底边和高度。

Public Function GetSimplifyDouglasPeucker(points As List(Of PointD), _
    tolerance As Double) As List(Of PointD)
    If points Is Nothing OrElse points.Count < 3 Then Return points 'nothing to do
    
    Dim pointfirst As Integer = 0
    Dim pointlast As Integer = points.Count-1
    
    'first and last point cannot be the same
    While points(pointfirst).IsCoincident(points(pointlast))
        pointlast -= 1
    End While

    Dim nli As New List(Of Integer)
    nli.Add(pointfirst)
    nli.Add(pointlast)
    
    'loop DPT
    Dim farPoint, newFarPoint As Integer
    Dim busyI As Boolean = True
    
    While busyI
        
        newFarPoint = 0
        For i As Integer = 0 To nli.Count-2
            farPoint = iterateDPT(points,nli(i),nli(i+1),tolerance)
            If farPoint > newFarPoint Then
                newFarPoint = farPoint
                'we can skip the rest of the for next loop
                Exit For
            End If
        Next
        
        If newFarPoint > 0 Then
            nli.Add(newFarPoint)
            nli.Sort()
        Else
            busyI = False
        End If
        
    End While
    
    'process and output the results
    Dim nl As New List(Of PointD)
    nli.Sort()
    For Each idx As Integer In nli
        nl.Add(points(idx))
    Next
    Return nl
End Function

Private Function iterateDPT(points As List(Of PointD), _
firstpoint As Integer, lastpoint As Integer, tolerance As Double) As Integer
    Dim maxDist As Double = 0
    Dim idxFarthest As Integer = 0
    Dim distance As Double = 0
    
    For i As Integer = firstpoint To lastpoint
        distance = GetTriangleHeight(points(firstPoint),points(lastPoint),points(i))
        If distance > maxDist Then
            maxDist = distance
            idxFarthest = i
        End If
    Next
    
    'the iteration should stop if either the tolerance is not exceeded or
    'when there are no more far (outlier) points
    'this is indicated by the output being either 0 or a valid farpoint
    If maxDist > tolerance AndAlso idxFarthest <> 0 Then
        Return idxFarthest
    Else
        Return 0
    End If
End Function

Public Function GetTriangleArea(base1 As PointD, base2 As PointD, apex As PointD) As Double
    Try
        Return Math.Abs(.5*(base1.X*base2.Y + base2.X*apex.Y + _
        apex.X*base1.Y - base2.X*base1.Y - apex.X*base2.Y - base1.X*apex.Y))
    Catch
        Return Double.NaN
    End Try
End Function

Public Function GetTriangleBase(base1 As PointD, base2 As PointD, apex As PointD) As Double
    Try
        Return ((base1.X-base2.X)^2 + (base1.Y-base2.Y)^2)^0.5
    Catch
        Return Double.NaN
    End Try
End Function

Public Function GetTriangleHeight(base1 As PointD, base2 As PointD, apex As PointD) As Double
    Return GetTriangleArea(base1,base2,apex)/GetTriangleBase(base1,base2,apex)*2
End Function

以及辅助类 PointD 的代码

Public Class PointD
    
    Public Sub New()
    End Sub
    
    Public Sub New(nx As Double, ny As Double)
        X = nx
        Y = ny
    End Sub
    
    Public Sub New(p As Veet.Geometry.PointD)
        X = p.X
        Y = p.Y
    End Sub
    
    Public X As Double = 0
    Public Y As Double = 0
    
    Public Function ToStringXY() As String
        Return ToStringXY(""," ",False)
    End Function
    
    Public Function ToStringXY(fmt As String) As String
        Return ToStringXY(fmt," ",True)
    End Function
    
    Public Function ToStringXY(fmt As String, seperator As String, withLabel As Boolean) As String
        
        Dim xstr As String = X.ToString
        Dim ystr As String = Y.ToString
        
        If fmt <> "" Then
            xstr = X.ToString(fmt)
            ystr = Y.ToString(fmt)
        End If
        
        If withLabel Then
            xstr = "X=" & xstr
            ystr = "Y=" & ystr
        End If
        
        Return String.Format("{0}{1}{2}",xstr,seperator,ystr)
    End Function
    
    Public Shared Operator +(p1 As PointD, p2 As PointD) As PointD
        Return New PointD(p1.X+p2.X,p1.Y+p2.Y)
    End Operator
    
    Public Shared Operator +(p As PointD, d As Double) As PointD
        Return New PointD(p.X+d,p.Y+d)
    End Operator
    
    Public Shared Operator +(d As Double, p As PointD) As PointD
        Return p+d
    End Operator
    
    Public Shared Operator -(p1 As PointD, p2 As PointD) As PointD
        Return New PointD(p1.X-p2.X,p1.Y-p2.Y)
    End Operator
    
    Public Shared Operator -(p As PointD, d As Double) As PointD
        Return New PointD(p.X-d,p.Y-d)
    End Operator
    
    Public Shared Operator -(d As Double, p As PointD) As PointD
        Return p-d
    End Operator
    
    Public Shared Operator *(p1 As PointD, p2 As PointD) As PointD
        Return New PointD(p1.X*p2.X,p1.Y*p2.Y)
    End Operator
    
    Public Shared Operator *(p As PointD, d As Double) As PointD
        Return New PointD(p.X*d,p.Y*d)
    End Operator
    
    Public Shared Operator *(d As Double, p As PointD) As PointD
        Return p*d
    End Operator
    
    Public Shared Operator /(p1 As PointD, p2 As PointD) As PointD
        Return New PointD(p1.X/p2.X,p1.Y/p2.Y)
    End Operator
    
    Public Shared Operator /(p As PointD, d As Double) As PointD
        Return New PointD(p.X/d,p.Y/d)
    End Operator
    
    Public Shared Operator /(d As Double, p As PointD) As PointD
        Return New PointD(d/p.X,d/p.Y)
    End Operator
End Class

就这样了...

历史

• 编写于 2016 年 4 月