#!/bin/python3 import osmapi import shapely import pickle from os.path import exists OSM_KOGGENLAND_ID = 161930 PICKLE_FILE = 'polygon.pickle' def chunks(lst, n): """Yield successive n-sized chunks from lst.""" for i in range(0, len(lst), n): yield lst[i:i + n] def dump_polygon(): osm_api = osmapi.OsmApi() koggenland = osm_api.RelationGet(OSM_KOGGENLAND_ID) ways = filter(lambda x: x['type'] == 'way', koggenland['member']) shell_ids = [] for way in ways: way = osm_api.WayGet(way['ref']) node_ids = way['nd'] if not shell_ids: shell_ids.append(node_ids[0]) for node_id in node_ids[1:]: shell_ids.append(node_id) node_data = dict() for chunk in chunks(shell_ids, 100): node_data |= osm_api.NodesGet(chunk) shell = list() for shell_id in shell_ids: node = node_data[shell_id] shell.append((node['lat'], node['lon'])) koggenland_polygon = shapely.Polygon(shell) with open(PICKLE_FILE, 'wb') as f: pickle.dump(koggenland_polygon, f) def main(): if not exists(PICKLE_FILE): print('Using OSM API to create polygon...') dump_polygon() print('Dumped polygon to file.') with open(PICKLE_FILE, 'rb') as f: area = pickle.load(f) print('Loaded polygon from file.') point_out = shapely.Point(52.68811, 4.91267) point_in = shapely.Point(52.69168, 4.91324) print('out', area.contains(point_out)) print('in', area.contains(point_in)) if __name__ == '__main__': main()