Getting geometries of hurricane's 'cone of uncertainty' using shapely?

Shapely geometries have a convex_hull method.

Should be as simple as polygons.convex_hull, but it will work with any Shapely geometry.

A note on cyclones as a domain: you should use the input cyclone positions as input rather than an interpolated curve: weather forecasts are typically made for moments in time, often 3, 6 or 12 hours apart, and the position in-between is uncertain simply because it is left uncalculated. The convex hull (a special kind of alpha shape) will encompass the spaces in-between the foretasted locations, exactly like in your sample images.

Also be careful with the antimeridian...

Edit: on second thought, you probably want a concave hull, or else to generate convex hulls sequentially, starting with the first pair of error shapes, then with the i+1 and i+2 pair, until complete. Then you union this set of pair-wise convex hulls together. If you do a simple convex hull, then your overall shape will be, well, convex rather than somewhat concave. But a naive concave hull may well be too "tight" and cause intrusions into the path that you don't want.

To illustrate (pseudo-code):

shapes = [a, b, c, d] # Ordered list of shapely geometries
parts = []
for first, second in zip(shapes, shapes[1:]):
    parts.append(union(first, second).convex_hull)
union(*parts)

If you need a polygon like in the image below, replace the last lines of your code by the following:

#### firstly, import Polygon class ####
from shapely.geometry import MultiPolygon, Polygon
.
.
.
# make polygon
thetas = np.linspace(0, 2 * np.pi, 360)
polygon_x = x[:,None] + r[:,None] * np.sin(thetas)
polygon_y = y[:,None] + r[:,None] * np.cos(thetas)

# circles
ps = [Polygon(i) for i in np.dstack((polygon_x, polygon_y))]

# list of convex hulls of subsequent circles
n = range(len(ps)-1)
convex_hulls = [MultiPolygon([ps[i], ps[i+1]]).convex_hull for i in n]

# Final polygon
polygons = cascaded_union(convex_hulls)

Convex hulls:

enter image description here

Final result:

enter image description here


There is an implementation of VariableWidthBuffer in the JTS Lab here. It uses a union of circles around each line vertex and "prisms" around each line segment. That could be a basis for a Python implementation.

This will make it into JTS sometime soon. Then perhaps into GEOS, where it can be exposed by Shapely.

Tags:

Python

Shapely