My background is very different so this post made me sad. Doing tons of point-in-polygon tests is a common problem and should use spatial index. If we're testing N points against P polygons, the problem should scale O(N * log(P)) not O(N * P). The author was on the right track with bounding box tests, but that's more of a micro optimization. You want something like nested bounding boxes where each level contains 2 or more smaller ones, this is how you go from P to log(P) time - by never testing most of the polygons at all. I have my own preferred spatial index, but all of them tend toward this complexity improvement.