资讯动态

山东SHP包实战:省、市、县边界与路网数据处理及空间分析

发布时间:2026/10/14 13:22:29 来源:尧图企业网站定制
简介这份资源面向GIS从业者、城市规划与交通研究人员及地理信息相关专业学生提供2022年7月山东省的矢量地理数据可用于行政区划制图、路网分析、空间查询与叠加运算等场景。压缩包共33个文件约37.59MB以shp、dbf、prj、shx、xml等格式为主其中shp存储省、市、县三级行政边界及道路、铁路的几何形状dbf保存名称、等级、里程等属性prj定义坐标系统shx与xml分别承担索引和元数据说明各类型配合构成完整可用的GIS数据集。目前已有2018人学习下载。读者可借助ArcGIS或QGIS直接加载完成行政边界可视化、交通网络规划、灾害救援路线评估、物流选址等分析并基于属性表开展统计与制图为区域研究与决策提供可靠的空间数据支撑。1. 从一份 2022.07 的山东 SHP 包说起省、市、县边界加路网到底能干什么如果你做过山东相关的 GIS 项目大概率遇到过这种局面底图边界对不上、路网数据缺属性、铁路和公路混在一个图层里分不开。这份 2022.07 的山东省行政区划 SHP 包把省、市、县三级行政边界和道路网、公路网放在了一起解决的正是「底图从哪来、路网怎么配」这两个最耗时间的前置问题。它适合做区域分析、选址评估、交通可达性建模的从业者也适合刚接触 GIS 想拿真实数据练手的新手。我拿到这类包的第一反应不是直接打开而是先看坐标系和字段结构——这两项决定了后面能不能省掉返工。2. 拆包先看坐标系与字段SHP 能不能直接用全在这两步2.1 坐标系不确认后面全是玄学SHP 文件本身不带坐标系信息靠的是同目录下的.prj文件。很多人拿到包直接拖进软件就开始画图结果量出来的距离和实际差了几十倍这就是坐标系没认对。国内常见的行政边界数据坐标系无非几种WGS84 地理坐标系EPSG:4326、CGCS2000EPSG:4490、或者投影后的高斯克吕格带号坐标。2022 年的数据CGCS2000 是主流但也不排除部分来源用 WGS84。确认方法很简单用 GDAL 命令行看一眼# 查看 prj 文件内容确认坐标系定义 ogrinfo -al -so shandong_boundary.shp # 如果 prj 缺失或可疑直接读 prj 文本 cat shandong_boundary.prjogrinfo -al -so会输出图层的元信息包括几何类型、要素数量、字段列表和坐标系。如果.prj文件缺失输出里坐标系一栏会是空或者unknown。这时候别猜用 QGIS 打开后手动指定一个候选坐标系叠加一个已知坐标系的底图比如在线影像看能不能对上。对不上就换下一个候选直到吻合。提示如果数据是投影坐标单位是米而你要做的是跨省或全国尺度的分析建议先统一转成地理坐标系再做后续处理避免不同带号之间拼接错位。2.2 字段结构决定你能做什么分析行政区划图层一般包含NAME名称、CODE行政区划代码、LEVEL层级这几个核心字段。道路网图层则通常有NAME道路名称、TYPE道路类型如高速、国道、省道、LEVEL等级等。拿到数据后先用ogrinfo把字段列出来# 列出属性表字段和前几条记录 ogrinfo -al shandong_road.shp | head -50这一步的目的是确认三件事一是名称字段是否规范有没有乱码、空格、全半角混用二是道路类型字段的取值有哪些决定你能不能按高速/国道/省道分层渲染三是行政区划代码字段是否完整后面做空间关联时靠它。如果发现名称字段有乱码大概率是编码问题。SHP 的属性表默认用 UTF-8 或 GBK用ogrinfo加-oo ENCODINGGBK试试ogrinfo -al -oo ENCODINGGBK shandong_boundary.shp能正常显示中文就说明原数据是 GBK 编码后续在 QGIS 或 ArcGIS 里加载时也要对应设置编码否则中文会变成问号。3. 把边界和路网用起来从加载到空间分析的完整链路3.1 在 QGIS 里加载并分层渲染QGIS 对 SHP 的支持最省心直接把.shp拖进去就行。但要注意如果同目录下有.dbf、.shx、.prj要保证它们和.shp同名且在同一目录否则会报错或丢属性。加载后先做分层渲染。行政区划按LEVEL字段分类省级一个样式、市级一个、县级一个颜色由浅到深或由粗到细。道路网按TYPE字段分类高速用红色粗线、国道用橙色、省道用黄色这样一眼就能看出路网密度和等级分布。# 用 geopandas 快速查看图层的基本统计信息 import geopandas as gpd # 读取行政区划和路网 boundary gpd.read_file(shandong_boundary.shp) road gpd.read_file(shandong_road.shp) # 查看坐标系 print(边界坐标系:, boundary.crs) print(路网坐标系:, road.crs) # 查看字段和要素数量 print(边界字段:, boundary.columns.tolist()) print(边界要素数:, len(boundary)) print(路网字段:, road.columns.tolist()) print(路网要素数:, len(road)) # 按层级统计行政区划数量 if LEVEL in boundary.columns: print(boundary[LEVEL].value_counts())这段代码的作用是快速摸底。boundary.crs和road.crs如果不一致后面做空间分析会直接报错或结果错位。len(boundary)告诉你省、市、县各有多少个要素正常来说省级 1 个、市级 16 个左右、县级 130 多个。如果数量对不上说明数据可能不完整或者层级字段有问题。3.2 用空间连接把路网挂到行政区划上一个常见需求是统计每个市或每个县的路网密度。做法是用空间连接spatial join把每条道路归属到它所在的行政区划。# 空间连接把道路挂到所属行政区划 # 先确保两者坐标系一致 if boundary.crs ! road.crs: road road.to_crs(boundary.crs) # 用 intersects 做空间连接 joined gpd.sjoin(road, boundary, howleft, predicateintersects) # 按行政区划名称统计道路数量 road_count joined.groupby(NAME_left).size().reset_index(nameroad_count) print(road_count.sort_values(road_count, ascendingFalse).head(10))sjoin的predicate参数决定判断规则intersects是只要相交就算within是道路完全在行政区划内才算。对于跨界的道路用intersects更合理但会导致一条路被算进多个行政区。如果要做精确统计可以先按行政区划裁剪路网再统计。howleft保证所有道路都保留没有匹配到行政区划的会标记为 NaN这些通常是数据边界外的道路或者坐标系不一致导致的漏匹配。3.3 路网密度计算与可视化路网密度 道路总长度 / 区域面积。这个指标在交通规划和选址里很常用。# 计算每个行政区的路网密度 # 先投影到适合山东的投影坐标系单位用米 boundary_proj boundary.to_crs(EPSG:2333) # CGCS2000 高斯克吕格 road_proj road.to_crs(EPSG:2333) # 计算道路长度单位米 road_proj[length_m] road_proj.geometry.length # 空间连接后按行政区汇总长度 joined_proj gpd.sjoin(road_proj, boundary_proj, howleft, predicateintersects) length_sum joined_proj.groupby(NAME_left)[length_m].sum().reset_index() # 计算行政区面积单位平方公里 boundary_proj[area_km2] boundary_proj.geometry.area / 1e6 # 合并计算密度 density length_sum.merge( boundary_proj[[NAME, area_km2]], left_onNAME_left, right_onNAME, howleft ) density[density_km_per_km2] density[length_m] / 1000 / density[area_km2] print(density[[NAME_left, density_km_per_km2]].sort_values(density_km_per_km2, ascendingFalse))这里的关键是投影坐标系的选择。EPSG:2333 是 CGCS2000 下山东所在的高斯克吕格投影带单位是米算长度和面积都准确。如果你用的数据是 WGS84对应的投影带号会不同需要根据经度范围选。山东大致在东经 114° 到 122° 之间3 度带的话带号在 38 到 41 之间。注意直接用地理坐标系算长度和面积会得到以度为单位的结果没有实际意义。必须先投影。4. 避坑与排查这类 SHP 包最容易翻车的五个地方4.1 中文乱码属性表打开全是问号现象在 QGIS 或 ArcGIS 里打开 SHP属性表里的中文名称显示为乱码或问号。原因SHP 的.dbf文件编码和软件默认编码不一致。国内数据常见 GBK 编码而 QGIS 默认用 UTF-8 读取。解决在 QGIS 图层属性里找到「数据源编码」手动改成 GBK 或 GB18030。或者在加载时用-oo ENCODINGGBK参数。如果已经加载了删掉重新加载并指定编码。4.2 坐标系缺失距离量出来差了几十倍现象用测量工具量两个城市之间的距离结果和实际差了几十倍。原因.prj文件缺失或坐标系定义错误软件按默认坐标系处理导致单位混乱。解决先确认数据来源的坐标系。如果没有.prj用 QGIS 的「另存为」功能手动指定一个候选坐标系叠加在线底图验证。验证通过后再做后续分析。4.3 路网和边界对不上道路飘到海里去了现象把路网和行政区划叠加后发现部分道路明显偏移甚至飘到行政区划外面。原因两个图层的坐标系不一致或者其中一个图层本身有偏移。解决用gpd.read_file读取后检查crs属性不一致就用to_crs统一。如果统一后仍然偏移可能是数据本身的问题需要用控制点做几何纠正。4.4 空间连接结果为空sjoin 后一条都没匹配上现象执行sjoin后结果里所有行政区划字段都是 NaN。原因两个图层的坐标系不一致或者几何类型不兼容比如一个是点、一个是面。解决先检查crs是否一致不一致先to_crs。再检查几何类型boundary.geom_type和road.geom_type要匹配。如果一个是多边形一个是线intersects是可以用的但within不行。4.5 数据裁剪后属性丢失现象用行政区划裁剪路网后裁剪结果里道路名称、类型等属性全没了。原因裁剪操作默认只保留几何不保留属性或者属性字段在裁剪过程中被丢弃。解决用gpd.clip而不是手动做几何运算clip会保留属性。如果还是丢检查裁剪前后的字段列表是否一致。# 正确的裁剪方式保留属性 clipped gpd.clip(road, boundary) print(clipped.columns.tolist()) # 应该和原始 road 的字段一致5. 进阶技巧用这份数据做可达性分析和批量出图5.1 基于路网的最短路径分析有了路网数据可以做简单的可达性分析。思路是把路网转成图结构用 NetworkX 或 igraph 算最短路径。但 SHP 里的路网是线要素需要先构建拓扑关系。import networkx as nx from shapely.geometry import LineString # 构建图每条道路的起点和终点作为节点 G nx.Graph() for idx, row in road.iterrows(): line row.geometry if isinstance(line, LineString): start line.coords[0] end line.coords[-1] # 用坐标作为节点标识长度作为权重 G.add_edge(start, end, weightline.length, namerow.get(NAME, )) # 计算两个节点之间的最短路径 # 这里需要先找到离起点和终点最近的节点 def nearest_node(G, point): nodes list(G.nodes) return min(nodes, keylambda n: (n[0]-point[0])**2 (n[1]-point[1])**2) # 示例假设起点和终点坐标 start_node nearest_node(G, (117.0, 36.6)) end_node nearest_node(G, (119.0, 35.4)) path nx.shortest_path(G, start_node, end_node, weightweight) print(路径节点数:, len(path))这段代码的核心是nx.Graph()构建无向图add_edge把每条道路的起终点连起来weight用道路长度。nearest_node函数用来把任意坐标映射到最近的图节点。实际项目中路网的拓扑关系可能更复杂比如立交桥、断头路需要先做拓扑修复。5.2 批量出图每个市一张路网密度图如果需要给每个地级市单独出一张图可以用循环加matplotlib批量生成。import matplotlib.pyplot as plt # 按城市分组 cities boundary[boundary[LEVEL] 市][NAME].unique() for city in cities: fig, ax plt.subplots(figsize(10, 10)) # 画该市的行政区划 city_boundary boundary[boundary[NAME] city] city_boundary.plot(axax, facecolornone, edgecolorblack, linewidth1) # 画该市的路网 city_road gpd.clip(road, city_boundary) city_road.plot(axax, colorred, linewidth0.5) ax.set_title(f{city} 路网分布) ax.set_axis_off() plt.savefig(f{city}_road.png, dpi150, bbox_inchestight) plt.close()gpd.clip把路网裁剪到每个市的范围内plot分别画边界和路网savefig保存图片。dpi150保证出图清晰度bbox_inchestight去掉多余白边。5.3 一个我踩过的坑投影带号选错导致面积算错有一次做路网密度我随手选了一个投影坐标系结果算出来的面积比实际大了将近一倍。排查了半天才发现那个投影带的中央经线不在山东范围内导致投影变形过大。从那以后我每次做投影转换前都会先确认数据的经度范围再选对应的带号。山东的数据3 度带用 EPSG:2333 到 2336 之间6 度带用 EPSG:2327 到 2330 之间具体选哪个要看数据覆盖的经度范围。提示如果不确定选哪个带号可以用pyproj的CRS.from_dict动态构建一个以数据中心经线为中央经线的高斯克吕格投影这样变形最小。from pyproj import CRS # 以数据中心经线构建自定义投影 center_lon boundary.to_crs(EPSG:4326).geometry.centroid.x.mean() custom_crs CRS.from_dict({ proj: tmerc, lat_0: 0, lon_0: center_lon, k: 1, x_0: 500000, y_0: 0, ellps: GRS80, units: m }) boundary_custom boundary.to_crs(custom_crs) print(自定义投影面积:, boundary_custom.geometry.area.sum() / 1e6, 平方公里)这段代码用数据本身的中心经线作为中央经线避免了固定带号可能带来的变形问题。tmerc是横轴墨卡托投影lon_0设为中心经线x_0加 500000 是避免负坐标。算出来的面积和实际面积对比一下如果接近就说明投影选对了。希望帮到你。本文还有配套的精品资源点击获取

读完文章,也想定制专属网站?

尧图设计师 24 小时内与您沟通定制方案

免费获取报价 →
↑