Skip to content
{{ message }}

Instantly share code, notes, and snippets.

# josiahcarlson/alternate_tile.py

Created Sep 1, 2010
 ''' Copyright 2010 Josiah Carlson Released into the public domain alternate_tile.py This will generate the integer-grid tiles of a provided polygon, similar in fashion to Polygon.Util.tile(), though this uses a BSP-like algorithm to partition the space for faster overall execution on large polygons with many tiles. Discussed: http://dr-josiah.blogspot.com/2010/08/binary-space-partitions-and-you.html Requires the Polygon library: http://polygon.origo.ethz.ch/ The underlying GPC library requires a license for commercial users: http://www.cs.man.ac.uk/~toby/alan/software/ ''' from collections import defaultdict import math import Polygon def _find_split(count_dict): ''' When provided a dictionary of counts for the number of points inside each of the unit grid rows/columns, this function will return the best column choice so as to come closest to cutting the points in half. It will also return the score, lower being better. Returns: (cutoff, score) ''' # find the prefix sums tmp = {} for i in xrange(min(count_dict), max(count_dict)+1): tmp[i] = tmp.get(i-1, 0) + count_dict.get(i, 0) by_index = sorted(tmp.items()) # the target number of points midpoint = by_index[-1] // 2 # calculate how far off from the target number each choice would be totals = [] for i in xrange(1, len(by_index)): totals.append(abs(by_index[i-1] - midpoint)) # choose the best target number mi = min(totals) index = totals.index(mi) return by_index[index+1], totals[index] def _single_poly(polygon, dim, maxv): for poly in polygon: if max(pt[dim] for pt in poly) > maxv: return False return True def tile(polygon): ''' When provided with a Polygon(), this function will yield tiles of the original polygon on the integer grid. ''' _int = int _floor = math.floor work = [polygon] while work: # we'll use an explicit stack to ensure that degenerate polygons don't # blow the system recursion limit polygon = work.pop() # find out how many points are in each row/column of the grid xs = defaultdict(_int) ys = defaultdict(_int) for poly in polygon: for x,y in poly: xs[_int(_floor(x))] += 1 ys[_int(_floor(y))] += 1 # handle empty polygons gracefully if not xs: continue # handle top and right-edge border points mvx = max(max(x for x,y in poly) for poly in polygon) vx = _int(_floor(mvx)) if len(xs) > 1 and mvx == vx: xs[vx-1] += xs.pop(vx, 0) mvy = max(max(y for x,y in poly) for poly in polygon) vy = _int(_floor(mvy)) if len(ys) > 1 and mvy == vy: ys[vy-1] += ys.pop(vy, 0) # we've got a single grid, yield it if len(xs) == len(ys) == 1: yield polygon continue # find the split if len(xs) < 2: spx, countx = xs.items() countx *= 3 else: spx, countx = _find_split(xs) if len(ys) < 2: spy, county = ys.items() county *= 3 else: spy, county = _find_split(ys) # get the grid bounds for the split minx = min(xs) maxx = max(xs) miny = min(ys) maxy = max(ys) # actually split the polygon and put the results back on the work # stack if (countx < county and not _single_poly(polygon, 0, minx + 1.0)) or _single_poly(polygon, 1, miny + 1.0): work.append(polygon & Polygon.Polygon([(minx, miny), (minx, maxy+1), (spx, maxy+1), (spx, miny)])) work.append(polygon & Polygon.Polygon([(spx, miny), (spx, maxy+1), (maxx+1, maxy+1), (maxx+1, miny)])) else: work.append(polygon & Polygon.Polygon([(minx, miny), (minx, spy), (maxx+1, spy), (maxx+1, miny)])) work.append(polygon & Polygon.Polygon([(minx, spy), (minx, maxy+1), (maxx+1, maxy+1), (maxx+1, spy)])) # Always recurse on the smallest set, which is a trick to ensure that # the stack size is O(log n) . if work[-2].nPoints() < work[-1].nPoints(): work.append(work.pop(-2))
to join this conversation on GitHub. Already have an account? Sign in to comment