这一部分重点介绍如何使用 OGR 工具处理矢量数据，包括 ETL（提取-转换-加载）流程、合并矢量文件以及执行空间查询和地理处理。 来源说明 本文翻译自 Spatial Thoughts，原作者为 Ujaval Gandhi。

**原文链接：** https://courses.spatialthoughts.com/gdal-tools.html

**本文仅供学习交流，版权归原作者所有。**

3. OGR 矢量工具 (OGR Tools)

我们将学习如何使用 OGR 工具处理矢量数据。这套工具是 GDAL 软件包的一部分，遵循相同的约定。除了格式转换外，OGR 工具还支持运行空间 SQL 查询，使其在构建空间 ETL 管道方面非常强大。

观看模块视频 ↗

3.1 ETL 基础 (ETL Basics)

在本节中，我们将逐步了解如何构建提取-转换-加载 (ETL) 流程。示例工作流将向你展示如何：

读取 CSV 数据源
将其转换为点数据图层
为其分配 CRS（坐标参考系统）
提取子集
更改列的数据类型
将结果写入 GeoPackage

3.1.1 读取 CSV 数据源

你的数据包中有一个名为 worldcities.csv 的 CSV 文件。此文件包含有关世界主要城市的基本信息及其坐标。

让我们使用 ogrinfo 命令检查数据集。

ogrinfo worldcities.csv

程序可以成功打开并读取文件，但没有显示任何其他信息。我们可以使用 -al 选项实际读取文件中的所有行，并结合 -so 选项打印摘要。

ogrinfo -so -al worldcities.csv

3.1.2 将其转换为点数据图层

我们现在获得了一个摘要，其中包含文件中的要素总数以及列及其类型。这是一个普通的 CSV 文件。让我们使用 lng 和 lat 字段中提供的 X 和 Y 坐标将其转换为空间数据图层。-oo 选项允许我们指定特定于格式的选项。OGR CSV 驱动程序允许使用 X_POSSIBLE_NAMES 和 Y_POSSIBLE_NAMES 选项指定包含几何信息的列名列表或名称模式（例如 Lat*）。

ogrinfo -so -al worldcities.csv -oo X_POSSIBLE_NAMES=lng -oo Y_POSSIBLE_NAMES=lat

OGR 现在能够将此图层识别为点几何图层。让我们将其写入空间数据格式。我们可以使用 ogr2ogr 实用程序在不同格式之间进行转换。以下命令从 CSV 文件创建一个名为 worldcities.gpkg 的新 GeoPackage 文件。

Windows

ogr2ogr -f GPKG worldcities.gpkg worldcities.csv ^
  -oo X_POSSIBLE_NAMES=lng -oo Y_POSSIBLE_NAMES=lat

Mac/Linux

ogr2ogr -f GPKG worldcities.gpkg worldcities.csv \
  -oo X_POSSIBLE_NAMES=lng -oo Y_POSSIBLE_NAMES=lat

3.1.3 为其分配 CRS

我们可以在 GIS 软件中打开结果，它会显示城市图层。虽然点图层加载正常，但它缺少 CRS。我们可以使用 -a_srs 选项为生成的图层分配 CRS。

Windows

ogr2ogr -f GPKG worldcities.gpkg worldcities.csv ^
  -oo X_POSSIBLE_NAMES=lng -oo Y_POSSIBLE_NAMES=lat -a_srs EPSG:4326

Mac/Linux

ogr2ogr -f GPKG worldcities.gpkg worldcities.csv \
  -oo X_POSSIBLE_NAMES=lng -oo Y_POSSIBLE_NAMES=lat -a_srs EPSG:4326

3.1.4 提取子集

OGR 工具允许对数据源执行 SQL 查询。它支持 OGR SQL 语法中描述的 SQL 功能子集。选择要素子集的一种简单方法是使用 -where 选项。这允许你指定属性查询来过滤结果。在这里，我们修改命令以仅提取 country 列值为 India 的城市。

Windows

ogr2ogr -f GPKG mycities.gpkg worldcities.csv ^
  -oo X_POSSIBLE_NAMES=lng -oo Y_POSSIBLE_NAMES=lat -a_srs EPSG:4326 ^
  -where "country = 'India'"

Mac/Linux

ogr2ogr -f GPKG mycities.gpkg worldcities.csv \
  -oo X_POSSIBLE_NAMES=lng -oo Y_POSSIBLE_NAMES=lat -a_srs EPSG:4326 \
  -where "country = 'India'"

3.1.5 更改列的数据类型

你还可以对输入数据源执行任何 SQL 语句。这是一个非常强大的功能，允许你过滤、连接、转换和汇总输入数据。我们可以将其应用于我们的数据图层，将列类型从字符串转换为整数。输入 CSV 中的 population（人口）列类型为 字符串。由于此字段主要用于存储整数值，我们可以使用 CAST() SQL 函数将类型更改为整数。SELECT 语句还允许我们只选择写入输出的相关字段。

Windows

ogr2ogr -f GPKG mycities.gpkg worldcities.csv ^
  -oo X_POSSIBLE_NAMES=lng -oo Y_POSSIBLE_NAMES=lat -a_srs EPSG:4326 ^
  -sql "SELECT city, country, CAST(population AS integer) as population ^
    from worldcities where country = 'India'"

Mac/Linux

ogr2ogr -f GPKG mycities.gpkg worldcities.csv \
  -oo X_POSSIBLE_NAMES=lng -oo Y_POSSIBLE_NAMES=lat -a_srs EPSG:4326 \
  -sql "SELECT city, country, CAST(population AS integer) as population \
    from worldcities where country = 'India'"

3.1.6 重命名 GeoPackage 中的图层

最后，我们可以使用 -nln 选项重命名输出图层的名称。

Windows

ogr2ogr -f GPKG mycities.gpkg worldcities.csv ^
  -oo X_POSSIBLE_NAMES=lng -oo Y_POSSIBLE_NAMES=lat -a_srs EPSG:4326 ^
  -sql "SELECT city, country, CAST(population AS integer) as population ^
  from worldcities where country = 'India'" -nln mycities

Mac/Linux

ogr2ogr -f GPKG mycities.gpkg worldcities.csv \
  -oo X_POSSIBLE_NAMES=lng -oo Y_POSSIBLE_NAMES=lat -a_srs EPSG:4326 \
  -sql "SELECT city, country, CAST(population AS integer) as population \
  from worldcities where country = 'India'" -nln mycities

这展示了 OGR 命令行工具的强大功能。只需一行代码，我们就可以读取、过滤、转换并将结果写入 OGR 支持的任何输出格式。

练习 6

编写一个使用 ogr2ogr 执行以下操作的命令：

从 mycities.gpkg 文件中读取新创建的 mycities 图层。
将图层重投影到新的 CRS EPSG:7755 (WGS84 / India NSF LCC)。
将结果保存为 shapefile mycities.shp。

所有上述操作应在单个命令中完成。

3.2 合并矢量文件 (Merging Vector Files)

OGR 提供的另一个有用的实用程序是 ogrmerge.py。此命令可以将多个矢量图层合并为单个数据图层。它还有几个选项来控制数据集的合并方式——使其成为数据处理任务的便捷工具。

我们将处理包含地震位置的多个 GeoJSON 文件。你的数据包在 earthquakes/ 目录中包含多个 GEOJSON 文件。我们将把这 12 个文件合并为一个文件。切换到该目录：

cd earthquakes

我们将通过合并目录中的所有文件来创建一个 geopackage 文件。Windows 用户可能需要指定 ogrmerge.py 脚本的完整路径。

Windows

python %CONDA_PREFIX%\Scripts\ogrmerge.py -o earthquakes.gpkg *.geojson

Mac/Linux

ogrmerge.py -o earthquakes.gpkg *.geojson

结果将是一个包含 12 个图层的单个 geopackage——每个源文件对应一个图层。对于大多数应用程序，最好将源文件合并为单个图层。我们可以使用 -single 选项来指示我们希望输出为单个图层。我们还使用 -nln 选项指定合并图层的名称。由于 GeoPackage 数据集已经存在，我们需要指定 -overwrite_ds 用新内容覆盖文件。

Windows

python %CONDA_PREFIX%\Scripts\ogrmerge.py -o earthquakes.gpkg *.geojson ^
  -single -nln all_earthquakes -overwrite_ds

Mac/Linux

ogrmerge.py -o earthquakes.gpkg *.geojson \
  -single -nln all_earthquakes -overwrite_ds

现在，GeoPackage 中将有一个名为 all_earthquakes 的单个图层，其中包含该年度记录的所有地震。

提示： 有一个非常有用的选项 -src_layer_field_name，它可以向输出图层添加一个新字段，其中包含贡献该特定记录的输入文件的名称。

练习 7

创建一个新的 GeoPackage large_earthquakes.gpkg，仅包含震级大于 4.5 的地震。

使用 ogr2ogr 实用程序读取上一节中创建的 earthquakes.gpkg 文件。
使用 -where 选项编写查询以过滤 mag 字段值大于 >4.5 的记录。

3.3 地理处理和空间查询 (Geoprocessing and Spatial Queries)

到目前为止，我们已经使用了符合 OGR SQL 方言的查询。OGR 工具还支持 SQLite 方言。使用 SQLite 方言的主要优点是可以使用 Spatialite 提供的空间 SQL 函数。这使得你可以使用 OGR 工具进行广泛的空间查询应用。使用 SQLite 方言时的一个主要区别是你必须显式指定几何列。

让我们做一些地理处理，看看如何使用 OGR 工具进行空间查询。你的数据包有一个 geopackage 文件 spatial_query.gpkg。此 geopackage 包含 2 个图层：

metro_stations：墨尔本市地铁站位置的点图层。
bars_and_pubs：墨尔本市酒吧和酒馆位置的点图层。

让我们看看如何向此 geopackage 添加一个新图层，其中包含距离地铁站 500 米以内的所有酒吧和酒馆。

3.3.1 重投影矢量图层

源数据图层位于 CRS EPSG:4326 中。由于我们要运行以米为单位距离的空间查询，让我们将两个图层都重投影到投影 CRS。让我们使用 ogr2ogr 命令重投影两个输入图层并将它们保存到同一个 geopackage 中。我们可以使用 -t_srs 选项指定 EPSG:7855 (GDA2020 / MGA Zone 55) 作为目标投影。-update 选项与 -nln 一起告诉 ogr2ogr 创建一个新图层并将其保存在同一个 geopackage 中。

Windows

ogr2ogr -t_srs EPSG:7855 spatial_query.gpkg spatial_query.gpkg ^
  bars_and_pubs -update -nln bars_and_pubs_reprojected
ogr2ogr -t_srs EPSG:7855 spatial_query.gpkg spatial_query.gpkg ^
  metro_stations -update -nln metro_stations_reprojected

Mac/Linux

ogr2ogr -t_srs EPSG:7855 spatial_query.gpkg spatial_query.gpkg \
  bars_and_pubs -update -nln bars_and_pubs_reprojected
ogr2ogr -t_srs EPSG:7855 spatial_query.gpkg spatial_query.gpkg \
  metro_stations -update -nln metro_stations_reprojected

3.3.2 创建缓冲区

让我们使用 500 米的距离对 metro_stations_reprojected 图层进行缓冲。这里我们使用 -sql 选项指定 SQL 查询。注意使用 Spatialite 引擎提供的 ST_Buffer() 函数。我们需要指定 -dialect SQLITE，因为查询使用了空间函数。

空间数据库函数通常带有 ST_ 前缀。ST 代表 Spatial Type（空间类型）。

Windows

ogr2ogr spatial_query.gpkg spatial_query.gpkg -update -nln metro_stations_buffer ^
  -sql "SELECT m.station, ST_Buffer(m.geom, 500) as geom FROM ^
    metro_stations_reprojected m" -dialect SQLITE

Mac/Linux

ogr2ogr spatial_query.gpkg spatial_query.gpkg -update -nln metro_stations_buffer \
  -sql "SELECT m.station, ST_Buffer(m.geom, 500) as geom FROM \
    metro_stations_reprojected m" -dialect SQLITE

Metro Station Buffers

此查询会产生单独的重叠缓冲区。我们可以使用 ST_Collect() 函数融合缓冲区。

Windows

ogr2ogr spatial_query.gpkg spatial_query.gpkg -update ^
  -nln metro_stations_buffer_dissolved -sql "SELECT ST_Union(d.geom) as geom ^
    FROM (SELECT ST_Collect(buffer(m.geom, 500)) as geom ^
    FROM metro_stations_reprojected m) as d" -dialect SQLITE

Mac/Linux

ogr2ogr spatial_query.gpkg spatial_query.gpkg -update \
  -nln metro_stations_buffer_dissolved -sql "SELECT ST_Union(d.geom) as geom \
  FROM (SELECT ST_Collect(buffer(m.geom, 500)) as geom \
  FROM metro_stations_reprojected m) as d" -dialect SQLITE

Dissolved Buffers

3.3.3 执行空间查询

现在我们有了一个包含地铁站周围缓冲区区域的多边形图层，我们可以使用 ogr_layer_algebra.py 工具选择落在缓冲区区域内的地铁站。此工具可以执行各种类型的空间查询，包括相交、差异（擦除）、联合等。我们将使用 Clip（裁剪）模式从输入图层中选择被次要图层要素裁剪的要素。

Windows

python %CONDA_PREFIX%\Scripts\ogr_layer_algebra.py Clip ^
  -input_ds spatial_query.gpkg -input_lyr bars_and_pubs_reprojected ^
  -method_ds spatial_query.gpkg -method_lyr metro_stations_buffer_dissolved ^
  -output_ds output.gpkg -output_lyr selected -nlt POINT

Mac/Linux

ogr_layer_algebra.py Clip \
  -input_ds spatial_query.gpkg -input_lyr bars_and_pubs_reprojected \
  -method_ds spatial_query.gpkg -method_lyr metro_stations_buffer_dissolved \
  -output_ds output.gpkg -output_lyr selected -nlt POINT

3.3.4 数据清洗 (Data Cleaning)

上一个查询完成了工作，你找到了缓冲区区域内的所有点。但你会注意到同一位置有许多重复点。经检查，你会发现源数据包含每个机构不同年份的多个要素。我们可以运行最终的 SQL 查询，通过为每个机构选择最新年份的要素来对数据进行去重。请注意，某些列名包含空格，因此我们将它们括在双引号中。

Windows

ogr2ogr -t_srs EPSG:4326 output.gpkg output.gpkg -update ^
  -nln cleaned -sql "SELECT t1.* FROM selected t1 ^
  JOIN (SELECT \"Property ID\", MAX(\"Census year\") as year FROM selected ^
  GROUP BY \"Property ID\") t2 ON t1.\"Property ID\" = t2.\"Property ID\" ^
  AND t1.\"Census year\" = t2.year"

Mac/Linux

ogr2ogr -t_srs EPSG:4326 output.gpkg output.gpkg -update \
  -nln cleaned -sql "SELECT t1.* FROM selected t1 \
  JOIN (SELECT \"Property ID\", MAX(\"Census year\") as year FROM selected \
  GROUP BY \"Property ID\") t2 ON t1.\"Property ID\" = t2.\"Property ID\" \
  AND t1.\"Census year\" = t2.year"