Performing polygon drilldown/overlay in PostGIS?
I think this works.
It is a windowing function, getting the difference between the intersection of each geometries intersection with the query box and the union of the preceding geometries.
The coalesce is needed since the union of the preceding geometries for the first geometry is null which gives a null result,instead of what is desired.
WITH a(geom, rank, val) AS
(
VALUES
('POLYGON((1 1, 1 5, 5 5, 5 1, 1 1))'::geometry,3,0.1),
('POLYGON((2 3, 2 8, 4 8, 5 3,2 3))'::geometry,1,0.2),
('POLYGON((4 4, 4 6, 6 6, 6 4,4 4))'::geometry,2,0.2)
)
,q AS
(
SELECT 'POLYGON((3 3, 3 4.5, 12 4.5,12 3, 3 3))'::geometry geom
)
SELECT
ST_Difference(
ST_Intersection(a.geom, q.geom),
COALESCE(ST_Union(a.geom)
OVER (ORDER BY rank ROWS BETWEEN UNBOUNDED PRECEDING and 1 PRECEDING),
'POLYGON EMPTY'::geometry)
) geom
FROM a,q
WHERE ST_Intersects(a.geom, q.geom);
I am not sure how it performs though. But since both ST_Union and ST_Intersection are marked immutable it might not be that bad.
A bit of a different approach to this. There is a caveat that I don't know how it will scale performance wise, but on an indexed table it should be ok. It performs about the same as Nicklas's query (a tad slower?), but the measurement on such a small sample is fraught.
It looks a lot uglier than Nicklas's query, but it does avoid recursion in the query.
WITH
a(geom, rank, val) AS
(
VALUES
('POLYGON((1 1, 1 5, 5 5, 5 1, 1 1))'::geometry,3,0.1),
('POLYGON((2 3, 2 8, 4 8, 5 3, 2 3))'::geometry,1,0.2),
('POLYGON((4 4, 4 6, 6 6, 6 4, 4 4))'::geometry,2,0.2)
),
polygonized AS (
-- get the basic building block polygons
SELECT (ST_Dump( -- 5. Dump out the polygons
ST_Polygonize(line)) -- 4. Polygonise the linework
).geom AS mypoly
FROM (
SELECT
ST_Node( -- 3. Node lines on intersection
ST_Union( -- 2. Union them for noding
ST_Boundary(geom) -- 1. Change to linestrings
)
)
AS line
FROM a
) line
),
overlayed AS (
-- Overlay with original polygons and get minimum rank.
-- Union and dissolve interior boundaries for like ranks.
SELECT (ST_Dump(ST_UnaryUnion(geom))).geom, rank
FROM (
-- union up the polygons by rank, unary union doesn't count as an aggregate function?
SELECT ST_Union(mypoly) geom, rank
FROM (
-- get the minimum rank for each of the polygons
SELECT p.mypoly, min(rank) rank
FROM polygonized p
INNER JOIN a ON ST_Contains(a.geom,ST_PointOnSurface(p.mypoly))
GROUP BY p.mypoly
) g
GROUP BY rank
) r
)
-- get the intersection of the query area with the overlayed polygons
SELECT ST_Intersection(o.geom,'POLYGON((3 3, 3 4.5, 12 4.5,12 3, 3 3))'::Geometry), rank
FROM overlayed o
WHERE ST_Intersects(o.geom,'POLYGON((3 3, 3 4.5, 12 4.5,12 3, 3 3))'::Geometry) and
-- Intersection can do funky things
GeometryType(ST_Intersection(o.geom,'POLYGON((3 3, 3 4.5, 12 4.5,12 3, 3 3))'::Geometry)) like '%POLYGON';