2.5D TINs in PostGIS
Posted by smathermather on December 18, 2013
(edited: changed TIN to TIN Z)
(edited again — function already exists as a flag in ST_DelaunayTriangles… .)
Uh, I wrote a function for nothin’…
As Regina points out in the commments, this functionality was rolled into 2.1 with a flag. Her example:
SELECT ST_AsText(ST_DelaunayTriangles( ST_Union(ST_GeomFromText( 'POLYGON((175 150, 20 40, 50 60, 125 100, 175 150))'), ST_Buffer(ST_GeomFromText(‘POINT(110 170)’), 20) ),0.001,2) ) AS dtriag;
For the record, my code is about 2% slower on a dataset that takes ~50minutes to triangulate .
Ever since ST_Delaunay made its way into PostGIS with GEOS 3.4, there’s been the promise of TINs from geometries in PostGIS. The promise is for 2.5D TINs. What are these, and why aren’t they 3D? It’s all about overhangs– imagine the edge of a building. A 3D TIN can handle that building edge with a true vertical wall, or even overhang. A 2.5D TIN is like a 2.5D raster– no overlaps allowed.
With those limitations in mind, you can have TINs today if you want:
-- A small function to convert ST_Delaunay (requires GEOS 3.4) to a 2.5D Tin -- Uses the hackerish approach of converting to text, and doing string replacement --- for format conversion. -- A Delaunay triangulation does not by itself formally make a TIN. To do this, -- we will require some modification to the output format. If we perform an ST_AsText -- on our new geometries we will see that they are POLYGON Z's wrapped in a -- GEOMETRYCOLLECTION. -- Thus there are two ways in which this is not a TIN-- the internal POLYGON Z are -- implicitly triangles, but explicitly POLYGON Zs. In addition, the external wrapper -- for the collection of triangles is a GEOMETRYCOLLECTION, not a TIN. -- Once that we have this geometry in text form, the easiest way to fix this is -- with string replacement to fix these two things, then convert back to binary -- We'll go from e.g. this: -- GEOMETRYCOLLECTION Z (POLYGON Z ((-1.14864 0.61002 -2.00108,-0.433998 -0.0288903 -2.13766, ... -- to this: -- TIN ( ((-1.14864 0.61002 -2.00108,-0.433998 -0.0288903 -2.13766, ... CREATE OR REPLACE FUNCTION AsTIN(geometry) RETURNS geometry AS $BODY$ WITH dt AS ( SELECT ST_AsText(ST_DelaunayTriangles(ST_Collect($1))) AS atext ), replacedt AS ( -- Remove polygon z designation, as TINs don't require it. SELECT replace(atext, 'POLYGON Z', '') as ttext FROM dt ), replacegc AS ( -- change leading declaration to TIN SELECT replace(ttext, 'GEOMETRYCOLLECTION Z', 'TIN Z') AS tintext from replacedt ), tingeom AS ( -- Aaaand convert back to binary. Voila! SELECT ST_GeomFromEWKT(tintext) AS the_geom FROM replacegc ) SELECT the_geom FROM tingeom $BODY$ LANGUAGE sql VOLATILE COST 100;