opengeodata.de

Centerline / Skeleton of Polygon with PostGIS

2015-09-10

Suppose you want to have the center line of a polygon. Further suppose you do not have access to proprietary means for this goal. PostGIS with SFCGAL comes to the rescue. SFCGAL will enable the ST_StraightSkeleton function in PostGIS and is currently available in PostGIS >2.1.x. User Zia posted a good How-To on SE. Once you are set with PostGIS and SFCGAL, you can go ahead using the following query:

with xxx as (
select objectid,
    -- dump MultiLineString into seperate parts
    (st_dump(ST_StraightSkeleton(geometry))).path[1],
    (st_dump(ST_StraightSkeleton(geometry))).geom as geometry
from table
)
select * from xxx
-- make sure the seperate parts which are within 1m of the exterior of the polygon do not get into the result
where not st_dwithin(xxx.geometry, st_exteriorring((select geometry from table)), 1)
-- get rid of some of the lose ends which do not touch any line
and ((st_touches(st_startpoint(xxx.geometry), (select st_union(geometry) from xxx)) AND st_touches(st_endpoint(xxx.geometry), (select st_union(geometry) from xxx))));

The operation will be quite costly, so better run it in pgsql2shp or ogr2ogr in order to write to a file rather than your DB application. The latter one would work like this:

ogr2ogr -f "ESRI Shapefile" shapefilename.shp PG:"host=localhost user=user dbname=db_name password=pass" -sql "the query"

After this you’ll need to clean up a bit. Or you can set a treshold for ST_Length and include it in the WHERE clause. It will not work perfectly but reasonably better than manual work in most occasions. Especially analysis on large polygons will benefit.