Module thunderlab.voronoi

Analyze Voronoi diagrams based on scipy.spatial.

Classes

  • class Voronoi: Compute and analyse Voronoi diagrams.

Functions

def main()
Expand source code
def main():        
    import scipy.stats as st
    import matplotlib.pyplot as plt

    # generate random points:
    rs = np.random.randint(0xffffffff)
    #rs = 1718315706  # n=10
    #rs = 803645102   # n=10
    #rs = 2751318392  # double infinite ridges at vertex 0 (n=10)
    #rs = 4226093154 out-hull fac < 1.0
    
    print(f'random seed: {rs}')
    np.random.seed(rs)
    n = 20    # number of points
    ar = 1.0  # aspect ratio of area in which the points should be placed
    points = np.random.rand(int(n//ar), 2)
    points = points[points[:, 1] < ar, :]
    
    # calculate Voronoi diagram:
    vor = Voronoi(points)
        
    # what we get is:
    print(f'dimension: {vor.ndim}')
    print(f'number of points: {vor.npoints}')
    print('')
    print('distances of nearest neighbors:')
    print(vor.nearest_distances)
    print('for each point its nearest neighbor:')
    print(vor.nearest_points)
    print('for each point all Voronoi distances:')
    print(vor.neighbor_distances)
    print('for each point all its neighbors:')
    print(vor.neighbor_points)
    print('for each ridge distances of neighbors:')
    print(vor.ridge_distances)
    print('corresponding neighbors enclosing ridges:')
    print(vor.ridge_points)
    print('')
    print('length of corresponding ridges:')
    print(vor.ridge_lengths())
    print('area of corresponding triangles:')
    print(vor.ridge_areas())
    for mode in ['inside', 'finite_inside', 'full', 'finite']:
        print(f'Voronoi area of each point ({mode}):')
        print(vor.areas(mode))
    print('Comparison of sum of Voronoi areas (inside < finite_inside < full < finit):')
    a = [f'{np.sum(vor.areas(mode)):.2f}' for mode in ['inside', 'finite_inside', 'full', 'finite']]
    print(' < '.join(a))
    print('Type of Voronoi area of each point:')
    print(vor.point_types())

    # plot Voronoi diagram:
    fig, ax = plt.subplots()
    ax.set_title('Voronoi')
    vor.fill_regions(ax, inside=True, colors=['red', 'green', 'blue', 'orange', 'cyan', 'magenta'], alpha=1.0, zorder=0)
    vor.fill_regions(ax, inside=False, colors=['red', 'green', 'blue', 'orange', 'cyan', 'magenta'], alpha=0.4, zorder=0)
    vor.fill_infinite_regions(ax, colors=['red', 'green', 'blue', 'orange', 'cyan', 'magenta'], alpha=0.1, zorder=0)
    vor.plot_distances(ax, color='red')
    #vor.plot_ridges(ax, c='k', lw=2)
    #vor.plot_infinite_ridges(ax, c='k', lw=2, linestyle='dashed')
    vor.plot_points(ax, text='p%d', c='c', s=100, zorder=10)
    vor.plot_vertices(ax, text='v%d', c='r', s=60)
    ex = 0.3
    ax.set_xlim(vor.min_bound[0]-ex, vor.max_bound[0]+ex)
    ax.set_ylim(vor.min_bound[1]-ex, vor.max_bound[1]+ex)
    ax.set_aspect('equal')

    # Convex hull:
    print('Convex hull:')
    print(f'Area of convex hull: {vor.hull_area():g}')
    print(f'Area of outer hull: {vor.outer_hull_area():g}')
    new_points = vor.random_points(poisson=True, mode='outer')

    # plot convex hull:
    fig, ax = plt.subplots()
    vor.fill_outer_hull(ax, color='black', alpha=0.1)
    vor.plot_outer_hull(ax, color='m', lw=2, label='outer hull')
    vor.fill_hull(ax, color='black', alpha=0.2)
    vor.plot_hull(ax, color='r', lw=2, label='convex hull')
    vor.plot_hull_center(ax, color='r', ms=16)
    vor.plot_hull_center(ax, color='k', ms=4)
    vor.plot_center(ax, color='c', ms=16)
    vor.plot_center(ax, color='k', ms=4)
    vor.plot_points(ax, text='p%d', c='c', s=100, zorder=10,
                    label='input points')
    ax.plot(new_points[:, 0], new_points[:, 1], 'ob', ms=10,
             label='random points')
    ax.set_xlim(vor.min_bound[0]-ex, vor.max_bound[0]+ex)
    ax.set_ylim(vor.min_bound[1]-ex, vor.max_bound[1]+ex)
    ax.set_aspect('equal')
    ax.legend()
        
    # bootstrap:
    def bootstrapped_nearest_distances(vor, n, poisson, mode, ax1, ax2, ax3, bins):
        ps = 'fixed'
        if poisson:
            ps = 'poisson'
        print(f'bootstrap {mode} {ps} ...')
        hist_distances = np.zeros((n, len(bins)-1))
        mean_distances = []
        cvs = []
        for j in range(n):
            points = vor.random_points(poisson=poisson, mode=mode)
            if len(points) < 4:
                continue
            try:
                bvor = Voronoi(points)
            except:
                continue
            hist_distances[j,:] = np.histogram(bvor.nearest_distances,
                                               bins=bins, density=True)[0]
            mean_distances.append(np.mean(bvor.nearest_distances))
            cvs.append(np.std(bvor.nearest_distances)/np.mean(bvor.nearest_distances))
        ax1.set_title(f'bootstrap {mode} {ps} ...')
        ax1.hist(mean_distances, 30, density=True, label='bootstrap')
        de = np.mean(vor.nearest_distances)
        sem = 2.1*0.26136*np.sqrt(vor.outer_hull_area())/vor.npoints # note factor 2.1!
        h = 0.4/sem
        x = np.linspace(de-4.0*sem, de+4.0*sem, 200)
        p = st.norm.pdf(x, loc=de, scale=sem)
        ax1.plot(x, p, 'r', lw=2, label='Gaussian')
        ax1.plot([de, de], [0.0, h], 'r', lw=4, label='observed a.n.n.')
        bmd = np.mean(mean_distances)
        ax1.plot([bmd, bmd], [0.0, h], 'g', lw=4, label='bootstrapped a.n.n.')
        ax1.set_xlabel('mean nearest-neighbor distance')
        ax1.legend()
        ax2.set_title(f'bootstrap {mode} {ps} ...')
        ax2.hist(cvs, 30, density=True, label='bootstrap')
        # observed cv:
        cvo = np.std(vor.nearest_distances)/np.mean(vor.nearest_distances)
        # significance from bootstrap:
        pb = 0.01*st.percentileofscore(cvs, cvo)
        if pb > 0.5:
            pb = 1.0-pb
        pb *= 2.0
        # significance from z-score:
        secv = 2.0*np.sqrt(1.0/np.pi-0.25)/np.sqrt(vor.npoints)
        zcv = (cvo - 2.0*0.26136)/secv
        pz = 2.0*(1.0 - st.norm.cdf(np.abs(zcv)))
        h = 0.4/secv
        ax2.plot([cvo, cvo], [0.0, h], 'r', lw=4, label=rf'observed CV $\alpha=${100.0*pz:.0f}%')
        cvb = np.mean(cvs)
        x = np.linspace(cvb-4.0*secv, cvb+4.0*secv, 200)
        p = st.norm.pdf(x, loc=cvb, scale=secv)
        ax2.plot(x, p, 'g', lw=2, label='Gaussian')
        ax2.plot([cvb, cvb], [0.0, h], 'g', lw=4, label=rf'bootstrapped CV $\alpha=${100.0*pb:.0f}%')
        ax2.set_xlabel('CV')
        ax2.legend()
        # nearest-neighbor distance histogram:
        ax3.set_title(f'bootstrap {mode} {ps} ...')
        ax3.fill_between(bins[:-1], *np.percentile(hist_distances, [2.5, 97.5], axis=0), facecolor='blue', alpha=0.5)
        ax3.plot(bins[:-1], np.median(hist_distances, axis=0), 'b', lw=2)
        hd = np.histogram(vor.nearest_distances, bins=bins, density=True)[0]
        ax3.plot(bins[:-1], hd, 'r', lw=2)
        ax3.set_xlim(0.0, 0.5)
        ax3.set_xlabel('nearest-neighbor distance')
        
    print(f'Mean distance: {np.mean(vor.nearest_distances):g}')
    print(f'CV: {np.std(vor.nearest_distances)/np.mean(vor.nearest_distances):g}')
    bins = np.linspace(0.0, 1.0, 30)
    nb = 300
    fig1, axs1 = plt.subplots(2, 2, sharex=True)
    fig2, axs2 = plt.subplots(2, 2, sharex=True)
    fig3, axs3 = plt.subplots(2, 2, sharex=True)
    for k, poisson in enumerate([False, True]):
        for j, mode in enumerate(['hull', 'outer']):
            bootstrapped_nearest_distances(vor, nb, poisson, mode,
                                           axs1[k,j], axs2[k,j], axs3[k,j],
                                           bins)
    for f in [fig1, fig2, fig3]:
        f.tight_layout()
    print('... done.')
    plt.tight_layout()
    plt.show()

Classes

class Voronoi (points, radius=None, qhull_options=None)
Expand source code
class Voronoi(object):
    
""" Parameters ---------- points: ndarray of floats, shape (npoints, 2) List of point coordinates. radius: float or None Radius for computing far points of infinite ridges. If None twice the maximum extent of the input points is used. qhull_options: string or None Options to be passed on to QHull. From the manual: - Qbb - scale last coordinate to [0,m] for Delaunay triangulations - Qc - keep coplanar points with nearest facet - Qx - exact pre-merges (skips coplanar and anglomaniacs facets) - Qz - add point-at-infinity to Delaunay triangulation - QJn - randomly joggle input in range [-n,n] - Qs - search all points for the initial simplex - Qz - add point-at-infinity to Voronoi diagram - QGn - Voronoi vertices if visible from point n, -n if not - QVn - Voronoi vertices for input point n, -n if not Default is: "Qbb Qc Qz Qx" for ndim > 4 and "Qbb Qc Qz" otherwise. Input points ------------ The points from which the Voronoi diagram is constructed. - `ndim`: int The dimension of the input points, i.e. number of coordinates. - `npoints`: int The number of input points. - `points`: 2-d ndarray List of input point coordinates. - `center`: ndarray of floats Center of mass of the input points (i.e. mean coordinates). - `min_bound`: ndarray of floats Lower left corner of the bounding-box of the input points. - `max_bound`: ndarray of floats Upper right corner of the bounding-box of the input points. Distances between input points ------------------------------ - `ridge_points`: 2-d ndarray of ints List of pairs of indices to `points` enclosing a Voronoi ridge. - `ridge_distances`: ndarray of floats For each ridge in `ridge_points` the distance between the enclosing points. - `neighbor_points`: list of ndarrays of ints For each point in `points` a list of indices of the Voronoi-cell's neighboring points. - `neighbor_distances`: list of ndarrays of floats For each point in `points` a list of distances to the Voronoi-cell's neighboring points, matching `neighbor_points`. - `nearest_points`: list of ints For each point in `points` the index of its nearest neighbor. - `nearest_distances`: ndarray of floats For each point in `points` the distance to its nearest neighbor. - `mean_nearest_distance`: float Unbiased estimate of the mean nearest-neighbor distance. This is the average nearest-neighbor distance corrected by one s.e.m, since points close to the border of the region have too large nearest-neighbor distances. Voronoi diagram --------------- - `vertices`: 2-d ndarray of floats List of vertex coordinates enclosing the Voronoi regions. - `regions`: list of list of ints List of lists of indices to vertices in `vertices` making up each Voronoi region. - `ridge_vertices`: list of list of ints List of lists of indices to `vertices` making up a Voronoi ridge. The pairs of `points` to poth sides of the ridge are listed in `ridge_points`. - `ridge_lengths()`: Length of Voronoi ridges between nearest neighbors. - `areas()`: The areas of the Voronoi regions for each input point in `points`. - `point_types()`: The type of Voronoi area (infinite, finite, inside) for each input point. - `inside_vertices`: list of ints Indices of `vertices` that are inside the convex hull of the input points. - `infinite_vertices`: list of ndarrays of floats For each ridge in `ridge_vertices` the coordinates of the far-point, if it is an infinite ridge. - `infinite_regions`: list of list of ints List of Voronoi regions with infinite ridges. If positive, the indices are indices to `vertices`. If negative they are indices into `infinite_vertices` (`-index-1`). Convex hull ----------- - `hull_points`: list of ints List of indices of the points in `points` making up the convex hull. - `hull_center`: ndarray of floats Center of mass of the points making up the convex hull. - `inside_vertices`: ndarray of boolean Indicates for each vertex in `vertices` whether it is inside the convex hull. - `in_hull()`: Test if points are within the convex hull of the input points. - `hull_area()`: The area contained in the convex hull. Outer hull ---------- The outer hull is constructed from the convex hull by expanding each point of the convex hull away from the center of mass of the convex hull such that randomly placing the same number of points in this area results on average in the same mean nearest-neighbor distance. This enlarged hull is needed for bootstrapping. Note that in some cases the outer hull can be *smaller* than the convex hull. - `outer_hull_points`: 2-d ndarray of floats Array of coordinates of the points of the outer hull. - `outer_min_bound`: ndarray of floats Lower left corner of the bounding-box of the outer hull. - `outer_max_bound`: ndarray of floats Upper right corner of the bounding-box of the outer hull. - `in_outer_hull()`: Test if points are within the outer hull. - `outer_hull_area()`: The area contained in the outer hull. Bootstrap Voronoi diagrams -------------------------- - `random_points()`: Generate random points within the area defined by the input points. Plotting the Voronoi diagram ---------------------------- - `plot_points()`: Plot and optionally annotate the input points of the Voronoi diagram. - `plot_center()`: Plot the center of mass of the input points. - `plot_vertices()`: Plot and optionally annotate the vertices of the Voronoi diagram. - `plot_distances()`: Plot lines connecting the neighbors in the Voronoi diagram. - `plot_ridges()`: Plot the finite ridges of the Voronoi diagram. - `plot_infinite_ridges(): Plot the infinite ridges of the Voronoi diagram. - `fill_regions()`: Fill each finite region of the Voronoi diagram with a color. - `fill_infinite_regions()`: Fill each infinite region of the Voronoi diagram with a color. Plotting the convex hull ------------------------ - `plot_hull()`: Plot the convex hull line containing the input points. - `fill_hull()`: Fill the convex hull. - `plot_hull_center()`: Plot the center of mass of the convex hull. - `plot_outer_hull()`: Plot the outer hull edge. - `fill_outer_hull()`: Fill the outer hull. Usage ----- Generate 20 random points in 2-D: ``` import numpy as np points = np.random.rand(20, 2) ``` Calculate the Voronoi diagram: ``` from thunderlab.voronoi import Voronoi vor = Voronoi(points) ``` Retrieve nearest-neighbor distances and compute Voronoi areas: ``` dists = vor.nearest_distances areas = vor.areas() ``` Generate random points drawn from a Poisson process with the same intensity and mean nearest neighbor distance as the data: ``` new_points = vor.random_points(poisson=True, mode='outer') ``` Plot Voronoi regions, distances, and input points: ``` import matplotlib.pyplot as plt vor.fill_regions(colors=['red', 'green', 'blue', 'orange', 'cyan'], alpha=0.3) vor.plot_distances(color='red') vor.plot_points(text='p%d', c='c', s=100) ``` """
def __init__(self, points, radius=None, qhull_options=None): self.vor = sp.Voronoi(points, furthest_site=False, incremental=False, qhull_options=qhull_options) if self.vor.ndim != 2: raise ValueError("Only 2D input points are supported.") self.hull = sp.Delaunay(points, furthest_site=False, incremental=False, qhull_options=qhull_options) self._compute_distances() self._compute_infinite_vertices() self._compute_hull(qhull_options) self.ndim = self.vor.ndim self.npoints = self.vor.npoints self.points = self.vor.points self.vertices = self.vor.vertices self.regions = self.vor.regions self.ridge_points = self.vor.ridge_points self.ridge_vertices = self.vor.ridge_vertices self.min_bound = self.vor.min_bound self.max_bound = self.vor.max_bound self.center = np.mean(self.points, axis=0) def _compute_distances(self): """Compute distances between points. """ # For each ridge the distance of the points enclosing the ridge: p1 = self.vor.points[self.vor.ridge_points[:,0]] p2 = self.vor.points[self.vor.ridge_points[:,1]] self.ridge_distances = sp.minkowski_distance(p1, p2) # For each point all its Voronoi distances: self.neighbor_points = [[] for k in range(len(self.vor.points))] self.neighbor_distances = [[] for k in range(len(self.vor.points))] for dist, points in zip(self.ridge_distances, self.vor.ridge_points): self.neighbor_points[points[0]].append(points[1]) self.neighbor_points[points[1]].append(points[0]) self.neighbor_distances[points[0]].append(dist) self.neighbor_distances[points[1]].append(dist) for k in range(len(self.neighbor_points)): inx = np.argsort(self.neighbor_distances[k]) self.neighbor_points[k] = np.array(self.neighbor_points[k])[inx] self.neighbor_distances[k] = np.array(self.neighbor_distances[k])[inx] # For each point the distance to its nearest neighbor: self.nearest_points = [] self.nearest_distances = np.zeros(len(self.neighbor_distances)) for k in range(len(self.neighbor_distances)): self.nearest_points.append(self.neighbor_points[k][0]) self.nearest_distances[k] = self.neighbor_distances[k][0] self.mean_nearest_distance = np.mean(self.nearest_distances) self.mean_nearest_distance *= 1.0-2.0*np.sqrt(1.0/np.pi-0.25)/np.sqrt(self.vor.npoints) def _compute_infinite_vertices(self, radius=None):
"""Compute far points of infinite ridges. Parameters ---------- radius: float or None Radius for computing far points of infinite ridges. If None twice the maximum extent of the input points is used. Note ---- Code inspired by http://stackoverflow.com/questions/20515554/colorize-voronoi-diagram """
# For each ridge, compute far point: center = self.vor.points.mean(axis=0) if radius is None: radius = np.max(np.ptp(2.0*self.vor.points, axis=0)) self.infinite_vertices = [] for points, vertices in zip(self.vor.ridge_points, self.vor.ridge_vertices): vertices = np.asarray(vertices) if np.all(vertices >= 0): self.infinite_vertices.append(np.array([])) else: i = vertices[vertices >= 0][0] # finite end Voronoi vertex t = self.vor.points[points[1]] - self.vor.points[points[0]] # tangent t /= np.linalg.norm(t) n = np.array([-t[1], t[0]]) # normal midpoint = self.vor.points[points].mean(axis=0) direction = np.sign(np.dot(midpoint - center, n)) * n far_point = self.vor.vertices[i] + direction * radius self.infinite_vertices.append(far_point) # Assemble list of infinite regions: # Indices to self.infinite_vertices are negative minus one. self.infinite_regions = [] for rvertices in self.vor.regions: if -1 in rvertices: new_rvertices = [] prev_vertex = rvertices[-1] # find index of data point enclosed by the region: ridge_points = [] for p, v in zip(self.vor.ridge_points, self.vor.ridge_vertices): if not -1 in v and set(v) <= set(rvertices): ridge_points.extend(p) region_point = None for rp in ridge_points: if ridge_points.count(rp) > 1: region_point = rp break # fill in far points for each region: for v_inx, v in enumerate(rvertices): if v >= 0: new_rvertices.append(v) else: for v1_inx, (points, vertices) in enumerate(zip(self.vor.ridge_points, self.vor.ridge_vertices)): if prev_vertex in vertices and -1 in vertices and \ (region_point is None or region_point in points): new_rvertices.append(-v1_inx-1) break next_vertex = rvertices[0] if v_inx+1 < len(rvertices): next_vertex = rvertices[v_inx+1] for v2_inx, (points, vertices) in enumerate(zip(self.vor.ridge_points, self.vor.ridge_vertices)): if next_vertex in vertices and -1 in vertices and \ (region_point is None or region_point in points) and \ new_rvertices[-1] != -v2_inx-1: new_rvertices.append(-v2_inx-1) break prev_vertex = v self.infinite_regions.append(new_rvertices) def _flatten_simplices(self, simplices):
"""Transforms list of simplex indices to list of vertex indices. In particular, transforms the Delaunay.convex_hull to a list of points of the hull that can then be easily plotted. Parameters ---------- simplices: 2-D array of ints List of pairs of indices of points forming each ridge of a polygon. Returns ------- indices: list of ints Indices of vertices of the polygon. """
if len(simplices) == 0: return [] indices = list(simplices[0]) simplices = np.delete(simplices, 0, 0) while len(simplices) > 0: for i, s in enumerate(simplices): if indices[-1] in s: if s[0] == indices[-1]: indices.append(s[1]) else: indices.append(s[0]) simplices = np.delete(simplices, i, 0) break return indices def _compute_hull(self, qhull_options): """Compute properties of the convex hull and set up the outer hull. """ self.inside_vertices = self.in_hull(self.vor.vertices) self.hull_points = self._flatten_simplices(self.hull.convex_hull) self.hull_center = np.mean(self.hull.points[self.hull_points], axis=0) # enlarge hull: # estimate area needed for a poisson process with same mean distance: new_area = 4.0 * self.mean_nearest_distance**2 * self.vor.npoints fac = np.sqrt(new_area/self.hull_area()) self.outer_hull_points = np.zeros((len(self.hull_points), self.vor.ndim)) for k, point in enumerate(self.hull.points[self.hull_points]): point -= self.hull_center point *= fac point += self.hull_center self.outer_hull_points[k] = point # compute outer hull: self.outer_hull = sp.Delaunay(self.outer_hull_points, furthest_site=False, incremental=False, qhull_options=qhull_options) self.outer_min_bound = np.min(self.outer_hull_points, axis=0) self.outer_max_bound = np.max(self.outer_hull_points, axis=0) def in_hull(self, p):
"""Test if points `p` are within convex hull of the input points. Parameters ---------- p: 2-d ndarray Array of points to be tested. Returns ------- inside: ndarray of booleans For each point in `p` whether it is inside the hull. """
inside = self.hull.find_simplex(p) >= 0 return inside def in_outer_hull(self, p):
"""Test if points `p` are within the outer hull. Parameters ---------- p: 2-d ndarray Array of points to be tested. Returns ------- inside: ndarray of booleans For each point in `p` whether it is inside the outer hull. """
inside = self.outer_hull.find_simplex(p) >= 0 return inside def point_types(self):
"""The type of the Voronoi regions for each input point. Returns ------- points: ndarray of ints Indicates the type of Voronoi region associated with each point in `points`: - 2: finite region with all vertices inside hull, - 1: finite region, - 0: infinite region. """
points = np.zeros(len(self.vor.points), dtype=int) + 2 for i in range(len(self.vor.points)): for j, rp in enumerate(self.vor.ridge_points): if i in rp: if not np.all(self.inside_vertices[self.vor.ridge_vertices[j]]) and\ points[i] > 0: points[i] = 1 if not np.all(np.array(self.vor.ridge_vertices[j]) >= 0) and\ points[i] > -1: points[i] = 0 return points def ridge_lengths(self):
"""Length of Voronoi ridges between nearest neighbors. May be used, for example, as a weigth for `ridge_distances`. Returns ------- distances: ndarray of floats The length of each ridge in `ridge_vertices`. np.inf if one vertex is at infinity. """
ridges = np.zeros(len(self.vor.ridge_vertices)) for k, p in enumerate(self.vor.ridge_vertices): if np.all(np.array(p)>=0): p1 = self.vor.vertices[p[0]] p2 = self.vor.vertices[p[1]] ridges[k] = sp.minkowski_distance(p1, p2) else: ridges[k] = np.inf return ridges def ridge_areas(self):
"""For each ridge the triangular area of the Voronoi region spanned by the center point and the ridge. Returns ------- areas: ndarray of floats For each ridge in `ridge_points` or `ridge_vertices` its corresponding triangular area. np.inf for infinite ridges. """
ridges = self.ridge_lengths() heights = 0.5*self.ridge_distances # area of a triangle: areas = 0.5*ridges*heights return areas def areas(self, mode='finite'):
"""The areas of the Voronoi regions for each input point. Parameters ---------- mode: string - 'inside': Calculate area of finite Voronoi regions whose vertices are all inside the hull, set all other to zero. - 'finite_inside': Calculate area of all Voronoi regions. Consider only areas of finite ridges whose vertices are all inside the hull. - 'full': Calculate area of finite Voronoi regions only, set all other to zero. - 'finite': Calculate area of all Voronoi regions. From infinite regions only areas contributed by finite ridges are considered. Returns ------- areas: array of floats For each point in `points` its corresponding area. """
ridge_areas = self.ridge_areas() areas = np.zeros(len(self.vor.points)) if mode == 'inside': for i in range(len(self.vor.points)): a = 0.0 for j, rp in enumerate(self.vor.ridge_points): if i in rp: if ridge_areas[j] != np.inf and \ np.all(self.inside_vertices[self.vor.ridge_vertices[j]]): a += ridge_areas[j] else: a = 0.0 break areas[i] = a elif mode == 'full': for i in range(len(self.vor.points)): a = 0.0 for j, rp in enumerate(self.vor.ridge_points): if i in rp: if ridge_areas[j] != np.inf: a += ridge_areas[j] else: a = 0.0 break areas[i] = a elif mode == 'finite_inside': for i in range(len(self.vor.points)): a = 0.0 for j, rp in enumerate(self.vor.ridge_points): if i in rp and ridge_areas[j] != np.inf and \ np.all(self.inside_vertices[self.vor.ridge_vertices[j]]): a += ridge_areas[j] areas[i] = a elif mode == 'finite': for i in range(len(self.vor.points)): a = 0.0 for j, rp in enumerate(self.vor.ridge_points): if i in rp and ridge_areas[j] != np.inf: a += ridge_areas[j] areas[i] = a else: print('') print(f'Voronoi.areas(): unknown value "{mode}" for the mode parameter:') print('Use one of the following values:') print(' inside: Finite Voronoi regions whose vertices are all inside the hull.') print(' finite_inside: Use all areas corresponding to finite ridges whose vertices are all inside the hull.') print(' full: Finite Voronoi regions only.') print(' finite: Use all areas corresponding to finite ridges.') print('') return areas def hull_area(self): """The area of the convex hull of the input points. Returns ------- area: float The area of the convex hull. """ # two sides of the simplex triangles: ab = self.hull.points[self.hull.simplices[:,0],:] - self.hull.points[self.hull.simplices[:,1],:] cb = self.hull.points[self.hull.simplices[:,2],:] - self.hull.points[self.hull.simplices[:,1],:] # area of each simplex is half of the absolute value of the cross product: area = 0.5*np.sum(np.abs(np.cross(ab, cb))) return area def outer_hull_area(self): """The area of the outer hull. Returns ------- area: float The area of the outer hull. """ # two sides of the simplex triangles: ab = self.outer_hull.points[self.outer_hull.simplices[:,0],:] - self.outer_hull.points[self.outer_hull.simplices[:,1],:] cb = self.outer_hull.points[self.outer_hull.simplices[:,2],:] - self.outer_hull.points[self.outer_hull.simplices[:,1],:] # area of each simplex is half of the absolute value of the cross product: area = 0.5*np.sum(np.abs(np.cross(ab, cb))) return area def random_points(self, n=None, poisson=True, mode='outer'):
"""Generate random points. Parameters ---------- n: int or None Number of random points to be generated. If None `n` is set to the number of input points. poisson: boolean If `True` then draw the number of points from a Poisson distribution with mean number of points given by `n`. mode: string - 'bbox': place points randomly in rectangular bounding box of the Voronoi diagram. - 'hull': place points randomly within convex hull of input data. - 'outer' place points randomly within outer hull. The mean nearest-neighbor distance between the generated points will be close to the observed one. Returns ------- points: 2-D array List of randomly generated points. """
# number of points: if n is None: n = self.npoints nn = n if poisson: nn = np.random.poisson(n) m = nn//2 if m < 5: m = 5 # get bounding box: if mode == 'outer': min_bound = self.outer_min_bound max_bound = self.outer_max_bound else: min_bound = self.min_bound max_bound = self.max_bound delta = np.max(max_bound - min_bound) points = np.zeros((0, self.ndim)) while len(points) < nn: # random points within bounding box: newpoints = np.random.rand(m, self.ndim) newpoints *= delta newpoints += min_bound if mode == 'outer': # only take the ones within outer hull: inside = self.in_outer_hull(newpoints) points = np.vstack((points, newpoints[inside])) elif mode == 'hull': # only take the ones within hull: inside = self.in_hull(newpoints) points = np.vstack((points, newpoints[inside])) elif mode == 'bbox': points = np.vstack((points, newpoints[np.all(newpoints<max_bound, axis=1),:])) else: print('') print(f'Voronoi.random_points(): unknown value "{mode}" for the mode parameter:') print('Use one of the following values:') print(' bbox: Place points within rectangular bounding box.') print(' hull: Place points inside the convex hull.') print(' outer: Place points inside the outer hull.') print('') return return points[:nn] def plot_points(self, ax, text=None, text_offs=(0, 0.05), text_align='center', **kwargs):
"""Plot and optionally annotate the input points of the Voronoi diagram. Parameters ---------- ax: matplotlib.Axes The axes to be used for plotting. text: string or None If not None the string that is placed at each point. A '%d' is replaced by the index of the point. text_offset: tuple of numbers The offset of the point labels. text_align: string The horizontal alignment of the point labels. **kwargs: dict Key-word arguments that are passed on to the matplotlib.scatter() function. """
ax.scatter(self.points[:,0], self.points[:,1], **kwargs) if text is not None: for i, p in enumerate(self.points): s = text if '%' in text: s = text % i ax.text(p[0]+text_offs[0], p[1]+text_offs[1], s, ha=text_align) def plot_center(self, ax, **kwargs):
"""Plot the center of mass of the input points. Parameters ---------- ax: matplotlib.Axes The axes to be used for plotting. **kwargs: dict Key-word arguments that are passed on to the `matplotlib.plot()` function. """
ax.plot(self.center[0], self.center[1], 'o', **kwargs) def plot_vertices(self, ax, text=None, text_offs=(0, 0.05), text_align='center', **kwargs):
"""Plot and optionally annotate the vertices of the Voronoi diagram. Parameters ---------- ax: matplotlib.Axes The axes to be used for plotting. text: string or None If not None the string that is placed at each vertex. A '%d' is replaced by the index of the vertex. text_offset: tuple of numbers The offset of the vertex labels. text_align: string The horizontal alignment of the vertex labels. **kwargs: dict Key-word arguments that are passed on to the matplotlib.scatter() function. """
ax.scatter(self.vertices[:,0], self.vertices[:,1], **kwargs) if text is not None: for i, p in enumerate(self.vertices): s = text if '%' in text: s = text % i ax.text(p[0]+text_offs[0], p[1]+text_offs[1], s, ha=text_align) def plot_distances(self, ax, **kwargs):
"""Plot lines connecting the nearest neighbors in the Voronoi diagram. Parameters ---------- ax: matplotlib.Axes The axes to be used for plotting. **kwargs: dict Key-word arguments that are passed on to the matplotlib.plot() function. """
for i, p in enumerate(self.ridge_points): ax.plot(self.points[p, 0], self.points[p, 1], **kwargs) def plot_ridges(self, ax, **kwargs):
"""Plot the finite ridges of the Voronoi diagram. Parameters ---------- ax: matplotlib.Axes or None The axes to be used for plotting. **kwargs: dict Key-word arguments that are passed on to the matplotlib.plot() function. """
for i, p in enumerate(self.ridge_vertices): if np.all(np.array(p)>=0): ax.plot(self.vertices[p, 0], self.vertices[p, 1], **kwargs) def plot_infinite_ridges(self, ax, **kwargs):
"""Plot the infinite ridges of the Voronoi diagram. Parameters ---------- ax: matplotlib.Axes The axes to be used for plotting. **kwargs: dict Key-word arguments that are passed on to the matplotlib.plot() function. """
for far_point, vertices in zip(self.infinite_vertices, self.vor.ridge_vertices): vertices = np.asarray(vertices) if not np.all(vertices >= 0): i = vertices[vertices >= 0][0] # finite end Voronoi vertex ax.plot([self.vor.vertices[i][0], far_point[0]], [self.vor.vertices[i][1], far_point[1]], **kwargs) def fill_regions(self, ax, inside=None, colors=None, **kwargs):
"""Fill each finite region of the Voronoi diagram with a color. Parameters ---------- ax: matplotlib.Axes The axes to be used for plotting. inside: boolean or None True: plot only finite regions with all vertices inside the hull False: plot only finite regions with at least one vertex outside the hull None: plot all finite regions colors: list of colors or None If not None then these colors are used in turn to fill the regions. **kwargs: dict Key-word arguments that are passed on to the matplotlib.fill() function. """
c = 0 for region in self.regions: if not -1 in region: polygon = self.vertices[region] if len(polygon) > 0: inside_hull = self.in_hull(polygon) if inside is None or (inside and all(inside_hull)) or (not inside and any(inside_hull)): c += 1 if colors is None: ax.fill(polygon[:, 0], polygon[:, 1], lw=0, **kwargs) else: ax.fill(polygon[:, 0], polygon[:, 1], color=colors[c % len(colors)], lw=0, **kwargs) def fill_infinite_regions(self, ax, colors=None, **kwargs):
"""Fill each infinite region of the Voronoi diagram with a color. Parameters ---------- ax: matplotlib.Axes The axes to be used for plotting. colors: list of colors or None If not None then these colors are used in turn to fill the regions. **kwargs: dict Key-word arguments that are passed on to the matplotlib.fill() function. """
c = 0 for region in self.infinite_regions: polygon = [] for p in region: if p >= 0: polygon.append(self.vertices[p]) else: polygon.append(self.infinite_vertices[-p-1]) if len(polygon) > 0: polygon = np.asarray(polygon) c += 1 if colors is None: ax.fill(polygon[:, 0], polygon[:, 1], lw=0, **kwargs) else: ax.fill(polygon[:, 0], polygon[:, 1], color=colors[c % len(colors)], lw=0, **kwargs) def plot_hull(self, ax, **kwargs):
"""Plot the hull line containing the input points. Parameters ---------- ax: matplotlib.Axes The axes to be used for plotting. **kwargs: dict Key-word arguments that are passed on to the matplotlib.plot() function. """
ax.plot(self.hull.points[self.hull_points, 0], self.hull.points[self.hull_points, 1], **kwargs) def fill_hull(self, ax, **kwargs):
"""Fill the hull containing the input points with a color. Parameters ---------- ax: matplotlib.Axes The axes to be used for plotting. **kwargs: dict Key-word arguments that are passed on to the matplotlib.fill() function. """
ax.fill(self.hull.points[self.hull_points, 0], self.hull.points[self.hull_points, 1], lw=0, **kwargs) def plot_hull_center(self, ax, **kwargs):
"""Plot the center of mass of the convex hull of the input points. Parameters ---------- ax: matplotlib.Axes The axes to be used for plotting. **kwargs: dict Key-word arguments that are passed on to the matplotlib.plot() function. """
ax.plot(self.hull_center[0], self.hull_center[1], 'o', **kwargs) def plot_outer_hull(self, ax, **kwargs):
"""Plot the hull line containing the input points and the vertices of the Voronoi diagram. Parameters ---------- ax: matplotlib.Axes or None The axes to be used for plotting. **kwargs: dict Key-word arguments that are passed on to the matplotlib.plot() function. """
ax.plot(self.outer_hull_points[:, 0], self.outer_hull_points[:, 1], **kwargs) def fill_outer_hull(self, ax, **kwargs):
"""Fill the hull containing the input points and the vertices of the Voronoi diagram. Parameters ---------- ax: matplotlib.Axes The axes to be used for plotting. **kwargs: dict Key-word arguments that are passed on to the matplotlib.fill() function. """
ax.fill(self.outer_hull_points[:, 0], self.outer_hull_points[:, 1], lw=0, **kwargs)

Parameters

points : ndarray of floats, shape (npoints, 2)
List of point coordinates.
radius : float or None
Radius for computing far points of infinite ridges. If None twice the maximum extent of the input points is used.
qhull_options : string or None

Options to be passed on to QHull. From the manual:

  • Qbb - scale last coordinate to [0,m] for Delaunay triangulations
  • Qc - keep coplanar points with nearest facet
  • Qx - exact pre-merges (skips coplanar and anglomaniacs facets)
  • Qz - add point-at-infinity to Delaunay triangulation
  • QJn - randomly joggle input in range [-n,n]
  • Qs - search all points for the initial simplex
  • Qz - add point-at-infinity to Voronoi diagram
  • QGn - Voronoi vertices if visible from point n, -n if not
  • QVn - Voronoi vertices for input point n, -n if not

Default is: "Qbb Qc Qz Qx" for ndim > 4 and "Qbb Qc Qz" otherwise.

Input Points

The points from which the Voronoi diagram is constructed.

  • ndim: int
    The dimension of the input points, i.e. number of coordinates.
  • npoints: int
    The number of input points.
  • points: 2-d ndarray
    List of input point coordinates.
  • center: ndarray of floats
    Center of mass of the input points (i.e. mean coordinates).
  • min_bound: ndarray of floats
    Lower left corner of the bounding-box of the input points.
  • max_bound: ndarray of floats
    Upper right corner of the bounding-box of the input points.

Distances Between Input Points

  • ridge_points: 2-d ndarray of ints
    List of pairs of indices to points enclosing a Voronoi ridge.
  • ridge_distances: ndarray of floats
    For each ridge in ridge_points the distance between the enclosing points.
  • neighbor_points: list of ndarrays of ints
    For each point in points a list of indices of the Voronoi-cell's neighboring points.
  • neighbor_distances: list of ndarrays of floats
    For each point in points a list of distances to the Voronoi-cell's neighboring points, matching neighbor_points.
  • nearest_points: list of ints
    For each point in points the index of its nearest neighbor.
  • nearest_distances: ndarray of floats
    For each point in points the distance to its nearest neighbor.
  • mean_nearest_distance: float
    Unbiased estimate of the mean nearest-neighbor distance. This is the average nearest-neighbor distance corrected by one s.e.m, since points close to the border of the region have too large nearest-neighbor distances.

Voronoi Diagram

  • vertices: 2-d ndarray of floats
    List of vertex coordinates enclosing the Voronoi regions.
  • regions: list of list of ints
    List of lists of indices to vertices in vertices making up each Voronoi region.
  • ridge_vertices: list of list of ints
    List of lists of indices to vertices making up a Voronoi ridge. The pairs of points to poth sides of the ridge are listed in ridge_points.
  • ridge_lengths(): Length of Voronoi ridges between nearest neighbors.
  • areas(): The areas of the Voronoi regions for each input point in points.
  • point_types(): The type of Voronoi area (infinite, finite, inside) for each input point.
  • inside_vertices: list of ints
    Indices of vertices that are inside the convex hull of the input points.
  • infinite_vertices: list of ndarrays of floats
    For each ridge in ridge_vertices the coordinates of the far-point, if it is an infinite ridge.
  • infinite_regions: list of list of ints
    List of Voronoi regions with infinite ridges. If positive, the indices are indices to vertices. If negative they are indices into infinite_vertices (-index-1).

Convex Hull

  • hull_points: list of ints
    List of indices of the points in points making up the convex hull.
  • hull_center: ndarray of floats
    Center of mass of the points making up the convex hull.
  • inside_vertices: ndarray of boolean
    Indicates for each vertex in vertices whether it is inside the convex hull.
  • in_hull(): Test if points are within the convex hull of the input points.
  • hull_area(): The area contained in the convex hull.

Outer Hull

The outer hull is constructed from the convex hull by expanding each point of the convex hull away from the center of mass of the convex hull such that randomly placing the same number of points in this area results on average in the same mean nearest-neighbor distance. This enlarged hull is needed for bootstrapping. Note that in some cases the outer hull can be smaller than the convex hull.

  • outer_hull_points: 2-d ndarray of floats
    Array of coordinates of the points of the outer hull.
  • outer_min_bound: ndarray of floats
    Lower left corner of the bounding-box of the outer hull.
  • outer_max_bound: ndarray of floats
    Upper right corner of the bounding-box of the outer hull.
  • in_outer_hull(): Test if points are within the outer hull.
  • outer_hull_area(): The area contained in the outer hull.

Bootstrap Voronoi Diagrams

  • random_points(): Generate random points within the area defined by the input points.

Plotting The Voronoi Diagram

  • plot_points(): Plot and optionally annotate the input points of the Voronoi diagram.
  • plot_center(): Plot the center of mass of the input points.
  • plot_vertices(): Plot and optionally annotate the vertices of the Voronoi diagram.
  • plot_distances(): Plot lines connecting the neighbors in the Voronoi diagram.
  • plot_ridges(): Plot the finite ridges of the Voronoi diagram.
  • `plot_infinite_ridges(): Plot the infinite ridges of the Voronoi diagram.
  • fill_regions(): Fill each finite region of the Voronoi diagram with a color.
  • fill_infinite_regions(): Fill each infinite region of the Voronoi diagram with a color.

Plotting The Convex Hull

  • plot_hull(): Plot the convex hull line containing the input points.
  • fill_hull(): Fill the convex hull.
  • plot_hull_center(): Plot the center of mass of the convex hull.
  • plot_outer_hull(): Plot the outer hull edge.
  • fill_outer_hull(): Fill the outer hull.

Usage

Generate 20 random points in 2-D:

import numpy as np
points = np.random.rand(20, 2)

Calculate the Voronoi diagram:

from thunderlab.voronoi import Voronoi
vor = Voronoi(points)

Retrieve nearest-neighbor distances and compute Voronoi areas:

dists = vor.nearest_distances
areas = vor.areas()

Generate random points drawn from a Poisson process with the same intensity and mean nearest neighbor distance as the data:

new_points = vor.random_points(poisson=True, mode='outer')

Plot Voronoi regions, distances, and input points:

import matplotlib.pyplot as plt
vor.fill_regions(colors=['red', 'green', 'blue', 'orange', 'cyan'], alpha=0.3)
vor.plot_distances(color='red')
vor.plot_points(text='p%d', c='c', s=100)

Methods

def in_hull(self, p)
Expand source code
def in_hull(self, p):
    
"""Test if points `p` are within convex hull of the input points. Parameters ---------- p: 2-d ndarray Array of points to be tested. Returns ------- inside: ndarray of booleans For each point in `p` whether it is inside the hull. """
inside = self.hull.find_simplex(p) >= 0 return inside

Test if points p are within convex hull of the input points.

Parameters

p : 2-d ndarray
Array of points to be tested.

Returns

inside : ndarray of booleans
For each point in p whether it is inside the hull.
def in_outer_hull(self, p)
Expand source code
def in_outer_hull(self, p):
    
"""Test if points `p` are within the outer hull. Parameters ---------- p: 2-d ndarray Array of points to be tested. Returns ------- inside: ndarray of booleans For each point in `p` whether it is inside the outer hull. """
inside = self.outer_hull.find_simplex(p) >= 0 return inside

Test if points p are within the outer hull.

Parameters

p : 2-d ndarray
Array of points to be tested.

Returns

inside : ndarray of booleans
For each point in p whether it is inside the outer hull.
def point_types(self)
Expand source code
def point_types(self):
    
"""The type of the Voronoi regions for each input point. Returns ------- points: ndarray of ints Indicates the type of Voronoi region associated with each point in `points`: - 2: finite region with all vertices inside hull, - 1: finite region, - 0: infinite region. """
points = np.zeros(len(self.vor.points), dtype=int) + 2 for i in range(len(self.vor.points)): for j, rp in enumerate(self.vor.ridge_points): if i in rp: if not np.all(self.inside_vertices[self.vor.ridge_vertices[j]]) and\ points[i] > 0: points[i] = 1 if not np.all(np.array(self.vor.ridge_vertices[j]) >= 0) and\ points[i] > -1: points[i] = 0 return points

The type of the Voronoi regions for each input point.

Returns

points : ndarray of ints

Indicates the type of Voronoi region associated with each point in points:

  • 2: finite region with all vertices inside hull,
  • 1: finite region,
  • 0: infinite region.
def ridge_lengths(self)
Expand source code
def ridge_lengths(self):
    
"""Length of Voronoi ridges between nearest neighbors. May be used, for example, as a weigth for `ridge_distances`. Returns ------- distances: ndarray of floats The length of each ridge in `ridge_vertices`. np.inf if one vertex is at infinity. """
ridges = np.zeros(len(self.vor.ridge_vertices)) for k, p in enumerate(self.vor.ridge_vertices): if np.all(np.array(p)>=0): p1 = self.vor.vertices[p[0]] p2 = self.vor.vertices[p[1]] ridges[k] = sp.minkowski_distance(p1, p2) else: ridges[k] = np.inf return ridges

Length of Voronoi ridges between nearest neighbors.

May be used, for example, as a weigth for ridge_distances.

Returns

distances : ndarray of floats
The length of each ridge in ridge_vertices. np.inf if one vertex is at infinity.
def ridge_areas(self)
Expand source code
def ridge_areas(self):
    
"""For each ridge the triangular area of the Voronoi region spanned by the center point and the ridge. Returns ------- areas: ndarray of floats For each ridge in `ridge_points` or `ridge_vertices` its corresponding triangular area. np.inf for infinite ridges. """
ridges = self.ridge_lengths() heights = 0.5*self.ridge_distances # area of a triangle: areas = 0.5*ridges*heights return areas

For each ridge the triangular area of the Voronoi region spanned by the center point and the ridge.

Returns

areas : ndarray of floats
For each ridge in ridge_points or ridge_vertices its corresponding triangular area. np.inf for infinite ridges.
def areas(self, mode='finite')
Expand source code
def areas(self, mode='finite'):
    
"""The areas of the Voronoi regions for each input point. Parameters ---------- mode: string - 'inside': Calculate area of finite Voronoi regions whose vertices are all inside the hull, set all other to zero. - 'finite_inside': Calculate area of all Voronoi regions. Consider only areas of finite ridges whose vertices are all inside the hull. - 'full': Calculate area of finite Voronoi regions only, set all other to zero. - 'finite': Calculate area of all Voronoi regions. From infinite regions only areas contributed by finite ridges are considered. Returns ------- areas: array of floats For each point in `points` its corresponding area. """
ridge_areas = self.ridge_areas() areas = np.zeros(len(self.vor.points)) if mode == 'inside': for i in range(len(self.vor.points)): a = 0.0 for j, rp in enumerate(self.vor.ridge_points): if i in rp: if ridge_areas[j] != np.inf and \ np.all(self.inside_vertices[self.vor.ridge_vertices[j]]): a += ridge_areas[j] else: a = 0.0 break areas[i] = a elif mode == 'full': for i in range(len(self.vor.points)): a = 0.0 for j, rp in enumerate(self.vor.ridge_points): if i in rp: if ridge_areas[j] != np.inf: a += ridge_areas[j] else: a = 0.0 break areas[i] = a elif mode == 'finite_inside': for i in range(len(self.vor.points)): a = 0.0 for j, rp in enumerate(self.vor.ridge_points): if i in rp and ridge_areas[j] != np.inf and \ np.all(self.inside_vertices[self.vor.ridge_vertices[j]]): a += ridge_areas[j] areas[i] = a elif mode == 'finite': for i in range(len(self.vor.points)): a = 0.0 for j, rp in enumerate(self.vor.ridge_points): if i in rp and ridge_areas[j] != np.inf: a += ridge_areas[j] areas[i] = a else: print('') print(f'Voronoi.areas(): unknown value "{mode}" for the mode parameter:') print('Use one of the following values:') print(' inside: Finite Voronoi regions whose vertices are all inside the hull.') print(' finite_inside: Use all areas corresponding to finite ridges whose vertices are all inside the hull.') print(' full: Finite Voronoi regions only.') print(' finite: Use all areas corresponding to finite ridges.') print('') return areas

The areas of the Voronoi regions for each input point.

Parameters

mode : string
  • 'inside': Calculate area of finite Voronoi regions whose vertices are all inside the hull, set all other to zero.
  • 'finite_inside': Calculate area of all Voronoi regions. Consider only areas of finite ridges whose vertices are all inside the hull.
  • 'full': Calculate area of finite Voronoi regions only, set all other to zero.
  • 'finite': Calculate area of all Voronoi regions. From infinite regions only areas contributed by finite ridges are considered.

Returns

areas : array of floats
For each point in points its corresponding area.
def hull_area(self)
Expand source code
def hull_area(self):
    """The area of the convex hull of the input points.
    
    Returns
    -------
    area: float
        The area of the convex hull.
    """
    # two sides of the simplex triangles:
    ab = self.hull.points[self.hull.simplices[:,0],:] - self.hull.points[self.hull.simplices[:,1],:]
    cb = self.hull.points[self.hull.simplices[:,2],:] - self.hull.points[self.hull.simplices[:,1],:]
    # area of each simplex is half of the absolute value of the cross product:
    area = 0.5*np.sum(np.abs(np.cross(ab, cb)))
    return area

The area of the convex hull of the input points.

Returns

area : float
The area of the convex hull.
def outer_hull_area(self)
Expand source code
def outer_hull_area(self):
    """The area of the outer hull.
    
    Returns
    -------
    area: float
        The area of the outer hull.
    """
    # two sides of the simplex triangles:
    ab = self.outer_hull.points[self.outer_hull.simplices[:,0],:] - self.outer_hull.points[self.outer_hull.simplices[:,1],:]
    cb = self.outer_hull.points[self.outer_hull.simplices[:,2],:] - self.outer_hull.points[self.outer_hull.simplices[:,1],:]
    # area of each simplex is half of the absolute value of the cross product:
    area = 0.5*np.sum(np.abs(np.cross(ab, cb)))
    return area

The area of the outer hull.

Returns

area : float
The area of the outer hull.
def random_points(self, n=None, poisson=True, mode='outer')
Expand source code
def random_points(self, n=None, poisson=True, mode='outer'):
    
"""Generate random points. Parameters ---------- n: int or None Number of random points to be generated. If None `n` is set to the number of input points. poisson: boolean If `True` then draw the number of points from a Poisson distribution with mean number of points given by `n`. mode: string - 'bbox': place points randomly in rectangular bounding box of the Voronoi diagram. - 'hull': place points randomly within convex hull of input data. - 'outer' place points randomly within outer hull. The mean nearest-neighbor distance between the generated points will be close to the observed one. Returns ------- points: 2-D array List of randomly generated points. """
# number of points: if n is None: n = self.npoints nn = n if poisson: nn = np.random.poisson(n) m = nn//2 if m < 5: m = 5 # get bounding box: if mode == 'outer': min_bound = self.outer_min_bound max_bound = self.outer_max_bound else: min_bound = self.min_bound max_bound = self.max_bound delta = np.max(max_bound - min_bound) points = np.zeros((0, self.ndim)) while len(points) < nn: # random points within bounding box: newpoints = np.random.rand(m, self.ndim) newpoints *= delta newpoints += min_bound if mode == 'outer': # only take the ones within outer hull: inside = self.in_outer_hull(newpoints) points = np.vstack((points, newpoints[inside])) elif mode == 'hull': # only take the ones within hull: inside = self.in_hull(newpoints) points = np.vstack((points, newpoints[inside])) elif mode == 'bbox': points = np.vstack((points, newpoints[np.all(newpoints<max_bound, axis=1),:])) else: print('') print(f'Voronoi.random_points(): unknown value "{mode}" for the mode parameter:') print('Use one of the following values:') print(' bbox: Place points within rectangular bounding box.') print(' hull: Place points inside the convex hull.') print(' outer: Place points inside the outer hull.') print('') return return points[:nn]

Generate random points.

Parameters

n : int or None
Number of random points to be generated. If None n is set to the number of input points.
poisson : boolean
If True then draw the number of points from a Poisson distribution with mean number of points given by n.
mode : string
  • 'bbox': place points randomly in rectangular bounding box of the Voronoi diagram.
  • 'hull': place points randomly within convex hull of input data.
  • 'outer' place points randomly within outer hull. The mean nearest-neighbor distance between the generated points will be close to the observed one.

Returns

points : 2-D array
List of randomly generated points.
def plot_points(self, ax, text=None, text_offs=(0, 0.05), text_align='center', **kwargs)
Expand source code
def plot_points(self, ax, text=None, text_offs=(0, 0.05),
                text_align='center', **kwargs):
    
"""Plot and optionally annotate the input points of the Voronoi diagram. Parameters ---------- ax: matplotlib.Axes The axes to be used for plotting. text: string or None If not None the string that is placed at each point. A '%d' is replaced by the index of the point. text_offset: tuple of numbers The offset of the point labels. text_align: string The horizontal alignment of the point labels. **kwargs: dict Key-word arguments that are passed on to the matplotlib.scatter() function. """
ax.scatter(self.points[:,0], self.points[:,1], **kwargs) if text is not None: for i, p in enumerate(self.points): s = text if '%' in text: s = text % i ax.text(p[0]+text_offs[0], p[1]+text_offs[1], s, ha=text_align)

Plot and optionally annotate the input points of the Voronoi diagram.

Parameters

ax : matplotlib.Axes
The axes to be used for plotting.
text : string or None
If not None the string that is placed at each point. A '%d' is replaced by the index of the point.
text_offset : tuple of numbers
The offset of the point labels.
text_align : string
The horizontal alignment of the point labels.
**kwargs : dict
Key-word arguments that are passed on to the matplotlib.scatter() function.
def plot_center(self, ax, **kwargs)
Expand source code
def plot_center(self, ax, **kwargs):
    
"""Plot the center of mass of the input points. Parameters ---------- ax: matplotlib.Axes The axes to be used for plotting. **kwargs: dict Key-word arguments that are passed on to the `matplotlib.plot()` function. """
ax.plot(self.center[0], self.center[1], 'o', **kwargs)

Plot the center of mass of the input points.

Parameters

ax : matplotlib.Axes
The axes to be used for plotting.
**kwargs : dict
Key-word arguments that are passed on to the matplotlib.plot() function.
def plot_vertices(self, ax, text=None, text_offs=(0, 0.05), text_align='center', **kwargs)
Expand source code
def plot_vertices(self, ax, text=None, text_offs=(0, 0.05),
                  text_align='center', **kwargs):
    
"""Plot and optionally annotate the vertices of the Voronoi diagram. Parameters ---------- ax: matplotlib.Axes The axes to be used for plotting. text: string or None If not None the string that is placed at each vertex. A '%d' is replaced by the index of the vertex. text_offset: tuple of numbers The offset of the vertex labels. text_align: string The horizontal alignment of the vertex labels. **kwargs: dict Key-word arguments that are passed on to the matplotlib.scatter() function. """
ax.scatter(self.vertices[:,0], self.vertices[:,1], **kwargs) if text is not None: for i, p in enumerate(self.vertices): s = text if '%' in text: s = text % i ax.text(p[0]+text_offs[0], p[1]+text_offs[1], s, ha=text_align)

Plot and optionally annotate the vertices of the Voronoi diagram.

Parameters

ax : matplotlib.Axes
The axes to be used for plotting.
text : string or None
If not None the string that is placed at each vertex. A '%d' is replaced by the index of the vertex.
text_offset : tuple of numbers
The offset of the vertex labels.
text_align : string
The horizontal alignment of the vertex labels.
**kwargs : dict
Key-word arguments that are passed on to the matplotlib.scatter() function.
def plot_distances(self, ax, **kwargs)
Expand source code
def plot_distances(self, ax, **kwargs):
    
"""Plot lines connecting the nearest neighbors in the Voronoi diagram. Parameters ---------- ax: matplotlib.Axes The axes to be used for plotting. **kwargs: dict Key-word arguments that are passed on to the matplotlib.plot() function. """
for i, p in enumerate(self.ridge_points): ax.plot(self.points[p, 0], self.points[p, 1], **kwargs)

Plot lines connecting the nearest neighbors in the Voronoi diagram.

Parameters

ax : matplotlib.Axes
The axes to be used for plotting.
**kwargs : dict
Key-word arguments that are passed on to the matplotlib.plot() function.
def plot_ridges(self, ax, **kwargs)
Expand source code
def plot_ridges(self, ax, **kwargs):
    
"""Plot the finite ridges of the Voronoi diagram. Parameters ---------- ax: matplotlib.Axes or None The axes to be used for plotting. **kwargs: dict Key-word arguments that are passed on to the matplotlib.plot() function. """
for i, p in enumerate(self.ridge_vertices): if np.all(np.array(p)>=0): ax.plot(self.vertices[p, 0], self.vertices[p, 1], **kwargs)

Plot the finite ridges of the Voronoi diagram.

Parameters

ax : matplotlib.Axes or None
The axes to be used for plotting.
**kwargs : dict
Key-word arguments that are passed on to the matplotlib.plot() function.
def plot_infinite_ridges(self, ax, **kwargs)
Expand source code
def plot_infinite_ridges(self, ax, **kwargs):
    
"""Plot the infinite ridges of the Voronoi diagram. Parameters ---------- ax: matplotlib.Axes The axes to be used for plotting. **kwargs: dict Key-word arguments that are passed on to the matplotlib.plot() function. """
for far_point, vertices in zip(self.infinite_vertices, self.vor.ridge_vertices): vertices = np.asarray(vertices) if not np.all(vertices >= 0): i = vertices[vertices >= 0][0] # finite end Voronoi vertex ax.plot([self.vor.vertices[i][0], far_point[0]], [self.vor.vertices[i][1], far_point[1]], **kwargs)

Plot the infinite ridges of the Voronoi diagram.

Parameters

ax : matplotlib.Axes
The axes to be used for plotting.
**kwargs : dict
Key-word arguments that are passed on to the matplotlib.plot() function.
def fill_regions(self, ax, inside=None, colors=None, **kwargs)
Expand source code
def fill_regions(self, ax, inside=None, colors=None, **kwargs):
    
"""Fill each finite region of the Voronoi diagram with a color. Parameters ---------- ax: matplotlib.Axes The axes to be used for plotting. inside: boolean or None True: plot only finite regions with all vertices inside the hull False: plot only finite regions with at least one vertex outside the hull None: plot all finite regions colors: list of colors or None If not None then these colors are used in turn to fill the regions. **kwargs: dict Key-word arguments that are passed on to the matplotlib.fill() function. """
c = 0 for region in self.regions: if not -1 in region: polygon = self.vertices[region] if len(polygon) > 0: inside_hull = self.in_hull(polygon) if inside is None or (inside and all(inside_hull)) or (not inside and any(inside_hull)): c += 1 if colors is None: ax.fill(polygon[:, 0], polygon[:, 1], lw=0, **kwargs) else: ax.fill(polygon[:, 0], polygon[:, 1], color=colors[c % len(colors)], lw=0, **kwargs)

Fill each finite region of the Voronoi diagram with a color.

Parameters

ax : matplotlib.Axes
The axes to be used for plotting.
inside : boolean or None
True: plot only finite regions with all vertices inside the hull False: plot only finite regions with at least one vertex outside the hull None: plot all finite regions
colors : list of colors or None
If not None then these colors are used in turn to fill the regions.
**kwargs : dict
Key-word arguments that are passed on to the matplotlib.fill() function.
def fill_infinite_regions(self, ax, colors=None, **kwargs)
Expand source code
def fill_infinite_regions(self, ax, colors=None, **kwargs):
    
"""Fill each infinite region of the Voronoi diagram with a color. Parameters ---------- ax: matplotlib.Axes The axes to be used for plotting. colors: list of colors or None If not None then these colors are used in turn to fill the regions. **kwargs: dict Key-word arguments that are passed on to the matplotlib.fill() function. """
c = 0 for region in self.infinite_regions: polygon = [] for p in region: if p >= 0: polygon.append(self.vertices[p]) else: polygon.append(self.infinite_vertices[-p-1]) if len(polygon) > 0: polygon = np.asarray(polygon) c += 1 if colors is None: ax.fill(polygon[:, 0], polygon[:, 1], lw=0, **kwargs) else: ax.fill(polygon[:, 0], polygon[:, 1], color=colors[c % len(colors)], lw=0, **kwargs)

Fill each infinite region of the Voronoi diagram with a color.

Parameters

ax : matplotlib.Axes
The axes to be used for plotting.
colors : list of colors or None
If not None then these colors are used in turn to fill the regions.
**kwargs : dict
Key-word arguments that are passed on to the matplotlib.fill() function.
def plot_hull(self, ax, **kwargs)
Expand source code
def plot_hull(self, ax, **kwargs):
    
"""Plot the hull line containing the input points. Parameters ---------- ax: matplotlib.Axes The axes to be used for plotting. **kwargs: dict Key-word arguments that are passed on to the matplotlib.plot() function. """
ax.plot(self.hull.points[self.hull_points, 0], self.hull.points[self.hull_points, 1], **kwargs)

Plot the hull line containing the input points.

Parameters

ax : matplotlib.Axes
The axes to be used for plotting.
**kwargs : dict
Key-word arguments that are passed on to the matplotlib.plot() function.
def fill_hull(self, ax, **kwargs)
Expand source code
def fill_hull(self, ax, **kwargs):
    
"""Fill the hull containing the input points with a color. Parameters ---------- ax: matplotlib.Axes The axes to be used for plotting. **kwargs: dict Key-word arguments that are passed on to the matplotlib.fill() function. """
ax.fill(self.hull.points[self.hull_points, 0], self.hull.points[self.hull_points, 1], lw=0, **kwargs)

Fill the hull containing the input points with a color.

Parameters

ax : matplotlib.Axes
The axes to be used for plotting.
**kwargs : dict
Key-word arguments that are passed on to the matplotlib.fill() function.
def plot_hull_center(self, ax, **kwargs)
Expand source code
def plot_hull_center(self, ax, **kwargs):
    
"""Plot the center of mass of the convex hull of the input points. Parameters ---------- ax: matplotlib.Axes The axes to be used for plotting. **kwargs: dict Key-word arguments that are passed on to the matplotlib.plot() function. """
ax.plot(self.hull_center[0], self.hull_center[1], 'o', **kwargs)

Plot the center of mass of the convex hull of the input points.

Parameters

ax : matplotlib.Axes
The axes to be used for plotting.
**kwargs : dict
Key-word arguments that are passed on to the matplotlib.plot() function.
def plot_outer_hull(self, ax, **kwargs)
Expand source code
def plot_outer_hull(self, ax, **kwargs):
    
"""Plot the hull line containing the input points and the vertices of the Voronoi diagram. Parameters ---------- ax: matplotlib.Axes or None The axes to be used for plotting. **kwargs: dict Key-word arguments that are passed on to the matplotlib.plot() function. """
ax.plot(self.outer_hull_points[:, 0], self.outer_hull_points[:, 1], **kwargs)

Plot the hull line containing the input points and the vertices of the Voronoi diagram.

Parameters

ax : matplotlib.Axes or None
The axes to be used for plotting.
**kwargs : dict
Key-word arguments that are passed on to the matplotlib.plot() function.
def fill_outer_hull(self, ax, **kwargs)
Expand source code
def fill_outer_hull(self, ax, **kwargs):
    
"""Fill the hull containing the input points and the vertices of the Voronoi diagram. Parameters ---------- ax: matplotlib.Axes The axes to be used for plotting. **kwargs: dict Key-word arguments that are passed on to the matplotlib.fill() function. """
ax.fill(self.outer_hull_points[:, 0], self.outer_hull_points[:, 1], lw=0, **kwargs)

Fill the hull containing the input points and the vertices of the Voronoi diagram.

Parameters

ax : matplotlib.Axes
The axes to be used for plotting.
**kwargs : dict
Key-word arguments that are passed on to the matplotlib.fill() function.