Python地圖柵格化範例
2022-10-13 21:04:24
Python地圖柵格化範例
引言
shapefile是GIS中的一種非常重要的資料型別，由ESRI開發的空間資料開放格式，目前該資料格式已經成為了GIS領域的開放標準。目前絕大多數開源以及收費的GIS軟體都支援該資料型別。事實上，shapefile檔案指的一種檔案儲存的方法，實際上該種檔案是由多個檔案組成的。組成一個shapefile有三種檔案必不可少, '.shp','.shx','.dbf'檔案。
geopandas對shapefile提供了很好的讀取與寫出支援。geopandas庫允許對幾何型別進行空間操作，其中DataFrame結構相當於GIS資料中的一張屬性表，使得可以直接操作向量資料屬性表，使得python中操作地理資料更加方便。本範例通過geopandas實現對地理資料的操作。
開發準備
由於geopandas庫的安裝需要一些前提庫，因此需要先安裝一些庫
pip install pipwin
pipwin install gdal
pipwin install fiona
pip install geopandas
實測以上方法可以成功在windows下安裝（注：如果在Anaconda下安裝geopandas更為方便）
資料準備
該資料是一段GPS掃描資料，包含經緯度。
程式碼範例
環境引入
import geopandas as gp
import matplotlib.pyplot as plt
from shapely import geometry
import math
GPS資料處理
lake_original_path = 'data.txt'
lake_original_data = ''
lake_points = []
# 讀取檔案
with open(lake_original_path) as f:
lake_original_data = f.read()
# 處理經緯度座標 並以Point的形式新增到list中
for xy in lake_original_data.split(';'):
x, _, y = xy.partition(',')
x = float(x.strip()) / 100
y = float(y.strip()) / 100
lake_points.append(geometry.Point(y, x))
建立要素
# 建立線狀要素
lake_line = geometry.LineString(lake_points)
# crs指定座標系
lake_ = gp.GeoSeries(lake_line, crs='EPSG:4326')
# 儲存shp檔案
lake_.to_file("boundary.shp", driver='ESRI Shapefile', encoding='utf-8')
# 記錄邊界條件 用於構建柵格
x_min, y_min, x_max, y_max = lake_line.bounds[:4]
# 繪圖
lake_.plot()
plt.show()
構建柵格
# 柵格大小
GRID_WIDTH = 0.009 * 2 / 100
grid_rows_num = int(math.ceil((y_max - y_min) / float(GRID_WIDTH)))
grid_columns_num = int(math.ceil((x_max - x_min) / float(GRID_WIDTH)))
grids = []
for r in range(grid_rows_num):
for c in range(grid_columns_num):
grid_4coords = []
# 左上角
x_lt = x_min + c * GRID_WIDTH
y_lt = y_max - r * GRID_WIDTH
# 右上角
x_rt = x_lt + GRID_WIDTH
y_rt = y_lt
# 左下角
x_lb = x_lt
y_lb = y_lt - GRID_WIDTH
# 右下角
x_rb = x_rt
y_rb = y_lb
# 兩個三角形拼接一個柵格
grid_4coords.append(geometry.Point(x_lt,y_lt))
grid_4coords.append(geometry.Point(x_rt,y_rt))
grid_4coords.append(geometry.Point(x_rb,y_rb))
grid_4coords.append(geometry.Point(x_lb,y_lb))
grid_4coords.append(geometry.Point(x_lt,y_lt))
# 建立一個網格
grids.append(geometry.LineString(grid_4coords))
grid_ = gp.GeoSeries(grids)
grid_.to_file('E:\\just\\海韻湖智慧技術實驗場\\data\\grids.shp',driver='ESRI Shapefile', encoding='utf-8')
grid_.plot()
plt.show()
要素疊加
# 要素疊加
elements = [lake_line]
elements += grids
elements_ = gp.GeoSeries(elements)
elements_.to_file('elements.shp', driver='ESRI Shapefile', encoding='utf-8')
elements_.plot()
plt.show()
參考連結
python-geopandas讀取、建立shapefile檔案、geopandas學習教學