landlab.utils.watershed module#
Functions to work with watersheds of model grids.
- get_watershed_mask(grid, outlet_id)[source]#
Get the watershed of an outlet returned as a boolean array.
- Parameters:
grid (RasterModelGrid) – A landlab RasterModelGrid.
outlet_id (integer) – The id of the outlet node.
- Returns:
watershed_mask – True elements of this array correspond to nodes with flow that is received by the outlet. The length of the array is equal to the grid number of nodes.
- Return type:
boolean ndarray
Examples
>>> import numpy as np >>> from landlab import RasterModelGrid >>> from landlab.components import FlowAccumulator >>> from landlab.utils import get_watershed_mask
>>> rmg = RasterModelGrid((7, 7)) >>> z = np.array([ ... -9999., -9999., -9999., -9999., -9999., -9999., -9999., ... -9999., 26., 0., 30., 32., 34., -9999., ... -9999., 28., 1., 25., 28., 32., -9999., ... -9999., 30., 3., 3., 11., 34., -9999., ... -9999., 32., 11., 25., 18., 38., -9999., ... -9999., 34., 32., 34., 36., 40., -9999., ... -9999., -9999., -9999., -9999., -9999., -9999., -9999.]) >>> rmg.at_node['topographic__elevation'] = z
Only the bottom boundary is set to open. >>> rmg.set_closed_boundaries_at_grid_edges(True, True, True, False) >>> rmg.set_fixed_value_boundaries_at_grid_edges(False, False, False, True)
Route flow. >>> fr = FlowAccumulator(rmg, flow_director=’D8’) >>> fr.run_one_step()
>>> get_watershed_mask(rmg, 2).reshape(rmg.shape) array([[False, False, True, False, False, False, False], [False, False, True, False, False, False, False], [False, True, True, True, True, True, False], [False, True, True, True, True, True, False], [False, True, True, True, True, True, False], [False, True, True, True, True, True, False], [False, False, False, False, False, False, False]], dtype=bool)
- get_watershed_masks(grid)[source]#
Assign the watershed outlet id to all nodes in the grid.
- Parameters:
grid (RasterModelGrid) – A landlab RasterModelGrid.
- Returns:
watershed_masks – The length of the array is equal to the grid number of nodes. Values of this array are the watershed ids. The value of a watershed id is the node id of the watershed outlet.
- Return type:
integer ndarray
Examples
>>> import numpy as np >>> from landlab import RasterModelGrid >>> from landlab.components import FlowAccumulator >>> from landlab.utils import get_watershed_masks
Create a grid with a node spacing of 200 meter. >>> rmg = RasterModelGrid((7, 7), xy_spacing=200) >>> z = np.array([ … -9999., -9999., -9999., -9999., -9999., -9999., -9999., … -9999., 26., 0., 26., 30., 34., -9999., … -9999., 28., 1., 28., 5., 32., -9999., … -9999., 30., 3., 30., 10., 34., -9999., … -9999., 32., 11., 32., 15., 38., -9999., … -9999., 34., 32., 34., 36., 40., -9999., … -9999., -9999., -9999., -9999., -9999., -9999., -9999.]) >>> rmg.at_node[‘topographic__elevation’] = z >>> rmg.set_closed_boundaries_at_grid_edges(True, True, True, False)
Route flow.
>>> fr = FlowAccumulator(rmg, flow_director='D8') >>> fr.run_one_step()
Assign mask.
>>> mask = get_watershed_masks(rmg) >>> mask.reshape(rmg.shape) array([[ 0, 1, 2, 3, 4, 5, 6], [ 7, 1, 2, 3, 4, 5, 13], [14, 2, 2, 2, 18, 18, 20], [21, 2, 2, 2, 18, 18, 27], [28, 2, 2, 2, 18, 18, 34], [35, 2, 2, 2, 18, 18, 41], [42, 43, 44, 45, 46, 47, 48]])
- get_watershed_masks_with_area_threshold(grid, critical_area)[source]#
Get masks of all of the watersheds with a minimum drainage area size.
- Parameters:
grid (RasterModelGrid) – A landlab RasterModelGrid.
critical_area (integer or float) – The minimum drainage area of the watersheds to identify.
- Returns:
watershed_masks – The length of the array is equal to the grid number of nodes. Values of this array are the watershed ids. The value of a watershed id is the node id of the watershed outlet. Nodes with a value of -1 have only downstream nodes with drainage areas below critical_area.
- Return type:
integer ndarray
Examples
>>> import numpy as np >>> from landlab import RasterModelGrid >>> from landlab.components import FlowAccumulator >>> from landlab.utils import get_watershed_masks_with_area_threshold
Create a grid with a node spacing of 200 meter.
>>> rmg = RasterModelGrid((7, 7), xy_spacing=200) >>> z = np.array([ ... -9999., -9999., -9999., -9999., -9999., -9999., -9999., ... -9999., 26., 0., 26., 30., 34., -9999., ... -9999., 28., 1., 28., 5., 32., -9999., ... -9999., 30., 3., 30., 10., 34., -9999., ... -9999., 32., 11., 32., 15., 38., -9999., ... -9999., 34., 32., 34., 36., 40., -9999., ... -9999., -9999., -9999., -9999., -9999., -9999., -9999.]) >>> rmg.at_node['topographic__elevation'] = z >>> rmg.set_closed_boundaries_at_grid_edges(True, True, True, False)
Route flow.
>>> fr = FlowAccumulator(rmg, flow_director='D8') >>> fr.run_one_step()
Get the masks of watersheds greater than or equal to 80,000 square-meters.
>>> critical_area = 80000 >>> mask = get_watershed_masks_with_area_threshold(rmg, critical_area)
Verify that all mask null nodes have a drainage area below critical area.
>>> null_nodes = np.where(mask == -1)[0] >>> A = rmg.at_node['drainage_area'][null_nodes] >>> below_critical_area_nodes = A < critical_area >>> np.all(below_critical_area_nodes) True
- get_watershed_nodes(grid, outlet_id)[source]#
Get the watershed of an outlet returned as a list of nodes.
- Parameters:
grid (RasterModelGrid) – A landlab RasterModelGrid.
outlet_id (integer) – The id of the outlet node.
- Returns:
watershed_nodes – The ids of the nodes that flow to the node with the id, outlet_id.
- Return type:
integer ndarray
Examples
>>> import numpy as np >>> from landlab import RasterModelGrid >>> from landlab.components import FlowAccumulator >>> from landlab.utils import get_watershed_nodes
>>> rmg = RasterModelGrid((7, 7)) >>> z = np.array([ ... -9999., -9999., -9999., -9999., -9999., -9999., -9999., ... -9999., 26., 0., 30., 32., 34., -9999., ... -9999., 28., 1., 25., 28., 32., -9999., ... -9999., 30., 3., 3., 11., 34., -9999., ... -9999., 32., 11., 25., 18., 38., -9999., ... -9999., 34., 32., 34., 36., 40., -9999., ... -9999., -9999., -9999., -9999., -9999., -9999., -9999.])
>>> rmg.at_node['topographic__elevation'] = z >>> rmg.set_watershed_boundary_condition_outlet_id(2, z, ... nodata_value=-9999.)
Route flow. >>> fr = FlowAccumulator(rmg, flow_director=’D8’) >>> fr.run_one_step()
Get the nodes of two watersheds. >>> mainstem_watershed_nodes = get_watershed_nodes(rmg, 2) >>> tributary_watershed_nodes = get_watershed_nodes(rmg, 24)
Given the watershed boundary conditions, the number of mainstem watershed nodes should be equal to the number of core nodes plus the outlet node. >>> len(mainstem_watershed_nodes) == rmg.number_of_core_nodes + 1 True
- get_watershed_outlet(grid, source_node_id)[source]#
Get the downstream-most node (the outlet) of the source node.
- Parameters:
grid (RasterModelGrid) – A landlab RasterModelGrid.
source_node_id (integer) – The id of the node in which to identify its outlet.
- Returns:
outlet_node – The id of the node that is the downstream-most node (the outlet) of the source node.
- Return type:
integer
Examples
>>> import numpy as np >>> from landlab import RasterModelGrid >>> from landlab.components import FlowAccumulator >>> from landlab.utils import get_watershed_outlet
>>> rmg = RasterModelGrid((7, 7)) >>> z = np.array([ ... -9999., -9999., -9999., -9999., -9999., -9999., -9999., ... -9999., 26., 0., 30., 32., 34., -9999., ... -9999., 28., 1., 25., 28., 32., -9999., ... -9999., 30., 3., 3., 11., 34., -9999., ... -9999., 32., 11., 25., 18., 38., -9999., ... -9999., 34., 32., 34., 36., 40., -9999., ... -9999., -9999., -9999., -9999., -9999., -9999., -9999.])
>>> rmg.at_node['topographic__elevation'] = z >>> imposed_outlet = 2 >>> rmg.set_watershed_boundary_condition_outlet_id(imposed_outlet, z, ... nodata_value=-9999.)
Route flow. >>> fr = FlowAccumulator(rmg, flow_director=’D8’) >>> fr.run_one_step()
Get the grid watershed outlet. >>> determined_outlet = get_watershed_outlet(rmg, 40) >>> determined_outlet == imposed_outlet True