from a set of points I built the Voronoi tessellation using scipy:
from scipy.spatial import Voronoi
vor = Voronoi(points)
Now I would like to build a Polygon in Shapely from the regions the Voronoi algorithm created. The problem is that the Polygon class requires a list of counter-clockwise vertices. Although I know how to order these vertices, I can't solve the problem because often this is my result:
(overlapping polygon). This is the code (ONE RANDOM EXAMPLE):
def order_vertices(l):
mlat = sum(x[0] for x in l) / len(l)
mlng = sum(x[1] for x in l) / len(l)
# https://stackoverflow.com/questions/1709283/how-can-i-sort-a-coordinate-list-for-a-rectangle-counterclockwise
def algo(x):
return (math.atan2(x[0] - mlat, x[1] - mlng) + 2 * math.pi) % 2*math.pi
l.sort(key=algo)
return l
a = np.asarray(order_vertices([(9.258054711746084, 45.486245994138976),
(9.239284166975443, 45.46805963143515),
(9.271640747003861, 45.48987234571072),
(9.25828782103321, 45.44377372506324),
(9.253993275176263, 45.44484395950612),
(9.250114174032936, 45.48417979682819)]))
plt.plot(a[:,0], a[:,1])
How can I solve this problem?
If you're just after a collection of polygons you don't need to pre-order the point to build them.
The scipy.spatial.Voronoi
object has a ridge_vertices
attribute containing indices of vertices forming the lines of the Voronoi ridge. If the index is -1
then the ridge goes to infinity.
First start with some random points to build the Voronoi object.
import numpy as np
from scipy.spatial import Voronoi, voronoi_plot_2d
import shapely.geometry
import shapely.ops
points = np.random.random((10, 2))
vor = Voronoi(points)
voronoi_plot_2d(vor)
You can use this to build a collection of Shapely LineString objects.
lines = [
shapely.geometry.LineString(vor.vertices[line])
for line in vor.ridge_vertices
if -1 not in line
]
The shapely.ops
module has a polygonize
that returns a generator for Shapely Polygon objects.
for poly in shapely.ops.polygonize(lines):
#do something with each polygon
Or if you wanted a single polygon formed from the region enclosed by the Voronoi tesselation you can use the Shapely unary_union
method:
shapely.ops.unary_union(list(shapely.ops.polygonize(lines)))