简介：辽宁省一级流域、二级流域矢量图层shp数据是一套面向GISer与区域规划人员的ArcMap矢量资源，覆盖黑龙江流域、辽河流域、海河流域及海岸线等一级分区，并细分至松花江水系、辽河干流水系、大凌河水系及辽东沿海诸河系、辽东半岛诸河系、鸭绿江水系、滦河水系等二级水系，可直接用于水系制图、流域边界提取与空间分析等场景。压缩包共15个文件，以shp主文件配合shx、dbf、prj、sbn、sbx等辅助文件形式组织，另含两个xml元数据文件和一个mxd工程文档，整体体积仅9.43MB，便于快速下载与应用。已有52人学习下载。依据这套数据，使用者可跳过数字化与拓扑处理环节，直接获得辽宁省分级的流域矢量底图，适用于教学演示、课题研究与水利规划辅助制图。
1. 辽宁省流域矢量数据：一张能直接进生产线的省级水系底图
做防汛一张图或者环评预审的时候，最耗时间的往往不是模型本身，而是手里的底图对不对。辽宁省一级流域、二级流域矢量图层shp数据，就是把辽河、鸭绿江、浑河这些大水系按水资源分区的一级、二级层级提前切好边界的矢量底图，文件格式就是最常见的shp。它解决三件事：水资源统计口径统一、空间裁切有现成边界、字段规整能直接入库。适合做水利信息化、GIS开发、环评预审的从业者。一个反直觉的结论：真正让这套数据翻车的不是缺线少面，而是坐标系和拓扑在静默出错。
2. 流域分级的地图密码：一级二级划分逻辑与shp属性字段
2.1 分水岭决定边界：一级与二级流域是怎么切出来的
流域划分的核心依据是分水岭。降水和径流顺着地形往低处汇，相邻两条山脊线之间的汇水区域就是一个流域。一级流域对应的是全国尺度的大水系格局，辽宁境内主要涉及辽河流域和松花江流域的边缘部分；二级流域则是在一级范围内按支流水系继续细分，比如辽河干流、浑河、太子河、大辽河、鸭绿江等这些能叫得上名字的水系。一级流域的边界通常跨省，二级流域的边界基本落在省内，这也是为什么省级项目拿二级流域当基础工作单元最顺手。
这套数据的制作流程一般是先基于数字高程模型自动提取汇水区，再用水文年鉴和实测水系线人工修编。自动提取的结果会有碎面、河道穿越分水岭等问题，人工修编要花大量时间。所以你在拿到shp时要注意：流域边界是一个完整的闭合多边形，多边形外是相邻流域，多边形的边严格贴着分水岭。如果发现边界和DEM山脊线有明显的交叉错位，那大概率是数据源或投影出了问题。
对实际项目而言，一级流域适合做宏观对标和汇报图，二级流域适合做具体的水资源量测算、面源污染负荷分摊和防汛责任区划分。两类数据叠加使用时，要留意两者不是严格的父子嵌套关系，后文避坑章节会专门说。
2.2 手把手读出属性表：编码、名称、面积字段
shp的本质是几何加属性表，属性表藏在dbase文件里。一套规范的流域shp，属性表至少包含这几个字段：要素唯一编号、一级流域名称、二级流域名称、流域编码、面积、周长，部分版本还带干流长度或河网密度。
流域编码是整张表里最值得研究的东西。常见的编码规则参考水利行业的水资源分区方案，用层次化编码表示级别关系：前几位代表一级流域，中间位代表二级流域，后面可能是三级或更细的编码。比如某个二级流域编码是A03，A暗示它属于某个一级大流域，03是它在二级序列里的序号。字段里如果同时存在数字编码和中文名称，日常筛选推荐用数字编码做关联，用中文名称做展示，因为中文名称在不同版本数据里可能有「浑河」「浑河水系」这样的写法差异，而编码相对稳定。
属性表里还会有一个AREA字段，但这里必须泼盆冷水：很多版本里的AREA是按要素原有坐标系直接算出来的，带不带投影、用什么投影算，结果都不一样。拿到数据后不要直接用这个字段做统计口径，应该用GIS工具在正确的投影坐标系下重算面积。这属于必做步骤，不是可选项。
读取shp属性表时还要注意字符编码。国内生产的shp，属性表可能是GBK，也可能是UTF-8，用错编码轻则中文乱码，重则读取中断。后面第4章会给出具体处理方式。
2.3 坐标系是真正的门槛：地理坐标与投影坐标怎么选
辽宁省一级流域、二级流域shp数据面世时，坐标系最常见的两种：地理坐标系用经纬度表示位置，通常基于CGCS2000或西安80；投影坐标系把球面展平到平面，单位是米，适合算面积和长度。做可视化展示时地理坐标系够用，叠在线底图上也不容易发现异常；一旦进入面积统计、缓冲区和overlay分析，就必须使用投影坐标系。
辽宁东西跨度大，横跨两个3度高斯-克吕格投影带，中央经线分别是120°E和123°E。选投影带的原则很简单：数据主要落在哪个带就用哪个带。辽西的朝阳、葫芦岛一带更接近120°E，沈阳、大连、抚顺等中东部地区用123°E更合适。全省统一出图或统计时，我一般会用兰伯特等积投影或阿尔伯斯等积投影做一个全省统一工作坐标，保证面积口径一致。
顶级的坑在于：很多数据发布时.prj文件是丢的，软件读不到坐标系信息，默认当成WGS84，而数据实际是西安80或CGCS2000，两种坐标系在辽宁地区平面位置差几十米到上百米。这个误差叠加到底图上不明显，但和周边省份数据或实地测绘数据一叠，立刻露馅。判断方法是用已知位置的河流控制点去做交叉验证，别靠感觉。
3. 拿到数据先做三件事：坐标系、拓扑与属性校验的完整流程
3.1 坐标系统一的参考命令：geopandas 读取与重投影
拿到shp后第一步不是打开看颜色，而是确认坐标系。用Python的geopandas可以在几分钟内完成批量的读取、检查和转换。以下是我处理辽宁省流域数据时的常规操作。
import geopandas as gpd gdf = gpd.read_file('liaoning_basin_level12.shp', encoding='utf-8') print(gdf.crs) # 看看有没有坐标系信息 # crs为空时，先根据元数据补上坐标系，再转换 if gdf.crs is None: gdf = gdf.set_crs('EPSG:4490') # CGCS2000地理坐标系 # 统一到CGCS2000经纬度坐标系 gdf = gdf.to_crs('EPSG:4490') print(gdf.head())这段逻辑的关键在于：先判断再转换，不要一上来就to_crs。如果crs本身就是空的，直接to_crs会被geopandas报错；如果原始数据是西安80却被当成CGCS2000，转换后整个边界会横向偏移，后期所有分析结果全部错位。EPSG:4490是CGCS2000地理坐标系的标准代码，我用它作为中间统一坐标系，之后要算面积再转到投影坐标系。
参数说明：encoding='utf-8'解决属性表中文乱码问题；如果读取报错，尝试把编码改成gbk，这是国内shp最常见的两种属性编码。gdf.crs输出类似<Geographic 2D CRS: EPSG:4490>，看到这个就说明坐标系信息完整。
如果确认数据是CGCS2000且需要转到投影坐标，可以这样：
# 辽宁中东部按中央经线123E的3度带转投影，单位米 gdf_proj = gdf.to_crs('+proj=tmerc +lat_0=0 +lon_0=123 +k=1 +x_0=500000 +y_0=0 +ellps=GRS80 +units=m +no_defs') print(gdf_proj.geometry.area.max() / 1e6, 'km²') # 用平面坐标直接算面积这个Proj4字符串指定了横轴墨卡托投影、中央经线123°E、克拉索夫斯基椭球被GRS80替代的CGCS2000基准。实际项目中如果你不确定该用哪一个带，就去看数据里河流的经度范围取中值，然后就近选带。
3.2 拓扑检查与修复：make_valid 和 buffer 的双保险
流域shp最容易出的拓扑问题是面与面之间的边界不闭合、自相交和微小缝隙。这些几何错误在画面上几乎看不见，但一跑overlay、求交、裁切就会触发GEOS拓扑异常，严重时整个分析直接中断。我一般用shapely的make_valid做批量修复。
from shapely.validation import make_valid # 检查无效几何 invalid_mask = ~gdf.geometry.is_valid print(f"无效几何数量: {invalid_mask.sum()}") # 全部几何过一遍修复 gdf['geometry'] = gdf.geometry.apply(make_valid) # 面要素之间如果存在细微缝隙，用buffer(0)焊一下 gdf['geometry'] = gdf.geometry.buffer(0)make_valid会把自相交的面拆成多个有效几何，把退化线段清理掉。它处理不了的是两个相邻流域共享边界上的「线缝」——相邻面的边界如果没精确重合，流域之间会有很细的空隙。此时buffer(0)的作用是不改变几何形状的前提下强制闭合拓扑。注意不要用正缓冲值，正缓冲会扩大边界，造成相邻流域之间的重叠。
逻辑说明：这一步在坐标系转换完成之后做，因为投影转换本身就可能产生新的几何自相交，先转坐标再修拓扑才是正确顺序。修复完成后建议重新导出shp作为工作版本，不破坏原始下载文件，以便随时回溯。
3.3 属性完整性核对：空值、重复编码与边界版本
数据能画出来不代表属性干净。shp的属性表经常出现编码重复、名称字段空值、面积字段与几何长度不一致的问题。下面这段代码用来快速验收属性质量。
# 空值检查 print(gdf.isnull().sum()) # 重复编码检查 dup = gdf[gdf.duplicated(subset=['流域编码'], keep=False)] print(f"重复编码记录: {len(dup)}") # 面积重算：投影坐标系下算真正的km² gdf_proj['calc_area_km2'] = gdf_proj.geometry.area / 1e6 # 对比属性表自带面积和重算面积 gdf_proj['area_diff_pct'] = ( (gdf_proj['calc_area_km2'] - gdf_proj['AREA']) / gdf_proj['AREA'] * 100 ) print(gdf_proj['area_diff_pct'].describe())逻辑说明：duplicated(subset=['流域编码'])找出编码重复的要素，这类问题多出现在多个县区级数据拼接成的省级图层中，拼接时县级各自编号，没有做全局唯一化。calc_area_km2除以1e6是因为投影坐标单位是米，面要素的面积单位是平方米，转成平方公里方便对比。
参数说明：如果area_diff_pct的绝对值大于5%，说明自带的AREA字段坐标系和当前投影有明显差异，后续统计绝对不能直接引用原字段。我把这一步叫「属性体检」，体检不过关的数据入库就是埋雷。空值字段如果集中在次要字段如备注、来源，可以留空；如果集中在流域名称或编码，必须回查原始数据源补录。
4. 把流域shp推进业务管线：属性筛选、dwg转shp与叠加统计
4.1 按水系名称筛选：快速提取特定流域的工作面
业务上经常会问「浑河流域的覆盖范围是哪些」「太子河流域内的行政区有哪些」，这些需求本质都是属性筛选加空间查询。用geopandas做筛选并导出，代码很直接。
huntai = gdf[gdf['二级名称'].str.contains('浑河')] print(huntai['二级名称'].unique()) # 导出为独立shp huntai.to_file('huntai_basin.shp', encoding='utf-8') # 把属性导出为txt/csv给外部系统用 huntai[['流域编码', '二级名称', 'calc_area_km2']].to_csv( 'huntai_basin_attr.txt', sep='|', index=False )逻辑说明：str.contains('浑河')是模糊匹配，能匹配到「浑河」「浑河水系」「浑河上游」等名称，比全等匹配更实用。导出txt时用|做分隔符是行业惯例，因为流域名称里可能带逗号或顿号，逗号分隔容易串列。
参数说明：如果模糊匹配误伤太多不相干流域，改用str.startswith('浑')做前缀匹配；如果要精确取某个二级流域，直接用流域编码 == '具体值'。导出shp时enconding同样要注意，下游如果是在ArcGIS里打开，用utf-8没问题；如果下游是老的MapGIS或自研系统，建议导出为gbk。
这一步做完，你就有了一个目标流域的独立工作面。下一步所有叠加分析都在这个面上做，不用每次全图层扫描。
4.2 dwg转shp的常见路径：CAD线稿变成能用的面
实际项目的可研阶段经常只提供CAD图纸，水系是dwg里的多段线，没有属性、没有闭合。把这类线稿转成能和流域shp叠加分析的面要素，是每个做水利信息化的人都要面对的活儿。dwg是闭源格式，常见做法是先批量转成dxf再用Python解析，或者直接用GIS桌面软件打开后另存为shp。下面演示用ezdxf读取dxf并生成shp的过程。
import ezdxf import geopandas as gpd from shapely.geometry import Polygon doc = ezdxf.readfile('water_line.dxf') polys = [] for e in doc.modelspace().query('LWPOLYLINE'): pts = [(p.dxf.location.x, p.dxf.location.y) for p in e.get_points()] # 只保留闭合且至少4个点的多段线作为候选水面 if len(pts) >= 4 and e.closed: polys.append(Polygon(pts)) gdf_cad = gpd.GeoDataFrame(geometry=polys, crs='EPSG:4490') gdf_cad.to_file('cad_water_polygon.shp', encoding='utf-8')逻辑说明：get_points()返回多段线的顶点坐标，e.closed判断是否闭合。CAD制图不规范时经常出现「看起来闭合但实际首尾点不重合」的线，这类线不会被e.closed识别，导出后面积会缺失。更稳妥的判断指标是首点与尾点距离小于容差。
参数说明：LWPOLYLINE是CAD轻量多段线类型，如果是旧版POLYLINE，查询时要写成POLYLINE或直接LINE。如果你的CAD水系是带高程的3D多段线，查询对象还要过滤掉z轴变化过大的线，否则生成的面会扭曲。dwg转shp后最关键的一步：把生成的面和流域shp叠加，人工核对有没有跨越流域边界的对象——如果有，说明CAD里那条河画错了位置或坐标系对不上，不要急着入库。
4.3 叠加行政区统计：面积分摊算法要这样做
流域和行政区从来不是互相包含的关系。一个流域跨多个县，一个县也同时跨多个流域。最常见的需求是算「某流域在各县内的面积」，做法是空间求交后按面积分摊。用gpd.overlay实现如下。
# admin是县级行政区shp，gdf_proj是投影坐标系下的流域数据 intersect = gpd.overlay(gdf_proj, admin, how='intersection') # 每个相交碎片的面积，单位km² intersect['part_area_km2'] = intersect.geometry.area / 1e6 # 按二级流域+县名汇总 result = intersect.groupby(['二级名称', '县名'])['part_area_km2'].sum().reset_index() print(result.head()) # 再按流域汇总，得到每个流域的总面积 basin_total = intersect.groupby('二级名称')['part_area_km2'].sum().reset_index()逻辑说明：overlay要求输入的两个GeoDataFrame坐标系一致，这里提前用投影坐标是关键。求交后每个县域内的流域碎片单独成一条记录，用groupby按流域加县聚合，就能得到「流域—县」的二维面积表，这是水环境容量计算和面源污染负荷分摊的常见输入格式。
参数说明：how='intersection'是取两图层相交的部分。如果想同时保留流域内未覆盖到行政区的面积，改用how='identity'或先做一个difference再合并，但多数项目只需要交集。另外很多分析会把研究区切成规则渔网shp，再做同样的overlay求交，可以先把两个图层合并成一个统一网格框架，再做聚合统计。渔网分割shp的做法在网格化污染源排查里很常见，思路和上面完全一样。
5. 避坑记录：这几个环节最容易让流域shp翻车
5.1 排查顺序与工具链
数据出问题不要靠肉眼扫图，要有固定排查顺序。我自己的工具链是：先用Python批处理完成坐标系、拓扑、属性三项体检，再用QGIS做目检，最后用桌面GIS的拓扑检查器跑一遍面重叠和缝隙检测。一套流程下来不超过十分钟，但能避免后面几天的返工。
QGIS的拓扑检查器能高亮所有几何错误，但它的问题是无法告诉你错误发生在哪个环节。所以我建议先跑命令再开可视化：命令行负责批量定位，QGIS负责确认和修编。遇到不确定的坐标系，用一个位于辽宁省中部的已知水文站坐标做锚点，检查shp里对应点的经纬度，这个方法比任何参数文档都可靠。
5.2 五个具体踩坑场景
场景一：坐标系静默偏移。现象是把流域shp叠加到在线底图，河流和影像对不上，整体偏移几十到上百米。原因是.prj文件缺失，软件默认按WGS84加载，而数据实际是西安80或CGCS2000。解决：先读.prj，缺失就用已知坐标的河流控制点反推，不要猜。用控制点做空间校正后再使用。
场景二：overlay报错GEOSGeometry拓扑异常。现象是gpd.overlay跑一半突然报TopologyException: side location conflict。原因是相邻流域的共享边界有微小自相交或缝隙，面积越大越容易触发。解决：先执行make_valid，再执行buffer(0)，然后重新导出shp再做overlay。这三个步骤缺一不可，跳过任何一个都可能复发。
场景三：面积和预期差一大截。现象是同一块流域，图上看起来挺大，统计出来面积却明显偏小或偏大。原因是用经纬度坐标系直接算面积，单位是度，换算到平方公里时系数选错。解决：先转投影坐标系再算面积；全省统口时用等积投影。属性表里的自带AREA字段仅供参考，绝不能直接出报告。
场景四：属性表中文全是乱码。现象是QGIS打开shp属性表，中文变成了æµæ²³一类的乱码。原因是shp属性文件是GBK编码，读成UTF-8了。解决：读取时指定encoding='gbk'；QGIS里在图层属性里修改数据源编码。这个坑几乎每个从老GIS平台拿到的数据都会遇到，不值得花时间研究为什么，直接转码就完了。
场景五：一级与二级图层嵌套对不上。现象是二级流域面积汇总后比对应一级流域面积少10%以上。原因是两个图层可能分别来自不同年份、不同精度的数据源，部分支流边界在修编后发生了变化。解决：做层级一致性校验，把二级流域按所属一级编码汇总后与一级图层求差，找到差异区域后以较新版本为准，或者统一重新裁切。我在辽宁一个地级市项目里就因为这个返工过一次，两个版本都是正规来源，但修编年份差了三年。
6. 进阶玩法：用流域边界裁栅格、向三维场景输送shp
流域shp不只是画图用的，它更值钱的地方在于作为空间分析的边界条件。最常见的进阶操作是用流域边界去裁切降雨、蒸散发、DEM这样的栅格数据，从而快速得到「这个流域平均下了多少雨」「这个流域高程分布如何」这类结论。用rasterio的mask接口可以一次性完成。
import rasterio from rasterio.mask import mask as rio_mask with rasterio.open('rainfall.tif') as src: # 取某个流域的外边界做裁切 geom = huntai_proj.geometry.iloc[0] out_img, out_transform = rio_mask( src, [geom], crop=True, nodata=0 ) print(out_img.shape, '像素数')逻辑说明：rio_mask要求栅格和矢量处于同一坐标系，裁切前最好确认两者都转到了投影坐标。nodata=0把裁切范围外的区域填为0，后续做像元统计时可以直接过滤。算出每个像元的面积乘上降雨值再求和，就是该流域的降雨总量。
三维场景方向，很多人问shp能不能直接转3dtiles。shp是二维矢量格式，不能直接用于三维场景，常见做法是把shp先转成GeoJSON，再转成带高程的glTF模型，最后切片成3dtiles。这里有个血泪经验：流域边界在二维平面看着很规整，一加地形高程，边界会被地形起伏拉扯，切片后视觉效果很糟。我的做法是在三维发布前先对边界做低通平滑，丢掉过于细碎的折点，保证贴地形时不会出现锯齿状的撕裂。
做这套数据做久了，最大的体会是坐标系和拓扑永远要比业务逻辑先到位。有一次我把二级流域shp的经纬度面积直接当成平方公里报给了规划院，整个区域算小了近三成，被对方退回重算，当场脸都丢尽了。自那以后，我接手任何shp的第一个动作永远是把crs打出来看一眼，这个习惯救了我很多次。希望帮到你。
本文还有配套的精品资源，点击获取
