鐵人賽2019 Day23 資料品質：幾何檢查(續篇)

Nov 7, 2018

data analysis

1. gap and sliver polygon
2. 檢查方法
3. 修正
4. 後記

今天繼續GIS的幾何檢查的部分，主要處理的是polygon的接邊問題，這是一種常見的topology的問題。

大綱:

gap and sliver polygon
檢查方法
修正

gap and sliver polygon

兩個相接的polygon可能因為繪製的錯誤，或是浮點數的位數問題，造成邊界些微的不吻合(shp檔沒有儲存topology)，可能有兩種情況，第一種是中間有縫隙(gap)，第二種則是中間有一些重疊(sliver)。

檢查方法

為了測試，我自行畫了一個範例：

總共有三個Polygon，紅色部分highlight起來凸顯邊界問題，另外一個Polygon則是分開的

要檢查是否有因邊界不吻合所造成gap或sliver，
最簡單的方法，可以將polygon組為multyPolygon，利用shapely的is_valid檢查：

A valid MultiPolygon may not collect any overlapping polygons. — Shapely

用Geopandas讀自己畫的Polygon

1
2
3

import geopandas as gpd
from shapely.geometry import MultiPolygon, Polygon
test=gpd.read_file('data/test.shp',encoding='utf-8')

將polygon丟進去multypolygon

from shapely.geometry import MultiPolygon, Polygon
data=[]
for i in range(len(test.geometry)):
    data.append(test['geometry'][i])
mp = MultiPolygon(data)
mp_obj = gpd.GeoSeries(mp)
mp_obj.is_valid

回傳:False
由於很明顯有不一致的現象，所以回傳是False

修正

確定這分資料有邊界的問題後，接下來是修正
昨天提到，修正最好還是在GIS軟體中編修介面中操作，
但我們可以試著將要修正的區塊找出來

首先可以使用unary_union將聯集的部分合併

1 2	union = test.geometry.unary_union union

結果由三個polygon被合併為兩個
1號polygon是被合併的部分，我們試著找出邊界不一致的區塊
可以例用difference跟intersection的運算找出gap或sliver

找gap:

1
2
3

# gap
gap = list(union)[1]  # 
gap=gap.difference(test['geometry'][2]).difference(test['geometry'][1])

找sliver

1 2	# sliver sliver=test['geometry'][2].intersection(test['geometry'][1])

把成果疊一起示意：

ax=mp_obj.plot()
gap_obj = gpd.GeoSeries(gap)
sliver_obj = gpd.GeoSeries(sliver)
ax=gap_obj.plot(ax=ax,color='yellow')
sliver_obj.plot(ax=ax,color='red')