[postgis-users] PostGIS transforms st_dwithin to && on geography column without preserving the custom spatial reference system

87 views
Skip to first unread message

Michael Loveridge

unread,
Jun 1, 2021, 3:18:40 PM6/1/21
to postgi...@lists.osgeo.org
When I use st_dwithin on a geography column with a GIST index, the
explain plan shows that the narrowing condition on the index scan is
changed into a && bounding box condition. With this different condition,
the result set is missing rows that are included if a st_dwithin is
performed using a full table scan when there is no index.

The condition in the explain plan without an index is:
st_dwithin(
geog,
'010100002010A4000000000000006066400000000000000000'::geography,
'3'::double precision,
false
)

But the explain plan indicates that it is changed into the following
when the column has a GIST index:
geog && _st_expand(
'010100002010A4000000000000006066400000000000000000'::geography,
'3'::double precision
)


The geography column has a custom spatial reference system that is not
being preserved in the transformation into a && condition. The custom
spatial reference system is a sphere where the radius is such that one
degree is equal to one meter of great-circle distance. So that a point
of (0 0) is one distance away from (0 1). But it seems that when using
st_dwithin, the search is being done in srid 4326 instead of my custom one.

The setup code is below:

-- Insert the custom spatial reference system into the db.
insert into spatial_ref_sys values (42000, 'customsrs', 1,
'GEOGCS[
"Normal Sphere (r=57.2957795)",
DATUM["unknown",
SPHEROID["Sphere",57.29577951308232087679,0]
],
PRIMEM["Greenwich",0],
CS[ellipsoidal,2],
AXIS["latitude",north],
AXIS["longitude",east],
UNIT["degree",0.0174532925199433]
]', '+proj=longlat +ellps=sphere +R=57.29577951308232087679
+no_defs');

-- Create a new table with a geography column in the new spatial
reference system.
CREATE TABLE geographytest(gid serial PRIMARY KEY, geog geography(POINT,
42000));

-- Insert some data around the dateline
insert into geographytest (gid, geog) values
(1, 'srid=42000;POINT(179 0)'),
(2, 'srid=42000;POINT(178 0)'),
(3, 'srid=42000;POINT(-179 0)'),
(4, 'srid=42000;POINT(-179 90)'),
(5, 'srid=42000;POINT(0 0)');


-- Select all points within a distance of 3 from POINT(179 0).
-- The expected 3 points are returned, with st_distances of 0, 1 and 2.
select
gid,
st_distance(
geog,
'srid=42000;POINT(179 0)',
false
),
st_dwithin(
geog,
'srid=42000;POINT(179 0)',
3,
false
)
from
geographytest
where
st_dwithin(
geog,
st_geogfromtext('srid=42000;POINT(179 0)'),
3,
false
);

-- Create a GIST index on our geography column
CREATE INDEX geographytestindex ON geographytest USING gist (geog);
VACUUM analyze geographytest (geog);

-- Now select again using the same query.
-- Now only one result is returned, the row that has POINT(179 0).
-- The explain plan indicates that the index scan is using
-- Index Cond: (geog &&
_st_expand('010100002010A4000000000000006066400000000000000000'::geography,
'3'::double precision))
-- when performing the select.
select
gid,
st_distance(
geog,
'srid=42000;POINT(179 0)',
false
),
st_dwithin(
geog,
'srid=42000;POINT(179 0)',
3,
false
)
from
geographytest
where
st_dwithin(
geog,
st_geogfromtext('srid=42000;POINT(179 0)'),
3,
false
);

Increasing the st_dwithin distance to 230000 instead of 3 returns all
three of the expected rows, since that is their distance in srid 4326.

How can I get PostGIS to use my custom spatial reference system in the
query when an index is present on the column?
_______________________________________________
postgis-users mailing list
postgi...@lists.osgeo.org
https://lists.osgeo.org/mailman/listinfo/postgis-users
Reply all
Reply to author
Forward
0 new messages