From b087b018a11f4bbff6eb2a7e3f6f30cc8d0a03c9 Mon Sep 17 00:00:00 2001
From: Tom Lane <tgl@sss.pgh.pa.us>
Date: Tue, 23 Jun 2009 16:25:02 +0000
Subject: [PATCH] Fix an ancient error in dist_ps (distance from point to line
segment), which a number of other geometric operators also depend on. It
miscalculated the slope of the perpendicular to the given line segment
anytime that slope was other than 0, infinite, or +/-1. In some cases the
error would be masked because the true closest point on the line segment was
one of its endpoints rather than the intersection point, but in other cases
it could give an arbitrarily bad answer. Per bug #4872 from Nick Roosevelt.
Bug goes clear back to Berkeley days, so patch all supported branches.
Make a couple of cosmetic adjustments while at it.
---
src/backend/utils/adt/geo_ops.c | 29 ++++++++++++++---------------
1 file changed, 14 insertions(+), 15 deletions(-)
diff --git a/src/backend/utils/adt/geo_ops.c b/src/backend/utils/adt/geo_ops.c
index c2444c38e88..24157a53c59 100644
--- a/src/backend/utils/adt/geo_ops.c
+++ b/src/backend/utils/adt/geo_ops.c
@@ -8,7 +8,7 @@
*
*
* IDENTIFICATION
- * $PostgreSQL: pgsql/src/backend/utils/adt/geo_ops.c,v 1.101 2009/01/01 17:23:49 momjian Exp $
+ * $PostgreSQL: pgsql/src/backend/utils/adt/geo_ops.c,v 1.102 2009/06/23 16:25:02 tgl Exp $
*
*-------------------------------------------------------------------------
*/
@@ -2449,18 +2449,16 @@ dist_ps_internal(Point *pt, LSEG *lseg)
tmpdist;
Point *ip;
-/*
- * Construct a line perpendicular to the input segment
- * and through the input point
- */
+ /*
+ * Construct a line perpendicular to the input segment
+ * and through the input point
+ */
if (lseg->p[1].x == lseg->p[0].x)
m = 0;
else if (lseg->p[1].y == lseg->p[0].y)
- { /* slope is infinite */
- m = (double) DBL_MAX;
- }
+ m = (double) DBL_MAX; /* slope is infinite */
else
- m = ((lseg->p[0].y - lseg->p[1].y) / (lseg->p[1].x - lseg->p[0].x));
+ m = (lseg->p[0].x - lseg->p[1].x) / (lseg->p[1].y - lseg->p[0].y);
ln = line_construct_pm(pt, m);
#ifdef GEODEBUG
@@ -2469,13 +2467,14 @@ dist_ps_internal(Point *pt, LSEG *lseg)
#endif
/*
- * Calculate distance to the line segment or to the endpoints of the
- * segment.
+ * Calculate distance to the line segment or to the nearest endpoint of
+ * the segment.
*/
/* intersection is on the line segment? */
if ((ip = interpt_sl(lseg, ln)) != NULL)
{
+ /* yes, so use distance to the intersection point */
result = point_dt(pt, ip);
#ifdef GEODEBUG
printf("dist_ps- distance is %f to intersection point is (%f,%f)\n",
@@ -2484,7 +2483,7 @@ dist_ps_internal(Point *pt, LSEG *lseg)
}
else
{
- /* intersection is not on line segment */
+ /* no, so use distance to the nearer endpoint */
result = point_dt(pt, &lseg->p[0]);
tmpdist = point_dt(pt, &lseg->p[1]);
if (tmpdist < result)
@@ -3790,7 +3789,7 @@ poly_contain(PG_FUNCTION_ARGS)
{
if (point_inside(&(polyb->p[i]), polya->npts, &(polya->p[0])) == 0)
{
-#if GEODEBUG
+#ifdef GEODEBUG
printf("poly_contain- point (%f,%f) not in polygon\n", polyb->p[i].x, polyb->p[i].y);
#endif
result = false;
@@ -3803,7 +3802,7 @@ poly_contain(PG_FUNCTION_ARGS)
{
if (point_inside(&(polya->p[i]), polyb->npts, &(polyb->p[0])) == 1)
{
-#if GEODEBUG
+#ifdef GEODEBUG
printf("poly_contain- point (%f,%f) in polygon\n", polya->p[i].x, polya->p[i].y);
#endif
result = false;
@@ -3814,7 +3813,7 @@ poly_contain(PG_FUNCTION_ARGS)
}
else
{
-#if GEODEBUG
+#ifdef GEODEBUG
printf("poly_contain- bound box ((%f,%f),(%f,%f)) not inside ((%f,%f),(%f,%f))\n",
polyb->boundbox.low.x, polyb->boundbox.low.y, polyb->boundbox.high.x, polyb->boundbox.high.y,
polya->boundbox.low.x, polya->boundbox.low.y, polya->boundbox.high.x, polya->boundbox.high.y);
--
GitLab