企业官网建设流程全解析

简介：这套面向GIS从业者、城市规划与交通研究人员及地理信息专业学生的空间数据包，聚焦四川省2022年7月的行政区划与交通网络底图需求。包内共33个文件，以shp、dbf、shx、prj为主，辅以xml元数据及sbn/sbx空间索引文件，压缩包约66.64MB，可在ArcGIS、QGIS等平台直接加载。数据覆盖省、市、县三级行政边界，属性表记录行政级别与区域代码，投影信息明确坐标系统；道路网区分高速公路、主干道和次要道路，铁路网包含高速铁路与普通线路及车站等属性，便于开展行政区域对比、人口与经济差异分析、交通流量评估、选址与服务覆盖范围测算。已有2637人学习，适合用于专题制图、空间查询、交通规划模型构建及灾害响应等场景，为理解四川省地理空间格局、识别交通瓶颈和优化规划提供较完整的基础数据支撑。

1. 从一份 2022.07 省级 SHP 说起：三级区划与路网能拼出什么

一份 2022.07 的省级 SHP 资源，把省、市、县三级行政区划面层和道路网、公路网线层放在同一个包里，做 GIS 落图时最费时间的往往不是加载，而是坐标、编码和拓扑对齐。它适合做县域统计、路网密度面板、专题图底图，也适合拿来做空间连接和叠加分析的练习数据。下面用“本省”指代这份 SHP 覆盖的省域，重点拆三件事：图层怎么体检、三级区划怎么挂接、道路网和公路网怎么落到县。铁路网线层如果也在包里，处理思路相近，但等级字段和道路网不是一套，后面会单独说。新手可以跟着代码跑一遍，熟手可以直接看参数和踩坑段。

2. SHP 拆包与图层体检：坐标系、字段、几何三项先过一遍

拿到 SHP 压缩包，很多人第一反应是拖进 QGIS 或 ArcGIS 看颜色，结果一算面积就翻车。血泪经验是：先做体检，再谈出图。体检只看三项——坐标系能不能认、字段类型稳不稳、几何有没有空值和自相交。三项过了，后面的空间连接和长度统计才有意义。这个资源是 2022.07 的行政区划加路网，年份不算新，但做县级统计足够；真正要小心的是不同图层可能来自不同来源，坐标系不一定统一。

2.1 压缩包里的目录结构与图层命名

SHP 不是单个文件，而是一组同名文件。常见的是.shp存几何、.shx存索引、.dbf存属性、.prj存坐标系、.cpg存编码。少一个.shx或.dbf，GDAL 可能还能读，但属性或空间索引会出问题；少.prj就只能靠猜坐标系。拿到包后先解压到一个纯英文路径，别放在带空格和中文的深层目录里，某些旧版驱动对路径很敏感。

我一般先列目录，不急着打开属性表：

# 解压后进入 shp 目录，先看文件配对情况 find ./shp -maxdepth 2 -type f | sort # 统计每个图层的同名文件是否齐全 for f in ./shp/*.shp; do base="${f%.shp}" for ext in shp shx dbf prj cpg; do [ -f "${base}.${ext}" ] || echo "缺少: ${base}.${ext}" done done

这段脚本只做一件事：按 SHP 的伴生文件逐个检查。base="${f%.shp}"去掉后缀拿到同名基底，[ -f ... ]判断文件是否存在。缺.prj不一定会报错，但后面to_crs会直接失败；缺.cpg则中文属性可能变成乱码。文件名以实际包为准，代码里用./shp只是示例目录。

图层命名通常能看出层级，比如省、市、县三级面层，道路网和公路网线层，可能还有铁路网。不要靠文件名猜字段，下一步用ogrinfo一次性看全。

2.2 用 GDAL/ogr 做一次图层体检

ogrinfo是 GDAL 自带的命令行工具，不依赖桌面软件，适合批量巡检。先看单个图层的摘要，再决定要不要读全量：

# -so 只输出摘要，-al 输出所有图层；单个文件也可以直接指定 ogrinfo -so -al ./shp/admin_county.shp # 批量看每个 SHP 的几何类型、要素数、坐标系和字段 for f in ./shp/*.shp; do echo "===== $f =====" ogrinfo -so -al "$f" | sed -n '1,80p' done

-so是 summary only，不读几何，速度快；-al表示所有图层。输出里重点看三行：Geometry:是 Polygon 还是 MultiPolygon，Feature Count:是不是 0，Layer SRS WKT:有没有坐标系。如果Layer SRS WKT是空，.prj大概率缺失或损坏。sed -n '1,80p'只是截断输出，避免刷屏，实际排查可以去掉。

字段部分要看类型。行政区划的编码字段如果是字符串，前导零能保住；如果变成整数，01会变成1，后面按编码挂接就会对不上。道路网的等级字段如果是数字，先确认它是不是编码值，别直接当排序用。长度字段如果原本有，只能参考，重新投影后要自己算。

2.3 坐标系判定与统一到投影坐标系的选择

坐标系是这类资源最容易埋雷的地方。行政区划面层常见的是地理坐标系，单位是度；道路网线层可能也是地理坐标系。用地理坐标系直接算面积，得到的是“平方度”，数值看起来像玄学，和平方公里差好几个数量级。正确做法是先确认 CRS，再统一转到投影坐标系。

用 GeoPandas 读前几行，看crs和字段：

import geopandas as gpd from pathlib import Path for p in sorted(Path("shp").glob("*.shp")): gdf = gpd.read_file(p, rows=5) print("图层:", p.name) print(" 几何类型:", gdf.geom_type.unique().tolist()) print(" 坐标系:", gdf.crs) print(" 字段:", list(gdf.columns))

rows=5只读前 5 行，用来快速看结构，不用于统计。gdf.crs如果是None，说明没有.prj或驱动没识别；这时不要直接to_crs，先找原始坐标系说明。geom_type.unique()能看出面层是不是混了Polygon和MultiPolygon，线层是不是混了LineString和MultiLineString。混着不是错，但后续line_merge、overlay要按类型处理。

投影怎么选？本省经度跨度不小，常见做法是按 3 度带高斯-克吕格投影，或者用 UTM 分带。做县级面积和路网长度统计，选一个覆盖研究区的等角投影即可，不必追求和官方图件完全一致。用 UTM 时按经度选带：

# 示例：东经 102°–108° 范围用 UTM 48N；更西的区域改用 47N crs_metric = "EPSG:32648" county = gpd.read_file("shp/admin_county.shp").to_crs(crs_metric) county["area_km2"] = county.geometry.area / 1e6 print(county[["county_code", "area_km2"]].head())

EPSG:32648是 WGS84 / UTM 48N，单位米；geometry.area得到平方米，除以1e6得平方公里。如果研究区跨两个带，不要硬用一个带，按县拆分或改用 Albers 等面积投影。参数不是固定的，关键是“先投影、再算长度面积”，这条顺序不能反。

提示：.prj缺失时，不要用“看起来像”的坐标系硬套。先用控制点或已知边界核对，再写入 CRS；否则后面所有面积和长度都会带着系统误差。

3. 省市县三级挂接：编码、面积与路网归属的落地做法

三级行政区划面层最常见的用法是“县挂市、市挂省”，然后按县汇总路网。听起来简单，实际有两个难点：一是层级字段不一定规整，二是空间关系不一定干净。2022.07 的资源做区划统计够用，但边界简化、飞地、跨县线状地物都会影响挂接结果。我的做法是先把属性挂接和空间挂接分开做，再互相校验，别一上来就用sjoin一把梭。

3.1 三级边界的层级字段与编码规则

省、市、县三层各自有名称和编码字段。常见编码是 6 位县级、4 位市级、2 位省级，但字段名可能是code、adcode、xzqhdm之类，长度也不一定统一。先看字段值和长度分布，不要假设编码都是字符串：

import geopandas as gpd for name in ["admin_province", "admin_city", "admin_county"]: gdf = gpd.read_file(f"shp/{name}.shp") print("=====", name, len(gdf)) for col in gdf.columns: if col == "geometry": continue sample = gdf[col].dropna().astype(str).head(3).tolist() print(f" {col}: dtype={gdf[col].dtype}, sample={sample}")

重点看编码字段的dtype和样本值。如果dtype是int64，前导零已经丢了，需要从名称或上级编码重新补；如果样本里出现"000000"、空字符串、"暂无"，说明有占位值，挂接前要过滤。名称字段可能有同名乡镇或同名区县，不能拿名称当唯一键，编码才是主键。

层级关系不要靠名称包含来判断，比如“某市”下面不一定所有县都带市名。常见做法是给县层补一个city_code，给市层补一个province_code，然后用编码做左连接。如果资源本身没有上级编码字段，就用空间关系反推，但反推结果要人工抽查。

3.2 用 Python + GeoPandas 做层级挂接与面积核算

先做属性挂接，再做面积核算。属性挂接用merge，空间挂接用sjoin，两者结果对不上时以编码为准，空间结果只做校验。下面是一段可抄的流程：

import geopandas as gpd import pandas as pd crs_metric = "EPSG:32648" prov = gpd.read_file("shp/admin_province.shp").to_crs(crs_metric) city = gpd.read_file("shp/admin_city.shp").to_crs(crs_metric) county = gpd.read_file("shp/admin_county.shp").to_crs(crs_metric) # 1) 面积：投影后算，单位平方公里 county["area_km2"] = county.geometry.area / 1e6 # 2) 县挂市：优先用属性编码；没有上级编码时用空间 within if "city_code" in county.columns: county_city = county.merge( city[["city_code", "city_name"]], on="city_code", how="left", suffixes=("", "_city") ) else: county_city = gpd.sjoin( county, city[["city_code", "city_name", "geometry"]], how="left", predicate="within" ) # 3) 检查一对多：一个县编码出现多条记录就是异常 dup = county_city.groupby("county_code").size() print("一对多县数:", (dup > 1).sum()) print(dup[dup > 1].head())

to_crs必须先做，geometry.area才是平方米。merge的how="left"保留所有县，挂不上的city_name为空，方便后面查。sjoin的predicate="within"要求县几何完全落在市几何内，边界重叠或简化不一致时会漏；这时可以改用intersects再按面积占比取最大，但不要默认intersects就是对的。groupby("county_code").size()是检查一对多的低成本办法，结果大于 1 就要看是编码重复还是空间重叠。

面积核算还要注意 MultiPolygon。geometry.area会把多部件面积相加，一般没问题；但如果面层里有洞，面积会自动扣掉。要核对总面积，可以和省级面层面积对比，差 1% 以内通常能接受，差太多就是投影或几何问题。

3.3 面层与线层的叠加：路网落在哪个县

路网落县有两个口径：一个是“道路经过哪些县”，用intersects，跨县道路会出现在多个县；另一个是“每个县境内的道路长度”，必须先裁剪再汇总。做密度面板要用后者，做连通性分析可以用前者。下面先给裁剪口径：

road = gpd.read_file("shp/road_network.shp").to_crs(crs_metric) county_sel = county[["county_code", "county_name", "geometry"]].copy() # overlay 会把跨县道路切碎，每个碎片只属于一个县 road_clip = gpd.overlay(road, county_sel, how="intersection") road_clip["len_km"] = road_clip.geometry.length / 1000 summary = ( road_clip.groupby(["county_code", "county_name"], as_index=False)["len_km"] .sum() .sort_values("len_km", ascending=False) ) print(summary.head(10))

overlay(how="intersection")的输出保留两个图层的属性，几何是相交部分；跨县道路会被切成多段，每段只落一个县，这样汇总不会重复。geometry.length在投影坐标系下单位是米，除以 1000 得公里。groupby用县编码和县名一起分组，避免同名县混淆。overlay对大数据慢，跑之前先裁研究区，或者给县层建空间索引。

如果包里有公路网和道路网两套线层，不要混在一起直接加总。道路网可能包含城市道路、乡村道路，公路网偏干线，等级字段和重复路段都不一样。铁路网线层如果也在，处理方式相同，但铁路的等级、复线、电气化字段和公路不是一套，统计口径要单独说明。常见做法是给每个线层加一个net_type字段，再合并汇总，避免把铁路里程算进公路密度。

注意：sjoin的intersects会把跨县道路复制到多个县，直接sum会重复计算里程。做县级里程统计时，用overlay裁剪，或者先按道路 ID 去重再汇总。

4. 道路网与公路网落图：线层清洗、分级符号化与里程汇总

路网线层比行政区划面层更容易出问题。面层最多是边界简化，线层会有悬挂点、重复线段、断头路、等级字段缺失。把道路网和公路网直接叠在区划底图上，视觉上很漂亮，一做里程统计就露馅。这个资源的价值在于“区划加路网”能直接支撑县级统计，但线层必须先清洗再符号化，否则图面好看、数字对不上。

4.1 线层几何清洗：打断、去重、悬挂点

线层清洗的目标不是把路网修成完美拓扑，而是让长度统计可重复。第一步去空几何和重复几何，第二步按道路等级或道路 ID 合并连续线段，第三步检查悬挂点。下面是一段基础清洗：

import geopandas as gpd from shapely.geometry import LineString, MultiLineString road = gpd.read_file("shp/road_network.shp") # 去空几何、去重复几何 road = road[road.geometry.notna() & ~road.geometry.is_empty].copy() road = road.drop_duplicates(subset=["geometry"]).copy() # 按道路等级合并，再尝试把多段线接起来 road["grade"] = road["grade"].fillna("unknown").astype(str) merged = road.dissolve(by="grade", as_index=False).explode(index_parts=False) def merge_line(g): if isinstance(g, MultiLineString): return g.line_merge() if isinstance(g, LineString): return g return g merged["geometry"] = merged.geometry.apply(merge_line) print(merged[["grade"]].assign(n=1).groupby("grade").count())

drop_duplicates(subset=["geometry"])只按几何去重，属性不同的重复线也会被删掉；如果重复线承载不同属性，改成按道路 ID 加几何去重。dissolve(by="grade")把同等级线合并成多段线，explode拆回单段，line_merge尝试把端点相接的线段接成一条。line_merge不会跨不相接的线段工作，所以悬挂点还在，下一步单独查。

悬挂点检查可以用端点计数：一条线有两个端点，如果某个端点只出现一次，大概率是断头路或边界裁断。下面只做统计，不强行修复：

from collections import Counter from shapely.geometry import Point endpoints = [] for geom in road.geometry: if geom.geom_type == "LineString": coords = list(geom.coords) endpoints.extend([Point(coords[0]), Point(coords[-1])]) elif geom.geom_type == "MultiLineString": for part in geom.geoms: coords = list(part.coords) endpoints.extend([Point(coords[0]), Point(coords[-1])]) counter = Counter((round(p.x, 6), round(p.y, 6)) for p in endpoints) hanging = [k for k, v in counter.items() if v == 1] print("端点总数:", len(endpoints), "只出现一次的端点:", len(hanging))

坐标四舍五入到 6 位是為了避免浮点误差导致同一点被当成两个点；这个精度在投影坐标系下约 0.1 米，够用。只出现一次的端点不一定是错误，行政边界裁断、数据分层都会造成断头。要做连通性分析再修，做里程统计可以先不修。

4.2 道路等级字段映射与符号化参数

道路等级字段常见的是数字编码或中文名称，直接拿来符号化会乱。先做一张映射表，把原始值归一到“高速、国道、省道、县道、乡道、城市道路、其他”几类。字段值以实际数据为准，下面只是示例：

映射不要写死在符号化里，先在属性表加一列grade_norm：

def normalize_grade(v): s = str(v).strip().upper() if s in {"高速", "G", "GS"}: return "高速" if s.startswith("G") or s in {"国道"}: return "国道" if s.startswith("S") or s in {"省道"}: return "省道" if s.startswith("X") or s in {"县道"}: return "县道" if s.startswith("Y") or s in {"乡道"}: return "乡道" if s in {"城市道路", "城区道路"}: return "城市道路" return "其他" road["grade_norm"] = road["grade"].apply(normalize_grade)

这个函数只是示例，实际字段可能是"G108"、"国道"、"高速"混用，需要按样本调整。str(v).strip().upper()先统一大小写和空格，避免" g108"漏掉。归一化后按等级设线宽，总图才能看出骨架；不要把所有线设成同一宽度，那样道路网和公路网会糊成一片。

4.3 按县汇总里程：投影选择与长度口径

里程汇总最怕两个问题：坐标系没投影、跨县重复计算。前面已经讲了投影和裁剪，这里补长度口径。道路长度有“几何长度”和“运营里程”两种，SHP 里能算的是几何长度；如果原始属性里有里程字段，只能做参考，不能直接替代。按县汇总时，先裁剪再算长度，顺序不能反。

import geopandas as gpd crs_metric = "EPSG:32648" county = gpd.read_file("shp/admin_county.shp").to_crs(crs_metric) road = gpd.read_file("shp/road_network.shp").to_crs(crs_metric) # 先裁到县界，再算长度 road_clip = gpd.overlay(road, county[["county_code", "county_name", "geometry"]], how="intersection") road_clip["len_km"] = road_clip.geometry.length / 1000 # 按县和等级汇总，方便后面做分级设色 by_county_grade = ( road_clip.groupby(["county_code", "county_name", "grade_norm"], as_index=False)["len_km"] .sum() ) print(by_county_grade.head())

overlay之后每条线只属于一个县，length不会重复。groupby同时保留县和等级，后面可以画堆叠图或做等级占比。如果只想要县总里程，去掉grade_norm再sum一次。crs_metric换成其他投影时，长度会有细微差异，做跨年对比要固定同一个投影，不要每年换带。

公路网和道路网如果分属两个线层，建议分别汇总再合并，不要先concat再裁剪。因为两个线层的重复路段、等级字段、几何精度可能不同，先合并会把重复里程带进来。常见做法是给每个线层加net_type，各自裁剪汇总，最后按县编码merge，并在报告里注明口径。

提示：长度统计前先确认线层有没有 Z 值或 M 值。geometry.length只看 XY，Z 值不影响；但如果几何是 3D 且投影不对，长度仍按平面算，不要误以为带了高程就是真实里程。

5. 避坑与排查：坐标系、编码、空几何这五处最容易翻车

这份资源用起来顺不顺，很大程度取决于前 30 分钟有没有做检查。下面五条是我在区划加路网数据上踩过的坑，每条按“现象、原因、解决”写，能复现的尽量给命令或代码。不是每个包都会全中，但遇到了能快速定位。

5.1 现象：面积算出来差了一个数量级

现象很直观：同一个县，QGIS 里量出来几千平方公里，Python 里算出来几万甚至几十万。原因通常只有一个——面层还是地理坐标系，geometry.area的单位是“平方度”。度不是长度单位，平方度更不是面积单位，数值没有可比性。另一种情况是投影带选错，本省西部和东部硬用同一个 UTM 带，边缘变形大，面积偏差也会上去。

解决分两步：先看gdf.crs是不是地理坐标系，is_geographic为True就说明还没投影；再用to_crs转到覆盖研究区的投影坐标系，然后算面积。跨带时不要一个 EPSG 用到底，按县判断经度范围，或者改用 Albers 等面积投影。算完拿省级面层总面积做一次对数，差 1% 以内算正常，差 10% 以上先查投影和几何。

import geopandas as gpd gdf = gpd.read_file("shp/admin_county.shp") print("是否地理坐标系:", gdf.crs.is_geographic if gdf.crs else "无 CRS") gdf_m = gdf.to_crs("EPSG:32648") gdf_m["area_km2"] = gdf_m.geometry.area / 1e6 print(gdf_m["area_km2"].describe())

describe()能快速看面积分布，出现 0 或异常大值就要查空几何和自相交。to_crs前如果crs是None，先补.prj或手动set_crs，不要直接转。

5.2 现象：市县挂接出现一对多

现象是sjoin之后一个县编码对应多条市记录，或者merge之后行数变多。原因有三种：县层编码重复、市层边界重叠、县界跨市界。边界简化数据里，飞地或跨界片区会让within失败，改用intersects又会把边界线附近的县挂到两个市。名称挂接更危险，同名区县直接撞车。

解决顺序是先用编码merge，再用空间关系校验。编码挂接后检查county_code是否唯一；空间校验时用within，对挂不上的县单独看几何，不要直接改成intersects了事。确实需要intersects的，按相交面积占比取最大，并记录占比，低于阈值的标记为待人工确认。

# 编码挂接后检查重复 dup = county_city.groupby("county_code").size() print("重复县编码:", dup[dup > 1].index.tolist()) # 空间挂接后检查一个县是否落到多个市 sj = gpd.sjoin(county, city[["city_code", "geometry"]], how="left", predicate="within") cnt = sj.groupby("county_code").size() print("空间一对多:", cnt[cnt > 1].head())

sjoin输出的索引是原表索引，groupby("county_code")前先确认字段还在。边界重叠严重时，先做make_valid修几何，再挂接。

5.3 现象：路网在县界处断成两截

现象是地图上路网在县界处突然断开，或者按县汇总后总里程比全省总里程多一截。原因不是路真的断了，而是线层按县界裁过，或者原始数据本身在边界处被打断。line_merge只能接端点完全重合的线，坐标差一点就接不上；overlay裁剪后每条线只保留县内部分，跨县道路自然分成两段。

解决要看目的：做连通性分析，先按道路 ID 或名称合并，再做line_merge，合并前把端点坐标四舍五入到统一精度；做县级里程统计，接受裁剪后的分段，但汇总时按县编码分组，不要按道路 ID 去重。检查悬挂点时，端点只出现一次的数量可以作为参考，不要追求零悬挂，行政边界数据很难做到。

# 按道路 ID 合并同一条路的多段线，再尝试接续 road["road_id"] = road["road_id"].fillna(road.index.astype(str)) dissolved = road.dissolve(by="road_id", as_index=False).explode(index_parts=False) dissolved["geometry"] = dissolved.geometry.apply( lambda g: g.line_merge() if g.geom_type == "MultiLineString" else g )

dissolve按road_id合并，explode拆回单段，line_merge再接。没有road_id时只能用名称加等级近似，误差会大一些。

5.4 现象：中文属性乱码

现象是属性表里县名、道路名变成问号或乱码，ogrinfo输出也不正常。原因是 DBF 的编码和读取端不一致，常见是 GBK 与 UTF-8 混用，.cpg缺失或写错。GDAL 有时能自动猜，猜错就乱码；GeoPandas 读进来直接是乱码字符串，后面按名称匹配全废。

解决先看.cpg内容，再用encoding参数显式读。GDAL 的SHAPE_ENCODING环境变量也能起作用。读进来后不要直接改.dbf，先复制一份，确认编码再写回。中文名称只用于显示，挂接尽量用编码，能避开大部分乱码问题。

# 先试 UTF-8，乱码再试 GBK；以实际 .cpg 为准 gdf = gpd.read_file("shp/admin_county.shp", encoding="GBK") print(gdf["county_name"].head())

encoding="GBK"不是万能，UTF-8 数据用 GBK 读也会乱。先看.cpg里写的是UTF-8还是GBK，再用对应编码。写回时用to_file(..., encoding="UTF-8")，并确认目标软件支持。

5.5 现象：导出 GeoJSON 后字段丢失

现象是 SHP 导出 GeoJSON 后，某些字段没了或字段名被截断。原因是 DBF 字段名有长度限制，GeoJSON 虽然支持长字段名，但中间转换可能被截断；另外中文字段名、特殊字符、空值类型也可能导致驱动写入失败。导出 GeoJSON 后坐标顺序变成[x, y]，如果下游按[lat, lon]读，位置会反。

解决是导出前把字段名改成英文短名，保留一个映射表；能用 GPKG 就用 GPKG，字段名和类型更稳。必须 GeoJSON 时，显式指定driver="GeoJSON"，导出后抽查字段数量和坐标顺序。

# 导出前重命名字段，保留映射 rename_map = {"county_code": "cty_code", "county_name": "cty_name"} out = county.rename(columns=rename_map) out.to_file("county.geojson", driver="GeoJSON", encoding="UTF-8")

rename_map按实际字段调整，导出后ogrinfo -so -al county.geojson对一遍字段数。坐标顺序在 GeoJSON 规范里是[x, y]，即[经度, 纬度]，下游库如果默认[lat, lon]要显式设置。

6. 进阶用法：把 SHP 做成县级路网密度面板

基础挂接跑通后，最有复用价值的产物不是一张图，而是一张县级面板表：每个县一行，包含面积、道路里程、公路里程、密度。后面换年份、换路网图层，只要字段对齐就能重复跑。下面这段代码把区划、道路网、公路网合成一张 CSV，字段名按实际数据改：

import geopandas as gpd import pandas as pd crs_metric = "EPSG:32648" county = gpd.read_file("shp/admin_county.shp").to_crs(crs_metric) county["area_km2"] = county.geometry.area / 1e6 def road_stats(path, prefix): road = gpd.read_file(path).to_crs(crs_metric) clip = gpd.overlay(road, county[["county_code", "geometry"]], how="intersection") clip[f"{prefix}_km"] = clip.geometry.length / 1000 return clip.groupby("county_code", as_index=False)[f"{prefix}_km"].sum() road_df = road_stats("shp/road_network.shp", "road") highway_df = road_stats("shp/highway_network.shp", "highway") panel = county[["county_code", "county_name", "area_km2"]].merge(road_df, on="county_code", how="left") panel = panel.merge(highway_df, on="county_code", how="left").fillna({"road_km": 0, "highway_km": 0}) panel["road_density"] = panel["road_km"] / panel["area_km2"] panel["highway_density"] = panel["highway_km"] / panel["area_km2"] panel.to_csv("county_road_panel.csv", index=False, encoding="utf-8-sig")

road_stats把“裁剪、算长度、按县汇总”封成函数，换图层只改路径和前缀。fillna(0)处理没有路网经过的县，避免密度出现空值。utf-8-sig让 CSV 在常见表格软件里打开不乱码。面板字段建议固定为：县编码、县名、面积、道路里程、公路里程、道路密度、公路密度，后面做分级设色或时序对比都够用。

这张表跑出来后，我一般先做两件事：一是按密度排序看尾部，密度为 0 的县是漏裁还是真没有路；二是把总面积加总，和省级面层对一遍，差太多就回查投影。铁路网如果也要进面板，建议单独加railway_km，不要和公路混在一个密度字段里，口径不同，混在一起后面解释不清。

从那以后我每次拿到新的 SHP，都先跑一遍坐标系、编码、空几何三项检查，再谈出图和统计；这套顺序帮我省掉了很多返工。希望帮到你。

本文还有配套的精品资源，点击获取