pythongisshapefilefiona

Check if point is inside boundaries


I have a list of points describing the boundaries of Spain. I want to be able to tell whether a pair of lat,lon is within these boundaries. I have tried the following:

import shapefile
import matplotlib.pyplot as plt
from shapely.geometry import MultiPoint, Point, Polygon
from shapely.geometry.polygon import Polygon
sf = shapefile.Reader(r"\ESP_adm0.shp")
shapes = sf.shapes()
lat = []; lon = []
for i in range(len(shapes[0].points)):
    lon.append(shapes[0].points[i][0]);lat.append(shapes[0].points[i][1])

I know I am retrieving the points, because I'm able to plot and get the desired results:

plt.plot(lon,lat,'.', ms=0.1)

(plot in the link below) plot result

I do the following to get the poitns into a polygon:

    coords = list(zip(lat,lon))
    spain_pol = Polygon(coords)

And then I use the contains function, always getting false.

    spain_pol.contains(Point(0,42))
    spain_pol.contains(Point(42,0))

These both return false. In fact I haven't been able to get a single point I've tried to return a True.

I have tried all sorts of things, and I think I must be missing something fundamental. Perhaps the facts Spain has islands and there's more than one polygon is the problem? I'm lost. Any help welcome.


Solution

  • Just in case someone else has the same issue. The following code worked perfectly:

    import fiona
    from shapely.geometry import MultiPoint, Point, Polygon,shape
    from shapely.geometry.polygon import Polygon
    
    multipol = fiona.open(r"C:\Users\Jordi\Downloads\ESP_adm_shp\ESP_adm0.shp")
    multi = next(iter(multipol))
    
    point = Point(0,42)
    point.within(shape(multi['geometry']))
    

    This returns a very welcome "True" :)