Colorize Voronoi Diagram

Question

Sto provando a colorare un diagramma di Voronoi creato usando scipy.spatial.Voronoi. Ecco il mio codice:

import numpy as np 
import matplotlib.pyplot as plt 
from scipy.spatial import Voronoi, voronoi_plot_2d 

# make up data points 
points = np.random.rand(15,2) 

# compute Voronoi tesselation 
vor = Voronoi(points) 

# plot 
voronoi_plot_2d(vor) 

# colorize 
for region in vor.regions: 
    if not -1 in region: 
     polygon = [vor.vertices[i] for i in region] 
     plt.fill(*zip(*polygon)) 

plt.show() 

L'immagine risultante:

Come potete vedere alcune delle regioni Voronoi al confine dell'immagine non sono colorati. Questo perché alcuni indici dei vertici di Voronoi per queste regioni sono impostati su -1, ad esempio per quei vertici al di fuori del diagramma di Voronoi. Secondo la documentazione:

limitata: (lista di lista di numeri interi, forma (nregions, *)) Indici dei vertici Voronoi formanti ciascuna regione Voronoi. -1 indica il vertice fuori dal diagramma di Voronoi.

Per colorare queste regioni così, ho cercato di rimuovere solo questi vertici "fuori" dal poligono, ma che non ha funzionato. Penso, ho bisogno di compilare alcuni punti al confine della regione dell'immagine, ma non riesco a capire come raggiungere questo ragionevole.

Qualcuno può aiutare?

Answer 1

La struttura dati Voronoi contiene tutte le informazioni necessarie per costruire posizioni per i "punti all'infinito". Qhull li riporta anche semplicemente come indici -1, quindi Scipy non li calcola per te.

https://gist.github.com/pv/8036995

import numpy as np 
import matplotlib.pyplot as plt 
from scipy.spatial import Voronoi 

def voronoi_finite_polygons_2d(vor, radius=None): 
    """ 
    Reconstruct infinite voronoi regions in a 2D diagram to finite 
    regions. 

    Parameters 
    ---------- 
    vor : Voronoi 
     Input diagram 
    radius : float, optional 
     Distance to 'points at infinity'. 

    Returns 
    ------- 
    regions : list of tuples 
     Indices of vertices in each revised Voronoi regions. 
    vertices : list of tuples 
     Coordinates for revised Voronoi vertices. Same as coordinates 
     of input vertices, with 'points at infinity' appended to the 
     end. 

    """ 

    if vor.points.shape[1] != 2: 
     raise ValueError("Requires 2D input") 

    new_regions = [] 
    new_vertices = vor.vertices.tolist() 

    center = vor.points.mean(axis=0) 
    if radius is None: 
     radius = vor.points.ptp().max() 

    # Construct a map containing all ridges for a given point 
    all_ridges = {} 
    for (p1, p2), (v1, v2) in zip(vor.ridge_points, vor.ridge_vertices): 
     all_ridges.setdefault(p1, []).append((p2, v1, v2)) 
     all_ridges.setdefault(p2, []).append((p1, v1, v2)) 

    # Reconstruct infinite regions 
    for p1, region in enumerate(vor.point_region): 
     vertices = vor.regions[region] 

     if all(v >= 0 for v in vertices): 
      # finite region 
      new_regions.append(vertices) 
      continue 

     # reconstruct a non-finite region 
     ridges = all_ridges[p1] 
     new_region = [v for v in vertices if v >= 0] 

     for p2, v1, v2 in ridges: 
      if v2 < 0: 
       v1, v2 = v2, v1 
      if v1 >= 0: 
       # finite ridge: already in the region 
       continue 

      # Compute the missing endpoint of an infinite ridge 

      t = vor.points[p2] - vor.points[p1] # tangent 
      t /= np.linalg.norm(t) 
      n = np.array([-t[1], t[0]]) # normal 

      midpoint = vor.points[[p1, p2]].mean(axis=0) 
      direction = np.sign(np.dot(midpoint - center, n)) * n 
      far_point = vor.vertices[v2] + direction * radius 

      new_region.append(len(new_vertices)) 
      new_vertices.append(far_point.tolist()) 

     # sort region counterclockwise 
     vs = np.asarray([new_vertices[v] for v in new_region]) 
     c = vs.mean(axis=0) 
     angles = np.arctan2(vs[:,1] - c[1], vs[:,0] - c[0]) 
     new_region = np.array(new_region)[np.argsort(angles)] 

     # finish 
     new_regions.append(new_region.tolist()) 

    return new_regions, np.asarray(new_vertices) 

# make up data points 
np.random.seed(1234) 
points = np.random.rand(15, 2) 

# compute Voronoi tesselation 
vor = Voronoi(points) 

# plot 
regions, vertices = voronoi_finite_polygons_2d(vor) 
print "--" 
print regions 
print "--" 
print vertices 

# colorize 
for region in regions: 
    polygon = vertices[region] 
    plt.fill(*zip(*polygon), alpha=0.4) 

plt.plot(points[:,0], points[:,1], 'ko') 
plt.xlim(vor.min_bound[0] - 0.1, vor.max_bound[0] + 0.1) 
plt.ylim(vor.min_bound[1] - 0.1, vor.max_bound[1] + 0.1) 

plt.show() 

Answer 2

Non penso che ci siano sufficienti informazioni dai dati disponibili nella struttura di vor per capirlo senza fare almeno una parte del calcolo di voronoi. Poiché questo è il caso, qui ci sono le parti rilevanti della funzione voronoi_plot_2d originale che dovresti essere in grado di usare per estrarre i punti che si intersecano con vor.max_bound o vor.min_bound che sono gli angoli in basso a sinistra e in alto a destra del diagramma in ordina le altre coordinate per i tuoi poligoni.

for simplex in vor.ridge_vertices: 
    simplex = np.asarray(simplex) 
    if np.all(simplex >= 0): 
     ax.plot(vor.vertices[simplex,0], vor.vertices[simplex,1], 'k-') 

ptp_bound = vor.points.ptp(axis=0) 
center = vor.points.mean(axis=0) 
for pointidx, simplex in zip(vor.ridge_points, vor.ridge_vertices): 
    simplex = np.asarray(simplex) 
    if np.any(simplex < 0): 
     i = simplex[simplex >= 0][0] # finite end Voronoi vertex 

     t = vor.points[pointidx[1]] - vor.points[pointidx[0]] # tangent 
     t /= np.linalg.norm(t) 
     n = np.array([-t[1], t[0]]) # normal 

     midpoint = vor.points[pointidx].mean(axis=0) 
     direction = np.sign(np.dot(midpoint - center, n)) * n 
     far_point = vor.vertices[i] + direction * ptp_bound.max() 

     ax.plot([vor.vertices[i,0], far_point[0]], 
       [vor.vertices[i,1], far_point[1]], 'k--')