I'm using PyOsmium. Given a node on OpenStreetMap, for example a village or a city like Chicago, how should I derive information such as the This probably involves using asked 27 Jul '21, 12:40  fsaler |
One Answer:
You will save yourself a lot of trouble if you can afford to run a local "Nominatim" installation and then simply use Nominatim's reverse geocoding ability. Alternatively, try "Photon" for the same purpose; there are ready-to-go Photon data dumps available for download, obviating the need for a database and an import process. Otherwise, the likely easiest path is to pre-filter a planet file with the command-line osmium utility (only keep stuff tagged boundary=administrative and with an admin_level of interest to you), then load the entirety of that file into Pyosmium, construct geometries, use the Pyosmium/Shapely interoperability to stuff the geometries into a Shapely R-Tree. Once that is done, you can efficiently query which administrative bounding boxes a point is in; if this yields more than one result (eg half the points in New Mexico are also in the Texas bounding box), you need to do proper (and more costly) point-in-polygon checks for the resulting polygons but that is also easy with Shapely. A third option is to do essentially what I wrote above but employ a PostGIS database to do it, i.e. load admin boundaries into PostGIS with osm2pgsql, then run SQL queries to find which polygons a point is in. answered 27 Jul '21, 13:09  Frederik Ramm ♦ edited 27 Jul '21, 17:18  alester |
I'm going to try the Photon route, it sounds smartest. Where can I find the data dumps you mention? I need one for the whole planet, but I don't need data on ways or most kinds of features, only inhabited places (hamlet or bigger).