[postgis-users] Re: computing unions of intersecting polygons
mbdavis at refractions.net
Wed Jan 14 16:32:31 PST 2009
Reid Priedhorsky wrote:
> On 01/13/2009 11:44 AM, Martin Davis wrote:
>> The way I'd look at doing this is to compute a relation which
>> contains pairs (poly_id1, poly_id2) for all polygons which
>> intersect. Then compute the transitive closure of this relation to
>> find all connected sets of polygons. (Until the recursive CTE
>> support emerges in PostgresQL 8.4, you'll probably need to do this
>> procedurally in pgplsql). Label each group with the smallest id in
>> that group. Then you can use an aggregate query over the transitive
>> closure table to compute unions and counts.
> Thanks, Martin!
> I ended up doing it procedurally, though without the intermediate
> table you suggest, because in our case it was OK to muck with the
> table itself (wh_raw) to keep track of progress!
> Below is my solution! It runs in about 2 minutes, compared to the
> st_union solution which takes about 6.5 hours!
>> /* This function compacts the table wh_raw by the transitive closure
>> intersection: for each set of polygons in wh_raw that intersect on
>> (transitively), remove each member of that set and insert the
>> union! This
>> solution is procedural to work around PostGIS's limitations on this
>> operation! */
>> create function my_union() returns void
>> as $$
>> id_i int;
>> count_i int;
>> geometry_i geometry;
>> count_u int;
>> union_u geometry;
>> -- Fetch row from the table which is not done!
>> select id, count_, geometry
>> into id_i, count_i, geometry_i
>> from wh_raw where not done limit 1;
>> -- If no such rows, algorithm is complete!
>> exit when not found;
>> -- Fetch the union of all rows that intersect that row!
>> select sum(count_),
>> into count_u, union_u
>> from wh_raw
>> where not done and ST_Intersects(geometry_i, geometry);
>> if (count_i = count_u) then
>> -- There's only one; we're done with that set!
>> update wh_raw set done = True where id = id_i;
>> -- Update the initial row with the union and remove others!
>> update wh_raw
>> set count_ = count_u, geometry = union_u
>> where id = id_i;
>> delete from wh_raw
>> where id != id_i and ST_Intersects(geometry_i, geometry);
>> end if;
>> end loop;
>> $$ language plpgsql;
>> select my_union();
>> drop function my_union();
>> \qecho -- (sanity check: this should give no results!)
>> select a.id, b.id
>> from wh_raw a join wh_raw b on (a.id != b.id
>> and ST_Intersects(a.geometry,
> postgis-users mailing list
> postgis-users at postgis.refractions.net
Senior Technical Architect
Refractions Research, Inc.
More information about the postgis-users