Izak Joubert
Izak Joubert

Reputation: 1001

Intersection of a Shapely polygon with a Matplotlib wedge

In this scenario I am plotting matplotlib.patches.Wedge objects and also buffered shapely.geometry.LineString objects. I need to compute the overlapping areas of these two objects. However, the Wedge is a matplotlib.wedges object and cannot be used with Shapely's .intersection() method.

How can I do this? Here is some code:

from shapely.geometry import LineString
from matplotlib.patches import Wedge
from matplotlib import pyplot as plt
from descartes.patch import PolygonPatch

width = 5
radius = 1
rich = 1

circle_patch = Wedge((0, 0), radius+3,
                     0, 360, 3)

fig, ax = plt.subplots()

ax.add_patch(circle_patch)

ax.plot(0, 0, 'xr')
plt.autoscale()

coords = [
    [0, 0],
    [0, 1],
    [0, 2],
    [1, 2],
    [2, 2]
]

stick = LineString(coords)

stick_patch = PolygonPatch(stick.buffer(0.5))

ax.add_patch(stick_patch)

x, y = stick.xy
ax.plot(x, y, 'r-', zorder=1)

plt.show()

area = stick.buffer(0.5).intersection(circle_patch).area

enter image description here

P.S. It has to be a ring shape, not a circle

Upvotes: 1

Views: 2072

Answers (2)

Georgy
Georgy

Reputation: 13717

The simplest solution would be not to work with Matplotlib patches and construct the wedge-polygon with Shapely in the first place:

import matplotlib.pyplot as plt
from descartes.patch import PolygonPatch
from shapely.geometry import LineString, Point

outer_circle = Point(0, 0).buffer(4)
inner_circle = Point(0, 0).buffer(1)
wedge = outer_circle.difference(inner_circle)

stick = LineString([(0, 0), (0, 2), (2, 2)])
buffered_stick = stick.buffer(0.5)

intersection = buffered_stick.intersection(wedge)

wedge_patch = PolygonPatch(wedge)
stick_patch = PolygonPatch(buffered_stick, alpha=0.5, hatch='/')
intersection_patch = PolygonPatch(intersection, alpha=0.5, hatch='.')

fig, ax = plt.subplots()
ax.add_patch(wedge_patch)
ax.add_patch(stick_patch)
ax.add_patch(intersection_patch)
plt.autoscale()

enter image description here


If, for some reason, this is not possible, and you have to work with the Matplotlib's Wedge, then I can think of two ways to get its intersection area with Shapely's polygon. In both of them, I convert the patches to Shapely polygons first. You probably can't get intersection area using only Matplotlib.

1) Using .get_path() method on the Matplotlib's patch from which you can extract vertices as a NumPy array and convert it to a Shapely polygon using asPolygon:

import matplotlib.pyplot as plt
from descartes.patch import PolygonPatch
from matplotlib.patches import Wedge
from shapely.geometry import asPolygon, LineString

wedge_patch = Wedge(center=(0, 0), 
                    r=4,
                    theta1=0, 
                    theta2=360, 
                    width=3)
stick = LineString([(0, 0), (0, 2), (2, 2)])
buffered_stick = stick.buffer(0.5)

wedge_path = wedge_patch.get_path()
wedge_polygon = asPolygon(wedge_path.vertices).buffer(0)
intersection = buffered_stick.intersection(wedge_polygon)

stick_patch = PolygonPatch(buffered_stick, alpha=0.5, hatch='/')
intersection_patch = PolygonPatch(intersection, alpha=0.5, hatch='.')

fig, ax = plt.subplots()
ax.add_patch(wedge_patch)
ax.add_patch(stick_patch)
ax.add_patch(intersection_patch)
plt.autoscale()

Note the buffer(0) which I apply to the wedge polygon. This is a common trick in Shapely to make a valid polygon out of an invalid. In your answer you do something similar when removing zeros from ring_coords.

2) By accessing Wedge attributes: center, r and width, and using them to recreate a polygon:

import matplotlib.pyplot as plt
from descartes.patch import PolygonPatch
from matplotlib.patches import Wedge
from shapely.geometry import LineString, Point

wedge_patch = Wedge(center=(0, 0), 
                    r=4,
                    theta1=0, 
                    theta2=360, 
                    width=3)
stick = LineString([(0, 0), (0, 2), (2, 2)])
buffered_stick = stick.buffer(0.5)

outer_circle = Point(wedge_patch.center).buffer(wedge_patch.r)
inner_circle = Point(wedge_patch.center).buffer(wedge_patch.r - wedge_patch.width)
wedge_polygon = outer_circle.difference(inner_circle)
intersection = buffered_stick.intersection(wedge_polygon)

stick_patch = PolygonPatch(buffered_stick, alpha=0.5, hatch='/')
intersection_patch = PolygonPatch(intersection, alpha=0.5, hatch='.')

fig, ax = plt.subplots()
ax.add_patch(wedge_patch)
ax.add_patch(stick_patch)
ax.add_patch(intersection_patch)
plt.autoscale()

All solutions give the same visual output.

And all methods give roughly the same area:

>>> intersection.area
3.3774012986988513  # 1st case
3.3823210603713694  # 2nd case and the original without Matplotlib

Upvotes: 1

Izak Joubert
Izak Joubert

Reputation: 1001

Figured it out. There is a ._path.vertices member of the matplotlib.patches class which gives you the array of coordinates of the wedge object which you can then use with Shapely's LinearRing class to create a Shapely object like so:

from shapely.geometry import LineString, LinearRing
from matplotlib.patches import Wedge

width = 5
radius = 1
rich = 1

circle_patch = Wedge((0, 0), radius,
                     0, 360,)

ring_coords = circle_patch._path.vertices
ring_coords = ring_coords[(ring_coords[:, 0] != 0) & (ring_coords[:, 1] != 0)]

ring = LinearRing(ring_coords)

It does however need manipulation of the coordinate array which I don't think is the most robust method but it will do for me. Also the ring is not entirely smooth but I am sure one could do some smoothing of the coordinate array with some or other Numpy or Scipy function. ring

EDIT: To create the single wedge line one must remove the width member of the wedge. This can however be re-incorporated later using Shapely's buffer() function.

Upvotes: 1

Related Questions