Calculer et dessiner des diagrammes Boronoi bornés (fermés) en Python

introduction

Les bibliothèques qui calculent et dessinent des diagrammes de Boronoi en Python incluent scipy.spatial.Voronoi et scipy.spatial.voronoi_plot_2d. Pour plus de détails à ce sujet, consultez les sites suivants.

• Documentation officielle de scipy.spatial.Voronoi
• Documentation officielle de scipy.spatial.voronoi_plot_2d
• scipy Voronoi

Cependant, puisque les diagrammes de Boronoi créés par ces bibliothèques sont calculés dans un espace illimité (non fermé), il n'est pas possible d'obtenir des informations sur le côté extérieur de Boronoi ou sur la région de Boronoi (polygone). Par conséquent, cette fois, je présenterai une méthode pour diviser une zone fermée en bore.

la mise en oeuvre

Pour la partie qui délimite la zone de Boronoi, je me suis référé à cet article. Le flux général est le suivant.

1. Ajoutez trois points mères factices pour relier toutes les régions de Boronoi.
2. Calculez le diagramme de Boronoi avec scipy.spatial.Voronoi.
3. Calculez la partie commune entre la zone à diviser et chaque zone de Boronoi avec «galbé».
4. Dessinez le polygone obtenu en 3.

import numpy as np
import matplotlib.pyplot as plt
from matplotlib.collections import PolyCollection
from scipy.spatial import Voronoi, voronoi_plot_2d
from shapely.geometry import Polygon

def bounded_voronoi(bnd, pnts):
    """
Une fonction qui calcule et dessine un diagramme Boronoi borné.
    """
    
    #Ajout de 3 points mères factices pour que la région de Boronoi de tous les points mères soit limitée
    gn_pnts = np.concatenate([pnts, np.array([[100, 100], [100, -100], [-100, 0]])])
    
    #Calcul du diagramme de boreoï
    vor = Voronoi(gn_pnts)
    
    #Polygone pour la zone à diviser
    bnd_poly = Polygon(bnd)
    
    #Liste pour stocker chaque zone de Boronoi
    vor_polys = []
    
    #Répétez pour les points mères non factices
    for i in range(len(gn_pnts) - 3):
        
        #Zone de Boronoi qui ne considère pas les espaces clos
        vor_poly = [vor.vertices[v] for v in vor.regions[vor.point_region[i]]]
        #Calculer la partie commune de la zone de Boronoi pour la zone à diviser
        i_cell = bnd_poly.intersection(Polygon(vor_poly))
        
        #Stocke les coordonnées des sommets de la région de Boronoi en tenant compte de l'espace fermé
        vor_polys.append(list(i_cell.exterior.coords[:-1]))
    
    
    #Dessiner un diagramme de Boronoi
    fig = plt.figure(figsize=(7, 6))
    ax = fig.add_subplot(111)
    
    #Point mère
    ax.scatter(pnts[:,0], pnts[:,1])
    
    #Région de Boronoi
    poly_vor = PolyCollection(vor_polys, edgecolor="black",
                              facecolors="None", linewidth = 1.0)
    ax.add_collection(poly_vor)
    
    xmin = np.min(bnd[:,0])
    xmax = np.max(bnd[:,0])
    ymin = np.min(bnd[:,1])
    ymax = np.max(bnd[:,1])
    
    ax.set_xlim(xmin-0.1, xmax+0.1)
    ax.set_ylim(ymin-0.1, ymax+0.1)
    ax.set_aspect('equal')
    
    plt.show()
    
    return vor_polys

Divisez le carré de l'unité

La division du carré de l'unité à l'aide de la fonction ci-dessus donne ce qui suit.

#Zone de Boronoi à diviser
bnd = np.array([[0, 0], [1, 0], [1, 1], [0, 1]])

#Nombre de points mères
n = 30
#Coordonnées du point mère
pnts = np.random.rand(n, 2)

#Calcul et dessin du diagramme de Boronoi
vor_polys = bounded_voronoi(bnd, pnts)

Diviser un polygone convexe général

Les polygones convexes généraux qui ne sont pas des carrés peuvent être divisés de la même manière.

from scipy.spatial import ConvexHull

def points_in_convex_polygon(bnd, n):
    """
Une fonction qui génère aléatoirement n points à l'intérieur d'un polygone convexe.
    """
    
    #Créer une matrice qui représente les limites d'une région
    bndhull = ConvexHull(bnd)
    bndTmp = bndhull.equations
    bndMat = np.matrix(bndTmp)
    Abnd = np.array(bndMat[:,0:2])
    bbnd = np.array(bndMat[:,2])
    
    #Rectangle entourant la zone
    xmin = np.min(bnd[:,0])
    xmax = np.max(bnd[:,0])
    ymin = np.min(bnd[:,1])
    ymax = np.max(bnd[:,1])
    
    #Pour la répétition
    i = 0
    pnts = []
    
    while i < n:
        
        #Générer des points
        pnt = np.random.rand(2)
        pnt[0] = xmin + (xmax - xmin) * pnt[0]
        pnt[1] = ymin + (ymax - ymin) * pnt[1]
        
        #Si le point est à l'intérieur d'un polygone convexe
        if (np.round(np.dot(Abnd,pnt.transpose()),7) <= np.round(-bbnd.transpose(),7)).all():

            pnts.append(pnt.tolist())
            i += 1

    return np.array(pnts)


#Zone de Boronoi à diviser
bnd = np.array([[0.1,0.4],[0.3,0.2],[0.8,0.3],[0.9,0.6],[0.7,0.7],[0.4,0.7],[0.2,0.6]])

#Nombre de points mères
n = 10
#Coordonnées du point mère
pnts = points_in_convex_polygon(bnd, n)

#Calcul et dessin du diagramme de Boronoi
vor_polys = bounded_voronoi(bnd, pnts)