PostGIS 测试点是否与缩放的多边形相交有时会得到错误的答案

PostGIS testing if point intersects scaled polygon sometimes gets wrong answer

我需要一个算法来检查一个点是否落在一个小公差范围内的多边形内(意味着该点完全在多边形内或足够靠近边界)。为了实现这一点,我计划使用 st_scale 函数将多边形扩展一个因子,然后使用 st_within 检查点。请注意,我不能使用 st_buffer,因为 st_buffer 需要缓冲区的绝对半径,但我需要一个相对半径(例如,确定该点是在多边形内还是在已缩放的多边形内增加了 5%)。

为此,我编写了以下测试代码。我取一个多边形,将其放大 5%,然后测试原始多边形的质心是否落在放大后的多边形内:

SELECT ST_within(
  st_centroid(polygon),
  st_scale(
    polygon,
    'POINT(1.05 1.05)'::geometry,
    st_centroid(polygon)
));

但是,有时此函数会无缘无故地失败(这意味着多边形的质心不在按比例放大的多边形内),我不知道为什么。这感觉像是有问题,但我无法想象我在 postGIS 中发现了一个错误。

通过测试用例:

-- drawing a 1x1 box at the origin of the x-y plane
SELECT ST_Within(
  st_centroid('MULTIPOLYGON(((0 0, 1 0, 1 1, 0 1, 0 0)))'::geometry),
  st_scale(
    'MULTIPOLYGON(((0 0, 1 0, 1 1, 0 1, 0 0)))'::geometry,
    'POINT(1.05 1.05)'::geometry,
    st_centroid('MULTIPOLYGON(((0 0, 1 0, 1 1, 0 1, 0 0)))'::geometry))
);

-- box in the first quadrant at an origin of (10, 10)
SELECT ST_Within(
  st_centroid('MULTIPOLYGON(((10 10, 11 10, 11 11, 10 11, 10 10)))'::geometry),
  st_scale(
    'MULTIPOLYGON(((10 10, 11 10, 11 11, 10 11, 10 10)))'::geometry,
    'POINT(1.05 1.05)'::geometry,
    st_centroid('MULTIPOLYGON(((10 10, 11 10, 11 11, 10 11, 10 10)))'::geometry))
);

失败的测试用例:

-- 1x1 box with the bottom left corner at (98, 28)
SELECT ST_within(
  st_centroid('MULTIPOLYGON(((98 28,99 28,99 29,98 29,98 28)))'::geometry),
  st_scale(
    'MULTIPOLYGON(((98 28,99 28,99 29,98 29,98 28)))'::geometry,
    'POINT(1.05 1.05)'::geometry,
    st_centroid('MULTIPOLYGON(((98 28,99 28,99 29,98 29,98 28)))'::geometry))
);

不过,有趣的是,将缩放后的形状转换为 WKT,然后再转换回几何体确实有效。

-- 1x1 box with the bottom left corner at (98, 28)
SELECT ST_within(
  st_centroid('MULTIPOLYGON(((98 28,99 28,99 29,98 29,98 28)))'::geometry),
  st_geomfromewkt(st_asewkt(st_scale(
    'MULTIPOLYGON(((98 28,99 28,99 29,98 29,98 28)))'::geometry,
    'POINT(1.05 1.05)'::geometry,
    st_centroid('MULTIPOLYGON(((98 28,99 28,99 29,98 29,98 28)))'::geometry))))
);

所以,是的,这是怎么回事????我在 PostGIS 中发现了错误吗?

我认为您在旧版本的 postgis 中发现了一个错误。我使用 docker 和一些服务器来测试 postgres 和 postgis 的几个版本。使用 3.1.4 版代码可以正常工作。

过去我曾在一些极端情况下遇到过 postgis 给出了意想不到的结果。这些与 GEOS 有关,根据手册是实现 ST_Within 功能的库,但在这种情况下,GEOS 库在所有测试中都是相同的。

根据您转换为 ewkt 并返回的工作查询显示,这应该与已更正的 postgis 中的内部操作有关。我试过了,但没有找到与此相关的错误。

我 post 我做过的两项测试:

postgres 12, postgis 3.0.0

存在错误:

test=# select postgis_full_version();
                                                                                              postgis_full_version                                                                                              
----------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------
 POSTGIS="3.0.0 r17983" [EXTENSION] PGSQL="120" GEOS="3.7.1-CAPI-1.11.1 27a5e771" PROJ="Rel. 5.2.0, September 15th, 2018" LIBXML="2.9.4" LIBJSON="0.12.1" LIBPROTOBUF="1.3.1" WAGYU="0.4.3 (Internal)" TOPOLOGY
(1 row)

test=# SELECT ST_within(
  st_centroid('MULTIPOLYGON(((98 28,99 28,99 29,98 29,98 28)))'::geometry),
  st_scale(
    'MULTIPOLYGON(((98 28,99 28,99 29,98 29,98 28)))'::geometry,
    'POINT(1.05 1.05)'::geometry,
    st_centroid('MULTIPOLYGON(((98 28,99 28,99 29,98 29,98 28)))'::geometry))
);
 st_within 
-----------
 f
(1 row)

postgres 13, postgis 3.1.4

代码按预期运行:

test=# select postgis_full_version();
                                                                                              postgis_full_version                                                                                               
-----------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------
 POSTGIS="3.1.4 ded6c34" [EXTENSION] PGSQL="130" GEOS="3.7.1-CAPI-1.11.1 27a5e771" PROJ="Rel. 5.2.0, September 15th, 2018" LIBXML="2.9.4" LIBJSON="0.12.1" LIBPROTOBUF="1.3.1" WAGYU="0.5.0 (Internal)" TOPOLOGY
(1 row)

test=# SELECT ST_within(
  st_centroid('MULTIPOLYGON(((98 28,99 28,99 29,98 29,98 28)))'::geometry),
  st_scale(
    'MULTIPOLYGON(((98 28,99 28,99 29,98 29,98 28)))'::geometry,
    'POINT(1.05 1.05)'::geometry,
    st_centroid('MULTIPOLYGON(((98 28,99 28,99 29,98 29,98 28)))'::geometry))
);
 st_within 
-----------
 t
(1 row)