Get precise polygon coordinates for the outline of a solid shape, in order

Question

I want to create a map (via the HTML "map" and "area" tags) with several unusually-shaped areas. The shapes are detailed enough that I don't want to write out all the coordinates by hand, but I do want the map areas to be as precise as possible, so I want to generate them automatically.

I've converted each of my map areas to plain white-on-black images to make them easy to process, such as this:

import matplotlib.pyplot as plt
import numpy as np

# flatten image to a 2D array of 1s and 0s
img = plt.imread('map-area-1.png')
flat_img = np.argwhere(img[:,:,1] != 0)

I have tried using a shapely.geometry.Polygon:

poly = Polygon(flat_img)
ext = poly.exterior.coords

but the len of the supposed exterior outline is as many pixels as are in the image, including the filled inside.

I have also made a version of the image with a one-pixel-thick white outline, with black inside and outside, and attempted to work with the array in the same way. Unfortunately, I think the coordinates are out of order in this case.

I have tried sorting the coordinates by polar angle, but since the shape is unusual (the borders of it have protrusions and "dents", for example, although it is all one singular solid shape) this does not successfully sort it.

Attempting to use unsorted outline coordinates in the "area" tag makes a chunk of the inside of the shape interactive, but large portions are missing. I am guessing that it draws the shape in the order of the coordinates, which causes it to zig-zag between points across the shape and create gaps in the interactive portion.

Ideally, it would also be nice to have the minimum number of coordinates (i.e. anywhere there is a straight line on the shape, you only need one coordinate on either end of it, not every single coordinate along the line) but if it needs to be every single one, that's fine. Ultimately, I just need the coordinates in order of the shape's outline.

Answer 1

I suggest you use OpenCV to tackle that problem. Use cv2.findContours to extract the outline of the shape, afterwards approximate the contour using cv2.approxPolyDP. You can control the number of points by the epsilon parameter (in the example below controlled by the global variable PRECISION). Setting the precision to a very low value or completely removing the cv2.approxPolyDP call, you get the minimum number of points for the exact map in pixel resolution, without having repeated coordinates along straight lines (as you would using np.argwhere).

import cv2

PRECISION = 0.01  # Precision for contour approximation

# Get contours on the thresholded image
image       = cv2.imread('image.png', cv2.IMREAD_GRAYSCALE)
_, binary   = cv2.threshold(image, 127, 255, cv2.THRESH_BINARY)
contours, _ = cv2.findContours(binary, cv2.RETR_EXTERNAL, cv2.CHAIN_APPROX_SIMPLE)

contour = max(contours, key=cv2.contourArea)  # If there are multiple contours, select the largest one or...
contour = contours[0]                         # ...just take the first contour found

# Optionally simplify the contour to reduce the number of points
epsilon = PRECISION * cv2.arcLength(contour, True)
contour = cv2.approxPolyDP(contour, epsilon, True)

coordinates = contour[:, 0, :].tolist()
print("Ordered Coordinates:", coordinates)

Output:

Ordered Coordinates: [[12, 55], [13, 58], [51, 57], [63, 52], [71, 61], [77, 62], [66, 48], [68, 39], [58, 33], [58, 20], [46, 19], [44, 16], [35, 20], [39, 24], [35, 27], [39, 33], [33, 36], [31, 44], [21, 46]]

Answer 2

It's very inefficient, I'm afraid, but you can find a list of boundary pixels in order by:

• find an initial boundary point (e.g. where black and white pixels are adjacent)
• search anticlockwise around the last point until you find a new boundary point
• circumnavigate the shape until you get back to the start.

Potentially there would be a large number of pixels in the boundary.

import matplotlib.pyplot as plt
import numpy as np
istep = [ 1, 1, 0, -1, -1, -1,  0,  1 ]
jstep = [ 0, 1, 1,  1,  0, -1, -1, -1 ]

img = plt.imread('map-area-1.png')
bw = img[:,:,1] != 0
Ni, Nj = bw.shape

ival = []
jval = []
start_found = False
dir = 4         
# Find a starting point (left-right is black-white)
for i in range( 1, Ni ):
    if start_found: break
    for j in range( Nj ):
        if ( not bw[i-1,j] ) and bw[i,j]:
            ival.append( i )
            jval.append( j )
            start_found = True
            break

while True:
    if ival[-1] == ival[0] and jval[-1] == jval[0] and len( ival ) > 1: break  # back to beginning
    for d in range( dir + 1, dir + 8 ):                                        # scan directions from last
        dd = d % 8
        ip, jp  = ival[-1] + istep[dd], jval[-1] + jstep[dd]
        if ip < 0 or ip >= Ni or jp < 0 or jp >= Nj: continue                  # outside domain
        if bw[ip,jp]:
            ival.append( ip )
            jval.append( jp )
            dir = ( dd + 4 ) % 8
            break

# Fix the orientation
x = jval.copy()
y = [ Ni - 1 - i for i in ival ]

plt.scatter( x, y )
plt.show()