[postgis-users] Topology: Simplifying Coastlines with Islands

[Sandro Santilli]

On Mon, May 28, 2012 at 11:56:13AM +0200, Michiel J. van Heek wrote:

> Sandro,
> 
> Based on your suggestions, I changed the script in three ways (see
> below for the result):
> 
> 1. The number of iterations is now a parameter and running down the
> tolerance uses linear steps for each iteration (0.5, 0.45, 0.4,
> 0.35, etc.). The division by 2 solution (0.5, 0.25, 0.125, etc.) put
> a heavy weight on the lower tolerances, whereas it is often
> worthwile to also try higher tolerances (between initial tolerance
> and half of the initial tolerance).
> 
> 2. If an iteration fails because of "not simple" the same iteration
> is run again, but with "PreserveTopology" now.
> 
> 3. When all iterations fail, the function returns NULL now, and does
> not throw an exception. It is now possible to find "unwilling" edges
> by running:
> 
> SELECT *
> FROM (
>     SELECT
>         edge_id,
>         SimplifyEdgeGeom('countries_topology_4', edge_id, 0.5) AS tolerance
>     FROM countries_topology_4.edge
> ) AS foo
> WHERE (tolerance IS NULL);
> 
> Also, I removed the STABLE attribute, because the funtion *does*
> make changes to the database.

Good work ! Maybe make the "curve not simple" handling more tolerant of
eventual changes in the exception text (only look for "not simple"?).

How do your country boundaries look now ? :)

--strk;

  ,------o-. 
  |   __/  |    Delivering high quality PostGIS 2.0 !
  |  / 2.0 |    http://strk.keybit.net - http://vizzuality.com
  `-o------'


> CREATE OR REPLACE FUNCTION SimplifyEdgeGeom(topology_param text,
> edge_id_param integer, tolerance_param double precision,
> num_iterations_param integer DEFAULT 10) RETURNS double precision AS
> $BD$
> DECLARE
>     i integer := 1;
>     xi text := '';
>     tolerance_var double precision := tolerance_param;
>     sql_var text;
> BEGIN
>     WHILE (i <= num_iterations_param) LOOP
>         sql_var := 'SELECT topology.ST_ChangeEdgeGeom(' ||
> quote_literal(topology_param) || ', ' || edge_id_param
>           || ', ST_Simplify' || xi || '(geom, ' || tolerance_var || ')) FROM '
>           || quote_ident(topology_param) || '.edge WHERE edge_id = '
> || edge_id_param;
>         BEGIN
>             RAISE DEBUG 'Running %', sql_var;
>             EXECUTE sql_var;
> 
>             RETURN tolerance_var;
>         EXCEPTION
>             WHEN OTHERS THEN
>                 RAISE WARNING 'Simplification of edge % failed in
> iteration %_%, with tolerance %: %',
>                     edge_id_param, i, xi, tolerance_var, SQLERRM;
> 
>                 IF (SQLERRM = 'SQL/MM Spatial exception - curve not
> simple') AND (xi = '') THEN
>                     xi := 'PreserveTopology';
>                 ELSE
>                     i := i + 1;
>                     xi := '';
>                     tolerance_var := tolerance_var -
> (tolerance_param / num_iterations_param);
>                 END IF;
>             END;
>         -- END EXCEPTION
>     END LOOP;
> 
>     RETURN NULL;
> END;
> $BD$ LANGUAGE 'plpgsql' STRICT;