常见的矢量数据格式有Shapefile、GeoJSON、CSV，及文件数据库gdb和空间数据库PostGIS，不论是何种格式的数据或如何存储，一旦打开数据源、获取矢量图层后（详情参考OGR操作矢量数据的类结构图），对数据的操作都一样。下面对矢量数据的读写进行详细的介绍。

    一、打开不同的矢量文件

    1、定义打开数据源的函数，并遍历所有的图层，输出他们的名字和图层

'''
    输出数据源中的图层
    参数：fn 数据源的路径
         is_write 打开数据源的模式，0 表示只读模式，1 表示读写模式
'''
def print_layers(fn, is_write):
    ds = ogr.Open(fn, is_write)
    if ds is None:
        #raise OSError('Could not open %s', fn)
        raise OSError('Could not open {}'.format(fn))
    for i in range(ds.GetLayerCount()):
        lyr = ds.GetLayer(i)
        print('{0}:{1}'.format(i, lyr.GetName()))

    2、打开PostGIS数据源

    利用上面打开数据源并输出图层的方法，打开本地的Post GIS数据库 postgis_24_sample，输出图层如下：

#打开postgis数据库
    print_layers('PG:user=postgres password=postgis dbname=postgis_24_sample', 0)

0:tiger.county
1:tiger.state
2:tiger.place
3:tiger.cousub
4:tiger.edges
5:tiger.addrfeat
6:tiger.faces
7:tiger.zcta5
8:tiger.tract
9:tiger.tabblock
10:tiger.bg

    3、打开文件夹数据源（shapefiles和CSV)

    打开shapefiles文件夹，输出的结果如下所示：

#打开shapefile文件夹
    shp_fn = r'E:\HZZ_GXDP_Spatial_Data'
    print_layers(shp_fn, 0)

0:CommunityRiver
1:CountyRiver
2:StreetRiver
3:Video
4:WRY

    二、读取矢量数据

        读取矢量数据的步骤包括：

        1、打开数据源；

        2、获取图层；

        3、获取图层中的要素；

        4、获取要素的属性和几何信息。

        下面是读取矢量数据的代码片段：

import sys
from osgeo import ogr
import ospybook as pb

fn = r'D:\soft\geoprocessing-with-python\china_basic_map'
ds = ogr.Open(fn, 0)
if ds is None:
    sys.exit('Could not open {0}.'.format(fn))

'''索引获取数据源中的图层方式'''
lyr = ds.GetLayer(0)
print(lyr.GetName())

i = 0
for fea in lyr:
    pt = fea.geometry()
    x = pt.GetX()
    y = pt.GetY()
    '''字段名称 获取属性值'''
    code = fea.GetField('GBCODE')
    '''对象方式获取属性值'''
    length = fea.LENGTH
    '''对象方式获取属性值2'''
    lpoly_ = fea['LPOLY_']
    '''索引方式获取属性值'''
    fnode_ = fea.GetField(0)
    '''获取特定类型的属性值'''
    str_len = fea.GetFieldAsString('LENGTH')
    print(code, length, fnode_, lpoly_, str_len, x, y)
    i += 1
    if i == 20:
        break

'''名称获取数据源中的图层'''
lyr2 = ds.GetLayer('国家')
print(lyr2.GetName())

    三、获取数据元数据

        矢量数据的元数据，包括数据范围、几何类型、空间参考及图层对象的概要信息等，下面的代码片段展示如何获取这些元数据信息：

import sys
from osgeo import  ogr

fn = r'D:\soft\geoprocessing-with-python\china_basic_map'
ds = ogr.Open(fn, 0)
if ds is None:
    sys.exit('Could not open {0}.'.format(fn))

lyr = ds.GetLayer(0)
'''获取图层范围'''
extent = lyr.GetExtent()
print(extent)

'''(73.44696044921875, 135.08583068847656, 3.408477306365967, 53.557926177978516)'''
'''(min_x, max_x, min_y, max_y)'''

'''获取图层集合类型，返回值为数字'''
'''1:点， 2：线， 3：面'''
geom_type = lyr.GetGeomType()
print(geom_type) #2
print(geom_type == ogr.wkbPoint) #False
print(geom_type == ogr.wkbLineString) #True
print(geom_type == ogr.wkbPolygon) #False

'''通过要素获取几何类型'''
fea = lyr.GetFeature(0)
print(fea.geometry().GetGeometryType()) #2
print(fea.geometry().GetGeometryName()) #LINESTRING

'''获取图层的空间参考'''
print(lyr.GetSpatialRef())
# GEOGCS["GCS_WGS_1984",
#     DATUM["WGS_1984",
#         SPHEROID["WGS_84",6378137.0,298.257223563]],
#     PRIMEM["Greenwich",0.0],
#     UNIT["Degree",0.0174532925199433],
#     AUTHORITY["EPSG","4326"]]

print(lyr.schema)

'''获取图层属性名称及数据类型'''
for field in lyr.schema:
    print(field.name, field.GetTypeName())
# FNODE_ Integer64
# TNODE_ Integer64
# LPOLY_ Integer64
# RPOLY_ Integer64
# LENGTH Real
# BOU2_4M_ Integer64
# BOU2_4M_ID Integer64
# GBCODE Integer

    四、矢量数据写入

        1、矢量数据写入的思路如下：

        （1）以读写模式打开（创建）数据源；

        （2）获取待添加要素的图层（或新创建图层以添加要素）；

        （3）创建空要素，并为几何对象和属性赋值；

        （4）将创建的要素插入图层中。

        下面是利用一个图层中的要素创建一个新图层的详细代码：

# -*- coding:utf-8 -*-
import sys
from osgeo import ogr

'''创建一个图层（根据一个图层中的要素创建）'''

fn = r'D:\soft\geoprocessing-with-python\china_basic_map'
'''以读写模式打开数据源'''
ds = ogr.Open(fn, 1)
if ds is None:
    sys.exit('Could not open {0}.'.format(fn))

in_lyr = ds.GetLayer('省会城市')

# lyr = ds.GetLayer('国家')
#
# from ospybook.vectorplotter import VectorPlotter
# vp = VectorPlotter(False)
# vp.plot(lyr, fill=False)
# vp.draw()

'''若存在同名图层，则先删除再创键'''
if ds.GetLayer('capital_city'):
    ds.DeleteLayer('capital_city')
out_lyr = ds.CreateLayer('capital_city',
                         in_lyr.GetSpatialRef(),
                         ogr.wkbPoint)

'''根据元图层字段对象创建新图层的字段'''
out_lyr.CreateFields(in_lyr.schema)

'''根据图层定义创建一个空的要素'''
out_defn = out_lyr.GetLayerDefn()
out_fea = ogr.Feature(out_defn)
for in_fea in in_lyr:
    geom = in_fea.geometry()
    out_fea.SetGeometry(geom)
    for i in range(in_fea.GetFieldCount()):
        value = in_fea.GetField(i)
        if i != 5:
            out_fea.SetField(i, value)
        else:
            print(value)
    out_lyr.CreateFeature(out_fea)

        2、创建新的数据源

        创建新数据源的关键在使用正确的驱动程序，每种驱动程序只处理操作一种类型的矢量数据。有两种方式获取驱动程序：（1）从一个已经打开的数据集中获取，这将允许创建一个新的数据源，它和已存在的数据源矢量数据格式一样；（2）使用OGR中的GetDriverByName函数，传递给它驱动程序的简称，驱动程序的名称在OGR网站有介绍。下面是代码片段：

''''''''''创建新的数据源'''''''''
import sys
from osgeo import ogr

fn = r'D:\soft\geoprocessing-with-python\china_basic_map'

'''获取数据驱动的第一种方式：根据已有数据源获取'''
ds = ogr.Open(fn, 0)
if ds is None:
    sys.exit('Could not open {0}.'.format(fn))

driver = ds.GetDriver()
print(driver.GetName())

'''获取数据源的第二种方式：GetDriverByName'''
json_driver = ogr.GetDriverByName('GeoJSON')
print(json_driver.GetName())

'''创建GeoJson数据源，数据源路径应到具体文件名'''
json_fn = r'D:\soft\geoprocessing-with-python\china_basic_map\json_fn.json'
json_ds = json_driver.CreateDataSource(json_fn)
if json_ds is None:
    sys.exit('Could not create {0}'.format(json_fn))
print(json_ds)

        3、新建属性字段

        要将一个字段添加到图层中，需要一个包含字段名称、数据类型、字段和精度等重要信息的FieldDefn对象，下面是新建属性字段的代码片段：

'''UseExceptions'''
ogr.UseExceptions()
json_fn = r'D:\soft\geoprocessing-with-python\china_basic_map\json_fn4.json'
json_driver = ogr.GetDriverByName('GeoJSON')
print('start')

try:
    json_ds = json_driver.CreateDataSource(json_fn)
    '''新建属性数据字段'''
    lyr = json_ds.CreateLayer('layer')
    coord_fld = ogr.FieldDefn('X', ogr.OFTReal)
    coord_fld.SetWidth(8)
    coord_fld.SetPrecision(3)
    lyr.CreateField(coord_fld)
    coord_fld.SetName('Y')
    lyr.CreateField(coord_fld)

    #新建要素
    fea = ogr.Feature(lyr.GetLayerDefn())
    fea.SetField('X', 12.34)
    fea.SetFID(1)
    lyr.CreateFeature(fea)
    #结果同步至硬盘（保存到文件中）
    json_ds.SyncToDisk()
except RuntimeError as e:
    print(e)

print('end')

    五、更新现有数据

        在进行矢量数据处理时，有时需要更新现有数据，而不是创建一个全新的数据集。是否能够更新以及支持的编辑操作，取决于数据格式。常见的更新数据操作包括 改变图层定义（属性字段）、要素添加、更新和删除等。下面分别进行说明：

        1、改变图层定义

        包括修改属性字段的定义、新增属性字段及删除属性字段，具体代码片段如下：

import sys
from osgeo import ogr
'''''''更新现有数据'''''''
fn = r'D:\soft\geoprocessing-with-python\china_basic_map'
ds = ogr.Open(fn, 1)
if ds is None:
    sys.exit('Could not open {0}.'.format(fn))

'''更改图层定义(新增、删除、修改字段)'''
lyr = ds.GetLayer(0)
lyrdefn = lyr.GetLayerDefn()
i = lyrdefn.GetFieldIndex('GBCODE')
'''获取字段的数据类型'''
fld_type = lyrdefn.GetFieldDefn(i).GetType()
print(fld_type)
'''定义新的属性字段'''
fld_defn = ogr.FieldDefn('GBCODECODE', fld_type)
fld_defn2 = ogr.FieldDefn('NEWFIELD4', ogr.OFTInteger)
'''使用ALTER_NAME_FLAG更改属性字段名称，必须保证字段的数据类型一致'''
lyr.AlterFieldDefn(i, fld_defn, ogr.ALTER_NAME_FLAG)
'''新建属性字段'''
lyr.CreateField(fld_defn2)
'''删除字段(参数为字段索引值)'''
lyr.DeleteField(lyr.FindFieldIndex('NEWFIELD', 0))
lyr.DeleteField(lyrdefn.GetFieldIndex('NEWFIELD2'))

        2、要素添加、更新和删除

'''增加、删除和更新要素'''
fn = r'D:\soft\geoprocessing-with-python\china_basic_map'
ds = ogr.Open(fn, 1)
if ds is None:
    sys.exit('Could not open {0}.'.format(fn))

#更新要素（新增属性字段）
lyr = ds.GetLayer(0)
lyr.CreateField(ogr.FieldDefn('NEWFIELD', ogr.OFTInteger))
n = 1
for fea in lyr:
    fea.SetField('NEWFIELD', n*n)
    lyr.SetFeature(fea)
    n += 1

#删除要素
for fea in lyr:
    print(fea.GetField('NEWFIELD'))
    if fea.GetField('NEWFIELD') == 3179089:
        lyr.DeleteFeature(fea.GetFID())

print(lyr.GetFeatureCount())

        对数据进行更新操作后，一般需要对数据库进行压缩、重组或重新计算数据的空间范围等，如下：

ds.ExecuteSQL('REPACK ' + lyr.GetName()) #压缩数据库
ds.ExecuteSQL('VACUUM') #重组功能

ds.ExecuteSQL('RECOMPUTE EXTENT ON ' + lyr.GetName()) #重新计算图层的空间范围
print(lyr.GetExtent())

转载自：https://blog.csdn.net/shaxiaozilove/article/details/79421513