Saltar al contenido

Obtener las coordenadas de un polígono acotado de las celdas de Voronoi

Siéntete en la libertad de compartir nuestra página y códigos con tus amigos, ayúdanos a hacer crecer nuestra comunidad.

Solución:

Dado un cuadro delimitador rectangular, mi primera idea fue definir una especie de operación de intersección entre este cuadro delimitador y el diagrama de Voronoï producido por scipy.spatial.Voronoi. Una idea no necesariamente genial, ya que requiere codificar una gran cantidad de funciones básicas de geometría computacional.

Sin embargo, aquí está la segunda idea (¿truco?) que me vino a la mente: los algoritmos para calcular el diagrama de Voronoï de un conjunto de n puntos en el plano tienen una complejidad temporal de O(n ln(n)). ¿Qué pasa con la adición de puntos para restringir las celdas de Voronoï de los puntos iniciales para que se encuentren en el cuadro delimitador?

Solución para un diagrama de Voronoï acotado

Una imagen vale un gran discurso:

ingrese la descripción de la imagen aquí

¿Qué hice aquí? ¡Es bastante simple! Los puntos iniciales (en azul) se encuentran en [0.0, 1.0] x [0.0, 1.0]. Luego obtengo los puntos (en azul) a la izquierda (es decir, [-1.0, 0.0] x [0.0, 1.0]) por una simetría de reflexión según x = 0.0 (borde izquierdo del cuadro delimitador). Con simetrías de reflexión según x = 1.0, y = 0.0 y y = 1.0 (otros bordes del cuadro delimitador), obtengo todos los puntos (en azul) que necesito para hacer el trabajo.

entonces corro scipy.spatial.Voronoi. La imagen anterior muestra el diagrama de Voronoï resultante (yo uso scipy.spatial.voronoi_plot_2d).

¿Qué hacer a continuación? Simplemente filtre puntos, bordes o caras de acuerdo con el cuadro delimitador. Y obtenga el centroide de cada cara de acuerdo con la conocida fórmula para calcular el centroide de un polígono. Aquí hay una imagen del resultado (los centroides están en rojo):

ingrese la descripción de la imagen aquí

Un poco de diversión antes de mostrar el código

¡Genial! Parece funcionar. ¿Qué pasa si después de una iteración trato de volver a ejecutar el algoritmo en los centroides (en rojo) en lugar de los puntos iniciales (en azul)? ¿Qué pasa si lo intento una y otra vez?

Paso 2

ingrese la descripción de la imagen aquí

Paso 10

ingrese la descripción de la imagen aquí

Paso 25

ingrese la descripción de la imagen aquí

¡Frio! Las células de Voronoï tienden a minimizar su energía

Aquí está el código

import matplotlib.pyplot as pl
import numpy as np
import scipy as sp
import scipy.spatial
import sys

eps = sys.float_info.epsilon

n_towers = 100
towers = np.random.rand(n_towers, 2)
bounding_box = np.array([0., 1., 0., 1.]) # [x_min, x_max, y_min, y_max]

def in_box(towers, bounding_box):
    return np.logical_and(np.logical_and(bounding_box[0] <= towers[:, 0],
                                         towers[:, 0] <= bounding_box[1]),
                          np.logical_and(bounding_box[2] <= towers[:, 1],
                                         towers[:, 1] <= bounding_box[3]))


def voronoi(towers, bounding_box):
    # Select towers inside the bounding box
    i = in_box(towers, bounding_box)
    # Mirror points
    points_center = towers[i, :]
    points_left = np.copy(points_center)
    points_left[:, 0] = bounding_box[0] - (points_left[:, 0] - bounding_box[0])
    points_right = np.copy(points_center)
    points_right[:, 0] = bounding_box[1] + (bounding_box[1] - points_right[:, 0])
    points_down = np.copy(points_center)
    points_down[:, 1] = bounding_box[2] - (points_down[:, 1] - bounding_box[2])
    points_up = np.copy(points_center)
    points_up[:, 1] = bounding_box[3] + (bounding_box[3] - points_up[:, 1])
    points = np.append(points_center,
                       np.append(np.append(points_left,
                                           points_right,
                                           axis=0),
                                 np.append(points_down,
                                           points_up,
                                           axis=0),
                                 axis=0),
                       axis=0)
    # Compute Voronoi
    vor = sp.spatial.Voronoi(points)
    # Filter regions
    regions = []
    for region in vor.regions:
        flag = True
        for index in region:
            if index == -1:
                flag = False
                break
            else:
                x = vor.vertices[index, 0]
                y = vor.vertices[index, 1]
                if not(bounding_box[0] - eps <= x and x <= bounding_box[1] + eps and
                       bounding_box[2] - eps <= y and y <= bounding_box[3] + eps):
                    flag = False
                    break
        if region != [] and flag:
            regions.append(region)
    vor.filtered_points = points_center
    vor.filtered_regions = regions
    return vor

def centroid_region(vertices):
    # Polygon's signed area
    A = 0
    # Centroid's x
    C_x = 0
    # Centroid's y
    C_y = 0
    for i in range(0, len(vertices) - 1):
        s = (vertices[i, 0] * vertices[i + 1, 1] - vertices[i + 1, 0] * vertices[i, 1])
        A = A + s
        C_x = C_x + (vertices[i, 0] + vertices[i + 1, 0]) * s
        C_y = C_y + (vertices[i, 1] + vertices[i + 1, 1]) * s
    A = 0.5 * A
    C_x = (1.0 / (6.0 * A)) * C_x
    C_y = (1.0 / (6.0 * A)) * C_y
    return np.array([[C_x, C_y]])

vor = voronoi(towers, bounding_box)

fig = pl.figure()
ax = fig.gca()
# Plot initial points
ax.plot(vor.filtered_points[:, 0], vor.filtered_points[:, 1], 'b.')
# Plot ridges points
for region in vor.filtered_regions:
    vertices = vor.vertices[region, :]
    ax.plot(vertices[:, 0], vertices[:, 1], 'go')
# Plot ridges
for region in vor.filtered_regions:
    vertices = vor.vertices[region + [region[0]], :]
    ax.plot(vertices[:, 0], vertices[:, 1], 'k-')
# Compute and plot centroids
centroids = []
for region in vor.filtered_regions:
    vertices = vor.vertices[region + [region[0]], :]
    centroid = centroid_region(vertices)
    centroids.append(list(centroid[0, :]))
    ax.plot(centroid[:, 0], centroid[:, 1], 'r.')

ax.set_xlim([-0.1, 1.1])
ax.set_ylim([-0.1, 1.1])
pl.savefig("bounded_voronoi.png")

sp.spatial.voronoi_plot_2d(vor)
pl.savefig("voronoi.png")

Si conservas alguna duda y capacidad de ascender nuestro artículo puedes escribir una anotación y con deseo lo estudiaremos.

¡Haz clic para puntuar esta entrada!
(Votos: 0 Promedio: 0)



Utiliza Nuestro Buscador

Deja una respuesta

Tu dirección de correo electrónico no será publicada. Los campos obligatorios están marcados con *