--- /dev/null	Thu Jan 01 00:00:00 1970 +0000
+++ b/schema/updates/1423/02.adapt_func.sql	Wed Mar 11 17:11:23 2020 +0100
@@ -0,0 +1,105 @@
+CREATE OR REPLACE FUNCTION ISRSrange_axis(
+    stretch isrsrange,
+    tolerance float
+    -- in m, up to which linestrings will be connected at their boundary
+) RETURNS geometry
+AS $$
+DECLARE z int;
+DECLARE result_geom geometry;
+BEGIN
+    -- Find best matchting UTM zone
+    z = best_utm(stretch);
+
+    CREATE TEMP TABLE axis AS
+        SELECT row_number() OVER () AS id,
+                geom AS wtwaxs,
+                ST_Boundary(geom) AS bdr
+            FROM waterway.waterway_axis,
+                ST_Dump(ST_Transform(wtwaxs::geometry, z));
+    CREATE INDEX axs_bdr ON axis USING GiST (bdr);
+    ANALYZE axis;
+
+    WITH RECURSIVE
+        -- In order to guarantee the following ST_Covers to work,
+        -- snap distance mark coordinates to axis
+        points0 AS (
+            SELECT ST_ClosestPoint(
+                    wtwaxs,
+                    ST_Transform(geom, z)) AS geom
+                FROM ST_Dump(ISRSrange_points(stretch)), (
+                    SELECT ST_Collect(wtwaxs) AS wtwaxs
+                        FROM axis) AS ax),
+        -- Ensure two distinct points on axis have been found
+        points AS (
+            SELECT geom
+                FROM points0
+                WHERE 2 = (SELECT count(DISTINCT geom) FROM points0)),
+        axis_snapped AS (
+            -- Iteratively connect non-contiguous axis chunks
+            -- to find the contiguous axis on which given distance marks lie
+            (SELECT ARRAY[id] AS ids, wtwaxs
+                FROM axis, points
+                WHERE ST_Intersects(
+                    ST_Buffer(axis.wtwaxs, 0.0001), points.geom)
+                FETCH FIRST ROW ONLY)
+            UNION
+            -- Connect endpoint of next linestring with closest
+            -- endpoint of merged linestring until a contiguous
+            -- linestring connecting both distance marks is build up
+            (SELECT refids || id,
+                    ST_LineMerge(ST_Collect(ARRAY(
+                        -- Linestring build up so far
+                        SELECT refgeom
+                        UNION
+                        -- Fill eventual gap
+                        SELECT ST_ShortestLine(
+                                ST_Boundary(refgeom), bdr)
+                        UNION
+                        -- Linestring to be added
+                        SELECT geom)))
+                FROM axis_snapped AS axis_snapped (refids, refgeom),
+                    axis AS axis (id, geom, bdr),
+                    (SELECT ST_Collect(points.geom) AS pts
+                        FROM points) AS points
+                WHERE id <> ALL(refids)
+                    AND ST_DWithin(
+                        ST_Boundary(refgeom), bdr, tolerance)
+                    AND NOT ST_Covers(ST_Buffer(refgeom, 0.0001), points.pts)
+                ORDER BY ST_Boundary(refgeom) <-> bdr
+                FETCH FIRST ROW ONLY)),
+        axis_segment AS (
+            -- Fetch end result from snapping
+            SELECT wtwaxs AS line
+                FROM axis_snapped,
+                    (SELECT ST_Collect(points.geom) AS pts
+                        FROM points) AS points
+                -- Return end result only if both distance marks were connected
+                WHERE ST_Covers(ST_Buffer(wtwaxs, 0.0001), points.pts))
+        -- Use linear referencing to clip axis between distance marks.
+        -- Simplification is used to work-around the problem, that
+        -- ST_LineSubstring might generate very small line segments at an
+        -- end of the resulting linestring, that significantly differ from
+        -- the direction of the input linestring due to finite precision
+        -- of the calculation. The generated small segment of the
+        -- resulting line would lead e.g. to unexpected results in an area
+        -- generated by ISRSrange_area().
+        SELECT ST_SimplifyPreserveTopology(ST_LineSubstring(
+                    axis_segment.line, min(fractions.f), max(fractions.f)),
+                0.0001) AS line
+        INTO STRICT result_geom
+        FROM axis_segment, LATERAL (
+            SELECT ST_LineLocatePoint(axis_segment.line, points.geom) AS f
+                FROM points) AS fractions
+        GROUP BY axis_segment.line;
+
+    -- Drop temporary table to avoid side effects on PostgreSQL's MVCC,
+    -- because otherwise subsequent invocations of the function will not see
+    -- changes on the underlying waterway.waterway_axis that might have
+    -- occured.
+    DROP TABLE axis;
+
+    RETURN result_geom;
+END;
+    $$
+    LANGUAGE plpgsql
+    PARALLEL RESTRICTED;