Skip to content
Open
Show file tree
Hide file tree
Changes from 10 commits
Commits
File filter

Filter by extension

Filter by extension

Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
90 changes: 75 additions & 15 deletions libpysal/graph/_triangulation.py
Original file line number Diff line number Diff line change
@@ -1,8 +1,10 @@
import warnings
from functools import wraps

import geopandas
import numpy
import pandas
import shapely
from scipy import sparse, spatial

from libpysal.cg import voronoi_frames
Expand Down Expand Up @@ -342,20 +344,27 @@ def _relative_neighborhood(coordinates, coplanar):
return heads_ix, tails_ix, coplanar


@_validate_coplanar
def _voronoi(coordinates, coplanar, clip="bounding_box", rook=True):
def _voronoi(
geoms,
ids=None,
coplanar="raise",
clip="bounding_box",
rook=True,
seed=None,
**kwargs,
):
"""
Compute contiguity weights according to a clipped
Voronoi diagram.

Parameters
---------
coordinates : numpy.ndarray, geopandas.GeoSeries, geopandas.GeoDataFrame
geometries containing locations to compute the delaunay triangulation. If
a geopandas object with Point geoemtry is provided, the .geometry attribute
is used. If a numpy.ndarray with shapely geoemtry is used, then the
coordinates are extracted and used. If a numpy.ndarray of a shape (2,n) is
used, it is assumed to contain x, y coordinates.
geoms : numpy.ndarray, geopandas.GeoSeries, geopandas.GeoDataFrame
Geometries to compute the Voronoi diagram. Accepts Point, Polygon,
and MultiPolygon geometries. For non-point geometries, the boundaries
are discretised into points and dissolved back using
``voronoi_frames()``. If a numpy.ndarray of shape (2, n) is provided,
it is assumed to contain x, y coordinates.
ids : numpy.narray (default: None)
ids to use for each sample in coordinates. Generally, construction functions
that are accessed via Graph.build_kernel() will set this automatically from
Expand Down Expand Up @@ -398,6 +407,10 @@ def _voronoi(coordinates, coplanar, clip="bounding_box", rook=True):
An integer value used to ensure that the pseudorandom number generator provides
the same value over replications. By default, no seed is used, so results
will be random every time. This is only used if coplanar='jitter'.
**kwargs: Additional keyword arguments passed to ``voronoi_frames()``.
Supports ``segment`` (float) to control boundary discretisation
and ``shrink`` (float) to shrink geometries before tessellation.
Only used when ``method="voronoi"``.

Notes
-----
Expand All @@ -411,21 +424,68 @@ def _voronoi(coordinates, coplanar, clip="bounding_box", rook=True):
delaunay triangulations in many applied contexts and
generally will remove "long" links in the delaunay graph.
"""
if coplanar == "raise":
unique = numpy.unique(coordinates, axis=0)
if unique.shape != coordinates.shape:
if isinstance(geoms, (geopandas.GeoSeries, geopandas.GeoDataFrame)):
geoms = geoms.geometry
if ids is None:
ids = numpy.asarray(geoms.index)

if ids is None:
ids = numpy.arange(len(geoms))

if coplanar == "jitter":
coords = shapely.get_coordinates(geoms)
coords = _jitter_geoms(coords, seed=seed)
geoms = geopandas.GeoSeries(
shapely.points(coords),
index=geoms.index if hasattr(geoms, "index") else None,

Copy link
Copy Markdown
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

This is going to get very wrong when non-point geometry is given.

Please ensure that you write test suite alongside your contributions. Many of the issues you'd be able to catch yourself.

)

voronoi_kwargs = {k: v for k, v in kwargs.items() if k in ("segment", "shrink")}

if coplanar == "raise" and not isinstance(
geoms, (geopandas.GeoSeries, geopandas.GeoDataFrame)
):
unique = numpy.unique(geoms, axis=0)
if unique.shape != geoms.shape:
raise CoplanarError(
f"There are {len(unique)} unique locations in "
f"the dataset, but {len(coordinates)} observations. This means there "
f"the dataset, but {len(geoms)} observations. This means there "
"are multiple points in the same location, which is undefined "
"for this graph type. To address this issue, consider setting "
"`coplanar='clique'` or consult the documentation about "
"coplanar points."
)
cells = voronoi_frames(coordinates, clip=clip, return_input=False, as_gdf=False)
heads_ix, tails_ix, weights = _vertex_set_intersection(cells, rook=rook)

return heads_ix, tails_ix, numpy.array([])
elif coplanar == "raise" and isinstance(
geoms, (geopandas.GeoSeries, geopandas.GeoDataFrame)
):
if geoms.geom_type.isin(["Point"]).all():
coords = shapely.get_coordinates(geoms)
unique = numpy.unique(coords, axis=0)
if len(unique) != len(geoms):
raise CoplanarError(
f"There are {len(unique)} unique locations in "
f"the dataset, but {len(geoms)} observations. This means there "
"are multiple points in the same location, which is undefined "
"for this graph type. To address this issue, consider setting "
"`coplanar='clique'` or consult the documentation about "
"coplanar points."
)

cells = voronoi_frames(
geoms, clip=clip, return_input=False, as_gdf=False, **voronoi_kwargs
)
heads_ix, tails_ix, _ = _vertex_set_intersection(cells, rook=rook)
n = len(ids)
mask = (heads_ix < n) & (tails_ix < n)
heads_ix = heads_ix[mask]
tails_ix = tails_ix[mask]

heads = ids[heads_ix]
tails = ids[tails_ix]
weights = numpy.ones(len(heads_ix), dtype=numpy.int8)

return heads, tails, weights
Comment on lines +486 to +489

Copy link
Copy Markdown
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Suggested change
heads = heads_ix
tails = tails_ix
weights = numpy.ones(len(heads_ix), dtype=numpy.int8)
return heads, tails, weights
weights = numpy.ones(len(heads_ix), dtype=numpy.int8)
return heads_ix, tails_ix, weights



#### utilities
Expand Down
5 changes: 5 additions & 0 deletions libpysal/graph/base.py
Original file line number Diff line number Diff line change
Expand Up @@ -1377,6 +1377,7 @@ def build_triangulation(
coplanar="raise",
taper=True,
decay=False,
**kwargs,
):
"""Generate Graph from geometry based on triangulation

Expand Down Expand Up @@ -1442,6 +1443,9 @@ def build_triangulation(
or negative) at some very large (possibly infinite) distance.
Otherwise, kernel functions are treated as proper
volume-preserving probability distributions.
**kwargs: Additional keyword arguments passed to ``voronoi_frames()`` when
``method="voronoi"``. Supports ``segment`` (float) and ``shrink``
(float). Ignored for other methods.
Comment on lines +1453 to +1454

Copy link
Copy Markdown
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Only for non-point inputs.


Returns
-------
Expand Down Expand Up @@ -1527,6 +1531,7 @@ def build_triangulation(
coplanar=coplanar,
decay=decay,
taper=taper,
**kwargs,
)
else:
raise ValueError(
Expand Down