Effectively forcing the planner to do the thing you want might help. In this case, sub-setting the polygon table prior to executing the spatial join with the points table. You might be able to outwit the planner using "WITH" syntax:
WITH polys AS (
SELECT * FROM area
WHERE area.id in IN(28,29,30,31,32,33,25,26,27,18,19,20,21,12,13,14,15,16,17,34,35,1,2,3,4,5,6,22,23,24,7,8,9,10,11)
)
SELECT *
FROM polys JOIN data ON ST_Intersects(data.point, polys.polygon)
WHERE datasetid IN(3)
AND (readingdatetime BETWEEN '1948-01-01' AND '2012-11-19');
The trouble with trying to play these games is that you are coding into your statement the assumption "my polygon list will always be more selective than my other query portions". Which might not be true for all parameterizations of your query, or for all applications of a particular query over a heterogeneously distributed dataset.
But it might work.
UPDATE: This goes even further down the dangerous road of assuming you know the selectivity of your clauses beforehand, this time we also take the attribute selection on the point table out and do it separately before the spatial join:
WITH polys AS (
SELECT * FROM area
WHERE area.id in IN(28,29,30,31,32,33,25,26,27,18,19,20,21,12,13,14,15,16,17,34,35,1,2,3,4,5,6,22,23,24,7,8,9,10,11)
),
WITH points AS (
SELECT * FROM data
WHERE datasetid IN(3)
AND (readingdatetime BETWEEN '1948-01-01' AND '2012-11-19')
)
SELECT *
FROM polys JOIN points ON ST_Intersects(points, polys.polygon);
Your MultiPolygon is invalid. You can test it with ogrinfo:
ogrinfo -ro -dialect sqlite -sql "select IsValid(geometry) from
OGRGeoJSON" multipolygon.json
INFO: Open of `multipolygon.json'
using driver `GeoJSON' successful.
GEOS warning: Self-intersection at or near point -115.907104 39.162109000000001
Layer name: SELECT
Geometry: None
Feature Count: 1
Layer SRS WKT:
(unknown)
isvalid(geometry): Integer (0.0)
OGRFeature(SELECT):0
isvalid(geometry) (Integer) = 0
You can try to correct it in PostGIS with ST_MakeValid http://postgis.net/docs/ST_MakeValid.html. However, it is not simple
Step 1 Correct the faulty geometry:
create table valid as select st_makevalid(geometry) as geometry from multipolygon;
Step 2 Test:
select * from valid where ST_Intersects(ST_PointFromText('POINT( -116.024551 38.485773 )',4326), geometry);
Result: error
ERROR: Relate Operation called with a LWGEOMCOLLECTION type. This is unsupported.
HINT: Change argument 2: 'GEOMETRYCOLLECTION(POLYGON((-115.000793 38.677307,-115.000793 38.499878,-115....'
Pity that ST_Intersects does not handle geometrycollections. So you must get rid of the collection. First part of the collection is the polygon so you can take just that.
Step 3 Keep just the polygon part
update valid set geometry=ST_GeometryN(geometry,1);
Step 4 Test again
select * from valid where ST_Intersects(ST_PointFromText('POINT( -116.024551 38.485773 )',4326), geometry);
Result: Success!
Best Answer
This return true as it should
So my best quess is that st_covers has bug or feature that it does not work when polygon fills all quartes of globe.
EDIT:
Played with this (PostGIS 2.0)
So it splits polygon to its quarters and casts geometry to geography and then compares point to polygon. It has its own problems It Uses Geometry to split because split doesn't support grography, but if you split it from meridian and equator it shouldn't be problem because lines are "straight" ( i may be wrong , but i'm fairly sure about it )
This solution is nice because you can convert all your polygons to multipolygons, or if you really want you can do it run time.