polygon 'contains' and other operations on geometries not supported

Question

I'm quite new to postgis and rgeo. I suspect I may be tackling things in the wrong way but I have been a bit surprised to find out a few operations, in particular contains & within, aren't possible on spherical based objects.

I have a bunch of geographically distributed objects which I would like to group together based on things like postcode. For each of of these groupings I have a boundary and I would like to check if an object is inside that boundary and also check if one boundary is inside another. I'm using rails and this is the migration used to set up my collection model

class CreateGeoCollectionDefinition < ActiveRecord::Migration[5.0]
  def change
     create_table :geo_collection_definitions do |t|
       t.string :name
       t.string :geo_place_id
       t.string :geo_place_types, array: true, default: []
       t.st_polygon :boundary, geographic: true
       t.jsonb :boundary_json
       t.st_point :latlng, geographic: true
     end
   end
end

Currently the boundary is coming from a google reverse geocode lookup. The NorthEast and SouthWest bounding box coordinates are passed into this method on object creation

GEO_FACTORY = RGeo::Geographic.spherical_factory(srid: 4326)

def self.createBoundary(pointOne, pointTwo)
  point1 = GEO_FACTORY.point(pointOne['lat'], pointOne['lng'])
  point2 = GEO_FACTORY.point(pointTwo['lat'], pointTwo['lng'])
  boundingBox = RGeo::Cartesian::BoundingBox.create_from_points(point1,    point2).to_geometry
  boundingBox
end

I have written a few of specs to check that everything is behaving in the way that I expect. The simple distance based tests all pass as expected but there are issues with the ones used to test the boundary functionality. I run into the following

 # ------------------
 # --- Caused by: ---
 # PG::UndefinedFunction:
 #   ERROR:  function st_contains(geography, geography) does not exist
 #   LINE 1: ...COUNT(*) FROM "geo_collection_definitions" WHERE (ST_Contain...
 #                                                                ^
 #   HINT:  No function matches the given name and argument types. You might need to add explicit type casts.

When I try and run a query like (I know its a bit of a silly one)

 GeoCollectionDefinition.where("ST_Contains(boundary, boundary)")

or if I try using the rgeo object directly

it "should be possible to test is points belong in GeoCollection.boundary" do
  factory = RGeo::Geographic.spherical_factory(srid: 4326)
  externalPoint = factory.point(EmptyGeocodeLatLag['lng'], EmptyGeocodeLatLag['lat'])
  expect(someplace_def.boundary.contains?(someplace_def.latlng)).to be_truthy
  expect(someplace_def.boundary.contains?(externalPoint)).to be_falsy
end

I get

  RGeo::Error::UnsupportedOperation:
   Method Geometry#contains? not defined.

Digging around I found this rgeo issue and other evidence that these operations are just not supported for spherical factory based objects

I'm just wondering;

1. Is this definitely the case
2. Modelling the boundary and the location of objects in the way described in my migration seems to make sense to me but I'm guessing I'm wrong about that?
3. How should I go about finding out if a bunch of points are inside a polygon using postgis, rgeo and the active record adaptor?
4. Is it possible to check if one polygon is inside another?

EDIT:

I followed up on tilt's suggestion and ran some queries directly on my db. First, just to verify my set up

SELECT PostGIS_full_version();

 NOTICE:  Function postgis_topology_scripts_installed() not found. Is topology support enabled and topology.sql installed?                                  postgis_full_version

 POSTGIS="2.1.7 r13414" GEOS="3.5.0-CAPI-1.9.0 r4084" PROJ="Rel.   4.9.2, 08 September 2015" GDAL="GDAL 1.11.5, released 2016/07/01" LIBXML="2.9.2" LIBJSON="UNKNOWN" RASTER

I also ran some queries;

  SELECT name  FROM geo_collection_definitions WHERE st_contains(latlng, boundary);
 ERROR:  function st_contains(geography, geography) does not exist
 LINE 1: ...ELECT name  FROM geo_collection_definitions WHERE   st_contain...

So guess contains just won't work on geographies. I dug around and I found st_covers

  SELECT name  FROM geo_collection_definitions WHERE st_covers(boundary, latlng);
  name
  ------
  (0 rows)

 SELECT name  FROM geo_collection_definitions WHERE st_covers(boundary, ST_GeomFromText('POINT(12.9549709 55.5563043)', 4326));
  name
  ------
  (0 rows)

This is really surprising as latlng is the centre point of boundary. I'm pretty confused and sure I'm doing something pretty stupid. Any help would be greatly appreciated

Answer 1

I recommend using ST_DWithin, which is well supported for PostGIS' geography type. For the radius parameter, you can use 0 or maybe 10 (i.e. if your data has an accuracy of 10 m).

There are not any plans to make ST_Contains or ST_Within available for geography type.

Answer 2

Mike and Tilt's responses helped me get to the bottom of this and it's a bit of an embarrassing one but just in case any of this stuff might help someone else down the line...

ST_DWithin is definitely the way to go for my particular problem and it was great to get a shout out about this. I changed some of spec queries to things like this

  collectionDefs = GeoCollectionDefinition.where("ST_DWithin(boundary, '#{someplace.latlng}', 10)")
  expect(collectionDefs.count).to eq(2)
  collectionDefs = GeoCollectionDefinition.where("ST_DWithin(boundary, 'Point(0.0 0.0)', 10)")
  expect(collectionDefs.count).to eq(0) 

The query ran fine which was great but I still wasn't getting the results I was expecting. This is the embarrassing bit

Following Tilt's advice I printed out the polygon

POLYGON((55.4965351 12.8894595,55.6445967 12.8894595,55.6445967 13.151087,55.4965351 13.151087,55.4965351 12.8894595)) and the centre point is POINT(13.0108705 55.5790534)

and quickly spotted that coordinates where the wrong way round. The damage was being done in the createBoundary method I outlined above. I was creating the boundary definition points and passing lat value in place of the lng and vice versa. I updated the code to the following

 def self.createBoundary(pointOne, pointTwo)
    point1 = GEO_FACTORY.point(pointOne['lng'], pointOne['lat'])
    point2 = GEO_FACTORY.point(pointTwo['lng'], pointTwo['lat'])
    boundingBox =    RGeo::Cartesian::BoundingBox.create_from_points(point1, point2).to_geometry
   boundingBox
 end

I'll hang my head in shame :)

Thanks a million Tilt and Mike. I'm not sure I would have found the ST_DWithin and the pro postgist debugging tips definitely saved me many hours of frustration