1

When executing an ST_Intersection on 2 lines, the resulting point is offset. In the attached image, the offset is 0.000027m. I have tried ST_Snap to force it, but this has not changed the points position.

red line is the extended geom, black line is the original geom, and the point is the result of ST_Intersection of the two

    WITH line1 AS (SELECT ST_GeomFromText('LINESTRING(470874.945140537 6000126.53834916,470825.026548551 6000129.39039651,470813.970641131 6000131.67770864,470798.339625488 6000130.91512638,470770.509325729 6000128.62757298,470716.373201037 6000122.90909733,470592.358279565 6000112.18733429,470584.082725689 6000111.4719529,470471.616772575 6000098.31957132,470442.261380633 6000108.99375658,470402.993405271 6000132.6307069,470339.326451592 6000140.06487686,470261.553382752 6000143.87732586,470226.86058284 6000141.9710548,470158.999734862 6000139.4932543,470075.507824541 6000123.86262835,470042.340010929 6000112.04367615,470005.359698272 6000095.84140767,469987.441344206 6000088.9789069,469981.323382904 6000088.34605401,469976.38532627 6000087.83546738,469928.191247079 6000074.51890886)') as the_geom), 
     line2 AS (SELECT ST_GeomFromText('LINESTRING(470835.701330449 6000121.00308247,470835.612988186 6000125.16357655,470834.757333107 6000165.42524755,470834.681958091 6000168.96974537,470841.821713877 6000270.84345442,470838.892222296 6000668.7623786,470838.348652514 6000835.0334808)') as the_geom)
SELECT st_intersects(ST_Intersection(line1.the_geom, line2.the_geom),line1.the_geom) 
FROM line1
JOIN line2 ON true
1
  • That's floating-point math for you (both in the computation and the display). Ask for the distance of the point from the lines.
    – CL.
    Apr 23, 2020 at 10:38

1 Answer 1

1

After conferring with users of postgis-users mailing list, specifically Paul Ramsey and Sandro Santilli, this is a known issue with PostGIS. Sandro metioned:

  1. Take a grid of square cells.
  2. Draw line from upper-left to lower-right corner of a cell.
  3. Draw another line from upper-right to lower-left corner.
  4. Verify the 2 lines intersect
  5. Draw the intersection point WARNING: you can only draw at grid intersection points!

The above is what happens on your computer too. Only the grid has cells of variable size depending on their distance from the 0,0 coordinate.

Paul provided a solution to this problem:

WITH lines AS (
SELECT
'LINESTRING(470874.945140537 6000126.53834916,470825.026548551 6000129.39039651,470813.970641131 6000131.67770864,470798.339625488 6000130.91512638,470770.509325729 6000128.62757298,470716.373201037 6000122.90909733,470592.358279565 6000112.18733429,470584.082725689 6000111.4719529,470471.616772575 6000098.31957132,470442.261380633 6000108.99375658,470402.993405271 6000132.6307069,470339.326451592 6000140.06487686,470261.553382752 6000143.87732586,470226.86058284 6000141.9710548,470158.999734862 6000139.4932543,470075.507824541 6000123.86262835,470042.340010929 6000112.04367615,470005.359698272 6000095.84140767,469987.441344206 6000088.9789069,469981.323382904 6000088.34605401,469976.38532627 6000087.83546738,469928.191247079 6000074.51890886)'::geometry AS l1,
'LINESTRING(470835.701330449 6000121.00308247,470835.612988186 6000125.16357655,470834.757333107 6000165.42524755,470834.681958091 6000168.96974537,470841.821713877 6000270.84345442,470838.892222296 6000668.7623786,470838.348652514 6000835.0334808)'::geometry AS l2
),
outputs AS (
        SELECT
                ST_Snap(l1, ST_Intersection(l1, l2), 0.00001) AS l1,
                ST_Snap(l2, ST_Intersection(l1, l2), 0.00001) AS l2,
                ST_Intersection(l1, l2) AS p
        FROM lines
)
SELECT ST_Intersects(p, l1) AS p_l1, ST_Intersects(p, l2) as p_l2
FROM outputs;

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.