[postgis-users] Problems with ST_Within in an island of a multipolygon
Larry Reeder
lnreeder at gmail.com
Fri Mar 18 07:43:02 PDT 2011
I have a layer called "foo" of about 5000 polygons. One of the
polygons is a multipolygon with island shells. ST_Within isn't
working as expected for points within that island shell.
Here's my table:
foo_spatial=# \d foo
Table "public.foo"
Column | Type | Modifiers
----------+----------+---------------------------------------------------------
gid | integer | not null default nextval('foo_gid_seq'::regclass)
the_geom | geometry |
Indexes:
"foo_pkey" PRIMARY KEY, btree (gid)
"foo_the_geom_key" gist (the_geom)
Check constraints:
"enforce_dims_the_geom" CHECK (st_ndims(the_geom) = 2)
"enforce_geotype_the_geom" CHECK (geometrytype(the_geom) =
'MULTIPOLYGON'::text OR the_geom IS NULL)
"enforce_srid_the_geom" CHECK (st_srid(the_geom) = (-1))
ST_Within returns true when asked if a point inside the island belongs
to the polygon that contains the island:
foo_spatial=# select ST_Within(ST_GeomFromText('POINT (-115.04252
36.05137)', -1), the_geom) from foo where gid = 4933;
st_within
-----------
t
(1 row)
However, when I ask if any polygons contain this point, I don't get any back:
foo_spatial=# SELECT gid from foo where
ST_Within(ST_GeomFromText('POINT (-115.04252 36.05137)', -1),
the_geom);
gid
-----
(0 rows)
Does this make sense? I know PostGIS can handle the general problem
of ST_Within for island sub-polygons because when I extract this
multipolygon into a separate table by itself, the queries work as
expected. There seems to be something happening with other polygons
in the layer that is breaking this. Initially I did have some bad
polygons in the layer, but I used the cool script here:
http://linfiniti.com/2010/07/cleaning-geometries-inside-postgis/.
That fixed my bad polygons, and now all polygons are valid:
foo_spatial=# select count(*) from foo where ST_IsValid(the_geom) = 'f';
count
-------
0
(1 row)
But still ST_Within isn't working as expected. Here's my version info:
foo_spatial=# select version();
version
---------------------------------------------------------------------------------------------------------------------
PostgreSQL 8.4.7 on x86_64-pc-linux-gnu, compiled by GCC gcc-4.4.real
(Ubuntu/Linaro 4.4.4-14ubuntu5) 4.4.5, 64-bit
(1 row)
foo_spatial=# SELECT postgis_full_version();
postgis_full_version
-------------------------------------------------------------------------------------------------------
POSTGIS="1.5.1" GEOS="3.2.0-CAPI-1.6.0" PROJ="Rel. 4.7.1, 23
September 2009" LIBXML="2.7.6" USE_STATS
(1 row)
dash_spatial=# select count(*) from foo where ST_IsValid(the_geom) = 'f';
count
-------
0
(1 row)
Anyone have any ideas?
Thanks............ Larry
More information about the postgis-users
mailing list