[1]:
import geopandas
import numpy
import matplotlib.pyplot as plt
import geoplanar
import libpysal
from shapely.geometry import Polygon
Planar Enforcement Violation: non-planar enforced edges#
[2]:
c1 = [[0,0], [0, 10], [10, 10], [10, 0], [0, 0]]
p1 = Polygon(c1)
c2 = [[10, 2], [10, 8], [20, 8], [20, 2], [10, 2]]
p2 = Polygon(c2)
gdf = geopandas.GeoDataFrame(geometry=[p1, p2])
base = gdf.plot(edgecolor='k', facecolor="none",alpha=0.5)
c1 = numpy.array(c1)
c2 = numpy.array(c2)
_ = base.scatter(c1[:,0], c1[:,1])
_ =base.scatter(c2[:,0], c2[:,1])
The two polygons are visually contiguous, but are not planar enforced as the right edge of the left polygon and the left edge of right polygon share no vertices. This will result in the two polygons not being Queen neighbors, since a necessary (and sufficient) condition for the latter is at least one shared vertex.
[3]:
w = libpysal.weights.Queen.from_dataframe(gdf)
/Users/serge/miniconda3/envs/edu_concordance/lib/python3.9/site-packages/libpysal/weights/weights.py:172: UserWarning: The weights matrix is not fully connected:
There are 2 disconnected components.
There are 2 islands with ids: 0, 1.
warnings.warn(message)
Detecting nonplanar edges#
geoplanar can detect and report nonplanar edges:
[4]:
geoplanar.non_planar_edges(gdf)
[4]:
defaultdict(set, {0: {1}})
Fixing nonplanar edges#
[5]:
geoplanar.is_planar_enforced(gdf)
[5]:
False
[6]:
gdf1 = geoplanar.fix_npe_edges(gdf)
/Users/serge/miniconda3/envs/edu_concordance/lib/python3.9/site-packages/libpysal/weights/weights.py:172: UserWarning: The weights matrix is not fully connected:
There are 2 disconnected components.
There are 2 islands with ids: 0, 1.
warnings.warn(message)
[7]:
geoplanar.non_planar_edges(gdf1)
[7]:
defaultdict(set, {})
[8]:
w1 = libpysal.weights.Queen.from_dataframe(gdf1)
w1.neighbors
[8]:
{0: [1], 1: [0]}
[9]:
geoplanar.is_planar_enforced(gdf1)
[9]:
True
Planar Enforcement Violation: Overlapping and non-planar enforced edges#
[10]:
from shapely.geometry import Polygon
[11]:
t1 = Polygon([[0,0],[10,10], [20,0]])
[12]:
t1
[12]:
[13]:
b1 = Polygon([[5,5], [20,5], [20,-10], [0,-10]])
b1
[13]:
[14]:
gdf = geopandas.GeoDataFrame(geometry=[t1,b1])
[15]:
gdf.plot(edgecolor='k',facecolor="none",alpha=0.5) # non planar enforcement
[15]:
<AxesSubplot:>
The two features overlap and would appear to share vertices, but they in fact do not share vertices. Again, because this violates planar enforcement, this results in two polygons not being Queen neighbors:
[16]:
import libpysal
w = libpysal.weights.Queen.from_dataframe(gdf)
/Users/serge/miniconda3/envs/edu_concordance/lib/python3.9/site-packages/libpysal/weights/weights.py:172: UserWarning: The weights matrix is not fully connected:
There are 2 disconnected components.
There are 2 islands with ids: 0, 1.
warnings.warn(message)
Detecting nonplanar edges#
geoplanar will use a failed contiguity check as part of a check for nonplanar enforced edges in the polygons of a geoseries:
[17]:
geoplanar.non_planar_edges(gdf)
[17]:
defaultdict(set, {0: {1}})
Correcting nonplanar edges#
[18]:
gdf_fixed = geoplanar.fix_npe_edges(gdf)
geoplanar.non_planar_edges(gdf_fixed)
/Users/serge/miniconda3/envs/edu_concordance/lib/python3.9/site-packages/libpysal/weights/weights.py:172: UserWarning: The weights matrix is not fully connected:
There are 2 disconnected components.
There are 2 islands with ids: 0, 1.
warnings.warn(message)
[18]:
defaultdict(set, {})
Default is to work on a copy#
[19]:
geoplanar.non_planar_edges(gdf)
/Users/serge/miniconda3/envs/edu_concordance/lib/python3.9/site-packages/libpysal/weights/weights.py:172: UserWarning: The weights matrix is not fully connected:
There are 2 disconnected components.
There are 2 islands with ids: 0, 1.
warnings.warn(message)
[19]:
defaultdict(set, {0: {1}})
[20]:
geoplanar.fix_npe_edges(gdf, inplace=True)
[20]:
| geometry | |
|---|---|
| 0 | POLYGON ((0.00000 0.00000, 5.00000 5.00000, 10... |
| 1 | POLYGON ((5.00000 5.00000, 15.00000 5.00000, 2... |
[21]:
geoplanar.non_planar_edges(gdf)
[21]:
defaultdict(set, {})
[22]:
w = libpysal.weights.Queen.from_dataframe(gdf)
w.neighbors
[22]:
{0: [1], 1: [0]}
Handle nonplanar edges in multi polygon case#
[23]:
from shapely.geometry import MultiPolygon
[24]:
c1 = [[0,0], [0, 10], [10, 10], [10, 0], [0, 0]]
p1 = Polygon(c1)
c2 = [[10, 2], [10, 8], [20, 8], [20, 2], [10, 2]]
p2 = Polygon(c2)
p3 = Polygon([ [21, 2], [21, 4], [23,3] ])
p2 = MultiPolygon([Polygon(c2), p3])
gdf = geopandas.GeoDataFrame(geometry=[p1, p2])
base = gdf.plot(edgecolor='k', facecolor="none",alpha=0.5)
c1 = numpy.array(c1)
c2 = numpy.array(c2)
_ = base.scatter(c1[:,0], c1[:,1])
_ =base.scatter(c2[:,0], c2[:,1])
[25]:
w = libpysal.weights.Queen.from_dataframe(gdf)
/Users/serge/miniconda3/envs/edu_concordance/lib/python3.9/site-packages/libpysal/weights/_contW_lists.py:31: ShapelyDeprecationWarning: Iteration over multi-part geometries is deprecated and will be removed in Shapely 2.0. Use the `geoms` property to access the constituent parts of a multi-part geometry.
return list(it.chain(*(_get_boundary_points(part.boundary) for part in shape)))
/Users/serge/miniconda3/envs/edu_concordance/lib/python3.9/site-packages/libpysal/weights/weights.py:172: UserWarning: The weights matrix is not fully connected:
There are 2 disconnected components.
There are 2 islands with ids: 0, 1.
warnings.warn(message)
[26]:
geoplanar.non_planar_edges(gdf)
[26]:
defaultdict(set, {0: {1}})
[27]:
geoplanar.is_planar_enforced(gdf)
[27]:
False
[28]:
gdf1 = geoplanar.fix_npe_edges(gdf)
/Users/serge/miniconda3/envs/edu_concordance/lib/python3.9/site-packages/libpysal/weights/_contW_lists.py:31: ShapelyDeprecationWarning: Iteration over multi-part geometries is deprecated and will be removed in Shapely 2.0. Use the `geoms` property to access the constituent parts of a multi-part geometry.
return list(it.chain(*(_get_boundary_points(part.boundary) for part in shape)))
/Users/serge/miniconda3/envs/edu_concordance/lib/python3.9/site-packages/libpysal/weights/weights.py:172: UserWarning: The weights matrix is not fully connected:
There are 2 disconnected components.
There are 2 islands with ids: 0, 1.
warnings.warn(message)
[29]:
geoplanar.is_planar_enforced(gdf1)
/Users/serge/miniconda3/envs/edu_concordance/lib/python3.9/site-packages/libpysal/weights/_contW_lists.py:31: ShapelyDeprecationWarning: Iteration over multi-part geometries is deprecated and will be removed in Shapely 2.0. Use the `geoms` property to access the constituent parts of a multi-part geometry.
return list(it.chain(*(_get_boundary_points(part.boundary) for part in shape)))
[29]:
True