3. 汇流累积阈值对河网密度的影响分析

为什么阈值是河网提取的“密度旋钮”

在基于DEM的水文分析中，汇流累积量（Flow Accumulation）表示每个栅格上游汇入的栅格数量。设定阈值 T 后，累积量 ≥ T 的栅格被判定为河道。阈值越小，被判为河道的栅格越多，河网越密；阈值越大，只有汇水面积足够大的沟道才保留，河网越稀疏。因此阈值直接控制三个可观测指标：河网密度（河道总长度 / 流域面积）、河源数量（河网线段起点数）以及流域出口位置（河网最下游节点与流域边界的吻合程度）。

本案例采用一套可复现的栅格计算流程：先由流量累积栅格按不同阈值二值化，再统计河道栅格数与河源数，最后用“河源密度—出口偏移”双指标确定合理阈值区间。下面给出可直接运行的核心计算代码（Python + NumPy，假设 flow_acc 为已由填洼、流向计算得到的流量累积二维数组，单位：上游栅格数）。

import numpy as np

def extract_channel(flow_acc, threshold):
    """按阈值提取河道二值栅格：1=河道，0=非河道"""
    return (flow_acc >= threshold).astype(np.uint8)

def count_sources(channel):
    """统计河源数：河道栅格中，上游无河道邻居的起点数量（8邻域近似）"""
    src = 0
    h, w = channel.shape
    for i in range(1, h - 1):
        for j in range(1, w - 1):
            if channel[i, j] == 1:
                neigh = channel[i-1:i+2, j-1:j+2].sum() - 1
                if neigh == 0:
                    src += 1
    return src

def evaluate_thresholds(flow_acc, cell_area, thresholds):
    rows = []
    for t in thresholds:
        ch = extract_channel(flow_acc, t)
        length = ch.sum() * np.sqrt(cell_area)   # 近似：每河道栅格折算一个格网边长
        density = length / (flow_acc.size * cell_area)
        rows.append((t, int(ch.sum()), count_sources(ch), round(density, 4)))
    return rows

# 示例：30 m 分辨率 DEM，阈值从 100 到 5000 栅格
cell_area = 30 * 30
thresholds = [100, 300, 500, 1000, 2000, 3000, 5000]
for r in evaluate_thresholds(flow_acc, cell_area, thresholds):
    print("T=%5d  河道栅格=%6d  河源数=%4d  河网密度=%.4f" % r)

示例输出（教学示例数据，非实测）：T=100 时河源数约 480、密度 0.021；T=1000 时河源数约 96、密度 0.008；T=5000 时河源数约 21、密度 0.003。可以看到河源数随阈值增大呈近似幂律衰减，而密度下降更平缓——这正是判断“合理阈值”的关键：当阈值继续增大而河源数骤减、主干河道形态开始断裂时，说明已越过合理上限。

如何确定研究区的合理阈值范围

单看密度曲线不足以定阈值，需要把“形态合理性”作为约束。推荐按以下判据筛选：

• 河源密度稳定段：绘制河源数—阈值曲线，取曲线由陡变缓的拐点作为下限，避免把坡面细沟误判为河道。
• 主干连通性：阈值增大到主干河道出现断裂（河道栅格不连续）时，该阈值即为上限。
• 出口吻合度：将提取河网最下游节点与已知流域出口或实测水系出口比对，偏移超过 1–2 个格网即视为过大。
• 与实测水系叠加：在合理区间内取 2–3 个候选阈值，分别与实测水系做叠加目视检查，选择形态最接近者。

结果判读与交付要点

完成多阈值对比后，应输出一张阈值—指标对照表，记录每个候选阈值的河道栅格数、河源数、河网密度与出口偏移量，并标注推荐阈值。推荐阈值需同时满足：河源数处于稳定段、主干无断裂、出口偏移在允许范围内。最终河网矢量应保留阈值参数与DEM分辨率作为元数据，便于复现与横向比较。