Proximity

Allocation

xrspatial.proximity.allocation(raster: xarray.core.dataarray.DataArray, x: str = 'x', y: str = 'y', target_values: list = [], max_distance: float = inf, distance_metric: str = 'EUCLIDEAN')[source]

Calculates, for all cells in the array, the downward slope direction Calculates, for all pixels in the input raster, the nearest source based on a set of target values and a distance metric.

This function attempts to produce the value of nearest feature of all pixels in the image to a set of pixels in the source image. The following options are used to define the behavior of the function. By default all non-zero pixels in raster.values will be considered as”target”, and all allocation will be computed in pixels.

Allocation supports NumPy backed, and Dask with NumPy backed xarray DataArray. The return values of allocation are of the same type as the input type. If input raster is a NumPy-backed DataArray, the result is NumPy-backed. If input raster is a Dask-backed DataArray, the result is Dask-backed.

allocation uses the same approach as proximity, which is ported from GDAL. A dynamic programming approach is used for identifying nearest target of a pixel from its surrounding neighborhood in a 3x3 window. The implementation for Dask-backed uses dask.map_overlap to compute allocation chunk by chunk by expanding the chunk’s borders to cover the max_distance.

Parameters
  • raster (xr.DataArray) – 2D array of target data.

  • x (str, default='x') – Name of x-coordinates.

  • y (str, default='y') – Name of y-coordinates.

  • target_values (list) – Target pixel values to measure the distance from. If this option is not provided, allocation will be computed from non-zero pixel values.

  • max_distance (float, default=np.inf) –

    The maximum distance to search. Proximity distances greater than this value will be set to NaN. Should be given in the same distance unit as input. For example, if input raster is in lat-lon and distances between points within the raster is calculated using Euclidean distance metric, max_distance should also be provided in lat-lon unit. If using Great Circle distance metric, and thus all distances is in km, max_distance should also be provided in kilometer unit.

    When scaling with Dask, whether the function scales well depends on the max_distance value. If max_distance is infinite by default, this function only works on a single machine. It should scale well, however, if max_distance is relatively small compared to the maximum possible distance in two arbitrary points in the input raster. Note that if max_distance is equal or larger than the max possible distance between 2 arbitrary points in the input raster, the input data array will be rechunked.

  • distance_metric (str, default='EUCLIDEAN') – The metric for calculating distance between 2 points. Valid distance metrics are: ‘EUCLIDEAN’, ‘GREAT_CIRCLE’, and ‘MANHATTAN’.

Returns

allocation_agg – 2D array of allocation values. All other input attributes are preserved.

Return type

xr.DataArray of same type as raster

References

Examples

>>> raster
<xarray.DataArray 'raster' (y: 5, x: 5)>
array([[0., 0., 0., 0., 0.],
       [0., 1., 0., 2., 0.],
       [0., 0., 3., 0., 0.],
       [0., 0., 0., 0., 0.],
       [1., 0., 0., 0., 0.]])
Coordinates:
  * y        (y) int64 4 3 2 1 0
  * x        (x) int64 0 1 2 3 4

>>> allocation_agg
<xarray.DataArray (y: 5, x: 5)>
array([[1., 1., 2., 2., 2.],
       [1., 1., 1., 2., 2.],
       [1., 1., 3., 2., 2.],
       [1., 3., 3., 3., 2.],
       [3., 3., 3., 3., 3.]])
Coordinates:
  * y        (y) int64 4 3 2 1 0
  * x        (x) int64 0 1 2 3 4

Direction

xrspatial.proximity.direction(raster: xarray.core.dataarray.DataArray, x: str = 'x', y: str = 'y', target_values: list = [], max_distance: float = inf, distance_metric: str = 'EUCLIDEAN')[source]

Calculates, for all cells in the array, the downward slope direction Calculates, for all pixels in the input raster, the direction to nearest source based on a set of target values and a distance metric.

This function attempts to calculate for each cell, the the direction, in degrees, to the nearest source. The output values are based on compass directions, where 90 is for the east, 180 for the south, 270 for the west, 360 for the north, and 0 for the source cell itself. The following options are used to define the behavior of the function. By default all non-zero pixels in raster.values will be considered as “target”, and all direction will be computed in pixels.

Direction support NumPy backed, and Dask with NumPy backed xarray DataArray. The return values of direction are of the same type as the input type. If input raster is a NumPy-backed DataArray, the result is NumPy-backed. If input raster is a Dask-backed DataArray, the result is Dask-backed.

Similar to proximity, the implementation for NumPy-backed is ported from GDAL, which uses a dynamic programming approach to identify nearest target of a pixel from its surrounding neighborhood in a 3x3 window The implementation for Dask-backed uses dask.map_overlap to compute proximity direction chunk by chunk by expanding the chunk’s borders to cover the max_distance.

Parameters
  • raster (xr.DataArray) – 2D array image with raster.shape = (height, width).

  • x (str, default='x') – Name of x-coordinates.

  • y (str, default='y') – Name of y-coordinates.

  • target_values (list) – Target pixel values to measure the distance from. If this option is not provided, proximity will be computed from non-zero pixel values.

  • max_distance (float, default=np.inf) –

    The maximum distance to search. Proximity distances greater than this value will be set to NaN. Should be given in the same distance unit as input. For example, if input raster is in lat-lon and distances between points within the raster is calculated using Euclidean distance metric, max_distance should also be provided in lat-lon unit. If using Great Circle distance metric, and thus all distances is in km, max_distance should also be provided in kilometer unit.

    When scaling with Dask, whether the function scales well depends on the max_distance value. If max_distance is infinite by default, this function only works on a single machine. It should scale well, however, if max_distance is relatively small compared to the maximum possible distance in two arbitrary points in the input raster. Note that if max_distance is equal or larger than the max possible distance between 2 arbitrary points in the input raster, the input data array will be rechunked.

  • distance_metric (str, default='EUCLIDEAN') – The metric for calculating distance between 2 points. Valid distance_metrics are: ‘EUCLIDEAN’, ‘GREAT_CIRCLE’, and ‘MANHATTAN’.

Returns

direction_agg – 2D array of direction values. All other input attributes are preserved.

Return type

xr.DataArray of same type as raster

References

Examples

… sourcecode:: python

>>> raster
<xarray.DataArray 'raster' (y: 5, x: 5)>
array([[0., 0., 0., 0., 0.],
       [0., 0., 0., 0., 0.],
       [0., 0., 1., 0., 0.],
       [0., 0., 0., 0., 0.],
       [1., 0., 0., 0., 0.]])
Coordinates:
  * y        (y) int64 4 3 2 1 0
  * x        (x) int64 0 1 2 3 4
>>> direction_agg
<xarray.DataArray (y: 5, x: 5)>
array([[ 45.      ,  26.56505 , 360.      , 333.43494 , 315.      ],
       [ 63.434948,  45.      , 360.      , 315.      , 296.56506 ],
       [ 90.      ,  90.      ,   0.      , 270.      , 270.      ],
       [360.      , 135.      , 180.      , 225.      , 243.43495 ],
       [  0.      , 270.      , 180.      , 206.56505 , 225.      ]],
      dtype=float32)
Coordinates:
  * y        (y) int64 4 3 2 1 0
  * x        (x) int64 0 1 2 3 4

Proximity

xrspatial.proximity.euclidean_distance(x1: float, x2: float, y1: float, y2: float) float[source]

Calculates Euclidean (straight line) distance between (x1, y1) and (x2, y2).

Parameters
  • x1 (float) – x-coordinate of the first point.

  • x2 (float) – x-coordinate of the second point.

  • y1 (float) – y-coordinate of the first point.

  • y2 (float) – y-coordinate of the second point.

Returns

distance – Euclidean distance between two points.

Return type

float

References

Examples

>>> # Imports
>>> from xrspatial import euclidean_distance
>>> point_a = (142.32, 23.23)
>>> point_b = (312.54, 432.01)

>>> # Calculate Euclidean Distance
>>> dist = euclidean_distance(
        point_a[0],
        point_b[0],
        point_a[1],
        point_b[1],
    )
>>> print(dist)
442.80462599209596
xrspatial.proximity.great_circle_distance(x1: float, x2: float, y1: float, y2: float, radius: float = 6378137) float[source]

Calculates great-circle (orthodromic/spherical) distance between (x1, y1) and (x2, y2), assuming each point is a longitude, latitude pair.

Parameters
  • x1 (float) – x-coordinate (latitude) between -180 and 180 of the first point.

  • x2 (float) – x-coordinate (latitude) between -180 and 180 of the second point.

  • y1 (float) – y-coordinate (longitude) between -90 and 90 of the first point.

  • y2 (float) – y-coordinate (longitude) between -90 and 90 of the second point.

  • radius (float, default=6378137) – Radius of sphere (earth).

Returns

distance – Great-Circle distance between two points.

Return type

float

References

Examples

>>> # Imports
>>> from xrspatial import great_circle_distance

>>> point_a = (123.2, 82.32)
>>> point_b = (178.0, 65.09)

>>> # Calculate Euclidean Distance
>>> dist = great_circle_distance(
        point_a[0],
        point_b[0],
        point_a[1],
        point_b[1],
    )
>>> print(dist)
2378290.489801402
xrspatial.proximity.manhattan_distance(x1: float, x2: float, y1: float, y2: float) float[source]

Calculates Manhattan distance (sum of distance in x and y directions) between (x1, y1) and (x2, y2).

Parameters
  • x1 (float) – x-coordinate of the first point.

  • x2 (float) – x-coordinate of the second point.

  • y1 (float) – y-coordinate of the first point.

  • y2 (float) – y-coordinate of the second point.

Returns

distance – Manhattan distance between two points.

Return type

float

References

Examples

>>> # Imports
>>> from xrspatial import manhattan_distance

>>> point_a = (142.32, 23.23)
>>> point_b = (312.54, 432.01)

>>> # Calculate Euclidean Distance
>>> dist = manhattan_distance(
        point_a[0],
        point_b[0],
        point_a[1],
        point_b[1],
    )
>>> print(dist)
196075.9368
xrspatial.proximity.proximity(raster: xarray.core.dataarray.DataArray, x: str = 'x', y: str = 'y', target_values: list = [], max_distance: float = inf, distance_metric: str = 'EUCLIDEAN') xarray.core.dataarray.DataArray[source]

Computes the proximity of all pixels in the image to a set of pixels in the source image based on a distance metric.

This function attempts to compute the proximity of all pixels in the image to a set of pixels in the source image. The following options are used to define the behavior of the function. By default all non-zero pixels in raster.values will be considered the “target”, and all proximities will be computed in pixels. Note that target pixels are set to the value corresponding to a distance of zero.

Proximity support NumPy backed, and Dask with NumPy backed xarray DataArray. The return values of proximity are of the same type as the input type. If input raster is a NumPy-backed DataArray, the result is NumPy-backed. If input raster is a Dask-backed DataArray, the result is Dask-backed.

The implementation for NumPy-backed is ported from GDAL, which uses a dynamic programming approach to identify nearest target of a pixel from its surrounding neighborhood in a 3x3 window. The implementation for Dask-backed uses dask.map_overlap to compute proximity chunk by chunk by expanding the chunk’s borders to cover the max_distance.

Parameters
  • raster (xr.DataArray) – 2D array image with raster.shape = (height, width).

  • x (str, default='x') – Name of x-coordinates.

  • y (str, default='y') – Name of y-coordinates.

  • target_values (list) – Target pixel values to measure the distance from. If this option is not provided, proximity will be computed from non-zero pixel values.

  • max_distance (float, default=np.inf) –

    The maximum distance to search. Proximity distances greater than this value will be set to NaN. Should be given in the same distance unit as input. For example, if input raster is in lat-lon and distances between points within the raster is calculated using Euclidean distance metric, max_distance should also be provided in lat-lon unit. If using Great Circle distance metric, and thus all distances is in km, max_distance should also be provided in kilometer unit.

    When scaling with Dask, whether the function scales well depends on the max_distance value. If max_distance is infinite by default, this function only works on a single machine. It should scale well, however, if max_distance is relatively small compared to the maximum possible distance in two arbitrary points in the input raster. Note that if max_distance is equal or larger than the max possible distance between 2 arbitrary points in the input raster, the input data array will be rechunked.

  • distance_metric (str, default='EUCLIDEAN') – The metric for calculating distance between 2 points. Valid distance metrics are: ‘EUCLIDEAN’, ‘GREAT_CIRCLE’, and ‘MANHATTAN’.

Returns

proximity_agg – 2D array of proximity values. All other input attributes are preserved.

Return type

xr.DataArray of same type as raster

References

Examples

>>> raster
<xarray.DataArray 'raster' (y: 5, x: 5)>
array([[0., 0., 0., 0., 0.],
       [0., 0., 0., 1., 0.],
       [0., 0., 0., 0., 0.],
       [0., 0., 0., 0., 0.],
       [1., 0., 0., 0., 0.]])
Coordinates:
  * y        (y) int64 4 3 2 1 0
  * x        (x) int64 0 1 2 3 4

>>> proximity_agg
<xarray.DataArray (y: 5, x: 5)>
array([[3.1622777, 2.236068 , 1.4142135, 1.       , 1.4142135],
       [3.       , 2.       , 1.       , 0.       , 1.       ],
       [3.1622777, 2.236068 , 1.4142135, 1.       , 1.4142135],
       [3.6055512, 2.828427 , 2.236068 , 2.       , 2.236068 ],
       [4.2426405, 3.6055512, 3.1622777, 3.       , 3.1622777]],
      dtype=float32)
Coordinates:
  * y        (y) int64 4 3 2 1 0
  * x        (x) int64 0 1 2 3 4