Refactor winding number

src/Meshes.jl

@@ -77,10 +77,11 @@ include("neighborsearch.jl")
 include("predicates.jl")
 
 # operations
+include("winding.jl")
 include("sideof.jl")
+include("orientation.jl")
 include("measures.jl")
 include("boundary.jl")
-include("orientation.jl")
 include("merging.jl")
 include("clipping.jl")
 include("intersections.jl")
@@ -199,7 +200,6 @@ export
   vertex,
   vertices,
   nvertices,
-  windingnumber,
   rings,
   segments,
   angles,
@@ -372,6 +372,9 @@ export
   iscollinear,
   iscoplanar,
 
+  # winding number
+  winding,
+
   # sideof
   sideof,
   SideType,
@@ -381,15 +384,6 @@ export
   LEFT,
   RIGHT,
 
-  # measures
-  measure,
-  area,
-  volume,
-  perimeter,
-
-  # boundary
-  boundary,
-
   # orientation
   OrientationMethod,
   WindingOrientation,
@@ -399,6 +393,15 @@ export
   CW,
   CCW,
 
+  # measures
+  measure,
+  area,
+  volume,
+  perimeter,
+
+  # boundary
+  boundary,
+
   # clipping
   ClippingMethod,
   SutherlandHodgman,

src/orientation.jl

@@ -62,7 +62,7 @@ function orientation(r::Ring{2,T}, ::WindingOrientation) where {T}
   # pick any segment
   x1, x2 = r.vertices[1:2]
   x̄ = center(Segment(x1, x2))
-  w = T(2π) * windingnumber(x̄, r) - ∠(x1, x̄, x2)
+  w = T(2π) * winding(x̄, r) - ∠(x1, x̄, x2)
   isapprox(w, T(π), atol=atol(T)) ? CCW : CW
 end
 

src/polytopes.jl

@@ -167,13 +167,6 @@ Return the outer and inner rings of the polygon.
 """
 function rings end
 
-"""
-    windingnumber(point, polygon)
-
-Winding number of `point` with respect to the `polygon`.
-"""
-function windingnumber end
-
 # implementations of Polygon
 include("polytopes/ngon.jl")
 include("polytopes/polyarea.jl")

src/polytopes/polyarea.jl

@@ -82,8 +82,6 @@ centroid(p::PolyArea) = centroid(first(p.rings))
 
 rings(p::PolyArea) = p.rings
 
-windingnumber(point::Point, p::PolyArea) = windingnumber(point, first(p.rings))
-
 function Base.unique!(p::PolyArea)
   foreach(unique!, p.rings)
   inds = findall(r -> nvertices(r) ≤ 2, p.rings)

src/polytopes/ring.jl

@@ -55,27 +55,6 @@ function Random.rand(rng::Random.AbstractRNG, ::Random.SamplerType{<:Ring{Dim,T}
   Ring(v)
 end
 
-"""
-    windingnumber(point, ring)
-
-Winding number of `point` with respect to the `ring`.
-The winding number is the total number of times that
-the ring travels counterclockwise around the point.
-See https://en.wikipedia.org/wiki/Winding_number.
-
-## References
-
-* Balbes, R. and Siegel, J. 1990. [A robust method for calculating
-  the simplicity and orientation of planar polygons]
-  (https://www.sciencedirect.com/science/article/abs/pii/0167839691900198)
-"""
-function windingnumber(p::Point{2,T}, r::Ring{2,T}) where {T}
-  v = r.vertices
-  n = length(v)
-  ∑ = sum(∠(v[i], p, v[i + 1]) for i in 1:n)
-  ∑ / T(2π)
-end
-
 """
     innerangles(ring)
 

src/sideof.jl

@@ -51,7 +51,7 @@ Determines on which side the `point` is in relation to the `ring`.
 Possible results are `IN` or `OUT` the `ring`.
 """
 function sideof(point::Point{2,T}, ring::Ring{2,T}) where {T}
-  w = windingnumber(point, ring)
+  w = winding(point, ring)
   ifelse(isapprox(w, zero(T), atol=atol(T)), OUT, IN)
 end
 

src/winding.jl

@@ -0,0 +1,21 @@
+# ------------------------------------------------------------------
+# Licensed under the MIT License. See LICENSE in the project root.
+# ------------------------------------------------------------------
+
+"""
+    winding(points, object)
+
+Generalized winding number of `points` with respect to the geometric `object`.
+
+## References
+
+* Barill et al. 2018. [Fast winding numbers for soups and clouds]
+  (https://dl.acm.org/doi/10.1145/3197517.3201337)
+"""
+function winding end
+
+function winding(p::Point{2,T}, r::Ring{2,T}) where {T}
+  v = vertices(r)
+  n = length(v)
+  sum(∠(v[i], p, v[i + 1]) for i in 1:n) / T(2π)
+end

test/polytopes.jl

@@ -129,13 +129,6 @@
     @test angles(c) ≈ [-atan(2), -π / 2, +π / 2, -π / 2, -π / 2, -(π - atan(2))]
     @test innerangles(c) ≈ [atan(2), π / 2, 3π / 2, π / 2, π / 2, π - atan(2)]
 
-    c = Ring(P2[(0, 0), (1, 0), (1, 1), (0, 1)])
-    @test windingnumber(P2(0.5, 0.5), c) ≈ 1
-    @test windingnumber(P2(0.5, 0.5), reverse(c)) ≈ -1
-    c = Ring(P2[(0, 0), (1, 0), (1, 1), (0, 1), (0, 0), (1, 0), (1, 1), (0, 1)])
-    @test windingnumber(P2(0.5, 0.5), c) ≈ 2
-    @test windingnumber(P2(0.5, 0.5), reverse(c)) ≈ -2
-
     c1 = Ring(P2[(0, 0), (1, 0), (1, 1), (0, 1)])
     c2 = Ring(vertices(c1))
     @test c1 == c2

test/runtests.jl

@@ -90,6 +90,7 @@ testfiles = [
   "neighborhoods.jl",
   "neighborsearch.jl",
   "predicates.jl",
+  "winding.jl",
   "sideof.jl",
   "orientation.jl",
   "merging.jl",

test/winding.jl

@@ -0,0 +1,8 @@
+@testset "winding" begin
+  c = Ring(P2[(0, 0), (1, 0), (1, 1), (0, 1)])
+  @test winding(P2(0.5, 0.5), c) ≈ 1
+  @test winding(P2(0.5, 0.5), reverse(c)) ≈ -1
+  c = Ring(P2[(0, 0), (1, 0), (1, 1), (0, 1), (0, 0), (1, 0), (1, 1), (0, 1)])
+  @test winding(P2(0.5, 0.5), c) ≈ 2
+  @test winding(P2(0.5, 0.5), reverse(c)) ≈ -2
+end