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.