鐵人賽2019 Day22 資料品質：幾何檢查

GIS生產的過程中，在向量幾何資料部分有一些常用的位向關係檢查，這些檢查的品質指標在ArcGIS或是QGIS中都有各自的模組，例如商用軟體ArcGIS的topology系列模組(e.g. Creating a topology)。

大綱：

• check geometry
• 自相交
• 重複的節點
• 多邊形小於三個節點或面積長度過小

check geometry

在QGIS部分，個人使用過check geometry這個模組，他可以解決大部分的需求，它包含了許多在處理GIS資料時需檢查的幾何特性，主要操作方法是在QGIS投入shp，其提供一個互動介面，可以在UI上修正(也支援自動修正)這些問題。

以下是一些檢查項目的示意圖：


取自QGIS Geometry Checker Plugin

有些問題很明顯是數化上的失誤或是演算法處理所造成的不合理現象，應加以排除。

今天就以QGIS check geometry這個plug-in提供的一些功能，這個plug-in的指標包含了主要的位相關係檢查，這幾個指標可作為資料ETL的品質指標，今天來看看如何在python做這些檢查。

找到問題後，建議回到QGIS手動編修，某些品質問題採用全自動的編修是危險的。

自相交

polygon的自相交應是不被允許的：

from shapely.geometry import Polygon
pp = Polygon([(0, 0), (1, 1), (1, -1), (0, 1)])

檢查自相交的方法，可以用shapely的is_valid

pp.is_valid

回傳：False

shapely: A valid Polygon may not possess any overlapping exterior or interior rings.

或是另一種方法
把節點轉成lineString,記得把最後一點補上以完成封閉

from shapely.geometry import LineString
ls=LineString(data[:]+data[0:1])
ls

使用unary_union，回傳會得到MultiLineString

from shapely.ops import unary_union
mls = unary_union(ls) 
mls.geom_type # MultiLineString'

使用polygonize，將其分開

from shapely.ops import polygonize 
for polygon in polygonize(mls): 
    print(polygon)

回傳

POLYGON ((0.333 0.333, 0 0, 0 1, 0.333 0.333))
POLYGON ((0.333 0.333 0.333 0.333, 1 1, 1 -1, 0.333 0.333 0.333 0.333))

偵測到自相交的多邊形後，可能有兩種處理方式，一種是刪除其中一個（數化錯誤），一種是分離，而上面已經提供了分離的方法。

若要刪除，可進一步將分離的兩部分計算其面積，太小的則可能是錯誤：

計算面積(shapely)

polygon.area

數化錯誤產生的自相交，例如

(乍看問題不大)

(放大來看)

重複的節點

在polyline(lineString)或polygon重複的節點當然是不被允許的，以下我們故意重複一個節點

data2=[(0, 0), (1, 1),(1,1),(2,3),(0,2)]
pp2 = Polygon(data2)
pp2

要檢查這個問題，可以簡單的遍歷節點確認
或是使用以下方法直接將重複的點刪除

print(list(pp2.simplify(0).exterior.coords))

回傳：[(0.0, 0.0), (1.0, 1.0), (2.0, 3.0), (0.0, 2.0), (0.0, 0.0)]

多邊形小於三個節點或面積長度過小

多邊形一定要至少三個節點，否則是錯誤的，而面積大小也會是過濾或品質檢查的指標。

用shapely new一個polygon，兩個節點是不允損的

data=[(0, 0),(0,2)]
pp = Polygon(data3)
pp ### A LinearRing must have at least 3 coordinate tuples

polyhon面積不一定是檢查的項目，但可以避免錯誤

data=[(0, 0),(1,1),(2,3),(0,2)]
pp2= Polygon(data)
pp.area # if smaller than xxx

而多邊形中有一些特別小的點，除了可以每條edge檢查外，也可以使用shapely的simplify簡化以減少不必要的儲存量
simplify有兩種方法，preserve_topology=True僅檢查距離，如果preserve_topology=False則會使用Douglas-Peucker演算法，可以更簡化資料但是可能會有一些變形

data=[(0, 0),(0.002, 0.006), (1, 1),(1,1),(2,3),(0,2)]
pp = Polygon(data)
print(list(pp.simplify(0.1, ).exterior.coords))

回傳:[(0.0, 0.0), (1.0, 1.0), (2.0, 3.0), (0.0, 2.0), (0.0, 0.0)]

如果調成0.002，則這點被保留

print(list(pp3.simplify(0.002).exterior.coords))

回傳：[(0.0, 0.0), (0.002, 0.006), (1.0, 1.0), (2.0, 3.0), (0.0, 2.0), (0.0, 0.0)]