fix a number of corner cases with interpolation splitting

Snapping a line to a point before splitting was meant to ensure
that the split point is really on the line. However, ST_Snap() does
not always behave well for this case. It may shorten the interpolation
line in some cases with the result that two points housenumbers
suddenly fall on the same point. It might also shorten the line down
to a single point which then makes ST_Split() crash.

Switch to a combination of ST_LineLocatePoint and ST_LineSubString
instead, which guarantees to keep the original geometry. Explicitly
handle the corner cases, where the split point falls on the beginning
or end of the line.

diff --git a/lib-sql/functions/interpolation.sql b/lib-sql/functions/interpolation.sql
index 9bb9102172b81169028ba302a52fa8d7b7676a0b..928d55c546294bb37b7da710473425442c524a59 100644 (file)
--- a/lib-sql/functions/interpolation.sql
+++ b/lib-sql/functions/interpolation.sql
@@ -164,7 +164,7 @@ DECLARE
   newend INTEGER;
   moddiff SMALLINT;
   linegeo GEOMETRY;
   newend INTEGER;
   moddiff SMALLINT;
   linegeo GEOMETRY;

-  splitline GEOMETRY;
+  splitpoint FLOAT;

   sectiongeo GEOMETRY;
   postcode TEXT;
   stepmod SMALLINT;
   sectiongeo GEOMETRY;
   postcode TEXT;
   stepmod SMALLINT;


@@ -223,15 +223,27 @@ BEGIN
         FROM placex, generate_series(1, array_upper(waynodes, 1)) nodeidpos
         WHERE osm_type = 'N' and osm_id = waynodes[nodeidpos]::BIGINT
               and address is not NULL and address ? 'housenumber'
         FROM placex, generate_series(1, array_upper(waynodes, 1)) nodeidpos
         WHERE osm_type = 'N' and osm_id = waynodes[nodeidpos]::BIGINT
               and address is not NULL and address ? 'housenumber'

+              and ST_Distance(NEW.linegeo, geometry) < 0.0005

         ORDER BY nodeidpos
     LOOP
       {% if debug %}RAISE WARNING 'processing point % (%)', nextnode.hnr, ST_AsText(nextnode.geometry);{% endif %}
       IF linegeo is null THEN
         linegeo := NEW.linegeo;
       ELSE
         ORDER BY nodeidpos
     LOOP
       {% if debug %}RAISE WARNING 'processing point % (%)', nextnode.hnr, ST_AsText(nextnode.geometry);{% endif %}
       IF linegeo is null THEN
         linegeo := NEW.linegeo;
       ELSE

-        splitline := ST_Split(ST_Snap(linegeo, nextnode.geometry, 0.0005), nextnode.geometry);
-        sectiongeo := ST_GeometryN(splitline, 1);
-        linegeo := ST_GeometryN(splitline, 2);
+        splitpoint := ST_LineLocatePoint(linegeo, nextnode.geometry);
+        IF splitpoint = 0 THEN
+          -- Corner case where the splitpoint falls on the first point
+          -- and thus would not return a geometry. Skip that section.
+          sectiongeo := NULL;
+        ELSEIF splitpoint = 1 THEN
+          -- Point is at the end of the line.
+          sectiongeo := linegeo;
+          linegeo := NULL;
+        ELSE
+          -- Split the line.
+          sectiongeo := ST_LineSubstring(linegeo, 0, splitpoint);
+          linegeo := ST_LineSubstring(linegeo, splitpoint, 1);
+        END IF;

       END IF;
 
       IF prevnode.hnr is not null
       END IF;
 
       IF prevnode.hnr is not null


@@ -239,6 +251,9 @@ BEGIN
          -- regularly mapped housenumbers.
          -- (Conveniently also fails if one of the house numbers is not a number.)
          and abs(prevnode.hnr - nextnode.hnr) > NEW.step
          -- regularly mapped housenumbers.
          -- (Conveniently also fails if one of the house numbers is not a number.)
          and abs(prevnode.hnr - nextnode.hnr) > NEW.step

+         -- If the interpolation geometry is broken or two nodes are at the
+         -- same place, then splitting might produce a point. Ignore that.
+         and ST_GeometryType(sectiongeo) = 'ST_LineString'

       THEN
         IF prevnode.hnr < nextnode.hnr THEN
           startnumber := prevnode.hnr;
       THEN
         IF prevnode.hnr < nextnode.hnr THEN
           startnumber := prevnode.hnr;


@@ -300,12 +315,12 @@ BEGIN
                   NEW.address, postcode,
                   NEW.country_code, NEW.geometry_sector, 0);
         END IF;
                   NEW.address, postcode,
                   NEW.country_code, NEW.geometry_sector, 0);
         END IF;

+      END IF;

 
 

-        -- early break if we are out of line string,
-        -- might happen when a line string loops back on itself
-        IF ST_GeometryType(linegeo) != 'ST_LineString' THEN
-            RETURN NEW;
-        END IF;
+      -- early break if we are out of line string,
+      -- might happen when a line string loops back on itself
+      IF linegeo is null or ST_GeometryType(linegeo) != 'ST_LineString' THEN
+          RETURN NEW;

       END IF;
 
       prevnode := nextnode;
       END IF;
 
       prevnode := nextnode;

diff --git a/test/bdd/db/import/interpolation.feature b/test/bdd/db/import/interpolation.feature
index 0b939f593b4649803985ac42cab6fd9c5f5c83e8..3c67d652e42fa59bc01a658b4b820b784c1891c1 100644 (file)
--- a/test/bdd/db/import/interpolation.feature
+++ b/test/bdd/db/import/interpolation.feature
@@ -29,14 +29,14 @@ Feature: Import of address interpolations
           | N2  | place | house  | 8       |
         And the places
           | osm | class | type   | addr+interpolation | geometry |
           | N2  | place | house  | 8       |
         And the places
           | osm | class | type   | addr+interpolation | geometry |

-          | W1  | place | houses | even    | 1,2 |
+          | W1  | place | houses | even               | 2,1      |

         And the ways
           | id | nodes |
           | 1  | 2,1 |
         When importing
         Then W1 expands to interpolation
           | start | end | geometry |
         And the ways
           | id | nodes |
           | 1  | 2,1 |
         When importing
         Then W1 expands to interpolation
           | start | end | geometry |

-          | 4     | 6   | 8,9 |
+          | 4     | 6   | 9,8      |

 
     Scenario: Simple odd two point interpolation
         Given the grid with origin 1,1
 
     Scenario: Simple odd two point interpolation
         Given the grid with origin 1,1


@@ -341,7 +341,7 @@ Feature: Import of address interpolations
         Then W1 expands to interpolation
           | start | end | geometry |
           | 4     | 4   | 144.963016 -37.762946 |
         Then W1 expands to interpolation
           | start | end | geometry |
           | 4     | 4   | 144.963016 -37.762946 |

-          | 8     | 8   | 144.963144 -37.7622237 |
+          | 8     | 8   | 144.96314407 -37.762223692 |

 
     Scenario: Place with missing address information
         Given the grid
 
     Scenario: Place with missing address information
         Given the grid


@@ -456,3 +456,69 @@ Feature: Import of address interpolations
           | foo   |
           | x     |
           | 12-2  |
           | foo   |
           | x     |
           | 12-2  |

+
+
+    Scenario: Interpolation line where points have been moved (Github #3022)
+        Given the 0.00001 grid
+         | 1 | | | | | | | | 2 | 3 | 9 | | | | | | | | 4 |
+        Given the places
+          | osm | class | type   | housenr | geometry |
+          | N1  | place | house  | 2       | 1 |
+          | N2  | place | house  | 18      | 3 |
+          | N3  | place | house  | 24      | 9 |
+          | N4  | place | house  | 42      | 4 |
+        And the places
+          | osm | class | type   | addr+interpolation | geometry |
+          | W1  | place | houses | even               | 1,2,3,4  |
+        And the ways
+          | id | nodes   |
+          | 1  | 1,2,3,4 |
+        When importing
+        Then W1 expands to interpolation
+          | start | end |
+          | 4     | 16  |
+          | 20    | 22  |
+          | 26    | 40  |
+
+
+    Scenario: Interpolation line with duplicated points
+        Given the grid
+          | 7 | 10 | 8 | 11 | 9 |
+        Given the places
+          | osm | class | type   | housenr | geometry |
+          | N1  | place | house  | 2       | 7 |
+          | N2  | place | house  | 6       | 8 |
+          | N3  | place | house  | 10      | 8 |
+          | N4  | place | house  | 14      | 9 |
+        And the places
+          | osm | class | type   | addr+interpolation | geometry |
+          | W1  | place | houses | even               | 7,8,8,9  |
+        And the ways
+          | id | nodes   |
+          | 1  | 1,2,3,4 |
+        When importing
+        Then W1 expands to interpolation
+          | start | end | geometry |
+          | 4     | 4   | 10       |
+          | 12    | 12  | 11       |
+
+
+    Scenario: Interpolaton line with broken way geometry (Github #2986)
+        Given the grid
+          | 1 | 8 | 10 | 11 | 9 | 2 | 3 | 4 |
+        Given the places
+          | osm | class | type   | housenr |
+          | N1  | place | house  | 2       |
+          | N2  | place | house  | 8       |
+          | N3  | place | house  | 12      |
+          | N4  | place | house  | 14      |
+        And the places
+          | osm | class | type   | addr+interpolation | geometry |
+          | W1  | place | houses | even               | 8,9      |
+        And the ways
+          | id | nodes       |
+          | 1  | 1,8,9,2,3,4 |
+        When importing
+        Then W1 expands to interpolation
+          | start | end | geometry |
+          | 4     | 6   | 10,11    |