5

I have a table of polygons (land parcels in a city) and I am trying to identify "chains" of 5 or more adjacent polygons. For example, if each letter is a polygon and each . is empty space, AB.CD.EFG.HIJKL.MNOP.QRSTU.VWXYZ then I want HIJKL, QRSTU and VWXYZ. In real life they aren't in a line like that of course.

The messy solution I used for a first attempt is this query:

SELECT p1.parcel_c as "Parcel 1", p2.parcel_c as "Parcel 2", p3.parcel_c as "Parcel 3", p4.parcel_c as "Parcel 4", p5.parcel_c as "Parcel 5", ST_Union(array[p1.geom, p2.geom, p3.geom, p4.geom, p5.geom]) as merged_geom FROM abandonded_vacant_lots_table p1 INNER JOIN abandonded_vacant_lots_table as p2 ON ST_TOUCHES(p1.geom, p2.geom) INNER JOIN abandonded_vacant_lots_table as p3 ON ST_TOUCHES(p2.geom, p3.geom) INNER JOIN abandonded_vacant_lots_table as p4 ON ST_TOUCHES(p3.geom, p4.geom) INNER JOIN abandonded_vacant_lots_table as p5 ON ST_TOUCHES(p4.geom, p5.geom) WHERE p1.parcel_c not in(p2.parcel_c, p3.parcel_c, p4.parcel_c, p5.parcel_c) and p2.parcel_c not in(p1.parcel_c, p3.parcel_c, p4.parcel_c, p5.parcel_c) and p3.parcel_c not in(p1.parcel_c, p2.parcel_c, p4.parcel_c, p5.parcel_c) and p4.parcel_c not in(p1.parcel_c, p2.parcel_c, p3.parcel_c, p5.parcel_c) and p5.parcel_c not in(p1.parcel_c, p2.parcel_c, p3.parcel_c, p4.parcel_c)

As you can see this is clumsy and inelegant and I suspect inefficient. It also returns a lot of almost duplicates when there is a chain of more than 5 parcels, eg if the parcels are ABCDEF it returns ABCDE and BCDEF when a better solution would return the longest chain, ABCDEF. How can I make this query more efficient, more elegant and in general more versatile?

4

I think I have a solution for postgis 2.0+:

First subquery to select the dump of unioned polygons to get only simple polygons and no multies and also get rid of internal edges of touching polygons:

with merged_patches as
(select (st_dump(st_union(geom))).geom as geom from abandonded_vacant_lots_table),

2nd subquery dumps the simple geoms from the table:

all_patches as 
    (select (st_dump(geom)).geom from abandonded_vacant_lots_table)

The query will return the count of adjacent polygons in patches and the merged patch geometry:

select count(*),merged_patches.geom 
from merged_patches,all_patches
where

Where the bbox overlap to speed things up:

    all_patches.geom && merged_patches.geom AND

The real condition is here: we're checking if the current patch is within a merged patch:

    ST_Within(all_patches.geom,merged_patches.geom)

Let's group these by the merged_patches:

group by merged_patches.geom

if you want to return merged_patches with more than 5 sub-patches:

having count(*)>=5;

Alltogether:

with
merged_patches as
    (select (st_dump(st_union(geom))).geom as geom from abandonded_vacant_lots_table),
all_patches as 
(select (st_dump(geom)).geom from abandonded_vacant_lots_table)
select count(*),merged_patches.geom 
from merged_patches,all_patches
where
    all_patches.geom && merged_patches.geom AND
    ST_Within(all_patches.geom,merged_patches.geom)
group by merged_patches.geom
having count(*)>=5;

I had a somewhat similar problem (Remove islands and completely surrounded polygons after polygonization with QGIS) that inspired this answer.

7

If your "chains" don't necessarily need to be a straight line (i.e, if they can fork), then you can tackle this with the new clustering functions in PostGIS 2.2. For example:

SELECT geom FROM 
  (SELECT unnest(ST_ClusterIntersecting(geom)) AS geom 
   FROM abandonded_vacant_lots_table) clusters 
WHERE ST_NumGeometries(clusters.geom) >= 5

In this query, ST_ClusterIntersecting returns an array of GeometryCollection entities, each reflecting a potentially forked "chain." The PostgreSQL function unnest turns the array entries into a recordset, forming the clusters subquery. Then, we're selecting only the chains having 5 or more geometries.

The easiest way to link these results back to your original table ("what parcel IDs are in the clusters") is to use a spatial join between the clusters and an interior point of your original geometries.

You may find ST_ClusterWithin to be more useful, especially if you have parcel polygons that almost touch but don't quite.

  • Wow those are awesomely useful functions! I don't have access to PostGIS 2.2 yet. As soon as 2.2 is available on Ubuntu I'll give it a shot. – Chris Hartley Oct 22 '15 at 18:51

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.