I have a set of 1000 polygons that overlap and I want to count the areas of overlap. I am using this blog post's instructions to create multilines from the 1000 polygons then use this line file to generate polygons for the overlapping areas and then count the overlaps. http://boundlessgeo.com/2014/10/postgis-training-creating-overlays/
This works fine for datasets of 100 or fewer polygons but hangs on anything more than that on the first step. It seems like the ST_Union is what's so slow, but ST_Collect seems not to work here because it generates a multiline with far fewer nodes and doesn't create non-overlapping polygons when used in the second step. It's too bad because ST_Collect is so fast - but doesn't give the same result. Anyone have ideas on how I can modify this process to work faster on a dataset of 1000 polygons?
Step 1 CREATE TABLE boundaries_polygons1000 AS SELECT ST_Union(ST_ExteriorRing(wkb_geometry)) AS geom FROM polygons1000;
Step 2 CREATE SEQUENCE polyseq_polygons1000; CREATE TABLE polys_polygons1000 AS SELECT nextval('polyseq_polygons1000') AS id, (ST_Dump(ST_Polygonize(geom))).geom AS geom FROM boundaries_polygons1000;
Step 3 ALTER TABLE polys_polygons1000 ADD COLUMN count INTEGER DEFAULT 0; UPDATE polys_polygons1000 set count = p.count FROM ( SELECT count(*) AS count, p.id AS id
FROM polys_polygons1000 p JOIN polygons1000 c ON ST_Contains(c.wkb_geometry, ST_PointOnSurface(p.geom)) GROUP BY p.id ) AS p WHERE p.id = polys_polygons1000.id;
أكثر...
This works fine for datasets of 100 or fewer polygons but hangs on anything more than that on the first step. It seems like the ST_Union is what's so slow, but ST_Collect seems not to work here because it generates a multiline with far fewer nodes and doesn't create non-overlapping polygons when used in the second step. It's too bad because ST_Collect is so fast - but doesn't give the same result. Anyone have ideas on how I can modify this process to work faster on a dataset of 1000 polygons?
Step 1 CREATE TABLE boundaries_polygons1000 AS SELECT ST_Union(ST_ExteriorRing(wkb_geometry)) AS geom FROM polygons1000;
Step 2 CREATE SEQUENCE polyseq_polygons1000; CREATE TABLE polys_polygons1000 AS SELECT nextval('polyseq_polygons1000') AS id, (ST_Dump(ST_Polygonize(geom))).geom AS geom FROM boundaries_polygons1000;
Step 3 ALTER TABLE polys_polygons1000 ADD COLUMN count INTEGER DEFAULT 0; UPDATE polys_polygons1000 set count = p.count FROM ( SELECT count(*) AS count, p.id AS id
FROM polys_polygons1000 p JOIN polygons1000 c ON ST_Contains(c.wkb_geometry, ST_PointOnSurface(p.geom)) GROUP BY p.id ) AS p WHERE p.id = polys_polygons1000.id;
أكثر...