全国乡镇街道边界数据下载前先查一遍 code:重复要素与环方向错乱的清理顺序
福建宁德福鼎市的乡镇级文件打开是 23 个要素,把 code 抽出来排一遍,只剩 17 个。六个乡镇各出现两次,其中沙埕镇的两份几何一件 123 个顶点、一件 2344 个,包围盒西边界相差 0.042 度。做全国省市区县区划边界数据下载的人,最容易在这一步翻车:文件能解析、能画出来、看着还是一张完整的县域图,可一旦拿去按乡镇做统计,人口和面积都会翻倍。
这类问题不挑地区。把站内 2901 个区县文件全部扫一遍,189 个文件存在重复 code,涉及 880 个代码、914 个冗余要素。它们不会让任何软件报错。
重复 code 的两种长相,处理方式完全不同
扫出来的重复不是同一种毛病,分清楚才谈得上修。
一种是同一实体两个版本。福鼎市这六个就属于此类:一版来自较早的乡镇级数据,几何粗、部件少;另一版来自较新的行政区划更新,几何细、部件多。两版的 NAME、code、PAC 完全一致,只有几何和 FID 不同。这种是纯粹的冗余,新版即正确版,可安全折叠。
另一种是编码本身撞车。宿州市萧县 35 个要素对应 24 个唯一代码;信阳市固始县 42 个要素对应 33 个。后者更多是历史代码残缺——原代码缺失或重复时被统一补了同一条占位值,覆盖了名称各异的多个乡级单元。这种不能靠"留一个"解决,得回到区划变更时间线逐个归位。
判断捷径是比对 FID。同一批次生成的记录 FID 连号,跨批次合并进来的往往差得很远:福鼎这六对分别是 2644/45067、2645/45068、2646/45072,前一个连号,后一个跳号。
重复要素有多值钱:福鼎的面积被抬高了 59%
光看要素个数不容易有体感,换成面积就直观了。
| 口径 | 要素数 | 球面面积合计 | 偏差 |
|---|---|---|---|
| 直接全量求和 | 23 | 1514.08 km² | 基准 +59.4% |
| 按 code 去重后 | 17 | 949.60 km² | 真实值 |
重复的六份几何不是完全重合的副本,而是互相错位的两个版本。管阳镇两版包围盒纬度下界差 0.098 度,沙埕镇差 0.116 度。把它们一起求和,等于把同一片地算了两次,还把两版之间的错位区域也算了进去。1514 与 950 之间那 564 平方公里,全是不存在的面积。
县级汇总通常不受影响,因为多数平台按市、区县两级取数。一旦把全国省市区县geojson数据下载下来落到乡镇做人口密度、设施覆盖或者报表拆分,这个 59% 会原样传导到结果里,而且找不到出处。
环方向比重复 code 更隐蔽
第二个坑不看面积看不出来:同一个文件里,外环的绕向可能一半顺时针、一半逆时针。
2901 个区县文件扫下来,外环全部一致的只有 2308 个——916 个全逆时针,1392 个全顺时针;剩下 593 个是混合的。数量本身不惊人,麻烦在于方向约定不统一。同一批数据里既有逆时针外环,也有顺时针外环,说明它不是单次导出生成,而是多次更新的结果叠加。
福鼎市 104 个环全部逆时针,霞浦县 24 个环也全部逆时针,福州市鼓楼区 10 个环却全部顺时针,漳州市东山县 31 个环同样全部顺时针。四份文件分别看都"自洽",放在一起就矛盾。
GeoJSON 规范里,外环逆时针、内环顺时针是通行约定。方向反过来,渲染引擎多数照画不误,可做差集、求交、溶解的时候就会把实心区当成空心区,把挖洞当成实体。全国省市区县区划边界数据下载拿回的多边形一旦进入空间运算环节,这类反绕向会直接产出错误结果。
顺带一个附带损耗:鼓楼区 10 个环里有 26 处相邻重复顶点,东山县 31 个环里有 87 处——同一个坐标连着写了两次。不致命,但会让所有基于线段数的统计偏高。
清理顺序:先折叠、再修向、最后才抽稀
顺序错了,返工成本很高。可行的一次性流水线是这样:
import json
from collections import defaultdict
def signed_area(ring):
"""鞋带公式(经纬度平面近似),正值 = 逆时针"""
s = 0.0
for i in range(len(ring) - 1):
(x1, y1), (x2, y2) = ring[i], ring[i + 1]
s += (x2 - x1) * (y2 + y1)
return s / 2
def clean(path):
fc = json.load(open(path, encoding="utf-8"))
# 1) 折叠重复 code:保留顶点更多的那份,几何更细
best = {}
for f in fc["features"]:
code = f["properties"].get("code")
if not code:
continue
n = sum(len(r) for poly in (
[f["geometry"]["coordinates"]] if f["geometry"]["type"] == "Polygon"
else f["geometry"]["coordinates"]) for r in poly)
if code not in best or n > best[code][0]:
best[code] = (n, f)
# 2) 逐环修正绕向:外环逆时针,内环顺时针
for _, f in best.values():
g = f["geometry"]
polys = [g["coordinates"]] if g["type"] == "Polygon" else g["coordinates"]
for poly in polys:
for idx, ring in enumerate(poly):
ccw = signed_area(ring) > 0
if (idx == 0 and not ccw) or (idx > 0 and ccw):
poly[idx] = ring[::-1]
return {"type": "FeatureCollection", "features": [f for _, f in best.values()]}
两个关键点:方向修正必须在折叠之后做,因为折叠时按顶点数择优,不同来源的方向可能不一致,先修会被无效版本覆盖。抽稀放在最后,一是抽稀会改变顶点数,影响不了择优判断了;二是抽稀之后环的走向有时会翻转,早修等于白修。
折叠这一步还有一处要留意:不能简单取第一条。福鼎的例子说明,先出现的那份往往是粗版。按顶点数取大值,能在全量数据里稳定拿到细版,不需要额外配置。
校核脚本跑多快,决定它会不会被用
这套检查的全部成本在解析。福鼎、霞浦、蕉城、福安四份合起来 2.1 MB,逐环算面积加去重判定耗时 220 毫秒左右;2901 个文件全扫约 90 秒,单机单线程。写进构建脚本挂在每次数据更新之后跑,成本可以忽略。
判断放行与否建议只看三条:重复 code 数为 0;外环绕向在本文件内一致;相邻重复顶点数低于总顶点数的 1%。三条达标再去抽稀、切瓦片、做空间索引,下游的返工基本可以避免。
这也是 GeoJSON下载 之后值得固定下来的一道工序。它不产出任何新数据,只是把"这份文件到底包含多少个真实单元"这件事确认清楚。确认过了,后面每一步的结论才有意义。
顺带一提,乡镇级数据本身还有一个绕不开的现实:名称里带"办事处"的全国有 838 个,"农场"402 个,"林场"234 个,"开发区"206 个,"监狱"27 个,"原种场"26 个。这些单元在区划代码体系里的层级归属并不统一,有些挂着正式乡镇码,有些用的是园区自编码。做全国村级geojson数据下载或者往下挂接村级数据时,先按 code 的 9 位结构判断层级,比按名称后缀猜要可靠。
如果手上这份是要拿去做汇报底图,想让某个市、某个区县单独成图,也可以去 ppt.html 页面按省市逐级下钻,挑一张可编辑的地图 PPT 模板,先预览确认版式再决定是否下载,省得自己从边界文件开始描。至于 GeoJSON 文件本身,本站的省、市、区县、乡镇各级都能先看预览图核对轮廓,确认无误再按次导出——预览不收费,导出与 PPT 文件是付费项目。