[postgis-users] Best way to combine a polygon

[Brent Wood]

--- Rob Tester <robtester at gmail.com> wrote:

> I have a table with many single ring polygons. Some of the polygons rings in
> the table may be contained by other polygons in the same table. What I would
> like to do is to assemble all the contained rings into a single polygon. But
> it gets a little more convoluted than that. Some of the rings contained by
> another poly may contain rings of other poly's too, so: 
> 
> Poly1Ring contains poly2ring, poly3ring and poly4ring
> 
> Poly2Ring contains poly3ring.
>  
> What is the best method to combine these polygons to effectively cut holes
> where the rings are contained?

Hi Rob,

I have an idea as to how I'd go about it, but it may not be the best way :-)


What you want to do, if I understand you correctly, is generate the geomunion
of all the polygons within a polygon, then get the difference between this
polygon & the geomunion.

I think it can be done in a single update statement, but I'd scribble out a
script to iterate through each polygon, coz I'm more comfortable with that way,
& because I've found such SQL's in one can use very large amounts of memory,
whereas breaking things down into iterated simple steps can be faster. You can
also check the results from each step before actually committing anything....

Perl would be faster than bash, but this is quick & easy for me...


So (in pseudo bash/SQL syntax) something along the lines of:



#!/bin/bash

psql -d db -Atc "select id from table;" > file

# now iterate through the id's to fix each one

while read id ; do
  # get the geomunion of all polygons withion the current one
  # but specifically exclude the current one

  UNION=`psql -d db -Atc "select geomunion(geom) from table
                          where geom && (select geom from table
                                         where id=$id)
                            and within(geom, (select geom from table
                                       where id=$id)
                                   and id != $id;"

  # now set the geometry for that polygon to itself less the geomunion

  psql -d db -c "update table set geom=difference(geom, $UNION)
                 where id=$id;"

done < file


Cheers,

  Brent Wood