Create A Network Grid from Raster Grid#

This notebook demonstrates how to extract a NetworkModelGrid from a RasterModelGrid using the ‘create_network_from_raster’ function

[ ]:

import matplotlib.pyplot as plt
import numpy as np

from landlab.grid.create_network import network_grid_from_raster
from landlab.io import esri_ascii
from landlab.plot.graph import plot_graph

Create a Raster Model Grid#

Now let’s create a raster model grid from the .asc data

[ ]:

with open("hugo_site.asc") as fp:
    grid = esri_ascii.load(fp, name="topographic__elevation", at="node")
z = grid.at_node["topographic__elevation"]

grid.status_at_node[grid.nodes_at_right_edge] = grid.BC_NODE_IS_FIXED_VALUE
grid.status_at_node[np.isclose(z, -9999.0)] = grid.BC_NODE_IS_CLOSED

Visualize Basin topography

[ ]:

grid.imshow(
    "topographic__elevation",
    plot_name="Basin topography",
    color_for_closed=None,
    colorbar_label="$z$ [m]",
)

Create a Network Model Grid#

Now let’s create a network model grid using the ‘create_network_from_raster’ function:

[ ]:

network_grid = network_grid_from_raster(
    grid,
    minimum_channel_threshold=12000.0,
    include=["drainage_area", "topographic__elevation"],
)

Let’s plot our network#

The nodes are numbered in landlab convention of lower-left to upper-right.

[ ]:

plot_graph(network_grid, at="node")
plt.title("Nodes");

As are the links,

[ ]:

plot_graph(network_grid, at="node,link", with_id="link")
plt.title("Links");

[ ]:

grid.imshow(
    "topographic__elevation",
    plot_name="Basin topography",
    color_for_closed=None,
    colorbar_label="$z$ [m]",
)
ax = plt.gca()
xlim, ylim = ax.get_xlim(), ax.get_ylim()

plot_graph(network_grid, at="link", with_id=False, axes=ax)

ax.set_xlim(xlim)
ax.set_ylim(ylim)

[ ]:

if len(network_grid.nodes_at_link) > (len(network_grid.nodes) - 1):
    print("there are duplicate links")
else:
    print("link count looks good!")

Reduce the numer of nodes in channel segments#

In the above, we have preserved every node in each channel segment. You can reduce the number of nodes in each segment by passing a reducer to the network_grid_from_raster function. A reducer is simply a function that accepts a list of nodes defining a segment. Although you can write your own, landlab provides a set of commonly used reducers.

[ ]:

from landlab.grid.create_network import (
    AtMostNodes,
    JustEndNodes,
    SpacingAtLeast,
    spacing_from_drainage_area,
)

As an example, the following will retain only the downstream and upstream node of each segment.

[ ]:

network_grid = network_grid_from_raster(
    grid,
    reducer=JustEndNodes(),
    minimum_channel_threshold=12000.0,
    include=["drainage_area", "topographic__elevation"],
)

plot_graph(network_grid, at="node,link")
plt.title("Nodes and Links");

Reduce to a maximum number of nodes per segment (including end nodes), use AtMostNodes,

[ ]:

network_grid = network_grid_from_raster(
    grid,
    reducer=AtMostNodes(count=4),
    minimum_channel_threshold=12000.0,
    include=["drainage_area", "topographic__elevation"],
)

plot_graph(network_grid, at="node,link")
plt.title("Nodes and Links");

Another option is to specify a spacing at each node. The spacing can be either a scalar, meaning the spacing is the same for every node, or can vary from node to node. An example of the latter is a spacing that depends on upstream area at each node.

[ ]:

spacing = spacing_from_drainage_area(grid.at_node["drainage_area"])

network_grid = network_grid_from_raster(
    grid,
    reducer=SpacingAtLeast(grid.xy_of_node, spacing),
    minimum_channel_threshold=12000.0,
    include=["drainage_area", "topographic__elevation"],
)

plot_graph(network_grid, at="node,link", with_id=False)
plt.title("Nodes and Links");

Generated by nbsphinx from a Jupyter notebook.