[PostGIS] #5898: Incorrect output from ST_LineSubstring(geography, float8, float8)

[PostGIS]

#5898: Incorrect output from ST_LineSubstring(geography, float8, float8)
--------------------------------+---------------------------
 Reporter:  atlight             |      Owner:  pramsey
     Type:  defect              |     Status:  new
 Priority:  medium              |  Milestone:  PostGIS 3.5.4
Component:  postgis             |    Version:  3.5.x
 Keywords:  linear referencing  |
--------------------------------+---------------------------
 Consider the following query:

 {{{
 WITH
   line AS (
     SELECT ST_CurveToLine(ST_GeomFromText(
       'CIRCULARSTRING(143 -35,143.1 -34.9,143.2 -35)', 4326)) AS geom
   ),
   line_third_vertex AS (
     SELECT ST_PointN(geom, 3) AS pt
     FROM line
   ),
   line_third_vertex_linref AS (
     SELECT ST_LineLocatePoint(geom::geography, ST_PointN(geom,
 3)::geography) AS lr
     FROM line
   ),
   line_third_vertex_via_substring AS (
     SELECT ST_EndPoint(ST_LineSubstring(geom::geography, 0, lr)::geometry)
 AS pt
     FROM line, line_third_vertex_linref
   )
 SELECT
   ST_AsText(line_third_vertex.pt) AS line_third_vertex_pt,
   ST_AsText(line_third_vertex_via_substring.pt) AS
 line_third_vertex_via_substring_pt,
   ST_Distance(
     line_third_vertex.pt::geography,
     line_third_vertex_via_substring.pt::geography
   ) AS distance
 FROM line_third_vertex, line_third_vertex_via_substring;
 }}}

 The query finds the linear reference of a vertex along a `LineString`
 using ST_LineLocatePoint and attempts to recreate that vertex by finding
 the endpoint of the ST_LineSubstring that ends at this linear reference.

 The two points returned by the query should obviously be the same.
 However, they are 272 m apart:

 {{{
 -[ RECORD 1
 ]----------------------+---------------------------------------------
 line_third_vertex_pt               | POINT(143.00048152733277
 -34.99019828596705)
 line_third_vertex_via_substring_pt | POINT(143.00030096290968
 -34.99264621161169)
 distance                           | 272.07377371
 }}}

 The issue appears restricted to ST_LineSubstring(**geography**, float8,
 float8) and does not occur in the geometry overload. Removing all
 occurrences of `::geography` from the query gives the following result:

 {{{
 -[ RECORD 1
 ]----------------------+---------------------------------------------
 line_third_vertex_pt               | POINT(143.00048152733277
 -34.99019828596705)
 line_third_vertex_via_substring_pt | POINT(143.00048152733277
 -34.99019828596705)
 distance                           | 0
 }}}

 We can verify that ST_LineLocatePoint(geography, geography) is working
 properly:

 {{{
 WITH
   line AS (
     SELECT ST_CurveToLine(ST_GeomFromText(
       'CIRCULARSTRING(143 -35,143.1 -34.9,143.2 -35)', 4326)) AS geom
   ),
   line_third_vertex_linref AS (
     SELECT ST_LineLocatePoint(geom::geography, ST_PointN(geom,
 3)::geography) AS lr
     FROM line
   ),
   line_first_three_vertices AS (
     SELECT ST_MakeLine(vtx) AS geom
         FROM (SELECT (ST_DumpPoints(geom)).geom AS vtx
       FROM line
       LIMIT 3)
   )
 SELECT lr,
 ST_Length(line_first_three_vertices.geom::geography)/ST_Length(line.geom::geography)
 AS true_lr
 FROM line, line_third_vertex_linref, line_first_three_vertices;

           lr          |       true_lr
 ----------------------+----------------------
  0.034174334732531814 | 0.034174334732531814
 }}}

 The fault must lie with ST_LineSubstring.

 Version: 3.5
-- 
Ticket URL: <https://trac.osgeo.org/postgis/ticket/5898>
PostGIS <http://trac.osgeo.org/postgis/>
The PostGIS Trac is used for bug, enhancement & task tracking, a user and developer wiki, and a view into the subversion code repository of PostGIS project.