Skip to content

P
Postgres FD Implementation
Commit a7dc63d9 authored by Tomas Vondra's avatar Tomas Vondra
Browse files
Options

Refactor geometric functions and operators 

The primary goal of this patch is to eliminate duplicate code and share
code between different geometric data types more often, to prepare the
ground for additional patches.  Until now the code reuse was limited,
probably because the simpler types (line and point) were implemented
after the more complex ones.

The changes are quite extensive and can be summarised as:

* Eliminate SQL-level function calls.
* Re-use more functions to implement others.
* Unify internal function names and signatures.
* Remove private functions from geo_decls.h.
* Replace should-not-happen checks with assertions.
* Add comments describe for various functions.
* Remove some unreachable code.
* Define delimiter symbols of line datatype like the other ones.
* Remove the GEODEBUG macro and printf() calls.
* Unify code style of a few oddly formatted lines.

While the goal was to cause minimal user-visible changes, it was not
possible to keep the original behavior in all cases - for example when
handling NaN values, or when reusing code makes the functions return
consistent results.

Author: Emre Hasegeli
Reviewed-by: Kyotaro Horiguchi, me

Discussion: https://www.postgresql.org/message-id/CAE2gYzxF7-5djV6-cEvqQu-fNsnt%3DEqbOURx7ZDg%2BVv6ZMTWbg%40mail.gmail.com
parent fa7d5b70
Show whitespace changes
Inline Side-by-side
Showing with 870 additions and 1026 deletions
src/backend/utils/adt/geo_ops.c
View file @ a7dc63d9
......@@ -38,24 +38,62 @@ enum path_delim
PATH_NONE, PATH_OPEN, PATH_CLOSED
};
/* Routines for points */
static inline void point_construct(Point *result, float8 x, float8 y);
static inline void point_add_point(Point *result, Point *pt1, Point *pt2);
static inline void point_sub_point(Point *result, Point *pt1, Point *pt2);
static inline void point_mul_point(Point *result, Point *pt1, Point *pt2);
static inline void point_div_point(Point *result, Point *pt1, Point *pt2);
static inline bool point_eq_point(Point *pt1, Point *pt2);
static inline float8 point_dt(Point *pt1, Point *pt2);
static inline float8 point_sl(Point *pt1, Point *pt2);
static int point_inside(Point *p, int npts, Point *plist);
/* Routines for lines */
static inline void line_construct(LINE *result, Point *pt, float8 m);
static inline float8 line_invsl(LINE *line);
static bool line_interpt_line(Point *result, LINE *l1, LINE *l2);
static bool line_contain_point(LINE *line, Point *point);
static float8 line_closept_point(Point *result, LINE *line, Point *pt);
/* Routines for line segments */
static inline void statlseg_construct(LSEG *lseg, Point *pt1, Point *pt2);
static inline float8 lseg_sl(LSEG *lseg);
static inline float8 lseg_invsl(LSEG *lseg);
static bool lseg_interpt_line(Point *result, LSEG *lseg, LINE *line);
static bool lseg_interpt_lseg(Point *result, LSEG *l1, LSEG *l2);
static int lseg_crossing(double x, double y, double px, double py);
static BOX *box_construct(double x1, double x2, double y1, double y2);
static BOX *box_fill(BOX *result, double x1, double x2, double y1, double y2);
static bool lseg_contain_point(LSEG *lseg, Point *point);
static float8 lseg_closept_point(Point *result, LSEG *lseg, Point *pt);
static float8 lseg_closept_line(Point *result, LSEG *lseg, LINE *line);
static float8 lseg_closept_lseg(Point *result, LSEG *l1, LSEG *l2);
/* Routines for boxes */
static inline void box_construct(BOX *result, Point *pt1, Point *pt2);
static void box_cn(Point *center, BOX *box);
static bool box_ov(BOX *box1, BOX *box2);
static double box_ar(BOX *box);
static double box_ht(BOX *box);
static double box_wd(BOX *box);
static bool box_contain_point(BOX *box, Point *point);
static bool box_contain_box(BOX *box1, BOX *box2);
static bool box_contain_lseg(BOX *box, LSEG *lseg);
static bool box_interpt_lseg(Point *result, BOX *box, LSEG *lseg);
static float8 box_closept_point(Point *result, BOX *box, Point *point);
static float8 box_closept_lseg(Point *result, BOX *box, LSEG *lseg);
/* Routines for circles */
static double circle_ar(CIRCLE *circle);
static CIRCLE *circle_copy(CIRCLE *circle);
static LINE *line_construct_pm(Point *pt, double m);
static void line_construct_pts(LINE *line, Point *pt1, Point *pt2);
static bool lseg_intersect_internal(LSEG *l1, LSEG *l2);
static double lseg_dt(LSEG *l1, LSEG *l2);
static bool on_ps_internal(Point *pt, LSEG *lseg);
/* Routines for polygons */
static void make_bound_box(POLYGON *poly);
static void poly_to_circle(CIRCLE *result, POLYGON *poly);
static bool lseg_inside_poly(Point *a, Point *b, POLYGON *poly, int start);
static bool poly_contain_poly(POLYGON *polya, POLYGON *polyb);
static bool plist_same(int npts, Point *p1, Point *p2);
static Point *point_construct(double x, double y);
static Point *point_copy(Point *pt);
static float8 dist_ppoly_internal(Point *pt, POLYGON *poly);
/* Routines for encoding and decoding */
static double single_decode(char *num, char **endptr_p,
const char *type_name, const char *orig_string);
static void single_encode(float8 x, StringInfo str);
......@@ -67,17 +105,6 @@ static void path_decode(char *str, bool opentype, int npts, Point *p,
bool *isopen, char **endptr_p,
const char *type_name, const char *orig_string);
static char *path_encode(enum path_delim path_delim, int npts, Point *pt);
static void statlseg_construct(LSEG *lseg, Point *pt1, Point *pt2);
static double box_ar(BOX *box);
static void box_cn(Point *center, BOX *box);
static Point *interpt_sl(LSEG *lseg, LINE *line);
static bool has_interpt_sl(LSEG *lseg, LINE *line);
static double dist_pl_internal(Point *pt, LINE *line);
static double dist_ps_internal(Point *pt, LSEG *lseg);
static Point *line_interpt_internal(LINE *l1, LINE *l2);
static bool lseg_inside_poly(Point *a, Point *b, POLYGON *poly, int start);
static Point *lseg_interpt_internal(LSEG *l1, LSEG *l2);
static double dist_ppoly_internal(Point *pt, POLYGON *poly);
/*
......@@ -93,6 +120,8 @@ static double dist_ppoly_internal(Point *pt, POLYGON *poly);
#define RDELIM_EP ']'
#define LDELIM_C '<'
#define RDELIM_C '>'
#define LDELIM_L '{'
#define RDELIM_L '}'
/*
......@@ -236,12 +265,9 @@ path_decode(char *str, bool opentype, int npts, Point *p,
p++;
}
while (isspace((unsigned char) *str))
str++;
while (depth > 0)
{
if ((*str == RDELIM)
|| ((*str == RDELIM_EP) && (*isopen) && (depth == 1)))
if (*str == RDELIM || (*str == RDELIM_EP && *isopen && depth == 1))
{
depth--;
str++;
......@@ -440,55 +466,29 @@ box_send(PG_FUNCTION_ARGS)
/* box_construct - fill in a new box.
*/
static BOX *
box_construct(double x1, double x2, double y1, double y2)
{
BOX *result = (BOX *) palloc(sizeof(BOX));
return box_fill(result, x1, x2, y1, y2);
}
/* box_fill - fill in a given box struct
*/
static BOX *
box_fill(BOX *result, double x1, double x2, double y1, double y2)
static inline void
box_construct(BOX *result, Point *pt1, Point *pt2)
{
if (x1 > x2)
if (pt1->x > pt2->x)
{
result->high.x = x1;
result->low.x = x2;
result->high.x = pt1->x;
result->low.x = pt2->x;
}
else
{
result->high.x = x2;
result->low.x = x1;
result->high.x = pt2->x;
result->low.x = pt1->x;
}
if (y1 > y2)
if (pt1->y > pt2->y)
{
result->high.y = y1;
result->low.y = y2;
result->high.y = pt1->y;
result->low.y = pt2->y;
}
else
{
result->high.y = y2;
result->low.y = y1;
result->high.y = pt2->y;
result->low.y = pt1->y;
}
return result;
}
/* box_copy - copy a box
*/
BOX *
box_copy(BOX *box)
{
BOX *result = (BOX *) palloc(sizeof(BOX));
memcpy((char *) result, (char *) box, sizeof(BOX));
return result;
}
......@@ -505,10 +505,8 @@ box_same(PG_FUNCTION_ARGS)
BOX *box1 = PG_GETARG_BOX_P(0);
BOX *box2 = PG_GETARG_BOX_P(1);
PG_RETURN_BOOL(FPeq(box1->high.x, box2->high.x) &&
FPeq(box1->low.x, box2->low.x) &&
FPeq(box1->high.y, box2->high.y) &&
FPeq(box1->low.y, box2->low.y));
PG_RETURN_BOOL(point_eq_point(&box1->high, &box2->high) &&
point_eq_point(&box1->low, &box2->low));
}
/* box_overlap - does box1 overlap box2?
......@@ -637,10 +635,7 @@ box_contained(PG_FUNCTION_ARGS)
BOX *box1 = PG_GETARG_BOX_P(0);
BOX *box2 = PG_GETARG_BOX_P(1);
PG_RETURN_BOOL(FPle(box1->high.x, box2->high.x) &&
FPge(box1->low.x, box2->low.x) &&
FPle(box1->high.y, box2->high.y) &&
FPge(box1->low.y, box2->low.y));
PG_RETURN_BOOL(box_contain_box(box2, box1));
}
/* box_contain - does box1 contain box2?
......@@ -651,10 +646,19 @@ box_contain(PG_FUNCTION_ARGS)
BOX *box1 = PG_GETARG_BOX_P(0);
BOX *box2 = PG_GETARG_BOX_P(1);
PG_RETURN_BOOL(FPge(box1->high.x, box2->high.x) &&
PG_RETURN_BOOL(box_contain_box(box1, box2));
}
/*
* Check whether the box is in the box or on its border
*/
static bool
box_contain_box(BOX *box1, BOX *box2)
{
return FPge(box1->high.x, box2->high.x) &&
FPle(box1->low.x, box2->low.x) &&
FPge(box1->high.y, box2->high.y) &&
FPle(box1->low.y, box2->low.y));
FPle(box1->low.y, box2->low.y);
}
......@@ -757,7 +761,7 @@ box_width(PG_FUNCTION_ARGS)
{
BOX *box = PG_GETARG_BOX_P(0);
PG_RETURN_FLOAT8(box->high.x - box->low.x);
PG_RETURN_FLOAT8(box_wd(box));
}
......@@ -769,7 +773,7 @@ box_height(PG_FUNCTION_ARGS)
{
BOX *box = PG_GETARG_BOX_P(0);
PG_RETURN_FLOAT8(box->high.y - box->low.y);
PG_RETURN_FLOAT8(box_ht(box));
}
......@@ -787,7 +791,7 @@ box_distance(PG_FUNCTION_ARGS)
box_cn(&a, box1);
box_cn(&b, box2);
PG_RETURN_FLOAT8(HYPOT(a.x - b.x, a.y - b.y));
PG_RETURN_FLOAT8(point_dt(&a, &b));
}
......@@ -905,7 +909,7 @@ line_decode(char *s, const char *str, LINE *line)
if (*s++ != DELIM)
return false;
line->C = single_decode(s, &s, "line", str);
if (*s++ != '}')
if (*s++ != RDELIM_L)
return false;
while (isspace((unsigned char) *s))
s++;
......@@ -926,7 +930,7 @@ line_in(PG_FUNCTION_ARGS)
s = str;
while (isspace((unsigned char) *s))
s++;
if (*s == '{')
if (*s == LDELIM_L)
{
if (!line_decode(s + 1, str, line))
ereport(ERROR,
......@@ -940,12 +944,12 @@ line_in(PG_FUNCTION_ARGS)
}
else
{
path_decode(s, true, 2, &(lseg.p[0]), &isopen, NULL, "line", str);
if (FPeq(lseg.p[0].x, lseg.p[1].x) && FPeq(lseg.p[0].y, lseg.p[1].y))
path_decode(s, true, 2, &lseg.p[0], &isopen, NULL, "line", str);
if (point_eq_point(&lseg.p[0], &lseg.p[1]))
ereport(ERROR,
(errcode(ERRCODE_INVALID_TEXT_REPRESENTATION),
errmsg("invalid line specification: must be two distinct points")));
line_construct_pts(line, &lseg.p[0], &lseg.p[1]);
line_construct(line, &lseg.p[0], lseg_sl(&lseg));
}
PG_RETURN_LINE_P(line);
......@@ -960,7 +964,8 @@ line_out(PG_FUNCTION_ARGS)
char *bstr = float8out_internal(line->B);
char *cstr = float8out_internal(line->C);
PG_RETURN_CSTRING(psprintf("{%s,%s,%s}", astr, bstr, cstr));
PG_RETURN_CSTRING(psprintf("%c%s%c%s%c%s%c", LDELIM_L, astr, DELIM, bstr,
DELIM, cstr, RDELIM_L));
}
/*
......@@ -1003,14 +1008,12 @@ line_send(PG_FUNCTION_ARGS)
* Internal form: Ax+By+C=0
*---------------------------------------------------------*/
/* line_construct_pm()
* point-slope
/*
* Fill already-allocated LINE struct from the point and the slope
*/
static LINE *
line_construct_pm(Point *pt, double m)
static inline void
line_construct(LINE *result, Point *pt, float8 m)
{
LINE *result = (LINE *) palloc(sizeof(LINE));
if (m == DBL_MAX)
{
/* vertical - use "x = C" */
......@@ -1025,50 +1028,6 @@ line_construct_pm(Point *pt, double m)
result->B = -1.0;
result->C = pt->y - m * pt->x;
}
return result;
}
/*
* Fill already-allocated LINE struct from two points on the line
*/
static void
line_construct_pts(LINE *line, Point *pt1, Point *pt2)
{
if (FPeq(pt1->x, pt2->x))
{ /* vertical */
/* use "x = C" */
line->A = -1;
line->B = 0;
line->C = pt1->x;
#ifdef GEODEBUG
printf("line_construct_pts- line is vertical\n");
#endif
}
else if (FPeq(pt1->y, pt2->y))
{ /* horizontal */
/* use "y = C" */
line->A = 0;
line->B = -1;
line->C = pt1->y;
#ifdef GEODEBUG
printf("line_construct_pts- line is horizontal\n");
#endif
}
else
{
/* use "mx - y + yinter = 0" */
line->A = (pt2->y - pt1->y) / (pt2->x - pt1->x);
line->B = -1.0;
line->C = pt1->y - line->A * pt1->x;
/* on some platforms, the preceding expression tends to produce -0 */
if (line->C == 0.0)
line->C = 0.0;
#ifdef GEODEBUG
printf("line_construct_pts- line is neither vertical nor horizontal (diffs x=%.*g, y=%.*g\n",
DBL_DIG, (pt2->x - pt1->x), DBL_DIG, (pt2->y - pt1->y));
#endif
}
}
/* line_construct_pp()
......@@ -1081,7 +1040,8 @@ line_construct_pp(PG_FUNCTION_ARGS)
Point *pt2 = PG_GETARG_POINT_P(1);
LINE *result = (LINE *) palloc(sizeof(LINE));
line_construct_pts(result, pt1, pt2);
line_construct(result, pt1, point_sl(pt1, pt2));
PG_RETURN_LINE_P(result);
}
......@@ -1096,9 +1056,7 @@ line_intersect(PG_FUNCTION_ARGS)
LINE *l1 = PG_GETARG_LINE_P(0);
LINE *l2 = PG_GETARG_LINE_P(1);
PG_RETURN_BOOL(!DatumGetBool(DirectFunctionCall2(line_parallel,
LinePGetDatum(l1),
LinePGetDatum(l2))));
PG_RETURN_BOOL(line_interpt_line(NULL, l1, l2));
}
Datum
......@@ -1107,10 +1065,7 @@ line_parallel(PG_FUNCTION_ARGS)
LINE *l1 = PG_GETARG_LINE_P(0);
LINE *l2 = PG_GETARG_LINE_P(1);
if (FPzero(l1->B))
PG_RETURN_BOOL(FPzero(l2->B));
PG_RETURN_BOOL(FPeq(l2->A, l1->A * (l2->B / l1->B)));
PG_RETURN_BOOL(!line_interpt_line(NULL, l1, l2));
}
Datum
......@@ -1169,6 +1124,20 @@ line_eq(PG_FUNCTION_ARGS)
* Line arithmetic routines.
*---------------------------------------------------------*/
/*
* Return inverse slope of the line
*/
static inline float8
line_invsl(LINE *line)
{
if (FPzero(line->A))
return DBL_MAX;
if (FPzero(line->B))
return 0.0;
return line->B / line->A;
}
/* line_distance()
* Distance between two lines.
*/
......@@ -1178,16 +1147,14 @@ line_distance(PG_FUNCTION_ARGS)
LINE *l1 = PG_GETARG_LINE_P(0);
LINE *l2 = PG_GETARG_LINE_P(1);
float8 result;
Point *tmp;
Point tmp;
if (!DatumGetBool(DirectFunctionCall2(line_parallel,
LinePGetDatum(l1),
LinePGetDatum(l2))))
if (line_interpt_line(NULL, l1, l2)) /* intersecting? */
PG_RETURN_FLOAT8(0.0);
if (FPzero(l1->B)) /* vertical? */
PG_RETURN_FLOAT8(fabs(l1->C - l2->C));
tmp = point_construct(0.0, l1->C);
result = dist_pl_internal(tmp, l2);
point_construct(&tmp, 0.0, l1->C);
result = line_closept_point(NULL, l2, &tmp);
PG_RETURN_FLOAT8(result);
}
......@@ -1201,9 +1168,9 @@ line_interpt(PG_FUNCTION_ARGS)
LINE *l2 = PG_GETARG_LINE_P(1);
Point *result;
result = line_interpt_internal(l1, l2);
result = (Point *) palloc(sizeof(Point));
if (result == NULL)
if (!line_interpt_line(result, l1, l2))
PG_RETURN_NULL();
PG_RETURN_POINT_P(result);
}
......@@ -1211,27 +1178,25 @@ line_interpt(PG_FUNCTION_ARGS)
/*
* Internal version of line_interpt
*
* returns a NULL pointer if no intersection point
* This returns true if two lines intersect (they do, if they are not
* parallel), false if they do not. This also sets the intersection point
* to *result, if it is not NULL.
*
* NOTE: If the lines are identical then we will find they are parallel
* and report "no intersection". This is a little weird, but since
* there's no *unique* intersection, maybe it's appropriate behavior.
*/
static Point *
line_interpt_internal(LINE *l1, LINE *l2)
static bool
line_interpt_line(Point *result, LINE *l1, LINE *l2)
{
Point *result;
double x,
y;
/*
* NOTE: if the lines are identical then we will find they are parallel
* and report "no intersection". This is a little weird, but since
* there's no *unique* intersection, maybe it's appropriate behavior.
*/
if (DatumGetBool(DirectFunctionCall2(line_parallel,
LinePGetDatum(l1),
LinePGetDatum(l2))))
return NULL;
if (FPzero(l1->B)) /* l1 vertical? */
{
if (FPzero(l2->B)) /* l2 vertical? */
return false;
x = l1->C;
y = (l2->A * x + l2->C);
}
......@@ -1242,18 +1207,17 @@ line_interpt_internal(LINE *l1, LINE *l2)
}
else
{
if (FPeq(l2->A, l1->A * (l2->B / l1->B)))
return false;
x = (l1->C - l2->C) / (l2->A - l1->A);
y = (l1->A * x + l1->C);
}
result = point_construct(x, y);
#ifdef GEODEBUG
printf("line_interpt- lines are A=%.*g, B=%.*g, C=%.*g, A=%.*g, B=%.*g, C=%.*g\n",
DBL_DIG, l1->A, DBL_DIG, l1->B, DBL_DIG, l1->C, DBL_DIG, l2->A, DBL_DIG, l2->B, DBL_DIG, l2->C);
printf("line_interpt- lines intersect at (%.*g,%.*g)\n", DBL_DIG, x, DBL_DIG, y);
#endif
if (result != NULL)
point_construct(result, x, y);
return result;
return true;
}
......@@ -1562,8 +1526,7 @@ path_inter(PG_FUNCTION_ARGS)
LSEG seg1,
seg2;
if (p1->npts <= 0 || p2->npts <= 0)
PG_RETURN_BOOL(false);
Assert(p1->npts > 0 && p2->npts > 0);
b1.high.x = b1.low.x = p1->p[0].x;
b1.high.y = b1.low.y = p1->p[0].y;
......@@ -1615,7 +1578,7 @@ path_inter(PG_FUNCTION_ARGS)
statlseg_construct(&seg1, &p1->p[iprev], &p1->p[i]);
statlseg_construct(&seg2, &p2->p[jprev], &p2->p[j]);
if (lseg_intersect_internal(&seg1, &seg2))
if (lseg_interpt_lseg(NULL, &seg1, &seg2))
PG_RETURN_BOOL(true);
}
}
......@@ -1670,9 +1633,7 @@ path_distance(PG_FUNCTION_ARGS)
statlseg_construct(&seg1, &p1->p[iprev], &p1->p[i]);
statlseg_construct(&seg2, &p2->p[jprev], &p2->p[j]);
tmp = DatumGetFloat8(DirectFunctionCall2(lseg_distance,
LsegPGetDatum(&seg1),
LsegPGetDatum(&seg2)));
tmp = lseg_closept_lseg(NULL, &seg1, &seg2);
if (!have_min || tmp < min)
{
min = tmp;
......@@ -1780,30 +1741,14 @@ point_send(PG_FUNCTION_ARGS)
}
static Point *
point_construct(double x, double y)
/*
* Initialize a point
*/
static inline void
point_construct(Point *result, float8 x, float8 y)
{
Point *result = (Point *) palloc(sizeof(Point));
result->x = x;
result->y = y;
return result;
}
static Point *
point_copy(Point *pt)
{
Point *result;
if (!PointerIsValid(pt))
return NULL;
result = (Point *) palloc(sizeof(Point));
result->x = pt->x;
result->y = pt->y;
return result;
}
......@@ -1876,7 +1821,7 @@ point_eq(PG_FUNCTION_ARGS)
Point *pt1 = PG_GETARG_POINT_P(0);
Point *pt2 = PG_GETARG_POINT_P(1);
PG_RETURN_BOOL(FPeq(pt1->x, pt2->x) && FPeq(pt1->y, pt2->y));
PG_RETURN_BOOL(point_eq_point(pt1, pt2));
}
Datum
......@@ -1885,9 +1830,16 @@ point_ne(PG_FUNCTION_ARGS)
Point *pt1 = PG_GETARG_POINT_P(0);
Point *pt2 = PG_GETARG_POINT_P(1);
PG_RETURN_BOOL(FPne(pt1->x, pt2->x) || FPne(pt1->y, pt2->y));
PG_RETURN_BOOL(!point_eq_point(pt1, pt2));
}
static inline bool
point_eq_point(Point *pt1, Point *pt2)
{
return FPeq(pt1->x, pt2->x) && FPeq(pt1->y, pt2->y);
}
/*----------------------------------------------------------
* "Arithmetic" operators on points.
*---------------------------------------------------------*/
......@@ -1898,16 +1850,12 @@ point_distance(PG_FUNCTION_ARGS)
Point *pt1 = PG_GETARG_POINT_P(0);
Point *pt2 = PG_GETARG_POINT_P(1);
PG_RETURN_FLOAT8(HYPOT(pt1->x - pt2->x, pt1->y - pt2->y));
PG_RETURN_FLOAT8(point_dt(pt1, pt2));
}
double
static inline float8
point_dt(Point *pt1, Point *pt2)
{
#ifdef GEODEBUG
printf("point_dt- segment (%f,%f),(%f,%f) length is %f\n",
pt1->x, pt1->y, pt2->x, pt2->y, HYPOT(pt1->x - pt2->x, pt1->y - pt2->y));
#endif
return HYPOT(pt1->x - pt2->x, pt1->y - pt2->y);
}
......@@ -1921,12 +1869,35 @@ point_slope(PG_FUNCTION_ARGS)
}
double
/*
* Return slope of two points
*
* Note that this function returns DBL_MAX when the points are the same.
*/
static inline float8
point_sl(Point *pt1, Point *pt2)
{
return (FPeq(pt1->x, pt2->x)
? (double) DBL_MAX
: (pt1->y - pt2->y) / (pt1->x - pt2->x));
if (FPeq(pt1->x, pt2->x))
return DBL_MAX;
if (FPeq(pt1->y, pt2->y))
return 0.0;
return (pt1->y - pt2->y) / (pt1->x - pt2->x);
}
/*
* Return inverse slope of two points
*
* Note that this function returns 0.0 when the points are the same.
*/
static inline float8
point_invsl(Point *pt1, Point *pt2)
{
if (FPeq(pt1->x, pt2->x))
return 0.0;
if (FPeq(pt1->y, pt2->y))
return DBL_MAX;
return (pt1->x - pt2->x) / (pt2->y - pt1->y);
}
......@@ -1952,7 +1923,7 @@ lseg_in(PG_FUNCTION_ARGS)
LSEG *lseg = (LSEG *) palloc(sizeof(LSEG));
bool isopen;
path_decode(str, true, 2, &(lseg->p[0]), &isopen, NULL, "lseg", str);
path_decode(str, true, 2, &lseg->p[0], &isopen, NULL, "lseg", str);
PG_RETURN_LSEG_P(lseg);
}
......@@ -1962,7 +1933,7 @@ lseg_out(PG_FUNCTION_ARGS)
{
LSEG *ls = PG_GETARG_LSEG_P(0);
PG_RETURN_CSTRING(path_encode(PATH_OPEN, 2, (Point *) &(ls->p[0])));
PG_RETURN_CSTRING(path_encode(PATH_OPEN, 2, &ls->p[0]));
}
/*
......@@ -2012,16 +1983,13 @@ lseg_construct(PG_FUNCTION_ARGS)
Point *pt2 = PG_GETARG_POINT_P(1);
LSEG *result = (LSEG *) palloc(sizeof(LSEG));
result->p[0].x = pt1->x;
result->p[0].y = pt1->y;
result->p[1].x = pt2->x;
result->p[1].y = pt2->y;
statlseg_construct(result, pt1, pt2);
PG_RETURN_LSEG_P(result);
}
/* like lseg_construct, but assume space already allocated */
static void
static inline void
statlseg_construct(LSEG *lseg, Point *pt1, Point *pt2)
{
lseg->p[0].x = pt1->x;
......@@ -2030,6 +1998,27 @@ statlseg_construct(LSEG *lseg, Point *pt1, Point *pt2)
lseg->p[1].y = pt2->y;
}
/*
* Return slope of the line segment
*/
static inline float8
lseg_sl(LSEG *lseg)
{
return point_sl(&lseg->p[0], &lseg->p[1]);
}
/*
* Return inverse slope of the line segment
*/
static inline float8
lseg_invsl(LSEG *lseg)
{
return point_invsl(&lseg->p[0], &lseg->p[1]);
}
Datum
lseg_length(PG_FUNCTION_ARGS)
{
......@@ -2052,25 +2041,9 @@ lseg_intersect(PG_FUNCTION_ARGS)
LSEG *l1 = PG_GETARG_LSEG_P(0);
LSEG *l2 = PG_GETARG_LSEG_P(1);
PG_RETURN_BOOL(lseg_intersect_internal(l1, l2));
PG_RETURN_BOOL(lseg_interpt_lseg(NULL, l1, l2));
}
static bool
lseg_intersect_internal(LSEG *l1, LSEG *l2)
{
LINE ln;
Point *interpt;
bool retval;
line_construct_pts(&ln, &l2->p[0], &l2->p[1]);
interpt = interpt_sl(l1, &ln);
if (interpt != NULL && on_ps_internal(interpt, l2))
retval = true; /* interpt on l1 and l2 */
else
retval = false;
return retval;
}
Datum
lseg_parallel(PG_FUNCTION_ARGS)
......@@ -2078,39 +2051,19 @@ lseg_parallel(PG_FUNCTION_ARGS)
LSEG *l1 = PG_GETARG_LSEG_P(0);
LSEG *l2 = PG_GETARG_LSEG_P(1);
PG_RETURN_BOOL(FPeq(point_sl(&l1->p[0], &l1->p[1]),
point_sl(&l2->p[0], &l2->p[1])));
PG_RETURN_BOOL(FPeq(lseg_sl(l1), lseg_sl(l2)));
}
/* lseg_perp()
/*
* Determine if two line segments are perpendicular.
*
* This code did not get the correct answer for
* '((0,0),(0,1))'::lseg ?-| '((0,0),(1,0))'::lseg
* So, modified it to check explicitly for slope of vertical line
* returned by point_sl() and the results seem better.
* - thomas 1998-01-31
*/
Datum
lseg_perp(PG_FUNCTION_ARGS)
{
LSEG *l1 = PG_GETARG_LSEG_P(0);
LSEG *l2 = PG_GETARG_LSEG_P(1);
double m1,
m2;
m1 = point_sl(&(l1->p[0]), &(l1->p[1]));
m2 = point_sl(&(l2->p[0]), &(l2->p[1]));
#ifdef GEODEBUG
printf("lseg_perp- slopes are %g and %g\n", m1, m2);
#endif
if (FPzero(m1))
PG_RETURN_BOOL(FPeq(m2, DBL_MAX));
else if (FPzero(m2))
PG_RETURN_BOOL(FPeq(m1, DBL_MAX));
PG_RETURN_BOOL(FPeq(m1 / m2, -1.0));
PG_RETURN_BOOL(FPeq(lseg_sl(l1), lseg_invsl(l2)));
}
Datum
......@@ -2136,10 +2089,8 @@ lseg_eq(PG_FUNCTION_ARGS)
LSEG *l1 = PG_GETARG_LSEG_P(0);
LSEG *l2 = PG_GETARG_LSEG_P(1);
PG_RETURN_BOOL(FPeq(l1->p[0].x, l2->p[0].x) &&
FPeq(l1->p[0].y, l2->p[0].y) &&
FPeq(l1->p[1].x, l2->p[1].x) &&
FPeq(l1->p[1].y, l2->p[1].y));
PG_RETURN_BOOL(point_eq_point(&l1->p[0], &l2->p[0]) &&
point_eq_point(&l1->p[1], &l2->p[1]));
}
Datum
......@@ -2148,10 +2099,8 @@ lseg_ne(PG_FUNCTION_ARGS)
LSEG *l1 = PG_GETARG_LSEG_P(0);
LSEG *l2 = PG_GETARG_LSEG_P(1);
PG_RETURN_BOOL(!FPeq(l1->p[0].x, l2->p[0].x) ||
!FPeq(l1->p[0].y, l2->p[0].y) ||
!FPeq(l1->p[1].x, l2->p[1].x) ||
!FPeq(l1->p[1].y, l2->p[1].y));
PG_RETURN_BOOL(!point_eq_point(&l1->p[0], &l2->p[0]) ||
!point_eq_point(&l1->p[1], &l2->p[1]));
}
Datum
......@@ -2210,33 +2159,7 @@ lseg_distance(PG_FUNCTION_ARGS)
LSEG *l1 = PG_GETARG_LSEG_P(0);
LSEG *l2 = PG_GETARG_LSEG_P(1);
PG_RETURN_FLOAT8(lseg_dt(l1, l2));
}
/* lseg_dt()
* Distance between two line segments.
* Must check both sets of endpoints to ensure minimum distance is found.
* - thomas 1998-02-01
*/
static double
lseg_dt(LSEG *l1, LSEG *l2)
{
double result,
d;
if (lseg_intersect_internal(l1, l2))
return 0.0;
d = dist_ps_internal(&l1->p[0], l2);
result = d;
d = dist_ps_internal(&l1->p[1], l2);
result = Min(result, d);
d = dist_ps_internal(&l2->p[0], l1);
result = Min(result, d);
d = dist_ps_internal(&l2->p[1], l1);
result = Min(result, d);
return result;
PG_RETURN_FLOAT8(lseg_closept_lseg(NULL, l1, l2));
}
......@@ -2254,57 +2177,38 @@ lseg_center(PG_FUNCTION_ARGS)
PG_RETURN_POINT_P(result);
}
static Point *
lseg_interpt_internal(LSEG *l1, LSEG *l2)
{
Point *result;
LINE tmp1,
tmp2;
/*
* Find the intersection of the appropriate lines, if any.
/*
* Find the intersection point of two segments (if any).
*
* This returns true if two line segments intersect, false if they do not.
* This also sets the intersection point to *result, if it is not NULL.
* This function is almost perfectly symmetric, even though it doesn't look
* like it. See lseg_interpt_line() for the other half of it.
*/
line_construct_pts(&tmp1, &l1->p[0], &l1->p[1]);
line_construct_pts(&tmp2, &l2->p[0], &l2->p[1]);
result = line_interpt_internal(&tmp1, &tmp2);
if (!PointerIsValid(result))
return NULL;
static bool
lseg_interpt_lseg(Point *result, LSEG *l1, LSEG *l2)
{
Point interpt;
LINE tmp;
/*
* If the line intersection point isn't within l1 (or equivalently l2),
* there is no valid segment intersection point at all.
*/
if (!on_ps_internal(result, l1) ||
!on_ps_internal(result, l2))
{
pfree(result);
return NULL;
}
line_construct(&tmp, &l2->p[0], lseg_sl(l2));
if (!lseg_interpt_line(&interpt, l1, &tmp))
return false;
/*
* If there is an intersection, then check explicitly for matching
* endpoints since there may be rounding effects with annoying lsb
* residue. - tgl 1997-07-09
* If the line intersection point isn't within l2, there is no valid
* segment intersection point at all.
*/
if ((FPeq(l1->p[0].x, l2->p[0].x) && FPeq(l1->p[0].y, l2->p[0].y)) ||
(FPeq(l1->p[0].x, l2->p[1].x) && FPeq(l1->p[0].y, l2->p[1].y)))
{
result->x = l1->p[0].x;
result->y = l1->p[0].y;
}
else if ((FPeq(l1->p[1].x, l2->p[0].x) && FPeq(l1->p[1].y, l2->p[0].y)) ||
(FPeq(l1->p[1].x, l2->p[1].x) && FPeq(l1->p[1].y, l2->p[1].y)))
{
result->x = l1->p[1].x;
result->y = l1->p[1].y;
}
if (!lseg_contain_point(l2, &interpt))
return false;
return result;
if (result != NULL)
*result = interpt;
return true;
}
/* lseg_interpt -
* Find the intersection point of two segments (if any).
*/
Datum
lseg_interpt(PG_FUNCTION_ARGS)
{
......@@ -2312,10 +2216,10 @@ lseg_interpt(PG_FUNCTION_ARGS)
LSEG *l2 = PG_GETARG_LSEG_P(1);
Point *result;
result = lseg_interpt_internal(l1, l2);
if (!PointerIsValid(result))
PG_RETURN_NULL();
result = (Point *) palloc(sizeof(Point));
if (!lseg_interpt_lseg(result, l1, l2))
PG_RETURN_NULL();
PG_RETURN_POINT_P(result);
}
......@@ -2340,15 +2244,9 @@ dist_pl(PG_FUNCTION_ARGS)
Point *pt = PG_GETARG_POINT_P(0);
LINE *line = PG_GETARG_LINE_P(1);
PG_RETURN_FLOAT8(dist_pl_internal(pt, line));
PG_RETURN_FLOAT8(line_closept_point(NULL, line, pt));
}
static double
dist_pl_internal(Point *pt, LINE *line)
{
return fabs((line->A * pt->x + line->B * pt->y + line->C) /
HYPOT(line->A, line->B));
}
/*
* Distance from a point to a lseg
......@@ -2359,60 +2257,7 @@ dist_ps(PG_FUNCTION_ARGS)
Point *pt = PG_GETARG_POINT_P(0);
LSEG *lseg = PG_GETARG_LSEG_P(1);
PG_RETURN_FLOAT8(dist_ps_internal(pt, lseg));
}
static double
dist_ps_internal(Point *pt, LSEG *lseg)
{
double m; /* slope of perp. */
LINE *ln;
double result,
tmpdist;
Point *ip;
/*
* 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)
m = (double) DBL_MAX; /* slope is infinite */
else
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
printf("dist_ps- line is A=%g B=%g C=%g from (point) slope (%f,%f) %g\n",
ln->A, ln->B, ln->C, pt->x, pt->y, m);
#endif
/*
* 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",
result, ip->x, ip->y);
#endif
}
else
{
/* 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)
result = tmpdist;
}
return result;
PG_RETURN_FLOAT8(lseg_closept_point(NULL, lseg, pt));
}
/*
......@@ -2429,21 +2274,10 @@ dist_ppath(PG_FUNCTION_ARGS)
int i;
LSEG lseg;
switch (path->npts)
{
case 0:
/* no points in path? then result is undefined... */
PG_RETURN_NULL();
case 1:
/* one point in path? then get distance between two points... */
result = point_dt(pt, &path->p[0]);
break;
default:
/* make sure the path makes sense... */
Assert(path->npts > 1);
Assert(path->npts > 0);
/*
* the distance from a point to a path is the smallest distance
* The distance from a point to a path is the smallest distance
* from the point to any of its constituent segments.
*/
for (i = 0; i < path->npts; i++)
......@@ -2456,19 +2290,18 @@ dist_ppath(PG_FUNCTION_ARGS)
{
if (!path->closed)
continue;
iprev = path->npts - 1; /* include the closure segment */
iprev = path->npts - 1; /* Include the closure segment */
}
statlseg_construct(&lseg, &path->p[iprev], &path->p[i]);
tmp = dist_ps_internal(pt, &lseg);
tmp = lseg_closept_point(NULL, &lseg, pt);
if (!have_min || tmp < result)
{
result = tmp;
have_min = true;
}
}
break;
}
PG_RETURN_FLOAT8(result);
}
......@@ -2480,15 +2313,8 @@ dist_pb(PG_FUNCTION_ARGS)
{
Point *pt = PG_GETARG_POINT_P(0);
BOX *box = PG_GETARG_BOX_P(1);
float8 result;
Point *near;
near = DatumGetPointP(DirectFunctionCall2(close_pb,
PointPGetDatum(pt),
BoxPGetDatum(box)));
result = point_dt(near, pt);
PG_RETURN_FLOAT8(result);
PG_RETURN_FLOAT8(box_closept_point(NULL, box, pt));
}
/*
......@@ -2502,12 +2328,12 @@ dist_sl(PG_FUNCTION_ARGS)
float8 result,
d2;
if (has_interpt_sl(lseg, line))
if (lseg_interpt_line(NULL, lseg, line))
result = 0.0;
else
{
result = dist_pl_internal(&lseg->p[0], line);
d2 = dist_pl_internal(&lseg->p[1], line);
result = line_closept_point(NULL, line, &lseg->p[0]);
d2 = line_closept_point(NULL, line, &lseg->p[1]);
/* XXX shouldn't we take the min not max? */
if (d2 > result)
result = d2;
......@@ -2524,17 +2350,8 @@ dist_sb(PG_FUNCTION_ARGS)
{
LSEG *lseg = PG_GETARG_LSEG_P(0);
BOX *box = PG_GETARG_BOX_P(1);
Point *tmp;
Datum result;
tmp = DatumGetPointP(DirectFunctionCall2(close_sb,
LsegPGetDatum(lseg),
BoxPGetDatum(box)));
result = DirectFunctionCall2(dist_pb,
PointPGetDatum(tmp),
BoxPGetDatum(box));
PG_RETURN_DATUM(result);
PG_RETURN_FLOAT8(box_closept_lseg(NULL, box, lseg));
}
/*
......@@ -2584,11 +2401,8 @@ dist_ppoly(PG_FUNCTION_ARGS)
{
Point *point = PG_GETARG_POINT_P(0);
POLYGON *poly = PG_GETARG_POLYGON_P(1);
float8 result;
result = dist_ppoly_internal(point, poly);
PG_RETURN_FLOAT8(result);
PG_RETURN_FLOAT8(dist_ppoly_internal(point, poly));
}
Datum
......@@ -2596,11 +2410,8 @@ dist_polyp(PG_FUNCTION_ARGS)
{
POLYGON *poly = PG_GETARG_POLYGON_P(0);
Point *point = PG_GETARG_POINT_P(1);
float8 result;
result = dist_ppoly_internal(point, poly);
PG_RETURN_FLOAT8(result);
PG_RETURN_FLOAT8(dist_ppoly_internal(point, poly));
}
static double
......@@ -2624,19 +2435,19 @@ dist_ppoly_internal(Point *pt, POLYGON *poly)
seg.p[0].y = poly->p[0].y;
seg.p[1].x = poly->p[poly->npts - 1].x;
seg.p[1].y = poly->p[poly->npts - 1].y;
result = dist_ps_internal(pt, &seg);
result = lseg_closept_point(NULL, &seg, pt);
#ifdef GEODEBUG
printf("dist_ppoly_internal- segment 0/n distance is %f\n", result);
#endif
/* check distances for other segments */
for (i = 0; (i < poly->npts - 1); i++)
for (i = 0; i < poly->npts - 1; i++)
{
seg.p[0].x = poly->p[i].x;
seg.p[0].y = poly->p[i].y;
seg.p[1].x = poly->p[i + 1].x;
seg.p[1].y = poly->p[i + 1].y;
d = dist_ps_internal(pt, &seg);
d = lseg_closept_point(NULL, &seg, pt);
#ifdef GEODEBUG
printf("dist_ppoly_internal- segment %d distance is %f\n", (i + 1), d);
#endif
......@@ -2655,49 +2466,51 @@ dist_ppoly_internal(Point *pt, POLYGON *poly)
* lines and boxes, since there are typically two.
*-------------------------------------------------------------------*/
/* Get intersection point of lseg and line; returns NULL if no intersection */
static Point *
interpt_sl(LSEG *lseg, LINE *line)
/*
* Check if the line segment intersects with the line
*
* This returns true if line segment intersects with line, false if they
* do not. This also sets the intersection point to *result, if it is not
* NULL.
*/
static bool
lseg_interpt_line(Point *result, LSEG *lseg, LINE *line)
{
Point interpt;
LINE tmp;
Point *p;
line_construct_pts(&tmp, &lseg->p[0], &lseg->p[1]);
p = line_interpt_internal(&tmp, line);
#ifdef GEODEBUG
printf("interpt_sl- segment is (%.*g %.*g) (%.*g %.*g)\n",
DBL_DIG, lseg->p[0].x, DBL_DIG, lseg->p[0].y, DBL_DIG, lseg->p[1].x, DBL_DIG, lseg->p[1].y);
printf("interpt_sl- segment becomes line A=%.*g B=%.*g C=%.*g\n",
DBL_DIG, tmp.A, DBL_DIG, tmp.B, DBL_DIG, tmp.C);
#endif
if (PointerIsValid(p))
{
#ifdef GEODEBUG
printf("interpt_sl- intersection point is (%.*g %.*g)\n", DBL_DIG, p->x, DBL_DIG, p->y);
#endif
if (on_ps_internal(p, lseg))
{
#ifdef GEODEBUG
printf("interpt_sl- intersection point is on segment\n");
#endif
}
else
p = NULL;
}
return p;
}
/*
* First, we promote the line segment to a line, because we know how
* to find the intersection point of two lines. If they don't have
* an intersection point, we are done.
*/
line_construct(&tmp, &lseg->p[0], lseg_sl(lseg));
if (!line_interpt_line(&interpt, &tmp, line))
return false;
/* variant: just indicate if intersection point exists */
static bool
has_interpt_sl(LSEG *lseg, LINE *line)
{
Point *tmp;
/*
* Then, we check whether the intersection point is actually on the line
* segment.
*/
if (!lseg_contain_point(lseg, &interpt))
return false;
if (result == NULL)
return true;
/*
* If there is an intersection, then check explicitly for matching
* endpoints since there may be rounding effects with annoying LSB
* residue.
*/
if (point_eq_point(&lseg->p[0], &interpt))
*result = lseg->p[0];
else if (point_eq_point(&lseg->p[1], &interpt))
*result = lseg->p[1];
else
*result = interpt;
tmp = interpt_sl(lseg, line);
if (tmp)
return true;
return false;
}
/*---------------------------------------------------------------------
......@@ -2705,277 +2518,217 @@ has_interpt_sl(LSEG *lseg, LINE *line)
* Point of closest proximity between objects.
*-------------------------------------------------------------------*/
/* close_pl -
/*
* The intersection point of a perpendicular of the line
* through the point.
*
* This sets the closest point to the *result if it is not NULL and returns
* the distance to the closest point.
*/
static float8
line_closept_point(Point *result, LINE *line, Point *point)
{
bool retval;
Point closept;
LINE tmp;
/* We drop a perpendicular to find the intersection point. */
line_construct(&tmp, point, line_invsl(line));
retval = line_interpt_line(&closept, line, &tmp);
Assert(retval); /* XXX: We need something better. */
if (result != NULL)
*result = closept;
/* Then we calculate the distance between the points. */
return point_dt(&closept, point);
}
Datum
close_pl(PG_FUNCTION_ARGS)
{
Point *pt = PG_GETARG_POINT_P(0);
LINE *line = PG_GETARG_LINE_P(1);
Point *result;
LINE *tmp;
double invm;
result = (Point *) palloc(sizeof(Point));
if (FPzero(line->B)) /* vertical? */
{
result->x = line->C;
result->y = pt->y;
PG_RETURN_POINT_P(result);
}
if (FPzero(line->A)) /* horizontal? */
{
result->x = pt->x;
result->y = line->C;
PG_RETURN_POINT_P(result);
}
/* drop a perpendicular and find the intersection point */
line_closept_point(result, line, pt);
/* invert and flip the sign on the slope to get a perpendicular */
invm = line->B / line->A;
tmp = line_construct_pm(pt, invm);
result = line_interpt_internal(tmp, line);
Assert(result != NULL);
PG_RETURN_POINT_P(result);
}
/* close_ps()
/*
* Closest point on line segment to specified point.
* Take the closest endpoint if the point is left, right,
* above, or below the segment, otherwise find the intersection
* point of the segment and its perpendicular through the point.
*
* Some tricky code here, relying on boolean expressions
* evaluating to only zero or one to use as an array index.
* bug fixes by gthaker@atl.lmco.com; May 1, 1998
* This sets the closest point to the *result if it is not NULL and returns
* the distance to the closest point.
*/
Datum
close_ps(PG_FUNCTION_ARGS)
static float8
lseg_closept_point(Point *result, LSEG *lseg, Point *pt)
{
Point *pt = PG_GETARG_POINT_P(0);
LSEG *lseg = PG_GETARG_LSEG_P(1);
Point *result = NULL;
LINE *tmp;
double invm;
int xh,
yh;
#ifdef GEODEBUG
printf("close_sp:pt->x %f pt->y %f\nlseg(0).x %f lseg(0).y %f lseg(1).x %f lseg(1).y %f\n",
pt->x, pt->y, lseg->p[0].x, lseg->p[0].y,
lseg->p[1].x, lseg->p[1].y);
#endif
Point closept;
LINE tmp;
/* xh (or yh) is the index of upper x( or y) end point of lseg */
/* !xh (or !yh) is the index of lower x( or y) end point of lseg */
xh = lseg->p[0].x < lseg->p[1].x;
yh = lseg->p[0].y < lseg->p[1].y;
/*
* To find the closest point, we draw a perpendicular line from the point
* to the line segment.
*/
line_construct(&tmp, pt, point_invsl(&lseg->p[0], &lseg->p[1]));
lseg_closept_line(&closept, lseg, &tmp);
if (FPeq(lseg->p[0].x, lseg->p[1].x)) /* vertical? */
{
#ifdef GEODEBUG
printf("close_ps- segment is vertical\n");
#endif
/* first check if point is below or above the entire lseg. */
if (pt->y < lseg->p[!yh].y)
result = point_copy(&lseg->p[!yh]); /* below the lseg */
else if (pt->y > lseg->p[yh].y)
result = point_copy(&lseg->p[yh]); /* above the lseg */
if (result != NULL)
PG_RETURN_POINT_P(result);
*result = closept;
/* point lines along (to left or right) of the vertical lseg. */
return point_dt(&closept, pt);
}
result = (Point *) palloc(sizeof(Point));
result->x = lseg->p[0].x;
result->y = pt->y;
PG_RETURN_POINT_P(result);
}
else if (FPeq(lseg->p[0].y, lseg->p[1].y)) /* horizontal? */
{
#ifdef GEODEBUG
printf("close_ps- segment is horizontal\n");
#endif
/* first check if point is left or right of the entire lseg. */
if (pt->x < lseg->p[!xh].x)
result = point_copy(&lseg->p[!xh]); /* left of the lseg */
else if (pt->x > lseg->p[xh].x)
result = point_copy(&lseg->p[xh]); /* right of the lseg */
if (result != NULL)
PG_RETURN_POINT_P(result);
Datum
close_ps(PG_FUNCTION_ARGS)
{
Point *pt = PG_GETARG_POINT_P(0);
LSEG *lseg = PG_GETARG_LSEG_P(1);
Point *result;
/* point lines along (at top or below) the horiz. lseg. */
result = (Point *) palloc(sizeof(Point));
result->x = pt->x;
result->y = lseg->p[0].y;
PG_RETURN_POINT_P(result);
}
/*
* vert. and horiz. cases are down, now check if the closest point is one
* of the end points or someplace on the lseg.
*/
invm = -1.0 / point_sl(&(lseg->p[0]), &(lseg->p[1]));
tmp = line_construct_pm(&lseg->p[!yh], invm); /* lower edge of the
* "band" */
if (pt->y < (tmp->A * pt->x + tmp->C))
{ /* we are below the lower edge */
result = point_copy(&lseg->p[!yh]); /* below the lseg, take lower end
* pt */
#ifdef GEODEBUG
printf("close_ps below: tmp A %f B %f C %f\n",
tmp->A, tmp->B, tmp->C);
#endif
PG_RETURN_POINT_P(result);
}
tmp = line_construct_pm(&lseg->p[yh], invm); /* upper edge of the
* "band" */
if (pt->y > (tmp->A * pt->x + tmp->C))
{ /* we are below the lower edge */
result = point_copy(&lseg->p[yh]); /* above the lseg, take higher end
* pt */
#ifdef GEODEBUG
printf("close_ps above: tmp A %f B %f C %f\n",
tmp->A, tmp->B, tmp->C);
#endif
PG_RETURN_POINT_P(result);
}
/*
* at this point the "normal" from point will hit lseg. The closest point
* will be somewhere on the lseg
*/
tmp = line_construct_pm(pt, invm);
#ifdef GEODEBUG
printf("close_ps- tmp A %f B %f C %f\n",
tmp->A, tmp->B, tmp->C);
#endif
result = interpt_sl(lseg, tmp);
/*
* ordinarily we should always find an intersection point, but that could
* fail in the presence of NaN coordinates, and perhaps even from simple
* roundoff issues. Return a SQL NULL if so.
*/
if (result == NULL)
if (isnan(lseg_closept_point(result, lseg, pt)))
PG_RETURN_NULL();
#ifdef GEODEBUG
printf("close_ps- result.x %f result.y %f\n", result->x, result->y);
#endif
PG_RETURN_POINT_P(result);
}
/* close_lseg()
* Closest point to l1 on l2.
/*
* Closest point on line segment to line segment
*
* This sets the closest point to the *result if it is not NULL and returns
* the distance to the closest point.
*/
Datum
close_lseg(PG_FUNCTION_ARGS)
static float8
lseg_closept_lseg(Point *result, LSEG *l1, LSEG *l2)
{
LSEG *l1 = PG_GETARG_LSEG_P(0);
LSEG *l2 = PG_GETARG_LSEG_P(1);
Point *result = NULL;
Point point;
double dist;
double d;
d = dist_ps_internal(&l1->p[0], l2);
d = lseg_closept_point(NULL, l1, &l2->p[0]);
dist = d;
memcpy(&point, &l1->p[0], sizeof(Point));
if (result != NULL)
*result = l2->p[0];
if ((d = dist_ps_internal(&l1->p[1], l2)) < dist)
d = lseg_closept_point(NULL, l1, &l2->p[1]);
if (d < dist)
{
dist = d;
memcpy(&point, &l1->p[1], sizeof(Point));
if (result != NULL)
*result = l2->p[1];
}
if (dist_ps_internal(&l2->p[0], l1) < dist)
{
result = DatumGetPointP(DirectFunctionCall2(close_ps,
PointPGetDatum(&l2->p[0]),
LsegPGetDatum(l1)));
memcpy(&point, result, sizeof(Point));
result = DatumGetPointP(DirectFunctionCall2(close_ps,
PointPGetDatum(&point),
LsegPGetDatum(l2)));
}
if (lseg_closept_point(&point, l2, &l1->p[0]) < dist)
d = lseg_closept_point(result, l1, &point);
if (dist_ps_internal(&l2->p[1], l1) < dist)
{
result = DatumGetPointP(DirectFunctionCall2(close_ps,
PointPGetDatum(&l2->p[1]),
LsegPGetDatum(l1)));
memcpy(&point, result, sizeof(Point));
result = DatumGetPointP(DirectFunctionCall2(close_ps,
PointPGetDatum(&point),
LsegPGetDatum(l2)));
}
if (lseg_closept_point(&point, l2, &l1->p[1]) < dist)
d = lseg_closept_point(result, l1, &point);
if (result == NULL)
result = point_copy(&point);
if (d < dist)
dist = d;
return dist;
}
Datum
close_lseg(PG_FUNCTION_ARGS)
{
LSEG *l1 = PG_GETARG_LSEG_P(0);
LSEG *l2 = PG_GETARG_LSEG_P(1);
Point *result;
result = (Point *) palloc(sizeof(Point));
lseg_closept_lseg(result, l2, l1);
PG_RETURN_POINT_P(result);
}
/* close_pb()
/*
* Closest point on or in box to specified point.
*
* This sets the closest point to the *result if it is not NULL and returns
* the distance to the closest point.
*/
Datum
close_pb(PG_FUNCTION_ARGS)
static float8
box_closept_point(Point *result, BOX *box, Point *pt)
{
Point *pt = PG_GETARG_POINT_P(0);
BOX *box = PG_GETARG_BOX_P(1);
LSEG lseg,
seg;
Point point;
LSEG lseg;
Point point,
closept;
double dist,
d;
if (DatumGetBool(DirectFunctionCall2(on_pb,
PointPGetDatum(pt),
BoxPGetDatum(box))))
PG_RETURN_POINT_P(pt);
if (box_contain_point(box, pt))
{
if (result != NULL)
*result = *pt;
return 0.0;
}
/* pairwise check lseg distances */
point.x = box->low.x;
point.y = box->high.y;
statlseg_construct(&lseg, &box->low, &point);
dist = dist_ps_internal(pt, &lseg);
dist = lseg_closept_point(result, &lseg, pt);
statlseg_construct(&seg, &box->high, &point);
if ((d = dist_ps_internal(pt, &seg)) < dist)
statlseg_construct(&lseg, &box->high, &point);
d = lseg_closept_point(&closept, &lseg, pt);
if (d < dist)
{
dist = d;
memcpy(&lseg, &seg, sizeof(lseg));
if (result != NULL)
*result = closept;
}
point.x = box->high.x;
point.y = box->low.y;
statlseg_construct(&seg, &box->low, &point);
if ((d = dist_ps_internal(pt, &seg)) < dist)
statlseg_construct(&lseg, &box->low, &point);
d = lseg_closept_point(&closept, &lseg, pt);
if (d < dist)
{
dist = d;
memcpy(&lseg, &seg, sizeof(lseg));
if (result != NULL)
*result = closept;
}
statlseg_construct(&seg, &box->high, &point);
if ((d = dist_ps_internal(pt, &seg)) < dist)
statlseg_construct(&lseg, &box->high, &point);
d = lseg_closept_point(&closept, &lseg, pt);
if (d < dist)
{
dist = d;
memcpy(&lseg, &seg, sizeof(lseg));
if (result != NULL)
*result = closept;
}
PG_RETURN_DATUM(DirectFunctionCall2(close_ps,
PointPGetDatum(pt),
LsegPGetDatum(&lseg)));
return dist;
}
Datum
close_pb(PG_FUNCTION_ARGS)
{
Point *pt = PG_GETARG_POINT_P(0);
BOX *box = PG_GETARG_BOX_P(1);
Point *result;
result = (Point *) palloc(sizeof(Point));
box_closept_point(result, box, pt);
PG_RETURN_POINT_P(result);
}
/* close_sl()
* Closest point on line to line segment.
*
......@@ -2995,16 +2748,17 @@ close_sl(PG_FUNCTION_ARGS)
float8 d1,
d2;
result = interpt_sl(lseg, line);
if (result)
result = (Point *) palloc(sizeof(Point));
if (lseg_interpt_line(result, lseg, line))
PG_RETURN_POINT_P(result);
d1 = dist_pl_internal(&lseg->p[0], line);
d2 = dist_pl_internal(&lseg->p[1], line);
d1 = line_closept_point(NULL, line, &lseg->p[0]);
d2 = line_closept_point(NULL, line, &lseg->p[1]);
if (d1 < d2)
result = point_copy(&lseg->p[0]);
*result = lseg->p[0];
else
result = point_copy(&lseg->p[1]);
*result = lseg->p[1];
PG_RETURN_POINT_P(result);
#endif
......@@ -3016,92 +2770,132 @@ close_sl(PG_FUNCTION_ARGS)
PG_RETURN_NULL();
}
/* close_ls()
/*
* Closest point on line segment to line.
*
* This sets the closest point to the *result if it is not NULL and returns
* the distance to the closest point.
*
* NOTE: When the lines are parallel, endpoints of one of the line segment
* are FPeq(), in presence of NaN or Infinitive coordinates, or perhaps =
* even because of simple roundoff issues, there may not be a single closest
* point. We are likely to set the result to the second endpoint in these
* cases.
*/
static float8
lseg_closept_line(Point *result, LSEG *lseg, LINE *line)
{
float8 dist1,
dist2;
if (lseg_interpt_line(result, lseg, line))
return 0.0;
dist1 = line_closept_point(NULL, line, &lseg->p[0]);
dist2 = line_closept_point(NULL, line, &lseg->p[1]);
if (dist1 < dist2)
{
if (result != NULL)
*result = lseg->p[0];
return dist1;
}
else
{
if (result != NULL)
*result = lseg->p[1];
return dist2;
}
}
Datum
close_ls(PG_FUNCTION_ARGS)
{
LINE *line = PG_GETARG_LINE_P(0);
LSEG *lseg = PG_GETARG_LSEG_P(1);
Point *result;
float8 d1,
d2;
result = interpt_sl(lseg, line);
if (result)
PG_RETURN_POINT_P(result);
result = (Point *) palloc(sizeof(Point));
d1 = dist_pl_internal(&lseg->p[0], line);
d2 = dist_pl_internal(&lseg->p[1], line);
if (d1 < d2)
result = point_copy(&lseg->p[0]);
else
result = point_copy(&lseg->p[1]);
lseg_closept_line(result, lseg, line);
PG_RETURN_POINT_P(result);
}
/* close_sb()
/*
* Closest point on or in box to line segment.
*
* This sets the closest point to the *result if it is not NULL and returns
* the distance to the closest point.
*/
Datum
close_sb(PG_FUNCTION_ARGS)
static float8
box_closept_lseg(Point *result, BOX *box, LSEG *lseg)
{
LSEG *lseg = PG_GETARG_LSEG_P(0);
BOX *box = PG_GETARG_BOX_P(1);
Point point;
LSEG bseg,
seg;
Point point,
closept;
LSEG bseg;
double dist,
d;
/* segment intersects box? then just return closest point to center */
if (DatumGetBool(DirectFunctionCall2(inter_sb,
LsegPGetDatum(lseg),
BoxPGetDatum(box))))
{
box_cn(&point, box);
PG_RETURN_DATUM(DirectFunctionCall2(close_ps,
PointPGetDatum(&point),
LsegPGetDatum(lseg)));
}
if (box_interpt_lseg(result, box, lseg))
return 0.0;
/* pairwise check lseg distances */
point.x = box->low.x;
point.y = box->high.y;
statlseg_construct(&bseg, &box->low, &point);
dist = lseg_dt(lseg, &bseg);
dist = lseg_closept_lseg(result, &bseg, lseg);
statlseg_construct(&seg, &box->high, &point);
if ((d = lseg_dt(lseg, &seg)) < dist)
statlseg_construct(&bseg, &box->high, &point);
d = lseg_closept_lseg(&closept, &bseg, lseg);
if (d < dist)
{
dist = d;
memcpy(&bseg, &seg, sizeof(bseg));
if (result != NULL)
*result = closept;
}
point.x = box->high.x;
point.y = box->low.y;
statlseg_construct(&seg, &box->low, &point);
if ((d = lseg_dt(lseg, &seg)) < dist)
statlseg_construct(&bseg, &box->low, &point);
d = lseg_closept_lseg(&closept, &bseg, lseg);
if (d < dist)
{
dist = d;
memcpy(&bseg, &seg, sizeof(bseg));
if (result != NULL)
*result = closept;
}
statlseg_construct(&seg, &box->high, &point);
if ((d = lseg_dt(lseg, &seg)) < dist)
statlseg_construct(&bseg, &box->high, &point);
d = lseg_closept_lseg(&closept, &bseg, lseg);
if (d < dist)
{
dist = d;
memcpy(&bseg, &seg, sizeof(bseg));
if (result != NULL)
*result = closept;
}
/* OK, we now have the closest line segment on the box boundary */
PG_RETURN_DATUM(DirectFunctionCall2(close_lseg,
LsegPGetDatum(lseg),
LsegPGetDatum(&bseg)));
return dist;
}
Datum
close_sb(PG_FUNCTION_ARGS)
{
LSEG *lseg = PG_GETARG_LSEG_P(0);
BOX *box = PG_GETARG_BOX_P(1);
Point *result;
result = (Point *) palloc(sizeof(Point));
box_closept_lseg(result, box, lseg);
PG_RETURN_POINT_P(result);
}
Datum
close_lb(PG_FUNCTION_ARGS)
{
......@@ -3123,37 +2917,55 @@ close_lb(PG_FUNCTION_ARGS)
* Whether one object lies completely within another.
*-------------------------------------------------------------------*/
/* on_pl -
/*
* Does the point satisfy the equation?
*/
static bool
line_contain_point(LINE *line, Point *point)
{
return FPzero(line->A * point->x + line->B * point->y + line->C);
}
Datum
on_pl(PG_FUNCTION_ARGS)
{
Point *pt = PG_GETARG_POINT_P(0);
LINE *line = PG_GETARG_LINE_P(1);
PG_RETURN_BOOL(FPzero(line->A * pt->x + line->B * pt->y + line->C));
PG_RETURN_BOOL(line_contain_point(line, pt));
}
/* on_ps -
/*
* Determine colinearity by detecting a triangle inequality.
* This algorithm seems to behave nicely even with lsb residues - tgl 1997-07-09
*/
static bool
lseg_contain_point(LSEG *lseg, Point *pt)
{
return FPeq(point_dt(pt, &lseg->p[0]) +
point_dt(pt, &lseg->p[1]),
point_dt(&lseg->p[0], &lseg->p[1]));
}
Datum
on_ps(PG_FUNCTION_ARGS)
{
Point *pt = PG_GETARG_POINT_P(0);
LSEG *lseg = PG_GETARG_LSEG_P(1);
PG_RETURN_BOOL(on_ps_internal(pt, lseg));
PG_RETURN_BOOL(lseg_contain_point(lseg, pt));
}
/*
* Check whether the point is in the box or on its border
*/
static bool
on_ps_internal(Point *pt, LSEG *lseg)
box_contain_point(BOX *box, Point *point)
{
return FPeq(point_dt(pt, &lseg->p[0]) + point_dt(pt, &lseg->p[1]),
point_dt(&lseg->p[0], &lseg->p[1]));
return box->high.x >= point->x && box->low.x <= point->x &&
box->high.y >= point->y && box->low.y <= point-> y;
}
Datum
......@@ -3162,8 +2974,7 @@ on_pb(PG_FUNCTION_ARGS)
Point *pt = PG_GETARG_POINT_P(0);
BOX *box = PG_GETARG_BOX_P(1);
PG_RETURN_BOOL(pt->x <= box->high.x && pt->x >= box->low.x &&
pt->y <= box->high.y && pt->y >= box->low.y);
PG_RETURN_BOOL(box_contain_point(box, pt));
}
Datum
......@@ -3172,8 +2983,7 @@ box_contain_pt(PG_FUNCTION_ARGS)
BOX *box = PG_GETARG_BOX_P(0);
Point *pt = PG_GETARG_POINT_P(1);
PG_RETURN_BOOL(pt->x <= box->high.x && pt->x >= box->low.x &&
pt->y <= box->high.y && pt->y >= box->low.y);
PG_RETURN_BOOL(box_contain_point(box, pt));
}
/* on_ppath -
......@@ -3217,18 +3027,33 @@ on_ppath(PG_FUNCTION_ARGS)
PG_RETURN_BOOL(point_inside(pt, path->npts, path->p) != 0);
}
/*
* Check whether the line segment is on the line or close enough
*
* It is, if both of its points are on the line or close enough.
*/
Datum
on_sl(PG_FUNCTION_ARGS)
{
LSEG *lseg = PG_GETARG_LSEG_P(0);
LINE *line = PG_GETARG_LINE_P(1);
PG_RETURN_BOOL(DatumGetBool(DirectFunctionCall2(on_pl,
PointPGetDatum(&lseg->p[0]),
LinePGetDatum(line))) &&
DatumGetBool(DirectFunctionCall2(on_pl,
PointPGetDatum(&lseg->p[1]),
LinePGetDatum(line))));
PG_RETURN_BOOL(line_contain_point(line, &lseg->p[0]) &&
line_contain_point(line, &lseg->p[1]));
}
/*
* Check whether the line segment is in the box or on its border
*
* It is, if both of its points are in the box or on its border.
*/
static bool
box_contain_lseg(BOX *box, LSEG *lseg)
{
return box_contain_point(box, &lseg->p[0]) &&
box_contain_point(box, &lseg->p[1]);
}
Datum
......@@ -3237,12 +3062,7 @@ on_sb(PG_FUNCTION_ARGS)
LSEG *lseg = PG_GETARG_LSEG_P(0);
BOX *box = PG_GETARG_BOX_P(1);
PG_RETURN_BOOL(DatumGetBool(DirectFunctionCall2(on_pb,
PointPGetDatum(&lseg->p[0]),
BoxPGetDatum(box))) &&
DatumGetBool(DirectFunctionCall2(on_pb,
PointPGetDatum(&lseg->p[1]),
BoxPGetDatum(box))));
PG_RETURN_BOOL(box_contain_lseg(box, lseg));
}
/*---------------------------------------------------------------------
......@@ -3256,24 +3076,28 @@ inter_sl(PG_FUNCTION_ARGS)
LSEG *lseg = PG_GETARG_LSEG_P(0);
LINE *line = PG_GETARG_LINE_P(1);
PG_RETURN_BOOL(has_interpt_sl(lseg, line));
PG_RETURN_BOOL(lseg_interpt_line(NULL, lseg, line));
}
/* inter_sb()
/*
* Do line segment and box intersect?
*
* Segment completely inside box counts as intersection.
* If you want only segments crossing box boundaries,
* try converting box to path first.
*
* This function also sets the *result to the closest point on the line
* segment to the center of the box when they overlap and the result is
* not NULL. It is somewhat arbitrary, but maybe the best we can do as
* there are typically two points they intersect.
*
* Optimize for non-intersection by checking for box intersection first.
* - thomas 1998-01-30
*/
Datum
inter_sb(PG_FUNCTION_ARGS)
static bool
box_interpt_lseg(Point *result, BOX *box, LSEG *lseg)
{
LSEG *lseg = PG_GETARG_LSEG_P(0);
BOX *box = PG_GETARG_BOX_P(1);
BOX lbox;
LSEG bseg;
Point point;
......@@ -3285,42 +3109,54 @@ inter_sb(PG_FUNCTION_ARGS)
/* nothing close to overlap? then not going to intersect */
if (!box_ov(&lbox, box))
PG_RETURN_BOOL(false);
return false;
if (result != NULL)
{
box_cn(&point, box);
lseg_closept_point(result, lseg, &point);
}
/* an endpoint of segment is inside box? then clearly intersects */
if (DatumGetBool(DirectFunctionCall2(on_pb,
PointPGetDatum(&lseg->p[0]),
BoxPGetDatum(box))) ||
DatumGetBool(DirectFunctionCall2(on_pb,
PointPGetDatum(&lseg->p[1]),
BoxPGetDatum(box))))
PG_RETURN_BOOL(true);
if (box_contain_point(box, &lseg->p[0]) ||
box_contain_point(box, &lseg->p[1]))
return true;
/* pairwise check lseg intersections */
point.x = box->low.x;
point.y = box->high.y;
statlseg_construct(&bseg, &box->low, &point);
if (lseg_intersect_internal(&bseg, lseg))
PG_RETURN_BOOL(true);
if (lseg_interpt_lseg(NULL, &bseg, lseg))
return true;
statlseg_construct(&bseg, &box->high, &point);
if (lseg_intersect_internal(&bseg, lseg))
PG_RETURN_BOOL(true);
if (lseg_interpt_lseg(NULL, &bseg, lseg))
return true;
point.x = box->high.x;
point.y = box->low.y;
statlseg_construct(&bseg, &box->low, &point);
if (lseg_intersect_internal(&bseg, lseg))
PG_RETURN_BOOL(true);
if (lseg_interpt_lseg(NULL, &bseg, lseg))
return true;
statlseg_construct(&bseg, &box->high, &point);
if (lseg_intersect_internal(&bseg, lseg))
PG_RETURN_BOOL(true);
if (lseg_interpt_lseg(NULL, &bseg, lseg))
return true;
/* if we dropped through, no two segs intersected */
PG_RETURN_BOOL(false);
return false;
}
Datum
inter_sb(PG_FUNCTION_ARGS)
{
LSEG *lseg = PG_GETARG_LSEG_P(0);
BOX *box = PG_GETARG_BOX_P(1);
PG_RETURN_BOOL(box_interpt_lseg(NULL, box, lseg));
}
/* inter_lb()
* Do line and box intersect?
*/
......@@ -3339,22 +3175,22 @@ inter_lb(PG_FUNCTION_ARGS)
p2.x = box->low.x;
p2.y = box->high.y;
statlseg_construct(&bseg, &p1, &p2);
if (has_interpt_sl(&bseg, line))
if (lseg_interpt_line(NULL, &bseg, line))
PG_RETURN_BOOL(true);
p1.x = box->high.x;
p1.y = box->high.y;
statlseg_construct(&bseg, &p1, &p2);
if (has_interpt_sl(&bseg, line))
if (lseg_interpt_line(NULL, &bseg, line))
PG_RETURN_BOOL(true);
p2.x = box->high.x;
p2.y = box->low.y;
statlseg_construct(&bseg, &p1, &p2);
if (has_interpt_sl(&bseg, line))
if (lseg_interpt_line(NULL, &bseg, line))
PG_RETURN_BOOL(true);
p1.x = box->low.x;
p1.y = box->low.y;
statlseg_construct(&bseg, &p1, &p2);
if (has_interpt_sl(&bseg, line))
if (lseg_interpt_line(NULL, &bseg, line))
PG_RETURN_BOOL(true);
/* if we dropped through, no intersection */
......@@ -3381,9 +3217,9 @@ make_bound_box(POLYGON *poly)
x2,
y2;
if (poly->npts > 0)
{
x2 = x1 = poly->p[0].x;
Assert(poly->npts > 0);
x1 = x2 = poly->p[0].x;
y2 = y1 = poly->p[0].y;
for (i = 1; i < poly->npts; i++)
{
......@@ -3397,12 +3233,10 @@ make_bound_box(POLYGON *poly)
y2 = poly->p[i].y;
}
box_fill(&(poly->boundbox), x1, x2, y1, y2);
}
else
ereport(ERROR,
(errcode(ERRCODE_INVALID_PARAMETER_VALUE),
errmsg("cannot create bounding box for empty polygon")));
poly->boundbox.low.x = x1;
poly->boundbox.high.x = x2;
poly->boundbox.low.y = y1;
poly->boundbox.high.y = y2;
}
/*------------------------------------------------------------------
......@@ -3746,9 +3580,10 @@ poly_overlap(PG_FUNCTION_ARGS)
POLYGON *polyb = PG_GETARG_POLYGON_P(1);
bool result;
Assert(polya->npts > 0 && polyb->npts > 0);
/* Quick check by bounding box */
result = (polya->npts > 0 && polyb->npts > 0 &&
box_ov(&polya->boundbox, &polyb->boundbox)) ? true : false;
result = box_ov(&polya->boundbox, &polyb->boundbox);
/*
* Brute-force algorithm - try to find intersected edges, if so then
......@@ -3766,7 +3601,7 @@ poly_overlap(PG_FUNCTION_ARGS)
sa.p[0] = polya->p[polya->npts - 1];
result = false;
for (ia = 0; ia < polya->npts && result == false; ia++)
for (ia = 0; ia < polya->npts && !result; ia++)
{
/* Second point of polya's edge is a current one */
sa.p[1] = polya->p[ia];
......@@ -3774,10 +3609,10 @@ poly_overlap(PG_FUNCTION_ARGS)
/* Init first of polyb's edge with last point */
sb.p[0] = polyb->p[polyb->npts - 1];
for (ib = 0; ib < polyb->npts && result == false; ib++)
for (ib = 0; ib < polyb->npts && !result; ib++)
{
sb.p[1] = polyb->p[ib];
result = lseg_intersect_internal(&sa, &sb);
result = lseg_interpt_lseg(NULL, &sa, &sb);
sb.p[0] = sb.p[1];
}
......@@ -3787,10 +3622,9 @@ poly_overlap(PG_FUNCTION_ARGS)
sa.p[0] = sa.p[1];
}
if (result == false)
if (!result)
{
result = (point_inside(polya->p, polyb->npts, polyb->p)
||
result = (point_inside(polya->p, polyb->npts, polyb->p) ||
point_inside(polyb->p, polya->npts, polya->p));
}
}
......@@ -3824,22 +3658,21 @@ touched_lseg_inside_poly(Point *a, Point *b, LSEG *s, POLYGON *poly, int start)
t.p[0] = *a;
t.p[1] = *b;
#define POINTEQ(pt1, pt2) (FPeq((pt1)->x, (pt2)->x) && FPeq((pt1)->y, (pt2)->y))
if (POINTEQ(a, s->p))
if (point_eq_point(a, s->p))
{
if (on_ps_internal(s->p + 1, &t))
if (lseg_contain_point(&t, s->p + 1))
return lseg_inside_poly(b, s->p + 1, poly, start);
}
else if (POINTEQ(a, s->p + 1))
else if (point_eq_point(a, s->p + 1))
{
if (on_ps_internal(s->p, &t))
if (lseg_contain_point(&t, s->p))
return lseg_inside_poly(b, s->p, poly, start);
}
else if (on_ps_internal(s->p, &t))
else if (lseg_contain_point(&t, s->p))
{
return lseg_inside_poly(b, s->p, poly, start);
}
else if (on_ps_internal(s->p + 1, &t))
else if (lseg_contain_point(&t, s->p + 1))
{
return lseg_inside_poly(b, s->p + 1, poly, start);
}
......@@ -3867,36 +3700,35 @@ lseg_inside_poly(Point *a, Point *b, POLYGON *poly, int start)
for (i = start; i < poly->npts && res; i++)
{
Point *interpt;
Point interpt;
CHECK_FOR_INTERRUPTS();
s.p[1] = poly->p[i];
if (on_ps_internal(t.p, &s))
if (lseg_contain_point(&s, t.p))
{
if (on_ps_internal(t.p + 1, &s))
if (lseg_contain_point(&s, t.p + 1))
return true; /* t is contained by s */
/* Y-cross */
res = touched_lseg_inside_poly(t.p, t.p + 1, &s, poly, i + 1);
}
else if (on_ps_internal(t.p + 1, &s))
else if (lseg_contain_point(&s, t.p + 1))
{
/* Y-cross */
res = touched_lseg_inside_poly(t.p + 1, t.p, &s, poly, i + 1);
}
else if ((interpt = lseg_interpt_internal(&t, &s)) != NULL)
else if (lseg_interpt_lseg(&interpt, &t, &s))
{
/*
* segments are X-crossing, go to check each subsegment
*/
intersection = true;
res = lseg_inside_poly(t.p, interpt, poly, i + 1);
res = lseg_inside_poly(t.p, &interpt, poly, i + 1);
if (res)
res = lseg_inside_poly(t.p + 1, interpt, poly, i + 1);
pfree(interpt);
res = lseg_inside_poly(t.p + 1, &interpt, poly, i + 1);
}
s.p[0] = s.p[1];
......@@ -3922,38 +3754,41 @@ lseg_inside_poly(Point *a, Point *b, POLYGON *poly, int start)
/*-----------------------------------------------------------------
* Determine if polygon A contains polygon B.
*-----------------------------------------------------------------*/
Datum
poly_contain(PG_FUNCTION_ARGS)
static bool
poly_contain_poly(POLYGON *polya, POLYGON *polyb)
{
POLYGON *polya = PG_GETARG_POLYGON_P(0);
POLYGON *polyb = PG_GETARG_POLYGON_P(1);
bool result;
int i;
LSEG s;
Assert(polya->npts > 0 && polyb->npts > 0);
/*
* Quick check to see if bounding box is contained.
*/
if (polya->npts > 0 && polyb->npts > 0 &&
DatumGetBool(DirectFunctionCall2(box_contain,
BoxPGetDatum(&polya->boundbox),
BoxPGetDatum(&polyb->boundbox))))
{
int i;
LSEG s;
if (!box_contain_box(&polya->boundbox, &polyb->boundbox))
return false;
s.p[0] = polyb->p[polyb->npts - 1];
result = true;
for (i = 0; i < polyb->npts && result; i++)
for (i = 0; i < polyb->npts; i++)
{
s.p[1] = polyb->p[i];
result = lseg_inside_poly(s.p, s.p + 1, polya, 0);
if (!lseg_inside_poly(s.p, s.p + 1, polya, 0))
return false;
s.p[0] = s.p[1];
}
}
else
{
result = false;
}
return true;
}
Datum
poly_contain(PG_FUNCTION_ARGS)
{
POLYGON *polya = PG_GETARG_POLYGON_P(0);
POLYGON *polyb = PG_GETARG_POLYGON_P(1);
bool result;
result = poly_contain_poly(polya, polyb);
/*
* Avoid leaking memory for toasted inputs ... needed for rtree indexes
......@@ -3971,11 +3806,20 @@ poly_contain(PG_FUNCTION_ARGS)
Datum
poly_contained(PG_FUNCTION_ARGS)
{
Datum polya = PG_GETARG_DATUM(0);
Datum polyb = PG_GETARG_DATUM(1);
POLYGON *polya = PG_GETARG_POLYGON_P(0);
POLYGON *polyb = PG_GETARG_POLYGON_P(1);
bool result;
/* Just switch the arguments and pass it off to poly_contain */
PG_RETURN_DATUM(DirectFunctionCall2(poly_contain, polyb, polya));
result = poly_contain_poly(polyb, polya);
/*
* Avoid leaking memory for toasted inputs ... needed for rtree indexes
*/
PG_FREE_IF_COPY(polya, 0);
PG_FREE_IF_COPY(polyb, 1);
PG_RETURN_BOOL(result);
}
......@@ -4025,8 +3869,22 @@ construct_point(PG_FUNCTION_ARGS)
{
float8 x = PG_GETARG_FLOAT8(0);
float8 y = PG_GETARG_FLOAT8(1);
Point *result;
result = (Point *) palloc(sizeof(Point));
point_construct(result, x, y);
PG_RETURN_POINT_P(point_construct(x, y));
PG_RETURN_POINT_P(result);
}
static inline void
point_add_point(Point *result, Point *pt1, Point *pt2)
{
point_construct(result,
pt1->x + pt2->x,
pt1->y + pt2->y);
}
Datum
......@@ -4038,12 +3896,20 @@ point_add(PG_FUNCTION_ARGS)
result = (Point *) palloc(sizeof(Point));
result->x = (p1->x + p2->x);
result->y = (p1->y + p2->y);
point_add_point(result, p1, p2);
PG_RETURN_POINT_P(result);
}
static inline void
point_sub_point(Point *result, Point *pt1, Point *pt2)
{
point_construct(result,
pt1->x - pt2->x,
pt1->y - pt2->y);
}
Datum
point_sub(PG_FUNCTION_ARGS)
{
......@@ -4053,12 +3919,20 @@ point_sub(PG_FUNCTION_ARGS)
result = (Point *) palloc(sizeof(Point));
result->x = (p1->x - p2->x);
result->y = (p1->y - p2->y);
point_sub_point(result, p1, p2);
PG_RETURN_POINT_P(result);
}
static inline void
point_mul_point(Point *result, Point *pt1, Point *pt2)
{
point_construct(result,
(pt1->x * pt2->x) - (pt1->y * pt2->y),
(pt1->x * pt2->y) + (pt1->y * pt2->x));
}
Datum
point_mul(PG_FUNCTION_ARGS)
{
......@@ -4068,31 +3942,39 @@ point_mul(PG_FUNCTION_ARGS)
result = (Point *) palloc(sizeof(Point));
result->x = (p1->x * p2->x) - (p1->y * p2->y);
result->y = (p1->x * p2->y) + (p1->y * p2->x);
point_mul_point(result, p1, p2);
PG_RETURN_POINT_P(result);
}
Datum
point_div(PG_FUNCTION_ARGS)
static inline void
point_div_point(Point *result, Point *pt1, Point *pt2)
{
Point *p1 = PG_GETARG_POINT_P(0);
Point *p2 = PG_GETARG_POINT_P(1);
Point *result;
double div;
result = (Point *) palloc(sizeof(Point));
div = (p2->x * p2->x) + (p2->y * p2->y);
div = (pt2->x * pt2->x) + (pt2->y * pt2->y);
if (div == 0.0)
ereport(ERROR,
(errcode(ERRCODE_DIVISION_BY_ZERO),
errmsg("division by zero")));
result->x = ((p1->x * p2->x) + (p1->y * p2->y)) / div;
result->y = ((p2->x * p1->y) - (p2->y * p1->x)) / div;
point_construct(result,
((pt1->x * pt2->x) + (pt1->y * pt2->y)) / div,
((pt2->x * pt1->y) - (pt2->y * pt1->x)) / div);
}
Datum
point_div(PG_FUNCTION_ARGS)
{
Point *p1 = PG_GETARG_POINT_P(0);
Point *p2 = PG_GETARG_POINT_P(1);
Point *result;
result = (Point *) palloc(sizeof(Point));
point_div_point(result, p1, p2);
PG_RETURN_POINT_P(result);
}
......@@ -4109,8 +3991,13 @@ points_box(PG_FUNCTION_ARGS)
{
Point *p1 = PG_GETARG_POINT_P(0);
Point *p2 = PG_GETARG_POINT_P(1);
BOX *result;
result = (BOX *) palloc(sizeof(BOX));
box_construct(result, p1, p2);
PG_RETURN_BOX_P(box_construct(p1->x, p2->x, p1->y, p2->y));
PG_RETURN_BOX_P(result);
}
Datum
......@@ -4118,11 +4005,14 @@ box_add(PG_FUNCTION_ARGS)
{
BOX *box = PG_GETARG_BOX_P(0);
Point *p = PG_GETARG_POINT_P(1);
BOX *result;
result = (BOX *) palloc(sizeof(BOX));
point_add_point(&result->high, &box->high, p);
point_add_point(&result->low, &box->low, p);
PG_RETURN_BOX_P(box_construct((box->high.x + p->x),
(box->low.x + p->x),
(box->high.y + p->y),
(box->low.y + p->y)));
PG_RETURN_BOX_P(result);
}
Datum
......@@ -4130,11 +4020,14 @@ box_sub(PG_FUNCTION_ARGS)
{
BOX *box = PG_GETARG_BOX_P(0);
Point *p = PG_GETARG_POINT_P(1);
BOX *result;
result = (BOX *) palloc(sizeof(BOX));
PG_RETURN_BOX_P(box_construct((box->high.x - p->x),
(box->low.x - p->x),
(box->high.y - p->y),
(box->low.y - p->y)));
point_sub_point(&result->high, &box->high, p);
point_sub_point(&result->low, &box->low, p);
PG_RETURN_BOX_P(result);
}
Datum
......@@ -4143,17 +4036,15 @@ box_mul(PG_FUNCTION_ARGS)
BOX *box = PG_GETARG_BOX_P(0);
Point *p = PG_GETARG_POINT_P(1);
BOX *result;
Point *high,
*low;
Point high,
low;
high = DatumGetPointP(DirectFunctionCall2(point_mul,
PointPGetDatum(&box->high),
PointPGetDatum(p)));
low = DatumGetPointP(DirectFunctionCall2(point_mul,
PointPGetDatum(&box->low),
PointPGetDatum(p)));
result = (BOX *) palloc(sizeof(BOX));
result = box_construct(high->x, low->x, high->y, low->y);
point_mul_point(&high, &box->high, p);
point_mul_point(&low, &box->low, p);
box_construct(result, &high, &low);
PG_RETURN_BOX_P(result);
}
......@@ -4164,17 +4055,15 @@ box_div(PG_FUNCTION_ARGS)
BOX *box = PG_GETARG_BOX_P(0);
Point *p = PG_GETARG_POINT_P(1);
BOX *result;
Point *high,
*low;
Point high,
low;
result = (BOX *) palloc(sizeof(BOX));
high = DatumGetPointP(DirectFunctionCall2(point_div,
PointPGetDatum(&box->high),
PointPGetDatum(p)));
low = DatumGetPointP(DirectFunctionCall2(point_div,
PointPGetDatum(&box->low),
PointPGetDatum(p)));
point_div_point(&high, &box->high, p);
point_div_point(&low, &box->low, p);
result = box_construct(high->x, low->x, high->y, low->y);
box_construct(result, &high, &low);
PG_RETURN_BOX_P(result);
}
......@@ -4284,10 +4173,7 @@ path_add_pt(PG_FUNCTION_ARGS)
int i;
for (i = 0; i < path->npts; i++)
{
path->p[i].x += point->x;
path->p[i].y += point->y;
}
point_add_point(&path->p[i], &path->p[i], point);
PG_RETURN_PATH_P(path);
}
......@@ -4300,10 +4186,7 @@ path_sub_pt(PG_FUNCTION_ARGS)
int i;
for (i = 0; i < path->npts; i++)
{
path->p[i].x -= point->x;
path->p[i].y -= point->y;
}
point_sub_point(&path->p[i], &path->p[i], point);
PG_RETURN_PATH_P(path);
}
......@@ -4316,17 +4199,10 @@ path_mul_pt(PG_FUNCTION_ARGS)
{
PATH *path = PG_GETARG_PATH_P_COPY(0);
Point *point = PG_GETARG_POINT_P(1);
Point *p;
int i;
for (i = 0; i < path->npts; i++)
{
p = DatumGetPointP(DirectFunctionCall2(point_mul,
PointPGetDatum(&path->p[i]),
PointPGetDatum(point)));
path->p[i].x = p->x;
path->p[i].y = p->y;
}
point_mul_point(&path->p[i], &path->p[i], point);
PG_RETURN_PATH_P(path);
}
......@@ -4336,17 +4212,10 @@ path_div_pt(PG_FUNCTION_ARGS)
{
PATH *path = PG_GETARG_PATH_P_COPY(0);
Point *point = PG_GETARG_POINT_P(1);
Point *p;
int i;
for (i = 0; i < path->npts; i++)
{
p = DatumGetPointP(DirectFunctionCall2(point_div,
PointPGetDatum(&path->p[i]),
PointPGetDatum(point)));
path->p[i].x = p->x;
path->p[i].y = p->y;
}
point_div_point(&path->p[i], &path->p[i], point);
PG_RETURN_PATH_P(path);
}
......@@ -4421,15 +4290,15 @@ Datum
poly_center(PG_FUNCTION_ARGS)
{
POLYGON *poly = PG_GETARG_POLYGON_P(0);
Datum result;
CIRCLE *circle;
Point *result;
CIRCLE circle;
circle = DatumGetCircleP(DirectFunctionCall1(poly_circle,
PolygonPGetDatum(poly)));
result = DirectFunctionCall1(circle_center,
CirclePGetDatum(circle));
result = (Point *) palloc(sizeof(Point));
poly_to_circle(&circle, poly);
*result = circle.center;
PG_RETURN_DATUM(result);
PG_RETURN_POINT_P(result);
}
......@@ -4439,10 +4308,8 @@ poly_box(PG_FUNCTION_ARGS)
POLYGON *poly = PG_GETARG_POLYGON_P(0);
BOX *box;
if (poly->npts < 1)
PG_RETURN_NULL();
box = box_copy(&poly->boundbox);
box = (BOX *) palloc(sizeof(BOX));
*box = poly->boundbox;
PG_RETURN_BOX_P(box);
}
......@@ -4474,8 +4341,7 @@ box_poly(PG_FUNCTION_ARGS)
poly->p[3].x = box->high.x;
poly->p[3].y = box->low.y;
box_fill(&poly->boundbox, box->high.x, box->low.x,
box->high.y, box->low.y);
box_construct(&poly->boundbox, &box->high, &box->low);
PG_RETURN_POLYGON_P(poly);
}
......@@ -4564,8 +4430,7 @@ circle_in(PG_FUNCTION_ARGS)
while (depth > 0)
{
if ((*s == RDELIM)
|| ((*s == RDELIM_C) && (depth == 1)))
if ((*s == RDELIM) || ((*s == RDELIM_C) && (depth == 1)))
{
depth--;
s++;
......@@ -4663,8 +4528,7 @@ circle_same(PG_FUNCTION_ARGS)
CIRCLE *circle2 = PG_GETARG_CIRCLE_P(1);
PG_RETURN_BOOL(FPeq(circle1->radius, circle2->radius) &&
FPeq(circle1->center.x, circle2->center.x) &&
FPeq(circle1->center.y, circle2->center.y));
point_eq_point(&circle1->center, &circle2->center));
}
/* circle_overlap - does circle1 overlap circle2?
......@@ -4737,7 +4601,8 @@ circle_contained(PG_FUNCTION_ARGS)
CIRCLE *circle1 = PG_GETARG_CIRCLE_P(0);
CIRCLE *circle2 = PG_GETARG_CIRCLE_P(1);
PG_RETURN_BOOL(FPle((point_dt(&circle1->center, &circle2->center) + circle1->radius), circle2->radius));
PG_RETURN_BOOL(FPle(point_dt(&circle1->center, &circle2->center),
circle2->radius - circle1->radius));
}
/* circle_contain - does circle1 contain circle2?
......@@ -4748,7 +4613,8 @@ circle_contain(PG_FUNCTION_ARGS)
CIRCLE *circle1 = PG_GETARG_CIRCLE_P(0);
CIRCLE *circle2 = PG_GETARG_CIRCLE_P(1);
PG_RETURN_BOOL(FPle((point_dt(&circle1->center, &circle2->center) + circle2->radius), circle1->radius));
PG_RETURN_BOOL(FPle(point_dt(&circle1->center, &circle2->center),
circle1->radius - circle2->radius));
}
......@@ -4865,20 +4731,6 @@ circle_ge(PG_FUNCTION_ARGS)
* "Arithmetic" operators on circles.
*---------------------------------------------------------*/
static CIRCLE *
circle_copy(CIRCLE *circle)
{
CIRCLE *result;
if (!PointerIsValid(circle))
return NULL;
result = (CIRCLE *) palloc(sizeof(CIRCLE));
memcpy((char *) result, (char *) circle, sizeof(CIRCLE));
return result;
}
/* circle_add_pt()
* Translation operator.
*/
......@@ -4889,10 +4741,10 @@ circle_add_pt(PG_FUNCTION_ARGS)
Point *point = PG_GETARG_POINT_P(1);
CIRCLE *result;
result = circle_copy(circle);
result = (CIRCLE *) palloc(sizeof(CIRCLE));
result->center.x += point->x;
result->center.y += point->y;
point_add_point(&result->center, &circle->center, point);
result->radius = circle->radius;
PG_RETURN_CIRCLE_P(result);
}
......@@ -4904,10 +4756,10 @@ circle_sub_pt(PG_FUNCTION_ARGS)
Point *point = PG_GETARG_POINT_P(1);
CIRCLE *result;
result = circle_copy(circle);
result = (CIRCLE *) palloc(sizeof(CIRCLE));
result->center.x -= point->x;
result->center.y -= point->y;
point_sub_point(&result->center, &circle->center, point);
result->radius = circle->radius;
PG_RETURN_CIRCLE_P(result);
}
......@@ -4922,16 +4774,11 @@ circle_mul_pt(PG_FUNCTION_ARGS)
CIRCLE *circle = PG_GETARG_CIRCLE_P(0);
Point *point = PG_GETARG_POINT_P(1);
CIRCLE *result;
Point *p;
result = circle_copy(circle);
result = (CIRCLE *) palloc(sizeof(CIRCLE));
p = DatumGetPointP(DirectFunctionCall2(point_mul,
PointPGetDatum(&circle->center),
PointPGetDatum(point)));
result->center.x = p->x;
result->center.y = p->y;
result->radius *= HYPOT(point->x, point->y);
point_mul_point(&result->center, &circle->center, point);
result->radius = circle->radius * HYPOT(point->x, point->y);
PG_RETURN_CIRCLE_P(result);
}
......@@ -4942,16 +4789,11 @@ circle_div_pt(PG_FUNCTION_ARGS)
CIRCLE *circle = PG_GETARG_CIRCLE_P(0);
Point *point = PG_GETARG_POINT_P(1);
CIRCLE *result;
Point *p;
result = circle_copy(circle);
result = (CIRCLE *) palloc(sizeof(CIRCLE));
p = DatumGetPointP(DirectFunctionCall2(point_div,
PointPGetDatum(&circle->center),
PointPGetDatum(point)));
result->center.x = p->x;
result->center.y = p->y;
result->radius /= HYPOT(point->x, point->y);
point_div_point(&result->center, &circle->center, point);
result->radius = circle->radius / HYPOT(point->x, point->y);
PG_RETURN_CIRCLE_P(result);
}
......@@ -5000,8 +4842,8 @@ circle_distance(PG_FUNCTION_ARGS)
CIRCLE *circle2 = PG_GETARG_CIRCLE_P(1);
float8 result;
result = point_dt(&circle1->center, &circle2->center)
- (circle1->radius + circle2->radius);
result = point_dt(&circle1->center, &circle2->center) -
(circle1->radius + circle2->radius);
if (result < 0)
result = 0;
PG_RETURN_FLOAT8(result);
......@@ -5197,42 +5039,46 @@ circle_poly(PG_FUNCTION_ARGS)
PG_RETURN_POLYGON_P(poly);
}
/* poly_circle - convert polygon to circle
/*
* Convert polygon to circle
*
* The result must be preallocated.
*
* XXX This algorithm should use weighted means of line segments
* rather than straight average values of points - tgl 97/01/21.
*/
Datum
poly_circle(PG_FUNCTION_ARGS)
static void
poly_to_circle(CIRCLE *result, POLYGON *poly)
{
POLYGON *poly = PG_GETARG_POLYGON_P(0);
CIRCLE *circle;
int i;
if (poly->npts < 1)
ereport(ERROR,
(errcode(ERRCODE_INVALID_PARAMETER_VALUE),
errmsg("cannot convert empty polygon to circle")));
Assert(poly->npts > 0);
circle = (CIRCLE *) palloc(sizeof(CIRCLE));
circle->center.x = 0;
circle->center.y = 0;
circle->radius = 0;
result->center.x = 0;
result->center.y = 0;
result->radius = 0;
for (i = 0; i < poly->npts; i++)
{
circle->center.x += poly->p[i].x;
circle->center.y += poly->p[i].y;
}
circle->center.x /= poly->npts;
circle->center.y /= poly->npts;
point_add_point(&result->center, &result->center, &poly->p[i]);
result->center.x /= poly->npts;
result->center.y /= poly->npts;
for (i = 0; i < poly->npts; i++)
circle->radius += point_dt(&poly->p[i], &circle->center);
circle->radius /= poly->npts;
result->radius += point_dt(&poly->p[i], &result->center);
result->radius /= poly->npts;
}
PG_RETURN_CIRCLE_P(circle);
Datum
poly_circle(PG_FUNCTION_ARGS)
{
POLYGON *poly = PG_GETARG_POLYGON_P(0);
CIRCLE *result;
result = (CIRCLE *) palloc(sizeof(CIRCLE));
poly_to_circle(result, poly);
PG_RETURN_CIRCLE_P(result);
}
......@@ -5268,8 +5114,7 @@ point_inside(Point *p, int npts, Point *plist)
int cross,
total_cross = 0;
if (npts <= 0)
return 0;
Assert(npts > 0);
/* compute first polygon point relative to single point */
x0 = plist[0].x - p->x;
......@@ -5378,8 +5223,7 @@ plist_same(int npts, Point *p1, Point *p2)
/* find match for first point */
for (i = 0; i < npts; i++)
{
if ((FPeq(p2[i].x, p1[0].x))
&& (FPeq(p2[i].y, p1[0].y)))
if (point_eq_point(&p2[i], &p1[0]))
{
/* match found? then look forward through remaining points */
......@@ -5387,8 +5231,7 @@ plist_same(int npts, Point *p1, Point *p2)
{
if (j >= npts)
j = 0;
if ((!FPeq(p2[j].x, p1[ii].x))
|| (!FPeq(p2[j].y, p1[ii].y)))
if (!point_eq_point(&p2[j], &p1[ii]))
{
#ifdef GEODEBUG
printf("plist_same- %d failed forward match with %d\n", j, ii);
......@@ -5407,8 +5250,7 @@ plist_same(int npts, Point *p1, Point *p2)
{
if (j < 0)
j = (npts - 1);
if ((!FPeq(p2[j].x, p1[ii].x))
|| (!FPeq(p2[j].y, p1[ii].y)))
if (!point_eq_point(&p2[j], &p1[ii]))
{
#ifdef GEODEBUG
printf("plist_same- %d failed reverse match with %d\n", j, ii);
......
src/backend/utils/adt/geo_spgist.c
View file @ a7dc63d9
......@@ -793,7 +793,8 @@ spg_poly_quad_compress(PG_FUNCTION_ARGS)
POLYGON *polygon = PG_GETARG_POLYGON_P(0);
BOX *box;
box = box_copy(&polygon->boundbox);
box = (BOX *) palloc(sizeof(BOX));
*box = polygon->boundbox;
PG_RETURN_BOX_P(box);
}
src/include/utils/geo_decls.h
View file @ a7dc63d9
......@@ -178,10 +178,6 @@ typedef struct
* in geo_ops.c
*/
/* private routines */
extern double point_dt(Point *pt1, Point *pt2);
extern double point_sl(Point *pt1, Point *pt2);
extern double pg_hypot(double x, double y);
extern BOX *box_copy(BOX *box);
#endif /* GEO_DECLS_H */
src/test/regress/regress.c
View file @ a7dc63d9
......@@ -176,8 +176,13 @@ pt_in_widget(PG_FUNCTION_ARGS)
{
Point *point = PG_GETARG_POINT_P(0);
WIDGET *widget = (WIDGET *) PG_GETARG_POINTER(1);
float8 distance;
PG_RETURN_BOOL(point_dt(point, &widget->center) < widget->radius);
distance = DatumGetFloat8(DirectFunctionCall2(point_distance,
PointPGetDatum(point),
PointPGetDatum(&widget->center)));
PG_RETURN_BOOL(distance < widget->radius);
}
PG_FUNCTION_INFO_V1(reverse_name);
......
Markdown is supported
0% or
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment