From 684515cbabcb37ff14a87b09f195de9c60e7bbd0 Mon Sep 17 00:00:00 2001 From: Greg Burd Date: Sat, 24 Oct 2026 14:32:11 -0400 Subject: [PATCH v11 1/8] Compute the distance between a point and a box exactly box_closept_point(), behind <->(point, box), <->(box, point) and ##(point, box), found the closest point by projecting the point onto each edge's line with the generic line-segment code. That code compares slopes and intersections with the EPSILON fuzz, so for a point near an edge's line it can pick the wrong edge or corner. The point (0.9999996, 2) is 1 from the box ((0,0),(1,1)), but the operator said 1.00000000000008, and (1e-300, -1e-300) is 1e-300 from ((2,2),(0,0)), but the operator said 1.41421356237e-300, the distance to the corner. The GiST and SP-GiST box opclasses compute the distance to a box their own way, exactly, and order a nearest-neighbor scan by that. So an index scan could return rows in an order the operator disagrees with: for the point above, a box 1.00000000000004 away came after the box the operator put at 1.00000000000008. A box's edges are parallel to the axes, so the closest point is the point clamped to the box, and the distance follows from that with no line arithmetic. Do that, and have GiST's and SP-GiST's distance to a box (computeDistance() for non-leaf entries and point_box_distance()) call the same function, so the index orders by the operator's own distance. SP-GiST's version now errors on overflow, as GiST's and the operator always did, and returns NaN when any coordinate is NaN. The geometry regression test changes in four rows, all corrections: the 1e-300 distance above, its closest point, and the closest points to (1e+300,Infinity) and (Infinity,1e+300), which are now the box's top right corner rather than its top left. This is independent of the rest of the series, which uses it to let the box opclasses return their distances. Author: Greg Burd Discussion: https://postgr.es/m/8s5lT8erXzBMugXJQ6Wginbp_gc2B4hmCYhF0Q0GpnuK87eCAOcKs5K31AWCwraCNFIQt42wwkow_MPPubHg485MWv-zDBY4NL_JX9yJEEg=@burd.me --- src/backend/access/gist/gistproc.c | 59 +------------ src/backend/access/spgist/spgproc.c | 28 ------ src/backend/utils/adt/geo_ops.c | 78 ++++++++-------- src/include/utils/geo_decls.h | 2 + src/test/regress/expected/box.out | 35 ++++++++ src/test/regress/expected/geometry.out | 118 ++++++++++++------------- src/test/regress/sql/box.sql | 24 +++++ 7 files changed, 158 insertions(+), 186 deletions(-) diff --git a/src/backend/access/gist/gistproc.c b/src/backend/access/gist/gistproc.c index f1044f49d6c..34772ffd9c0 100644 --- a/src/backend/access/gist/gistproc.c +++ b/src/backend/access/gist/gistproc.c @@ -1223,67 +1223,14 @@ gist_point_fetch(PG_FUNCTION_ARGS) static float8 computeDistance(bool isLeaf, BOX *box, Point *point) { - float8 result = 0.0; - if (isLeaf) { /* simple point to point distance */ - result = point_point_distance(point, &box->low); - } - else if (point->x <= box->high.x && point->x >= box->low.x && - point->y <= box->high.y && point->y >= box->low.y) - { - /* point inside the box */ - result = 0.0; - } - else if (point->x <= box->high.x && point->x >= box->low.x) - { - /* point is over or below box */ - Assert(box->low.y <= box->high.y); - if (point->y > box->high.y) - result = float8_mi(point->y, box->high.y); - else if (point->y < box->low.y) - result = float8_mi(box->low.y, point->y); - else - elog(ERROR, "inconsistent point values"); - } - else if (point->y <= box->high.y && point->y >= box->low.y) - { - /* point is to left or right of box */ - Assert(box->low.x <= box->high.x); - if (point->x > box->high.x) - result = float8_mi(point->x, box->high.x); - else if (point->x < box->low.x) - result = float8_mi(box->low.x, point->x); - else - elog(ERROR, "inconsistent point values"); - } - else - { - /* closest point will be a vertex */ - Point p; - float8 subresult; - - result = point_point_distance(point, &box->low); - - subresult = point_point_distance(point, &box->high); - if (result > subresult) - result = subresult; - - p.x = box->low.x; - p.y = box->high.y; - subresult = point_point_distance(point, &p); - if (result > subresult) - result = subresult; - - p.x = box->high.x; - p.y = box->low.y; - subresult = point_point_distance(point, &p); - if (result > subresult) - result = subresult; + return point_point_distance(point, &box->low); } - return result; + /* the operator's own distance, so that box_ops' distances are exact */ + return point_box_distance(point, box); } static bool diff --git a/src/backend/access/spgist/spgproc.c b/src/backend/access/spgist/spgproc.c index 78873f23855..01b2091d310 100644 --- a/src/backend/access/spgist/spgproc.c +++ b/src/backend/access/spgist/spgproc.c @@ -26,34 +26,6 @@ DatumGetFloat8(DirectFunctionCall2(point_distance, \ PointPGetDatum(p1), PointPGetDatum(p2))) -/* Point-box distance in the assumption that box is aligned by axis */ -static double -point_box_distance(Point *point, BOX *box) -{ - double dx, - dy; - - if (isnan(point->x) || isnan(box->low.x) || - isnan(point->y) || isnan(box->low.y)) - return get_float8_nan(); - - if (point->x < box->low.x) - dx = box->low.x - point->x; - else if (point->x > box->high.x) - dx = point->x - box->high.x; - else - dx = 0.0; - - if (point->y < box->low.y) - dy = box->low.y - point->y; - else if (point->y > box->high.y) - dy = point->y - box->high.y; - else - dy = 0.0; - - return hypot(dx, dy); -} - /* * Returns distances from given key to array of ordering scan keys. Leaf key * is expected to be point, non-leaf key is expected to be box. Scan key diff --git a/src/backend/utils/adt/geo_ops.c b/src/backend/utils/adt/geo_ops.c index 73324b91fe5..5d99579953f 100644 --- a/src/backend/utils/adt/geo_ops.c +++ b/src/backend/utils/adt/geo_ops.c @@ -2951,64 +2951,56 @@ close_lseg(PG_FUNCTION_ARGS) /* - * Closest point on or in box to specified point. + * The closest point of 'box' to 'pt': 'pt' itself clamped to the box. * - * If *result is not NULL, set it to the closest point on the box to the - * given point, and return the distance of the two points. + * Return the distance between them. If *result is not NULL, set it to the + * closest point. A box's edges are parallel to the axes, so this needs no + * line arithmetic, and the result is exact up to hypot()'s rounding. The + * GiST and SP-GiST support functions use it through point_box_distance(), + * so the distances an index scan orders by are the operator's. */ static float8 box_closept_point(Point *result, BOX *box, Point *pt) { - float8 dist, - d; - Point point, - closept; - LSEG lseg; + Point closept; - if (box_contain_point(box, pt)) + if (isnan(pt->x) || isnan(pt->y) || + isnan(box->low.x) || isnan(box->low.y) || + isnan(box->high.x) || isnan(box->high.y)) { if (result != NULL) *result = *pt; - - return 0.0; + return get_float8_nan(); } - /* pairwise check lseg distances */ - point.x = box->low.x; - point.y = box->high.y; - statlseg_construct(&lseg, &box->low, &point); - dist = lseg_closept_point(result, &lseg, pt); + if (pt->x < box->low.x) + closept.x = box->low.x; + else if (pt->x > box->high.x) + closept.x = box->high.x; + else + closept.x = pt->x; - statlseg_construct(&lseg, &box->high, &point); - d = lseg_closept_point(&closept, &lseg, pt); - if (float8_lt(d, dist)) - { - dist = d; - if (result != NULL) - *result = closept; - } + if (pt->y < box->low.y) + closept.y = box->low.y; + else if (pt->y > box->high.y) + closept.y = box->high.y; + else + closept.y = pt->y; - point.x = box->high.x; - point.y = box->low.y; - statlseg_construct(&lseg, &box->low, &point); - d = lseg_closept_point(&closept, &lseg, pt); - if (float8_lt(d, dist)) - { - dist = d; - if (result != NULL) - *result = closept; - } + if (result != NULL) + *result = closept; - statlseg_construct(&lseg, &box->high, &point); - d = lseg_closept_point(&closept, &lseg, pt); - if (float8_lt(d, dist)) - { - dist = d; - if (result != NULL) - *result = closept; - } + return point_dt(&closept, pt, NULL); +} - return dist; +/* + * Distance from 'pt' to the nearest point of 'box', for the index support + * functions; see box_closept_point(). + */ +float8 +point_box_distance(Point *pt, BOX *box) +{ + return box_closept_point(NULL, box, pt); } Datum diff --git a/src/include/utils/geo_decls.h b/src/include/utils/geo_decls.h index a7014cded51..5f99e1e2007 100644 --- a/src/include/utils/geo_decls.h +++ b/src/include/utils/geo_decls.h @@ -274,4 +274,6 @@ CirclePGetDatum(const CIRCLE *X) #define PG_GETARG_CIRCLE_P(n) DatumGetCircleP(PG_GETARG_DATUM(n)) #define PG_RETURN_CIRCLE_P(x) return CirclePGetDatum(x) +extern float8 point_box_distance(Point *pt, BOX *box); + #endif /* GEO_DECLS_H */ diff --git a/src/test/regress/expected/box.out b/src/test/regress/expected/box.out index 10760870ce7..9b01c19745c 100644 --- a/src/test/regress/expected/box.out +++ b/src/test/regress/expected/box.out @@ -641,6 +641,41 @@ WHERE seq.id IS NULL OR idx.id IS NULL; RESET enable_seqscan; RESET enable_indexscan; RESET enable_bitmapscan; +-- The distance from a point to a box is exact, so a kNN index scan, which +-- orders by the index's own distance, returns rows in the operator's order. +-- For the point (0.9999996,2), box 1 is 1 away and box 2 is 1.00000000000004 +-- away. (The operator used to give box 1 1.00000000000008, out of order.) +-- Both distances are exact on any platform: the point is above both boxes, +-- so each is hypot(0, dy) = |dy|. +CREATE TEMP TABLE box_knn (id int, b box); +INSERT INTO box_knn VALUES (1, '((0,0),(1,1))'), + (2, box(point(0.4999996, 0), point(1.4999996, 0.99999999999996))); +SET extra_float_digits = 0; +SET enable_seqscan = OFF; +SET enable_bitmapscan = OFF; +CREATE INDEX box_knn_gist ON box_knn USING gist (b); +SELECT id, b <-> point(0.9999996, 2) AS dist FROM box_knn + ORDER BY b <-> point(0.9999996, 2); + id | dist +----+------------------ + 1 | 1 + 2 | 1.00000000000004 +(2 rows) + +DROP INDEX box_knn_gist; +CREATE INDEX box_knn_spgist ON box_knn USING spgist (b); +SELECT id, b <-> point(0.9999996, 2) AS dist FROM box_knn + ORDER BY b <-> point(0.9999996, 2); + id | dist +----+------------------ + 1 | 1 + 2 | 1.00000000000004 +(2 rows) + +RESET extra_float_digits; +RESET enable_seqscan; +RESET enable_bitmapscan; +DROP TABLE box_knn; -- test non-error-throwing API for some core types SELECT pg_input_is_valid('200', 'box'); pg_input_is_valid diff --git a/src/test/regress/expected/geometry.out b/src/test/regress/expected/geometry.out index 1d168b21cbc..a0a3a204f85 100644 --- a/src/test/regress/expected/geometry.out +++ b/src/test/regress/expected/geometry.out @@ -682,58 +682,58 @@ SELECT p.f1, l.s, p.f1 <-> l.s AS dist_ps, l.s <-> p.f1 AS dist_sp FROM POINT_TB -- Distance to box SELECT p.f1, b.f1, p.f1 <-> b.f1 AS dist_pb, b.f1 <-> p.f1 AS dist_bp FROM POINT_TBL p, BOX_TBL b; - f1 | f1 | dist_pb | dist_bp --------------------+---------------------+--------------------+-------------------- - (0,0) | (2,2),(0,0) | 0 | 0 - (0,0) | (3,3),(1,1) | 1.41421356237 | 1.41421356237 - (0,0) | (-2,2),(-8,-10) | 2 | 2 - (0,0) | (2.5,3.5),(2.5,2.5) | 3.53553390593 | 3.53553390593 - (0,0) | (3,3),(3,3) | 4.24264068712 | 4.24264068712 - (-10,0) | (2,2),(0,0) | 10 | 10 - (-10,0) | (3,3),(1,1) | 11.0453610172 | 11.0453610172 - (-10,0) | (-2,2),(-8,-10) | 2 | 2 - (-10,0) | (2.5,3.5),(2.5,2.5) | 12.747548784 | 12.747548784 - (-10,0) | (3,3),(3,3) | 13.3416640641 | 13.3416640641 - (-3,4) | (2,2),(0,0) | 3.60555127546 | 3.60555127546 - (-3,4) | (3,3),(1,1) | 4.12310562562 | 4.12310562562 - (-3,4) | (-2,2),(-8,-10) | 2 | 2 - (-3,4) | (2.5,3.5),(2.5,2.5) | 5.52268050859 | 5.52268050859 - (-3,4) | (3,3),(3,3) | 6.0827625303 | 6.0827625303 - (5.1,34.5) | (2,2),(0,0) | 32.6475113906 | 32.6475113906 - (5.1,34.5) | (3,3),(1,1) | 31.5699223946 | 31.5699223946 - (5.1,34.5) | (-2,2),(-8,-10) | 33.2664996656 | 33.2664996656 - (5.1,34.5) | (2.5,3.5),(2.5,2.5) | 31.108841187 | 31.108841187 - (5.1,34.5) | (3,3),(3,3) | 31.5699223946 | 31.5699223946 - (-5,-12) | (2,2),(0,0) | 13 | 13 - (-5,-12) | (3,3),(1,1) | 14.3178210633 | 14.3178210633 - (-5,-12) | (-2,2),(-8,-10) | 2 | 2 - (-5,-12) | (2.5,3.5),(2.5,2.5) | 16.3248277173 | 16.3248277173 - (-5,-12) | (3,3),(3,3) | 17 | 17 - (1e-300,-1e-300) | (2,2),(0,0) | 1.41421356237e-300 | 1.41421356237e-300 - (1e-300,-1e-300) | (3,3),(1,1) | 1.41421356237 | 1.41421356237 - (1e-300,-1e-300) | (-2,2),(-8,-10) | 2 | 2 - (1e-300,-1e-300) | (2.5,3.5),(2.5,2.5) | 3.53553390593 | 3.53553390593 - (1e-300,-1e-300) | (3,3),(3,3) | 4.24264068712 | 4.24264068712 - (1e+300,Infinity) | (2,2),(0,0) | Infinity | Infinity - (1e+300,Infinity) | (3,3),(1,1) | Infinity | Infinity - (1e+300,Infinity) | (-2,2),(-8,-10) | Infinity | Infinity - (1e+300,Infinity) | (2.5,3.5),(2.5,2.5) | Infinity | Infinity - (1e+300,Infinity) | (3,3),(3,3) | Infinity | Infinity - (Infinity,1e+300) | (2,2),(0,0) | Infinity | Infinity - (Infinity,1e+300) | (3,3),(1,1) | Infinity | Infinity - (Infinity,1e+300) | (-2,2),(-8,-10) | Infinity | Infinity - (Infinity,1e+300) | (2.5,3.5),(2.5,2.5) | Infinity | Infinity - (Infinity,1e+300) | (3,3),(3,3) | Infinity | Infinity - (NaN,NaN) | (2,2),(0,0) | NaN | NaN - (NaN,NaN) | (3,3),(1,1) | NaN | NaN - (NaN,NaN) | (-2,2),(-8,-10) | NaN | NaN - (NaN,NaN) | (2.5,3.5),(2.5,2.5) | NaN | NaN - (NaN,NaN) | (3,3),(3,3) | NaN | NaN - (10,10) | (2,2),(0,0) | 11.313708499 | 11.313708499 - (10,10) | (3,3),(1,1) | 9.89949493661 | 9.89949493661 - (10,10) | (-2,2),(-8,-10) | 14.4222051019 | 14.4222051019 - (10,10) | (2.5,3.5),(2.5,2.5) | 9.92471662064 | 9.92471662064 - (10,10) | (3,3),(3,3) | 9.89949493661 | 9.89949493661 + f1 | f1 | dist_pb | dist_bp +-------------------+---------------------+---------------+--------------- + (0,0) | (2,2),(0,0) | 0 | 0 + (0,0) | (3,3),(1,1) | 1.41421356237 | 1.41421356237 + (0,0) | (-2,2),(-8,-10) | 2 | 2 + (0,0) | (2.5,3.5),(2.5,2.5) | 3.53553390593 | 3.53553390593 + (0,0) | (3,3),(3,3) | 4.24264068712 | 4.24264068712 + (-10,0) | (2,2),(0,0) | 10 | 10 + (-10,0) | (3,3),(1,1) | 11.0453610172 | 11.0453610172 + (-10,0) | (-2,2),(-8,-10) | 2 | 2 + (-10,0) | (2.5,3.5),(2.5,2.5) | 12.747548784 | 12.747548784 + (-10,0) | (3,3),(3,3) | 13.3416640641 | 13.3416640641 + (-3,4) | (2,2),(0,0) | 3.60555127546 | 3.60555127546 + (-3,4) | (3,3),(1,1) | 4.12310562562 | 4.12310562562 + (-3,4) | (-2,2),(-8,-10) | 2 | 2 + (-3,4) | (2.5,3.5),(2.5,2.5) | 5.52268050859 | 5.52268050859 + (-3,4) | (3,3),(3,3) | 6.0827625303 | 6.0827625303 + (5.1,34.5) | (2,2),(0,0) | 32.6475113906 | 32.6475113906 + (5.1,34.5) | (3,3),(1,1) | 31.5699223946 | 31.5699223946 + (5.1,34.5) | (-2,2),(-8,-10) | 33.2664996656 | 33.2664996656 + (5.1,34.5) | (2.5,3.5),(2.5,2.5) | 31.108841187 | 31.108841187 + (5.1,34.5) | (3,3),(3,3) | 31.5699223946 | 31.5699223946 + (-5,-12) | (2,2),(0,0) | 13 | 13 + (-5,-12) | (3,3),(1,1) | 14.3178210633 | 14.3178210633 + (-5,-12) | (-2,2),(-8,-10) | 2 | 2 + (-5,-12) | (2.5,3.5),(2.5,2.5) | 16.3248277173 | 16.3248277173 + (-5,-12) | (3,3),(3,3) | 17 | 17 + (1e-300,-1e-300) | (2,2),(0,0) | 1e-300 | 1e-300 + (1e-300,-1e-300) | (3,3),(1,1) | 1.41421356237 | 1.41421356237 + (1e-300,-1e-300) | (-2,2),(-8,-10) | 2 | 2 + (1e-300,-1e-300) | (2.5,3.5),(2.5,2.5) | 3.53553390593 | 3.53553390593 + (1e-300,-1e-300) | (3,3),(3,3) | 4.24264068712 | 4.24264068712 + (1e+300,Infinity) | (2,2),(0,0) | Infinity | Infinity + (1e+300,Infinity) | (3,3),(1,1) | Infinity | Infinity + (1e+300,Infinity) | (-2,2),(-8,-10) | Infinity | Infinity + (1e+300,Infinity) | (2.5,3.5),(2.5,2.5) | Infinity | Infinity + (1e+300,Infinity) | (3,3),(3,3) | Infinity | Infinity + (Infinity,1e+300) | (2,2),(0,0) | Infinity | Infinity + (Infinity,1e+300) | (3,3),(1,1) | Infinity | Infinity + (Infinity,1e+300) | (-2,2),(-8,-10) | Infinity | Infinity + (Infinity,1e+300) | (2.5,3.5),(2.5,2.5) | Infinity | Infinity + (Infinity,1e+300) | (3,3),(3,3) | Infinity | Infinity + (NaN,NaN) | (2,2),(0,0) | NaN | NaN + (NaN,NaN) | (3,3),(1,1) | NaN | NaN + (NaN,NaN) | (-2,2),(-8,-10) | NaN | NaN + (NaN,NaN) | (2.5,3.5),(2.5,2.5) | NaN | NaN + (NaN,NaN) | (3,3),(3,3) | NaN | NaN + (10,10) | (2,2),(0,0) | 11.313708499 | 11.313708499 + (10,10) | (3,3),(1,1) | 9.89949493661 | 9.89949493661 + (10,10) | (-2,2),(-8,-10) | 14.4222051019 | 14.4222051019 + (10,10) | (2.5,3.5),(2.5,2.5) | 9.92471662064 | 9.92471662064 + (10,10) | (3,3),(3,3) | 9.89949493661 | 9.89949493661 (50 rows) -- Distance to path @@ -1224,19 +1224,19 @@ SELECT p.f1, b.f1, p.f1 ## b.f1 FROM POINT_TBL p, BOX_TBL b; (-5,-12) | (-2,2),(-8,-10) | (-5,-10) (-5,-12) | (2.5,3.5),(2.5,2.5) | (2.5,2.5) (-5,-12) | (3,3),(3,3) | (3,3) - (1e-300,-1e-300) | (2,2),(0,0) | (0,0) + (1e-300,-1e-300) | (2,2),(0,0) | (1e-300,0) (1e-300,-1e-300) | (3,3),(1,1) | (1,1) (1e-300,-1e-300) | (-2,2),(-8,-10) | (-2,-1e-300) (1e-300,-1e-300) | (2.5,3.5),(2.5,2.5) | (2.5,2.5) (1e-300,-1e-300) | (3,3),(3,3) | (3,3) - (1e+300,Infinity) | (2,2),(0,0) | (0,2) - (1e+300,Infinity) | (3,3),(1,1) | (1,3) - (1e+300,Infinity) | (-2,2),(-8,-10) | (-8,2) + (1e+300,Infinity) | (2,2),(0,0) | (2,2) + (1e+300,Infinity) | (3,3),(1,1) | (3,3) + (1e+300,Infinity) | (-2,2),(-8,-10) | (-2,2) (1e+300,Infinity) | (2.5,3.5),(2.5,2.5) | (2.5,3.5) (1e+300,Infinity) | (3,3),(3,3) | (3,3) - (Infinity,1e+300) | (2,2),(0,0) | (0,2) - (Infinity,1e+300) | (3,3),(1,1) | (1,3) - (Infinity,1e+300) | (-2,2),(-8,-10) | (-8,2) + (Infinity,1e+300) | (2,2),(0,0) | (2,2) + (Infinity,1e+300) | (3,3),(1,1) | (3,3) + (Infinity,1e+300) | (-2,2),(-8,-10) | (-2,2) (Infinity,1e+300) | (2.5,3.5),(2.5,2.5) | (2.5,3.5) (Infinity,1e+300) | (3,3),(3,3) | (3,3) (NaN,NaN) | (2,2),(0,0) | diff --git a/src/test/regress/sql/box.sql b/src/test/regress/sql/box.sql index 2d0868e8898..13f89c18b6c 100644 --- a/src/test/regress/sql/box.sql +++ b/src/test/regress/sql/box.sql @@ -282,6 +282,30 @@ RESET enable_seqscan; RESET enable_indexscan; RESET enable_bitmapscan; +-- The distance from a point to a box is exact, so a kNN index scan, which +-- orders by the index's own distance, returns rows in the operator's order. +-- For the point (0.9999996,2), box 1 is 1 away and box 2 is 1.00000000000004 +-- away. (The operator used to give box 1 1.00000000000008, out of order.) +-- Both distances are exact on any platform: the point is above both boxes, +-- so each is hypot(0, dy) = |dy|. +CREATE TEMP TABLE box_knn (id int, b box); +INSERT INTO box_knn VALUES (1, '((0,0),(1,1))'), + (2, box(point(0.4999996, 0), point(1.4999996, 0.99999999999996))); +SET extra_float_digits = 0; +SET enable_seqscan = OFF; +SET enable_bitmapscan = OFF; +CREATE INDEX box_knn_gist ON box_knn USING gist (b); +SELECT id, b <-> point(0.9999996, 2) AS dist FROM box_knn + ORDER BY b <-> point(0.9999996, 2); +DROP INDEX box_knn_gist; +CREATE INDEX box_knn_spgist ON box_knn USING spgist (b); +SELECT id, b <-> point(0.9999996, 2) AS dist FROM box_knn + ORDER BY b <-> point(0.9999996, 2); +RESET extra_float_digits; +RESET enable_seqscan; +RESET enable_bitmapscan; +DROP TABLE box_knn; + -- test non-error-throwing API for some core types SELECT pg_input_is_valid('200', 'box'); SELECT * FROM pg_input_error_info('200', 'box'); -- 2.50.1