Postgis st_intersects fails on an obtained point on the line with the line itself

Question

I am trying to split a line from its intersection with another line. I know I can use st_split directly here but I wish to split a line from several positions where it intersects other lines. What I discover is the following:

I know that lines 1 and 31 intersect

select st_intersects(a.geom,b.geom) from dumped_roads as a,
       dumped_roads as b where a.gid=1 and b.gid=31;

Output:

st_intersects 
---------------
t
(1 row)

Also, st_intersection(a.geom,b.geom) returns an intersection point.

But when I do,

select st_intersects(a.geom, st_intersection(a.geom,b.geom))
from
  dumped_roads as a,dumped_roads as b where a.gid=1 and b.gid=31;

st_intersects 
---------------
f
(1 row)

I can't understand what the problem is. Although a.geom && st_intersection(a.geom,b.geom) returns true but that wouldn't help me in splitting the line with the point obtained from intersection.

If it is indeed how it is, is there any way out on breaking a line from an obtained point which physically lies on the line, but is not in the linestrings point set?

Answer 1

This is a robustness issue. ST_Intersection creates a new point that is probably not a vertex from either of the original LineStrings, but it probably is close to the line. A more robust test of intersection is to check if the distance is zero.

SELECT ST_AsText(ST_Intersection(a.geom, b.geom)),
       ST_Intersects(a.geom, ST_Intersection(a.geom, b.geom)),
       ST_Distance(a.geom, ST_Intersection(a.geom, b.geom))
FROM (select 'LINESTRING (156 258, 607 488)'::geometry as geom) as a,
     (select 'LINESTRING (105 418, 737 218)'::geometry as geom) as b;

-[ RECORD 1 ]-+-----------------------------------------
st_astext     | POINT(330.074206147054 346.773985396502)
st_intersects | f
st_distance   | 0

Answer 2

Ok, so I have been able to work out the task of splitting a line from its intersection with other lines. There is indeed a robustness issue with postgis in which the intersection point of two lines doesnt lie on either of the lines sometimes even though the distance of the intersection point from the line is 0. So I used st_line_locate_point to locate the point on the line closest to the intersection point (which is the intersection point itself since its distance is 0 from the line), and then I used st_line_substring to split the line.

I can explain it here using an example:

Suppose I want to split line joining (0,0) and (10,10) from its intersection with the line joining (0,10) and (10,0), I write the following query:

    select ST_ASTEXT( 
      st_line_substring( st_geomfromtext('LINESTRING(0 0, 10 10)'),0,
         st_line_locate_point( st_geomfromtext('LINESTRING(0 0, 10 10)'), 
            st_intersection( 
              st_geomfromtext('LINESTRING(0 0, 10 10)'), st_geomfromtext('LINESTRING(0 10, 10 0)')  
                       )                 
                   ) 
                )
      );

The output is correctly:

    st_astext      
    ---------------------
    LINESTRING(0 0,5 5)
    (1 row)

Answer 3

try

select st_intersects(a.geom, select st_intersection(a.geom,b.geom) from
    dumped_roads as a,dumped_roads as b where a.gid=1 and b.gid=31);