Opened 11 years ago

Closed 11 years ago

#1792 closed defect (invalid)

[raster] ST_Intersect returns false positives

Reported by: febus982 Owned by: pramsey
Priority: high Milestone: PostGIS 2.0.1
Component: postgis Version: 2.0.x
Keywords: Cc:


Software and DB Structure: PostgreSQL 9.1.3 GEOS 3.3.2 GDAL 1.9.0

Table 1 (vector): Single column (geom) - About 100.000 points geometries exported with GRASS 6.4.2 (v.out.ogr command).

Table 2 (raster): Column 1 (rast): raster column Column 2 (filename): filename

Raster is geotiff file imported with: /usr/lib/postgresql/9.1/bin/raster2pgsql -s 4326 -a -Y -F -t 20x20

View 1:

SELECT v.geom, r.filename, st_value(r.rast, v.geom) AS value

FROM vector v, raster r

WHERE r.filename = 'file.tif'::text AND st_intersects(r.rast, 1, v.geom);

I got multiple points in this view values (some are null values) They are all aligned on same X and/or same Y.

If i materialize this view in a table using: CREATE TABLE table AS SELECT…… (same command) I got a bunch of notices in console: st_intersects NOTICE: Attempting to get pixel value with out of range raster coordinates:

View 2:

SELECT v.geom, r.filename, st_value(r.rast, v.geom) AS value

FROM geoserver_layers.z1 v, raster.raster r

WHERE r.filename = '2012042611_me.tif'::text AND st_intersects(r.rast, 1, v.geom) AND st_value(r.rast, v.geom) IS NOT NULL;

I Tried filtering out NULL values by adding "st_value(r.rast, v.geom) IS NOT NULL" but I got some points missing instead of multiple points;

This happens in all forms of ST_Intersects st_intersects(v.geom, r.rast, 1) st_intersects(r.rast, 1, v.geom)

I Tried also changing tile size in raster2pgsql

Change History (2)

comment:1 by Bborie Park, 11 years ago

Summary: [raster] ST_Intersection returns false positives[raster] ST_Intersect returns false positives

Since you're using points, the point is being converted to a raster for ST_Intersects(raster, band, geometry). Hence the point becomes an area which may intersect with a raster that the point itself wouldn't.

As for ST_Intersects(geometry, raster, band), the raster is converted to polygons before running a geometry, geometry intersects test.

So by the looks of it, the problem is that in your query the ST_Intersects(raster, band, geometry) is operating by pixel (an "area") while ST_Value(raster, geometry) is essentially sticking a pin (your point) through the raster. The ST_Intersects returns true while ST_Value returns NULL.

You may want to use ST_Intersection(raster, band, geometry) instead of ST_Value.

comment:2 by febus982, 11 years ago

Resolution: invalid
Status: newclosed

I don't know the reason but i had to reboot and all seems working good now. Didn't have time to test with ST_Intersection.

Maybe some memory corruption :/

Note: See TracTickets for help on using tickets.