GDAL python教程基础篇(5)OGR空间滤波器
创始人
2025-05-28 19:59:07
0

1.空间过滤器

如果说按照属性筛选要素是带有数据库特征的话,那么,根据空间位置的筛选就是纯GIS了。在OGR中,使用了Spatial filters(空间过滤)这一术语表征这一功能。

OGR提供的空间过滤功能有两种,一种是SetSpatialFilter(geom)—过滤某一类型的Feature,如参数中的Polygon,效用就是选出Layer中的所有Polygon覆盖的要素(注意,只要相交即可,不必完全包含)。另外还有 SetSpatialFilterRect(, , , ) 参数,可输入四个坐标选中矩形内的要素。

1.1SetSpatialFilter

下面我们来做一个示例,首先这里有两份数据分别为全国的市级行政区划矢量文件和江苏省的行政区划。然后我们通过空间滤波器选择江苏省的下辖市,并将其写入单独的文件当中。
在这里插入图片描述
示例代码

from osgeo import ogrprovice_path = r'D:/work/全国行政区划/江苏省.shp'
city_path = r'D:/work/全国行政区划/84坐标行政区划/市84坐标.shp'driver = ogr.GetDriverByName('ESRI Shapefile')    # 驱动provice_ds = driver.Open(provice_path)    # 获取数据集
city_ds = driver.Open(city_path)provice_layer = provice_ds.GetLayer(0)    # 获取图层
city_layer = city_ds.GetLayer(0)print(city_layer.GetFeatureCount())   # 打印出全国一共有371个地级市provice_feat = provice_layer.GetNextFeature()    # 获取要素poly = provice_feat.GetGeometryRef()       # 获取几何
city_layer.SetSpatialFilter(poly)           # 空间滤波器获取包含在几何内的cityimport osout_shp = 'D:/work/全国行政区划/江苏省_市级行政区划.shp'       # 定义输出shp
if os.access(out_shp, os.F_OK):           # 判断是否已经存在该文件driver.DeleteDataSource(out_shp)          # 若存在该文件则删除
newds = driver. CreateDataSource (out_shp)      # 创建数据集
pt_layer = newds.CopyLayer (city_layer, '')        # 通过图层拷贝方法获取图层
newds.Destroy ()                                              # 销毁数据源,若不销毁无法写入磁盘

将生成的文件加载入ArcGIS Pro,其中黄色的图层即我们通过空间滤波器得到的江苏省的市级行政区划。这里由于
在这里插入图片描述

1.2 SetSpatialFilterRect

SetSpatialFilterRect(, , , )` 参数,可输入四个坐标选中矩形内的要素。

print(city_layer.GetFeatureCount())
>>>27
city_layer.SetSpatialFilter(poly)    #SetSpatialFilter(None)则是清空空间属性的过滤器,可查看输出图层要素的数目。
print(city_layer.GetFeatureCount())
>>>371
city_layer.SetSpatialFilterRect(118.055315, 29.877947, 121.431257, 31.4785669)
print(city_layer.GetFeatureCount())
>>>16

2.根据属性条件选择要素

根据属性选择数据库记录是数据库应用中必不可少的一项功能。地理数据也是如此,譬如筛选人口在百万以上的城市,面积在百万平方千米以上的国家等等。

Layer对象中有 SetAttributeFilter()方法,可根据属性将Layer中符合某一条件的Feature过滤出来。设定了Filter之后用GetNextFeature()方法依次筛选出符合条件的Feature:

>>> from osgeo import ogr
>>> import os
>>> shpfile =  r'D:/work/全国行政区划/84坐标行政区划/市84坐标.shp'
>>> ds = ogr.Open(shpfile)
>>> layer = ds.GetLayer(0)
>>> lyr_count = layer.GetFeatureCount()
>>> print(lyr_count)
371
>>> layer.SetAttributeFilter("FID = 0")
>>> lyr_count = layer.GetFeatureCount()
>>> print(lyr_count)
1

3.在OGR中使用SQL语句进行查询

属性与空间的筛选可以看作是OGR的内置功能,这两种功能可以解决大部分实际应用问题。但是也有查询条件复杂的情况。针对这种情况,OGR也提供了灵活的解决方案:支持SQL查询语句。

例如执行SQL查询语句ExecuteSQL(),凭借SQL的强大功能,可执行更复杂的任务。例如下面这段代码,是从东北地区的分县数据中选择出吉林省的县级行政单位(对应的Prov_ID为22),并且按行政代码()降序输出。

>>> from osgeo import ogr
>>> import os
>>> driver = ogr.GetDriverByName("ESRI Shapefile")
>>> world_shp = '/gdata/GSHHS_h.shp'
>>> world_ds = ogr.Open(world_shp)
>>> world_layer = world_ds.GetLayer()
>>> world_layer_name = world_layer.GetName()
>>> result = world_ds.ExecuteSQL("select * from %s " % (world_layer_name)) # )
>>> resultFeat = result.GetNextFeature ()
>>> out_shp = '/tmp/sql_res.shp'
>>> create_shp_by_layer(out_shp, result)
>>> world_ds.ReleaseResultSet(result)

上面使用的SQL语句与平常的SQL语句无太大区别,相同点是使用了SELECT语句与WHERE条件,不同点是在OGR中此语句会生成空间数据。

最后一句ReleaseResultSet()是将查询结果释放,在执行下一条SQL语句之前一定要先释放。为简便而言,可同样将查询的结果当成是数据来查看。

相关内容

热门资讯

监控摄像头接入GB28181平... 流程简介将监控摄像头的视频在网站和APP中直播,要解决的几个问题是:1&...
Windows10添加群晖磁盘... 在使用群晖NAS时,我们需要通过本地映射的方式把NAS映射成本地的一块磁盘使用。 通过...
protocol buffer... 目录 目录 什么是protocol buffer 1.protobuf 1.1安装  1.2使用...
Fluent中创建监测点 1 概述某些仿真问题,需要创建监测点,用于获取空间定点的数据࿰...
educoder数据结构与算法...                                                   ...
MySQL下载和安装(Wind... 前言:刚换了一台电脑,里面所有东西都需要重新配置,习惯了所...
MFC文件操作  MFC提供了一个文件操作的基类CFile,这个类提供了一个没有缓存的二进制格式的磁盘...
在Word、WPS中插入AxM... 引言 我最近需要写一些文章,在排版时发现AxMath插入的公式竟然会导致行间距异常&#...
有效的括号 一、题目 给定一个只包括 '(',')','{','}'...
【Ctfer训练计划】——(三... 作者名:Demo不是emo  主页面链接:主页传送门 创作初心ÿ...