Efficient Spatial Joins on Large Geometry Sets

We consider the following standard spatial-join problem: Given two sets of geometric objects in 2D (points, lines, polygonal areas, and collections of these), compute the spatial relations of all pairs of intersecting objects as a standard DE-9IM matrix. Most previous work focuses on one aspect of the problem, like the candidate generation, candidate reduction heuristics, efficient data structures, or parallelization. We provide a complete, fully functional, and carefully engineered implementation, as well as an extensive experimental evaluation of the relevance of various heuristics and of two variants for the exact geometry comparisons: our own implementation which preprocesses the geometries, and one using the GEOS library, which powers spatial joins in the widely used PostgreSQL+PostGIS. In particular, we find that the former speeds up spatial joins by more than an order of magnitude when complex geometries are involved. Our best approach can compute the full self join of the 1.4 billion geometries from OpenStreeMap in less than 3 hours on a commodity PC. This was out of reach for any existing implementation we tried. Our code and all the materials needed to reproduce our results are freely available on GitHub.