Python+GDAL面数据中心点提取为单独shape文件，并复制属性数据

最近需要使用GDAL实现shape面数据批量将中心点提取为点shape文件，并复制所有的属性信息。整体的逻辑不是非常复杂，但是在网上现成的代码也找不到。

实现思路

具体思路如下：

#mermaid-svg-8m1tffzStNM4t1hK {font-family:"trebuchet ms",verdana,arial,sans-serif;font-size:16px;fill:#333;}#mermaid-svg-8m1tffzStNM4t1hK .error-icon{fill:#552222;}#mermaid-svg-8m1tffzStNM4t1hK .error-text{fill:#552222;stroke:#552222;}#mermaid-svg-8m1tffzStNM4t1hK .edge-thickness-normal{stroke-width:2px;}#mermaid-svg-8m1tffzStNM4t1hK .edge-thickness-thick{stroke-width:3.5px;}#mermaid-svg-8m1tffzStNM4t1hK .edge-pattern-solid{stroke-dasharray:0;}#mermaid-svg-8m1tffzStNM4t1hK .edge-pattern-dashed{stroke-dasharray:3;}#mermaid-svg-8m1tffzStNM4t1hK .edge-pattern-dotted{stroke-dasharray:2;}#mermaid-svg-8m1tffzStNM4t1hK .marker{fill:#333333;stroke:#333333;}#mermaid-svg-8m1tffzStNM4t1hK .marker.cross{stroke:#333333;}#mermaid-svg-8m1tffzStNM4t1hK svg{font-family:"trebuchet ms",verdana,arial,sans-serif;font-size:16px;}#mermaid-svg-8m1tffzStNM4t1hK .label{font-family:"trebuchet ms",verdana,arial,sans-serif;color:#333;}#mermaid-svg-8m1tffzStNM4t1hK .cluster-label text{fill:#333;}#mermaid-svg-8m1tffzStNM4t1hK .cluster-label span{color:#333;}#mermaid-svg-8m1tffzStNM4t1hK .label text,#mermaid-svg-8m1tffzStNM4t1hK span{fill:#333;color:#333;}#mermaid-svg-8m1tffzStNM4t1hK .node rect,#mermaid-svg-8m1tffzStNM4t1hK .node circle,#mermaid-svg-8m1tffzStNM4t1hK .node ellipse,#mermaid-svg-8m1tffzStNM4t1hK .node polygon,#mermaid-svg-8m1tffzStNM4t1hK .node path{fill:#ECECFF;stroke:#9370DB;stroke-width:1px;}#mermaid-svg-8m1tffzStNM4t1hK .node .label{text-align:center;}#mermaid-svg-8m1tffzStNM4t1hK .node.clickable{cursor:pointer;}#mermaid-svg-8m1tffzStNM4t1hK .arrowheadPath{fill:#333333;}#mermaid-svg-8m1tffzStNM4t1hK .edgePath .path{stroke:#333333;stroke-width:2.0px;}#mermaid-svg-8m1tffzStNM4t1hK .flowchart-link{stroke:#333333;fill:none;}#mermaid-svg-8m1tffzStNM4t1hK .edgeLabel{background-color:#e8e8e8;text-align:center;}#mermaid-svg-8m1tffzStNM4t1hK .edgeLabel rect{opacity:0.5;background-color:#e8e8e8;fill:#e8e8e8;}#mermaid-svg-8m1tffzStNM4t1hK .cluster rect{fill:#ffffde;stroke:#aaaa33;stroke-width:1px;}#mermaid-svg-8m1tffzStNM4t1hK .cluster text{fill:#333;}#mermaid-svg-8m1tffzStNM4t1hK .cluster span{color:#333;}#mermaid-svg-8m1tffzStNM4t1hK div.mermaidTooltip{position:absolute;text-align:center;max-width:200px;padding:2px;font-family:"trebuchet ms",verdana,arial,sans-serif;font-size:12px;background:hsl(80, 100%, 96.2745098039%);border:1px solid #aaaa33;border-radius:2px;pointer-events:none;z-index:100;}#mermaid-svg-8m1tffzStNM4t1hK :root{--mermaid-font-family:"trebuchet ms",verdana,arial,sans-serif;}

面shp文件

坐标系信息提取

字段信息提取

创建点shp文件

遍历要素

字段值提取

中心点提取

写入新shp文件

实现代码

from osgeo import ogr
from osgeo import gdaldef area_to_point(in_path,out_path):# 解决中文字符问题# 为了支持中文路径，请添加下面这句代码gdal.SetConfigOption("GDAL_FILENAME_IS_UTF8", "YES")# 为了使属性表字段支持中文，请添加下面这句gdal.SetConfigOption("SHAPE_ENCODING", "CP936")# 注册所有的驱动ogr.RegisterAll()ds = ogr.Open(in_path, 0)if ds is None:  # 打开失败print("打开失败")# 获取该数据源中的图层个数，一般shp数据图层只有一个，如果是mdb、dxf等图层就会有多个layer_count = ds.GetLayerCount()if layer_count != 1:print("图层数量异常")driver = ogr.GetDriverByName("ESRI Shapefile")out_ds = driver.CreateDataSource(out_path)in_lyr = ds.GetLayer()out_lyr = out_ds.CreateLayer('Point', in_lyr.GetSpatialRef(), geom_type=ogr.wkbPoint)# 获取文件字段属性in_feature = in_lyr.GetNextFeature()# 复制全部的字段in_feature_defn = in_lyr.GetLayerDefn()field_names = []for i in range(0, in_feature_defn.GetFieldCount()):field_defn = in_feature_defn.GetFieldDefn(i)field_name = field_defn.namefield_names.append(field_name)# 输出layer创建字段out_lyr.CreateField(field_defn)# 获取输出图层属性表信息out_feature_defn = out_lyr.GetLayerDefn()feature = in_lyr.GetFeature(0)# 遍历要素while feature:new_feature = ogr.Feature(out_feature_defn)# 添加点要素geom = feature.GetGeometryRef()center = geom.Centroid()new_feature.SetGeometry(center)# 复制所有字段属性for field_name in field_names:area = feature.GetField(field_name)# 添加点的字段值new_feature.SetField(field_name, area)# 添加要素到图层out_lyr.CreateFeature(new_feature)feature = in_lyr.GetNextFeature()out_ds.Destroy()if __name__ == '__main__':shp_path = "E:/测试数据/DK4107262017.shp"out_path = "E:/测试数据/out.shp"area_to_point(shp_path,out_path)

提取效果