landlab.graph.sort package#

Subpackages#

Submodules#

landlab.graph.sort.sort module#

Sort the elements of a graph.

This module provides functions that sort the elements of a graph structure.

argsort_spokes_at_hub(spokes_at_hub, xy_of_hub, xy_of_spokes)[source]#
argsort_spokes_at_hub_on_graph(graph, spoke=None, at='node')[source]#

Order spokes clockwise around spokes.

Parameters:
  • graph (Graph-like) – A landlab graph.

  • spoke (str) – Name of graph elements that make the spokes.

  • at (str) – Namve of graph elements that make the hubs.

Returns:

Spokes ordered around each hub.

Return type:

ndarray of int, shape (n_spokes, n_hubs)

Examples

>>> import numpy as np
>>> from landlab.graph import UniformRectilinearGraph
>>> from landlab.graph.sort.sort import argsort_spokes_at_hub_on_graph
>>> graph = UniformRectilinearGraph((3, 3))
>>> argsort_spokes_at_hub_on_graph(graph, "link", at="node")
array([[ 0,  1,  2,  3],
       [ 4,  5,  6,  7],
       [ 9, 10,  8, 11],
       [12, 13, 15, 14],
       [16, 17, 18, 19],
       [21, 22, 23, 20],
       [24, 27, 25, 26],
       [28, 30, 31, 29],
       [34, 35, 32, 33]])
calc_angle_of_spoke(spokes_at_hub, xy_of_hub, xy_of_spoke, badval=None)[source]#
calc_angle_of_spoke_on_graph(graph, spoke=None, at='node', badval=None)[source]#

Calculate angles spokes make with a hub.

Parameters:
  • graph (Graph-like) – A landlab graph.

  • spoke (str) – Name of graph elements that make the spokes.

  • at (str) – Namve of graph elements that make the hubs.

  • badval (float or iterable of float, optional) – Value to insert for missing spokes. If an iterable, use items as bad values to use for different spokes.

Returns:

Angle that spoke elements make with each hub element.

Return type:

ndarray of float, shape (n_spokes, n_hubs)

Examples

>>> import numpy as np
>>> from landlab.graph import UniformRectilinearGraph
>>> from landlab.graph.sort.sort import calc_angle_of_spoke_on_graph
>>> graph = UniformRectilinearGraph((3, 3))
>>> np.rad2deg(
...     calc_angle_of_spoke_on_graph(graph, "link", at="node", badval=np.nan)
... )
array([[   0.,   90.,   nan,   nan],
       [   0.,   90.,  180.,   nan],
       [  nan,   90.,  180.,   nan],
       [   0.,   90.,   nan,  270.],
       [   0.,   90.,  180.,  270.],
       [  nan,   90.,  180.,  270.],
       [   0.,   nan,   nan,  270.],
       [   0.,   nan,  180.,  270.],
       [  nan,   nan,  180.,  270.]])
reindex_by_xy(graph)[source]#
reindex_nodes_by_xy(graph)[source]#
reindex_patches_by_xy(graph)[source]#
remap(src, mapping, out=None, inplace=False)[source]#

Remap elements in an id array.

Parameters:
  • src (ndarray of int) – Initial array of ids.

  • mapping (ndarray of int) – Mapping of ids.

  • out (ndarray of int, optional) – Buffer into which to place output.

  • inplace (bool, optional) – Mapping of values will include inplace.

Returns:

Array of remapped values.

Return type:

ndarray of int

Examples

>>> from landlab.graph.sort.sort import remap
>>> import numpy as np
>>> src = np.array([1, 2, 3, 4])
>>> mapping = np.array([-1, 10, 20, 30, 40])
>>> remap(src, mapping)
array([10, 20, 30, 40])
reverse_one_to_many(ids, min_counts=0)[source]#

Reverse a one-to-many mapping.

Parameters:

ids (ndarray of int, shape (M, N)) – Array of identifier mapping.

Returns:

Array of the reverse mapping.

Return type:

ndarray of int, shape (m, n)

Examples

>>> from landlab.graph.sort.sort import reverse_one_to_many
>>> ids = np.array([[1, 2, 3], [-1, -1, -1], [2, 3, -1]], dtype=int)
>>> reverse_one_to_many(ids)
array([[-1, -1],
       [ 0, -1],
       [ 0,  2],
       [ 0,  2]])
reverse_one_to_one(ids, minlength=None)[source]#

Reverse a one-to-one mapping.

Parameters:
  • ids (ndarray of int, shape (N, )) – Array of identifier mapping.

  • minlength (int, optional) – A minimum number of identifiers for the output array.

Returns:

Array of the reverse mapping.

Return type:

ndarray of int, shape (n, )

Examples

>>> from landlab.graph.sort.sort import reverse_one_to_one
>>> ids = np.array([-1, -1, 6, 3, -1, 2, 4, 1, -1, 5, 0], dtype=int)
>>> reverse_one_to_one(ids)
array([10,  7,  5,  3,  6,  9,  2])
sort_graph(nodes, links=None, patches=None)[source]#

Sort elements of a graph by x, then y.

Parameters:
  • nodes (tuple of ndarray) – Coordinate of nodes as (y, x).

  • links (ndarray of int, shape (n_links, 2), optional) – Indices into nodes array of link tail then head.

  • patches (array_like of array_like, optional) – Indices into links array of the links that define each patch.

Returns:

Sorted nodes, links, and patches.

Return type:

(nodes, links, patches)

Examples

o—o—o | | / | o—o—o

>>> from landlab.graph.sort import sort_graph
>>> import numpy as np
>>> x = np.array([1.0, 2.0, 2.0, 0.0, 1.0, 0.0])
>>> y = np.array([0.0, 0.0, 1.0, 0.0, 1.0, 1.0])

Sort a graph with just points - no links or patches.

>>> _ = sort_graph((y, x))
>>> y
array([ 0.,  0.,  0.,  1.,  1.,  1.])
>>> x
array([ 0.,  1.,  2.,  0.,  1.,  2.])

Sort the points and links of a graph.

>>> x = np.array([1.0, 2.0, 2.0, 0.0, 1.0, 0.0])
>>> y = np.array([0.0, 0.0, 1.0, 0.0, 1.0, 1.0])
>>> links = np.array(
...     [[3, 0], [0, 4], [4, 5], [5, 3], [0, 1], [1, 2], [2, 0], [2, 4]]
... )
>>> _ = sort_graph((y, x), links)
>>> links
array([[0, 1], [1, 2], [3, 0], [1, 4], [5, 1], [2, 5], [4, 3], [5, 4]])

Sort the points, links, and patches of a graph.

>>> x = np.array([1.0, 2.0, 2.0, 0.0, 1.0, 0.0])
>>> y = np.array([0.0, 0.0, 1.0, 0.0, 1.0, 1.0])
>>> links = np.array(
...     [[3, 0], [0, 4], [4, 5], [5, 3], [0, 1], [1, 2], [2, 0], [2, 4]]
... )
>>> patches = (np.array([1, 6, 7, 4, 5, 6, 0, 1, 2, 3]), np.array([0, 3, 6, 10]))
>>> _ = sort_graph((y, x), links, patches)
>>> patches[0]
array([1, 5, 4, 0, 3, 6, 2, 3, 4, 7])
>>> patches[1]
array([ 0,  3,  7, 10])

Sort links by their midpoint.

Parameters:
  • nodes_at_link (ndarray of int, shape (n_links, 2)) – Node for each link tail and head.

  • nodes (tuple of ndarray of float) – Node coordinates.

  • midpoint_of_link (ndarray of float, shape (n_links, 2), optional) – Buffer to store the link midpoints that were used for sorting.

Returns:

Array of indices that sort the links.

Return type:

ndarray of int

Examples

>>> from landlab.graph.sort import sort_nodes
>>> import numpy as np
>>> nodes = np.array([[0, 0, 0, 1, 1, 1], [0, 1, 2, 0, 1, 2]])
>>> links = np.array([[0, 1], [0, 3], [1, 2], [1, 4], [2, 5], [3, 4], [4, 5]])
>>> sort_links(links, nodes)
array([0, 2, 1, 3, 4, 5, 6])
>>> links
array([[0, 1],
       [1, 2],
       [0, 3],
       [1, 4],
       [2, 5],
       [3, 4],
       [4, 5]])

Reorder links around a patch to be counterclockwise.

Examples

>>> from landlab.graph.sort.sort import sort_links_at_patch
>>> import numpy as np
>>> xy_of_node = np.array(
...     [
...         [0.0, 0.0],
...         [0.0, 1.0],
...         [1.0, 1.0],
...         [1.0, 0.0],
...     ]
... )
>>> nodes_at_link = np.array([[0, 1], [1, 2], [2, 3], [3, 0]])
>>> links_at_patch = np.array([[0, 1, 3, 2]])
>>> sort_links_at_patch(links_at_patch, nodes_at_link, xy_of_node)
>>> links_at_patch
array([[2, 1, 0, 3]])
>>> xy_of_node = np.array(
...     [
...         [0.0, 0.0],
...         [0.0, 1.0],
...         [1.0, 1.0],
...         [1.0, 0.0],
...         [2.0, 0.0],
...     ]
... )
>>> nodes_at_link = np.array([[0, 1], [1, 2], [2, 3], [3, 0], [2, 4], [3, 4]])
>>> links_at_patch = np.array([[0, 1, 3, 2], [2, 4, 5, -1]])
>>> sort_links_at_patch(links_at_patch, nodes_at_link, xy_of_node)
>>> links_at_patch
array([[ 2,  1,  0,  3],
       [ 4,  2,  5, -1]])
>>> xy_of_node = np.array(
...     [
...         [0.0, 0.0],
...         [0.0, 1.0],
...         [1.0, 1.0],
...         [1.0, 0.0],
...         [2.0, 1.0],
...     ]
... )
>>> nodes_at_link = np.array([[0, 1], [1, 2], [2, 3], [3, 0], [2, 4], [3, 4]])
>>> links_at_patch = np.array([[0, 1, 3, 2], [2, -1, 4, 5]])
>>> sort_links_at_patch(links_at_patch, nodes_at_link, xy_of_node)
>>> links_at_patch
array([[ 2,  1,  0,  3],
       [ 4,  2,  5, -1]])
sort_nodes(nodes)[source]#

Sort nodes based on their position.

Parameters:

nodes (tuple of ndarray of float) – Coordinates of nodes as (y, x).

Returns:

Array of indices that sort the nodes.

Return type:

ndarray of int

Examples

>>> from landlab.graph.sort import sort_nodes
>>> import numpy as np
>>> x = np.array([0.0, 1.0, 2.0])
>>> y = np.array([0.5, 0.0, 1.0])
>>> sort_nodes((y, x))
array([1, 0, 2])
>>> x
array([ 1.,  0.,  2.])
>>> y
array([ 0. ,  0.5,  1. ])
sort_patches(links_at_patch, offset_to_patch, xy_of_link)[source]#

Sort patches by their centroid.

Parameters:
  • links_at_patch (ndarray of int) – Links that define patches.

  • offset_to_patch (ndarray of int) – Offsets into links_at_patch for each patch.

  • xy_of_link (ndarray of float, shape (n_links, 2)) – Midpoint coordinates for each link.

Examples

>>> from landlab.graph.sort import sort_patches
>>> import numpy as np
>>> links_at_patch = np.array([0, 1, 2, 3, 2, 4])
>>> offset_to_patch = np.array([0, 3, 6])
>>> xy_of_link = np.array(
...     [[0.0, 0.5], [0.5, 1.0], [0.5, 0.5], [0.5, 0.0], [1.0, 0.5]]
... )
>>> sort_patches(links_at_patch, offset_to_patch, xy_of_link)
array([1, 0])
>>> links_at_patch
array([3, 2, 4, 0, 1, 2])
>>> offset_to_patch
array([0, 3, 6])
sort_spokes_at_hub(spokes_at_hub, xy_of_hub, xy_of_spokes, inplace=False)[source]#
sort_spokes_at_hub_on_graph(graph, spoke=None, at='node', inplace=False)[source]#

Order spokes of a graph clockwise around spokes.

Parameters:
  • graph (Graph-like) – A landlab graph.

  • spoke (str) – Name of graph elements that make the spokes.

  • at ({'node', 'corner', 'link', 'face', 'patch', 'cell'}) – Namve of graph elements that make the hubs.

Returns:

Spokes ordered around each hub.

Return type:

ndarray of int, shape (n_spokes, n_hubs)

Examples

>>> import numpy as np
>>> from landlab.graph import UniformRectilinearGraph, TriGraph
>>> from landlab.graph.sort.sort import sort_spokes_at_hub_on_graph
>>> graph = UniformRectilinearGraph((3, 3))
>>> sort_spokes_at_hub_on_graph(graph, "link", at="node")
array([[ 0,  2, -1, -1],
       [ 1,  3,  0, -1],
       [ 4,  1, -1, -1],
       [ 5,  7,  2, -1],
       [ 6,  8,  5,  3],
       [ 9,  6,  4, -1],
       [10,  7, -1, -1],
       [11, 10,  8, -1],
       [11,  9, -1, -1]])
>>> graph = TriGraph((3, 3), node_layout="hex", sort=True)
>>> sort_spokes_at_hub_on_graph(graph, "patch", at="node")
array([[ 0,  2, -1, -1, -1, -1],
       [ 1,  3,  0, -1, -1, -1],
       [ 4,  1, -1, -1, -1, -1],
       [ 5,  2, -1, -1, -1, -1],
       [ 6,  8,  5,  2,  0,  3],
       [ 7,  9,  6,  3,  1,  4],
       [ 7,  4, -1, -1, -1, -1],
       [ 5,  8, -1, -1, -1, -1],
       [ 8,  6,  9, -1, -1, -1],
       [ 9,  7, -1, -1, -1, -1]])

Module contents#

reindex_by_xy(graph)[source]#
sort_graph(nodes, links=None, patches=None)[source]#

Sort elements of a graph by x, then y.

Parameters:
  • nodes (tuple of ndarray) – Coordinate of nodes as (y, x).

  • links (ndarray of int, shape (n_links, 2), optional) – Indices into nodes array of link tail then head.

  • patches (array_like of array_like, optional) – Indices into links array of the links that define each patch.

Returns:

Sorted nodes, links, and patches.

Return type:

(nodes, links, patches)

Examples

o—o—o | | / | o—o—o

>>> from landlab.graph.sort import sort_graph
>>> import numpy as np
>>> x = np.array([1.0, 2.0, 2.0, 0.0, 1.0, 0.0])
>>> y = np.array([0.0, 0.0, 1.0, 0.0, 1.0, 1.0])

Sort a graph with just points - no links or patches.

>>> _ = sort_graph((y, x))
>>> y
array([ 0.,  0.,  0.,  1.,  1.,  1.])
>>> x
array([ 0.,  1.,  2.,  0.,  1.,  2.])

Sort the points and links of a graph.

>>> x = np.array([1.0, 2.0, 2.0, 0.0, 1.0, 0.0])
>>> y = np.array([0.0, 0.0, 1.0, 0.0, 1.0, 1.0])
>>> links = np.array(
...     [[3, 0], [0, 4], [4, 5], [5, 3], [0, 1], [1, 2], [2, 0], [2, 4]]
... )
>>> _ = sort_graph((y, x), links)
>>> links
array([[0, 1], [1, 2], [3, 0], [1, 4], [5, 1], [2, 5], [4, 3], [5, 4]])

Sort the points, links, and patches of a graph.

>>> x = np.array([1.0, 2.0, 2.0, 0.0, 1.0, 0.0])
>>> y = np.array([0.0, 0.0, 1.0, 0.0, 1.0, 1.0])
>>> links = np.array(
...     [[3, 0], [0, 4], [4, 5], [5, 3], [0, 1], [1, 2], [2, 0], [2, 4]]
... )
>>> patches = (np.array([1, 6, 7, 4, 5, 6, 0, 1, 2, 3]), np.array([0, 3, 6, 10]))
>>> _ = sort_graph((y, x), links, patches)
>>> patches[0]
array([1, 5, 4, 0, 3, 6, 2, 3, 4, 7])
>>> patches[1]
array([ 0,  3,  7, 10])

Reorder links around a patch to be counterclockwise.

Examples

>>> from landlab.graph.sort.sort import sort_links_at_patch
>>> import numpy as np
>>> xy_of_node = np.array(
...     [
...         [0.0, 0.0],
...         [0.0, 1.0],
...         [1.0, 1.0],
...         [1.0, 0.0],
...     ]
... )
>>> nodes_at_link = np.array([[0, 1], [1, 2], [2, 3], [3, 0]])
>>> links_at_patch = np.array([[0, 1, 3, 2]])
>>> sort_links_at_patch(links_at_patch, nodes_at_link, xy_of_node)
>>> links_at_patch
array([[2, 1, 0, 3]])
>>> xy_of_node = np.array(
...     [
...         [0.0, 0.0],
...         [0.0, 1.0],
...         [1.0, 1.0],
...         [1.0, 0.0],
...         [2.0, 0.0],
...     ]
... )
>>> nodes_at_link = np.array([[0, 1], [1, 2], [2, 3], [3, 0], [2, 4], [3, 4]])
>>> links_at_patch = np.array([[0, 1, 3, 2], [2, 4, 5, -1]])
>>> sort_links_at_patch(links_at_patch, nodes_at_link, xy_of_node)
>>> links_at_patch
array([[ 2,  1,  0,  3],
       [ 4,  2,  5, -1]])
>>> xy_of_node = np.array(
...     [
...         [0.0, 0.0],
...         [0.0, 1.0],
...         [1.0, 1.0],
...         [1.0, 0.0],
...         [2.0, 1.0],
...     ]
... )
>>> nodes_at_link = np.array([[0, 1], [1, 2], [2, 3], [3, 0], [2, 4], [3, 4]])
>>> links_at_patch = np.array([[0, 1, 3, 2], [2, -1, 4, 5]])
>>> sort_links_at_patch(links_at_patch, nodes_at_link, xy_of_node)
>>> links_at_patch
array([[ 2,  1,  0,  3],
       [ 4,  2,  5, -1]])
sort_nodes(nodes)[source]#

Sort nodes based on their position.

Parameters:

nodes (tuple of ndarray of float) – Coordinates of nodes as (y, x).

Returns:

Array of indices that sort the nodes.

Return type:

ndarray of int

Examples

>>> from landlab.graph.sort import sort_nodes
>>> import numpy as np
>>> x = np.array([0.0, 1.0, 2.0])
>>> y = np.array([0.5, 0.0, 1.0])
>>> sort_nodes((y, x))
array([1, 0, 2])
>>> x
array([ 1.,  0.,  2.])
>>> y
array([ 0. ,  0.5,  1. ])
sort_patches(links_at_patch, offset_to_patch, xy_of_link)[source]#

Sort patches by their centroid.

Parameters:
  • links_at_patch (ndarray of int) – Links that define patches.

  • offset_to_patch (ndarray of int) – Offsets into links_at_patch for each patch.

  • xy_of_link (ndarray of float, shape (n_links, 2)) – Midpoint coordinates for each link.

Examples

>>> from landlab.graph.sort import sort_patches
>>> import numpy as np
>>> links_at_patch = np.array([0, 1, 2, 3, 2, 4])
>>> offset_to_patch = np.array([0, 3, 6])
>>> xy_of_link = np.array(
...     [[0.0, 0.5], [0.5, 1.0], [0.5, 0.5], [0.5, 0.0], [1.0, 0.5]]
... )
>>> sort_patches(links_at_patch, offset_to_patch, xy_of_link)
array([1, 0])
>>> links_at_patch
array([3, 2, 4, 0, 1, 2])
>>> offset_to_patch
array([0, 3, 6])