2026/6/20 6:57:20
网站建设
项目流程
中山营销型网站设计,新公司做网站多少钱,绍兴网站建设方案服务,永久免费手机建站平台前言
❝
在GIS开发中#xff0c;空间查询和属性查询都是常见的基础操作#xff0c;是每一个GISer都要掌握的必备技能。实现高效的数据查询功能可以提升用户体验#xff0c;提升数据可视化效率。在之前的文章中讲了如何使用GDAL或者ogr2ogr工具将txt以及csv文本数据转换为Shp…前言❝在GIS开发中空间查询和属性查询都是常见的基础操作是每一个GISer都要掌握的必备技能。实现高效的数据查询功能可以提升用户体验提升数据可视化效率。在之前的文章中讲了如何使用GDAL或者ogr2ogr工具将txt以及csv文本数据转换为Shp格式本篇教程在之前一系列文章的基础上讲解如何使用GDAL实现数据空间查询功能。GDAL 简介GDAL 下载安装GDAL 开发起步GDAL 实现 GIS 数据读取转换全如果你还没有看过建议从以上内容开始。1. 开发环境本文使用如下开发环境以供参考。时间2025年系统Windows 11Python3.11.7GDAL3.11.12. 空间查询在GDAL中有两个图层方法可以用于实现空间查询分别是SetSpatialFilter和SetSpatialFilterRect。此种查询方式为直接在源数据上操作返回结果为查询图层。参数-geom用于几何查询的几何对象SetSpatialFilter(geom)参数-minx最小x坐标-miny最小y坐标-maxx最大x坐标-maxy最大y坐标SetSpatialFilterRect(minx,miny,maxx,maxy)传入None值时重置查询。sourceLayer.SetSpatialFilter(None)在进行正式查询前和以前的文章一样首先获取数据驱动添加数据源获取图层数据。定义一个方法SpatialFilter用于实现空间查询该方法接受两个参数一个sourcePath传递源数据图层路径另一个selectPath用于定义查询图层路径。说明图层属性过滤参数-sourcePath待查询 Shp 文件路径-selectPath查询 Shp 文件路径def SpatialFilter(sourcePath,selectPath):在以下代码中完成图层数据的读取操作。# 注册所有驱动ogr.RegisterAll()# 添加数据驱动shpDriver ogr.GetDriverByName(ESRI Shapefile)checkFilePath(sourcePath,shpDriver)checkFilePath(selectPath,shpDriver)# 打开数据源sourceDs shpDriver.Open(sourcePath)selectDs shpDriver.Open(selectPath)hasDs sourceDs and selectDsifhasDs is None:print(数据源打开异常请检查路径)returnFalse# 获取图层sourceLayer sourceDs.GetLayer(0)selectLayer selectDs.GetLayer(0)下面对几何查询以及矩形查询两种实现方式进行介绍。2.1. 几何查询通过图层方法GetNextFeature获取第一个要素对象。# 获取查询几何对象queryFeat selectLayer.GetNextFeature()queryGeom queryFeat.GetGeometryRef()print(f查询要素id{queryFeat.GetField(Id)})print(f查询几何对象{queryGeom})将几何对象传入SetSpatialFilter方法进行空间过滤。# 空间过滤sourceLayer.SetSpatialFilter(queryGeom)查询完成之后传入None值结束查询。# 结束查询sourceLayer.SetSpatialFilter(None)以下为几何查询部分代码。# 获取查询几何对象queryFeat selectLayer.GetNextFeature()queryGeom queryFeat.GetGeometryRef()print(f查询要素id{queryFeat.GetField(Id)})print(f查询几何对象{queryGeom})# 获取要素数量featureCount sourceLayer.GetFeatureCount()print(f所有要素数量{featureCount})# 空间过滤sourceLayer.SetSpatialFilter(queryGeom)queryFeatCount sourceLayer.GetFeatureCount()print(f查询要素数量{queryFeatCount})# 结束查询sourceLayer.SetSpatialFilter(None)finalFeatCount sourceLayer.GetFeatureCount()print(f重置查询后要素数量{finalFeatCount})print(n~~~~~~~方式一结束几何查询~~~~~~~)若有兴趣还可以创建自定义几何对象进行空间查询。# 自定义Geometry查询对象customGeom ogr.Geometry(ogr.wkbPolygon)ringGeom ogr.Geometry(ogr.wkbLinearRing)ringGeom.AddPoint(102.884350,32.501570)ringGeom.AddPoint(105.025865,74.974949)ringGeom.AddPoint(50.417235,55.701314)ringGeom.AddPoint(50.417235,32.501570)ringGeom.CloseRings()customGeom.AddGeometry(ringGeom)# 获取要素数量featureCount sourceLayer.GetFeatureCount()print(f所有要素数量{featureCount})# 空间过滤sourceLayer.SetSpatialFilter(customGeom)queryFeatCount sourceLayer.GetFeatureCount()print(f查询要素数量{queryFeatCount})# 结束查询sourceLayer.SetSpatialFilter(None)finalFeatCount sourceLayer.GetFeatureCount()print(f重置查询后要素数量{finalFeatCount})print(n~~~~~~~方式一结束几何查询~~~~~~~)如下为几何查询输出结果如下为该要素在ArcGIS中查询显示结果。2.2. 矩形查询矩形查询主要调用SetSpatialFilterRect方法传入x、y坐标即可。# 获取要素数量featureCount2 sourceLayer.GetFeatureCount()print(f所有要素数量{featureCount2})# 空间过滤sourceLayer.SetSpatialFilterRect(-16.9994,23.1235,42.9111,49.5852,)queryFeatCount2 sourceLayer.GetFeatureCount()如下为矩形查询输出结果3. 注意事项在windows开发环境中同时安装GDAL与PostGIS其中投影库PROJ的环境变量指向PostGIS的安装路径在运行GDAL程序时涉及到要素、几何与投影操作时会导致异常。具体意思为GDAL不支持PostGIS插件中的投影库版本需要更换投影库或者升级版本。RuntimeError: PROJ: proj_identify: D:Program FilesPostgreSQL13sharecontribpostgis-3.5projproj.db contains DATABASE.LAYOUT.VERSION.MINOR 2 whereas a number 5 is expected. It comes from another PROJ installation.解决办法为修改PROJ的环境变量到GDAL支持的版本或者在GDAL程序开头添加以下代码os.environ[PROJ_LIB] rD:\Programs\Python\Python311\Libsite-packages\osgeo\data\proj