3

I have a strange st_difference behaviour. I’ve got two overlapping lines (each one has two points), but st_difference result is complete geometry of the first line.

The code to reproduce is here:

WITH
mylines AS
(
SELECT
 st_geomfromtext('LINESTRING(4377625.44 7399187.54,4377633.47 7398278.09)') as ro_geom,
 st_geomfromtext('LINESTRING(4377625.44 7399187.54,4377630.67957507 7398594.12386704)') as lk_geom
)

select
  1 as gid,
  st_difference(ro_geom, lk_geom) as diff_geom
from mylines

Visually, it looks like this: enter image description here

So the st_difference between red and blue lines returns a green line, which isn't expected. I'd expect the green line to be shorter and below the blue one (without overlapping it).

This is just a test example I've made for demonstration of the problem. In practice, I'm using st_difference for two huge datasets and every difference looks unexpected.

Any clues?

1 Answer 1

6

At the scale of your display, these lines appear to be coincident, but in reality they're not:

enter image description here

I've exaggerated the difference between the two lines using the "Magnify Topology" tool in JTS TestBuilder, which I highly recommend for looking into cases like this.

Depending on the goal of your analysis, it might be suitable to snap your vertices to a grid first:

ST_AsText(ST_Difference(ST_SnapToGrid(ro_geom, 10), ST_SnapToGrid(lk_geom, 10)))
6
  • Nice tool, didn't know about that. You beat me to the answer, but another less pretty way of showing the same thing would have been .... SELECT ST_intersection(ro_geom, lk_geom)) as diff_geom which, of course, just returns a single point. Commented May 4, 2016 at 15:09
  • Thanks @dbaston. The smaller line is made by st_locatebetween function from the main line. How come it's second node isn't located on the parent line?
    – mofoyoda
    Commented May 4, 2016 at 15:09
  • 1
    It's the cruel world of floating point numbers. Vertices that are constructed (from ST_Intersection, or ST_LocateBetween, or anything else) may not actually intersect the geometries that produced them.
    – dbaston
    Commented May 4, 2016 at 15:17
  • @dbaston So how can I snap the constructed gometry to match the intial line then?
    – mofoyoda
    Commented May 4, 2016 at 15:19
  • 2
    read up on ST_Snap. WITH mylines AS ( SELECT st_geomfromtext('LINESTRING(4377625.44 7399187.54,4377633.47 7398278.09)') as ro_geom, st_geomfromtext('LINESTRING(4377625.44 7399187.54,4377630.67957507 7398594.12386704)') as lk_geom ) select 1 as gid, st_astext(st_difference(st_snap(ro_geom, lk_geom, 0.001), st_snap(lk_geom, ro_geom, 0.001))) as diff_geom from mylines; Commented May 4, 2016 at 15:28

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge you have read our privacy policy.

Not the answer you're looking for? Browse other questions tagged or ask your own question.