在数字高程模型(DEM)处理中,执行填洼(Fill Sinks)是一项至关重要的预处理步骤,尤其在进行水文分析(如流域划分、水流路径提取、汇流累积量计算等)时不可或缺。
在DEM中,洼地是指被更高地形包围的局部低点,即没有出水口的负地形。
- 真实地貌:如喀斯特漏斗、湖泊盆地。
- 更常见的是:由遥感数据采集与插值过程中的误差导致。
大多数水文模拟算法(如D8等)假设水流应能从每个像元流向下游,最终汇入河流或海洋,不允许存在“死胡同”。
环境准备
准备PyQGIS运行环境,安装GRASS依赖。
修改变量
dem_part:修改为你的DEM图层名称
D:/output:修改为你的输出目录
代码
# -*- coding: utf-8 -*-
from qgis.core import QgsProject, QgsRasterLayer
import processing
import os
# 获取当前项目
project = QgsProject.instance()
layers = project.mapLayersByName('dem_part')
if not layers:
raise Exception("❌ 未找到名为 'dem_part' 的图层!")
input_dem = layers[0] # 获取 dem_part 图层
# 设置输出目录为 D:/output,并确保该路径存在
output_dir = r'D:/output'
if not os.path.exists(output_dir):
os.makedirs(output_dir) # 自动创建目录
# 构建输出文件路径
output_filled_dem = os.path.join(output_dir, 'filled_dem_grass.tif') # 填洼后DEM
output_flow_dir = os.path.join(output_dir, 'flow_direction_grass.tif') # 流向图
output_problem_areas = os.path.join(output_dir, 'problem_areas.tif') # 问题区域(可选)
# 使用 GRASS r.fill.dir 工具进行填洼
result = processing.run("grass7:r.fill.dir", {
'input': input_dem,
'output': output_filled_dem,
'direction': output_flow_dir,
'areas': output_problem_areas,
'format': 0, # 输出流向格式:0=grass(角度,逆时针从东开始)
'-f': False, # 不启用“仅查找”模式,正常执行填洼
'GRASS_REGION_PARAMETER': None, # 使用当前图层范围
'GRASS_REGION_CELLSIZE_PARAMETER': 0, # 使用原始分辨率
'GRASS_RASTER_FORMAT_OPT': '',
'GRASS_RASTER_FORMAT_META': ''
})
# 加载结果到地图
filled_raster = QgsRasterLayer(output_filled_dem, "填洼后DEM (r.fill.dir)")
flow_dir_raster = QgsRasterLayer(output_flow_dir, "流向图")
problem_raster = QgsRasterLayer(output_problem_areas, "问题区域")
if filled_raster.isValid():
project.addMapLayer(filled_raster)
else:
print("⚠️ 无法加载填洼后DEM图层")
if flow_dir_raster.isValid():
project.addMapLayer(flow_dir_raster)
else:
print("⚠️ 无法加载流向图图层")
if problem_raster.isValid():
stats = problem_raster.dataProvider().bandStatistics(1)
if stats.sum != 0: # 存在问题区域数据
project.addMapLayer(problem_raster)
else:
print("🗑️ 问题区域图为空,未添加到地图")
else:
print("⚠️ 无法加载问题区域图层")
print("✅ r.fill.dir 填洼完成!")
print(f"📁 输出目录: {output_dir}")
print(f"📊 填洼后DEM: {output_filled_dem}")
print(f"🧭 流向图: {output_flow_dir}")
print(f"⚠️ 问题区域: {output_problem_areas}")
更多算法问题,欢迎留言或联系我们。转载须注明出处。