0

Using Postgres 9.2

I have a strange issue. To simplify it:

I have some data with line points. The query in question is using st_intersects to determine if the line points overlap with a polygon. Both the line and the polygon are stored in a 3d representation, with the z-axis being 0. This is for geospatial data.

In this case, I have a line point where the start and end are the same value. Two records have seemingly the same value, the X,Y,Z components on the start and end points are the same. Comparing the two points using =~, they are equal. Using =, they are equal. Using st_equals, the result is false, but comparing the components that make up the lines, the values all seem to be equal, including comparing the binary representation by visual inspection.

When I do st_intersects(my_line, some_polygon), one record returns true, the other false, even though the value of the line for both records appear identical. I didn't create the original values, so I do not know how they were originally created. There is a vehicle associated with each record, and for whatever reason, one of the vehicles has this problem for several of its records.

If I change the function from st_intersects to the presumably more expensive st_3dintersects, they both return true as expected, and the problem goes away. The polygon being compared against is quite large, and this affects several records with different points, so its unlikely we're hitting some fringe rounding error of any kind. Using st_force2d doesn't work either.

Any ideas why I might be seeing the behavior that I'm seeing?

Here's the EWKT of the line, with the coordinates changed:

SRID=4326;LINESTRING(-85.6600021 30.7976979 0,-85.6600021 30.7976979 0)

Both records have this exact same value for ST_AsEWKT, and yet one of them returns false for st_intersects(my_line, the_poly) and the other returns true for st_intersects(my_line, the_poly). Even if I hard code the EWKT value I still see this discrepency:

ST_Intersects(
  ST_GeomFromEWKT('SRID=4326;LINESTRING(-85.6600021 30.7976979 0,-85.6600021 30.7976979 0)')),
  x.geom
)

It seems like it is always affecting the very first record in the result set and no other records. If I change everything else in my query, this always returns false for the first record and true for all the subsequent records.

Edit: More investigation, it appears that the linestring is not valid with both the start and end being the same value. Casting st_makevalid fixes it by making it a point. Apparently the invalid linestring is evaluated inconsistently.

7
  • Provide the ST_AsEWKT() for both of the things you're sending to ST_Equals. Dec 7, 2016 at 6:26
  • What do you mean when you say when I do st_intersects(my_line, some_polygon), one returns true, the other false What other? Dec 7, 2016 at 6:37
  • This is likely a better question for gis.stackexchange.com (you can flag it and ask it to be moved if you want) Dec 7, 2016 at 6:37
  • I did more research and found that the first record with this specific line string always returns false for st_intersects if its the first record returned, and true for all the subsequent records. This seems to indicate that it is a bug in gis? Dec 7, 2016 at 17:14
  • No idea what you're talking about no idea what the question is, no sample data. Dec 7, 2016 at 17:43

2 Answers 2

1

It's most likely there are sub-decimal differences with coordinates that you can only see with the WKB, which means there are small differences that you cannot see with the WKT formatting. Here is an example:

SELECT ST_AsEWKT(A) AS wkt_a, ST_AsEWKT(B) AS wkt_b,
    ST_AsEWKT(A) = ST_AsEWKT(B) AS wkt_are_equal,
    A::text = B::text AS wkb_are_equal,
    ST_Intersects(A, B), ST_Distance(A, B),
    ST_Distance(A, B) < 1e-12 AS pretty_much_intersect
FROM (
    SELECT
       '01010000A0E6100000A5B272793D6A55C07E96F8ED35CC3E400000000000000000'::geometry AS A,
       '01010000A0E6100000A5B272793D6A55C07F96F8ED35CC3E400000000000000000'::geometry AS B
) f;
-[ RECORD 1 ]---------+------------------------------------------
wkt_a                 | SRID=4326;POINT(-85.6600021 30.7976979 0)
wkt_b                 | SRID=4326;POINT(-85.6600021 30.7976979 0)
wkt_are_equal         | t
wkb_are_equal         | f
st_intersects         | f
st_distance           | 3.5527136788005e-015
pretty_much_intersect | t

So you can see that the WKT are equal, but the WKB are not. There is a tiny distance between the two, therefore ST_Intersects will return false, as these predicate functions require exact noding.

A more robust metric to find geometries that essentially intersect is shown by testing if the distance is within a small distance, as demonstrated by the last column. Another solution is to see ST_Snap.


Now just seeing the invalid geometries in the question, my answer is to not use invalid geometries!

Behaviour is reproduced here:

DROP TABLE IF EXISTS invalid;

CREATE TEMP TABLE invalid(id integer primary key, geom geometry);
INSERT INTO invalid(id, geom) VALUES
(1, 'LINESTRING(-85.6600021 30.7976979,-85.6600021 30.7976979)'),
(2, 'LINESTRING(-85.6600021 30.7976979,-85.6600021 30.7976979)'),
(3, 'LINESTRING(-85.6600021 30.7976979,-85.6600021 30.7976979)');

SELECT id, ST_Intersects(
  ST_GeomFromEWKT('LINESTRING(-85.6600021 30.7976979,-85.6600021 30.7976979)'), x.geom)
    FROM invalid x;

 id | st_intersects
----+---------------
  1 | f
  2 | t
  3 | t
(3 rows)
3
  • That's not the case because you can hard code the geometry and reproduce the issue by simply cross joining it with any table. The problem seems to be with how postgres handles an invalid line with both the start and end containing the same value. Dec 9, 2016 at 1:41
  • @JamesHutchison Ok, I see your edit about invalid geometries, and I can get PostGIS to show what you are trying to explain. My advice is to not use invalid geometries.
    – Mike T
    Dec 9, 2016 at 2:34
  • Yes, I don't have control over the data so I am using st_makevalid on the suspect geometry to fix it, which resolves the issue. It's a shame though that postgres has this behavior and does so silently. Dec 9, 2016 at 20:49
-2

I have no idea what you're talking about. Just provide data because you're not good at describing the problem.

\set linestring ST_GeomFromEWKT($$SRID=4326;LINESTRING(-85.6600021 30.7976979 0,-80.6600021 30.7976979 0)$$)
\set point ST_GeomFromEWKT($$SRID=4326;POINT(-85.6600021 30.7976979 0)$$)

SELECT ST_Intersects( :linestring, :point ) AS linestringPoint
  , ST_Intersects( :linestring, :linestring ) AS linestringLinestring

Your Answer

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

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