STRTree#

class STRtree(geoms: Iterable[shapely.geometry.base.BaseGeometry], node_capacity: int = 10)#

A query-only R-tree spatial index created using the Sort-Tile-Recursive (STR) [1] algorithm.

For two-dimensional spatial data. The tree is constructed directly at initialization. The tree is immutable and query-only, meaning that once created nodes cannot be added or removed.

All operations return indices of the input geometries. These indices can be used to index into anything associated with the input geometries, including the input geometries themselves, or custom items stored in another object of the same length as the geometries.

Any mixture of geometry types may be stored in the tree.

Parameters
geomssequence

A sequence of geometry objects.

node_capacityint, default 10

The maximum number of child nodes per parent node in the tree.

References

1

Leutenegger, Scott T.; Edgington, Jeffrey M.; Lopez, Mario A. (February 1997). “STR: A Simple and Efficient Algorithm for R-Tree Packing”. https://ia600900.us.archive.org/27/items/nasa_techdoc_19970016975/19970016975.pdf

property geometries#

Geometries stored in the tree in the order used to construct the tree.

The order of this array corresponds to the tree indices returned by other STRtree methods.

Do not attempt to modify items in the returned array.

Returns
ndarray of Geometry objects
nearest(geometry) Optional[Any]#

Return the index of the nearest geometry in the tree for each input geometry.

Note

‘nearest’ requires at least GEOS 3.6.0.

If there are multiple equidistant or intersected geometries in the tree, only a single result is returned for each input geometry, based on the order that tree geometries are visited; this order may be nondeterministic.

If any input geometry is None or empty, an error is raised.

Parameters
geometryGeometry or array_like

Input geometries to query the tree.

Returns
scalar or ndarray

Indices of geometries in tree. Return value will have the same shape as the input.

None is returned if this index is empty. This may change in version 2.0.

See also

query_nearest

returns all equidistant geometries, exclusive geometries, and optional distances

Examples

>>> from shapely.geometry import Point
>>> tree = STRtree([Point(i, i) for i in range(10)])

Query the tree for nearest using a scalar geometry:

>>> index = tree.nearest(Point(2.2, 2.2))
>>> index
2
>>> tree.geometries.take(index)
<POINT (2 2)>

Query the tree for nearest using an array of geometries:

>>> indices = tree.nearest([Point(2.2, 2.2), Point(4.4, 4.4)])
>>> indices.tolist()
[2, 4]
>>> tree.geometries.take(indices).tolist()
[<POINT (2 2)>, <POINT (4 4)>]

Nearest only return one object if there are multiple equidistant results:

>>> tree = STRtree ([Point(0, 0), Point(0, 0)])
>>> tree.nearest(Point(0, 0))
0
query(geometry, predicate=None, distance=None)#

Return the integer indices of all combinations of each input geometry and tree geometries where the extent of each input geometry intersects the extent of a tree geometry.

If the input geometry is a scalar, this returns an array of shape (n, ) with the indices of the matching tree geometries. If the input geometry is an array_like, this returns an array with shape (2,n) where the subarrays correspond to the indices of the input geometries and indices of the tree geometries associated with each. To generate an array of pairs of input geometry index and tree geometry index, simply transpose the result.

If a predicate is provided, the tree geometries are further filtered to those that meet the predicate when comparing the input geometry to the tree geometry: predicate(geom, tree_geometry)

The ‘dwithin’ predicate requires GEOS >= 3.10.

Any input geometry that is None or empty will never match geometries in the tree.

Parameters
geometryGeometry or array_like

Input geometries to query the tree and filter results using the optional predicate.

predicate{None, ‘intersects’, ‘within’, ‘contains’, ‘overlaps’, ‘crosses’,’touches’, ‘covers’, ‘covered_by’, ‘contains_properly’, ‘dwithin’}, optional

The predicate to use for testing geometries from the tree that are within the input geometry’s extent.

distancenumber or array_like, optional

Distances around each input geometry within which to query the tree for the ‘dwithin’ predicate. If array_like, shape must be broadcastable to shape of geometry. Required if predicate=’dwithin’.

Returns
ndarray with shape (n,) if geometry is a scalar

Contains tree geometry indices.

OR
ndarray with shape (2, n) if geometry is an array_like

The first subarray contains input geometry indices. The second subarray contains tree geometry indices.

Notes

In the context of a spatial join, input geometries are the “left” geometries that determine the order of the results, and tree geometries are “right” geometries that are joined against the left geometries. This effectively performs an inner join, where only those combinations of geometries that can be joined based on overlapping extents or optional predicate are returned.

Examples

>>> from shapely import box, Point
>>> import numpy as np
>>> points = [Point(0, 0), Point(1, 1), Point(2,2), Point(3, 3)]
>>> tree = STRtree(points)

Query the tree using a scalar geometry:

>>> indices = tree.query(box(0, 0, 1, 1))
>>> indices.tolist()
[0, 1]

Query using an array of geometries:

>>> boxes = np.array([box(0, 0, 1, 1), box(2, 2, 3, 3)])
>>> arr_indices = tree.query(boxes)
>>> arr_indices.tolist()
[[0, 0, 1, 1], [0, 1, 2, 3]]

Or transpose to get all pairs of input and tree indices:

>>> arr_indices.T.tolist()
[[0, 0], [0, 1], [1, 2], [1, 3]]

Retrieve the tree geometries by results of query:

>>> tree.geometries.take(indices).tolist()
[<POINT (0 0)>, <POINT (1 1)>]

Retrieve all pairs of input and tree geometries:

>>> np.array([boxes.take(arr_indices[0]),tree.geometries.take(arr_indices[1])]).T.tolist()
[[<POLYGON ((1 0, 1 1, 0 1, 0 0, 1 0))>, <POINT (0 0)>],
 [<POLYGON ((1 0, 1 1, 0 1, 0 0, 1 0))>, <POINT (1 1)>],
 [<POLYGON ((3 2, 3 3, 2 3, 2 2, 3 2))>, <POINT (2 2)>],
 [<POLYGON ((3 2, 3 3, 2 3, 2 2, 3 2))>, <POINT (3 3)>]]

Query using a predicate:

>>> tree = STRtree([box(0, 0, 0.5, 0.5), box(0.5, 0.5, 1, 1), box(1, 1, 2, 2)])
>>> tree.query(box(0, 0, 1, 1), predicate="contains").tolist()
[0, 1]
>>> tree.query(Point(0.75, 0.75), predicate="dwithin", distance=0.5).tolist()
[0, 1, 2]
>>> tree.query(boxes, predicate="contains").tolist()
[[0, 0], [0, 1]]
>>> tree.query(boxes, predicate="dwithin", distance=0.5).tolist()
[[0, 0, 0, 1], [0, 1, 2, 2]]

Retrieve custom items associated with tree geometries (records can be in whatever data structure so long as geometries and custom data can be extracted into arrays of the same length and order):

>>> records = [
...     {"geometry": Point(0, 0), "value": "A"},
...     {"geometry": Point(2, 2), "value": "B"}
... ]
>>> tree = STRtree([record["geometry"] for record in records])
>>> items = np.array([record["value"] for record in records])
>>> items.take(tree.query(box(0, 0, 1, 1))).tolist()
['A']
query_nearest(geometry, max_distance=None, return_distance=False, exclusive=False, all_matches=True)#

Returns the index of the nearest geometries in the tree for each input geometry.

Note

‘query_nearest’ requires at least GEOS 3.6.0.

If there are multiple equidistant or intersected geometries in tree and all_matches is True (the default), all matching tree geometries are returned; otherwise only the first matching tree geometry is returned. Tree indices are returned in the order they are visited for each input geometry and may not be in ascending index order; no meaningful order is implied.

The max_distance used to search for nearest items in the tree may have a significant impact on performance by reducing the number of input geometries that are evaluated for nearest items in the tree. Only those input geometries with at least one tree geometry within +/- max_distance beyond their envelope will be evaluated. However, using a large max_distance may have a negative performance impact because many tree geometries will be queried for each input geometry.

The distance, if returned, will be 0 for any intersected geometries in the tree.

Any geometry that is None or empty in the input geometries is omitted from the output.

Parameters
geometryGeometry or array_like

Input geometries to query the tree.

max_distancefloat, optional

Maximum distance within which to query for nearest items in tree. Must be greater than 0.

return_distancebool, default False

If True, will return distances in addition to indices.

exclusivebool, default False

If True, the nearest tree geometries that are equal to the input geometry will not be returned.

all_matchesbool, default True

If True, all equidistant and intersected geometries will be returned for each input geometry. If False, only the first nearest geometry will be returned.

Returns
tree indices or tuple of (tree indices, distances) if geometry is a scalar

indices is an ndarray of shape (n, ) and distances (if present) an ndarray of shape (n, )

OR
indices or tuple of (indices, distances)

indices is an ndarray of shape (2,n) and distances (if present) an ndarray of shape (n). The first subarray of indices contains input geometry indices. The second subarray of indices contains tree geometry indices.

See also

nearest

returns singular nearest geometry for each input

Examples

>>> import numpy as np
>>> from shapely import box, Point
>>> points = [Point(0, 0), Point(1, 1), Point(2,2), Point(3, 3)]
>>> tree = STRtree(points)

Find the nearest tree geometries to a scalar geometry:

>>> indices = tree.query_nearest(Point(0.25, 0.25))
>>> indices.tolist()
[0]

Retrieve the tree geometries by results of query:

>>> tree.geometries.take(indices).tolist()
[<POINT (0 0)>]

Find the nearest tree geometries to an array of geometries:

>>> query_points = np.array([Point(2.25, 2.25), Point(1, 1)])
>>> arr_indices = tree.query_nearest(query_points)
>>> arr_indices.tolist()
[[0, 1], [2, 1]]

Or transpose to get all pairs of input and tree indices:

>>> arr_indices.T.tolist()
[[0, 2], [1, 1]]

Retrieve all pairs of input and tree geometries:

>>> list(zip(query_points.take(arr_indices[0]), tree.geometries.take(arr_indices[1])))
[(<POINT (2.25 2.25)>, <POINT (2 2)>), (<POINT (1 1)>, <POINT (1 1)>)]

All intersecting geometries in the tree are returned by default:

>>> tree.query_nearest(box(1,1,3,3)).tolist()
[1, 2, 3]

Set all_matches to False to to return a single match per input geometry:

>>> tree.query_nearest(box(1,1,3,3), all_matches=False).tolist()
[1]

Return the distance to each nearest tree geometry:

>>> index, distance = tree.query_nearest(Point(0.5, 0.5), return_distance=True)
>>> index.tolist()
[0, 1]
>>> distance.round(4).tolist()
[0.7071, 0.7071]

Return the distance for each input and nearest tree geometry for an array of geometries:

>>> indices, distance = tree.query_nearest([Point(0.5, 0.5), Point(1, 1)], return_distance=True)
>>> indices.tolist()
[[0, 0, 1], [0, 1, 1]]
>>> distance.round(4).tolist()
[0.7071, 0.7071, 0.0]

Retrieve custom items associated with tree geometries (records can be in whatever data structure so long as geometries and custom data can be extracted into arrays of the same length and order):

>>> records = [
...     {"geometry": Point(0, 0), "value": "A"},
...     {"geometry": Point(2, 2), "value": "B"}
... ]
>>> tree = STRtree([record["geometry"] for record in records])
>>> items = np.array([record["value"] for record in records])
>>> items.take(tree.query_nearest(Point(0.5, 0.5))).tolist()
['A']