Python读写修改Shapefile

说明

当我们需要借助编程手段对 Shapefile 进行操作时，会想到用 ArcGIS Engine、ArcGIS Runtime SDK或开源的 GDAL。在没有一定基础的情况下，这些弄起来属实挺麻烦(配环境、写程序等)。当我们做些简单的 Shapefile 数据处理，有没有好办法呢？这里非常推荐使用 Python 语言，介绍三个常用的开发库：

ArcPy
什么是快… ArcPy，下面是文档地址。
https://desktop.arcgis.com/zh-cn/arcmap/10.5/analyze/arcpy/what-is-arcpy-.htm
简言之就是 ArcGIS 的配套工具，好用，但是有限定条件：使用它的前提是安装好 ArcGIS 的环境，另外无论是运行还是打包后执行，都将永久依赖 ArcGIS 。
pyshp
由纯 Python 打造的读、写 Shapefile 的轻量级开源工具，Github 地址。
https://github.com/GeospatialPython/pyshp
如果你只做 Shapefile 读写的话，用它是非常合适的。
GDAL for Python
嗷它还是 GDAL，但是显然 Python 版本更适合于小型项目(可能叫工具更合适)，它比 C++版本的好配太多，也更省去源码编译的过程。功能强大且开源，选它选它。

1. GDAL安装

Python环境就不说了，先到这里下载 GDAL‑<gdal版本>‑cp<python版本>‑<平台架构>.whl 安装包，
https://www.lfd.uci.edu/~gohlke/pythonlibs/#gdal
下载完成后安装

pip install xxx.whl

安装后试一下能否 import 进来包。

D:\test>python
Python 3.7.6 (default, Jan  8 2020, 20:23:39) [MSC v.1916 64 bit (AMD64)] :: Anaconda, Inc. on win32>>> from osgeo import gdal
>>>

没有问题。

2. GDAL读取shp

#!/usr/bin/python
# -*- coding: UTF-8 -*-
from osgeo import gdal
from osgeo import ogr
from osgeo import osrdef ReadShp(_filename):# 支持中文路径gdal.SetConfigOption("GDAL_FILENAME_IS_UTF8", "YES")# 支持中文编码gdal.SetConfigOption("SHAPE_ENCODING", "UTF-8")# 注册所有的驱动ogr.RegisterAll()# 打开数据ds = ogr.Open(_filename, 0)if ds == None:return ("打开文件失败！")# 获取数据源中的图层个数，shp数据图层只有一个，gdb、dxf会有多个# iLayerCount = ds.GetLayerCount()# 获取第一个图层oLayer = ds.GetLayerByIndex(0)if oLayer == None:return ("获取图层失败！")# 对图层进行初始化oLayer.ResetReading()# 输出图层中的要素个数num = oLayer.GetFeatureCount(0)result_list = []# 获取要素for i in range(0, num):ofeature = oLayer.GetFeature(i)# shp 文件的 NAME 字段fieldName = ofeature.GetFieldAsString('NAME')# geom = str(ofeature.GetGeometryRef())result_list.append(fieldName)ds.Destroy()del dsreturn result_listif __name__ == '__main__':shpFile = "TM_WORLD_BORDERS-0.3.shp"nameList = ReadShp(shpFile)for name in nameList:print (name)

3. GDAL创建shp

#!/usr/bin/python
# -*- coding: UTF-8 -*-from osgeo import gdal
from osgeo import ogr
from osgeo import osrdef NewShp(_strDriverName):# 支持中文路径gdal.SetConfigOption("GDAL_FILENAME_IS_UTF8", "YES")# 属性表字段支持中文gdal.SetConfigOption("SHAPE_ENCODING", "UTF-8")# 注册驱动ogr.RegisterAll()# 创建shp数据strDriverName = "ESRI Shapefile"oDriver = ogr.GetDriverByName(strDriverName)if oDriver == None:return ("驱动不可用："+strDriverName)# 创建数据源oDS = oDriver.CreateDataSource(_strDriverName)if oDS == None:return ("创建文件失败"+_strDriverName)# 创建一个多边形图层，指定坐标系为WGS84papszLCO = []geosrs = osr.SpatialReference()geosrs.SetWellKnownGeogCS("WGS84")# 线：ogr_type = ogr.wkbLineString# 点：ogr_type = ogr.wkbPointogr_type = ogr.wkbPolygon# 面的类型为Polygon，线的类型为Polyline，点的类型为PointoLayer = oDS.CreateLayer("Polygon", geosrs, ogr_type, papszLCO)if oLayer == None:return ("图层创建失败！")# 创建属性表# 创建id字段oId = ogr.FieldDefn("id", ogr.OFTInteger)oLayer.CreateField(oId, 1)# 创建name字段oName = ogr.FieldDefn("name", ogr.OFTString)oLayer.CreateField(oName, 1)oDefn = oLayer.GetLayerDefn()# 创建要素# 数据集# wkt_geom id namefeatures = ['test0;POLYGON((-1.58 0.53, -0.79 0.55, -0.79 -0.23, -1.57 -0.25, -1.58 0.53))', 'test1;POLYGON((-1.58 0.53, -0.79 0.55, -0.79 -0.23, -1.57 -0.25, -1.58 0.53))']for index, f in enumerate(features):oFeaturePolygon = ogr.Feature(oDefn)oFeaturePolygon.SetField("id",index)oFeaturePolygon.SetField("name",f.split(";")[0])geomPolygon = ogr.CreateGeometryFromWkt(f.split(";")[1])oFeaturePolygon.SetGeometry(geomPolygon)oLayer.CreateFeature(oFeaturePolygon)# 创建完成后，关闭进程oDS.Destroy()print ("数据集创建完成！")if __name__ == '__main__':shpFile = "TestCreateShp.shp"NewShp(shpFile)

4. GDAL修改shp

#!/usr/bin/python
# -*- coding: UTF-8 -*-
from osgeo import gdal
from osgeo import ogr
from osgeo import osr
def UpdateShp(_filename):# 支持中文路径gdal.SetConfigOption("GDAL_FILENAME_IS_UTF8", "YES")# 属性表字段支持中文gdal.SetConfigOption("SHAPE_ENCODING", "UTF-8")# 注册驱动ogr.RegisterAll()strDriverName = "ESRI Shapefile"oDriver = ogr.GetDriverByName(strDriverName)if oDriver == None:warninglist.append("驱动不可用："+strDriverName)return ("驱动不可用："+strDriverName)# 打开数据oDS = oDriver.Open(_filename, 1)oLayer = oDS.GetLayer()# 编辑属性字段for feature in oLayer: feature.SetField('NAME', "EditName")oLayer.SetFeature(feature)feature.Destroy()oDS.Destroy()print ("编辑完成！")if __name__ == '__main__':shpFile = "TM_WORLD_BORDERS-0.3.shp"UpdateShp(shpFile)

5. 参考

本文参考了

[1]. https://cloud.tencent.com/developer/article/1743325
[2]. https://blog.csdn.net/sinat_41310868/article/details/112722839
[3]. https://www.cnblogs.com/weihongli/p/7843451.html

6. 相关资源

本文用到的示例数据和上贴出的代码，没什么其他特别的东西。

https://download.csdn.net/download/ShyLoneGirl/16752858

.
.
.
.
.
.

桃花仙人种桃树，又摘桃花换酒钱_