Ytukyg

I am trying to make a shapely Polygon from a binary mask, but I always end up with an invalid Polygon. How can I make a valid polygon from an arbitrary binary mask? Below is an example using a circular mask. I suspect that it is because the points I get from the mask contour are out of order, which is apparent when I plot the points (see images below).

import matplotlib.pyplot as plt

import numpy as np

from shapely.geometry import Point, Polygon

from scipy.ndimage.morphology import binary_erosion

from skimage import draw



def get_circular_se(radius=2):



    N = (radius * 2) + 1

    se = np.zeros(shape=[N,N])

    for i in range(N):

        for j in range(N):

                se[i,j] = (i - N / 2)**2 + (j - N / 2)**2 <= radius**2

    se = np.array(se, dtype="uint8")

    return se



    return new_regions, np.asarray(new_vertices)



#generates a circular mask

side_len = 512

rad = 100

mask = np.zeros(shape=(side_len, side_len))

rr, cc = draw.circle(side_len/2, side_len/2, radius=rad, shape=mask.shape)

mask[rr, cc] = 1



#makes a polygon from the mask perimeter

se = get_circular_se(radius=1)

contour = mask - binary_erosion(mask, structure=se)

pixels_mask = np.array(np.where(contour==1)[::-1]).T

polygon = Polygon(pixels_mask)

print polygon.is_valid



>>False



#plots the results

fig, ax = plt.subplots()

ax.imshow(mask,cmap='Greys_r')

ax.plot(pixels_mask[:,0],pixels_mask[:,1],'b-',lw=0.5)

plt.tight_layout()

plt.show()

In fact I already found a solution that worked for me, but maybe someone has a better one. The problem was indeed that my points were out of order. Input coordinate order is crucial for making valid polygons. So, one just has to put the points in the right order first. Below is an example solution using a nearest neighbor approach with a KDTree, which I've already posted elsewhere for related problems.

from sklearn.neighbors import KDTree



def polygonize_by_nearest_neighbor(pp):

    """Takes a set of xy coordinates pp Numpy array(n,2) and reorders the array to make

    a polygon using a nearest neighbor approach.



    """



    # start with first index

    pp_new = np.zeros_like(pp)

    pp_new[0] = pp[0]

    p_current_idx = 0



    tree = KDTree(pp)



    for i in range(len(pp) - 1):



        nearest_dist, nearest_idx = tree.query([pp[p_current_idx]], k=4)  # k1 = identity

        nearest_idx = nearest_idx[0]



        # finds next nearest point along the contour and adds it

        for min_idx in nearest_idx[1:]:  # skip the first point (will be zero for same pixel)

            if not pp[min_idx].tolist() in pp_new.tolist():  # make sure it's not already in the list

                pp_new[i + 1] = pp[min_idx]

                p_current_idx = min_idx

                break



    pp_new[-1] = pp[0]

    return pp_new



pixels_mask_ordered = polygonize_by_nearest_neighbor(pixels_mask)

polygon = Polygon(pixels_mask_ordered)

print polygon.is_valid



>>True



#plots the results

fig, ax = plt.subplots()

ax.imshow(mask,cmap='Greys_r')

ax.plot(pixels_mask_ordered[:,0],pixels_mask_ordered[:,1],'b-',lw=2)

plt.tight_layout()

plt.show()