I'm trying to create a layer of all cities with 100km of the coast of Australia. The following is returning coastCities as an empty DataFrame.

I have performed the same operations in ArcGIS and it works fine, returning 135 records. I suspect it's the within() part that I'm getting wrong but I can't work out what's going on. cities is a GeoDataFrame with several thousand records with point geometry and Australia is a GeoDataFrame with 1 MultiPolygon row. The buffering works fine, I plotted that and it's good.

import geopandas as gpd
import matplotlib.pyplot as plt

Australia = gpd.read_file('australia_Main_Dissolve.shp')
cities = gpd.read_file('ne_10m_populated_places.shp')

AusCoast_buffer = Australia.boundary.buffer(100000)

coastCities = cities.loc[cities.within(AusCoast_buffer)]

Edit: corrected description of variables

1 Answer 1


I found the answer in this post.

By changing the lines after the file reads to

AusCoast_buffer = gpd.GeoDataFrame(Australia.boundary.buffer(100000), columns=['geometry'])

AusCoast_buffer_u = AusCoast_buffer.geometry.unary_union

AusCoast_cities = cities.loc[cities.within(AusCoast_buffer_u)]

I get what I want. It selects the correct records in cities. However AusCoast_buffer_u is of type shapely.geometry.multipolygon.MultiPolygon and if I try to convert it to a GeoDataFrame (so that it can be easily plotted) before doing the selection then within() goes back to not selecting any cities records. I really don't understand what is happening here, any explanation would be very welcome.

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge you have read our privacy policy.

Not the answer you're looking for? Browse other questions tagged or ask your own question.