I'm using the following procedure to calculate hexagonal polygon coordinates of a given radius for a square grid of a given extent (lower left --> upper right):
def calc_polygons(startx, starty, endx, endy, radius):
sl = (2 * radius) * math.tan(math.pi / 6)
# calculate coordinates of the hexagon points
p = sl * 0.5
b = sl * math.cos(math.radians(30))
w = b * 2
h = 2 * sl
origx = startx
origy = starty
# offsets for moving along and up rows
xoffset = b
yoffset = 3 * p
polygons = []
row = 1
counter = 0
while starty < endy:
if row % 2 == 0:
startx = origx + xoffset
else:
startx = origx
while startx < endx:
p1x = startx
p1y = starty + p
p2x = startx
p2y = starty + (3 * p)
p3x = startx + b
p3y = starty + h
p4x = startx + w
p4y = starty + (3 * p)
p5x = startx + w
p5y = starty + p
p6x = startx + b
p6y = starty
poly = [
(p1x, p1y),
(p2x, p2y),
(p3x, p3y),
(p4x, p4y),
(p5x, p5y),
(p6x, p6y),
(p1x, p1y)]
polygons.append(poly)
counter += 1
startx += w
starty += yoffset
row += 1
return polygons
This works well for polygons into the millions, but quickly slows down (and takes up very large amounts of memory) for large grids. I'm wondering whether there's a way to optimise this, possibly by zipping together numpy arrays of vertices that have been calculated based on the extents, and removing the loops altogether – my geometry isn't good enough for this, however, so any suggestions for improvements are welcome.
Decompose the problem into the regular square grids (not-contiguous). One list will contain all shifted hexes (i.e. the even rows) and the other will contain unshifted (straight) rows.
def calc_polygons_new(startx, starty, endx, endy, radius):
sl = (2 * radius) * math.tan(math.pi / 6)
# calculate coordinates of the hexagon points
p = sl * 0.5
b = sl * math.cos(math.radians(30))
w = b * 2
h = 2 * sl
# offsets for moving along and up rows
xoffset = b
yoffset = 3 * p
row = 1
shifted_xs = []
straight_xs = []
shifted_ys = []
straight_ys = []
while startx < endx:
xs = [startx, startx, startx + b, startx + w, startx + w, startx + b, startx]
straight_xs.append(xs)
shifted_xs.append([xoffset + x for x in xs])
startx += w
while starty < endy:
ys = [starty + p, starty + (3 * p), starty + h, starty + (3 * p), starty + p, starty, starty + p]
(straight_ys if row % 2 else shifted_ys).append(ys)
starty += yoffset
row += 1
polygons = [zip(xs, ys) for xs in shifted_xs for ys in shifted_ys] + [zip(xs, ys) for xs in straight_xs for ys in straight_ys]
return polygons
As you predicted, zipping results in much faster performance, especially for larger grids. On my laptop I saw 3x speedup when calculating 30 hexagon grid - 10x speed for 2900 hexagon grid.
>>> from timeit import Timer
>>> t_old = Timer('calc_polygons_orig(1, 1, 100, 100, 10)', 'from hexagons import calc_polygons_orig')
>>> t_new = Timer('calc_polygons_new(1, 1, 100, 100, 10)', 'from hexagons import calc_polygons_new')
>>> t_old.timeit(20000)
9.23395299911499
>>> t_new.timeit(20000)
3.12791109085083
>>> t_old_large = Timer('calc_polygons_orig(1, 1, 1000, 1000, 10)', 'from hexagons import calc_polygons_orig')
>>> t_new_large = Timer('calc_polygons_new(1, 1, 1000, 1000, 10)', 'from hexagons import calc_polygons_new')
>>> t_old_large.timeit(200)
9.09613299369812
>>> t_new_large.timeit(200)
0.7804560661315918
There might be an opportunity to create an iterator rather than the list to save memory. Depends on how your code is using the list of the polygons.