简介：中国三级流域数据集是一套面向地理信息系统应用与水资源研究的矢量地理数据，适合地理信息学习者、水文研究人员及相关从业者使用，可用于全国尺度水文分析与流域管理教学。资源内含全国三级流域的空间边界及其属性信息，每个流域拥有唯一编码，可用于洪灾风险评估、水资源分配、生态保护等场景；结合“重力与水文”的研究方向，还可辅助分析地下水补给区与排泄区，支撑可持续水资源管理。压缩包共8个文件，约4.09MB，主要包含 Shapefile 格式的空间数据（.shp）、属性数据表（.dbf）、坐标投影文件（.prj）、空间索引文件（.shx/.sbn/.sbx）以及元数据 XML，格式规范、配套齐全，可直接加载到 ArcGIS 或 QGIS 等平台开展流域查询、统计、制图与专题分析。目前已有4275人学习下载，对需要处理全国尺度水文数据或开展流域划分教学的用户具有实用价值。
1. 中国三级流域数据集：水文分析绕不开的那张底图资源
中国三级流域分区数据，是做水文分析时绕不开的那张底图。有一回我拿省界裁切径流数据核算生态流量，结果和水利公报差了将近四个百分点，查了三天才发现问题不在算法，在分区方式——省界和流域分水岭并不重合，山脊线不会顺着行政边界走。换上三级流域矢量数据后，同一套算法当场就对上了。这份资源就是一套全国三级流域分区边界数据，包含面状图斑、层级编码和面积属性，能直接用于水文统计、生态流量核算、水土流失评价，也能配合重力异常辅助判断边界合理性。适合做水资源、GIS 数据处理和环境评价的从业者，也适合想搞懂流域编码与空间边界关系的初学者。
2. 三级流域的分区逻辑：编码结构、层级关系与数据字段
拿到数据先别急着往地图上拖。三级流域分区不是随便画的多边形，它的编码体系里带了一整套水文拓扑关系，读懂了编码，后面做汇总、裁剪和关联数据都会省很多时间。这一章先把分区逻辑和属性字段讲清楚，再给一个能直接跑的加载检查流程。
2.1 一级到三级：流域划分的依据与边界特征
全国流域分区体系把国土按水循环的自然边界切成三层。一级流域对应大江大河干流及其入海水系，包括长江流域、黄河流域、珠江流域、松花江流域等，一共十个左右；二级流域在一级流域内按主要支流和大型湖泊控制区继续切分，比如长江流域里的若干主要支流分别对应各自的二级流域；三级流域则是二级流域内部进一步细分的汇水单元，通常以更小支流的集水区或水文测站控制区间为边界。
贯穿整个数据使用过程的一个原则是：边界是分水岭，不是行政区。山脊线决定降水往哪边流，不决定归哪个省管。所以三级流域边界和省界、县界会相交，同一个县域跨两三个三级流域是常态。做水资源量核算时，这种“跨界”恰恰保证了数据正确性——把每个站的径流按测站控制面积归并到对应三级流域，再逐级向上汇总，水量就能闭合。
三级流域的面积差异非常大。东部湿润地区一个三级流域可能只有几百平方千米，西部干旱区有的达到几万平方千米，这和河流发育密度与测绘条件都有关系。使用时要意识到，这套分区体系的重点是“控制单元”而不是均匀网格，后续统计应该按单元面积加权，而不是简单按图斑数量取平均。如果代码里出现直接对图斑数求平均的操作，通常都是方向搞反了。
另外值得留意的是流域分区与自然地理单元的对应关系。比如盆地地形被周边山地包围，盆地内部的支流水系构成一个相对独立的二级或三级单元；而平原地区地势平缓，分水岭不显著，边界往往借助人工渠道或测站控制区间划定。这两种情况在属性表里都能通过各级名称看出来，读到“控制区间”这类字眼时，说明这个单元的边界带有一定人为约定成分，不能完全按地形分水岭去理解。
2.2 属性表里的字段与编码规则
把 shp 文件加载进 GIS 后，属性表里最常见的字段是下面这几个：
| 字段名 | 示例值 | 说明 |
|---|---|---|
| code_level1 | F | 一级流域编码 |
| name_level1 | 长江流域 | 一级流域名称 |
| code_level2 | F02 | 二级流域编码 |
| name_level2 | 某支流流域 | 二级流域名称 |
| code_level3 | F0203 | 三级流域编码 |
| name_level3 | 某支流汇水区 | 三级流域名称 |
| area_km2 | 3564.78 | 椭球面积（平方千米） |
| length_km | 128.33 | 边界周长（千米） |
编码规则值得优先琢磨。F0203 不是随机编号：第一位 F 对应一级流域，前三位 F02 对应二级流域，最后两位是三级序号。也就是说，三级编码天然包含了它的一级、二级归属，用字符串前缀就能向上汇总。比如统计整个一级流域的某个指标，用 code_level1 == "F" 过滤；要统计某一支流控制区，用 code_level2.str.startswith("F02") 前缀判断就能实现，不需要再跑空间关联丢掉性能。这套前缀递进的设计在不同版本数据里可能有差异，有的版本用纯数字，有的用字母加数字，但“短编码是长编码的前缀”这个规律几乎不变。
还有一种常见现象：部分三级流域名称后面带括号，比如“某支流（含某某区间）”，表示这个单元合并了多个小型汇水区。遇到这类名称，查询时用字符串包含而不是精确匹配更稳妥。内陆河流域和部分闭合流域单独占用一级编码，比如西北诸河、西南诸河这类区域，降水少、蒸发强，三级流域划分的尺度明显偏大，面积字段比东部大一个数量级。对比东西部统计结果时，要格外留意分母不统一的问题，看到面积平均值差几十倍先别当数据错误，先看是不是跨越了一级流域边界。
2.3 打开文件先做三件事：投影检查、字段预览与图斑校验
我第一次用这类数据时直接拖进 GIS 就开始画图，结果图层和底图偏移了几百千米，查了半天发现 shp 里缺少投影定义文件，软件按默认坐标系去猜了。拿到数据先做三件套检查，能省掉后面很多玄学问题。下面这段代码用 GeoPandas 完成加载和基础体检：
import geopandas as gpd # 加载三级流域shp，显式指定CRS，避免默认猜测 gdf = gpd.read_file( "China_Third_Level_Basin.shp", encoding="utf-8", crs="EPSG:4326" ) # 第一件事：确认坐标系 print("CRS:", gdf.crs) # 第二件事：预览字段类型，检查空值 print(gdf.dtypes) print("空值统计：\n", gdf.isnull().sum()) # 第三件事：检查图斑数量与属性完整性 print("三级流域图斑数:", len(gdf)) print("一级流域列表:", gdf["name_level1"].unique()) # 面积合计与全国陆地面积做粗对照 total_area = gdf["area_km2"].sum() print("面积合计(万km²):", round(total_area / 10000, 2))这段代码做三件事：第一，打印 CRS 对象，确认数据坐标系，这类数据常见的是 WGS84 或 CGCS2000；第二，检查字段类型和空值，空值集中在面积字段时，说明原始数据清洗没做完整；第三，输出图斑总数和一级流域分布，面积合计与全国陆地面积对照能判断数据是否覆盖完整。如果数据自带投影定义，read_file 里的 crs 参数不会覆盖原有定义，但显式写出来能防止不带投影文件的 shp 被错误读取。
参数说明：EPSG:4326 是 WGS84 经纬度坐标系，适合全国尺度查看和与在线底图叠加；后续只要涉及面积、长度和坡度计算，我一般会先投影到 Albers 等积投影或对应区域的 UTM 投影，避免经纬度直接参与距离和面积运算带来误差。判断标准很简单：输出单位带“°”就只做定位，不做度量。这三步跑完，数据的基本健康状况就清楚了，接下来才能放心地做过滤、裁剪和叠加分析。
3. 在 GIS 与 Python 中过滤与叠加：从图斑选择到坡度分析
三级流域数据真正发挥价值是在和其他数据叠加的时候。单独看边界只是一个个多边形，把降水、径流、地形、重力异常挂到这些多边形上，才能回答“哪个单元缺水”“哪个区域水土流失严重”这类具体问题。这一章按三个递进场景给出可复现的操作。
3.1 按一级流域上色与分组统计：先掌握整体分布
拿到数据第一步我习惯先做一张全局图，但不用默认的单一色。按一级流域编码分组上色，能直观看到主要流域的空间分布，也能发现有没有完全不认识的图斑混进来。在 QGIS 里做上色很直接：打开图层属性，符号化选“分类”，按一级流域编码分类，颜色渐变选择一个跨度大的色带，避免相邻流域颜色分不清。保存时记得把图例标题改成“一级流域”，因为字段名对外部看图的人来说不够直观。
用 GeoPandas 做分组统计可以顺带验证层级关系：
import geopandas as gpd gdf = gpd.read_file("China_Third_Level_Basin.shp", encoding="utf-8") # 按二级流域统计图斑数和面积，观察层级分布 level2_stat = gdf.groupby(["code_level2", "name_level2"]).agg( 图斑数=("code_level3", "count"), 面积合计_km2=("area_km2", "sum") ).reset_index() print(level2_stat.sort_values("面积合计_km2", ascending=False).head(10))groupby 的两个字段分别是二级流域编码和名称，agg 里的 count 统计三级图斑数量，sum 汇总面积。输出结果里如果出现图斑数很多但面积很小的二级流域，说明这个区域分区粒度细、水系密；反过来图斑数少但面积巨大的二级流域，多半在西部干旱区。这个分布一眼就能看出数据的分区风格是否符合预期。如果你拿到的是 geojson 格式，read_file 同样能读，只是坐标系定义信息可能藏在文件头里，打印 gdf.crs 确认一下就好。
3.2 用编码前缀过滤目标区域：字符串筛选实战
做区域研究时经常只需要其中几个流域，比如只做某个二级流域或某一条三级流域的详细分析。字符串前缀过滤比空间裁剪更快，而且不会引入边界上的裁剪误差。
import geopandas as gpd gdf = gpd.read_file("China_Third_Level_Basin.shp", encoding="utf-8") # 先做一次编码清洗，防止混入空格 gdf["code_level2"] = gdf["code_level2"].str.strip() gdf["code_level3"] = gdf["code_level3"].str.strip() # 过滤出二级流域为F02的所有三级图斑 f02 = gdf[gdf["code_level2"] == "F02"] # 前缀匹配更宽松，适合编码层级不统一的情况 f02_alt = gdf[gdf["code_level2"].str.startswith("F02")] # 过滤单条三级流域：精确匹配 target = gdf[gdf["code_level3"] == "F0203"] print("F02 三级图斑数:", len(f02)) print(target[["code_level3", "name_level3", "area_km2"]])这里有两个容易踩的点。第一个是编码里可能混入空格，比如 "F02 " 和 "F02"，精确匹配会漏掉，所以先做一次 strip 清洗。第二个是 startswith 匹配的语义，startswith("F02")会把 F0200 到 F0299 都包含进来，如果你只想要 F02 这一个二级流域且编码体系里没有歧义，直接用==更严谨；如果编码体系里 F02 下面只有 F0201 到 F0209 这些三级码，startswith 也没问题。
过滤后的 GeoDataFrame 可以直接导出成新文件，方便在别的项目里复用：
f02.to_file("F02_basin.shp", encoding="utf-8")导出时注意编码参数。to_file 的 encoding 设成 utf-8 后，在部分旧版 GIS 软件里打开中文名称可能会乱码，遇到这种情况可以另存一份编码为 gb18030 的副本，交付给不同环境的合作方时都能打开。肉眼判断乱码的简单方法是看一级流域名称是不是显示成“鈥滃悗”这类字，不对就换编码重导出。
提示：过滤后如果发现图斑面积之和与预期相差很大，先看是不是漏掉了内部飞地或跨区域的碎片图斑，不要急着调整过滤条件。
3.3 叠加 DEM 做坡度与流向统计：把重力作用量化
流域分析绕不开地形。降落到地面的雨水在重力作用下沿坡度方向汇流，坡度越大，重力沿坡面的切向分量越大，径流速度越快。三级流域边界正好提供了统计单元，把 DEM 的坡度、流向信息按流域汇总，就能量化每个单元的重力驱动效应。
先做坡度分析：
import geopandas as gpd import rasterio import numpy as np from rasterio.mask import mask gdf = gpd.read_file("China_Third_Level_Basin.shp", encoding="utf-8") # 读取DEM并裁剪到某个三级流域 dem_path = "DEM_30m.tif" target = gdf[gdf["code_level3"] == "F0203"] with rasterio.open(dem_path) as src: out_image, out_transform = mask( src, [target.geometry.values[0]], crop=True, nodata=-9999 ) dem = out_image[0].astype("float32") dem[dem == -9999] = np.nan # 填充NaN参与梯度计算 dem_filled = np.where(np.isnan(dem), np.nanmean(dem), dem) # DEM分辨率 cell_size = 30 # 米 gy, gx = np.gradient(dem_filled, cell_size) slope = np.degrees(np.arctan(np.sqrt(gx**2 + gy**2))) print("平均坡度(度):", round(np.nanmean(slope), 2)) print("最大坡度(度):", round(np.nanmax(slope), 2))这段代码先把 DEM 裁剪到目标三级流域范围内，再用 numpy.gradient 计算南北和东西方向的梯度，合成坡度的正切值，arctan 转成角度后，slope 数组里每个像元都有一个坡度值。平均坡度和最大坡度是常用的重力指标：平均坡度高说明单元内地势起伏大，重力驱动的侵蚀潜力大；最大坡度可以用来判断是否发育陡崖、滑坡易发地形。
真正做水文分析时更专业的做法是调用 pysheds 或 whitebox 计算流向和汇流累积量，而不是用 np.gradient。np.gradient 适合快速看趋势，不适合严谨的汇流路径模拟。从流域管理角度看，坡度统计结果可以直接用来筛重点治理区域：平均坡度超过 25° 的三级流域，水土流失风险明显偏高，在生态流量核算和泥沙输移计算里需要单独标注。这个阈值不是硬标准，不同地区要结合土壤类型和植被覆盖调整。
4. 用三级流域做水文统计：面积加权、径流归并与重力参数化
这一章的主题是把水文数据挂到三级流域上。很多项目的核心目标一句话能说清：某个指标在每个三级流域里是多少。指标可能是年降水量、径流量、人口密度、化肥施用量，方法都一样——先空间关联，再按权重统计。本章用径流深和重力相关地形参数做例子，把完整流程跑一遍。
4.1 为什么统计单元选三级流域而不是省界
前面的翻车经历已经说明了省界的问题。再往深一层说，用三级流域做统计单元有三个好处。第一是水量闭合，三级流域边界就是分水岭，降水落在哪个单元就归哪个单元，不产生跨界的模糊地带；第二是尺度适配，气象站点、水文站的控制面积与三级流域的大小区间很接近，空间插值误差比直接插值到县级要小；第三是便于向上汇总，三级到二级到一级的编码前缀关系让任何一次聚合都只需要一条 groupby 语句。
需要注意区分两类统计口径：一类是按流域面积平均，适用于降水量、蒸散发这类面状指标；另一类是按测站控制范围归并，适用于径流量、输沙量这类点状观测。用面积平均去处理径流数据是常见误用，会把测站上游和下游的贡献混在一起。判断口径的简单方法：看原始数据的空间载体，如果是网格或面状分布按面积加权，如果是站点观测按控制区间归并。
4.2 把站点观测数据归并到三级流域：空间连接与泰森加权
用雨量站数据做示范。假设手头有一批雨量站的多年平均降水观测值，先做站点和三级流域的空间连接，再计算每个流域的面积加权降水：
import geopandas as gpd import pandas as pd gdf = gpd.read_file("China_Third_Level_Basin.shp", encoding="utf-8") # 站点表：经度、纬度、多年平均降水（mm） stations = pd.DataFrame({ "lon": [104.2, 104.8, 103.9], "lat": [31.1, 30.7, 31.5], "precip": [820, 765, 890], }) stations_gdf = gpd.GeoDataFrame( stations, geometry=gpd.points_from_xy(stations["lon"], stations["lat"]), crs="EPSG:4326", ) # 站点挂到三级流域 joined = gpd.sjoin(stations_gdf, gdf, how="left", predicate="within") print(joined[["precip", "code_level3", "name_level3"]])sjoin 的 predicate 参数用 "within"，确保只把落在流域内的站点挂上去。落在所有流域外的站点会保留在连接结果里但编码字段为 NaN，需要打印出来检查——有些站点正好落在边界线上，within 查询会失败，换成 "intersects" 可以兜底，但这样同一站点可能同时匹配两个流域，后续汇总时要做去重。
如果两个站落在同一个三级流域，需要按控制面积加权平均。这里的常见做法是用泰森多边形：把每个站点生成一个泰森多边形，代表它的控制区间，再与三级流域求交，按相交面积占比分配站点值。在 QGIS 里的流程是：站点图层生成泰森多边形，与三级流域求交，按面积比例加权汇总。Python 里用 shapely 的 voronoi_diagram 生成多边形再和流域做 intersection，逻辑比较绕，但结果更接近真实面雨量。
提示：简单平均只在站点密度很高、分布均匀时可用，正式成果建议至少做泰森加权，否则山区降水点的权重会被严重高估。
4.3 坡度、流向与汇流累积：重力在流域分析里的三张表
重力在地表水文学里的作用最终体现在三张派生数据上：坡度、流向、汇流累积量。坡度决定势能转化为动能的速度，流向决定径流往哪个方向走，汇流累积量决定某个像元上方汇集了多少格子水。用 pysheds 完整实现这三步：
from pysheds.grid import Grid import geopandas as gpd gdf = gpd.read_file("China_Third_Level_Basin.shp", encoding="utf-8") target = gdf[gdf["code_level3"] == "F0203"] # 从DEM创建网格，读取数据 grid = Grid.from_raster("DEM_clip.tif") dem = grid.read_raster("DEM_clip.tif") # 填洼，确保水流不反向 flooded_dem = grid.fill_depressions(dem) # 计算流向（D8方向编码） dirmap = (64, 128, 1, 2, 4, 8, 16, 32) flowdir = grid.flowdir(flooded_dem, dirmap=dirmap) # 汇流累积量 acc = grid.accumulation(flowdir, dirmap=dirmap) # 提取目标流域内的累积量并进行统计 acc_clip = grid.clip(acc, gdf=target) print("平均汇流累积量:", round(float(acc_clip.mean()), 2))pysheds 的核心调用顺序是：from_raster 创建网格对象，fill_depressions 填洼，flowdir 计算 D8 流向，accumulation 生成汇流累积量。dirmap 定义了流向编码到邻域方向的映射，标准 D8 编码是 64、128、1、2、4、8、16、32，对应东南、南、西南、西、西北、北、东北、东这八个方向，顺序不能写错。flowdir 输出的是一个整数网格，每个值代表该像元水流去的方向。
参数说明：fill_depressions 默认把洼地填到与最低溢出口齐平，如果 DEM 里存在真实湖泊洼地，需要先用掩膜把湖泊区保护起来，否则湖盆会被当成噪声填掉。汇流累积量的单位是栅格数，乘以像元面积再除以汇水面积，可以得到径流模数的空间分布。这套三图层组合起来能回答很多实际问题：抓取流向网格和降雨网格可以追踪污染物迁移路径；找汇流累积量超过某个阈值的像元可以自动提取河网；把累积量和坡度相乘可以粗略估算泥沙运移的功率指标。三级流域边界在这里的意义是把计算结果的“归属”钉死，每一个统计值都能回溯到具体的汇水单元。
5. 常见问题与避坑：坐标系偏移、边界缝隙与编码对照
数据本身再干净，落到不同人的手里也会撞上各式各样的坑。下面五条是我和身边同事实际踩过的，按出现的频率排了序，每一条都按现象、原因、解决写清楚。
5.1 图层加载后整体偏移，与底图错位几百米到几千米
现象：把三级流域 shp 拖进 GIS，底图是卫星影像，流域边界整体朝某个方向平移，有时偏几百米，高纬度地区甚至偏几千米。
原因：shp 缺少投影定义文件，软件按默认的经纬度坐标系读取，而底图用的是 Web 墨卡托投影，两者之间的变换没被识别，几何体直接被当作经纬度坐标绘制。
解决：先确认原始数据坐标系。大多数流域数据是 WGS84 或 CGCS2000，在 GIS 里右键图层设置坐标系，选对应 EPSG 编码；如果确实不确定，根据数据分布范围试，先试 4326 和 4490，再和底图对照。定位正确后另存一份带投影定义的新文件，避免下次打开再猜。
5.2 相邻流域之间有缝隙或重叠
现象：导出图斑后，两个相邻三级流域边界没有完全贴合，放大到高倍率能看到一条细缝；有的数据里出现小面积重叠，叠加分析时属性被重复计算。
原因：原始数据经历过多次拓扑编辑，或者不同片区拼接时没有做顶点捕捉。缝隙导致面积求和漏掉一部分，重叠则导致重复统计，两者都会污染汇总结果。
解决：在 QGIS 里用修复几何图形工具跑一遍，拓扑容差先按 0.0001 度试，不要一上来就设到 0.001，否则会把真实细小支流边界抹掉。处理完重新计算面积字段，再做面积求和校验。如果问题图斑集中在某些省域交界处，说明原始建库时拼接精度不一致，用 GRASS 的 v.clean 按 break 和 rmdupl 处理比手动编辑更效率。
5.3 三级编码和最新资料对不上
现象：按某个三级流域代码去查最新水资源资料里的数据，名称或范围对不上，周边流域编码也和之前的记录不一致。
原因：流域分区方案有过多次修订，部分三级流域被拆分或合并，旧版本数据使用了早期编码方案。
解决：以最新版分区方案为准做对照。如果手头没有对照表，用流域名称字符串匹配校验——编码会变，但主要支流的名称不会大变。我的习惯是保留一份新旧编码对照表，在项目交付报告里注明数据版本年份，避免后续合作者拿到数据后对不上号。这个坑最容易发生在跨年份项目上，旧数据只跑过一次但第二年接着用的情况特别多。
5.4 面积字段大量为0或明显偏小
现象：面积字段里不少图斑为 0，还有些偏小到个位数，用面积做权重时结果严重失真。
原因：原始面积在未投影的经纬度坐标系下直接计算，单位是平方度，数值完全不正确，或者建库时面积字段没有回填。
解决：用投影后的数据重算面积。全国范围用 Albers 等积投影，中央经线 105°E，标准纬线 25°N 和 47°N：
import geopandas as gpd gdf = gpd.read_file("China_Third_Level_Basin.shp", encoding="utf-8") # 全国范围投影到Albers等积投影，再计算真实面积 albers_crs = "+proj=aea +lat_1=25 +lat_2=47 +lat_0=0 +lon_0=105 +datum=WGS84" gdf_proj = gdf.to_crs(albers_crs) gdf["area_km2_recalc"] = gdf_proj.geometry.area / 1_000_000 print(gdf[["code_level3", "area_km2", "area_km2_recalc"]].head())数据原始坐标系如果是 CGCS2000，把 datum 参数换成对应的定义。重算后把新面积写回业务表，后续所有面积相关计算都用新字段，不再碰原始面积列。
5.5 与DEM自动提取的流域边界不一致
现象：用 DEM 和流向算法自动提取出来的小流域边界，与这份三级流域边界有出入，平原区域尤其明显，两边差出好几千米。
原因：官方三级流域边界兼顾了水文测站控制区间和部分人工约定，纯地形提取的边界以数字高程模型的分水岭为准，两者在平原地区没有明显山脊可依赖，分歧几乎必然存在。
解决：不要试图让两套边界完全一致。做水资源行政管理用官方三级流域，做精细水文建模用 DEM 提取结果。如果必须统一，用填洼后的 DEM 重新模拟分水岭，设定一个集水面积阈值，提取子流域后与三级流域做包含关系检查，找出差异最大的区域重点复核。我在一个跨流域项目中就是这么操作的，最终在双方边界差异区补了一项人工踏勘，才把争议消掉。
这五条里其实可以按问题性质分成两类：前两条是数据本身的物理伤，后三条是使用方式问题。排查时先定位是哪一类，能少走很多弯路。物理伤找源数据重新处理，方式问题调整自己的流程就行，不要盲目改数据。
6. 数据质量进阶验证：用面积守恒与排水密度反推边界合理性
拿到任何一份流域数据，我都会先花二十分钟做一次体检，不跑模型，只做三个快速检查：面积守恒、排水密度、重力异常梯度对照。这一步能挡住大部分低质量数据，避免把错误带进正式建模。
第一个检查是面积守恒。把所有三级流域面积加起来，应该接近全国陆地面积。偏差在 2% 以内算正常，超过 5% 就要怀疑有图斑缺失或属性面积没回填，用第 5 章讲过的 Albers 投影重算字段直接做。
第二个检查是排水密度。每个三级流域统计河网总长度除以面积，得到排水密度。这个值有鲜明的区域特征：湿润山区偏高，能到 0.5 到 1.2 km/km²；干旱区偏低，常低于 0.2 km/km²。如果湿润山区某三级流域排水密度只有 0.05，大概率是河网数据在这儿有空缺，不是自然现象。
import geopandas as gpd gdf = gpd.read_file("China_Third_Level_Basin.shp", encoding="utf-8") # 河网线数据，只保留在三级流域内部的河流 rivers = gpd.read_file("rivers.shp", encoding="utf-8") gdf_proj = gdf.to_crs("+proj=aea +lat_1=25 +lat_2=47 +lat_0=0 +lon_0=105 +datum=WGS84") rivers_proj = rivers.to_crs(gdf_proj.crs) # 空间连接，找出每条河归属的三级流域 joined = gpd.sjoin(rivers_proj, gdf_proj, predicate="within", how="inner") # 按三级流域汇总河网长度（km） river_len_km = joined.groupby("code_level3").geometry.length.sum() / 1000 area_km2 = gdf_proj.set_index("code_level3").geometry.area / 1_000_000 # 排水密度 = 河网长度 / 流域面积，单位 km/km² density = (river_len_km / area_km2).dropna() print("排水密度异常偏低的流域:\n", density.nsmallest(10))第三个检查更有意思：布格重力异常梯度与流域边界的对应。布格重力异常反映地壳浅层密度变化，流域分水岭往往与地形突变带或构造边界有关，这些位置常出现明显的重力梯度带。把重力异常梯度模长加载进 GIS，叠加三级流域边界，目视检查有多少边界正好落在梯度带上。重合率高的数据边界可信度更高；全部不重合的，我就要怀疑这套边界经过过度平滑，失去了地形控制意义。
三个检查加起来不到半小时，挡住的返工却不少。从那以后，我每次拿到新的流域数据都强制走一遍面积守恒和排水密度，再决定要不要继续用。重力梯度对照不常做，但做一次能明显增强对边界的信任感。希望这套快速体检帮你在下一个项目里少走弯路。
本文还有配套的精品资源，点击获取
