pythonimage-processingscipygeometryvoronoi

Area Voronoi Diagram for arbitrary shapes using scipy.spatial.Voronoi


Instead of computing the Point Voronoi Diagram, we need to compute the Area Voronoi Diagram. A point Voronoi Diagram is not a viable option due to the rather complex shapes of the figures involved.

The Area Voronoi Diagram is a suitable option and can be approximated by an adapted version of the Point Voronoi Diagram algorithm.

The approximation involves three steps. Instead of calculating the Voronoi Diagram from the center coordinates of the objects, it uses points on the contour of each object. Besides the desired effects this will also lead to an unnecessary partitioning of the objects, originating from points that lie on the same contour of an object, that has to be corrected.

  1. Evaluate sets of points on the contours of the figures
  2. Compute the Point Voronoi Diagram of those points
  3. delete Voronoi edges generated from points on the same figure

While 1. and 2. are straight forward using scipy.spatial.Voronoi we wonder how can be determined which edges were generated from points on the same figure without writing our own Voronoi Diagram implementation?

import numpy as np
import matplotlib.pyplot as plt

from scipy.spatial import Voronoi, voronoi_plot_2d
from skimage.filters import threshold_otsu
from skimage.measure import find_contours

def create_sample_img():
    object_selectors = []
    object_selectors.append((slice(50, 70), slice(100, 120)))
    object_selectors.append((slice(30, 80), slice(20, 60)))
    object_selectors.append((slice(10, 35), slice(80, 105)))
    object_selectors.append((slice(10, 90), slice(140, 190)))

    image = np.zeros((100,200))
    for selector in object_selectors:
        image[selector] = 1

    return image

def imshow_contours(image, contours, ax):
    ax.imshow(image, alpha=0.75)
    for contour in contours:
        ax.plot(contour[:, 1], contour[:, 0], linewidth=2)
    return ax

def area_voronoi(binary_image, ):
    downsampling = 3
    contour_index_mask = [] # to be used to see from which contour a point originates form
    contour_coordinates = []

    contours_grouped = find_contours(binary_image, 0)

    for contour_index, contours in enumerate(contours_grouped):
        contour_coordinates.extend(contours[::downsampling])
        contour_index_mask.extend([contour_index, ]*len(contours[::downsampling]))

    contour_coordinates = [[coord[1], coord[0]] for coord in contour_coordinates]

    vor = Voronoi(contour_coordinates)

    # TODO remove edges that originate from the same contours

    fig, ax = plt.subplots()
    imshow_contours(image, contours_grouped, ax)
    voronoi_plot_2d(vor, ax=ax)
    plt.show()


image = create_sample_img()
binary_image = image > threshold_otsu(image)
area_voronoi(binary_image)

Solution

  • Note that, when computing the Voronoi diagram of the points, you should have information about which two points each edge separates.

    Next, you should compute the labeling (connected component analysis) of the original image you extracted the points from. Now, for each edge, you can look up which two points it separates, and in the labeled image you can see which connected component each of these points belongs to. If the connected components are equal, delete the edge.

    From the connected component image you can also immediately derive the area of each connected component. This gives you enough information for the next step, which is applying a certain threshold on the distance between the two points and on the area ratio of the two connected components, to decide whether to keep or delete the edge.


    The scipy.spatial.Voronoi object has an attributes vor.ridge_points (a NumPy array) that tells you which two points are separated by the edge. So you can look up the coordinates in the input image for the two points, and thus look up the connected component each point belongs to.

    There is also an attribute vor.ridge_vertices (a list), which you will need to know how to draw the edge later. This list and the array above need to be kept in sync, when you delete an edge, you must delete both the row in the ridge_points array, and the corresponding element in the ridge_vertices list. Though deleting is expensive, so maybe keep a logical array indicating which edges to keep and which to discard?