1

I have three points:

  • Point A
  • Point SW
  • Point NE

I want to know if Point A resides within the bounds determined by points SW (for south-west coordinates) and NE (for north-east coordinates).

I've so far tried both ST_Intersects and && to determine this, which works fine for small regions, but fails whenever the bounding area determined by SW and NE points covers more than half the globe. In this case it seems to instead use a different, smaller, bounding area using the same four points.

As an example, if I'd want to cover everything from the SW point of South Africa to the NE point of Ireland it would instead catch points where the the NW point is in Ireland and the the SE point is in Madagascar. That is, China would not considered to be within the bounds.

Here are some SQL queries to illustrate this:

select 
ST_GeographyFromText('SRID=4326;POINT(0 0)') && ST_GeographyFromText('SRID=4326;POLYGON((180 -80,  -80 -80, -80 80, 180 80, 180 -80))') as expectTrue_1, 
ST_Intersects(ST_GeographyFromText('SRID=4326;POINT(0 0)'), ST_GeographyFromText('SRID=4326;POLYGON((180 -80,  -80 -80, -80 80, 180 80, 180 -80))')) as expectTrue_2, 
ST_GeographyFromText('SRID=4326;POINT(0 0)') && ST_GeographyFromText('SRID=4326;POLYGON((180 -80,  180 80, -80 80, -80 -80, 180 -80))') as expectTrue_3, 
ST_Intersects(ST_GeographyFromText('SRID=4326;POINT(0 0)'), ST_GeographyFromText('SRID=4326;POLYGON((180 -80,  180 80, -80 80, -80 -80, 180 -80))')) as expectTrue_5, 
ST_GeographyFromText('SRID=4326;POINT(-90 0)') && ST_GeographyFromText('SRID=4326;POLYGON((180 -80,  -80 -80, -80 80, 180 80, 180 -80))') as expectFalse_1, 
ST_GeographyFromText('SRID=4326;POINT(0 0)') && ST_GeographyFromText('SRID=4326;POLYGON((1 -1,  -1 -1, -1 1, 1 1, 1 -1))') as expectTrue_5;

I'd want for this to return t t t t f t but it returns f f f f t t.

I've also tried ST_Intersects(ST_MakeBox2D()) as per the following SO post: https://stackoverflow.com/questions/24212355/query-by-coordinates-takes-too-long-options-to-optimize/24218861?iemail=1&noredirect=1#24218861

I'm using PostGIS 2.1.3 with the Geography data type.

What functions or approaches should I use to have PostGIS treat my more-than-half-globe-spanning polygons as more-than-half-globe-spanning ones instead of 'optimizing' them to be less?

I understand I could just split the bounds in, say, two parts, and check for each, but there must be a cleaner solution.

|improve this question
2

I'm guessing that you probably figure the connectors between your points are "straight lines" and would run parallel to the edges of a mercator map of the world. If you are using the geography type, that won't happen, the edges will be great circles and will run in directions you'd never guess.

http://boundlessgeo.com/2012/04/the-earth-is-not-flat-volume-2/

In any event, no edge in geography can be longer than 180d, since any two points define a great circle, the system has to decide what portion of the circle the points bound, and chooses the shortest one (the one less than 180d). If you want longer edges you have to add intermediate points that fall in the places you'd expect.

But remember, you aren't going to get "straight lines" out of this, you're going to get great circles.

If you want a straight line solution, work in geometry, and recognize that you're going to have to put in your own logic to handle dateline splitting.

(Dateline splitting is not hard to do, if your input box is not defined as "two points" but rather as a "lower left point" and "upper right point". Then you can easily test for cases in which the "lower left" point appears to be to the right of the "upper right point" and say "aha! this is actually a dateline crossing case!" and construct two query polygons with the dateline as the bounding edge.)

|improve this answer
  • That's a good point. I had assumed that the bounding box edges will run parallel to Mercator map edges (as I'm trying to retrieve locations that the user can "see" on a Google Map). Which Geometry SRID is best suited to store global data for this purpose? – John M Jun 18 '14 at 16:51
  • 4326 is the simplest one to use. You can always do casts out to geography for things like measuring distances and so on in your database if you need to for application logic. – Paul Ramsey Jun 18 '14 at 19:04
  • OK, I implemented it this way - by switching to geometry and handling the date-line spanning ranges as a special case. Thanks! – John M Jun 19 '14 at 16:03
1

Apply ST_Shift_Longitude() to both geometry variables (the point and the BBOX). Like this:

select 
ST_Shift_Longitude(ST_GeographyFromText('SRID=4326;POINT(0 0)')) && ST_Shift_Longitude(ST_GeographyFromText('SRID=4326;POLYGON((180 -80,  -80 -80, -80 80, 180 80, 180 -80))')) as expectTrue_1,
|improve this answer
  • 1
    Actually, ST_Shift_Longitude() won't work for you in this case... however, you can add 180deg to ALL of your longitudes and the && math would then work correctly: SELECT ST_GeomFromEWKT('SRID=4326;POINT(180 0)') && ST_GeomFromEWKT('SRID=4326;POLYGON((360 -80, 100 -80, 100 80, 360 80, 360 -80))') as expectTrue_1 – mtn.biker Jun 18 '14 at 13:28
  • Not critical here but if the aim is to get exact results it is better to use Intersects than && which gives false positives because it is only comparing the bounding boxes. From postgis.net/docs/reference.html: "&& — Returns TRUE if A's 2D bounding box intersects B's 2D bounding box" – user30184 Jun 18 '14 at 16:08
  • Actually instead of ST_GeographyFromText('SRID=4326;POINT(0 0)') (as Point A) I have to use a column of a largish table in a database. I could either store it already with 180 added or do something like ST_MakePoint(ST_Y(geometry(location)) + 180, ST_X(geometry(location))) - or is there an easier way? – John M Jun 18 '14 at 16:47
  • 1
    I'd advise you to store your point data as actual locations, then do your +180 when you perform the intersection with the bounding box. However, I'm not totally confident that you'd receive the benefit of the spatial index if doing it this way? I would test it both ways to see if your performance improves. – mtn.biker Jun 18 '14 at 17:32
  • 1
    Also, Paul brings up a good point to consider with the great circle bounding box while working with geography data type. I would, however, suggest that your user/client may not be so sophisticated and will concentrate on simply whether the point longitude value is less than (or greater than) the bounding box longitude value. You will need to consider what is your definition of the correct answer (actually correct vs user's perception)? P.S. SRID=4326 is a fine coordinate system for storing world-wide geometry data. – mtn.biker Jun 18 '14 at 17:39

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

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