diff --git a/docs/src/algorithms/clipping.md b/docs/src/algorithms/clipping.md
index c2655598f..f67932e7b 100644
--- a/docs/src/algorithms/clipping.md
+++ b/docs/src/algorithms/clipping.md
@@ -28,6 +28,30 @@ other = Box((0,1), (3,7))
 # clipped polygon
 clipped = clip(poly, other, SutherlandHodgmanClipping())
 
+viz(poly)
+viz!(other, color = :black, alpha = 0.2)
+viz!(boundary(clipped), color = :red, segmentsize = 3)
+Mke.current_figure()
+```
+
+## Weiler-Atherton
+
+```@docs
+WeilerAthertonClipping
+```
+
+```@example clipping
+# polygon to clip
+outer = Ring((8, 0), (4, 8), (2, 8), (-2, 0), (0, 0), (1, 2), (5, 2), (6, 0))
+inner = Ring((4, 4), (2, 4), (3, 6))
+poly = PolyArea([outer, inner])
+
+# clipping geometry can be a polygon with holes
+other = poly |> Rotate(pi) |> Translate(8.0, 8.0)
+
+# clipped polygon
+clipped = clip(poly, other, WeilerAthertonClipping())
+
 viz(poly)
 viz!(other, color = :black, alpha = 0.2)
 viz!(boundary(clipped), color = :red, segmentsize = 3)
diff --git a/src/Meshes.jl b/src/Meshes.jl
index 2d2064416..24630bf7f 100644
--- a/src/Meshes.jl
+++ b/src/Meshes.jl
@@ -434,6 +434,7 @@ export
   # clipping
   ClippingMethod,
   SutherlandHodgmanClipping,
+  WeilerAthertonClipping,
   clip,
 
   # intersections
diff --git a/src/clipping.jl b/src/clipping.jl
index aa9e1560d..0aab27182 100644
--- a/src/clipping.jl
+++ b/src/clipping.jl
@@ -21,3 +21,4 @@ function clip end
 # ----------------
 
 include("clipping/sutherlandhodgman.jl")
+include("clipping/weileratherton.jl")
diff --git a/src/clipping/weileratherton.jl b/src/clipping/weileratherton.jl
new file mode 100644
index 000000000..a566cbb13
--- /dev/null
+++ b/src/clipping/weileratherton.jl
@@ -0,0 +1,408 @@
+# ------------------------------------------------------------------
+# Licensed under the MIT License. See LICENSE in the project root.
+# ------------------------------------------------------------------
+
+"""
+    WeilerAthertonClipping()
+
+The Wieler-Atherton algorithm for clipping polygons.
+
+## References
+
+* Weiler, K., & Atherton, P. 1977. [Hidden surface removal using polygon area sorting]
+  (https://dl.acm.org/doi/pdf/10.1145/563858.563896)
+"""
+struct WeilerAthertonClipping <: ClippingMethod end
+
+# VertexType has three different types of vertices used for deciding how to navigate the data
+# structure when collecting the polygon rings after clipping.
+abstract type VertexType end
+abstract type Normal <: VertexType end
+abstract type Entering <: VertexType end
+abstract type Exiting <: VertexType end
+
+# Data structure for clipping the polygons. Fields left and right are used depending on the
+# VertexType. If Normal, left points to the following vertex of the original ring, and right
+# to the next  intersection vertex on the edge between point and left.point, or the following
+# original ring vertex if no intersections on the edge. For Entering and Exiting types, left
+# poitns to the following vertex on the clipping ring and right to the following vertex on one
+# of the clipped rings.
+mutable struct RingVertex{VT<:VertexType,M<:Manifold,C<:CRS}
+  point::Point{M,C}
+  left::RingVertex{VT,M,C} where {VT<:VertexType}
+  right::RingVertex{VT,M,C} where {VT<:VertexType}
+
+  function RingVertex{VT,M,C}(point::Point{M,C}) where {VT<:VertexType,M<:Manifold,C<:CRS}
+    v = new(point)
+    v.left = v
+    v.right = v
+    return v
+  end
+end
+
+RingVertex{VT}(point::Point{M,C}) where {VT<:VertexType,M<:Manifold,C<:CRS} =
+  RingVertex{VT,manifold(point),typeof(point.coords)}(point)
+
+function appendvertices!(v1::RingVertex{Normal}, v2::RingVertex{Normal})
+  v2.left = v1.left
+  v1.left = v2
+  v2.right = v1.right
+  v1.right = v2
+end
+
+function gettraversalring(ring::Ring)
+  vs = vertices(ring)
+  start = RingVertex{Normal}(vs[1])
+
+  for v in vs[2:end]
+    new = RingVertex{Normal}(v)
+    appendvertices!(start, new)
+    start = new
+  end
+
+  return start.left
+end
+
+isnormal(_::RingVertex{Normal}) = true
+isnormal(_::RingVertex) = false
+
+function insertintersection!(head::RingVertex, intersection::RingVertex, side::Symbol)
+  tail = head.left
+  vertex = head
+
+  new = measure(Segment(head.point, intersection.point))
+  while true
+    os = isnormal(vertex) ? :right : side
+    current = measure(Segment(head.point, getfield(vertex, os).point))
+    if (new < current) || (getfield(vertex, os) == tail)
+      next = getfield(vertex, os)
+      setfield!(intersection, side, next)
+      setfield!(vertex, os, intersection)
+      break
+    end
+    vertex = getfield(vertex, os)
+  end
+end
+
+# Traversing and selecting following vertices for collecting of the rings after clipping depends
+# on the vertex type. Helper nextvertex follows the correct path.
+nextvertex(v::RingVertex{Normal}) = v.right
+nextvertex(v::RingVertex{Entering}) = v.right
+nextvertex(v::RingVertex{Exiting}) = v.left
+
+function collectclipped(entering::Vector{RingVertex{Entering}})
+  rings::Vector{Ring} = []
+  visited::Vector{RingVertex} = []
+  for i in eachindex(entering)
+    if entering[i] in visited
+      continue
+    end
+
+    ring::Vector{RingVertex} = []
+    vertex = entering[i]
+    while !(vertex in ring)
+      if vertex in visited
+        break
+      end
+
+      push!(ring, vertex)
+      push!(visited, vertex)
+      vertex = nextvertex(vertex)
+    end
+
+    # Remove duplicates.
+    newring::Vector{RingVertex} = [ring[1]]
+    for i in 2:length(ring)
+      if !isapprox(ring[i].point, newring[end].point)
+        push!(newring, ring[i])
+      end
+    end
+    ring = newring
+
+    # Polygon might start several vertices after the first collected.
+    # This generally happens when there are overlapping edges that lead
+    # to several entering vertices without exiting in between. Then, the
+    # actual polygon is found by discarding the extra vertices before the
+    # proper loop.
+    k = findfirst(x -> ring[end].point == x.point, ring[1:(end - 1)])
+    if !isnothing(k)
+      ring = ring[(k + 1):end]
+    end
+
+    if length(ring) > 2
+      push!(rings, Ring([r.point for r in ring]...))
+    end
+  end
+  return rings
+end
+
+function decide_direction(clipped₁, clipped₂, clipping₁, clipping₂)
+  # The input segments should all share one common vertex.
+  # Form vectors from the common vertex to other vertices.
+  a = minimum(clipped₁) - maximum(clipped₁)
+  b = maximum(clipped₂) - minimum(clipped₂)
+  c = minimum(clipping₁) - maximum(clipping₁)
+  d = maximum(clipping₂) - minimum(clipping₂)
+
+  α = mod(atan(a[2], a[1]), 2pi) # Angle formed by the clipped₁ segment and a horizontal line.
+  β = mod(atan(b[2], b[1]) - α, 2pi) # Angle between clipped₁ and clipped₂.
+  γ = mod(atan(c[2], c[1]) - α, 2pi) # Angle between clipped₁ and clipping₁.
+  δ = mod(atan(d[2], d[1]) - α, 2pi) # Angle between clipped₁ and clipping₂.
+
+  if isapprox(γ, 0.0) || isapprox(γ, 2pi)
+    if δ < β < 2pi
+      return Entering
+    else
+      return Exiting
+    end
+  elseif isapprox(δ, 0.0) || isapprox(δ, 2pi)
+    if γ < β < 2pi
+      return Exiting
+    else
+      return Entering
+    end
+  end
+  if γ < δ && (γ < β || isapprox(β, γ)) && (isapprox(β, δ) || β < δ)
+    return Exiting
+  elseif δ < γ && (δ < β || isapprox(β, δ)) && (isapprox(β, γ) || β < γ)
+    return Entering
+  end
+  return nothing
+end
+
+function insertintersections!(clipping, clipped, point, vtype, entering)
+  intersection = RingVertex{vtype}(point)
+  insertintersection!(clipping, intersection, :left)
+  insertintersection!(clipped, intersection, :right)
+
+  if vtype == Entering
+    push!(entering, intersection)
+  end
+end
+
+function clip(poly::Polygon, ring::Ring, _::WeilerAthertonClipping)
+  polyrings = rings(poly)
+
+  # Convert the subject polygon rings and the clipping ring to the RingVertex data structure.
+  clippedrings = [gettraversalring(r) for r in polyrings]
+  startclipping = gettraversalring(ring)
+
+  # For keeping track of intersected rings, as the non-intersected ones need additional
+  # processing at the end.
+  intersected = zeros(Bool, length(clippedrings))
+
+  # For collecting all entering vertices, as they are used as starting points for collection
+  # of the output rings.
+  entering::Vector{RingVertex{Entering}} = []
+
+  clipping = startclipping
+  while true
+    # Three consecutive clipping vertices are used to properly identify intersection vertex type
+    # for corner cases, so the clipping segment is constructed with the two following vertices
+    # after the current one.
+    clippingsegment = Segment(clipping.left.point, clipping.left.left.point)
+
+    for (k, startclipped) in enumerate(clippedrings)
+      clipped = startclipped
+      while true
+        # Like for the clipping, the clipped also uses three consecutive vertices.
+        clippedsegment = Segment(clipped.left.point, clipped.left.left.point)
+
+        vertex_type = nothing
+        I = intersection(clippingsegment, clippedsegment)
+
+        # First try to handle Crossing, EdgeTouching and CornerTouching intersections as they might
+        # add only a single intersection.
+
+        if type(I) == Crossing
+          point = get(I)
+          cl = Line(clipping.left.point, clipping.left.left.point)
+          vertex_type = sideof(clipped.left.left.point, cl) == LEFT ? Entering : Exiting
+        elseif type(I) == EdgeTouching
+          point = get(I)
+          if isapprox(point, clipped.left.point)
+            # When intersection is at the shared vertex of two edges of the clipped ring,
+            # then split the interscting edge of the clipping ring at the intersection point.
+            vertex_type = decide_direction(
+              Segment(clipped.point, clipped.left.point),
+              clippedsegment,
+              Segment(clipping.left.point, point),
+              Segment(point, clipping.left.left.point)
+            )
+          elseif isapprox(point, clipping.left.point)
+            # When intersection is at the shared vertex of two edges of the clipping ring,
+            # then split the interscting edge of the clipped ring at the intersection point.
+            vertex_type = decide_direction(
+              Segment(clipped.left.point, point),
+              Segment(point, clipped.left.left.point),
+              Segment(clipping.point, clipping.left.point),
+              clippingsegment
+            )
+          end
+        elseif type(I) == CornerTouching
+          # When intersection is at the shared vertices of both the clipping and the clipped rings.
+          point = get(I)
+          if isapprox(point, clipped.left.point) && isapprox(point, clipping.left.point)
+            # Only applies if the intersection coincides with the middles of the currently observed
+            # vertices for both clipping and clipped rings.
+            vertex_type = decide_direction(
+              Segment(clipped.point, point),
+              Segment(point, clipped.left.left.point),
+              Segment(clipping.point, point),
+              Segment(point, clipping.left.left.point)
+            )
+          end
+        end
+        if !isnothing(vertex_type)
+          insertintersections!(clipping.left, clipped.left, point, vertex_type, entering)
+          intersected[k] = true
+        end
+
+        # Overlapping intersections might add up to two intersections, so handle separately.
+
+        if type(I) == Overlapping
+          for point in extrema(get(I))
+            # For both ends of the intersecting segment, check if it coincides with the middle of
+            # the obeserved vertices for clipped and clipping rings. If it does, attempt adding a
+            # point.
+
+            if isapprox(point, clipped.left.point)
+              vertex_type = decide_direction(
+                Segment(clipped.point, point),
+                Segment(point, clipped.left.left.point),
+                Segment(clipping.left.point, point),
+                Segment(point, clipping.left.left.point)
+              )
+              if !isnothing(vertex_type)
+                insertintersections!(clipping.left, clipped.left, point, vertex_type, entering)
+                intersected[k] = true
+              end
+            end
+
+            if isapprox(point, clipping.left.point)
+              vertex_type = decide_direction(
+                Segment(clipped.left.point, point),
+                Segment(point, clipped.left.left.point),
+                Segment(clipping.point, point),
+                Segment(point, clipping.left.left.point)
+              )
+              if !isnothing(vertex_type)
+                insertintersections!(clipping.left, clipped.left, point, vertex_type, entering)
+                intersected[k] = true
+              end
+            end
+          end
+        end
+
+        clipped = clipped.left
+        if clipped.left == startclipped.left
+          break
+        end
+      end
+    end
+
+    clipping = clipping.left
+    if clipping.left == startclipping.left
+      break
+    end
+  end
+
+  collectedrings::typeof(polyrings) = collectclipped(entering)
+
+  # When no interesction have been registered with any of the rings, the clipping ring either
+  # encircles everything or is completely contained in the clipped polygon.
+  if all(isequal(false), intersected)
+    o = orientation(ring)
+    contained = [v ∈ PolyArea(ring) for v in vertices(polyrings[1])]
+    if (o == CCW && all(contained))
+      return poly
+    end
+    if (o == CW && !all(contained))
+      push!(collectedrings, polyrings[1])
+      intersected[1] = true
+      for poly_ring in polyrings[2:end]
+        if !all([v ∈ PolyArea(ring) for v in vertices(poly_ring)])
+          push!(collectedrings, poly_ring)
+        end
+      end
+      if any([v ∈ PolyArea(polyrings[1]) for v in vertices(ring)])
+        push!(collectedrings, ring)
+      end
+      return PolyArea(collectedrings)
+    end
+    if all([v ∈ PolyArea(polyrings[1]) for v in vertices(ring)])
+      push!(collectedrings, ring)
+    end
+  end
+
+  # Handle all the non-intersected rings. Add them if they are contained in the clipping ring.
+  polys::Vector{PolyArea} = []
+  for r in collectedrings
+    newpolyrings = [r]
+    valid = true
+    for k in eachindex(intersected)
+      if !intersected[k]
+        if orientation(ring) == CCW
+          # Check if majority of vertices are inside the clipping ring.
+          ins = count(isequal(true), [v ∈ PolyArea(ring) for v in vertices(polyrings[k])])
+          if ins > (length(vertices(polyrings[k])) - ins)
+            if orientation(polyrings[k]) == CCW
+              pushfirst!(newpolyrings, polyrings[k])
+            else
+              push!(newpolyrings, polyrings[k])
+            end
+            intersected[k] = true
+          elseif vertices(ring)[1] ∈ PolyArea(polyrings[k]) && k > 1
+            # If the ring is contained within one of the inner rings, invalidate it.
+            valid = false
+          end
+        else
+          # Check if majority of vertices are outside the clipping ring.
+          ins = count(isequal(true), [v ∈ PolyArea(ring) for v in vertices(polyrings[k])])
+          if ins < (length(vertices(polyrings[k])) - ins)
+            if orientation(polyrings[k]) == CCW
+              pushfirst!(newpolyrings, polyrings[k])
+            else
+              push!(newpolyrings, polyrings[k])
+            end
+            intersected[k] = true
+          end
+        end
+      end
+    end
+    if valid
+      push!(polys, PolyArea(newpolyrings))
+    end
+  end
+
+  n = length(polys)
+  out = n == 0 ? nothing : (n == 1 ? polys[1] : GeometrySet(polys))
+  return out
+end
+
+function clip(poly::Polygon, other::Geometry, method::WeilerAthertonClipping)
+  b = boundary(other)
+  if typeof(b) <: Multi
+    for r in parent(b)
+      poly = clip(poly, r, method)
+      if isnothing(poly)
+        return nothing
+      end
+    end
+    return poly
+  else
+    return clip(poly, b, method)
+  end
+end
+
+function clip(gs::GeometrySet, other::Geometry, method::WeilerAthertonClipping)
+  polys::Vector{PolyArea} = []
+  for i in 1:nelements(gs)
+    poly = clip(element(gs, i), other, method)
+    if !isnothing(poly)
+      push!(polys, poly)
+    end
+  end
+  n = length(polys)
+  return n == 0 ? nothing : (n == 1 ? polys[1] : GeometrySet(polys))
+end
diff --git a/src/intersections/polygons.jl b/src/intersections/polygons.jl
index 908c9fd13..8ea00d714 100644
--- a/src/intersections/polygons.jl
+++ b/src/intersections/polygons.jl
@@ -3,14 +3,7 @@
 # ------------------------------------------------------------------
 
 function intersection(f, poly₁::Polygon, poly₂::Polygon)
-  # TODO: use Weiler-Atherton or other more general clipping method
-  clipped = if isconvex(poly₂)
-    clip(poly₁, poly₂, SutherlandHodgmanClipping())
-  elseif isconvex(poly₁)
-    clip(poly₂, poly₁, SutherlandHodgmanClipping())
-  else
-    throw(ErrorException("intersection not implemented between two non-convex polygons"))
-  end
+  clipped = clip(poly₁, poly₂, WeilerAthertonClipping())
 
   if isnothing(clipped)
     @IT NotIntersecting nothing f
diff --git a/test/clipping.jl b/test/clipping.jl
index 117f8a8ce..66992bc5c 100644
--- a/test/clipping.jl
+++ b/test/clipping.jl
@@ -74,3 +74,164 @@
   @test all(vertices(crings[1]) .≈ [cart(6, 4), cart(6, 5), cart(1, 2.5), cart(1, 1), cart(5, 2)])
   @test all(vertices(crings[2]) .≈ [cart(3.0, 3.0), cart(3.0, 3.5), cart(10 / 3, 11 / 3), cart(4.0, 3.0)])
 end
+
+@testitem "WeilerAtherton" setup = [Setup] begin
+  # triangle
+  poly = Triangle(cart(6, 2), cart(3, 5), cart(0, 2))
+  other = Quadrangle(cart(5, 0), cart(5, 4), cart(0, 4), cart(0, 0))
+  clipped = clip(poly, other, WeilerAthertonClipping())
+  @test issimple(clipped)
+  @test all(vertices(clipped) .≈ [cart(2, 4), cart(0, 2), cart(5, 2), cart(5, 3), cart(4, 4)])
+
+  # octagon
+  poly = Octagon(cart(8, -2), cart(8, 5), cart(2, 5), cart(4, 3), cart(6, 3), cart(4, 1), cart(2, 1), cart(2, -2))
+  other = Quadrangle(cart(5, 0), cart(5, 4), cart(0, 4), cart(0, 0))
+  clipped = clip(poly, other, WeilerAthertonClipping())
+  @test length(clipped) == 2
+  @test issimple(clipped[1])
+  @test issimple(clipped[2])
+  out = vertices.(clipped)
+  @test all(out[1] .≈ [cart(3, 4), cart(4, 3), cart(5, 3), cart(5, 4)])
+  @test all(out[2] .≈ [cart(5, 2), cart(4, 1), cart(2, 1), cart(2, 0), cart(5, 0)])
+
+  # inside
+  poly = Quadrangle(cart(1, 0), cart(1, 1), cart(0, 1), cart(0, 0))
+  other = Quadrangle(cart(5, 0), cart(5, 4), cart(0, 4), cart(0, 0))
+  clipped = clip(poly, other, WeilerAthertonClipping())
+  @test issimple(clipped)
+  @test all(vertices(clipped) .≈ [cart(1, 0), cart(1, 1), cart(0, 1), cart(0, 0)])
+
+  # outside
+  poly = Quadrangle(cart(7, 6), cart(7, 7), cart(6, 7), cart(6, 6))
+  other = Quadrangle(cart(5, 0), cart(5, 4), cart(0, 4), cart(0, 0))
+  clipped = clip(poly, other, WeilerAthertonClipping())
+  @test isnothing(clipped)
+
+  # surrounded
+  poly = Hexagon(cart(0, 2), cart(-2, 2), cart(-2, 0), cart(0, -2), cart(2, -2), cart(2, 0))
+  other = Hexagon(cart(1, 0), cart(0, 1), cart(-1, 1), cart(-1, 0), cart(0, -1), cart(1, -1))
+  clipped = clip(poly, other, WeilerAthertonClipping())
+  @test issimple(clipped)
+  @test all(vertices(clipped) .≈ [cart(1, 0), cart(0, 1), cart(-1, 1), cart(-1, 0), cart(0, -1), cart(1, -1)])
+
+  # PolyArea with box
+  outer = Ring(cart(8, 0), cart(4, 8), cart(2, 8), cart(-2, 0), cart(0, 0), cart(1, 2), cart(5, 2), cart(6, 0))
+  inner = Ring(cart(4, 4), cart(2, 4), cart(3, 6))
+  poly = PolyArea([outer, inner])
+  other = Box(cart(0, 1), cart(3, 7))
+  clipped = clip(poly, other, WeilerAthertonClipping())
+  crings = rings(clipped)
+  @test issimple(clipped)
+  @test all(
+    vertices(crings[1]) .≈ [
+      cart(3.0, 4.0),
+      cart(2.0, 4.0),
+      cart(3.0, 6.0),
+      cart(3.0, 7.0),
+      cart(1.5, 7.0),
+      cart(0.0, 4.0),
+      cart(0.0, 1.0),
+      cart(0.5, 1.0),
+      cart(1.0, 2.0),
+      cart(3.0, 2.0)
+    ]
+  )
+
+  # PolyArea with outer ring outside and inner ring inside
+  outer = Ring(cart(8, 0), cart(2, 6), cart(-4, 0))
+  inner = Ring(cart(1, 3), cart(3, 3), cart(3, 1), cart(1, 1))
+  poly = PolyArea([outer, inner])
+  other = Quadrangle(cart(4, 4), cart(0, 4), cart(0, 0), cart(4, 0))
+  clipped = clip(poly, other, WeilerAthertonClipping())
+  @test !issimple(clipped)
+  crings = rings(clipped)
+  @test all(vertices(crings[1]) .≈ [cart(4, 0), cart(4, 4), cart(0, 4), cart(0, 0)])
+  @test all(vertices(crings[2]) .≈ [cart(1, 3), cart(3, 3), cart(3, 1), cart(1, 1)])
+
+  # PolyArea with one inner ring inside `other` and another inner ring outside `other`
+  outer = Ring(cart(6, 4), cart(6, 7), cart(1, 6), cart(1, 1), cart(5, 2))
+  inner₁ = Ring(cart(3, 3), cart(3, 4), cart(4, 3))
+  inner₂ = Ring(cart(2, 5), cart(2, 6), cart(3, 5))
+  poly = PolyArea([outer, inner₁, inner₂])
+  other = PolyArea(Ring(cart(6, 1), cart(7, 2), cart(6, 5), cart(0, 2), cart(1, 1)))
+  clipped = clip(poly, other, WeilerAthertonClipping())
+  crings = rings(clipped)
+  @test issimple(clipped)
+  @test length(crings) == 1
+  @test all(
+    vertices(crings[1]) .≈ [
+      cart(1, 2.5),
+      cart(1, 1),
+      cart(5, 2),
+      cart(6, 4),
+      cart(6, 5),
+      cart(10 / 3, 11 / 3),
+      cart(4.0, 3.0),
+      cart(3.0, 3.0),
+      cart(3.0, 3.5)
+    ]
+  )
+
+  # Two polygons with holes.
+  outer = Ring((8, 0), (4, 8), (2, 8), (-2, 0), (0, 0), (1, 2), (5, 2), (6, 0))
+  inner = Ring((4, 4), (2, 4), (3, 6))
+  poly = PolyArea([outer, inner])
+  other = poly |> Rotate(pi) |> Translate(8.0, 8.0)
+  clipped = clip(poly, other, WeilerAthertonClipping())
+  crings = rings(clipped)
+  @test all(
+    vertices(crings[1]) .≈
+    [cart(7, 2), cart(5, 6), cart(3, 6), cart(2, 8), cart(1, 6), cart(3, 2), cart(5, 2), cart(6, 0)]
+  )
+  @test all(vertices(crings[2]) .≈ [cart(4, 4), cart(2, 4), cart(3, 6)])
+  @test all(vertices(crings[3]) .≈ [cart(4, 4), cart(6, 4), cart(5, 2)])
+
+  # Overlapping clipping polygon with the edges of the hole of the subject polygon.
+  triangle = Triangle((0.0, 0.0), (1.0, 0.0), (0.5, 1.0))
+  rectangle = Quadrangle((0.0, 0.0), (1.0, 0.0), (1.0, 1.0), (0.0, 1.0))
+  poly = PolyArea(rings(rectangle)[1], reverse(rings(triangle |> Scale(0.5) |> Translate((0.3, 0.3)))[1]))
+  other = triangle |> Scale(0.9) |> Translate(0.1, -0.1)
+  clipped = clip(poly, other, WeilerAthertonClipping())
+  crings = rings(clipped)
+  @test length(crings) == 1
+  @test all(vertices(crings[1]) .≈ [cart(0.8, 0.3), cart(0.3, 0.3), cart(0.15, 0.0), cart(0.95, 0.0)])
+
+  # Intersection only at a corner of two polygons with holes.
+  other = poly |> Rotate(pi)
+  clipped = clip(poly, other, WeilerAthertonClipping())
+  @test isnothing(clipped)
+
+  # Proper intersection of two polygons with holes.
+  other = poly |> Rotate(pi) |> Translate(1.0, 1.0)
+  clipped = clip(poly, other, WeilerAthertonClipping())
+  crings = rings(clipped)
+  @test length(crings) == 2
+  @test all(vertices(crings[1]) .≈ [cart(0.0, 1.0), cart(0.0, 0.0), cart(1.0, 0.0), cart(1.0, 1.0), cart(0.0, 1.0)])
+  @test all(
+    vertices(crings[2]) .≈ [
+      cart(0.4, 0.3),
+      cart(0.3, 0.3),
+      cart(0.35, 0.4),
+      cart(0.2, 0.7),
+      cart(0.5, 0.7),
+      cart(0.55, 0.8),
+      cart(0.6, 0.7),
+      cart(0.7, 0.7),
+      cart(0.65, 0.6),
+      cart(0.8, 0.3),
+      cart(0.5, 0.3),
+      cart(0.45, 0.2)
+    ]
+  )
+
+  # Clipping a GeometrySet.
+  poly = Quadrangle(cart(0.0, 0.0), cart(1.0, 0.0), cart(1.0, 1.0), cart(0.0, 1.0))
+  polyset = GeometrySet([poly, poly |> Translate(2.0, 0.0)])
+  other = Quadrangle(cart(0.5, 0.25), cart(2.5, 0.25), cart(2.5, 0.75), cart(0.5, 0.75))
+  clipped = clip(polyset, other, WeilerAthertonClipping())
+  @test typeof(clipped) <: GeometrySet
+  crings = rings.(clipped)
+  @test length(crings) == 2
+  @test all(vertices(crings[1][1]) .≈ [cart(1.0, 0.25), cart(1.0, 0.75), cart(0.5, 0.75), cart(0.5, 0.25)])
+  @test all(vertices(crings[2][1]) .≈ [cart(2.0, 0.75), cart(2.0, 0.25), cart(2.5, 0.25), cart(2.5, 0.75)])
+end
diff --git a/test/discretization.jl b/test/discretization.jl
index ebd728bac..d5342470a 100644
--- a/test/discretization.jl
+++ b/test/discretization.jl
@@ -398,7 +398,7 @@ end
   @test eltype(mesh) <: Ngon
   @test nvertices.(mesh) ⊆ [3, 4]
 
-  poly = PolyArea(cart.([(0, 0), (0, 1), (1, 2), (2, 1), (2, 0)]))
+  poly = PolyArea(cart.([(0, 0), (2, 0), (2, 1), (1, 2), (0, 1)]))
   mesh = discretize(poly, RegularDiscretization(50))
   @test mesh isa Meshes.SubGrid
   grid = parent(mesh)
diff --git a/test/intersections.jl b/test/intersections.jl
index edcc1c43f..4b1eee270 100644
--- a/test/intersections.jl
+++ b/test/intersections.jl
@@ -1149,22 +1149,22 @@ end
   poly = Triangle(cart(6, 2), cart(3, 5), cart(0, 2))
   other = Quadrangle(cart(5, 0), cart(5, 4), cart(0, 4), cart(0, 0))
   @test intersection(poly, other) |> type == Intersecting
-  @test all(vertices(poly ∩ other) .≈ [cart(5, 3), cart(4, 4), cart(2, 4), cart(0, 2), cart(5, 2)])
+  @test all(vertices(poly ∩ other) .≈ [cart(2, 4), cart(0, 2), cart(5, 2), cart(5, 3), cart(4, 4)])
 
   # octagon
   poly = Octagon(cart(8, -2), cart(8, 5), cart(2, 5), cart(4, 3), cart(6, 3), cart(4, 1), cart(2, 1), cart(2, -2))
   other = Quadrangle(cart(5, 0), cart(5, 4), cart(0, 4), cart(0, 0))
   @test intersection(poly, other) |> type == Intersecting
-  @test all(
-    vertices(poly ∩ other) .≈
-    [cart(3, 4), cart(4, 3), cart(5, 3), cart(5, 2), cart(4, 1), cart(2, 1), cart(2, 0), cart(5, 0), cart(5, 4)]
-  )
+  out = vertices.(poly ∩ other)
+  @test length(out) == 2
+  @test all(out[1] .≈ [cart(3, 4), cart(4, 3), cart(5, 3), cart(5, 4)])
+  @test all(out[2] .≈ [cart(5, 2), cart(4, 1), cart(2, 1), cart(2, 0), cart(5, 0)])
 
   # inside
   poly = Quadrangle(cart(1, 0), cart(1, 1), cart(0, 1), cart(0, 0))
   other = Quadrangle(cart(5, 0), cart(5, 4), cart(0, 4), cart(0, 0))
   @test intersection(poly, other) |> type == Intersecting
-  @test all(vertices(poly ∩ other) .≈ vertices(poly))
+  @test all(vertices(poly ∩ other) .≈ [cart(1, 0), cart(1, 1), cart(0, 1), cart(0, 0)])
 
   # outside
   poly = Quadrangle(cart(7, 6), cart(7, 7), cart(6, 7), cart(6, 6))
@@ -1176,7 +1176,7 @@ end
   quad = Quadrangle(cart(0, 0), cart(0.1, 0.0), cart(0.1, 0.1), cart(0.0, 0.1))
   poly = PolyArea(cart(0, 0), cart(2, 0), cart(1, 1), cart(1, 0.5))
   @test intersection(quad, poly) |> type == Intersecting
-  @test all(vertices(quad ∩ poly) .≈ [cart(0, 0), cart(0.1, 0), cart(0.1, 0.05)])
+  @test all(vertices(quad ∩ poly) .≈ [cart(0.1, 0.05), cart(0, 0), cart(0.1, 0)])
 end
 
 @testitem "Domain intersection" setup = [Setup] begin
diff --git a/test/viewing.jl b/test/viewing.jl
index 91c856316..f674bbda9 100644
--- a/test/viewing.jl
+++ b/test/viewing.jl
@@ -91,7 +91,7 @@
   @test linds[10, 6] ∈ indices(grid, multi)
 
   # clipping
-  tri = Triangle(cart(-4, 10), cart(5, 19), cart(5, 1))
+  tri = Triangle(cart(-4, 10), cart(5, 1), cart(5, 19))
   grid = cartgrid(20, 20)
   linds = LinearIndices(size(grid))
   @test linds[3, 10] ∈ indices(grid, tri)