Accessing Eevation Map via Cartopy

I’m trying to overlay an elevation map over a cartopy map but had no success in finding one. Is there potentially a way to create an elevation map using cartopy?

I tried browsing cartopy and Google to find potential digital elevation maps but was unable to work out a solution. I attempted to use a stock image, but the resolution was too low.

Maybe contextily could work if you just need some shaded relief. What sort of scale are you plotting at, and what kind of data are you trying to overlay on top of the elevation map, are they points, lines or polygons?

If you need proper DEMs and aren’t tied to Cartopy, I can also recommend PyGMT which has a pygmt.datasets.load_earth_relief function that downloads the SRTM15+V2.6 grid for a region of interest (see examples here). That function will return an xarray.DataArray that can then be plotted with grdimage. Happy to show some example code if you provide a bounding box region (min/max longitude/latitude coordinates).

Hello Wei,
Thanks for your response! I’m plotting a small scale map, focusing on a region in Northern California. I’m also overlaying weather station locations on top of the map, which are just points.
I’m not tied to Cartopy and can give DEM a try.


The longitude range is [-125, -123]
Tha latitude range is [38.8, 41]

Ok, I’ll just show a PyGMT example because that’s what I’m most familiar with. Here’s the code:

import numpy as np
import pygmt

# Load DEM for Northern California, with resolution of 30 arc minutes
region = [-125, -123, 38.8, 41]  # minx, maxx, miny, maxy
dem_grid = pygmt.datasets.load_earth_relief(resolution="30s", region=region)

# Generate random weather station points
rng = np.random.default_rng()
x_coordinates = rng.uniform(low=-125, high=-123, size=30)
y_coordinates = rng.uniform(low=38.8, high=41, size=30)

# %%
# Initialize PyGMT figure
fig = pygmt.Figure()
# Basemap with plot title (+t), x-axis and y-axis labels
fig.basemap(
    region=region,  # list [minx, maxx, miny, maxy]
    frame=["WSne+tNor-Cal GHCNd Stations", "xaf+lLongitude", "yaf+lLatitude"],
)
# Plot DEM with a colormap picked from
# https://docs.generic-mapping-tools.org/6.5/reference/cpts.html
fig.grdimage(grid=dem_grid, cmap="gmt/terra")
fig.colorbar(
    frame=["xaf+lElevation", "y+lm"], position="JCR"  # colorbar at Center-Right
)
# Plot rivers on land using darkblue color
fig.coast(
    rivers="a/0.5p,darkblue,solid"  # plot all rivers and canals (a) with 0.5p darkblue line
)
# Plot weather station points
fig.plot(x=x_coordinates, y=y_coordinates, style="c0.2c", fill="black")
# Save figure to PNG file
fig.savefig(fname="stations.png", dpi=300)
fig.show()

produces a map like this:

This tutorial might also be a good reference on plotting DEMs. The fig.plot() call also supports passing in data from a geopandas.GeoDataFrame or other formats, I wasn’t sure what your weather station points looked like so just used some random data. Do let me know though if you need help tweaking some of the map elements, or you can visit the GMT forum for similar examples like https://forum.generic-mapping-tools.org/t/map-of-seismic-stations-using-pygmt/828.

2 Likes

Hello Wei,
Thank you so much for these useful resources. I’m trying to use your code right now, but am having issues. Currently in the process of debugging.

GMTCLibNotFoundError: Error loading GMT shared library at ‘gmt.dll’.
Could not find module ‘gmt.dll’ (or one of its dependencies). Try using the full path with constructor syntax.
Error loading GMT shared library at ‘gmt_w64.dll’.
Could not find module ‘gmt_w64.dll’ (or one of its dependencies). Try using the full path with constructor syntax.
Error loading GMT shared library at ‘gmt_w32.dll’.
Could not find module ‘gmt_w32.dll’ (or one of its dependencies). Try using the full path with constructor syntax.

Ah yes, I’m guessing you did pip install pygmt instead of conda install pygmt? The error you’re getting is because of a missing GMT C library (PyGMT is a Python wrapper around the GMT C lib, so both need to be installed).

I’d recommend following the instructions at Installing — PyGMT to install pygmt from conda-forge. If you have Anaconda or miniconda installed, the command would be something like conda install -c conda-forge pygmt. Alternatively, you could try the Binder instance at GitHub - GenericMappingTools/try-gmt: Try GMT, PyGMT, and GMT.jl online! All in one place! if you’re after just a one-off map, and your data is not too big to upload.

Hi Wei,
Thank you so much for your suggestions! I was just able to get the map to work.

1 Like