landlab

ChiFinder: Perform chi-index analysis for gridded topography

class ChiFinder(grid, reference_concavity=0.5, min_drainage_area=1000000.0, reference_area=1.0, use_true_dx=False, noclobber=True, **kwds)[source]

Bases: landlab.core.model_component.Component

This component calculates chi indices, sensu Perron & Royden, 2013, for a Landlab landscape.

Examples

>>> import numpy as np
>>> from landlab import RasterModelGrid, CLOSED_BOUNDARY
>>> from landlab.components import FlowAccumulator, FastscapeEroder
>>> from landlab.components import ChiFinder
>>> mg = RasterModelGrid((3, 4))
>>> for nodes in (mg.nodes_at_right_edge, mg.nodes_at_bottom_edge,
...               mg.nodes_at_top_edge):
...     mg.status_at_node[nodes] = CLOSED_BOUNDARY
>>> _ = mg.add_field('node', 'topographic__elevation', mg.node_x)
>>> fr = FlowAccumulator(mg, flow_director='D8')
>>> cf = ChiFinder(mg,
...     min_drainage_area=1.,
...     reference_concavity=1.)
>>> fr.run_one_step()
>>> cf.calculate_chi()
>>> mg.at_node['channel__chi_index'].reshape(mg.shape)[1, :]
array([ 0.5,  1. ,  2. ,  0. ])
>>> mg2 = RasterModelGrid((5, 5), xy_spacing=100.)
>>> for nodes in (mg2.nodes_at_right_edge, mg2.nodes_at_bottom_edge,
...               mg2.nodes_at_top_edge):
...     mg2.status_at_node[nodes] = CLOSED_BOUNDARY
>>> _ = mg2.add_zeros('node', 'topographic__elevation')
>>> mg2.at_node['topographic__elevation'][mg2.core_nodes] = mg2.node_x[
...     mg2.core_nodes]/1000.
>>> np.random.seed(0)
>>> mg2.at_node['topographic__elevation'][
...     mg2.core_nodes] += np.random.rand(mg2.number_of_core_nodes)
>>> fr2 = FlowAccumulator(mg2, flow_director='D8')
>>> sp2 = FastscapeEroder(mg2, K_sp=0.01)
>>> cf2 = ChiFinder(
...     mg2,
...     min_drainage_area=0.,
...     reference_concavity=0.5)
>>> for i in range(10):
...     mg2.at_node['topographic__elevation'][mg2.core_nodes] += 10.
...     fr2.run_one_step()
...     sp2.run_one_step(1000.)
>>> fr2.run_one_step()
>>> cf2.calculate_chi()
>>> mg2.at_node['channel__chi_index'].reshape(
...     mg2.shape)  
array([[ 0.        ,  0.        ,  0.        ,  0.        ,  0.        ],
       [ 0.77219416,  1.54438833,  2.63643578,  2.61419437,  0.        ],
       [ 1.09204746,  2.18409492,  1.52214691,  2.61419437,  0.        ],
       [ 0.44582651,  0.89165302,  1.66384718,  2.75589464,  0.        ],
       [ 0.        ,  0.        ,  0.        ,  0.        ,  0.        ]])
>>> cf2.calculate_chi(min_drainage_area=20000., use_true_dx=True,
...                   reference_area=mg2.at_node['drainage_area'].max())
>>> cf2.chi_indices.reshape(mg2.shape)  
array([[   0. ,   0.        ,   0.        ,   0. ,   0. ],
       [   0. , 173.20508076,   0.        ,   0. ,   0. ],
       [   0. ,   0.        , 270.71067812,   0. ,   0. ],
       [   0. , 100.        , 236.60254038,   0. ,   0. ],
       [   0. ,   0.        ,   0.        ,   0. ,   0. ]])
>>> cf2.hillslope_mask.reshape(mg2.shape)
array([[ True,  True,  True,  True,  True],
       [False, False,  True,  True,  True],
       [ True,  True, False,  True,  True],
       [False, False, False,  True,  True],
       [ True,  True,  True,  True,  True]], dtype=bool)
Parameters
  • grid (RasterModelGrid) – A landlab RasterModelGrid.

  • reference_concavity (float) – The reference concavity to use in the calculation.

  • min_drainage_area (float (m**2)) – The drainage area down to which to calculate chi.

  • reference_area (float or None (m**2)) – If None, will default to the mean core cell area on the grid. Else, provide a value to use. Essentially becomes a prefactor on the value of chi.

  • use_true_dx (bool (default False)) – If True, integration to give chi is performed using each value of node spacing along the channel (which can lead to a quantization effect, and is not preferred by Taylor & Royden). If False, the mean value of node spacing along the all channels is assumed everywhere.

  • noclobber (bool (default True)) – Raise an exception if adding an already existing field.

__init__(grid, reference_concavity=0.5, min_drainage_area=1000000.0, reference_area=1.0, use_true_dx=False, noclobber=True, **kwds)[source]
Parameters
  • grid (RasterModelGrid) – A landlab RasterModelGrid.

  • reference_concavity (float) – The reference concavity to use in the calculation.

  • min_drainage_area (float (m**2)) – The drainage area down to which to calculate chi.

  • reference_area (float or None (m**2)) – If None, will default to the mean core cell area on the grid. Else, provide a value to use. Essentially becomes a prefactor on the value of chi.

  • use_true_dx (bool (default False)) – If True, integration to give chi is performed using each value of node spacing along the channel (which can lead to a quantization effect, and is not preferred by Taylor & Royden). If False, the mean value of node spacing along the all channels is assumed everywhere.

  • noclobber (bool (default True)) – Raise an exception if adding an already existing field.

best_fit_chi_elevation_gradient_and_intercept(ch_nodes=None)[source]

Returns least squares best fit for a straight line through a chi plot.

Parameters

ch_nodes (array of ints or None) – Nodes at which to consider chi and elevation values. If None, will use all nodes in grid with area greater than the component min_drainage_area.

Returns

coeffs – A len-2 array containing the m then z0, where z = z0 + m * chi.

Return type

array(gradient, intercept)

Examples

>>> import numpy as np
>>> from landlab import RasterModelGrid, CLOSED_BOUNDARY
>>> from landlab.components import FlowAccumulator
>>> from landlab.components import ChiFinder
>>> mg = RasterModelGrid((3, 4))
>>> for nodes in (mg.nodes_at_right_edge, mg.nodes_at_bottom_edge,
...               mg.nodes_at_top_edge):
...     mg.status_at_node[nodes] = CLOSED_BOUNDARY
>>> z = mg.add_field('node', 'topographic__elevation',
...                  mg.node_x.copy())
>>> z[4:8] = np.array([0.5, 1., 2., 0.])
>>> fr = FlowAccumulator(mg, flow_director='D8')
>>> cf = ChiFinder(
...     mg,
...     min_drainage_area=1.,
...     reference_concavity=1.)
>>> fr.run_one_step()
>>> cf.calculate_chi()
>>> mg.at_node['channel__chi_index'].reshape(mg.shape)[1, :]
array([ 0.5,  1. ,  2. ,  0. ])
>>> coeffs = cf.best_fit_chi_elevation_gradient_and_intercept()
>>> np.allclose(np.array([1., 0.]), coeffs)
True
calculate_chi(**kwds)[source]

This is the main method. Call it to calculate local chi indices at all points with drainage areas greater than min_drainage_area.

This “run” method can optionally take the same parameter set as provided at instantiation. If they are provided, they will override the existing values from instantiation.

Chi of any node without a defined value is reported as 0. These nodes are also identified in the mask retrieved with hillslope_mask.

property chi_indices

Return the array of channel steepness indices.

Nodes not in the channel receive zeros.

create_chi_plot(channel_heads=None, label_axes=True, symbol='kx', plot_line=False, line_symbol='r-')[source]

Plots a “chi plot” (chi vs elevation for points in channel network).

If channel_heads is provided, only the channel nodes downstream of the provided points (and with area > min_drainage_area) will be plotted.

Parameters
  • channel_heads (int, list or array of ints, or None) – Node IDs of channel heads to from which plot downstream.

  • label_axes (bool) – If True, labels the axes as “Chi” and “Elevation (m)”.

  • symbol (str) – A matplotlib-style string for the style to use for the points.

  • plot_line (bool) – If True, will plot a linear best fit line through the data cloud.

  • line_symbol (str) – A matplotlib-style string for the style to use for the line, if plot_line.

property hillslope_mask

Return a boolean array, False where steepness indices exist.

integrate_chi_avg_dx(valid_upstr_order, chi_integrand, chi_array, mean_dx)[source]

Calculates chi at each channel node by summing chi_integrand.

This method assumes a uniform, mean spacing between nodes. Method is deliberately split out for potential cythonization at a later stage.

Parameters
  • valid_upstr_order (array of ints) – nodes in the channel network in upstream order.

  • chi_integrand (array of floats) – The value (A0/A)**concavity, in upstream order.

  • chi_array (array of floats) – Array in which to store chi.

  • mean_dx (float) – The mean node spacing in the network.

Examples

>>> import numpy as np
>>> from landlab import RasterModelGrid, CLOSED_BOUNDARY
>>> from landlab.components import FlowAccumulator
>>> from landlab.components import ChiFinder
>>> mg = RasterModelGrid((5, 4))
>>> for nodes in (mg.nodes_at_right_edge, mg.nodes_at_bottom_edge,
...               mg.nodes_at_top_edge):
...     mg.status_at_node[nodes] = CLOSED_BOUNDARY
>>> z = mg.node_x.copy()
>>> z[[5, 13]] = z[6]  # guard nodes
>>> _ = mg.add_field('node', 'topographic__elevation', z)
>>> fr = FlowAccumulator(mg, flow_director='D8')
>>> cf = ChiFinder(mg)
>>> fr.run_one_step()
>>> ch_nodes = np.array([4, 8, 12, 5, 9, 13, 6, 10, 14])
>>> ch_integrand = 3.*np.ones(9, dtype=float)  # to make calc clearer
>>> chi_array = np.zeros(mg.number_of_nodes, dtype=float)
>>> cf.integrate_chi_avg_dx(ch_nodes, ch_integrand, chi_array, 0.5)
>>> chi_array.reshape(mg.shape)
array([[ 0. ,  0. ,  0. ,  0. ],
       [ 1.5,  3. ,  4.5,  0. ],
       [ 1.5,  3. ,  4.5,  0. ],
       [ 1.5,  3. ,  4.5,  0. ],
       [ 0. ,  0. ,  0. ,  0. ]])
integrate_chi_each_dx(valid_upstr_order, chi_integrand_at_nodes, chi_array)[source]

Calculates chi at each channel node by summing chi_integrand*dx.

This method accounts explicitly for spacing between each node. Method is deliberately split out for potential cythonization at a later stage. Uses a trapezium integration method.

Parameters
  • valid_upstr_order (array of ints) – nodes in the channel network in upstream order.

  • chi_integrand_at_nodes (array of floats) – The value (A0/A)**concavity, in node order.

  • chi_array (array of floats) – Array in which to store chi.

Examples

>>> import numpy as np
>>> from landlab import RasterModelGrid, CLOSED_BOUNDARY
>>> from landlab.components import FlowAccumulator
>>> from landlab.components import ChiFinder
>>> mg = RasterModelGrid((5, 4), xy_spacing=3.)
>>> for nodes in (mg.nodes_at_right_edge, mg.nodes_at_bottom_edge,
...               mg.nodes_at_top_edge):
...     mg.status_at_node[nodes] = CLOSED_BOUNDARY
>>> z = mg.node_x.copy()
>>> z[[5, 13]] = z[6]  # guard nodes
>>> _ = mg.add_field('node', 'topographic__elevation', z)
>>> fr = FlowAccumulator(mg, flow_director='D8')
>>> cf = ChiFinder(mg)
>>> fr.run_one_step()
>>> ch_nodes = np.array([4, 8, 12, 5, 9, 13, 6, 10, 14])
>>> ch_integrand = 2.*np.ones(mg.number_of_nodes,
...                           dtype=float)  # to make calc clearer
>>> chi_array = np.zeros(mg.number_of_nodes, dtype=float)
>>> cf.integrate_chi_each_dx(ch_nodes, ch_integrand, chi_array)
>>> chi_array.reshape(mg.shape)
array([[  0.        ,   0.        ,   0.        ,   0.        ],
       [  0.        ,   6.        ,  14.48528137,   0.        ],
       [  0.        ,   6.        ,  12.        ,   0.        ],
       [  0.        ,   6.        ,  14.48528137,   0.        ],
       [  0.        ,   0.        ,   0.        ,   0.        ]])
>>> from landlab.components import FastscapeEroder
>>> mg2 = RasterModelGrid((5, 5), xy_spacing=100.)
>>> for nodes in (mg2.nodes_at_right_edge, mg2.nodes_at_bottom_edge,
...               mg2.nodes_at_top_edge):
...     mg2.status_at_node[nodes] = CLOSED_BOUNDARY
>>> _ = mg2.add_zeros('node', 'topographic__elevation')
>>> mg2.at_node['topographic__elevation'][mg2.core_nodes] = mg2.node_x[
...     mg2.core_nodes]/1000.
>>> np.random.seed(0)
>>> mg2.at_node['topographic__elevation'][
...     mg2.core_nodes] += np.random.rand(mg2.number_of_core_nodes)
>>> fr2 = FlowAccumulator(mg2, flow_director='D8')
>>> sp2 = FastscapeEroder(mg2, K_sp=0.01)
>>> cf2 = ChiFinder(mg2, min_drainage_area=1., reference_concavity=0.5,
...                 use_true_dx=True)
>>> for i in range(10):
...     mg2.at_node['topographic__elevation'][mg2.core_nodes] += 10.
...     fr2.run_one_step()
...     sp2.run_one_step(1000.)
>>> fr2.run_one_step()
>>> output_array = np.zeros(25, dtype=float)
>>> cf2.integrate_chi_each_dx(mg2.at_node['flow__upstream_node_order'],
...                           np.ones(25, dtype=float),
...                           output_array)
>>> output_array.reshape(mg2.shape)
array([[   0. ,    0. ,    0.        ,    0.        ,    0. ],
       [   0. ,  100. ,  200.        ,  382.84271247,    0. ],
       [   0. ,  100. ,  241.42135624,  341.42135624,    0. ],
       [   0. ,  100. ,  200.        ,  300.        ,    0. ],
       [   0. ,    0. ,    0.        ,    0.        ,    0. ]])
property masked_chi_indices

Returns a masked array version of the ‘channel__chi_index’ field. This enables easier plotting of the values with landlab.imshow_grid_at_node or similar.

Examples

Make a topographic map with an overlay of chi values:

>>> from landlab import imshow_grid_at_node
>>> from landlab import RasterModelGrid, CLOSED_BOUNDARY
>>> from landlab.components import FlowAccumulator, FastscapeEroder
>>> from landlab.components import ChiFinder
>>> mg = RasterModelGrid((5, 5), xy_spacing=100.)
>>> for nodes in (mg.nodes_at_right_edge, mg.nodes_at_bottom_edge,
...               mg.nodes_at_top_edge):
...     mg.status_at_node[nodes] = CLOSED_BOUNDARY
>>> _ = mg.add_zeros('node', 'topographic__elevation')
>>> mg.at_node['topographic__elevation'][mg.core_nodes] = mg.node_x[
...     mg.core_nodes]/1000.
>>> np.random.seed(0)
>>> mg.at_node['topographic__elevation'][
...     mg.core_nodes] += np.random.rand(mg.number_of_core_nodes)
>>> fr = FlowAccumulator(mg, flow_director='D8')
>>> sp = FastscapeEroder(mg, K_sp=0.01)
>>> cf = ChiFinder(mg, min_drainage_area=20000.)
>>> for i in range(10):
...     mg.at_node['topographic__elevation'][mg.core_nodes] += 10.
...     fr.run_one_step()
...     sp.run_one_step(1000.)
>>> fr.run_one_step()
>>> cf.calculate_chi()
>>> imshow_grid_at_node(mg, 'topographic__elevation',
...                     allow_colorbar=False)
>>> imshow_grid_at_node(mg, cf.masked_chi_indices,
...                     color_for_closed=None, cmap='winter')
mean_channel_node_spacing(ch_nodes)[source]

Calculates the mean spacing between all adjacent channel nodes.

Parameters

ch_nodes (array of ints) – The nodes within the defined channel network.

Returns

mean_spacing – The mean spacing between all nodes in the network.

Return type

float (m)

Examples

>>> import numpy as np
>>> from landlab import RasterModelGrid, CLOSED_BOUNDARY
>>> from landlab.components import FlowAccumulator
>>> from landlab.components import ChiFinder
>>> mg = RasterModelGrid((5, 4), xy_spacing=2.)
>>> for nodes in (mg.nodes_at_right_edge, mg.nodes_at_bottom_edge,
...               mg.nodes_at_top_edge):
...     mg.status_at_node[nodes] = CLOSED_BOUNDARY
>>> z = mg.node_x.copy()
>>> z[[5, 13]] = z[6]  # guard nodes
>>> _ = mg.add_field('node', 'topographic__elevation', z)
>>> fr = FlowAccumulator(mg, flow_director='D8')
>>> cf = ChiFinder(mg)
>>> fr.run_one_step()
>>> ch_nodes = np.array([4, 8, 12, 5, 9, 13, 6, 10, 14])
>>> cf.mean_channel_node_spacing(ch_nodes)
2.2761423749153966
nodes_downstream_of_channel_head(channel_head)[source]

Find and return an array with nodes downstream of channel_head.

Parameters

channel_head (int) – Node ID of channel head from which to get downstream nodes.

Examples

>>> import numpy as np
>>> from landlab import RasterModelGrid, CLOSED_BOUNDARY
>>> from landlab.components import FlowAccumulator, ChiFinder
>>> mg = RasterModelGrid((3, 4))
>>> for nodes in (mg.nodes_at_right_edge, mg.nodes_at_bottom_edge,
...               mg.nodes_at_top_edge):
...     mg.status_at_node[nodes] = CLOSED_BOUNDARY
>>> z = mg.add_field('node', 'topographic__elevation',
...                  mg.node_x.copy())
>>> z[4:8] = np.array([0.5, 1., 2., 0.])
>>> fr = FlowAccumulator(mg, flow_director='D8')
>>> fr.run_one_step()
>>> mg.at_node['flow__receiver_node']
array([ 0,  1,  2,  3,  4,  4,  5,  7,  8,  9, 10, 11])
>>> cf = ChiFinder(mg, min_drainage_area=0., reference_concavity=1.)
>>> cf.calculate_chi()
>>> cf.nodes_downstream_of_channel_head(6)
[6, 5, 4]