PostGIS简化英国行政区划边界出现几何损坏、缝隙重叠如何解决

行政边界无重叠无碎片简化方案 #

现有方案问题根因

你当前使用的SQL在处理时存在两个核心缺陷：

1. 仅提取外边界做合并简化，会丢失多边形内部的洞边界，同时st_linemerge遇到未完全衔接的线串时会生成断链，后续st_polygonize会生成无效碎面
2. 0.5的面积匹配阈值在小面积行政区较多的场景下容易出现属性错配，进一步加剧几何和属性不对应的问题

方案1：无需拓扑扩展的通用简化方案

适配英国行政区（坐标系推荐用EPSG:27700英国国家格网）场景的优化SQL如下：

WITH 
-- 提取所有多边形的内外边界，建议给原始表加gid主键，比name匹配更准确
all_rings AS (
    SELECT 
        p.gid,
        p.name,
        st_exteriorRing((st_dumpRings(p.geom)).geom) AS ring_geom
    FROM new_bdline_county p
),
-- 合并所有边界后统一简化，保证共享边简化结果完全一致
simplified_edges AS (
    SELECT 
        (st_dump(st_simplifyPreserveTopology(st_union(ring_geom), 100))).geom AS edge_geom
    FROM all_rings
),
-- 用简化后的边重新构建面
polygonized AS (
    SELECT (st_dump(st_polygonize(edge_geom))).geom AS geom
    FROM simplified_edges
)
-- 精准匹配原始面和新生成面的属性
SELECT 
    a.name,
    p.geom
FROM polygonized p
JOIN new_bdline_county a 
    ON st_contains(a.geom, st_pointonsurface(p.geom)) -- *用面内点匹配比面积比更准确，不会错配小面*
WHERE st_isValid(p.geom);

参数调整说明：

• 简化容差：如果用EPSG:27700坐标系，单位为米，100代表简化后边界偏移不超过100米，可根据你的精度需求调整
• 若原始行政区存在重名，务必用主键gid做匹配，避免属性错配

方案2：PostGIS拓扑扩展正确用法

之前用拓扑扩展出问题大概率是没有设置拓扑容差，导致共享边没有被正确合并，正确步骤如下：

-- 1. 启用拓扑扩展
CREATE EXTENSION IF NOT EXISTS postgis_topology;
-- 2. 创建拓扑，参数为拓扑名、坐标系SRID、拓扑容差（单位同坐标系，这里设0.1米用于合并精度误差导致的共享边缝隙）
SELECT CreateTopology('uk_admin_topo', 27700, 0.1);
-- 3. 给原始表加拓扑列
SELECT AddTopoGeometryColumn('uk_admin_topo', 'public', 'new_bdline_county', 'topo_geom', 'MULTIPOLYGON');
-- 4. 把原始几何转换成拓扑几何，自动拆分共享边
UPDATE new_bdline_county SET topo_geom = toTopoGeom(geom, 'uk_admin_topo', 1);
-- 5. 统一简化所有拓扑边
UPDATE uk_admin_topo.edge_data SET geom = st_simplifyPreserveTopology(geom, 100);
-- 6. 从简化后的拓扑生成最终简化面
CREATE TABLE simplified_uk_admin AS
SELECT name, topo_geom::geometry AS geom FROM new_bdline_county;

后置校验步骤

简化完成后执行以下检查，确认无问题：

• 检查无效几何：SELECT * FROM simplified_uk_admin WHERE NOT st_isValid(geom);
• 检查重叠面：SELECT a.name, b.name FROM simplified_uk_admin a JOIN simplified_uk_admin b ON a.gid < b.gid AND st_overlaps(a.geom, b.geom);
• 检查面积损失：SELECT name, st_area(geom)/st_area(original_geom) AS area_ratio FROM simplified_uk_admin WHERE st_area(geom)/st_area(original_geom) < 0.9;

内容的提问来源于stack exchange，提问作者gcj