# How to get a Multipolygon object from Overpass QL?

I'm new to OSM/Overpass QL and there's a simple problem that I can't solve. I've been able to figure out the a query that will return the maritime boundary of a country. For example:

(rel[boundary="administrative"][name="Puerto Rico"][admin_level="4"];>;);way(r);(._;>;);out;


This returns the boundary as a set of ways. Using overpy I've been able to get the returned ways and display them using shapely/matplotlib, and I'm confident they have the data I need.

What I can't figure out is how can I reliably convert the returned set of ways to a MultiPolygon object, ideally a shapely.geometry.MyltiPolygon.

The above query is only an example that produces a fairly simple set of ways, but there are other queries that return much more complex results, such as

(rel[boundary="administrative"]["name:en"="French Polynesia"][admin_level="3"];>;);way(r);(._;>;);out;


I understand that Overpass doesn't handle (multi)polygons, but I assume there must be a simple way to achieve this using python/GeoPandas.shapely. I've tried some manual approaches, but they are error-prone and I'd be extremely surprised if there are no implemented librariessolutions to achieve this. Any ideas?

• Extracting such amounts of data makes only sense if you run your own Overpass instance! For coastlines you should take a look at openstreetmapdata.com/data/coast – mmd Oct 24 '17 at 12:58
• Thanks for the comment @mmd. I agree that's a lot of data if you do the whole world, but if you only care about a few countries, I think is a valid approach. That said, thanks for the reference to the coast data! – jorgeh Oct 24 '17 at 16:47

I was able to solve my problem using the function ways2poly function I wrote in this gist.

I also added code to run the Overpass queries (using overpy), in case anyone needs it.

I convert OSM ways to polygons using Shapely's linemerge, unary_union, and polygonize.

Like so:

import overpass
import shapely.geometry as geometry
from shapely.ops import linemerge, unary_union, polygonize

api = overpy.Overpass()

query = """ YOUR QUERY HERE """
response = api.query(query)

lss = [] #convert ways to linstrings

for ii_w,way in enumerate(response.ways):
ls_coords = []

for node in way.nodes:
ls_coords.append((node.lon,node.lat)) # create a list of node coordinates

lss.append(geometry.LineString(ls_coords)) # create a LineString from coords

merged = linemerge([*lss]) # merge LineStrings
borders = unary_union(merged) # linestrings to a MultiLineString
polygons = list(polygonize(borders))


Running your query and plotting in matplotlib: