Create local topography maps from different DEM sources with OGGM

Create local topography maps from different DEM sources with OGGM#

There are a number of datasets available out-of-the box in OGGM. This Notebook will show you how to download the original sources and create the local glacier centered map. It is also possible to use your own DEM data in OGGM.

Set-up#

import rioxarray as rioxr
import matplotlib.pyplot as plt
import oggm
from oggm import cfg, utils, workflow, tasks, graphics
from oggm.core import gis
cfg.initialize(logging_level='WARNING')
cfg.PARAMS['border'] = 10

2024-08-25 21:08:38: oggm.cfg: Reading default parameters from the OGGM `params.cfg` configuration file.
2024-08-25 21:08:38: oggm.cfg: Multiprocessing switched OFF according to the parameter file.
2024-08-25 21:08:38: oggm.cfg: Multiprocessing: using all available processors (N=4)
2024-08-25 21:08:39: oggm.cfg: PARAMS['border'] changed from `80` to `10`.

RGI outlines#

We use the RGI outlines to identify the necessary DEM tiles. If you haven’t downloaded the RGI files yet, this will also download them. Feel free to use your desired RGI-ID here, otherwise let’s use the Hintereisferner glacier as an example.

entity = utils.get_rgi_glacier_entities(['RGI60-11.00897'])
entity.plot();

../../_images/2438363810202ff8748ec15a5f65e2d4f7e364d334fe21bc351a529c2efd3f1d.png

Choose a DEM source (e.g. SRTM)#

If not specifying anything, OGGM will use it’s default settings, i.e. NASADEM for mid- and low-latitudes (60°S-60°N). However, this needs registration at NASA Earthdata (see “Register” below). Here, we choose the SRTM source as example DEM (no registration necessary).

# Let's make a working directory for this DEM 
cfg.PATHS['working_dir'] = utils.gettempdir('default', reset=True)
gdir = workflow.init_glacier_directories(entity)[0]
tasks.define_glacier_region(gdir, source='SRTM')
# if not yet installed, you have to install rasterio

2024-08-25 21:08:44: oggm.workflow: Execute entity tasks [GlacierDirectory] on 1 glaciers

You can access the (reprojected and interpolated) DEM file in the working directory:

dem_path = gdir.get_filepath('dem')
dem_path

'/tmp/OGGM/default/per_glacier/RGI60-11/RGI60-11.00/RGI60-11.00897/dem.tif'

It is a geotiff file. Xarray can open them thanks to rasterio:

da = rioxr.open_rasterio(dem_path)
f, ax = plt.subplots()
da.plot(cmap='terrain', ax=ax);
# Add the outlines
gdir.read_shapefile('outlines').plot(ax=ax, color='none', edgecolor='black');

../../_images/36774181971dc74e7e9cd6c0ff6403ffe36fa317419d96679d0e5112cb40c320.png

The source of the DEM is documented in the directory itself:

with open(gdir.get_filepath('dem_source'), 'r') as f:
    print(f.read())

SRTM V4

Original resolution: 3 arcsec (~90m)
Date of acquisition: February 2000
Date range: 2000-2000

# Data Citation

We kindly ask  any users to  cite this data  in any published  material
produced using this data, and if possible link web pages to the CIAT-CSI
SRTM website (http://srtm.csi.cgiar.org).

Citations should be made as follows:

Jarvis A., H.I. Reuter, A.  Nelson, E. Guevara, 2008, Hole-filled seamless SRTM
data V4, International  Centre for Tropical  Agriculture (CIAT), available
from http://srtm.csi.cgiar.org,


# Data files

srtm_39_03.tif

OGGM is neither the owner nor the distributer of these datasets! OGGM only provides tools to access it. It is your responsibility as the data user to read the individual usage requirements and cite and acknowledge the original data sources accordingly.

OGGM provided datasets#

At the moment OGGM is able to download and process the following DEM sources:

for src in utils.DEM_SOURCES:
    print('{:<10}: {}'.format(src, gis.DEM_SOURCE_INFO[src].split('\n')[0]))

GIMP      : GIMP V1.1
ARCTICDEM : ARCTICDEM V3
RAMP      : RAMP V2
TANDEM    : TanDEM-X 90m
AW3D30    : AW3D30 V1804
MAPZEN    : MAPZEN TERRAIN TILES
DEM3      : VIEWFINDER PANORAMAS DEMs
ASTER     : ASTGTMV3
SRTM      : SRTM V4
REMA      : REMA V1.1
ALASKA    : ALASKA V3
COPDEM30  : Copernicus DEM GLO-30
COPDEM90  : Copernicus DEM GLO-90
NASADEM   : NASADEM Merged DEM Global

Register for online datasets#

The default DEM source for low and mid-latitudes (60°S-60°N), NASADEM, requires a user account to download data, so you need to register at NASA Earthdata. There are other DEM sources where a registration is necessary; for ASTGTMV3 at NASA Earthdata, for TanDEM-X at DLR, and for COPDEM at spacedata.copernicus.eu/.

After that you can use the command line functionality oggm_netrc_credentials to store your user credentials in a local ~/.netrc file. Your user credentials are only stored locally and are only used by the download function for authentification with the original DEM source. Credentials are not needed if you use the RGI-TOPO data (see below).

Use pre-processed DEMs from RGI-TOPO#

The RGI-TOPO dataset is an RGI-provided dataset in beta release. These data are available for everyone, and were created with OGGM. Of course you can easily use these data in OGGM as well:

# use NASADEM, the default DEM for low and mid-latitudes in OGGM, you can also change to e.g. 'COPDEM'
from oggm.shop import rgitopo
cfg.PATHS['working_dir'] = utils.gettempdir('rgitopo', reset=True)
gdir = rgitopo.init_glacier_directories_from_rgitopo(['RGI60-11.00897'], dem_source='NASADEM')[0]
graphics.plot_domain(gdir)

2024-08-25 21:08:45: oggm.workflow: init_glacier_directories from prepro level 1 on 1 glaciers.
2024-08-25 21:08:45: oggm.workflow: Execute entity tasks [gdir_from_prepro] on 1 glaciers
2024-08-25 21:08:45: oggm.workflow: Execute entity tasks [select_dem_from_dir] on 1 glaciers

../../_images/effef083b4d56da861f6b5e6997039851a66f8235302ec7259a2cc7a0535c321.png

Use another DEM source#

Using RGI-TOPO DEMs is by far the easiest since all data is prepared for you and ready to use. But if you really want, you can go back to the original data sources:

# Let's make a working directory for this DEM 
cfg.PATHS['working_dir'] = utils.gettempdir('alternative')
try:
    gdir = workflow.init_glacier_directories(entity)[0]
    tasks.define_glacier_region(gdir, source='DEM3')
except oggm.exceptions.InvalidDEMError as err:
    print(err)

2024-08-25 21:08:46: oggm.workflow: Execute entity tasks [GlacierDirectory] on 1 glaciers

Let’s check that the source text is updated as well:

with open(gdir.get_filepath('dem_source'), 'r') as f:
    print(f.read())

VIEWFINDER PANORAMAS DEMs

Original resolution: 3 arcsec (~90m)
Date of acquisition: Depending on original source

These data are downloaded from Jonathan de Ferranti's webpage:
http://viewfinderpanoramas.org/dem3.html

# Data Citation

There is no recommended citation for these data. Please refer to the website
above in case of doubt.

# Data files

L32.tif

f, ax = plt.subplots()
da_dem3 = rioxr.open_rasterio(gdir.get_filepath('dem'))
da_dem3.plot(cmap='terrain', ax=ax);
gdir.read_shapefile('outlines').plot(ax=ax, color='none', edgecolor='black');

../../_images/fe39fecf07b90399b205a009ca24f5b330393797b6b2adcaf2767bd1130349e1.png

There might not be much difference a first sight, but by subtracting them the difference become clear:

f, ax = plt.subplots()
(da_dem3 - da).plot(ax=ax);
plt.title('DEM3 - SRTM');
gdir.read_shapefile('outlines').plot(ax=ax, color='none', edgecolor='black');

../../_images/309cec19c28860925b2067b6abe2c9190c8f5da83fe9eb8aeed1f8387d8cc60d.png

Regional DEMs / DEM availability#

Of course not all sources are available for every glacier as some of the DEMs are regional only. If we for example try the GIMP DEM, which is a Greenland specific DEM, it will not work for glaciers outside that region:

# Let's make a working directory for this DEM 
cfg.PATHS['working_dir'] = utils.gettempdir('gimp', reset=True)
try:
    gdir = workflow.init_glacier_directories(entity)[0]
    tasks.define_glacier_region(gdir, source='GIMP')
except oggm.exceptions.InvalidWorkflowError as err:
    print(err)

Source: GIMP is not available for the glacier RGI60-11.00897 with border 10

2024-08-25 21:08:47: oggm.workflow: Execute entity tasks [GlacierDirectory] on 1 glaciers
2024-08-25 21:08:48: oggm.core.gis: InvalidWorkflowError occurred during task define_glacier_region on RGI60-11.00897: Source: GIMP is not available for the glacier RGI60-11.00897 with border 10

User provided DEM#

Users should be able to use any DEM file which can be opened by rasterio (i.e. geotiff). Here, we use a subset SRTM file shipped with OGGM as an example:

custom_dem_path = utils.get_demo_file('hef_srtm.tif')
custom_dem_path

'/github/home/.oggm/oggm-sample-data-e6da53aea8438a06d7702769e73fe5df8ab0d669/test-files/hef_srtm.tif'

We tell OGGM to use it by changing the entry in the RGI table and by giving the path to the file:

cfg.PATHS['dem_file'] = custom_dem_path

cfg.PATHS['working_dir'] = utils.gettempdir('user', reset=True)
gdir = workflow.init_glacier_directories(entity)[0]
tasks.define_glacier_region(gdir, source='USER')

2024-08-25 21:08:48: oggm.workflow: Execute entity tasks [GlacierDirectory] on 1 glaciers

Now the user provided DEM is used:

f, ax = plt.subplots()
da_user = rioxr.open_rasterio(gdir.get_filepath('dem'))
da_user.plot(cmap='terrain', ax=ax);
gdir.read_shapefile('outlines').plot(ax=ax, color='none', edgecolor='black');

../../_images/7c80c237115b0366eadb2a738394b97e6ec4db676ca43bae7ef73f7ac002a00f.png

The border value, or how to chose the size of the topographic map#

It is possible to specify the extent of the local topographic map. All maps are centered on the glacier and the size of the map is determined in grid points around the glacier. The number of grid points that was used in this example are 10 in order to save storage. But depending on your study you might need a larger topographic map.

OGGM’s pre-processed directories come in 4 border sizes: 10, 40, 80 and 160. But if you process the topography yourself you can chose every value.

# print the currently used number of gridpoints around a glacier
cfg.PARAMS['border']

cfg.PARAMS['border'] = 1

2024-08-25 21:08:48: oggm.cfg: PARAMS['border'] changed from `10` to `1`.

# Let's make a working directory for this DEM 
cfg.PATHS['working_dir'] = utils.gettempdir('border1')
gdir = workflow.init_glacier_directories(entity)[0]
tasks.define_glacier_region(gdir)
da = rioxr.open_rasterio(gdir.get_filepath('dem'))
f, ax = plt.subplots()
da.plot(cmap='terrain', ax=ax);
# Add the outlines
gdir.read_shapefile('outlines').plot(ax=ax, color='none', edgecolor='black');

2024-08-25 21:08:48: oggm.workflow: Execute entity tasks [GlacierDirectory] on 1 glaciers

../../_images/136ae47b1dbc14fee064eb36ea613bb16c977574a33021eb6e9cea1f437a4f75.png

cfg.PARAMS['border'] = 100

2024-08-25 21:08:49: oggm.cfg: PARAMS['border'] changed from `1` to `100`.

# Let's make a working directory for this DEM 
cfg.PATHS['working_dir'] = utils.gettempdir('border100')
gdir = workflow.init_glacier_directories(entity)[0]
tasks.define_glacier_region(gdir)
da = rioxr.open_rasterio(gdir.get_filepath('dem'))
f, ax = plt.subplots()
da.plot(cmap='terrain', ax=ax);
# Add the outlines
gdir.read_shapefile('outlines').plot(ax=ax, color='none', edgecolor='black');

2024-08-25 21:08:49: oggm.workflow: Execute entity tasks [GlacierDirectory] on 1 glaciers

../../_images/1a2d4920ee223b15399d906b52125fea476e7b2bb28f8f39537aebfa9e296898.png

What’s next?#

return to the OGGM documentation
back to the table of contents