全国矢量地图SHP文件实战:Python处理行政、路网、河网与水系 简介这份全国矢量地图SHP文件包面向GIS从业者、地理科研人员及城市规划学习者提供覆盖全国范围的行政边界、路网、河网与水系空间数据可用于地图制作、空间分析与规划决策等场景。压缩包共113个文件约10.49MB以shp、shx、dbf、prj、sbn、sbx等Shapefile标准组件为主分别承载几何对象、空间索引、属性表与坐标系统定义另含少量adf、xml、mxd等辅助文件便于在ArcGIS、QGIS等软件中直接导入与可视化。资源涵盖省级、县级行政区划多边形地级市与县级驻地位置点以及各级河流水系和公路铁路线要素属性中带有名称、行政级别、里程等信息可支撑区域分布研究、交通网络分析与洪水预警等工作。目前已有5828人学习下载适合需要完整全国基础地理底图的中高级GIS用户参考使用。1. 全国矢量地图 SHP 文件到底能干什么从一次区划底图翻车说起做空间分析的人迟早会撞上一个场景手头有一份业务数据带经纬度想按行政区汇总或者想看看点位落在哪条河、哪条路附近。这时候你打开常用地图服务发现它只给你瓦片不给你边界你想下载行政边界要么精度不够要么字段残缺要么干脆是加密格式。全国矢量地图 SHP 文件这个方向解决的正是这件事——把行政边界、路网、河网、水系这几类基础地理要素以 Shapefile 这种最通用的矢量格式落到本地让你能离线、可复现地做裁剪、叠加、统计和出图。它适合三类人做国土空间、交通、水利、环保方向数据分析的从业者需要给业务系统配底图、又不想依赖在线服务的后端工程师以及做论文或课题、需要一套干净矢量底图的研究人员。SHP 是 Esri 定义的老牌矢量格式一个图层由 .shp、.shx、.dbf 三个必需文件加若干可选文件组成几乎所有 GIS 软件和 Python 生态都能直接读。全国范围的数据量不小行政边界到县级通常几十万面路网和水系动辄百万级线要素所以拿到 .rar 之后第一件事不是急着打开而是先想清楚坐标系、字段和精度这三件事否则后面每一步都在还债。2. 拿到压缩包先别急着解压SHP 文件结构与坐标系选型2.1 一个 SHP 图层到底由哪些文件组成很多人第一次拿到 SHP 会懵解压出来一堆同名不同后缀的文件删错一个整个图层就打不开。先把这套结构记牢后面排错全靠它。文件后缀是否必需作用缺失后果.shp必需存储几何形状图层无法读取.shx必需几何索引加速定位部分软件报错或读取极慢.dbf必需属性表存字段数据只剩图形没有属性.prj强烈建议坐标系定义WKT坐标含义丢失叠加错位.cpg可选属性表字符编码中文属性乱码.sbn/.sbx可选空间索引仅影响查询性能血泪经验.prj 丢了是最坑的因为文件还能打开图形也画得出来但你不确定它是 WGS84 经纬度还是某种投影坐标一旦拿去和别的数据叠加偏移可能几百米甚至几公里而且这种错误很隐蔽出图看着差不多一量距离全错。2.2 地理坐标系和投影坐标系全国数据该选哪个这是选型的核心。全国矢量地图常见两种坐标状态地理坐标系GCS单位是度常见 CGCS2000 或 WGS84。优点是通用、便于存储和交换缺点是直接算面积、算长度会失真因为一度经度在不同纬度对应的实际距离不一样。投影坐标系PCS单位是米常见各种高斯-克吕格 3 度带或 6 度带以及 Albers 等面积投影。优点是算面积、算距离准确缺点是全国跨多个带单带投影边缘变形大。我的建议是存储用地理坐标系分析前按需投影。全国尺度的面积统计用 Albers 等面积投影中央经线 105°E双标准纬线常见 25°N 和 47°N比高斯分带更合适因为它保证全国范围内面积变形可控。而做局部路网长度、河网长度计算时用对应 3 度带高斯投影精度更高。用 Python 检查一个图层的坐标系最省事的办法是读 .prj 或直接用 GeoPandasimport geopandas as gpd # 读取一个 SHP 图层GeoPandas 会自动解析 .prj gdf gpd.read_file(boundary_county.shp) # 看坐标系定义确认是地理坐标还是投影坐标 print(CRS:, gdf.crs) # 看几何类型和要素数量判断数据规模 print(几何类型:, gdf.geom_type.unique()) print(要素数量:, len(gdf)) # 看属性字段确认有没有行政区代码、名称这类关键字段 print(字段列表:, list(gdf.columns))这段代码的逻辑是先确认坐标系决定后面要不要投影再确认几何类型面/线/点决定能做什么分析最后看字段决定能不能按行政区关联。参数上没什么可调的重点是gdf.crs的输出——如果显示None说明 .prj 缺失你必须自己判断坐标系并手动指定否则后续所有空间操作都不可信。2.3 解压与目录组织别把四类数据混在一起全国矢量地图通常包含行政、路网、河网、水系四类每类可能又按层级或区域分多个图层。解压后建议按下面的结构组织避免文件名冲突和路径混乱# 假设压缩包解压到 data_raw 目录 mkdir -p data/{admin,road,river,water} # 行政边界按层级放命名体现层级 mv data_raw/boundary_province.* data/admin/ mv data_raw/boundary_city.* data/admin/ mv data_raw/boundary_county.* data/admin/ # 路网、河网、水系各自独立目录 mv data_raw/road_*.shp* data/road/ mv data_raw/river_*.shp* data/river/ mv data_raw/water_*.shp* data/water/注意mv用了通配符.*是为了把 .shp/.shx/.dbf/.prj 等同名文件一起搬走。SHP 的坑就在于它是多文件一体只搬 .shp 等于把图层搬废了。目录分好后后面写批处理脚本时路径清晰也方便按类别做不同的投影和简化处理。3. 用 Python 把行政、路网、河网、水系跑通一遍3.1 环境准备与依赖版本工欲善其事先把环境弄干净。SHP 读取依赖 GDAL而 GDAL 的安装是新手最容易翻车的地方。推荐用 conda 或 mamba 装比 pip 省心得多。# 用 conda 创建独立环境避免和系统 GDAL 冲突 conda create -n gis python3.10 -y conda activate gis # 核心依赖geopandas 读矢量shapely 做几何运算pyproj 管投影 conda install -c conda-forge geopandas shapely pyproj fiona rtree -y # 可选matplotlib 出图pandas 做属性统计 conda install -c conda-forge matplotlib pandas -y参数说明-c conda-forge指定社区源GDAL 相关包在这里版本最全、依赖最干净。rtree是空间索引库做叠加分析比如判断点是否落在某行政区时如果没有它速度会慢到无法接受。如果你坚持用 pip务必先单独装好 GDAL 的 wheel再装 geopandas否则大概率卡在编译报错上。3.2 读取四类图层并统一坐标系拿到数据后第一步不是分析是统一坐标系。四类数据来源可能不同坐标系未必一致直接叠加必然错位。import geopandas as gpd # 统一目标坐标系CGCS2000 地理坐标便于存储和交换 TARGET_CRS EPSG:4490 def load_and_reproject(path, target_crsTARGET_CRS): 读取 SHP 并统一到目标坐标系 gdf gpd.read_file(path) # 如果原数据没有坐标系定义这里会报错需要先手动指定 if gdf.crs is None: raise ValueError(f{path} 缺少坐标系定义请先确认原始 CRS) # 只有坐标系不同才转换相同则跳过省时间 if gdf.crs ! target_crs: gdf gdf.to_crs(target_crs) return gdf # 分别读取四类数据 admin load_and_reproject(data/admin/boundary_county.shp) road load_and_reproject(data/road/road_national.shp) river load_and_reproject(data/river/river_main.shp) water load_and_reproject(data/water/water_polygon.shp) # 打印各自范围确认覆盖全国 for name, gdf in [(行政, admin), (路网, road), (河网, river), (水系, water)]: print(name, 范围:, gdf.total_bounds)逻辑说明to_crs做的是坐标转换不是简单改标签它会真正重算每个顶点的坐标。total_bounds返回[minx, miny, maxx, maxy]全国数据大致落在经度 73~135、纬度 3~54 这个范围地理坐标下。如果某个图层的范围明显偏离要么坐标系判断错了要么数据本身有问题必须先查清楚再往下走。3.3 按行政区裁剪路网与河网这是最常见的需求我只关心某个省或某个市范围内的路网、河网。用 GeoPandas 的clip一步到位。# 假设要裁剪出某个省的范围先按名称筛出该省 province admin[admin[name] 某省] # 用省界裁剪路网只保留落在省内的部分 road_clipped gpd.clip(road, province) # 裁剪河网同理 river_clipped gpd.clip(river, province) # 裁剪后要素数量通常会减少打印对比 print(路网裁剪前:, len(road), 裁剪后:, len(road_clipped)) print(河网裁剪前:, len(river), 裁剪后:, len(river_clipped)) # 保存结果注意 keep_geom_type 保持几何类型一致 road_clipped.to_file(data/out/road_in_province.shp, encodingutf-8) river_clipped.to_file(data/out/river_in_province.shp, encodingutf-8)参数说明gpd.clip会把跨越边界的线在边界处切断只保留内部部分这正是我们要的。to_file的encodingutf-8很关键不写的话中文属性在部分软件里会乱码。注意裁剪是计算密集型操作全国路网裁剪单省要素百万级时可能要几分钟耐心等别以为卡死了。3.4 统计每个行政区的路网密度有了裁剪结果就能做业务统计了。路网密度 区内道路总长度 / 区面积是交通分析的常见指标。# 先投影到等面积投影否则长度和面积都不可信 # Albers 等面积投影适合全国尺度 albers projaea lat_125 lat_247 lat_00 lon_0105 x_00 y_00 datumCGCS2000 unitsm admin_m admin.to_crs(albers) road_m road.to_crs(albers) # 计算每个行政区的面积平方米转平方公里 admin_m[area_km2] admin_m.geometry.area / 1e6 # 空间连接把每条道路归属到它所在的行政区 road_join gpd.sjoin(road_m, admin_m, predicatewithin, howleft) # 按行政区汇总道路长度 road_join[length_km] road_join.geometry.length / 1000 density road_join.groupby(name)[length_km].sum().reset_index() # 合并面积算密度 density density.merge(admin_m[[name, area_km2]], onname) density[density] density[length_km] / density[area_km2] print(density.sort_values(density, ascendingFalse).head(10))逻辑说明sjoin的predicatewithin表示道路完全落在行政区内才归属跨界的道路会被切分后分别归属所以前面裁剪那步其实可以省但先裁后连更直观。groupby按行政区名汇总长度再除以面积得到密度。这里有个坑如果行政区名有重名比如多个某城区groupby 会把它们合并所以生产环境应该用行政区代码而不是名称做关联键。4. 避坑指南SHP 处理里最容易翻车的五件事4.1 中文属性乱码字段名变成问号现象打开 SHP属性表里中文全是乱码或者字段名直接显示成????。原因.dbf 文件的字符编码没有声明或者声明了但和实际编码不符。老数据常见 GBK新数据多为 UTF-8而 .cpg 文件缺失时软件只能猜。解决先确认编码再补 .cpg 或用参数指定。读取时显式指定编码# 尝试用 GBK 读取如果乱码再换 UTF-8 gdf gpd.read_file(data/admin/boundary_county.shp, encodinggbk)如果读进来还是乱码用encodingutf-8再试。确定编码后写一个同名 .cpg 文件内容就一行编码名后续软件就能自动识别。4.2 叠加分析结果整体偏移几百米现象把路网和行政边界叠在一起发现道路整体偏出边界或者点位落在河里。原因两个图层坐标系不一致或者其中一个 .prj 缺失被误判。地理坐标和投影坐标混用是最常见的。解决叠加前强制统一坐标系并且用gdf.crs逐个确认不要假设。如果某个图层 crs 是 None先查清它的真实坐标系再手动set_crs绝不能随便猜一个。4.3 全国数据读取内存爆掉现象读全国路网时进程被系统杀掉或者报 MemoryError。原因SHP 格式本身不支持高效的空间索引和分块读取百万级要素一次性载入内存加上 GeoPandas 的几何对象开销几个 G 内存瞬间见底。解决分块处理或者先按区域筛选再读。用 Fiona 逐要素读取或者用bbox参数只读感兴趣范围# 只读取某个经纬度范围内的要素大幅降低内存占用 gdf gpd.read_file(data/road/road_national.shp, bbox(116.0, 39.0, 117.0, 40.0))bbox是(minx, miny, maxx, maxy)只返回与该范围相交的要素。做局部分析时这是最有效的省内存手段。4.4 裁剪后几何类型变了面变成线现象裁剪水系面数据后结果里混进了线要素后续面积统计报错。原因clip在边界处切割时极窄的面可能退化成线或者原始数据本身就混了几何类型。解决裁剪后过滤几何类型只保留需要的# 只保留面要素 water_poly water_clipped[water_clipped.geom_type Polygon]同时检查原始数据如果 SHP 里混了多种几何类型SHP 规范其实不允许但劣质数据常见最好先拆分再处理。4.5 属性关联时字段类型不匹配关联全空现象用行政区代码关联两张表结果全匹配不上关联结果为空。原因一边代码是字符串110000另一边是整数110000或者一边有前导零一边没有。解决关联前统一类型字符串转字符串并去掉首尾空格# 统一转成字符串并去空格避免类型和空格导致的匹配失败 admin[code] admin[code].astype(str).str.strip() stats[code] stats[code].astype(str).str.strip() merged admin.merge(stats, oncode, howleft)这个坑极其隐蔽因为数据看着都对就是关联不上排查时优先怀疑类型和空格。5. 进阶把 SHP 转成更适合分析的格式以及一个提速技巧SHP 有个绕不开的硬伤单文件 2GB 上限字段名最长 10 个字符不支持存储 NULL 和日期时间类型全国级数据用起来处处受限。所以真正做分析时我一般会把 SHP 转成 GeoPackage 或 Parquet前者是 OGC 标准、单文件、支持长字段名和索引后者读取速度极快适合反复迭代。# SHP 转 GeoPackage保留坐标系和属性突破 2GB 限制 gdf.to_file(data/out/admin.gpkg, layercounty, driverGPKG) # 转 Parquet读取速度比 SHP 快一个数量级适合频繁分析 gdf.to_parquet(data/out/admin.parquet) # 下次直接读 Parquet省去解析 SHP 的开销 gdf_fast gpd.read_parquet(data/out/admin.parquet)参数说明driverGPKG指定输出格式layer是图层名一个 gpkg 文件可以装多个图层比一堆 SHP 文件清爽得多。Parquet 不支持空间索引的持久化但读取快适合中间结果缓存。再给一个提速技巧做空间连接sjoin之前先给两个图层建空间索引。GeoPandas 在 sjoin 时会自动建但如果你反复用同一个图层做连接手动建一次能省不少时间# 手动建空间索引反复做空间查询时提速明显 admin_m.sindex # 首次访问即构建后续查询复用 road_m.sindex # 之后再做 sjoin 或 within 查询会快很多 result gpd.sjoin(road_m, admin_m, predicatewithin, howleft)最后说个我自己的习惯拿到任何一份全国矢量数据先做三件事——确认坐标系、确认字段、抽样看几何是否合法用gdf.is_valid检查。这三步花不了五分钟但能挡掉后面百分之八十的翻车。数据这东西脏起来是没有下限的早发现早处理比分析到一半才发现问题强得多。希望帮到你。本文还有配套的精品资源点击获取