My dataset currently contains about 15.000 buildings which partially touch each other. The objective is that those which do touch should belong to a cluster. This "basically" works but the main issue I am facing is the speed which is taking hours and if the dataset of polygons is too large it feels like the recursion is entering an endless loop.
-- generating the clusters table which should feature the row number and an array of building ids which make up the cluster
CREATE TABLE clusters (
rownum INTEGER NOT NULL,
ret INTEGER[]
);
-- recursion to find the polygons which are very close to each other by buffering them and finding checking if they intersect
WITH RECURSIVE
source (rownum, geom, ret) AS (
SELECT row_number() OVER (ORDER BY id ASC), ST_Multi(geom), ARRAY[id::numeric::integer] FROM buildings),
result (rownum, geom, ret, incroci) AS (
SELECT rownum, geom, ret, 1 FROM source
UNION ALL
SELECT s.rownum, s.geom, (r.ret || s.ret), (r.incroci + 1)
FROM source AS s
INNER JOIN result as r
ON s.rownum > r.rownum
AND ST_Intersects(ST_Buffer(s.geom,0.2), ST_Buffer(r.geom,0.2))
)
-- ultimately insert them to the clusters table if they contain more than 1 building
INSERT INTO clusters(rownum, ret)
SELECT row_number() OVER (ORDER BY rownum ASC), ret FROM result
WHERE array_length(ret, 1) > 1;
-- we create an index on the clusters table array column which will be needed later
CREATE INDEX clusters_idx_ret_with_intarray ON clusters USING GIN(ret gin__int_ops);
CREATE INDEX clusters_idx_rownum ON clusters (rownum);
So basically this script is recursively finding polygons which intersect with each other and afterwards inserts the results into a new table. Currently I am using a subset of the 15.000 buildings of about 1.429 buildings which takes approximately 3 minutes (however if I use a subset of 5.000 buildings it runs for hours and never stops).
The result will look something like this:
gis=# SELECT * FROM clusters ORDER BY ret;
rownum | ret
--------+--------------------------------
8 | {3,2262}
37 | {3,2262,3456}
209 | {3,2262,3456,5827}
451 | {3,2262,3456,5827,7653}
456 | {3,2262,3456,7653}
529 | {3,2262,8144}
570 | {3,2262,8144,8437}
682 | {3,2262,8144,8913}
1089 | {3,2262,8144,8913,10436}
1086 | {3,2262,8144,10436}
566 | {3,2262,8437}
38 | {3,3456}
207 | {3,3456,5827}
457 | {3,3456,5827,7653}
455 | {3,3456,7653}
208 | {3,5827}
453 | {3,5827,7653}
955 | {9,9988}
1372 | {9,9988,11689}
647 | {12,8740}
811 | {12,8740,9422}
9 | {31,2264}
1 | {33,71}
1136 | {33,71,10663}
1137 | {33,10663}
954 | {62,9988}
1373 | {62,9988,11689}
1135 | {71,10663}
28 | {91,3345}
297 | {95,6545}
541 | {95,6545,8246}
617 | {95,6545,8246,8597}
615 | {95,6545,8597}
375 | {101,7135}
7 | {118,2105}
101 | {118,2105,4445}
130 | {118,2105,4806}
316 | {118,2105,6683}
...
The obvious issue now is that the result will contain many subsets of arrays considering other rows (I guess this is got to do with the rownum comparator?).
To overcome this I added a second script which will remove subsets, e.g.
DELETE FROM clusters
WHERE ret NOT IN (
SELECT ret FROM clusters AS a
WHERE NOT EXISTS (
SELECT *
FROM clusters AS b
WHERE a.ret <@ b.ret
AND a.ret <> b.ret
)
);
Ultimately this will leave me with a table without subsets but it will yet have rows with common array members such as
570 | {3,2262,8144,8437}
682 | {3,2262,8144,8913}
...
..
And to "merge" these I figured to introduce a second recursion script similarly to https://stackoverflow.com/questions/46485699/postgresql-group-by-array-elements-in-common but this is absolutely slow. It will take 30 minutes for this small dataset and if I use a dataset with approx 5.000 buildings it will never stop running after 5 hours or more. The script in particular I am using for the array merges is
WITH RECURSIVE
cte(rownum, ret) AS (
SELECT rownum, ret FROM clusters
UNION
SELECT cs.rownum, array_merge(c.ret, cs.ret)
FROM cte AS c
JOIN clusters AS cs
ON c.rownum <> cs.rownum
AND c.ret && cs.ret
)
SELECT count(DISTINCT ret) FROM cte;
-- keeping it simple but here I join with buildings again and group them by cluster and union them...
Where array_merge is a helper function:
CREATE OR REPLACE FUNCTION public.array_merge(arr1 anyarray, arr2 anyarray)
RETURNS anyarray
LANGUAGE SQL IMMUTABLE
AS $function$
SELECT array_agg(DISTINCT elem ORDER BY elem)
FROM (
SELECT unnest(arr1) elem
UNION
SELECT unnest(arr2)
) s
$function$;
I have the feeling that I am missing something in the first recursion script and possibly the second is not even needed? Can I anyhow simplify this pipeline?
Some further notes:
the buildings have a geometry index added
to the generated clusters table I add a gin__int_ops index on the building_ids array column to speed up the && operator
WITH RECURSIVE
functionality: if you add aLIMIT <n>
to the query that fetches rows from the actual recursive term, PostgreSQL will only execute the iteration<n>
times. Good to debug.