9

I want to know given a latitude and longitude if a coordinate is land or sea

According to https://gis.stackexchange.com/questions/235133/checking-if-a-geocoordinate-point-is-land-or-ocean

from mpl_toolkits.basemap import Basemap
bm = Basemap()   # default: projection='cyl'
print bm.is_land(99.675, 13.104)  #True
print bm.is_land(100.539, 13.104)  #False

The problem is that basemap is deprecated. how di perform this with cartopy?

David Michael Gang
  • 7,107
  • 8
  • 53
  • 98

1 Answers1

14

A question which deals with point containment testing of country geometries using cartopy can be found at Polygon containment test in matplotlib artist.

Cartopy has the tools to achieve this, but there is no built-in method such as "is_land". Instead, you need to get hold of the appropriate geometry data, and query that using standard shapely predicates.

import cartopy.io.shapereader as shpreader
import shapely.geometry as sgeom
from shapely.ops import unary_union
from shapely.prepared import prep

land_shp_fname = shpreader.natural_earth(resolution='50m',
                                       category='physical', name='land')

land_geom = unary_union(list(shpreader.Reader(land_shp_fname).geometries()))
land = prep(land_geom)

def is_land(x, y):
    return land.contains(sgeom.Point(x, y))

This gives the expected results for two sample points:

>>> print(is_land(0, 0))
False
>>> print(is_land(0, 10))
True

If you have access to it, fiona will make this easier (and snappier):

import fiona
import cartopy.io.shapereader as shpreader
import shapely.geometry as sgeom
from shapely.prepared import prep

geoms = fiona.open(
            shpreader.natural_earth(resolution='50m',
                                    category='physical', name='land'))

land_geom = sgeom.MultiPolygon([sgeom.shape(geom['geometry'])
                                for geom in geoms])

land = prep(land_geom)

Finally, I produced (back in 2011) the shapely.vectorized functionality to speed up this kind of operation when testing many points at the same time. The code is available as a gist at https://gist.github.com/pelson/9785576, and produces the following proof-of-concept for testing land containment for the UK:

uk_containment

Another tool you may be interested in reading about is geopandas, as this kind of containment testing is one of its core capabilities.

pelson
  • 21,252
  • 4
  • 92
  • 99
  • I had a heck of a time getting all of the cartopy dependencies installed. I might add a little blurb warning users of your answer that it's not just a straight-forward pip install. – rocksNwaves Oct 16 '20 at 16:16
  • Now I am discovering those dependencies have a setup that is non-trivial, even after following the sparse directions on the GEOS and cartopy web-pages. Then after installation of all requirements, imports error out unless you install even more dependencies that are not listed anywhere except on this github issue: https://github.com/SciTools/cartopy/issues/1239 – rocksNwaves Oct 16 '20 at 16:41