Suppose we wish to understand if a given X/Y pair falls within a set of non-intersecting arbitrary polygons. For example, we are trying to see if a given lat,long falls within a Michigan minor civil division. There are 1500+ such divisions each one of which could be a polygon with hundreds or thousands of points. The following is a practical solution that is fast enough to return to a user.
- Break up Michigan into
$N$ box oriented tiles. We can easily calculate ifx, ywill only fall within a tile based on the value(x % WIDTH, y % HEIGHT). Select$N$ such that most (e.g. ~60%?) tiles are fully contained within a polygon, and the vast majority (e.g. 95%) are in no more than 4 polygons. I imagine this would happen easy at$N$ less than 100,000. - Preprocess the polygon shape files using Shapely to determine a map