PostGIS多边形的递归交集问题

13

我正在尝试在一个空间表中的所有多边形之间执行递归交集,并获取每个交集的(多)多边形以及有关它们的每个交集的信息。

下面是一张(非比例)图片来解释它: Example

假设表格中有A,B,C三个正方形。我希望输出A,B,C,A+B,A+C,B+C,A+B+C七个多边形,并且需要知道A+BAB的交集等等。

目前,我已经有了一个查询可以执行交集,但它没有“剪掉”原始多边形的相交部分。例如:

Polygon A should be      A - (A+B) - (A+C) - (A+B+C)
Polygon A+C should be    A+C - (A+B+C)

这是我现在得到的 AA+C 多边形结果的图片:

当前错误结果

以下是一个测试脚本,使用图片中的正方形作为数据。从 area 列可以看出,缺少某些递归的 ST_Difference,但我无法想出如何实现。欢迎任何想法。

-- Create a test table
CREATE TABLE test (
    name text PRIMARY KEY,
    geom geometry(POLYGON)
);

-- Insert test data
INSERT INTO test (name, geom) VALUES 
    ('A', ST_GeomFromText('POLYGON((1 2, 1 6, 5 6, 5 2, 1 2))')),
    ('B', ST_GeomFromText('POLYGON((0 0, 0 4, 4 4, 4 0, 0 0))')),
    ('C', ST_GeomFromText('POLYGON((2 0, 2 4, 6 4, 6 0, 2 0))'));


-- Query    
WITH RECURSIVE 
source (rownum, geom, ret) AS (
    SELECT row_number() OVER (ORDER BY name ASC), ST_Multi(geom), ARRAY[name] FROM test 
),
r (rownum, geom, ret, incroci) AS (
    SELECT rownum, geom, ret, 0 FROM source 
    UNION ALL
    SELECT s.rownum, ST_CollectionExtract(ST_Intersection(s.geom, r.geom), 3), (r.ret || s.ret), (r.incroci + 1) 
        FROM source AS s INNER JOIN r ON s.rownum > r.rownum AND ST_Intersects(s.geom, r.geom) AND ST_Area(ST_Intersection(s.geom, r.geom)) > 0.5
),
result (geom, ret) AS (
    SELECT ST_Union(geom) AS geom, ret FROM r GROUP BY ret
)
SELECT geom, ST_Area(geom) AS area, ret FROM result ORDER BY ret

当然,在这个特定的例子中,窗口函数并不是必需的,但这段代码是我真实情况的简化版本,它还在一些其他方面做了更多的事情。

我正在使用PostgreSQL 9.2和PostGIS 2.0。

1个回答

8

ST_DIFFRENCE不必递归,因为您已经拥有所有的多边形,所以从每个geom减去包含该ret但不等于该ret的其他geom的并集。这样可以实现,您应该按照以下方式进行操作:

    WITH RECURSIVE 
source (rownum, geom, ret) AS (
    SELECT row_number() OVER (ORDER BY name ASC), ST_Multi(geom), ARRAY[name] FROM test 
),
r (rownum, geom, ret, incroci) AS (
    SELECT rownum, geom, ret, 0 FROM source 
    UNION ALL
    SELECT s.rownum, ST_CollectionExtract(ST_Intersection(s.geom, r.geom), 3), (r.ret || s.ret), (r.incroci + 1) 
        FROM source AS s INNER JOIN r ON s.rownum > r.rownum AND ST_Intersects(s.geom, r.geom) AND ST_Area(ST_Intersection(s.geom, r.geom)) > 0.5
),
result (geom, ret) AS (
    SELECT ST_Difference(ST_Union(r.geom),q.geom) AS geom, r.ret FROM r JOIN (SELECT r.ret,ST_UNION(COALESCE(r2.geom,ST_GeomFromText('POLYGON EMPTY'))) as geom FROM r LEFT JOIN r AS r2 ON r.ret<@r2.ret AND r.ret!=r2.ret GROUP BY r.ret) AS q on r.ret=q.ret GROUP BY r.ret,q.geom
)
SELECT geom, ST_Area(geom) AS area, ret FROM result ORDER BY ret

1
谢谢,使用自连接确实非常聪明! - Eggplant

网页内容由stack overflow 提供, 点击上面的
可以查看英文原文,
原文链接