Gradient calculation functions#

Calculate gradients on a raster grid.

Gradient calculators for raster grids#

calc_grad_at_link(node_values[, out])

Calculate gradients in node_values at links.

calc_grad_across_cell_faces(node_values, ...)

grid.calc_grad_across_cell_faces(node_values, [cell_ids], out=None)

calc_grad_across_cell_corners(node_values, ...)

grid.calc_grad_across_cell_corners(node_values, [cell_ids], out=None)

calc_aspect_at_cell_subtriangles(grid, elevs='topographic__elevation', subtriangle_unit_normals=None, unit='degrees')[source]#

Get tuple of arrays of aspect of each of the eight cell subtriangles.

Aspect is returned as radians clockwise of north, unless input parameter units is set to ‘degrees’.

If subtriangle_unit_normals is provided the aspect will be calculated from these data.

If it is not, it will be derived from elevation data at the nodes, which can either be a string referring to a grid field (default: ‘topographic__elevation’), or an nnodes-long numpy array of the values themselves.

Parameters
  • elevs (str or array (optional)) – Node field name or node array of elevations. If subtriangle_unit_normals is not provided, must be set, but unused otherwise.

  • subtriangle_unit_normals (tuple of 8 (ncels, 3) arrays (optional)) – The unit normal vectors for the eight subtriangles of each cell, if already known. Order is from north of east, counter clockwise to south of east (East North East, North North East, North North West, West North West, West South West, South South West, South South East, East South East).

  • unit ({'degrees', 'radians'}) – Controls the unit that the aspect is returned as.

Returns

each a length num-cells array Len-8 tuple of the aspect of each of the eight cell subtriangles. Aspect is returned as angle clockwise of north. Units are given as radians unless input parameter units is set to ‘degrees’. Order is from north of east, counter clockwise to south of east (East North East, North North East, North North West, West North West, West South West, South South West, South South East, East South East).

Return type

(a_ENE, a_NNE, a_NNW, a_WNW, a_WSW, a_SSW, a_SSE, a_ESE)

Examples

>>> import numpy as np
>>> from landlab import RasterModelGrid
>>> mg = RasterModelGrid((3, 3))
>>> z = np.array([1., 0., 1., 0., 0., 0., 1., 0., 1.])
>>> eight_tris = mg.calc_unit_normals_at_cell_subtriangles(z)
>>> A = mg.calc_aspect_at_cell_subtriangles(z, eight_tris)
>>> A0 = mg.calc_aspect_at_cell_subtriangles(z)
>>> np.allclose(A, A0)
True
>>> type(A) is tuple
True
>>> len(A)
8
>>> len(A[0]) == mg.number_of_cells
True
>>> A0  
(array([ 180.]), array([ 270.]), array([ 90.]), array([ 180.]),
 array([ 0.]), array([ 90.]), array([ 270.]), array([ 0.]))
calc_diff_at_d8(node_values, out=None)[source]#

Calculate differences of node values over links and diagonals.

Calculates the difference in quantity node_values at each link in the grid.

Parameters
  • node_values (ndarray or field name) – Values at grid nodes.

  • out (ndarray, optional) – Buffer to hold the result.

Returns

Differences across links.

Return type

ndarray

Examples

>>> import numpy as np
>>> from landlab import RasterModelGrid
>>> grid = RasterModelGrid((3, 4), xy_spacing=(4, 3))
>>> z = [
...     [60.0, 60.0, 60.0, 60.0],
...     [60.0, 60.0,  0.0,  0.0],
...     [60.0,  0.0,  0.0,  0.0],
... ]
>>> grid.calc_diff_at_d8(z)
array([  0.,   0.,   0.,   0.,   0., -60., -60.,   0., -60.,   0.,   0.,
       -60.,   0.,   0., -60.,   0.,   0.,   0.,   0., -60.,   0., -60.,
       -60., -60.,   0., -60.,   0.,   0.,   0.])
calc_diff_at_diagonal(node_values, out=None)[source]#

Calculate differences of node values over diagonals.

Calculates the difference in quantity node_values at each link in the grid.

Parameters
  • node_values (ndarray or field name) – Values at grid nodes.

  • out (ndarray, optional) – Buffer to hold the result.

Returns

Differences across links.

Return type

ndarray

Examples

>>> import numpy as np
>>> from landlab import RasterModelGrid
>>> grid = RasterModelGrid((3, 4), xy_spacing=(4, 3))
>>> z = [
...     [5.0, 5.0, 5.0, 5.0],
...     [5.0, 5.0, 0.0, 0.0],
...     [5.0, 0.0, 0.0, 0.0],
... ]
>>> grid.calc_diff_at_diagonal(z)
array([ 0.,  0., -5.,  0., -5., -5., -5.,  0., -5.,  0.,  0.,  0.])
calc_grad_across_cell_corners(node_values, *args, **kwds)[source]#

grid.calc_grad_across_cell_corners(node_values, [cell_ids], out=None)

Get gradients to diagonally opposite nodes.

Calculate gradient of the value field provided by node_values to the values at diagonally opposite nodes. The returned gradients are ordered as upper-right, upper-left, lower-left and lower-right.

Parameters
  • node_values (array_like or field name) – Quantity to take the gradient of defined at each node.

  • cell_ids (array_like, optional) – If provided, cell ids to measure gradients. Otherwise, find gradients for all cells.

  • out (array_like, optional) – Alternative output array in which to place the result. Must be of the same shape and buffer length as the expected output.

Returns

Gradients to each diagonal node.

Return type

(N, 4) ndarray

Examples

Create a grid with two cells.

>>> from landlab import RasterModelGrid
>>> grid = RasterModelGrid((3, 4))
>>> x = np.array([1., 0., 0., 1.,
...               0., 0., 1., 1.,
...               3., 3., 3., 3.])

A decrease in quantity to a diagonal node is a negative gradient.

>>> from math import sqrt
>>> grid.calc_grad_across_cell_corners(x) * sqrt(2.)
array([[ 3.,  3.,  1.,  0.],
       [ 2.,  2., -1.,  0.]])
>>> grid = RasterModelGrid((3, 4), xy_spacing=(4, 3))
>>> grid.calc_grad_across_cell_corners(x)
array([[ 0.6,  0.6,  0.2,  0. ],
       [ 0.4,  0.4, -0.2,  0. ]])
calc_grad_across_cell_faces(node_values, *args, **kwds)[source]#

grid.calc_grad_across_cell_faces(node_values, [cell_ids], out=None)

Get gradients across the faces of a cell.

Calculate gradient of the value field provided by node_values across each of the faces of the cells of a grid. The returned gradients are ordered as right, top, left, and bottom.

Note that the returned gradients are masked to exclude neighbor nodes which are closed. Beneath the mask is the value -1.

Parameters
  • node_values (array_like or field name) – Quantity to take the gradient of defined at each node.

  • cell_ids (array_like, optional) – If provided, cell ids to measure gradients. Otherwise, find gradients for all cells.

  • out (array_like, optional) – Alternative output array in which to place the result. Must be of the same shape and buffer length as the expected output.

Returns

Gradients for each face of the cell.

Return type

(N, 4) Masked ndarray

Examples

Create a grid with two cells.

>>> from landlab import RasterModelGrid
>>> grid = RasterModelGrid((3, 4))
>>> x = np.array([0., 0., 0., 0.,
...               0., 0., 1., 1.,
...               3., 3., 3., 3.])

A decrease in quantity across a face is a negative gradient.

>>> grid.calc_grad_across_cell_faces(x) 
masked_array(data =
 [[ 1.  3.  0.  0.]
 [ 0.  2. -1. -1.]],
             mask =
 False,
       fill_value = 1e+20)
>>> grid = RasterModelGrid((3, 4), xy_spacing=(1, 2))
>>> grid.calc_grad_across_cell_faces(x) 
masked_array(data =
 [[ 1.   1.5  0.   0. ]
 [ 0.   1.  -1.  -0.5]],
              mask =
 False,
       fill_value = 1e+20)

grid.calc_grad_along_node_links(node_values, [cell_ids], out=None)

Get gradients along links touching a node.

Calculate gradient of the value field provided by node_values across each of the faces of the nodes of a grid. The returned gradients are ordered as right, top, left, and bottom. All returned values follow our standard sign convention, where a link pointing N or E and increasing in value is positive, a link pointing S or W and increasing in value is negative.

Note that the returned gradients are masked to exclude neighbor nodes which are closed. Beneath the mask is the value -1.

Parameters
  • node_values (array_like or field name) – Quantity to take the gradient of defined at each node.

  • node_ids (array_like, optional) – If provided, node ids to measure gradients. Otherwise, find gradients for all nodes.

  • out (array_like, optional) – Alternative output array in which to place the result. Must be of the same shape and buffer length as the expected output.

Returns

Gradients for each link of the node. Ordering is E,N,W,S.

Return type

(N, 4) Masked ndarray

Examples

Create a grid with nine nodes.

>>> from landlab import RasterModelGrid
>>> grid = RasterModelGrid((3, 3))
>>> x = np.array([0., 0., 0.,
...               0., 1., 2.,
...               2., 2., 2.])

A decrease in quantity across a face is a negative gradient.

>>> grid.calc_grad_along_node_links(x) 
masked_array(data =
 [[-- -- -- --]
 [-- 1.0 -- --]
 [-- -- -- --]
 [1.0 -- -- --]
 [1.0 1.0 1.0 1.0]
 [-- -- 1.0 --]
 [-- -- -- --]
 [-- -- -- 1.0]
 [-- -- -- --]],
             mask =
 [[ True  True  True  True]
 [ True False  True  True]
 [ True  True  True  True]
 [False  True  True  True]
 [False False False False]
 [ True  True False  True]
 [ True  True  True  True]
 [ True  True  True False]
 [ True  True  True  True]],
       fill_value = 1e+20)
>>> grid = RasterModelGrid((3, 3), xy_spacing=(4, 2))
>>> grid.calc_grad_along_node_links(x) 
masked_array(data =
 [[-- -- -- --]
 [-- 0.5 -- --]
 [-- -- -- --]
 [0.25 -- -- --]
 [0.25 0.5 0.25 0.5]
 [-- -- 0.25 --]
 [-- -- -- --]
 [-- -- -- 0.5]
 [-- -- -- --]],
             mask =
 [[ True  True  True  True]
 [ True False  True  True]
 [ True  True  True  True]
 [False  True  True  True]
 [False False False False]
 [ True  True False  True]
 [ True  True  True  True]
 [ True  True  True False]
 [ True  True  True  True]],
       fill_value = 1e+20)
calc_grad_at_d8(grid, node_values, out=None)[source]#

Calculate gradients over all diagonals and links.

Parameters
  • node_values (array_like or field name) – Values at nodes.

  • out (ndarray, optional) – Buffer to hold result. If None, create a new array.

Examples

>>> from landlab import RasterModelGrid
>>> import numpy as np
>>> grid = RasterModelGrid((3, 4), xy_spacing=(4, 3))
>>> z = [
...     [60.0, 60.0, 60.0, 60.0],
...     [60.0, 60.0,  0.0,  0.0],
...     [60.0,  0.0,  0.0,  0.0],
... ]
>>> grid.calc_grad_at_d8(z)
array([  0.,   0.,   0.,   0.,   0., -20., -20.,   0., -15.,   0.,   0.,
       -20.,   0.,   0., -15.,   0.,   0.,   0.,   0., -12.,   0., -12.,
       -12., -12.,   0., -12.,   0.,   0.,   0.])
calc_grad_at_diagonal(grid, node_values, out=None)[source]#

Calculate gradients over all diagonals.

Parameters
  • node_values (array_like or field name) – Values at nodes.

  • out (ndarray, optional) – Buffer to hold result. If None, create a new array.

Examples

>>> from landlab import RasterModelGrid
>>> import numpy as np
>>> grid = RasterModelGrid((3, 4), xy_spacing=(4, 3))
>>> z = [
...     [5.0, 5.0, 5.0, 5.0],
...     [5.0, 5.0, 0.0, 0.0],
...     [5.0, 0.0, 0.0, 0.0],
... ]
>>> grid.calc_grad_at_diagonal(z)
array([ 0.,  0., -1.,  0., -1., -1., -1.,  0., -1.,  0.,  0.,  0.])

Calculate gradients in node_values at links.

Parameters
  • node_values (array_like or field name) – Values at nodes.

  • out (ndarray, optional) – Buffer to hold result. If None, create a new array.

Returns

Gradients of the nodes values for each link.

Return type

ndarray

Examples

>>> from landlab import RasterModelGrid
>>> grid = RasterModelGrid((3, 3))
>>> node_values = [0., 0., 0.,
...                1., 3., 1.,
...                2., 2., 2.]
>>> grid.calc_grad_at_link(node_values)
array([ 0.,  0.,  1.,  3.,  1.,  2., -2.,  1., -1.,  1.,  0.,  0.])
>>> out = np.empty(grid.number_of_links, dtype=float)
>>> rtn = grid.calc_grad_at_link(node_values, out=out)
>>> rtn is out
True
>>> out
array([ 0.,  0.,  1.,  3.,  1.,  2., -2.,  1., -1.,  1.,  0.,  0.])
>>> grid = RasterModelGrid((3, 3), xy_spacing=(2, 1))
>>> grid.calc_grad_at_link(node_values)
array([ 0.,  0.,  1.,  3.,  1.,  1., -1.,  1., -1.,  1.,  0.,  0.])
>>> _ = grid.add_field("elevation", node_values, at="node")
>>> grid.calc_grad_at_link('elevation')
array([ 0.,  0.,  1.,  3.,  1.,  1., -1.,  1., -1.,  1.,  0.,  0.])
calc_grad_at_patch(grid, elevs='topographic__elevation', ignore_closed_nodes=True, subtriangle_unit_normals=None, slope_magnitude=None)[source]#

Calculate the components of the gradient of each raster patch.

Returns the mean gradient of the four possible patch subtriangles, in radians.

If ignore_closed_nodes is True, closed nodes do not affect gradient calculations. If more than one closed node is present in a patch, the patch gradients in both x and y directions are set to zero.

Parameters
  • elevs (str or ndarray, optional) – Field name or array of node values.

  • ignore_closed_nodes (bool) – If True, do not incorporate values at closed nodes into the calc.

  • subtriangle_unit_normals (tuple of 4 (npatches, 3) arrays (optional)) – The unit normal vectors for the four subtriangles of each patch, if already known. Order is TR, TL, BL, BR.

  • slope_magnitude (array with size num_patches (optional)) – The mean slope of each patch, if already known. Units must be the same as provided here!

Returns

gradient_tuple – Len-2 tuple of arrays giving components of gradient in the x and y directions, in the units of radians.

Return type

(x_component_at_patch, y_component_at_patch)

Examples

>>> import numpy as np
>>> from landlab import RasterModelGrid
>>> mg = RasterModelGrid((4, 5))
>>> z = mg.node_y
>>> (x_grad, y_grad) = mg.calc_grad_at_patch(elevs=z)
>>> np.allclose(y_grad, np.pi/4.)
True
>>> np.allclose(x_grad, 0.)
True
>>> z = mg.node_x.copy()
>>> mg.set_closed_boundaries_at_grid_edges(True, True, True, True)
>>> mg.status_at_node[11] = mg.BC_NODE_IS_CLOSED
>>> mg.status_at_node[[9, 2]] = mg.BC_NODE_IS_FIXED_VALUE
>>> z[11] = 100.  # this should get ignored now
>>> z[9] = 2.  # this should be felt by patch 7 only
>>> z[2] = 1.  # should be felt by patches 1 and 2
>>> xgrad, ygrad = mg.calc_grad_at_patch(
...     elevs=z, ignore_closed_nodes=True)
>>> (xgrad.reshape((3, 4)) * 4./np.pi)[1, 1:]
array([ 1.,  1., -1.])
>>> np.allclose(ygrad[1:3], xgrad[1:3])
True
calc_slope_at_cell_subtriangles(grid, elevs='topographic__elevation', subtriangle_unit_normals=None)[source]#

Calculate the slope (positive magnitude of gradient) at each of the eight cell subtriangles.

Parameters
  • elevs (str or ndarray, optional) – Field name or array of node values.

  • subtriangle_unit_normals (tuple of 8 (ncells, 3) arrays (optional)) – The unit normal vectors for the eight subtriangles of each cell, if already known. Order is from north of east, counter clockwise to south of east (East North East, North North East, North North West, West North West, West South West, South South West, South South East, East South East).

Returns

each a length num-cells array Len-8 tuple of the slopes (positive gradient magnitude) of each of the eight cell subtriangles, in radians. Order is from north of east, counter clockwise to south of east (East North East, North North East, North North West, West North West, West South West, South South West, South South East, East South East).

Return type

(s_ENE, s_NNE, s_NNW, s_WNW, s_WSW, s_SSW, s_SSE, s_ESE)

Examples

>>> import numpy as np
>>> from landlab import RasterModelGrid
>>> mg = RasterModelGrid((3, 3))
>>> z = np.array([np.sqrt(3.), 0., 4./3.,
...               0., 0., 0.,
...               1., 0., 1./np.sqrt(3.)])
>>> eight_tris = mg.calc_unit_normals_at_cell_subtriangles(z)
>>> S = mg.calc_slope_at_cell_subtriangles(z, eight_tris)
>>> S0 = mg.calc_slope_at_cell_subtriangles(z)
>>> np.allclose(S, S0)
True
>>> type(S) is tuple
True
>>> len(S)
8
>>> len(S[0]) == mg.number_of_cells
True
>>> np.allclose(S[0], S[1])
True
>>> np.allclose(S[2], S[3])
True
>>> np.allclose(S[4], S[5])
True
>>> np.allclose(S[6], S[7])
True
>>> np.allclose(np.rad2deg(S[0])[0], 30.)
True
>>> np.allclose(np.rad2deg(S[2])[0], 45.)
True
>>> np.allclose(np.rad2deg(S[4])[0], 60.)
True
>>> np.allclose(np.cos(S[6])[0], 3./5.)
True
calc_slope_at_node(grid, elevs='topographic__elevation', method='patch_mean', ignore_closed_nodes=True, return_components=False)[source]#

Array of slopes at nodes, averaged over neighboring patches.

Produces a value for node slope (i.e., mean gradient magnitude) at each node in a manner analogous to a GIS-style slope map. If method==’patch_mean’, it averages the gradient on each of the patches surrounding the node; if method==’Horn’, it returns the resolved slope direction. Directional information can still be returned through use of the return_components keyword. All values are returned in radians, including the components; take the tan to recover the rise/run.

Note that under these definitions, it is not always true that:

mag, cmp = mg.calc_slope_at_node(z)
mag**2 == cmp[0]**2 + cmp[1]**2  # only if method=='Horn'

If ignore_closed_nodes is False, all proximal elevation values will be used in the calculation. If True, only unclosed nodes are used.

This is a verion of this code specialized for a raster. It subdivides the four square patches around each node into subtriangles, in order to ensure more correct solutions that incorporate equally weighted information from all surrounding nodes on rough surfaces.

Parameters
  • elevs (str or ndarray, optional) – Field name or array of node values.

  • method ({'patch_mean', 'Horn'}) – Controls the slope algorithm. Current options are ‘patch_mean’, which takes the mean slope of each pf the four neighboring square patches, and ‘Horn’, which is the standard ArcGIS slope algorithm. These produce very similar solutions; the Horn method gives a vector mean and the patch_mean gives a scalar mean.

  • ignore_closed_nodes (bool) – If True, do not incorporate values at closed nodes into the calc.

  • return_components (bool) – If True, return a tuple, (array_of_magnitude, (array_of_slope_x_radians, array_of_slope_y_radians)). If false, return an array of floats of the slope magnitude.

Returns

If return_components, returns (array_of_magnitude, (array_of_slope_x_radians, array_of_slope_y_radians)). If not return_components, returns an array of slope magnitudes.

Return type

float array or length-2 tuple of float arrays

Examples

>>> import numpy as np
>>> from landlab import RadialModelGrid, RasterModelGrid
>>> mg = RasterModelGrid((5, 5))
>>> z = mg.node_x
>>> slopes = mg.calc_slope_at_node(elevs=z)
>>> np.allclose(slopes, np.pi / 4.)
True
>>> mg = RasterModelGrid((4, 5), xy_spacing=2.)
>>> z = - mg.node_y
>>> slope_mag, cmp = mg.calc_slope_at_node(elevs=z,
...                                        return_components=True)
>>> np.allclose(slope_mag, np.pi / 4.)
True
>>> np.allclose(cmp[0], 0.)
True
>>> np.allclose(cmp[1], - np.pi / 4.)
True
>>> mg = RasterModelGrid((4, 4))
>>> z = mg.node_x ** 2 + mg.node_y ** 2
>>> slopes, cmp = mg.calc_slope_at_node(z, return_components=True)
>>> slopes
array([ 0.95531662,  1.10991779,  1.32082849,  1.37713803,  1.10991779,
        1.20591837,  1.3454815 ,  1.38904403,  1.32082849,  1.3454815 ,
        1.39288142,  1.41562833,  1.37713803,  1.38904403,  1.41562833,
        1.43030663])
>>> np.allclose(cmp[0].reshape((4, 4))[:, 0],
...             cmp[1].reshape((4, 4))[0, :])  # test radial symmetry
True
calc_slope_at_patch(grid, elevs='topographic__elevation', ignore_closed_nodes=True, subtriangle_unit_normals=None)[source]#

Calculate the slope (positive magnitude of gradient) at raster patches.

Returns the mean of the slopes of the four possible patch subtriangles.

If ignore_closed_nodes is True, closed nodes do not affect slope calculations. If more than one closed node is present in a patch, the patch slope is set to zero.

Parameters
  • elevs (str or ndarray, optional) – Field name or array of node values.

  • ignore_closed_nodes (bool) – If True, do not incorporate values at closed nodes into the calc.

  • subtriangle_unit_normals (tuple of 4 (npatches, 3) arrays (optional)) – The unit normal vectors for the four subtriangles of each patch, if already known. Order is TR, TL, BL, BR.

Returns

slopes_at_patch – The slope (positive gradient magnitude) of each patch, in radians.

Return type

n_patches-long array

Examples

>>> import numpy as np
>>> from landlab import RasterModelGrid
>>> mg = RasterModelGrid((4, 5))
>>> z = mg.node_x
>>> S = mg.calc_slope_at_patch(elevs=z)
>>> S.size == mg.number_of_patches
True
>>> np.allclose(S, np.pi/4.)
True
>>> z = mg.node_y**2
>>> mg.calc_slope_at_patch(elevs=z).reshape((3, 4))
array([[ 0.78539816,  0.78539816,  0.78539816,  0.78539816],
       [ 1.24904577,  1.24904577,  1.24904577,  1.24904577],
       [ 1.37340077,  1.37340077,  1.37340077,  1.37340077]])
>>> z = mg.node_x.copy()
>>> mg.set_closed_boundaries_at_grid_edges(True, True, True, True)
>>> mg.status_at_node[11] = mg.BC_NODE_IS_CLOSED
>>> mg.status_at_node[9] = mg.BC_NODE_IS_FIXED_VALUE
>>> z[11] = 100.  # this should get ignored now
>>> z[9] = 2.  # this should be felt by patch 7 only
>>> mg.calc_slope_at_patch(elevs=z, ignore_closed_nodes=True).reshape(
...     (3, 4)) * 4./np.pi
array([[ 0.,  0.,  0.,  0.],
       [ 0.,  1.,  1.,  1.],
       [ 0.,  0.,  0.,  0.]])
calc_unit_normals_at_cell_subtriangles(grid, elevs='topographic__elevation')[source]#

Calculate unit normals on a cell.

Calculate the eight unit normal vectors <a, b, c> to the eight subtriangles of a four-cornered (raster) cell.

Parameters

elevs (str or ndarray, optional) – Field name or array of node values.

Returns

each a num-cells x length-3 array Len-8 tuple of the eight unit normal vectors <a, b, c> for the eight subtriangles in the cell. Order is from north of east, counter clockwise to south of east (East North East, North North East, North North West, West North West, West South West, South South West, South South East, East South East).

Return type

(n_ENE, n_NNE, n_NNW, n_WNW, n_WSW, n_SSW, n_SSE, n_ESE)

Examples

>>> import numpy as np
>>> from landlab import RasterModelGrid
>>> mg = RasterModelGrid((3, 3))
>>> z = mg.node_x ** 2
>>> eight_tris = mg.calc_unit_normals_at_cell_subtriangles(z)
>>> type(eight_tris) is tuple
True
>>> len(eight_tris)
8
>>> eight_tris[0].shape == (mg.number_of_cells, 3)
True
>>> eight_tris 
(array([[-0.9486833 ,  0.        ,  0.31622777]]),
 array([[-0.9486833 ,  0.        ,  0.31622777]]),
 array([[-0.70710678,  0.        ,  0.70710678]]),
 array([[-0.70710678,  0.        ,  0.70710678]]),
 array([[-0.70710678,  0.        ,  0.70710678]]),
 array([[-0.70710678,  0.        ,  0.70710678]]),
 array([[-0.9486833 ,  0.        ,  0.31622777]]),
 array([[-0.9486833 ,  0.        ,  0.31622777]]))
calc_unit_normals_at_patch_subtriangles(grid, elevs='topographic__elevation')[source]#

Calculate unit normals on a patch.

Calculate the four unit normal vectors <a, b, c> to the four possible subtriangles of a four-cornered (raster) patch.

Parameters

elevs (str or ndarray, optional) – Field name or array of node values.

Returns

(n_TR, n_TL, n_BL, n_BR) – Len-4 tuple of the four unit normal vectors <a, b, c> for the four possible subtriangles in the patch. Order is (topright, topleft, bottomleft, bottomright).

Return type

each a num-patches x length-3 array

Examples

>>> import numpy as np
>>> from landlab import RasterModelGrid
>>> mg = RasterModelGrid((4, 5))
>>> z = mg.node_x ** 2
>>> four_tris = mg.calc_unit_normals_at_patch_subtriangles(z)
>>> type(four_tris) is tuple
True
>>> len(four_tris)
4
>>> np.allclose(four_tris[0], four_tris[1])
True
>>> np.allclose(four_tris[2], four_tris[3])
True
>>> np.allclose(four_tris[0], four_tris[2])
True
>>> np.allclose(np.square(four_tris[0]).sum(axis=1), 1.)
True
>>> four_tris[0]
array([[-0.70710678,  0.        ,  0.70710678],
       [-0.9486833 ,  0.        ,  0.31622777],
       [-0.98058068,  0.        ,  0.19611614],
       [-0.98994949,  0.        ,  0.14142136],
       [-0.70710678,  0.        ,  0.70710678],
       [-0.9486833 ,  0.        ,  0.31622777],
       [-0.98058068,  0.        ,  0.19611614],
       [-0.98994949,  0.        ,  0.14142136],
       [-0.70710678,  0.        ,  0.70710678],
       [-0.9486833 ,  0.        ,  0.31622777],
       [-0.98058068,  0.        ,  0.19611614],
       [-0.98994949,  0.        ,  0.14142136]])