Glacier flowlines

Centerlines

Computing the centerlines is the first task to run after the initialisation of the local glacier directories and of the local topography.

Our algorithm is an implementation of the procedure described by Kienholz et al., (2014). Appart from some minor changes (mostly the choice of certain parameters), we stayed close to the original algorithm.

The relevant task is tasks.compute_centerlines():

In [1]: graphics.plot_centerlines(gdir)

The glacier has a major centerline (the longest one), and tributaries (in this case two ). The Centerline objects are stored as a list, the last one being the major one. Navigation between inflows (can be more than one) and outflow (only one or none) is facilitated by the inflows and flows_to attributes:

In [2]: fls = gdir.read_pickle('centerlines')

In [3]: fls[0]  # a Centerline object
Out[3]: <oggm.core.centerlines.Centerline at 0x7fc77058aef0>

# make sure the first flowline realy flows into the major one:
In [4]: assert fls[0].flows_to is fls[-1]

At this stage, the centerline coordinates are still defined on the original grid, and they are not considered as “flowlines” by OGGM. A rather simple task (tasks.initialize_flowlines()) converts them to flowlines which now have a regular coordinate spacing along the flowline (which they will keep for the rest of the workflow). The tail of the tributaries are cut according to a distance threshold rule:

In [5]: graphics.plot_centerlines(gdir, use_flowlines=True)

Downstream lines

For the glacier to be able to grow we need to determine the flowlines downstream of the current glacier geometry. This is done by the tasks.compute_downstream_line() task:

In [6]: graphics.plot_centerlines(gdir, use_flowlines=True, add_downstream=True)

The downsteam lines area also computed using a routing algorithm minimizing the distance to cover and upward slopes.

Catchment areas

Each flowline has it’s own “catchment area”. These areas are computed using similar flow routing methods as the one used for determining the flowlines. Their role is to attribute each glacier pixel to the right tributory (this will also influence the later computation of the glacier widths).

In [7]: tasks.catchment_area(gdir)

In [8]: graphics.plot_catchment_areas(gdir)

Flowline widths

Finally, the glacier widths are computed in two steps.

First, we compute the geometrical width at each grid point. The width is drawn from the intersection of a line normal to the flowline and either the glacier or the catchment outlines (when there are tributaries):

In [9]: tasks.catchment_width_geom(gdir)

In [10]: graphics.plot_catchment_width(gdir)

Then, these geometrical widths are corrected so that the altitude-area distribution of the “flowline-glacier” is as close as possible as the actual distribution of the glacier using its full 2D geometry. This job is done by the tasks.catchment_width_correction() task:

In [11]: tasks.catchment_width_correction(gdir)

In [12]: graphics.plot_catchment_width(gdir, corrected=True)

Note that a perfect distribution is not possible since the sample size is not the same between the “1.5D” and the 2D representation of the glacier. OGGM deals with this by iteratively search for an altidute bin size which ensures that both representations have at least one element for each bin.

Implementation details

Shared setup for these examples:

import os
import geopandas as gpd
import oggm
from oggm import cfg, tasks, graphics
from oggm.utils import get_demo_file

cfg.initialize()
cfg.set_intersects_db(get_demo_file('rgi_intersect_oetztal.shp'))
cfg.PATHS['dem_file'] = get_demo_file('hef_srtm.tif')

base_dir = os.path.join(os.path.expanduser('~'), 'OGGM_docs', 'Flowlines')
entity = gpd.read_file(get_demo_file('HEF_MajDivide.shp')).iloc[0]
gdir = oggm.GlacierDirectory(entity, base_dir=base_dir)

tasks.define_glacier_region(gdir, entity=entity)
tasks.glacier_masks(gdir)
tasks.compute_centerlines(gdir)
tasks.initialize_flowlines(gdir)
tasks.compute_downstream_line(gdir)