--- a/schema/isrs_functions.sql	Fri Mar 01 18:28:50 2019 +0100
+++ b/schema/isrs_functions.sql	Sun Mar 03 13:00:50 2019 +0100
@@ -12,19 +12,20 @@
 --  * Tom Gottfried <tom@intevation.de>
 --  * Sascha Wilde <wilde@intevation.de>
 
--- Clip an area to a stretch given by a pair of ISRS location codes.
+CREATE OR REPLACE FUNCTION best_utm(stretch isrsrange) RETURNS integer
+AS $$
+    SELECT best_utm(ST_Collect(geom::geometry))
+        FROM waterway.distance_marks_virtual
+        WHERE location_code IN (lower(stretch), upper(stretch))
+    $$
+    LANGUAGE sql
+    STABLE PARALLEL SAFE;
+
+-- Clip waterway axis to a stretch given by a pair of ISRS location codes.
 -- Uses the table waterway.distance_marks_virtual to map ISRS location codes
--- to their geo-location and the table waterway.waterway_axis to retrieve
--- perpendicular direction at these geo-locations.
--- Distance marks are assumed to be near the axis and the area passed as
--- argument is assumed to intersect with the axis
--- (use e.g. waterway area or fairway dimensions).
--- If a multipolygon is passed, the union of the polygons intersecting with the
--- relevant part of the axis is used for clipping.
-CREATE OR REPLACE FUNCTION ISRSrange_area(
-    stretch isrsrange,
-    area geometry
-) RETURNS geometry
+-- to their geo-location.
+-- Distance marks are assumed to be near the axis.
+CREATE OR REPLACE FUNCTION ISRSrange_axis(stretch isrsrange) RETURNS geometry
 AS $$
     WITH RECURSIVE
         -- Get coordinates of location codes
@@ -34,8 +35,7 @@
                     OR location_code = upper(stretch)),
         utm_zone AS (
             -- Find best matchting UTM zone
-            SELECT best_utm(ST_Collect(geom::geometry)) AS z
-                FROM points_geog),
+            SELECT best_utm(stretch) AS z),
         axis AS (
             SELECT id, ST_Transform(wtwaxs::geometry, z) AS wtwaxs
                 FROM waterway.waterway_axis, utm_zone),
@@ -84,23 +84,43 @@
             SELECT wtwaxs AS line
                 FROM axis_snapped
                 WHERE array_length(ids, 1) = (
-                    SELECT max(array_length(ids, 1)) FROM axis_snapped)),
+                    SELECT max(array_length(ids, 1)) FROM axis_snapped))
+        -- 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 leads to unexpected results of the buffer with
+        -- endcap=flat in the CTE below.
+        SELECT ST_SimplifyPreserveTopology(ST_LineSubstring(
+                    axis_segment.line, min(fractions.f), max(fractions.f)),
+                0.0001) AS line
+        FROM axis_segment, LATERAL (
+            SELECT ST_LineLocatePoint(axis_segment.line, points.geom) AS f
+                FROM points) AS fractions
+        GROUP BY axis_segment.line
+    $$
+    LANGUAGE sql
+    STABLE PARALLEL SAFE;
+
+-- Clip an area to a stretch given by a pair of ISRS location codes.
+-- Uses ISRSrange_axis() to retrieve the respective clipped axis used to find
+-- perpendicular direction at geo-locations of ISRS codes.
+-- The area passed as argument is assumed to intersect with the axis
+-- (use e.g. waterway area or fairway dimensions).
+-- If a multipolygon is passed, the union of the polygons intersecting with the
+-- relevant part of the axis is used for clipping.
+CREATE OR REPLACE FUNCTION ISRSrange_area(
+    stretch isrsrange,
+    area geometry
+) RETURNS geometry
+AS $$
+    WITH
         axis_substring AS (
-            -- 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 leads to unexpected results of the buffer with
-            -- endcap=flat in the CTE below.
-            SELECT ST_SimplifyPreserveTopology(ST_LineSubstring(
-                        axis_segment.line, min(fractions.f), max(fractions.f)),
-                    0.0001) AS line
-            FROM axis_segment, LATERAL (
-                SELECT ST_LineLocatePoint(axis_segment.line, points.geom) AS f
-                    FROM points) AS fractions
-            GROUP BY axis_segment.line),
+            SELECT ISRSrange_axis(stretch) AS line),
+        utm_zone AS (
+            SELECT best_utm(stretch) AS z),
         area_subset AS (
             -- In case area is a multipolygon, process the union of those
             -- polygons, which intersect with the axis. The union is to avoid