1

I have two tables, 'protein' and 'data_aggreg', in a PostGIS database. I am simply trying to query the points that intersect (all of the points in 'data_aggreg' share exactly the same geometry as a portion of the points in 'protein' table).

This is typically a simple operation. However, this time PostGIS and QGIS return nothing and ArcGIS returns the correct selection of points. The data in both QGIS and ArcGIS are sourced directly from the same PostGIS database via their respective methods of connecting to external DBs (ie, no shapefiles, csv's, etc are being used for this analysis). Both tables are in projected coordinate system EPSG 26912 and both tables are single point geometry.

When I run ST_Intersects in PostGIS on the two geom columns nothing is returned. When I run the 'select by location' tool in QGIS nothing is returned (I've also tried spatial join with same effect).

Syntax for ST_Intersects, which returns nothing:

SELECT protein.protein, data_aggreg.pro_2016
from production.protein, scratch.data_aggreg
where ST_Intersects(geom26912, geom);

However when I run 'select by location' in ArcGIS I get the expected result. That is, the subset of points in the 'protein' table that share the same geometry with 'data_aggreg' are correctly selected. Also, when I change my PostGIS 'where' query from ST_Intersects to geom1 = geom2 I also get the expected result.

Syntax when using '=' in PostGIS, which returns the expected result:

SELECT protein.protein, data_aggreg.pro_2016
from production.protein, scratch.data_aggreg
where geom26912 = geom;

It is my understanding that if any part of features share any space then ST_Intersects will return 'true'. So it follows that if geom1 = geom2 then ST_Intersects should return 'true'. But this is not the case.

So why the difference in results within PostGIS when using ST_Intersects vs '=', and why the difference among PostGIS when using ST_Intersects, QGIS 'select by location', and ArcGIS 'select by location'?

9
  • It may be helpful to post the syntax you used for postgres.
    – Aaron
    Nov 16, 2016 at 18:15
  • Both SELECT ST_Intersects('POINT(0 0)'::geometry, 'POINT(0 0)'::geometry); and SELECT 'POINT(0 0)'::geometry = 'POINT(0 0)'::geometry; return true in my tests ...
    – underdark
    Nov 16, 2016 at 18:18
  • Shouldn't you still use Equals when comparing geometries instead of "="?
    – user30184
    Nov 16, 2016 at 18:34
  • What if you add ST_SnapToGrid with a rather small tolerance to ST_Intersects?
    – user30184
    Nov 16, 2016 at 18:37
  • 1
    I have a feeling that the problem is in the main SQL but I must think a bit further.
    – user30184
    Nov 16, 2016 at 18:47

1 Answer 1

0

= in PostGIS means bounding box equality. There is long history of that and people have debated changing it. But that's where that stands.

So for anything but points (even points you can't quite trust because of box rounding), it does not mean what you think it means.

As someone mentioned use something like ST_Equals (note this doesn't care about the orientation geometries), or ST_OrderingEquals.

I think ST_OrderingEquals is more akin to what people think of when they think of equality in the database sense.

2
  • According to the documentation, this is what = means Returns TRUE if the coordinates and coordinate order geometry/geography A are the same as the coordinates and coordinate order of geometry/geography B. Coordinates are floats, so this is very rarely a useful or even accurate comparison.
    – jbalk
    Sep 13, 2022 at 0:01
  • I'll change my downvote to an upvote if you fix that. Possibly just outdated?
    – jbalk
    Sep 13, 2022 at 0:02

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.