20_计算栅格面积 QGIS3

原文链接: https://www.qgistutorials.com/en/docs/3/calculating_raster_area.html

计算栅格面积 (QGIS3)¶

许多应用需要对区域内的土地利用模式进行量化。土地利用土地覆盖 (LULC) 数据集以栅格文件的形式提供，其中每个像素都被分配了一个类别值。GIS 分析师通常需要基于此数据生成报告，计算给定区域内每个类别的面积。QGIS 自带了许多内置工具来计算和汇总栅格面积。

注意

历史上，计算栅格面积的建议方法是将栅格转换为矢量图层并使用矢量面积计算技术。这种方法计算量相当大且容易出错。推荐的方法是使用处理工具 栅格图层唯一值报告，它可以直接计算像素面积。如果你有一个包含许多多边形的图层，并且需要每个多边形的面积，你可以使用 分区直方图 工具获取每个类别的像素计数，然后乘以每个像素的面积。

任务概述¶

我们将使用一个包含 11 个土地覆盖类别的栅格图层，计算一个国家公园内每个类别的面积。我们还将对结果进行后处理，创建一个包含类别名称和面积的电子表格。

您将学到的其他技能¶

• 如何从 .qml 格式的样式文件为图层应用符号系统。
• 如何使用 CASE 语句编写包含多个 if/else 条件的表达式。
• 如何将 QGIS 表格导出为 Excel 电子表格。

获取数据¶

我们将下载以下数据集：

1. 世界保护区数据库 (WDPA)：我们将下载印度卡齐兰加国家公园边界的 Shapefile。
2. ESA WorldCover 2020：欧洲航天局 (ESA) WorldCover 10 米 2020 产品提供了 2020 年全球 10 米分辨率的土地覆盖图。我们将下载覆盖我们感兴趣区域的图块。

公园边界¶

1. 访问 Protected Planet 网站，点击搜索工具箱。搜索 Kaziranga National Park。

1. 保护区矢量边界将作为搜索结果显示。点击它以查看 Kaziranga National Park 的页面。

1. 此页面将包含总面积、创建年份等附加信息。点击下载并点击 SHP，以Shapefile 格式下载数据。

1. 将出现两个下载选项。在非商业用途下点击继续。现在将下载一个包含国家公园边界的 zip 文件。

土地覆盖数据¶

1. 访问 ESA WorldCover 网站，点击 DATA ACCESS 菜单。

1. 向下滚动到 DATA DOWNLOAD 部分，点击链接打开 WorldCover 查看器。

1. 您需要创建一个免费账户才能下载数据。点击右上角的注册链接。按照说明创建新账户。

1. 创建账户后，使用凭据登录。本教程的感兴趣区域是卡齐兰加国家公园。缩放到印度东北地区。一旦您放大到足够程度，土地覆盖图像图块的边界框将开始出现。

1. 搜索并定位 N24E093 图块区域。

1. 选择该图块并点击 NEXT。

1. 点击 DOWNLOAD 下载一个包含土地覆盖信息（栅格格式）的 zip 文件。请确保选择下载 WorldCover 版本 1 的数据。

1. 我们还将下载 ESA 提供的符号系统文件。访问 ESA WorldCover 数据访问 页面，向下滚动到符号系统部分。点击 QGIS 下载 ESAWorldCover_ColorLegend.qml 文件，该文件可用于使用适当的颜色和类别标签来样式化栅格图层。

为方便起见，您可以直接通过以下链接下载数据集的副本：

• 卡齐兰加国家公园边界 WDPA_WDOECM_Apr2022_Public_10744_shp_0.zip
• ESA WorldCover 图块 N24E093 子集 terrascope_download_20220422_114733.zip
• ESA WorldCover QML 样式文件 ESAWorldCover_ColorLegend.qml

数据来源 [WDPA] [WorldCover]

操作步骤¶

1. 解压缩所有下载的文件。在浏览器中，找到包含边界文件 WDPA_WDOECM_Apr2022_Publicc_10744_shp-polygons.shp 的文件夹，并将其拖放到 QGIS 画布中。

1. 现在找到世界覆盖栅格图块 ESA_WorldCover_10m_2020_v100_N24_E093_Map.tif 并将其拖放到 QGIS 画布中。

1. 现在图层面板中已加载矢量边界和土地覆盖栅格图层。让我们将土地覆盖栅格裁剪到国家公园边界。转到 处理 ‣ 工具箱 以打开处理工具箱。搜索并找到 GDAL ‣ 栅格提取 ‣ 按掩膜图层裁剪栅格 算法。双击启动它。

1. 在"按掩膜图层裁剪栅格"对话框中，选择 ESA_WorldCover_10m_2020_v100_N24_E093_Map 图层作为输入图层，选择 WDPA_WDOECM_Apr2022_Publicc_10744_shp-polygons 图层作为掩膜图层。在"为输出波段指定无数据值"部分输入 -9999。

1. 现在打开"高级参数"部分，在"配置文件"中选择 High Compression。然后在"已裁剪(掩膜)"下，点击 ... 并选择 保存到文件…。输入文件名 worldcover_clipped.tif。点击 运行。

1. 现在 worldcover_clipped 图层将加载到 QGIS 画布中。右键单击 ESA_WorldCover_10m_2020_v100_N24_E093_Map 图层，选择 移除图层…

1. 我们的两个图层都使用地理坐标参考系统 EPSG:4326。此 CRS 的单位是度，不适合计算面积。我们必须首先将图层重新投影到投影坐标参考系统。对于此类区域分析，UTM 是投影 CRS 的良好选择。我们将把图层重新投影到本地 UTM 分带的 CRS。打开处理工具箱，搜索 矢量通用 ‣ 重投影图层 算法。双击启动它。

1. 在"重投影图层"对话框中，选择 WDPA_WDOECM_Apr2022_Publicc_10744_shp-polygons 图层作为输入图层，点击"目标 CRS"下的选择 CRS 按钮。

1. 我们的感兴趣区域位于 UTM 46N 分带。搜索 46 N 并选择 WGS 84 / UTM zone 46N CRS。

注意

要查找您所在区域的 UTM 分带，请参考 What UTM Zone am I in 网站。

1. 在"重投影后的"部分，点击 ... 并选择 保存到文件…。输入名称 nationalpark_reprojected.gpkg。点击 运行。

1. 现在 nationalpark_reprojected 图层将加载到画布中。右键单击 WDPA_WDOECM_Apr2022_Publicc_10744_shp-polygons 图层，选择 移除图层… 将其移除。现在我们将重新投影栅格图层。在处理工具箱中，搜索并找到 GDAL ‣ 栅格投影 ‣ 扭曲（重投影）。

1. 在"扭曲（重投影）"对话框中，选择 worldcover_clipped 作为输入图层，在"目标 CRS"中选择 WGS 84 / UTM zone 46N CRS。打开"高级参数"并在"配置文件"中选择 High Compression。

1. 现在在"重投影后的"下，点击 ... 并选择 保存到文件…。输入名称 worldcover_reprojected.tif。点击 运行。

1. 现在 worldcover_reprojected 图层将加载到画布中，移除 worldcover_clipped 图层。让我们将项目 CRS 设置为 UTM 分带。点击任何图层，选择 图层 CRS ‣ 从图层设置项目 CRS。

1. 现在项目 CRS 将更新。让我们根据 ESA WorldCover 数据集的类别名称和颜色设置栅格图层的符号系统。右键单击 worldcover_reprojected 图层，点击 属性…

1. 在"图层属性"对话框中，选择 符号系统。您会注意到图层颜色以黑白调显示。要解决此问题，点击 样式 ‣ 加载样式…。浏览并选择 ESAWorldCover_ColorLegend.qml 文件。

1. 现在您可以看到更新后的符号颜色和类别描述。点击 确定。

1. 在图层面板中展开 worldcover_reprojected 图层，查看带有正确类别描述的图例。

1. 现在让我们计算每个类别的面积。在处理工具箱中，搜索并找到 栅格图层唯一值报告 工具。双击打开它。

1. 在"栅格图层唯一值报告"对话框中，选择 worldcover_reprojected 作为输入图层。在"唯一值表"下点击 ... 并选择 保存到文件…。输入名称 class_areas.gpkg。点击 运行。

1. 现在 class_areas 图层将添加到图层面板。右键单击该图层，点击 打开属性表。m2 列包含每个类别的面积（单位：平方米）。

1. 让我们将面积转换为平方公里。在处理工具箱中，搜索并选择 矢量表 ‣ 字段计算器。

1. 在"字段计算器"对话框中，选择 class_areas 图层作为输入图层。输入字段名称为 area_sqkm。在结果字段类型中选择 浮点型 (real)。在表达式窗口中，输入以下表达式。这将把平方米转换为平方公里，并将结果四舍五入到 2 位小数。在"已计算"下点击 ... 并选择 保存到文件…。输入名称 class_area_sqkm.gpkg。点击 运行。

bash round("m2"/ 1e6, 2)

1. 现在 class_area_sqkm 图层将加载到画布中。打开属性表，查看新添加的 area_sqkm 列。您会注意到值列包含每个类别的编号。为了使结果更容易解释，我们还要为每个类别编号添加描述。类别描述可在 ESA 产品用户手册 中找到。

1. 打开字段计算器，在输入图层中选择 class_area_sqkm 图层。输入字段名称为 landcover，在结果字段类型中选择 字符串 (string)。在表达式窗口中输入以下表达式。此表达式使用 CASE 语句根据多个条件分配值。在"已计算"下点击 ... 并选择 保存到文件…。输入名称 class_area_with_landcover.gpkg。点击 运行。

bash CASE WHEN "value" = 10 THEN 'Tree cover' WHEN "value" = 20 THEN 'Shrubland' WHEN "value" = 30 THEN 'Grassland' WHEN "value" = 40 THEN 'Cropland' WHEN "value" = 50 THEN 'Built-up' WHEN "value" = 60 THEN 'Bare / sparse vegetation' WHEN "value" = 70 THEN 'Snow and Ice' WHEN "value" = 80 THEN 'Permanent water bodies' WHEN "value" = 90 THEN 'Herbaceous wetland' WHEN "value" = 95 THEN 'Mangroves' WHEN "value" = 100 THEN 'Moss and lichen' END

1. 现在 class_area_with_landcover 图层将加载到画布中。打开属性表。landcover 列将包含与每个土地覆盖值对应的土地覆盖名称。

1. 让我们将此结果导出为 Excel 文件。在导出之前，我们还将整理表格并删除不需要的字段。在处理工具箱中，搜索并选择 矢量表 ‣ 重构字段。

1. 在"重构字段"对话框中，选择 class_area_with_landcover 图层作为输入图层。选中除 area_sqkm 和 landcover 之外的所有列，然后点击删除选定字段。

1. 您还可以使用"上移选定字段"/"下移选定字段"按钮更改表中字段的顺序。完成编辑后，点击"重构后的"旁边的 ... 按钮，选择 保存到文件…。选择 XLSX Files (*.xlsx) 作为格式。输入文件名 park_area_by_landcover.xlsx 并点击 保存。在"重构字段"对话框中，点击 运行 以应用更改。

1. 结果将是一个包含 landcover 和 area_sqkm 列的电子表格。

---

如果您想提供反馈或分享您对本教程的体验，请在下方留言。（需要 GitHub 账户）