基於Python Shapely的幾何集合操作

  • 2019 年 10 月 4 日
  • 筆記

shapely是基於笛卡爾坐標的幾何對象操作和分析Python庫，底層基於GEOS和JTS庫。

shapely無法讀取和寫數據文件，但可以基於應用廣泛的一些格式和協議進行序列化(serialize)和去序列化(deserialize)操作。

shapely不關注數據格式和坐標系統，但shapely的整合性很強，可以和GIS之類的工具協同工作。這種黏性類似python。

安裝

基於構建的發行版

windows

  • conda install shapely
  • 基於 wheels 安裝 (http://www.lfd.uci.edu/~gohlke/pythonlibs/#shapely)

Mac OS和Linux

  • pip install shapely，如果需要針對向量化加速版本可通過pip install shapely[vectorized]安裝
  • 通過系統包管理器，比如 apt，yum 和 Homebrew等。
  • 也可以通過Canopy和Anaconda等Python發行版工具安裝，比如Anaconda，conda install shapely

基於源碼

當需要兼容基於GEOS的更多模組，或者想要使用不同的GEOS版本，可以基於源碼進行安裝：

pip install shapely --no-binaryshapely

如果使用自定義GEOS版本進行安裝時，可能需要指定geos-config程式的路徑，

GEOS_CONFIG = /path/to/geos-config pip install shapely

基本操作

  • 創建點 from shapely.geometry import Point point = Point(0, 0) # Point((0, 0)) point.area # 獲取點的面積 point.length # 獲取點的長度 point.bounds # 獲取點的邊界
  • 創建圓 In[20]: circle = Point(0, 0).buffer(10) # 創建以(0, 0)為圓心，10為半徑的圓 In[25]: circle.area # 獲取創建的圓的面積 Out[25]: 313.6548490545939 從上述結果可以看出，所創建的圓的面積小於pi r^2，這是因為buffer方法默認參數resolution為16，resolution 的值越大圓越完整。 In[30]: circle= Point(0, 0).buffer(10, resolution=1000) In[31]: circle.area Out[31]: 314.15913616617644 In[32]: circle= Point(0, 0).buffer(10, resolution=1000000) In[33]: circle.area Out[33]: 314.1592653588436
  • 創建多邊形 from shapely.geometry import Polygon polygon = Polygon([(0, 1), (0, 2), (0, 3), (1, 1), (1, 2), (1, 3), (0, 3)]) 注意： Polygon 函數僅能基於有序的點創建多邊形，且點的集合必須要是閉合的。使用MultiPoint 函數創建，並使用 convex_hull 方法創建多邊形。 from shapely.geometry import MultiPoint coords = [(0, 1), (1, 2), (1, 4), (2, 0), (3, 2)] # coords不一定要是閉合點集合 poly = MultiPoint(coords).convex_hull

集合操作

  • 判斷點是否在多邊形
In[50]: p1 = Point(24.952242, 60.1696017)  In[51]: p2 = Point(24.976567, 60.1612500)  In[52]: coords = [(24.950899, 60.169158), (24.953492, 60.169158), (24.953510, 60.170104), (24.950958, 60.169990)]  In[53]: poly = Polygon(coords)    In[54]: poly.contains(p1)  Out[54]: True    In[55]: p1.within(poly)  Out[55]: True    In[56]: poly.contains(p2)  Out[56]: False
  • 判斷多邊形的集合操作
In[57]: poly2 = Polygon([(23.2154, 59.1156), (24.83151, 59.41516), (25.11667, 60.311561), (24.16178, 60.13315)] )    In[58]: poly2.intersects(poly)  Out[58]: True    In[59]: poly2.contains(poly)  Out[59]: True    In[60]: poly.contains(poly2)  Out[60]: False
  • .contains：判斷polygon1是否包含polygon2
  • .intersects：判斷polygon1和polygon2是否重疊
  • .intersections ：返回兩個polygon重疊的部分

參考鏈接：

1. https://stackoverflow.com/questions/36399381/whats-the-fastest-way-of-checking-if-a-point-is-inside-a-polygon-in-python

2. https://gis.stackexchange.com/questions/90055/finding-if-two-polygons-intersect-in-python

3. https://automating-gis-processes.github.io/CSC18/lessons/L4/point-in-polygon.html

VirMach 便宜 VPS

QNews

QNews