以下程序包括两个功能:
基于shapefile建立缓冲区
把单个shapefile所有的要素都合并在一起
# -*- coding: utf-8 -*-
"""
Created on Tue Feb 28 11:38:39 2023
@author: Asus
"""
from pathlib import Path
from osgeo import ogr, osr
import geopandas as gpd
def buffer(inShp, fname, bdistance=10000):
jingkai_polygon = gpd.read_file(inShp)
jingkai_polygon = jingkai_polygon.to_crs("EPSG:3413")
#地理坐标系转投影坐标系,地理坐标意味着distance单位为度,而投影坐标单位为米
jingkai_polygon=jingkai_polygon.buffer(bdistance)
jingkai_polygon.to_file(fname,
driver='ESRI Shapefile',
encoding='utf-8')
def uni(shpPath, fname):
"""
:param shpPath: 输入的矢量路径
:param fname: 输出的矢量路径
:return:
"""
driver = ogr.GetDriverByName("ESRI Shapefile")
dataSource = driver.Open(shpPath, 1)
layer = dataSource.GetLayer()
# 新建DataSource,Layer
out_ds = driver.CreateDataSource(fname)
out_lyr = out_ds.CreateLayer(fname, layer.GetSpatialRef(), ogr.wkbPolygon)
def_feature = out_lyr.GetLayerDefn()
# 遍历原始的Shapefile文件给每个Geometry做Buffer操作
# current_union = layer[0].Clone()
print('the length of layer:', len(layer))
if len(layer) == 0:
return
for i, feature in enumerate(layer):
geometry = feature.GetGeometryRef()
if i == 0:
current_union = geometry.Clone()
current_union = current_union.Union(geometry).Clone()
if i == len(layer) - 1:
out_feature = ogr.Feature(def_feature)
out_feature.SetGeometry(current_union)
out_lyr.ResetReading()
out_lyr.CreateFeature(out_feature)
if __name__ == '__main__':
inShp = 'G:\mask\simplified-land-polygons-complete-3857\simplified_land_polygons.shp'
fname = 'G:\mask\simple\simple_land_10kmBuffer2.shp'
finame='G:\mask\simple\simple_land_10kmBuffer2_merge.shp'
# buffer(inShp, fname)
uni(fname,finame)