GDAL读写矢量文件——Python
在Python中使用OGR时，先要导入OGR库，如果需要对中文的支持，还需要导入GDAL库，具体代码如下。Python创建的shp结果如图1所示。
图1 Python创建矢量结果
-
#-*- coding: cp936 -*- -
try: -
from osgeo import gdal -
from osgeo import ogr -
exceptImportError: -
import gdal -
import ogr
1．读取矢量
-
#-*- coding: cp936 -*- -
try: -
from osgeo import gdal -
from osgeo import ogr -
exceptImportError: -
import gdal -
import ogr -
defReadVectorFile(): -
# 为了支持中文路径，请添加下面这句代码 -
gdal.SetConfigOption("GDAL_FILENAME_IS_UTF8","NO") -
# 为了使属性表字段支持中文，请添加下面这句 -
gdal.SetConfigOption("SHAPE_ENCODING","") -
strVectorFile ="E:\\Datum\\GDALCsTest\\Debug\\beijing.shp" -
# 注册所有的驱动 -
ogr.RegisterAll() -
#打开数据 -
ds = ogr.Open(strVectorFile, 0) -
if ds == None: -
print("打开文件【%s】失败！", strVectorFile) -
return -
print("打开文件【%s】成功！", strVectorFile) -
# 获取该数据源中的图层个数，一般shp数据图层只有一个，如果是mdb、dxf等图层就会有多个 -
iLayerCount = ds.GetLayerCount() -
# 获取第一个图层 -
oLayer = ds.GetLayerByIndex(0) -
if oLayer == None: -
print("获取第%d个图层失败！\n", 0) -
return -
# 对图层进行初始化，如果对图层进行了过滤操作，执行这句后，之前的过滤全部清空 -
oLayer.ResetReading() -
# 通过属性表的SQL语句对图层中的要素进行筛选，这部分详细参考SQL查询章节内容 -
oLayer.SetAttributeFilter("\"NAME99\"LIKE \"北京市市辖区\"") -
# 通过指定的几何对象对图层中的要素进行筛选 -
#oLayer.SetSpatialFilter() -
# 通过指定的四至范围对图层中的要素进行筛选 -
#oLayer.SetSpatialFilterRect() -
# 获取图层中的属性表表头并输出 -
print("属性表结构信息：") -
oDefn = oLayer.GetLayerDefn() -
iFieldCount = oDefn.GetFieldCount() -
for iAttr in range(iFieldCount): -
oField =oDefn.GetFieldDefn(iAttr) -
print( "%s: %s(%d.%d)" % ( \ -
oField.GetNameRef(),\ -
oField.GetFieldTypeName(oField.GetType() ), \ -
oField.GetWidth(),\ -
oField.GetPrecision())) -
# 输出图层中的要素个数 -
print("要素个数 = %d", oLayer.GetFeatureCount(0)) -
oFeature = oLayer.GetNextFeature() -
# 下面开始遍历图层中的要素 -
while oFeature is not None: -
print("当前处理第%d个: \n属性值：", oFeature.GetFID()) -
# 获取要素中的属性表内容 -
for iField inrange(iFieldCount): -
oFieldDefn =oDefn.GetFieldDefn(iField) -
line = " %s (%s) = " % ( \ -
oFieldDefn.GetNameRef(),\ -
ogr.GetFieldTypeName(oFieldDefn.GetType())) -
ifoFeature.IsFieldSet( iField ): -
line = line+ "%s" % (oFeature.GetFieldAsString( iField ) ) -
else: -
line = line+ "(null)" -
print(line) -
# 获取要素中的几何体 -
oGeometry =oFeature.GetGeometryRef() -
# 为了演示，只输出一个要素信息 -
break -
print("数据集关闭！")
执行上面的代码，如果不设置属性过滤，输出内容如图3‑9上半部分所示，如过设置了属性过滤，输出内容如图3‑9下半部分所示。（Python输出的中文转为编码了）。
图2 OGR库使用Python读取矢量示例
2．写入矢量
在使用Python创建矢量图形的时候，使用的WKT格式的字符串来进行创建。也可以使用其他的方式进行创建。代码如下，写出来的矢量图形和属性表如图1所示。
-
#-*- coding: cp936 -*- -
try: -
from osgeo import gdal -
from osgeo import ogr -
exceptImportError: -
import gdal -
import ogr -
defWriteVectorFile(): -
# 为了支持中文路径，请添加下面这句代码 -
gdal.SetConfigOption("GDAL_FILENAME_IS_UTF8","NO") -
# 为了使属性表字段支持中文，请添加下面这句 -
gdal.SetConfigOption("SHAPE_ENCODING","") -
strVectorFile ="E:\\TestPolygon.shp" -
# 注册所有的驱动 -
ogr.RegisterAll() -
# 创建数据，这里以创建ESRI的shp文件为例 -
strDriverName = "ESRIShapefile" -
oDriver =ogr.GetDriverByName(strDriverName) -
if oDriver == None: -
print("%s 驱动不可用！\n", strDriverName) -
return -
# 创建数据源 -
oDS =oDriver.CreateDataSource(strVectorFile) -
if oDS == None: -
print("创建文件【%s】失败！", strVectorFile) -
return -
# 创建图层，创建一个多边形图层，这里没有指定空间参考，如果需要的话，需要在这里进行指定 -
papszLCO = [] -
oLayer =oDS.CreateLayer("TestPolygon", None, ogr.wkbPolygon, papszLCO) -
if oLayer == None: -
print("图层创建失败！\n") -
return -
# 下面创建属性表 -
# 先创建一个叫FieldID的整型属性 -
oFieldID =ogr.FieldDefn("FieldID", ogr.OFTInteger) -
oLayer.CreateField(oFieldID, 1) -
# 再创建一个叫FeatureName的字符型属性，字符长度为50 -
oFieldName =ogr.FieldDefn("FieldName", ogr.OFTString) -
oFieldName.SetWidth(100) -
oLayer.CreateField(oFieldName, 1) -
oDefn = oLayer.GetLayerDefn() -
# 创建三角形要素 -
oFeatureTriangle = ogr.Feature(oDefn) -
oFeatureTriangle.SetField(0, 0) -
oFeatureTriangle.SetField(1, "三角形") -
geomTriangle =ogr.CreateGeometryFromWkt("POLYGON ((0 0,20 0,10 15,0 0))") -
oFeatureTriangle.SetGeometry(geomTriangle) -
oLayer.CreateFeature(oFeatureTriangle) -
# 创建矩形要素 -
oFeatureRectangle = ogr.Feature(oDefn) -
oFeatureRectangle.SetField(0, 1) -
oFeatureRectangle.SetField(1, "矩形") -
geomRectangle =ogr.CreateGeometryFromWkt("POLYGON ((30 0,60 0,60 30,30 30,30 0))") -
oFeatureRectangle.SetGeometry(geomRectangle) -
oLayer.CreateFeature(oFeatureRectangle) -
# 创建五角形要素 -
oFeaturePentagon = ogr.Feature(oDefn) -
oFeaturePentagon.SetField(0, 2) -
oFeaturePentagon.SetField(1, "五角形") -
geomPentagon =ogr.CreateGeometryFromWkt("POLYGON ((70 0,85 0,90 15,80 30,65 15,700))") -
oFeaturePentagon.SetGeometry(geomPentagon) -
oLayer.CreateFeature(oFeaturePentagon) -
oDS.Destroy() -
print("数据集创建完成！\n")
3．矢量数据管理
-
defVectorDelete(strVectorFile): -
# 注册所有的驱动 -
ogr.RegisterAll() -
oDriver = None -
#打开矢量 -
oDS = ogr.Open(strVectorFile, 0) -
if oDS == None: -
os.remove(strVectorFile) -
return -
oDriver = oDS.GetDriver() -
if oDriver == None: -
os.remove(strVectorFile) -
return -
ifoDriver.DeleteDataSource(strVectorFile) == ogr.OGRERR_NONE: -
return -
else: -
os.remove(strVectorFile) -
defVectorRename(strOldFile, strNewFile): -
# 注册所有的驱动 -
ogr.RegisterAll() -
oDriver = None -
#打开矢量 -
oDS = Ogr.Open(strOldFile, 0) -
if oDS == None : -
os.rename(strOldFile,strNewFile) -
return -
oDriver = oDS.GetDriver() -
if oDriver == None: -
os.rename(strOldFile,strNewFile) -
return -
oDDS = oDriver.CopyDataSource(oDS,strNewFile, None) -
if oDDS == None: -
os.rename(strOldFile,strNewFile) -
if oDriver.DeleteDataSource(strOldFile)== ogr.OGRERR_NONE: -
return -
else : -
os.rename(strOldFile,strNewFile)