What is the fastest way to find the "visual" center of an irregularly shaped polygon?

Question

I have found a very good solution to this from MapBox called Polylabel. The full source is available on their Github too.

Essentially it tries to find the visual centre of the polygon as T Austin said.

Certain details suggest this may be a practical solution:

Unfortunately, calculating [the ideal solution ] is both complex and slow. The published solutions to the problem require either Constrained Delaunay Triangulation or computing a straight skeleton as preprocessing steps — both of which are slow and error-prone.

For our use case, we don’t need an exact solution — we’re willing to trade some precision to get more speed. When we’re placing a label on a map, it’s more important for it to be computed in milliseconds than to be mathematically perfect.

A quick note about usage though. The source code works great for Javascript out of the box however if you intend on using this with a "normal" polygon then you should wrap it in an empty array as the functions here take GeoJSONPolygons rather than normal polygons i.e.

var myPolygon = [[x1, y1], [x2, y2], [x3, y3]];
var center = polylabel([myPolygon]);

If you can convert the polygon into a binary image, then you can use the foundation that exists in the field of image processing, e.g.: A Fast Skeleton Algorithm on Block Represented Binary Images.

But this is not really reasonable in the general case, because of discretization errors and extra work.

However, maybe you find these useful:

• Straight skeleton of a simple polygon
• Determining the Skeleton of a Simple Polygon in (Almost) Linear Time

EDIT: Maybe you want to look for the point that is the center of the largest circle contained in the polygon. It is not necessarily always in the observed centre, but most of the time would probably give the expected result, and only in slightly pathological cases something that is totally off.

How about:

If the centroid of the polygon is inside the polygon then use it, else:

1) Extend a line from the centroid through the polygon dividing the polygon into two halves of equal area

2) The "visual center" is the point half way between the nearest point where the line touches the perimeter and the next point cutting the perimeter in the direction going away from the centroid

Here are a couple of pictures to illustrate it:

Have you looked into using the centroid formula?

http://en.wikipedia.org/wiki/Centroid

http://en.wikipedia.org/wiki/K-means_algorithm

Here are four different approaches I have tried.

1. cv2 based center of mass (get_center_of_mass)
2. shapely based representative point (get_representative_point)
3. cv2 + skimage.skeleton based center of mass of the skeletonized shape (get_skeleton_center_of_mass)
4. scipy based furthest distance to border (get_furthest_point_from_edge)

import numpy as np
import cv2
from shapely.geometry import Polygon
from skimage.morphology import skeletonize, medial_axis
from scipy.ndimage.morphology import distance_transform_edt
import matplotlib.pyplot as plt
H, W = 300, 300

def get_random_contour():
    xs = np.random.randint(0, W, 4)
    ys = np.random.randint(0, H, 4)
    cnt = np.array([[x,y] for x,y in zip(xs,ys)])
    mask = draw_contour_on_mask((H,W), cnt)
    cnt, _ = cv2.findContours(mask, 1, 2)
    cnt = cnt[0]
    return cnt

def draw_contour_on_mask(size, cnt):
    mask = np.zeros(size, dtype='uint8')
    mask = cv2.drawContours(mask, [cnt], -1, 255, -1)
    return mask

def get_center_of_mass(cnt):
    M = cv2.moments(cnt)
    cx = int(M['m10']/M['m00'])
    cy = int(M['m01']/M['m00'])
    return cx, cy

def get_representative_point(cnt):
    poly = Polygon(cnt.squeeze())
    cx = poly.representative_point().x
    cy = poly.representative_point().y
    return cx, cy

def get_skeleton_center_of_mass(cnt):
    mask = draw_contour_on_mask((H,W), cnt)
    skel = medial_axis(mask//255).astype(np.uint8) #<- medial_axis wants binary masks with value 0 and 1
    skel_cnt,_ = cv2.findContours(skel,1,2)
    skel_cnt = skel_cnt[0]
    M = cv2.moments(skel_cnt) 
    if(M["m00"]==0): # this is a line
        cx = int(np.mean(skel_cnt[...,0]))
        cy = int(np.mean(skel_cnt[...,1]))
    else:
        cx = int(M['m10']/M['m00'])
        cy = int(M['m01']/M['m00'])
    return cx, cy

def get_furthest_point_from_edge(cnt):
    mask = draw_contour_on_mask((H,W), cnt)
    d = distance_transform_edt(mask)
    cy, cx = np.unravel_index(d.argmax(), d.shape)
    return cx, cy

Here's my analysis on the topic:

• get_center_of_mass is fastest but, as mentioned in this thread, the center of mass can be located outside the shape for non-convex shapes.
• get_representative_point is also fast but the identified point, although always guaranteed to stay inside the shape (or with minor edits even multiple disconnected shapes!), does not have much if anything to do with center of the object
• get_skeleton_center_of_mass returns a perceptually nice center point, but is slow and requires a logic for disconnected shapes
• get_furthest_point_from_edge is relatively fast, generalizes easily to disconnected shapes and the center point is visually pleasing

rows = 4
cols = 4
markers = ['x', '+', "*", "o"]
colors = ['r','b','g','orange']
functions = [get_center_of_mass, get_representative_point, get_skeleton_center_of_mass, get_furthest_point_from_edge]

plt.figure(figsize=(2*cols, 2*rows, ))
for i in range(rows*cols): 
    cnt = get_random_contour()
    mask = draw_contour_on_mask((H,W), cnt)
    
    plt.subplot(cols,rows, i+1)
    plt.imshow(mask, cmap='gray')
    for c, m, f in zip(colors, markers, functions):
        l = f.__name__
        cx, cy = f(cnt)
        plt.scatter(cx, cy, c=c, s=100, label=l, marker=m, alpha=0.7)

plt.tight_layout()    
plt.legend(loc=3)
plt.show()

Here's how the algorithms, run on 100 random examples, compare in speed:

N_EXAMPLES = 100
cnts = [get_random_contour() for _ in range(N_EXAMPLES)]

%%time
_ = [get_center_of_mass(cnt) for cnt in cnts]
CPU times: user 7.07 ms, sys: 155 µs, total: 7.23 ms
Wall time: 5.75 ms

%%time
_ = [get_representative_point(cnt) for cnt in cnts]
CPU times: user 23.6 ms, sys: 7.84 ms, total: 31.4 ms
Wall time: 28.5 ms

%%time
_ = [get_skeleton_center_of_mass(cnt) for cnt in cnts]
CPU times: user 5.56 s, sys: 3.31 ms, total: 5.56 s
Wall time: 5.55 s

%%time
_ = [get_furthest_point_from_edge(cnt) for cnt in cnts]
CPU times: user 486 ms, sys: 53 µs, total: 486 ms
Wall time: 485 ms

Compute the centre position (x,y) of each edge of the polygon. You can do this by finding the difference between the positions of the ends of each edge. Take the average of each centre in each dimension. This will be the centre of the polygon.

Centroid method has already been suggested multiple times. I think this is an excellent resource that describes the process (and many other useful tricks with polygons) very intuitively:

http://paulbourke.net/geometry/polygonmesh/centroid.pdf

Also, for placing a simple UI label, it might be sufficient to just calculate the bounding box of the polygon (a rectangle defined by the lowest and highest x and y coordinates of any vertex in the polygon), and getting its center at:

{
    x = min_x + (max_x - min_x)/2,
    y = min_y + (max_y - min_y)/2
}

This is a bit faster than calculating the centroid, which might be significant for a real-time or embedded application.

Also notice, that if your polygons are static (they don't change form), you could optimize by saving the result of the BB center / center of mass calculation (relative to e.g. the first vertex of the polygon) to the data structure of the polygon.