[postgis-tickets] r17773 - Support 2D TINs in ST_3DIntersection
Darafei
komzpa at gmail.com
Fri Aug 23 02:56:38 PDT 2019
Author: komzpa
Date: 2019-08-23 02:56:37 -0700 (Fri, 23 Aug 2019)
New Revision: 17773
Modified:
trunk/NEWS
trunk/liblwgeom/measures.c
trunk/liblwgeom/measures.h
trunk/regress/core/measures.sql
trunk/regress/core/measures_expected
Log:
Support 2D TINs in ST_3DIntersection
Brings TIN into dist2d calculation. Removes homebrew atan2.
Closes #4328
Closes https://github.com/postgis/postgis/pull/461
Modified: trunk/NEWS
===================================================================
--- trunk/NEWS 2019-08-23 06:37:52 UTC (rev 17772)
+++ trunk/NEWS 2019-08-23 09:56:37 UTC (rev 17773)
@@ -81,8 +81,7 @@
- #4258, #4278 Enhance postgis functions ST_Area, ST_Distance, ST_Intersection,
ST_Difference, ST_Union, ST_Intersects, ST_3DIntersects, ST_3DDistance
to handle types (Solid PolyhedralSurface, TINS) that were only available in SFCGAL.
- Remove postgis.backend switch.
- (Darafei Praliaskouski)
+ Remove postgis.backend switch. (Darafei Praliaskouski)
- #4230, SP-GiST and GiST support for ND box operators overlaps, contains,
within, equals (Esteban Zimányi and Arthur Lesuisse from Université
Libre de Bruxelles (ULB), Darafei Praliaskouski)
@@ -206,6 +205,7 @@
ST_AsGeoJSON sets SRID in JSON output if it differs from 4326.
(Darafei Praliaskouski)
- #3979, postgis_sfcgal_noop() round trip function (Lucas C. Villa Real)
+ - #4328, ST_3DIntersects for 2D TINs. (Darafei Praliaskouski)
* Fixes *
- #4342, Move deprecated functions into legacy.sql file
Modified: trunk/liblwgeom/measures.c
===================================================================
--- trunk/liblwgeom/measures.c 2019-08-23 06:37:52 UTC (rev 17772)
+++ trunk/liblwgeom/measures.c 2019-08-23 09:56:37 UTC (rev 17773)
@@ -21,10 +21,10 @@
* Copyright 2001-2006 Refractions Research Inc.
* Copyright 2010 Nicklas Avén
* Copyright 2012 Paul Ramsey
+ * Copyright 2019 Darafei Praliaskouski
*
**********************************************************************/
-
#include <string.h>
#include <stdlib.h>
@@ -31,8 +31,6 @@
#include "measures.h"
#include "lwgeom_log.h"
-
-
/*------------------------------------------------------------------------------------------------------------
Initializing functions
The functions starting the distance-calculation processses
@@ -41,28 +39,27 @@
LWGEOM *
lwgeom_closest_line(const LWGEOM *lw1, const LWGEOM *lw2)
{
- return lw_dist2d_distanceline(lw1, lw2, lw1->srid, DIST_MIN);
+ return lw_dist2d_distanceline(lw1, lw2, lw1->srid, DIST_MIN);
}
LWGEOM *
lwgeom_furthest_line(const LWGEOM *lw1, const LWGEOM *lw2)
{
- return lw_dist2d_distanceline(lw1, lw2, lw1->srid, DIST_MAX);
+ return lw_dist2d_distanceline(lw1, lw2, lw1->srid, DIST_MAX);
}
LWGEOM *
lwgeom_closest_point(const LWGEOM *lw1, const LWGEOM *lw2)
{
- return lw_dist2d_distancepoint(lw1, lw2, lw1->srid, DIST_MIN);
+ return lw_dist2d_distancepoint(lw1, lw2, lw1->srid, DIST_MIN);
}
LWGEOM *
lwgeom_furthest_point(const LWGEOM *lw1, const LWGEOM *lw2)
{
- return lw_dist2d_distancepoint(lw1, lw2, lw1->srid, DIST_MAX);
+ return lw_dist2d_distancepoint(lw1, lw2, lw1->srid, DIST_MAX);
}
-
void
lw_dist2d_distpts_init(DISTPTS *dl, int mode)
{
@@ -71,7 +68,7 @@
dl->p2.x = dl->p2.y = 0.0;
dl->mode = mode;
dl->tolerance = 0.0;
- if ( mode == DIST_MIN )
+ if (mode == DIST_MIN)
dl->distance = FLT_MAX;
else
dl->distance = -1 * FLT_MAX;
@@ -83,9 +80,9 @@
LWGEOM *
lw_dist2d_distanceline(const LWGEOM *lw1, const LWGEOM *lw2, int32_t srid, int mode)
{
- double x1,x2,y1,y2;
+ double x1, x2, y1, y2;
- double initdistance = ( mode == DIST_MIN ? FLT_MAX : -1.0);
+ double initdistance = (mode == DIST_MIN ? FLT_MAX : -1.0);
DISTPTS thedl;
LWPOINT *lwpoints[2];
LWGEOM *result;
@@ -96,7 +93,7 @@
LWDEBUG(2, "lw_dist2d_distanceline is called");
- if (!lw_dist2d_comp( lw1,lw2,&thedl))
+ if (!lw_dist2d_comp(lw1, lw2, &thedl))
{
/*should never get here. all cases ought to be error handled earlier*/
lwerror("Some unspecified error.");
@@ -111,10 +108,10 @@
}
else
{
- x1=thedl.p1.x;
- y1=thedl.p1.y;
- x2=thedl.p2.x;
- y2=thedl.p2.y;
+ x1 = thedl.p1.x;
+ y1 = thedl.p1.y;
+ x2 = thedl.p2.x;
+ y2 = thedl.p2.y;
lwpoints[0] = lwpoint_make2d(srid, x1, y1);
lwpoints[1] = lwpoint_make2d(srid, x2, y2);
@@ -130,18 +127,18 @@
LWGEOM *
lw_dist2d_distancepoint(const LWGEOM *lw1, const LWGEOM *lw2, int32_t srid, int mode)
{
- double x,y;
+ double x, y;
DISTPTS thedl;
double initdistance = FLT_MAX;
LWGEOM *result;
thedl.mode = mode;
- thedl.distance= initdistance;
+ thedl.distance = initdistance;
thedl.tolerance = 0;
LWDEBUG(2, "lw_dist2d_distancepoint is called");
- if (!lw_dist2d_comp( lw1,lw2,&thedl))
+ if (!lw_dist2d_comp(lw1, lw2, &thedl))
{
/*should never get here. all cases ought to be error handled earlier*/
lwerror("Some unspecified error.");
@@ -154,14 +151,13 @@
}
else
{
- x=thedl.p1.x;
- y=thedl.p1.y;
+ x = thedl.p1.x;
+ y = thedl.p1.y;
result = (LWGEOM *)lwpoint_make2d(srid, x, y);
}
return result;
}
-
/**
Function initializing max distance calculation
*/
@@ -170,7 +166,7 @@
{
LWDEBUG(2, "lwgeom_maxdistance2d is called");
- return lwgeom_maxdistance2d_tolerance( lw1, lw2, 0.0 );
+ return lwgeom_maxdistance2d_tolerance(lw1, lw2, 0.0);
}
/**
@@ -184,12 +180,11 @@
DISTPTS thedl;
LWDEBUG(2, "lwgeom_maxdistance2d_tolerance is called");
thedl.mode = DIST_MAX;
- thedl.distance= -1;
+ thedl.distance = -1;
thedl.tolerance = tolerance;
- if (lw_dist2d_comp( lw1,lw2,&thedl))
- {
+ if (lw_dist2d_comp(lw1, lw2, &thedl))
return thedl.distance;
- }
+
/*should never get here. all cases ought to be error handled earlier*/
lwerror("Some unspecified error.");
return -1;
@@ -201,7 +196,7 @@
double
lwgeom_mindistance2d(const LWGEOM *lw1, const LWGEOM *lw2)
{
- return lwgeom_mindistance2d_tolerance( lw1, lw2, 0.0 );
+ return lwgeom_mindistance2d_tolerance(lw1, lw2, 0.0);
}
/**
@@ -214,18 +209,15 @@
DISTPTS thedl;
LWDEBUG(2, "lwgeom_mindistance2d_tolerance is called");
thedl.mode = DIST_MIN;
- thedl.distance= FLT_MAX;
+ thedl.distance = FLT_MAX;
thedl.tolerance = tolerance;
- if (lw_dist2d_comp( lw1,lw2,&thedl))
- {
+ if (lw_dist2d_comp(lw1, lw2, &thedl))
return thedl.distance;
- }
/*should never get here. all cases ought to be error handled earlier*/
lwerror("Some unspecified error.");
return FLT_MAX;
}
-
/*------------------------------------------------------------------------------------------------------------
End of Initializing functions
--------------------------------------------------------------------------------------------------------------*/
@@ -241,10 +233,8 @@
bboxes we will use anyway.
*/
int
-lw_dist2d_comp(const LWGEOM *lw1,const LWGEOM *lw2, DISTPTS *dl)
+lw_dist2d_comp(const LWGEOM *lw1, const LWGEOM *lw2, DISTPTS *dl)
{
- LWDEBUG(2, "lw_dist2d_comp is called");
-
return lw_dist2d_recursive(lw1, lw2, dl);
}
@@ -251,9 +241,10 @@
static int
lw_dist2d_is_collection(const LWGEOM *g)
{
-
+ /* Differs from lwgeom_is_collection by not treating CURVEPOLYGON as collection */
switch (g->type)
{
+ case TINTYPE:
case MULTIPOINTTYPE:
case MULTILINETYPE:
case MULTIPOLYGONTYPE:
@@ -273,11 +264,12 @@
/**
This is a recursive function delivering every possible combination of subgeometries
*/
-int lw_dist2d_recursive(const LWGEOM *lwg1, const LWGEOM *lwg2, DISTPTS *dl)
+int
+lw_dist2d_recursive(const LWGEOM *lwg1, const LWGEOM *lwg2, DISTPTS *dl)
{
int i, j;
- int n1=1;
- int n2=1;
+ int n1 = 1;
+ int n2 = 1;
LWGEOM *g1 = NULL;
LWGEOM *g2 = NULL;
LWCOLLECTION *c1 = NULL;
@@ -298,66 +290,63 @@
n2 = c2->ngeoms;
}
- for ( i = 0; i < n1; i++ )
+ for (i = 0; i < n1; i++)
{
if (lw_dist2d_is_collection(lwg1))
- {
g1 = c1->geoms[i];
- }
else
- {
- g1 = (LWGEOM*)lwg1;
- }
+ g1 = (LWGEOM *)lwg1;
- if (lwgeom_is_empty(g1)) return LW_TRUE;
+ if (lwgeom_is_empty(g1))
+ return LW_TRUE;
if (lw_dist2d_is_collection(g1))
{
LWDEBUG(3, "Found collection inside first geometry collection, recursing");
- if (!lw_dist2d_recursive(g1, lwg2, dl)) return LW_FALSE;
+ if (!lw_dist2d_recursive(g1, lwg2, dl))
+ return LW_FALSE;
continue;
}
- for ( j = 0; j < n2; j++ )
+ for (j = 0; j < n2; j++)
{
if (lw_dist2d_is_collection(lwg2))
- {
g2 = c2->geoms[j];
- }
else
- {
- g2 = (LWGEOM*)lwg2;
- }
+ g2 = (LWGEOM *)lwg2;
+
if (lw_dist2d_is_collection(g2))
{
LWDEBUG(3, "Found collection inside second geometry collection, recursing");
- if (!lw_dist2d_recursive(g1, g2, dl)) return LW_FALSE;
+ if (!lw_dist2d_recursive(g1, g2, dl))
+ return LW_FALSE;
continue;
}
- if ( ! g1->bbox )
- {
+ if (!g1->bbox)
lwgeom_add_bbox(g1);
- }
- if ( ! g2->bbox )
- {
+
+ if (!g2->bbox)
lwgeom_add_bbox(g2);
- }
- /*If one of geometries is empty, return. True here only means continue searching. False would have stopped the process*/
- if (lwgeom_is_empty(g1)||lwgeom_is_empty(g2)) return LW_TRUE;
+ /* If one of geometries is empty, return. True here only means continue searching. False would
+ * have stopped the process*/
+ if (lwgeom_is_empty(g1) || lwgeom_is_empty(g2))
+ return LW_TRUE;
- if ( (dl->mode != DIST_MAX) &&
- (! lw_dist2d_check_overlap(g1, g2)) &&
- (g1->type == LINETYPE || g1->type == POLYGONTYPE) &&
- (g2->type == LINETYPE || g2->type == POLYGONTYPE) )
+ if ((dl->mode != DIST_MAX) && (!lw_dist2d_check_overlap(g1, g2)) &&
+ (g1->type == LINETYPE || g1->type == POLYGONTYPE || g1->type == TRIANGLETYPE) &&
+ (g2->type == LINETYPE || g2->type == POLYGONTYPE || g2->type == TRIANGLETYPE))
{
- if (!lw_dist2d_distribute_fast(g1, g2, dl)) return LW_FALSE;
+ if (!lw_dist2d_distribute_fast(g1, g2, dl))
+ return LW_FALSE;
}
else
{
- if (!lw_dist2d_distribute_bruteforce(g1, g2, dl)) return LW_FALSE;
- if (dl->distance<=dl->tolerance && dl->mode == DIST_MIN) return LW_TRUE; /*just a check if the answer is already given*/
+ if (!lw_dist2d_distribute_bruteforce(g1, g2, dl))
+ return LW_FALSE;
+ if (dl->distance <= dl->tolerance && dl->mode == DIST_MIN)
+ return LW_TRUE; /*just a check if the answer is already given*/
}
}
}
@@ -364,151 +353,179 @@
return LW_TRUE;
}
-
int
-lw_dist2d_distribute_bruteforce(const LWGEOM *lwg1,const LWGEOM *lwg2, DISTPTS *dl)
+lw_dist2d_distribute_bruteforce(const LWGEOM *lwg1, const LWGEOM *lwg2, DISTPTS *dl)
{
- int t1 = lwg1->type;
- int t2 = lwg2->type;
+ int t1 = lwg1->type;
+ int t2 = lwg2->type;
- switch ( t1 )
+ switch (t1)
{
+ case POINTTYPE:
+ {
+ dl->twisted = 1;
+ switch (t2)
+ {
case POINTTYPE:
+ return lw_dist2d_point_point((LWPOINT *)lwg1, (LWPOINT *)lwg2, dl);
+ case LINETYPE:
+ return lw_dist2d_point_line((LWPOINT *)lwg1, (LWLINE *)lwg2, dl);
+ case TRIANGLETYPE:
+ return lw_dist2d_point_tri((LWPOINT *)lwg1, (LWTRIANGLE *)lwg2, dl);
+ case POLYGONTYPE:
+ return lw_dist2d_point_poly((LWPOINT *)lwg1, (LWPOLY *)lwg2, dl);
+ case CIRCSTRINGTYPE:
+ return lw_dist2d_point_circstring((LWPOINT *)lwg1, (LWCIRCSTRING *)lwg2, dl);
+ case CURVEPOLYTYPE:
+ return lw_dist2d_point_curvepoly((LWPOINT *)lwg1, (LWCURVEPOLY *)lwg2, dl);
+ default:
+ lwerror("%s: Unsupported geometry type: %s", __func__, lwtype_name(t2));
+ return LW_FALSE;
+ }
+ }
+ case LINETYPE:
+ {
+ dl->twisted = 1;
+ switch (t2)
{
- dl->twisted = 1;
- switch ( t2 )
- {
- case POINTTYPE:
- return lw_dist2d_point_point((LWPOINT *)lwg1, (LWPOINT *)lwg2, dl);
- case LINETYPE:
- return lw_dist2d_point_line((LWPOINT *)lwg1, (LWLINE *)lwg2, dl);
- case POLYGONTYPE:
- return lw_dist2d_point_poly((LWPOINT *)lwg1, (LWPOLY *)lwg2, dl);
- case CIRCSTRINGTYPE:
- return lw_dist2d_point_circstring((LWPOINT *)lwg1, (LWCIRCSTRING *)lwg2, dl);
- case CURVEPOLYTYPE:
- return lw_dist2d_point_curvepoly((LWPOINT *)lwg1, (LWCURVEPOLY *)lwg2, dl);
- default:
- lwerror("%s: Unsupported geometry type: %s", __func__, lwtype_name(t2));
- return LW_FALSE;
- }
+ case POINTTYPE:
+ dl->twisted = -1;
+ return lw_dist2d_point_line((LWPOINT *)lwg2, (LWLINE *)lwg1, dl);
+ case LINETYPE:
+ return lw_dist2d_line_line((LWLINE *)lwg1, (LWLINE *)lwg2, dl);
+ case TRIANGLETYPE:
+ return lw_dist2d_line_tri((LWLINE *)lwg1, (LWTRIANGLE *)lwg2, dl);
+ case POLYGONTYPE:
+ return lw_dist2d_line_poly((LWLINE *)lwg1, (LWPOLY *)lwg2, dl);
+ case CIRCSTRINGTYPE:
+ return lw_dist2d_line_circstring((LWLINE *)lwg1, (LWCIRCSTRING *)lwg2, dl);
+ case CURVEPOLYTYPE:
+ return lw_dist2d_line_curvepoly((LWLINE *)lwg1, (LWCURVEPOLY *)lwg2, dl);
+ default:
+ lwerror("%s: Unsupported geometry type: %s", __func__, lwtype_name(t2));
+ return LW_FALSE;
}
+ }
+ case TRIANGLETYPE:
+ {
+ dl->twisted = 1;
+ switch (t2)
+ {
+ case POINTTYPE:
+ dl->twisted = -1;
+ return lw_dist2d_point_tri((LWPOINT *)lwg2, (LWTRIANGLE *)lwg1, dl);
case LINETYPE:
+ dl->twisted = -1;
+ return lw_dist2d_line_tri((LWLINE *)lwg2, (LWTRIANGLE *)lwg1, dl);
+ case TRIANGLETYPE:
+ return lw_dist2d_tri_tri((LWTRIANGLE *)lwg1, (LWTRIANGLE *)lwg2, dl);
+ case POLYGONTYPE:
+ return lw_dist2d_tri_poly((LWTRIANGLE *)lwg1, (LWPOLY *)lwg2, dl);
+ case CIRCSTRINGTYPE:
+ return lw_dist2d_tri_circstring((LWTRIANGLE *)lwg1, (LWCIRCSTRING *)lwg2, dl);
+ case CURVEPOLYTYPE:
+ return lw_dist2d_tri_curvepoly((LWTRIANGLE *)lwg1, (LWCURVEPOLY *)lwg2, dl);
+ default:
+ lwerror("%s: Unsupported geometry type: %s", __func__, lwtype_name(t2));
+ return LW_FALSE;
+ }
+ }
+ case CIRCSTRINGTYPE:
+ {
+ dl->twisted = 1;
+ switch (t2)
{
- dl->twisted = 1;
- switch ( t2 )
- {
- case POINTTYPE:
- dl->twisted=(-1);
- return lw_dist2d_point_line((LWPOINT *)lwg2, (LWLINE *)lwg1, dl);
- case LINETYPE:
- return lw_dist2d_line_line((LWLINE *)lwg1, (LWLINE *)lwg2, dl);
- case POLYGONTYPE:
- return lw_dist2d_line_poly((LWLINE *)lwg1, (LWPOLY *)lwg2, dl);
- case CIRCSTRINGTYPE:
- return lw_dist2d_line_circstring((LWLINE *)lwg1, (LWCIRCSTRING *)lwg2, dl);
- case CURVEPOLYTYPE:
- return lw_dist2d_line_curvepoly((LWLINE *)lwg1, (LWCURVEPOLY *)lwg2, dl);
- default:
- lwerror("%s: Unsupported geometry type: %s", __func__, lwtype_name(t2));
- return LW_FALSE;
- }
+ case POINTTYPE:
+ dl->twisted = -1;
+ return lw_dist2d_point_circstring((LWPOINT *)lwg2, (LWCIRCSTRING *)lwg1, dl);
+ case LINETYPE:
+ dl->twisted = -1;
+ return lw_dist2d_line_circstring((LWLINE *)lwg2, (LWCIRCSTRING *)lwg1, dl);
+ case TRIANGLETYPE:
+ dl->twisted = -1;
+ return lw_dist2d_tri_circstring((LWTRIANGLE *)lwg2, (LWCIRCSTRING *)lwg1, dl);
+ case POLYGONTYPE:
+ return lw_dist2d_circstring_poly((LWCIRCSTRING *)lwg1, (LWPOLY *)lwg2, dl);
+ case CIRCSTRINGTYPE:
+ return lw_dist2d_circstring_circstring((LWCIRCSTRING *)lwg1, (LWCIRCSTRING *)lwg2, dl);
+ case CURVEPOLYTYPE:
+ return lw_dist2d_circstring_curvepoly((LWCIRCSTRING *)lwg1, (LWCURVEPOLY *)lwg2, dl);
+ default:
+ lwerror("%s: Unsupported geometry type: %s", __func__, lwtype_name(t2));
+ return LW_FALSE;
}
+ }
+ case POLYGONTYPE:
+ {
+ dl->twisted = -1;
+ switch (t2)
+ {
+ case POINTTYPE:
+ return lw_dist2d_point_poly((LWPOINT *)lwg2, (LWPOLY *)lwg1, dl);
+ case LINETYPE:
+ return lw_dist2d_line_poly((LWLINE *)lwg2, (LWPOLY *)lwg1, dl);
+ case TRIANGLETYPE:
+ return lw_dist2d_tri_poly((LWTRIANGLE *)lwg2, (LWPOLY *)lwg1, dl);
case CIRCSTRINGTYPE:
- {
+ return lw_dist2d_circstring_poly((LWCIRCSTRING *)lwg2, (LWPOLY *)lwg1, dl);
+ case POLYGONTYPE:
dl->twisted = 1;
- switch ( t2 )
- {
- case POINTTYPE:
- dl->twisted = -1;
- return lw_dist2d_point_circstring((LWPOINT *)lwg2, (LWCIRCSTRING *)lwg1, dl);
- case LINETYPE:
- dl->twisted = -1;
- return lw_dist2d_line_circstring((LWLINE *)lwg2, (LWCIRCSTRING *)lwg1, dl);
- case POLYGONTYPE:
- return lw_dist2d_circstring_poly((LWCIRCSTRING *)lwg1, (LWPOLY *)lwg2, dl);
- case CIRCSTRINGTYPE:
- return lw_dist2d_circstring_circstring((LWCIRCSTRING *)lwg1, (LWCIRCSTRING *)lwg2, dl);
- case CURVEPOLYTYPE:
- return lw_dist2d_circstring_curvepoly((LWCIRCSTRING *)lwg1, (LWCURVEPOLY *)lwg2, dl);
- default:
- lwerror("%s: Unsupported geometry type: %s", __func__, lwtype_name(t2));
- return LW_FALSE;
- }
+ return lw_dist2d_poly_poly((LWPOLY *)lwg1, (LWPOLY *)lwg2, dl);
+ case CURVEPOLYTYPE:
+ dl->twisted = 1;
+ return lw_dist2d_poly_curvepoly((LWPOLY *)lwg1, (LWCURVEPOLY *)lwg2, dl);
+ default:
+ lwerror("%s: Unsupported geometry type: %s", __func__, lwtype_name(t2));
+ return LW_FALSE;
}
+ }
+ case CURVEPOLYTYPE:
+ {
+ dl->twisted = -1;
+ switch (t2)
+ {
+ case POINTTYPE:
+ return lw_dist2d_point_curvepoly((LWPOINT *)lwg2, (LWCURVEPOLY *)lwg1, dl);
+ case LINETYPE:
+ return lw_dist2d_line_curvepoly((LWLINE *)lwg2, (LWCURVEPOLY *)lwg1, dl);
+ case TRIANGLETYPE:
+ return lw_dist2d_tri_curvepoly((LWTRIANGLE *)lwg2, (LWCURVEPOLY *)lwg1, dl);
case POLYGONTYPE:
- {
- dl->twisted = -1;
- switch ( t2 )
- {
- case POINTTYPE:
- return lw_dist2d_point_poly((LWPOINT *)lwg2, (LWPOLY *)lwg1, dl);
- case LINETYPE:
- return lw_dist2d_line_poly((LWLINE *)lwg2, (LWPOLY *)lwg1, dl);
- case CIRCSTRINGTYPE:
- return lw_dist2d_circstring_poly((LWCIRCSTRING *)lwg2, (LWPOLY *)lwg1, dl);
- case POLYGONTYPE:
- dl->twisted = 1;
- return lw_dist2d_poly_poly((LWPOLY *)lwg1, (LWPOLY *)lwg2, dl);
- case CURVEPOLYTYPE:
- dl->twisted = 1;
- return lw_dist2d_poly_curvepoly((LWPOLY *)lwg1, (LWCURVEPOLY *)lwg2, dl);
- default:
- lwerror("%s: Unsupported geometry type: %s", __func__, lwtype_name(t2));
- return LW_FALSE;
- }
- }
+ return lw_dist2d_poly_curvepoly((LWPOLY *)lwg2, (LWCURVEPOLY *)lwg1, dl);
+ case CIRCSTRINGTYPE:
+ return lw_dist2d_circstring_curvepoly((LWCIRCSTRING *)lwg2, (LWCURVEPOLY *)lwg1, dl);
case CURVEPOLYTYPE:
- {
- dl->twisted = (-1);
- switch ( t2 )
- {
- case POINTTYPE:
- return lw_dist2d_point_curvepoly((LWPOINT *)lwg2, (LWCURVEPOLY *)lwg1, dl);
- case LINETYPE:
- return lw_dist2d_line_curvepoly((LWLINE *)lwg2, (LWCURVEPOLY *)lwg1, dl);
- case POLYGONTYPE:
- return lw_dist2d_poly_curvepoly((LWPOLY *)lwg2, (LWCURVEPOLY *)lwg1, dl);
- case CIRCSTRINGTYPE:
- return lw_dist2d_circstring_curvepoly((LWCIRCSTRING *)lwg2, (LWCURVEPOLY *)lwg1, dl);
- case CURVEPOLYTYPE:
- dl->twisted = 1;
- return lw_dist2d_curvepoly_curvepoly((LWCURVEPOLY *)lwg1, (LWCURVEPOLY *)lwg2, dl);
- default:
- lwerror("%s: Unsupported geometry type: %s", __func__, lwtype_name(t2));
- return LW_FALSE;
- }
- }
+ dl->twisted = 1;
+ return lw_dist2d_curvepoly_curvepoly((LWCURVEPOLY *)lwg1, (LWCURVEPOLY *)lwg2, dl);
default:
- {
- lwerror("%s: Unsupported geometry type: %s", __func__, lwtype_name(t1));
+ lwerror("%s: Unsupported geometry type: %s", __func__, lwtype_name(t2));
return LW_FALSE;
}
}
+ default:
+ {
+ lwerror("%s: Unsupported geometry type: %s", __func__, lwtype_name(t1));
+ return LW_FALSE;
+ }
+ }
return LW_FALSE;
}
-
-
-
-/**
-
-We have to check for overlapping bboxes
-*/
+/* Check for overlapping bboxes */
int
-lw_dist2d_check_overlap(LWGEOM *lwg1,LWGEOM *lwg2)
+lw_dist2d_check_overlap(LWGEOM *lwg1, LWGEOM *lwg2)
{
LWDEBUG(2, "lw_dist2d_check_overlap is called");
- if ( ! lwg1->bbox )
+ if (!lwg1->bbox)
lwgeom_calculate_gbox(lwg1, lwg1->bbox);
- if ( ! lwg2->bbox )
+ if (!lwg2->bbox)
lwgeom_calculate_gbox(lwg2, lwg2->bbox);
- /*Check if the geometries intersect.
- */
- if ((lwg1->bbox->xmax<lwg2->bbox->xmin||lwg1->bbox->xmin>lwg2->bbox->xmax||lwg1->bbox->ymax<lwg2->bbox->ymin||lwg1->bbox->ymin>lwg2->bbox->ymax))
+ /* Check if the geometries intersect. */
+ if ((lwg1->bbox->xmax < lwg2->bbox->xmin || lwg1->bbox->xmin > lwg2->bbox->xmax ||
+ lwg1->bbox->ymax < lwg2->bbox->ymin || lwg1->bbox->ymin > lwg2->bbox->ymax))
{
LWDEBUG(3, "geometries bboxes did not overlap");
return LW_FALSE;
@@ -517,16 +534,13 @@
return LW_TRUE;
}
-/**
-
-Here the geometries are distributed for the new faster distance-calculations
-*/
+/** Geometries are distributed for the new faster distance-calculations */
int
lw_dist2d_distribute_fast(LWGEOM *lwg1, LWGEOM *lwg2, DISTPTS *dl)
{
POINTARRAY *pa1, *pa2;
- int type1 = lwg1->type;
- int type2 = lwg2->type;
+ int type1 = lwg1->type;
+ int type2 = lwg2->type;
LWDEBUGF(2, "lw_dist2d_distribute_fast is called with typ1=%d, type2=%d", lwg1->type, lwg2->type);
@@ -538,6 +552,9 @@
case POLYGONTYPE:
pa1 = ((LWPOLY *)lwg1)->rings[0];
break;
+ case TRIANGLETYPE:
+ pa1 = ((LWTRIANGLE *)lwg1)->points;
+ break;
default:
lwerror("Unsupported geometry1 type: %s", lwtype_name(type1));
return LW_FALSE;
@@ -550,11 +567,14 @@
case POLYGONTYPE:
pa2 = ((LWPOLY *)lwg2)->rings[0];
break;
+ case TRIANGLETYPE:
+ pa2 = ((LWTRIANGLE *)lwg2)->points;
+ break;
default:
lwerror("Unsupported geometry2 type: %s", lwtype_name(type1));
return LW_FALSE;
}
- dl->twisted=1;
+ dl->twisted = 1;
return lw_dist2d_fast_ptarray_ptarray(pa1, pa2, dl, lwg1->bbox, lwg2->bbox);
}
@@ -562,7 +582,6 @@
End of Preprocessing functions
--------------------------------------------------------------------------------------------------------------*/
-
/*------------------------------------------------------------------------------------------------------------
Brute force functions
The old way of calculating distances, now used for:
@@ -571,33 +590,42 @@
--------------------------------------------------------------------------------------------------------------*/
/**
-
point to point calculation
*/
int
lw_dist2d_point_point(LWPOINT *point1, LWPOINT *point2, DISTPTS *dl)
{
- const POINT2D *p1, *p2;
-
- p1 = getPoint2d_cp(point1->point, 0);
- p2 = getPoint2d_cp(point2->point, 0);
-
+ const POINT2D *p1 = getPoint2d_cp(point1->point, 0);
+ const POINT2D *p2 = getPoint2d_cp(point2->point, 0);
return lw_dist2d_pt_pt(p1, p2, dl);
}
/**
-
point to line calculation
*/
int
lw_dist2d_point_line(LWPOINT *point, LWLINE *line, DISTPTS *dl)
{
- const POINT2D *p;
- LWDEBUG(2, "lw_dist2d_point_line is called");
- p = getPoint2d_cp(point->point, 0);
+ const POINT2D *p = getPoint2d_cp(point->point, 0);
return lw_dist2d_pt_ptarray(p, line->points, dl);
}
int
+lw_dist2d_point_tri(LWPOINT *point, LWTRIANGLE *tri, DISTPTS *dl)
+{
+ const POINT2D *pt = getPoint2d_cp(point->point, 0);
+ /* Is point inside triangle? */
+ if (dl->mode == DIST_MIN && ptarray_contains_point(tri->points, pt) != LW_OUTSIDE)
+ {
+ dl->distance = 0.0;
+ dl->p1.x = dl->p2.x = pt->x;
+ dl->p1.y = dl->p2.y = pt->y;
+ return LW_TRUE;
+ }
+
+ return lw_dist2d_pt_ptarray(pt, tri->points, dl);
+}
+
+int
lw_dist2d_point_circstring(LWPOINT *point, LWCIRCSTRING *circ, DISTPTS *dl)
{
const POINT2D *p;
@@ -613,24 +641,15 @@
int
lw_dist2d_point_poly(LWPOINT *point, LWPOLY *poly, DISTPTS *dl)
{
- const POINT2D *p;
- uint32_t i;
+ const POINT2D *p = getPoint2d_cp(point->point, 0);
- LWDEBUG(2, "lw_dist2d_point_poly called");
-
- p = getPoint2d_cp(point->point, 0);
-
+ /* Max distance? Check only outer ring.*/
if (dl->mode == DIST_MAX)
- {
- LWDEBUG(3, "looking for maxdistance");
return lw_dist2d_pt_ptarray(p, poly->rings[0], dl);
- }
+
/* Return distance to outer ring if not inside it */
- if ( ptarray_contains_point(poly->rings[0], p) == LW_OUTSIDE )
- {
- LWDEBUG(3, "first point not inside outer-ring");
+ if (ptarray_contains_point(poly->rings[0], p) == LW_OUTSIDE)
return lw_dist2d_pt_ptarray(p, poly->rings[0], dl);
- }
/*
* Inside the outer ring.
@@ -638,84 +657,69 @@
* see if its inside. If not, distance==0.
* Otherwise, distance = pt to ring distance
*/
- for ( i = 1; i < poly->nrings; i++)
- {
- /* Inside a hole. Distance = pt -> ring */
- if ( ptarray_contains_point(poly->rings[i], p) != LW_OUTSIDE )
- {
- LWDEBUG(3, " inside an hole");
+ for (uint32_t i = 1; i < poly->nrings; i++)
+ if (ptarray_contains_point(poly->rings[i], p) != LW_OUTSIDE)
return lw_dist2d_pt_ptarray(p, poly->rings[i], dl);
- }
- }
- LWDEBUG(3, " inside the polygon");
- if (dl->mode == DIST_MIN)
- {
- dl->distance = 0.0;
- dl->p1.x = dl->p2.x = p->x;
- dl->p1.y = dl->p2.y = p->y;
- }
- return LW_TRUE; /* Is inside the polygon */
+ /* Is inside the polygon */
+ dl->distance = 0.0;
+ dl->p1.x = dl->p2.x = p->x;
+ dl->p1.y = dl->p2.y = p->y;
+ return LW_TRUE;
}
int
lw_dist2d_point_curvepoly(LWPOINT *point, LWCURVEPOLY *poly, DISTPTS *dl)
{
- const POINT2D *p;
- uint32_t i;
+ const POINT2D *p = getPoint2d_cp(point->point, 0);
- p = getPoint2d_cp(point->point, 0);
-
if (dl->mode == DIST_MAX)
lwerror("lw_dist2d_point_curvepoly cannot calculate max distance");
/* Return distance to outer ring if not inside it */
- if ( lwgeom_contains_point(poly->rings[0], p) == LW_OUTSIDE )
- {
- return lw_dist2d_recursive((LWGEOM*)point, poly->rings[0], dl);
- }
+ if (lwgeom_contains_point(poly->rings[0], p) == LW_OUTSIDE)
+ return lw_dist2d_recursive((LWGEOM *)point, poly->rings[0], dl);
- /*
- * Inside the outer ring.
- * Scan though each of the inner rings looking to
- * see if its inside. If not, distance==0.
- * Otherwise, distance = pt to ring distance
+ /* Inside the outer ring.
+ * Scan though each of the inner rings looking to see if its inside. If not, distance==0.
+ * Otherwise, distance = pt to ring distance.
*/
- for ( i = 1; i < poly->nrings; i++)
- {
- /* Inside a hole. Distance = pt -> ring */
- if ( lwgeom_contains_point(poly->rings[i], p) != LW_OUTSIDE )
- {
- LWDEBUG(3, " inside a hole");
- return lw_dist2d_recursive((LWGEOM*)point, poly->rings[i], dl);
- }
- }
+ for (uint32_t i = 1; i < poly->nrings; i++)
+ if (lwgeom_contains_point(poly->rings[i], p) == LW_INSIDE)
+ return lw_dist2d_recursive((LWGEOM *)point, poly->rings[i], dl);
- LWDEBUG(3, " inside the polygon");
- if (dl->mode == DIST_MIN)
- {
- dl->distance = 0.0;
- dl->p1.x = dl->p2.x = p->x;
- dl->p1.y = dl->p2.y = p->y;
- }
-
- return LW_TRUE; /* Is inside the polygon */
+ /* Is inside the polygon */
+ dl->distance = 0.0;
+ dl->p1.x = dl->p2.x = p->x;
+ dl->p1.y = dl->p2.y = p->y;
+ return LW_TRUE;
}
-/**
-
-line to line calculation
-*/
int
lw_dist2d_line_line(LWLINE *line1, LWLINE *line2, DISTPTS *dl)
{
POINTARRAY *pa1 = line1->points;
POINTARRAY *pa2 = line2->points;
- LWDEBUG(2, "lw_dist2d_line_line is called");
return lw_dist2d_ptarray_ptarray(pa1, pa2, dl);
}
int
+lw_dist2d_line_tri(LWLINE *line, LWTRIANGLE *tri, DISTPTS *dl)
+{
+ const POINT2D *pt = getPoint2d_cp(line->points, 0);
+ /* Is there a point inside triangle? */
+ if (dl->mode == DIST_MIN && ptarray_contains_point(tri->points, pt) != LW_OUTSIDE)
+ {
+ dl->distance = 0.0;
+ dl->p1.x = dl->p2.x = pt->x;
+ dl->p1.y = dl->p2.y = pt->y;
+ return LW_TRUE;
+ }
+
+ return lw_dist2d_ptarray_ptarray(line->points, tri->points, dl);
+}
+
+int
lw_dist2d_line_circstring(LWLINE *line1, LWCIRCSTRING *line2, DISTPTS *dl)
{
return lw_dist2d_ptarray_ptarrayarc(line1->points, line2->points, dl);
@@ -735,57 +739,63 @@
int
lw_dist2d_line_poly(LWLINE *line, LWPOLY *poly, DISTPTS *dl)
{
- const POINT2D *pt;
- uint32_t i;
+ POINTARRAY *pa = line->points;
+ const POINT2D *pt = getPoint2d_cp(pa, 0);
- LWDEBUGF(2, "lw_dist2d_line_poly called (%d rings)", poly->nrings);
+ /* Line has a pount outside poly. Check distance to outer ring only. */
+ if (ptarray_contains_point(poly->rings[0], pt) == LW_OUTSIDE)
+ return lw_dist2d_ptarray_ptarray(pa, poly->rings[0], dl);
- pt = getPoint2d_cp(line->points, 0);
- if ( ptarray_contains_point(poly->rings[0], pt) == LW_OUTSIDE )
+ for (uint32_t i = 1; i < poly->nrings; i++)
{
- return lw_dist2d_ptarray_ptarray(line->points, poly->rings[0], dl);
+ if (!lw_dist2d_ptarray_ptarray(pa, poly->rings[i], dl))
+ return LW_FALSE;
+
+ /* just a check if the answer is already given */
+ if (dl->distance <= dl->tolerance && dl->mode == DIST_MIN)
+ return LW_TRUE;
}
- for (i=1; i<poly->nrings; i++)
+ /* It's inside a hole, then the actual distance is the min ring distance */
+ for (uint32_t i = 1; i < poly->nrings; i++)
+ if (ptarray_contains_point(poly->rings[i], pt) != LW_OUTSIDE)
+ return LW_TRUE;
+
+ /* Not in hole, so inside polygon */
+ if (dl->mode == DIST_MIN)
{
- if (!lw_dist2d_ptarray_ptarray(line->points, poly->rings[i], dl)) return LW_FALSE;
-
- LWDEBUGF(3, " distance from ring %d: %f, mindist: %f",
- i, dl->distance, dl->tolerance);
- /* just a check if the answer is already given */
- if (dl->distance<=dl->tolerance && dl->mode == DIST_MIN) return LW_TRUE;
+ dl->distance = 0.0;
+ dl->p1.x = dl->p2.x = pt->x;
+ dl->p1.y = dl->p2.y = pt->y;
}
+ return LW_TRUE;
+}
- /*
- * No intersection, have to check if a point is
- * inside polygon
- */
- pt = getPoint2d_cp(line->points, 0);
+int
+lw_dist2d_line_curvepoly(LWLINE *line, LWCURVEPOLY *poly, DISTPTS *dl)
+{
+ const POINT2D *pt = getPoint2d_cp(line->points, 0);
- /*
- * Outside outer ring, so min distance to a ring
- * is the actual min distance
+ /* Line has a pount outside curvepoly. Check distance to outer ring only. */
+ if (lwgeom_contains_point(poly->rings[0], pt) == LW_OUTSIDE)
+ return lw_dist2d_recursive((LWGEOM *)line, poly->rings[0], dl);
- if ( ! pt_in_ring_2d(&pt, poly->rings[0]) )
+ for (uint32_t i = 1; i < poly->nrings; i++)
{
- return ;
- } */
+ if (!lw_dist2d_recursive((LWGEOM *)line, poly->rings[i], dl))
+ return LW_FALSE;
- /*
- * Its in the outer ring.
- * Have to check if its inside a hole
- */
- for (i=1; i<poly->nrings; i++)
- {
- if ( ptarray_contains_point(poly->rings[i], pt) != LW_OUTSIDE )
- {
- /*
- * Its inside a hole, then the actual
- * distance is the min ring distance
- */
+ /* just a check if the answer is already given */
+ if (dl->distance <= dl->tolerance && dl->mode == DIST_MIN)
return LW_TRUE;
- }
}
+
+ /* It's inside a hole, then distance is actual min distance */
+ for (uint32_t i = 1; i < poly->nrings; i++)
+ if (lwgeom_contains_point(poly->rings[i], pt) != LW_OUTSIDE)
+ return LW_TRUE;
+
+ /* Not in hole, so inside polygon */
if (dl->mode == DIST_MIN)
{
dl->distance = 0.0;
@@ -792,46 +802,167 @@
dl->p1.x = dl->p2.x = pt->x;
dl->p1.y = dl->p2.y = pt->y;
}
- return LW_TRUE; /* Not in hole, so inside polygon */
+ return LW_TRUE;
}
int
-lw_dist2d_line_curvepoly(LWLINE *line, LWCURVEPOLY *poly, DISTPTS *dl)
+lw_dist2d_tri_tri(LWTRIANGLE *tri1, LWTRIANGLE *tri2, DISTPTS *dl)
{
- const POINT2D *pt = getPoint2d_cp(line->points, 0);
- uint32_t i;
+ POINTARRAY *pa1 = tri1->points;
+ POINTARRAY *pa2 = tri2->points;
+ const POINT2D *pt = getPoint2d_cp(pa2, 0);
+ if (dl->mode == DIST_MIN && ptarray_contains_point(pa1, pt) != LW_OUTSIDE)
+ {
+ dl->distance = 0.0;
+ dl->p1.x = dl->p2.x = pt->x;
+ dl->p1.y = dl->p2.y = pt->y;
+ return LW_TRUE;
+ }
- if ( lwgeom_contains_point(poly->rings[0], pt) == LW_OUTSIDE )
+ pt = getPoint2d_cp(pa1, 0);
+ if (dl->mode == DIST_MIN && ptarray_contains_point(pa2, pt) != LW_OUTSIDE)
{
- return lw_dist2d_recursive((LWGEOM*)line, poly->rings[0], dl);
+ dl->distance = 0.0;
+ dl->p1.x = dl->p2.x = pt->x;
+ dl->p1.y = dl->p2.y = pt->y;
+ return LW_TRUE;
}
- for ( i = 1; i < poly->nrings; i++ )
+ return lw_dist2d_ptarray_ptarray(pa1, pa2, dl);
+}
+
+int
+lw_dist2d_tri_poly(LWTRIANGLE *tri, LWPOLY *poly, DISTPTS *dl)
+{
+ POINTARRAY *pa = tri->points;
+ const POINT2D *pt = getPoint2d_cp(pa, 0);
+
+ /* If we are looking for maxdistance, just check the outer rings.*/
+ if (dl->mode == DIST_MAX)
+ return lw_dist2d_ptarray_ptarray(pa, poly->rings[0], dl);
+
+ /* Triangle has a point outside poly. Check distance to outer ring only. */
+ if (ptarray_contains_point(poly->rings[0], pt) == LW_OUTSIDE)
{
- if ( ! lw_dist2d_recursive((LWGEOM*)line, poly->rings[i], dl) )
+ if (!lw_dist2d_ptarray_ptarray(pa, poly->rings[0], dl))
return LW_FALSE;
- if ( dl->distance<=dl->tolerance && dl->mode == DIST_MIN )
+ /* just a check if the answer is already given */
+ if (dl->distance <= dl->tolerance)
return LW_TRUE;
+
+ /* Maybe poly is inside triangle? */
+ const POINT2D *pt2 = getPoint2d_cp(poly->rings[0], 0);
+ if (ptarray_contains_point(pa, pt2) != LW_OUTSIDE)
+ {
+ dl->distance = 0.0;
+ dl->p1.x = dl->p2.x = pt2->x;
+ dl->p1.y = dl->p2.y = pt2->y;
+ return LW_TRUE;
+ }
}
- for ( i=1; i < poly->nrings; i++ )
+ for (uint32_t i = 1; i < poly->nrings; i++)
{
- if ( lwgeom_contains_point(poly->rings[i],pt) != LW_OUTSIDE )
+ if (!lw_dist2d_ptarray_ptarray(pa, poly->rings[i], dl))
+ return LW_FALSE;
+
+ /* just a check if the answer is already given */
+ if (dl->distance <= dl->tolerance && dl->mode == DIST_MIN)
+ return LW_TRUE;
+ }
+
+ /* It's inside a hole, then the actual distance is the min ring distance */
+ for (uint32_t i = 1; i < poly->nrings; i++)
+ if (ptarray_contains_point(poly->rings[i], pt) != LW_OUTSIDE)
+ return LW_TRUE;
+
+ /* Not in hole, so inside polygon */
+ dl->distance = 0.0;
+ dl->p1.x = dl->p2.x = pt->x;
+ dl->p1.y = dl->p2.y = pt->y;
+ return LW_TRUE;
+}
+static const POINT2D *
+lw_curvering_getfirstpoint2d_cp(LWGEOM *geom)
+{
+ switch (geom->type)
+ {
+ case LINETYPE:
+ return getPoint2d_cp(((LWLINE *)geom)->points, 0);
+ case CIRCSTRINGTYPE:
+ return getPoint2d_cp(((LWCIRCSTRING *)geom)->points, 0);
+ case COMPOUNDTYPE:
+ {
+ LWCOMPOUND *comp = (LWCOMPOUND *)geom;
+ LWLINE *line = (LWLINE *)(comp->geoms[0]);
+ return getPoint2d_cp(line->points, 0);
+ }
+ default:
+ lwerror("lw_curvering_getfirstpoint2d_cp: unknown type");
+ }
+ return NULL;
+}
+
+int
+lw_dist2d_tri_curvepoly(LWTRIANGLE *tri, LWCURVEPOLY *poly, DISTPTS *dl)
+{
+ const POINT2D *pt = getPoint2d_cp(tri->points, 0);
+
+ /* If we are looking for maxdistance, just check the outer rings.*/
+ if (dl->mode == DIST_MAX)
+ return lw_dist2d_recursive((LWGEOM *)tri, poly->rings[0], dl);
+
+ /* Line has a pount outside curvepoly. Check distance to outer ring only. */
+ if (lwgeom_contains_point(poly->rings[0], pt) == LW_OUTSIDE)
+ {
+ if (lw_dist2d_recursive((LWGEOM *)tri, poly->rings[0], dl))
+ return LW_TRUE;
+ /* Maybe poly is inside triangle? */
+ if (lwgeom_contains_point((LWGEOM *)tri, lw_curvering_getfirstpoint2d_cp(poly->rings[0])) != LW_OUTSIDE)
{
- /* Its inside a hole, then the actual */
+ dl->distance = 0.0;
+ dl->p1.x = dl->p2.x = pt->x;
+ dl->p1.y = dl->p2.y = pt->y;
return LW_TRUE;
}
}
- if (dl->mode == DIST_MIN)
+ for (uint32_t i = 1; i < poly->nrings; i++)
{
+ if (!lw_dist2d_recursive((LWGEOM *)tri, poly->rings[i], dl))
+ return LW_FALSE;
+
+ /* just a check if the answer is already given */
+ if (dl->distance <= dl->tolerance && dl->mode == DIST_MIN)
+ return LW_TRUE;
+ }
+
+ /* It's inside a hole, then distance is actual min distance */
+ for (uint32_t i = 1; i < poly->nrings; i++)
+ if (lwgeom_contains_point(poly->rings[i], pt) != LW_OUTSIDE)
+ return LW_TRUE;
+
+ /* Not in hole, so inside polygon */
+ dl->distance = 0.0;
+ dl->p1.x = dl->p2.x = pt->x;
+ dl->p1.y = dl->p2.y = pt->y;
+ return LW_TRUE;
+}
+
+int
+lw_dist2d_tri_circstring(LWTRIANGLE *tri, LWCIRCSTRING *line, DISTPTS *dl)
+{
+ const POINT2D *pt = lw_curvering_getfirstpoint2d_cp((LWGEOM *)line);
+ if (ptarray_contains_point(tri->points, pt) != LW_OUTSIDE && dl->mode == DIST_MIN)
+ {
dl->distance = 0.0;
dl->p1.x = dl->p2.x = pt->x;
dl->p1.y = dl->p2.y = pt->y;
+ return LW_TRUE;
}
- return LW_TRUE; /* Not in hole, so inside polygon */
+ return lw_dist2d_ptarray_ptarrayarc(tri->points, line->points, dl);
}
/**
@@ -840,62 +971,49 @@
2 check if poly1 has first point outside poly2 and vice versa, if so, just check outer rings
3 check if first point of poly2 is in a hole of poly1. If so check outer ring of poly2 against that hole of poly1
4 check if first point of poly1 is in a hole of poly2. If so check outer ring of poly1 against that hole of poly2
-5 If we have come all the way here we know that the first point of one of them is inside the other ones outer ring and not in holes so we check which one is inside.
+5 If we have come all the way here we know that the first point of one of them is inside the other ones outer ring
+and not in holes so we check which one is inside.
*/
int
lw_dist2d_poly_poly(LWPOLY *poly1, LWPOLY *poly2, DISTPTS *dl)
{
-
const POINT2D *pt;
- uint32_t i;
LWDEBUG(2, "lw_dist2d_poly_poly called");
/*1 if we are looking for maxdistance, just check the outer rings.*/
if (dl->mode == DIST_MAX)
- {
return lw_dist2d_ptarray_ptarray(poly1->rings[0], poly2->rings[0], dl);
- }
-
/* 2 check if poly1 has first point outside poly2 and vice versa, if so, just check outer rings
- here it would be possible to handle the information about which one is inside which one and only search for the smaller ones in the bigger ones holes.*/
+ here it would be possible to handle the information about which one is inside which one and only search for the
+ smaller ones in the bigger ones holes.*/
pt = getPoint2d_cp(poly1->rings[0], 0);
- if ( ptarray_contains_point(poly2->rings[0], pt) == LW_OUTSIDE )
+ if (ptarray_contains_point(poly2->rings[0], pt) == LW_OUTSIDE)
{
pt = getPoint2d_cp(poly2->rings[0], 0);
- if ( ptarray_contains_point(poly1->rings[0], pt) == LW_OUTSIDE )
- {
+ if (ptarray_contains_point(poly1->rings[0], pt) == LW_OUTSIDE)
return lw_dist2d_ptarray_ptarray(poly1->rings[0], poly2->rings[0], dl);
- }
}
- /*3 check if first point of poly2 is in a hole of poly1. If so check outer ring of poly2 against that hole of poly1*/
+ /*3 check if first point of poly2 is in a hole of poly1. If so check outer ring of poly2 against that hole
+ * of poly1*/
pt = getPoint2d_cp(poly2->rings[0], 0);
- for (i=1; i<poly1->nrings; i++)
- {
- /* Inside a hole */
- if ( ptarray_contains_point(poly1->rings[i], pt) != LW_OUTSIDE )
- {
+ for (uint32_t i = 1; i < poly1->nrings; i++)
+ if (ptarray_contains_point(poly1->rings[i], pt) != LW_OUTSIDE)
return lw_dist2d_ptarray_ptarray(poly1->rings[i], poly2->rings[0], dl);
- }
- }
- /*4 check if first point of poly1 is in a hole of poly2. If so check outer ring of poly1 against that hole of poly2*/
+ /*4 check if first point of poly1 is in a hole of poly2. If so check outer ring of poly1 against that hole
+ * of poly2*/
pt = getPoint2d_cp(poly1->rings[0], 0);
- for (i=1; i<poly2->nrings; i++)
- {
- /* Inside a hole */
- if ( ptarray_contains_point(poly2->rings[i], pt) != LW_OUTSIDE )
- {
+ for (uint32_t i = 1; i < poly2->nrings; i++)
+ if (ptarray_contains_point(poly2->rings[i], pt) != LW_OUTSIDE)
return lw_dist2d_ptarray_ptarray(poly1->rings[0], poly2->rings[i], dl);
- }
- }
-
- /*5 If we have come all the way here we know that the first point of one of them is inside the other ones outer ring and not in holes so we check wich one is inside.*/
+ /*5 If we have come all the way here we know that the first point of one of them is inside the other ones
+ * outer ring and not in holes so we check wich one is inside.*/
pt = getPoint2d_cp(poly1->rings[0], 0);
- if ( ptarray_contains_point(poly2->rings[0], pt) != LW_OUTSIDE )
+ if (ptarray_contains_point(poly2->rings[0], pt) != LW_OUTSIDE)
{
dl->distance = 0.0;
dl->p1.x = dl->p2.x = pt->x;
@@ -904,7 +1022,7 @@
}
pt = getPoint2d_cp(poly2->rings[0], 0);
- if ( ptarray_contains_point(poly1->rings[0], pt) != LW_OUTSIDE )
+ if (ptarray_contains_point(poly1->rings[0], pt) != LW_OUTSIDE)
{
dl->distance = 0.0;
dl->p1.x = dl->p2.x = pt->x;
@@ -912,7 +1030,6 @@
return LW_TRUE;
}
-
lwerror("Unspecified error in function lw_dist2d_poly_poly");
return LW_FALSE;
}
@@ -922,7 +1039,7 @@
{
LWCURVEPOLY *curvepoly1 = lwcurvepoly_construct_from_lwpoly(poly1);
int rv = lw_dist2d_curvepoly_curvepoly(curvepoly1, curvepoly2, dl);
- lwgeom_free((LWGEOM*)curvepoly1);
+ lwgeom_free((LWGEOM *)curvepoly1);
return rv;
}
@@ -930,16 +1047,15 @@
lw_dist2d_circstring_poly(LWCIRCSTRING *circ, LWPOLY *poly, DISTPTS *dl)
{
LWCURVEPOLY *curvepoly = lwcurvepoly_construct_from_lwpoly(poly);
- int rv = lw_dist2d_line_curvepoly((LWLINE*)circ, curvepoly, dl);
- lwgeom_free((LWGEOM*)curvepoly);
+ int rv = lw_dist2d_line_curvepoly((LWLINE *)circ, curvepoly, dl);
+ lwgeom_free((LWGEOM *)curvepoly);
return rv;
}
-
int
lw_dist2d_circstring_curvepoly(LWCIRCSTRING *circ, LWCURVEPOLY *poly, DISTPTS *dl)
{
- return lw_dist2d_line_curvepoly((LWLINE*)circ, poly, dl);
+ return lw_dist2d_line_curvepoly((LWLINE *)circ, poly, dl);
}
int
@@ -948,80 +1064,44 @@
return lw_dist2d_ptarrayarc_ptarrayarc(line1->points, line2->points, dl);
}
-static const POINT2D *
-lw_curvering_getfirstpoint2d_cp(LWGEOM *geom)
-{
- switch( geom->type )
- {
- case LINETYPE:
- return getPoint2d_cp(((LWLINE*)geom)->points, 0);
- case CIRCSTRINGTYPE:
- return getPoint2d_cp(((LWCIRCSTRING*)geom)->points, 0);
- case COMPOUNDTYPE:
- {
- LWCOMPOUND *comp = (LWCOMPOUND*)geom;
- LWLINE *line = (LWLINE*)(comp->geoms[0]);
- return getPoint2d_cp(line->points, 0);
- }
- default:
- lwerror("lw_curvering_getfirstpoint2d_cp: unknown type");
- }
- return NULL;
-}
-
int
lw_dist2d_curvepoly_curvepoly(LWCURVEPOLY *poly1, LWCURVEPOLY *poly2, DISTPTS *dl)
{
const POINT2D *pt;
- uint32_t i;
- LWDEBUG(2, "lw_dist2d_curvepoly_curvepoly called");
-
/*1 if we are looking for maxdistance, just check the outer rings.*/
if (dl->mode == DIST_MAX)
- {
- return lw_dist2d_recursive(poly1->rings[0], poly2->rings[0], dl);
- }
+ return lw_dist2d_recursive(poly1->rings[0], poly2->rings[0], dl);
-
/* 2 check if poly1 has first point outside poly2 and vice versa, if so, just check outer rings
- here it would be possible to handle the information about which one is inside which one and only search for the smaller ones in the bigger ones holes.*/
+ here it would be possible to handle the information about which one is inside which one and only search for the
+ smaller ones in the bigger ones holes.*/
pt = lw_curvering_getfirstpoint2d_cp(poly1->rings[0]);
- if ( lwgeom_contains_point(poly2->rings[0], pt) == LW_OUTSIDE )
+ if (lwgeom_contains_point(poly2->rings[0], pt) == LW_OUTSIDE)
{
pt = lw_curvering_getfirstpoint2d_cp(poly2->rings[0]);
- if ( lwgeom_contains_point(poly1->rings[0], pt) == LW_OUTSIDE )
- {
+ if (lwgeom_contains_point(poly1->rings[0], pt) == LW_OUTSIDE)
return lw_dist2d_recursive(poly1->rings[0], poly2->rings[0], dl);
- }
}
- /*3 check if first point of poly2 is in a hole of poly1. If so check outer ring of poly2 against that hole of poly1*/
+ /*3 check if first point of poly2 is in a hole of poly1. If so check outer ring of poly2 against that hole
+ * of poly1*/
pt = lw_curvering_getfirstpoint2d_cp(poly2->rings[0]);
- for (i = 1; i < poly1->nrings; i++)
- {
- /* Inside a hole */
- if ( lwgeom_contains_point(poly1->rings[i], pt) != LW_OUTSIDE )
- {
+ for (uint32_t i = 1; i < poly1->nrings; i++)
+ if (lwgeom_contains_point(poly1->rings[i], pt) != LW_OUTSIDE)
return lw_dist2d_recursive(poly1->rings[i], poly2->rings[0], dl);
- }
- }
- /*4 check if first point of poly1 is in a hole of poly2. If so check outer ring of poly1 against that hole of poly2*/
+ /*4 check if first point of poly1 is in a hole of poly2. If so check outer ring of poly1 against that hole
+ * of poly2*/
pt = lw_curvering_getfirstpoint2d_cp(poly1->rings[0]);
- for (i = 1; i < poly2->nrings; i++)
- {
- /* Inside a hole */
- if ( lwgeom_contains_point(poly2->rings[i], pt) != LW_OUTSIDE )
- {
- return lw_dist2d_recursive(poly1->rings[0], poly2->rings[i], dl);
- }
- }
+ for (uint32_t i = 1; i < poly2->nrings; i++)
+ if (lwgeom_contains_point(poly2->rings[i], pt) != LW_OUTSIDE)
+ return lw_dist2d_recursive(poly1->rings[0], poly2->rings[i], dl);
-
- /*5 If we have come all the way here we know that the first point of one of them is inside the other ones outer ring and not in holes so we check which one is inside.*/
+ /*5 If we have come all the way here we know that the first point of one of them is inside the other ones
+ * outer ring and not in holes so we check which one is inside.*/
pt = lw_curvering_getfirstpoint2d_cp(poly1->rings[0]);
- if ( lwgeom_contains_point(poly2->rings[0], pt) != LW_OUTSIDE )
+ if (lwgeom_contains_point(poly2->rings[0], pt) != LW_OUTSIDE)
{
dl->distance = 0.0;
dl->p1.x = dl->p2.x = pt->x;
@@ -1030,7 +1110,7 @@
}
pt = lw_curvering_getfirstpoint2d_cp(poly2->rings[0]);
- if ( lwgeom_contains_point(poly1->rings[0], pt) != LW_OUTSIDE )
+ if (lwgeom_contains_point(poly1->rings[0], pt) != LW_OUTSIDE)
{
dl->distance = 0.0;
dl->p1.x = dl->p2.x = pt->x;
@@ -1042,32 +1122,30 @@
return LW_FALSE;
}
-
-
/**
* search all the segments of pointarray to see which one is closest to p1
* Returns minimum distance between point and pointarray
*/
int
-lw_dist2d_pt_ptarray(const POINT2D *p, POINTARRAY *pa,DISTPTS *dl)
+lw_dist2d_pt_ptarray(const POINT2D *p, POINTARRAY *pa, DISTPTS *dl)
{
- uint32_t t;
const POINT2D *start, *end;
int twist = dl->twisted;
- LWDEBUG(2, "lw_dist2d_pt_ptarray is called");
-
start = getPoint2d_cp(pa, 0);
- if ( !lw_dist2d_pt_pt(p, start, dl) ) return LW_FALSE;
+ if (!lw_dist2d_pt_pt(p, start, dl))
+ return LW_FALSE;
- for (t=1; t<pa->npoints; t++)
+ for (uint32_t t = 1; t < pa->npoints; t++)
{
- dl->twisted=twist;
+ dl->twisted = twist;
end = getPoint2d_cp(pa, t);
- if (!lw_dist2d_pt_seg(p, start, end, dl)) return LW_FALSE;
+ if (!lw_dist2d_pt_seg(p, start, end, dl))
+ return LW_FALSE;
- if (dl->distance<=dl->tolerance && dl->mode == DIST_MIN) return LW_TRUE; /*just a check if the answer is already given*/
+ if (dl->distance <= dl->tolerance && dl->mode == DIST_MIN)
+ return LW_TRUE; /*just a check if the answer is already given*/
start = end;
}
@@ -1075,9 +1153,9 @@
}
/**
-* Search all the arcs of pointarray to see which one is closest to p1
-* Returns minimum distance between point and arc pointarray.
-*/
+ * Search all the arcs of pointarray to see which one is closest to p1
+ * Returns minimum distance between point and arc pointarray.
+ */
int
lw_dist2d_pt_ptarrayarc(const POINT2D *p, const POINTARRAY *pa, DISTPTS *dl)
{
@@ -1089,7 +1167,7 @@
LWDEBUG(2, "lw_dist2d_pt_ptarrayarc is called");
- if ( pa->npoints % 2 == 0 || pa->npoints < 3 )
+ if (pa->npoints % 2 == 0 || pa->npoints < 3)
{
lwerror("lw_dist2d_pt_ptarrayarc called with non-arc input");
return LW_FALSE;
@@ -1103,20 +1181,20 @@
A1 = getPoint2d_cp(pa, 0);
- if ( ! lw_dist2d_pt_pt(p, A1, dl) )
+ if (!lw_dist2d_pt_pt(p, A1, dl))
return LW_FALSE;
- for ( t=1; t<pa->npoints; t += 2 )
+ for (t = 1; t < pa->npoints; t += 2)
{
dl->twisted = twist;
A2 = getPoint2d_cp(pa, t);
- A3 = getPoint2d_cp(pa, t+1);
+ A3 = getPoint2d_cp(pa, t + 1);
- if ( lw_dist2d_pt_arc(p, A1, A2, A3, dl) == LW_FALSE )
+ if (lw_dist2d_pt_arc(p, A1, A2, A3, dl) == LW_FALSE)
return LW_FALSE;
- if ( dl->distance <= dl->tolerance && dl->mode == DIST_MIN )
- return LW_TRUE; /*just a check if the answer is already given*/
+ if (dl->distance <= dl->tolerance && dl->mode == DIST_MIN)
+ return LW_TRUE; /*just a check if the answer is already given*/
A1 = A3;
}
@@ -1124,34 +1202,30 @@
return LW_TRUE;
}
-
-
-
/**
-* test each segment of l1 against each segment of l2.
-*/
+ * test each segment of l1 against each segment of l2.
+ */
int
-lw_dist2d_ptarray_ptarray(POINTARRAY *l1, POINTARRAY *l2,DISTPTS *dl)
+lw_dist2d_ptarray_ptarray(POINTARRAY *l1, POINTARRAY *l2, DISTPTS *dl)
{
- uint32_t t,u;
- const POINT2D *start, *end;
- const POINT2D *start2, *end2;
+ uint32_t t, u;
+ const POINT2D *start, *end;
+ const POINT2D *start2, *end2;
int twist = dl->twisted;
- LWDEBUGF(2, "lw_dist2d_ptarray_ptarray called (points: %d-%d)",l1->npoints, l2->npoints);
+ LWDEBUGF(2, "lw_dist2d_ptarray_ptarray called (points: %d-%d)", l1->npoints, l2->npoints);
- if (dl->mode == DIST_MAX)/*If we are searching for maxdistance we go straight to point-point calculation since the maxdistance have to be between two vertexes*/
+ /* If we are searching for maxdistance we go straight to point-point calculation since the maxdistance have
+ * to be between two vertexes*/
+ if (dl->mode == DIST_MAX)
{
- for (t=0; t<l1->npoints; t++) /*for each segment in L1 */
+ for (t = 0; t < l1->npoints; t++) /*for each segment in L1 */
{
start = getPoint2d_cp(l1, t);
- for (u=0; u<l2->npoints; u++) /*for each segment in L2 */
+ for (u = 0; u < l2->npoints; u++) /*for each segment in L2 */
{
start2 = getPoint2d_cp(l2, u);
lw_dist2d_pt_pt(start, start2, dl);
- LWDEBUGF(4, "maxdist_ptarray_ptarray; seg %i * seg %i, dist = %g\n",t,u,dl->distance);
- LWDEBUGF(3, " seg%d-seg%d dist: %f, mindist: %f",
- t, u, dl->distance, dl->tolerance);
}
}
}
@@ -1158,19 +1232,17 @@
else
{
start = getPoint2d_cp(l1, 0);
- for (t=1; t<l1->npoints; t++) /*for each segment in L1 */
+ for (t = 1; t < l1->npoints; t++) /*for each segment in L1 */
{
end = getPoint2d_cp(l1, t);
start2 = getPoint2d_cp(l2, 0);
- for (u=1; u<l2->npoints; u++) /*for each segment in L2 */
+ for (u = 1; u < l2->npoints; u++) /*for each segment in L2 */
{
end2 = getPoint2d_cp(l2, u);
- dl->twisted=twist;
+ dl->twisted = twist;
lw_dist2d_seg_seg(start, end, start2, end2, dl);
- LWDEBUGF(4, "mindist_ptarray_ptarray; seg %i * seg %i, dist = %g\n",t,u,dl->distance);
- LWDEBUGF(3, " seg%d-seg%d dist: %f, mindist: %f",
- t, u, dl->distance, dl->tolerance);
- if (dl->distance<=dl->tolerance && dl->mode == DIST_MIN) return LW_TRUE; /*just a check if the answer is already given*/
+ if (dl->distance <= dl->tolerance && dl->mode == DIST_MIN)
+ return LW_TRUE; /*just a check if the answer is already given*/
start2 = end2;
}
start = end;
@@ -1180,8 +1252,8 @@
}
/**
-* Test each segment of pa against each arc of pb for distance.
-*/
+ * Test each segment of pa against each arc of pb for distance.
+ */
int
lw_dist2d_ptarray_ptarrayarc(const POINTARRAY *pa, const POINTARRAY *pb, DISTPTS *dl)
{
@@ -1193,15 +1265,15 @@
const POINT2D *B3;
int twist = dl->twisted;
- LWDEBUGF(2, "lw_dist2d_ptarray_ptarrayarc called (points: %d-%d)",pa->npoints, pb->npoints);
+ LWDEBUGF(2, "lw_dist2d_ptarray_ptarrayarc called (points: %d-%d)", pa->npoints, pb->npoints);
- if ( pb->npoints % 2 == 0 || pb->npoints < 3 )
+ if (pb->npoints % 2 == 0 || pb->npoints < 3)
{
lwerror("lw_dist2d_ptarray_ptarrayarc called with non-arc input");
return LW_FALSE;
}
- if ( dl->mode == DIST_MAX )
+ if (dl->mode == DIST_MAX)
{
lwerror("lw_dist2d_ptarray_ptarrayarc does not currently support DIST_MAX mode");
return LW_FALSE;
@@ -1209,20 +1281,20 @@
else
{
A1 = getPoint2d_cp(pa, 0);
- for ( t=1; t < pa->npoints; t++ ) /* For each segment in pa */
+ for (t = 1; t < pa->npoints; t++) /* For each segment in pa */
{
A2 = getPoint2d_cp(pa, t);
B1 = getPoint2d_cp(pb, 0);
- for ( u=1; u < pb->npoints; u += 2 ) /* For each arc in pb */
+ for (u = 1; u < pb->npoints; u += 2) /* For each arc in pb */
{
B2 = getPoint2d_cp(pb, u);
- B3 = getPoint2d_cp(pb, u+1);
+ B3 = getPoint2d_cp(pb, u + 1);
dl->twisted = twist;
lw_dist2d_seg_arc(A1, A2, B1, B2, B3, dl);
/* If we've found a distance within tolerance, we're done */
- if ( dl->distance <= dl->tolerance && dl->mode == DIST_MIN )
+ if (dl->distance <= dl->tolerance && dl->mode == DIST_MIN)
return LW_TRUE;
B1 = B3;
@@ -1234,8 +1306,8 @@
}
/**
-* Test each arc of pa against each arc of pb for distance.
-*/
+ * Test each arc of pa against each arc of pb for distance.
+ */
int
lw_dist2d_ptarrayarc_ptarrayarc(const POINTARRAY *pa, const POINTARRAY *pb, DISTPTS *dl)
{
@@ -1248,7 +1320,7 @@
const POINT2D *B3;
int twist = dl->twisted;
- LWDEBUGF(2, "lw_dist2d_ptarrayarc_ptarrayarc called (points: %d-%d)",pa->npoints, pb->npoints);
+ LWDEBUGF(2, "lw_dist2d_ptarrayarc_ptarrayarc called (points: %d-%d)", pa->npoints, pb->npoints);
if (dl->mode == DIST_MAX)
{
@@ -1258,21 +1330,21 @@
else
{
A1 = getPoint2d_cp(pa, 0);
- for ( t=1; t < pa->npoints; t += 2 ) /* For each segment in pa */
+ for (t = 1; t < pa->npoints; t += 2) /* For each segment in pa */
{
A2 = getPoint2d_cp(pa, t);
- A3 = getPoint2d_cp(pa, t+1);
+ A3 = getPoint2d_cp(pa, t + 1);
B1 = getPoint2d_cp(pb, 0);
- for ( u=1; u < pb->npoints; u += 2 ) /* For each arc in pb */
+ for (u = 1; u < pb->npoints; u += 2) /* For each arc in pb */
{
B2 = getPoint2d_cp(pb, u);
- B3 = getPoint2d_cp(pb, u+1);
+ B3 = getPoint2d_cp(pb, u + 1);
dl->twisted = twist;
lw_dist2d_arc_arc(A1, A2, A3, B1, B2, B3, dl);
/* If we've found a distance within tolerance, we're done */
- if ( dl->distance <= dl->tolerance && dl->mode == DIST_MIN )
+ if (dl->distance <= dl->tolerance && dl->mode == DIST_MIN)
return LW_TRUE;
B1 = B3;
@@ -1284,25 +1356,30 @@
}
/**
-* Calculate the shortest distance between an arc and an edge.
-* Line/circle approach from http://stackoverflow.com/questions/1073336/circle-line-collision-detection
-*/
+ * Calculate the shortest distance between an arc and an edge.
+ * Line/circle approach from http://stackoverflow.com/questions/1073336/circle-line-collision-detection
+ */
int
-lw_dist2d_seg_arc(const POINT2D *A1, const POINT2D *A2, const POINT2D *B1, const POINT2D *B2, const POINT2D *B3, DISTPTS *dl)
+lw_dist2d_seg_arc(const POINT2D *A1,
+ const POINT2D *A2,
+ const POINT2D *B1,
+ const POINT2D *B2,
+ const POINT2D *B3,
+ DISTPTS *dl)
{
- POINT2D C; /* center of arc circle */
+ POINT2D C; /* center of arc circle */
double radius_C; /* radius of arc circle */
- POINT2D D; /* point on A closest to C */
+ POINT2D D; /* point on A closest to C */
double dist_C_D; /* distance from C to D */
int pt_in_arc, pt_in_seg;
DISTPTS dltmp;
/* Bail out on crazy modes */
- if ( dl->mode < 0 )
+ if (dl->mode < 0)
lwerror("lw_dist2d_seg_arc does not support maxdistance mode");
/* What if the "arc" is a point? */
- if ( lw_arc_is_pt(B1, B2, B3) )
+ if (lw_arc_is_pt(B1, B2, B3))
return lw_dist2d_pt_seg(B1, A1, A2, dl);
/* Calculate center and radius of the circle. */
@@ -1309,12 +1386,12 @@
radius_C = lw_arc_center(B1, B2, B3, &C);
/* This "arc" is actually a line (B2 is collinear with B1,B3) */
- if ( radius_C < 0.0 )
+ if (radius_C < 0.0)
return lw_dist2d_seg_seg(A1, A2, B1, B3, dl);
/* Calculate distance between the line and circle center */
lw_dist2d_distpts_init(&dltmp, DIST_MIN);
- if ( lw_dist2d_pt_seg(&C, A1, A2, &dltmp) == LW_FALSE )
+ if (lw_dist2d_pt_seg(&C, A1, A2, &dltmp) == LW_FALSE)
lwerror("lw_dist2d_pt_seg failed in lw_dist2d_seg_arc");
D = dltmp.p1;
@@ -1323,28 +1400,27 @@
/* Line intersects circle, maybe arc intersects edge? */
/* If so, that's the closest point. */
/* If not, the closest point is one of the end points of A */
- if ( dist_C_D < radius_C )
+ if (dist_C_D < radius_C)
{
- double length_A; /* length of the segment A */
- POINT2D E, F; /* points of intersection of edge A and circle(B) */
+ double length_A; /* length of the segment A */
+ POINT2D E, F; /* points of intersection of edge A and circle(B) */
double dist_D_EF; /* distance from D to E or F (same distance both ways) */
- dist_D_EF = sqrt(radius_C*radius_C - dist_C_D*dist_C_D);
- length_A = sqrt((A2->x-A1->x)*(A2->x-A1->x)+(A2->y-A1->y)*(A2->y-A1->y));
+ dist_D_EF = sqrt(radius_C * radius_C - dist_C_D * dist_C_D);
+ length_A = sqrt((A2->x - A1->x) * (A2->x - A1->x) + (A2->y - A1->y) * (A2->y - A1->y));
/* Point of intersection E */
- E.x = D.x - (A2->x-A1->x) * dist_D_EF / length_A;
- E.y = D.y - (A2->y-A1->y) * dist_D_EF / length_A;
+ E.x = D.x - (A2->x - A1->x) * dist_D_EF / length_A;
+ E.y = D.y - (A2->y - A1->y) * dist_D_EF / length_A;
/* Point of intersection F */
- F.x = D.x + (A2->x-A1->x) * dist_D_EF / length_A;
- F.y = D.y + (A2->y-A1->y) * dist_D_EF / length_A;
+ F.x = D.x + (A2->x - A1->x) * dist_D_EF / length_A;
+ F.y = D.y + (A2->y - A1->y) * dist_D_EF / length_A;
-
/* If E is within A and within B then it's an intersection point */
pt_in_arc = lw_pt_in_arc(&E, B1, B2, B3);
pt_in_seg = lw_pt_in_seg(&E, A1, A2);
- if ( pt_in_arc && pt_in_seg )
+ if (pt_in_arc && pt_in_seg)
{
dl->distance = 0.0;
dl->p1 = E;
@@ -1356,7 +1432,7 @@
pt_in_arc = lw_pt_in_arc(&F, B1, B2, B3);
pt_in_seg = lw_pt_in_seg(&F, A1, A2);
- if ( pt_in_arc && pt_in_seg )
+ if (pt_in_arc && pt_in_seg)
{
dl->distance = 0.0;
dl->p1 = F;
@@ -1368,7 +1444,7 @@
/* Line grazes circle, maybe arc intersects edge? */
/* If so, grazing point is the closest point. */
/* If not, the closest point is one of the end points of A */
- else if ( dist_C_D == radius_C )
+ else if (dist_C_D == radius_C)
{
/* Closest point D is also the point of grazing */
pt_in_arc = lw_pt_in_arc(&D, B1, B2, B3);
@@ -1375,7 +1451,7 @@
pt_in_seg = lw_pt_in_seg(&D, A1, A2);
/* Is D contained in both A and B? */
- if ( pt_in_arc && pt_in_seg )
+ if (pt_in_arc && pt_in_seg)
{
dl->distance = 0.0;
dl->p1 = D;
@@ -1389,16 +1465,15 @@
else
{
POINT2D G; /* Point on circle closest to A */
- G.x = C.x + (D.x-C.x) * radius_C / dist_C_D;
- G.y = C.y + (D.y-C.y) * radius_C / dist_C_D;
+ G.x = C.x + (D.x - C.x) * radius_C / dist_C_D;
+ G.y = C.y + (D.y - C.y) * radius_C / dist_C_D;
pt_in_arc = lw_pt_in_arc(&G, B1, B2, B3);
pt_in_seg = lw_pt_in_seg(&D, A1, A2);
/* Closest point is on the interior of A and B */
- if ( pt_in_arc && pt_in_seg )
+ if (pt_in_arc && pt_in_seg)
return lw_dist2d_pt_pt(&D, &G, dl);
-
}
/* Now we test the many combinations of end points with either */
@@ -1435,17 +1510,17 @@
}
int
-lw_dist2d_pt_arc(const POINT2D* P, const POINT2D* A1, const POINT2D* A2, const POINT2D* A3, DISTPTS* dl)
+lw_dist2d_pt_arc(const POINT2D *P, const POINT2D *A1, const POINT2D *A2, const POINT2D *A3, DISTPTS *dl)
{
double radius_A, d;
POINT2D C; /* center of circle defined by arc A */
POINT2D X; /* point circle(A) where line from C to P crosses */
- if ( dl->mode < 0 )
+ if (dl->mode < 0)
lwerror("lw_dist2d_pt_arc does not support maxdistance mode");
/* What if the arc is a point? */
- if ( lw_arc_is_pt(A1, A2, A3) )
+ if (lw_arc_is_pt(A1, A2, A3))
return lw_dist2d_pt_pt(P, A1, dl);
/* Calculate centers and radii of circles. */
@@ -1452,7 +1527,7 @@
radius_A = lw_arc_center(A1, A2, A3, &C);
/* This "arc" is actually a line (A2 is colinear with A1,A3) */
- if ( radius_A < 0.0 )
+ if (radius_A < 0.0)
return lw_dist2d_pt_seg(P, A1, A3, dl);
/* Distance from point to center */
@@ -1459,7 +1534,7 @@
d = distance2d_pt_pt(&C, P);
/* P is the center of the circle */
- if ( FP_EQUALS(d, 0.0) )
+ if (FP_EQUALS(d, 0.0))
{
dl->distance = radius_A;
dl->p1 = *A1;
@@ -1472,7 +1547,7 @@
X.y = C.y + (P->y - C.y) * radius_A / d;
/* Is crossing point inside the arc? Or arc is actually circle? */
- if ( p2d_same(A1, A3) || lw_pt_in_arc(&X, A1, A2, A3) )
+ if (p2d_same(A1, A3) || lw_pt_in_arc(&X, A1, A2, A3))
{
lw_dist2d_pt_pt(P, &X, dl);
}
@@ -1486,33 +1561,42 @@
}
/* Auxiliary function to calculate the distance between 2 concentric arcs*/
-int lw_dist2d_arc_arc_concentric( const POINT2D *A1, const POINT2D *A2,
- const POINT2D *A3, double radius_A,
- const POINT2D *B1, const POINT2D *B2,
- const POINT2D *B3, double radius_B,
- const POINT2D *CENTER, DISTPTS *dl);
+int lw_dist2d_arc_arc_concentric(const POINT2D *A1,
+ const POINT2D *A2,
+ const POINT2D *A3,
+ double radius_A,
+ const POINT2D *B1,
+ const POINT2D *B2,
+ const POINT2D *B3,
+ double radius_B,
+ const POINT2D *CENTER,
+ DISTPTS *dl);
int
-lw_dist2d_arc_arc(const POINT2D *A1, const POINT2D *A2, const POINT2D *A3,
- const POINT2D *B1, const POINT2D *B2, const POINT2D *B3,
- DISTPTS *dl)
+lw_dist2d_arc_arc(const POINT2D *A1,
+ const POINT2D *A2,
+ const POINT2D *A3,
+ const POINT2D *B1,
+ const POINT2D *B2,
+ const POINT2D *B3,
+ DISTPTS *dl)
{
- POINT2D CA, CB; /* Center points of arcs A and B */
+ POINT2D CA, CB; /* Center points of arcs A and B */
double radius_A, radius_B, d; /* Radii of arcs A and B */
- POINT2D D; /* Mid-point between the centers CA and CB */
+ POINT2D D; /* Mid-point between the centers CA and CB */
int pt_in_arc_A, pt_in_arc_B; /* Test whether potential intersection point is within the arc */
- if ( dl->mode != DIST_MIN )
+ if (dl->mode != DIST_MIN)
lwerror("lw_dist2d_arc_arc only supports mindistance");
/* TODO: Handle case where arc is closed circle (A1 = A3) */
/* What if one or both of our "arcs" is actually a point? */
- if ( lw_arc_is_pt(B1, B2, B3) && lw_arc_is_pt(A1, A2, A3) )
+ if (lw_arc_is_pt(B1, B2, B3) && lw_arc_is_pt(A1, A2, A3))
return lw_dist2d_pt_pt(B1, A1, dl);
- else if ( lw_arc_is_pt(B1, B2, B3) )
+ else if (lw_arc_is_pt(B1, B2, B3))
return lw_dist2d_pt_arc(B1, A1, A2, A3, dl);
- else if ( lw_arc_is_pt(A1, A2, A3) )
+ else if (lw_arc_is_pt(A1, A2, A3))
return lw_dist2d_pt_arc(A1, B1, B2, B3, dl);
/* Calculate centers and radii of circles. */
@@ -1520,15 +1604,15 @@
radius_B = lw_arc_center(B1, B2, B3, &CB);
/* Two co-linear arcs?!? That's two segments. */
- if ( radius_A < 0 && radius_B < 0 )
+ if (radius_A < 0 && radius_B < 0)
return lw_dist2d_seg_seg(A1, A3, B1, B3, dl);
/* A is co-linear, delegate to lw_dist_seg_arc here. */
- if ( radius_A < 0 )
+ if (radius_A < 0)
return lw_dist2d_seg_arc(A1, A3, B1, B2, B3, dl);
/* B is co-linear, delegate to lw_dist_seg_arc here. */
- if ( radius_B < 0 )
+ if (radius_B < 0)
return lw_dist2d_seg_arc(B1, B3, A1, A2, A3, dl);
/* Center-center distance */
@@ -1535,26 +1619,34 @@
d = distance2d_pt_pt(&CA, &CB);
/* Concentric arcs */
- if ( FP_EQUALS(d, 0.0) )
- return lw_dist2d_arc_arc_concentric(A1, A2, A3, radius_A,
- B1, B2, B3, radius_B,
- &CA, dl);
+ if (FP_EQUALS(d, 0.0))
+ return lw_dist2d_arc_arc_concentric(A1, A2, A3, radius_A, B1, B2, B3, radius_B, &CA, dl);
/* Make sure that arc "A" has the bigger radius */
- if ( radius_B > radius_A )
+ if (radius_B > radius_A)
{
const POINT2D *tmp;
POINT2D TP; /* Temporary point P */
double td;
- tmp = B1; B1 = A1; A1 = tmp;
- tmp = B2; B2 = A2; A2 = tmp;
- tmp = B3; B3 = A3; A3 = tmp;
- TP = CB; CB = CA; CA = TP;
- td = radius_B; radius_B = radius_A; radius_A = td;
+ tmp = B1;
+ B1 = A1;
+ A1 = tmp;
+ tmp = B2;
+ B2 = A2;
+ A2 = tmp;
+ tmp = B3;
+ B3 = A3;
+ A3 = tmp;
+ TP = CB;
+ CB = CA;
+ CA = TP;
+ td = radius_B;
+ radius_B = radius_A;
+ radius_A = td;
}
/* Circles touch at a point. Is that point within the arcs? */
- if ( d == (radius_A + radius_B) )
+ if (d == (radius_A + radius_B))
{
D.x = CA.x + (CB.x - CA.x) * radius_A / d;
D.y = CA.y + (CB.y - CA.y) * radius_A / d;
@@ -1563,7 +1655,7 @@
pt_in_arc_B = lw_pt_in_arc(&D, B1, B2, B3);
/* Arcs do touch at D, return it */
- if ( pt_in_arc_A && pt_in_arc_B )
+ if (pt_in_arc_A && pt_in_arc_B)
{
dl->distance = 0.0;
dl->p1 = D;
@@ -1573,7 +1665,7 @@
}
/* Disjoint or contained circles don't intersect. Closest point may be on */
/* the line joining CA to CB. */
- else if ( d > (radius_A + radius_B) /* Disjoint */ || d < (radius_A - radius_B) /* Contained */ )
+ else if (d > (radius_A + radius_B) /* Disjoint */ || d < (radius_A - radius_B) /* Contained */)
{
POINT2D XA, XB; /* Points where the line from CA to CB cross their circle bounds */
@@ -1590,7 +1682,7 @@
/* If the nearest points are both within the arcs, that's our answer */
/* the shortest distance is at the nearest points */
- if ( pt_in_arc_A && pt_in_arc_B )
+ if (pt_in_arc_A && pt_in_arc_B)
{
return lw_dist2d_pt_pt(&XA, &XB, dl);
}
@@ -1597,13 +1689,13 @@
}
/* Circles cross at two points, are either of those points in both arcs? */
/* http://paulbourke.net/geometry/2circle/ */
- else if ( d < (radius_A + radius_B) )
+ else if (d < (radius_A + radius_B))
{
POINT2D E, F; /* Points where circle(A) and circle(B) cross */
/* Distance from CA to D */
- double a = (radius_A*radius_A - radius_B*radius_B + d*d) / (2*d);
+ double a = (radius_A * radius_A - radius_B * radius_B + d * d) / (2 * d);
/* Distance from D to E or F */
- double h = sqrt(radius_A*radius_A - a*a);
+ double h = sqrt(radius_A * radius_A - a * a);
/* Location of D */
D.x = CA.x + (CB.x - CA.x) * a / d;
@@ -1617,7 +1709,7 @@
pt_in_arc_A = lw_pt_in_arc(&E, A1, A2, A3);
pt_in_arc_B = lw_pt_in_arc(&E, B1, B2, B3);
- if ( pt_in_arc_A && pt_in_arc_B )
+ if (pt_in_arc_A && pt_in_arc_B)
{
dl->p1 = dl->p2 = E;
dl->distance = 0.0;
@@ -1632,7 +1724,7 @@
pt_in_arc_A = lw_pt_in_arc(&F, A1, A2, A3);
pt_in_arc_B = lw_pt_in_arc(&F, B1, B2, B3);
- if ( pt_in_arc_A && pt_in_arc_B )
+ if (pt_in_arc_A && pt_in_arc_B)
{
dl->p1 = dl->p2 = F;
dl->distance = 0.0;
@@ -1655,7 +1747,7 @@
}
/* Closest point is in the arc B, but not in the arc A, so */
/* one of the A end points must be the closest. */
- else if ( pt_in_arc_B && ! pt_in_arc_A )
+ else if (pt_in_arc_B && !pt_in_arc_A)
{
lw_dist2d_pt_arc(A1, B1, B2, B3, dl);
lw_dist2d_pt_arc(A3, B1, B2, B3, dl);
@@ -1675,11 +1767,16 @@
}
int
-lw_dist2d_arc_arc_concentric( const POINT2D *A1, const POINT2D *A2,
- const POINT2D *A3, double radius_A,
- const POINT2D *B1, const POINT2D *B2,
- const POINT2D *B3, double radius_B,
- const POINT2D *CENTER, DISTPTS *dl)
+lw_dist2d_arc_arc_concentric(const POINT2D *A1,
+ const POINT2D *A2,
+ const POINT2D *A3,
+ double radius_A,
+ const POINT2D *B1,
+ const POINT2D *B2,
+ const POINT2D *B3,
+ double radius_B,
+ const POINT2D *CENTER,
+ DISTPTS *dl)
{
int seg_size;
double dist_sqr, shortest_sqr;
@@ -1819,33 +1916,31 @@
int
lw_dist2d_seg_seg(const POINT2D *A, const POINT2D *B, const POINT2D *C, const POINT2D *D, DISTPTS *dl)
{
- double s_top, s_bot,s;
- double r_top, r_bot,r;
+ double s_top, s_bot, s;
+ double r_top, r_bot, r;
- LWDEBUGF(2, "lw_dist2d_seg_seg [%g,%g]->[%g,%g] by [%g,%g]->[%g,%g]",
- A->x,A->y,B->x,B->y, C->x,C->y, D->x, D->y);
-
/*A and B are the same point */
- if ( ( A->x == B->x) && (A->y == B->y) )
+ if ((A->x == B->x) && (A->y == B->y))
{
- return lw_dist2d_pt_seg(A,C,D,dl);
+ return lw_dist2d_pt_seg(A, C, D, dl);
}
+
/*U and V are the same point */
-
- if ( ( C->x == D->x) && (C->y == D->y) )
+ if ((C->x == D->x) && (C->y == D->y))
{
- dl->twisted= ((dl->twisted) * (-1));
- return lw_dist2d_pt_seg(D,A,B,dl);
+ dl->twisted = ((dl->twisted) * (-1));
+ return lw_dist2d_pt_seg(D, A, B, dl);
}
+
/* AB and CD are line segments */
/* from comp.graphics.algo
Solving the above for r and s yields
(Ay-Cy)(Dx-Cx)-(Ax-Cx)(Dy-Cy)
- r = ----------------------------- (eqn 1)
+ r = ----------------------------- (eqn 1)
(Bx-Ax)(Dy-Cy)-(By-Ay)(Dx-Cx)
- (Ay-Cy)(Bx-Ax)-(Ax-Cx)(By-Ay)
+ (Ay-Cy)(Bx-Ax)-(Ax-Cx)(By-Ay)
s = ----------------------------- (eqn 2)
(Bx-Ax)(Dy-Cy)-(By-Ay)(Dx-Cx)
Let P be the position vector of the intersection point, then
@@ -1859,52 +1954,52 @@
If the numerator in eqn 1 is also zero, AB & CD are collinear.
*/
- r_top = (A->y-C->y)*(D->x-C->x) - (A->x-C->x)*(D->y-C->y);
- r_bot = (B->x-A->x)*(D->y-C->y) - (B->y-A->y)*(D->x-C->x);
+ r_top = (A->y - C->y) * (D->x - C->x) - (A->x - C->x) * (D->y - C->y);
+ r_bot = (B->x - A->x) * (D->y - C->y) - (B->y - A->y) * (D->x - C->x);
- s_top = (A->y-C->y)*(B->x-A->x) - (A->x-C->x)*(B->y-A->y);
- s_bot = (B->x-A->x)*(D->y-C->y) - (B->y-A->y)*(D->x-C->x);
+ s_top = (A->y - C->y) * (B->x - A->x) - (A->x - C->x) * (B->y - A->y);
+ s_bot = (B->x - A->x) * (D->y - C->y) - (B->y - A->y) * (D->x - C->x);
- if ( (r_bot==0) || (s_bot == 0) )
+ if ((r_bot == 0) || (s_bot == 0))
{
- if ((lw_dist2d_pt_seg(A,C,D,dl)) && (lw_dist2d_pt_seg(B,C,D,dl)))
+ if ((lw_dist2d_pt_seg(A, C, D, dl)) && (lw_dist2d_pt_seg(B, C, D, dl)))
{
- dl->twisted= ((dl->twisted) * (-1)); /*here we change the order of inputted geometries and that we notice by changing sign on dl->twisted*/
- return ((lw_dist2d_pt_seg(C,A,B,dl)) && (lw_dist2d_pt_seg(D,A,B,dl))); /*if all is successful we return true*/
+ /* change the order of inputted geometries and that we notice by changing sign on dl->twisted*/
+ dl->twisted *= -1;
+ return ((lw_dist2d_pt_seg(C, A, B, dl)) && (lw_dist2d_pt_seg(D, A, B, dl)));
}
else
- {
return LW_FALSE; /* if any of the calls to lw_dist2d_pt_seg goes wrong we return false*/
- }
}
- s = s_top/s_bot;
- r= r_top/r_bot;
+ s = s_top / s_bot;
+ r = r_top / r_bot;
- if (((r<0) || (r>1) || (s<0) || (s>1)) || (dl->mode == DIST_MAX))
+ if (((r < 0) || (r > 1) || (s < 0) || (s > 1)) || (dl->mode == DIST_MAX))
{
- if ((lw_dist2d_pt_seg(A,C,D,dl)) && (lw_dist2d_pt_seg(B,C,D,dl)))
+ if ((lw_dist2d_pt_seg(A, C, D, dl)) && (lw_dist2d_pt_seg(B, C, D, dl)))
{
- dl->twisted= ((dl->twisted) * (-1)); /*here we change the order of inputted geometries and that we notice by changing sign on dl->twisted*/
- return ((lw_dist2d_pt_seg(C,A,B,dl)) && (lw_dist2d_pt_seg(D,A,B,dl))); /*if all is successful we return true*/
+ /* change the order of inputted geometries and that we notice by changing sign on dl->twisted*/
+ dl->twisted *= -1;
+ return ((lw_dist2d_pt_seg(C, A, B, dl)) && (lw_dist2d_pt_seg(D, A, B, dl)));
}
else
- {
return LW_FALSE; /* if any of the calls to lw_dist2d_pt_seg goes wrong we return false*/
- }
}
else
{
- if (dl->mode == DIST_MIN) /*If there is intersection we identify the intersection point and return it but only if we are looking for mindistance*/
+ /* If there is intersection we identify the intersection point and return it but only if we are looking
+ * for mindistance */
+ if (dl->mode == DIST_MIN)
{
POINT2D theP;
- if (((A->x==C->x)&&(A->y==C->y))||((A->x==D->x)&&(A->y==D->y)))
+ if (((A->x == C->x) && (A->y == C->y)) || ((A->x == D->x) && (A->y == D->y)))
{
theP.x = A->x;
theP.y = A->y;
}
- else if (((B->x==C->x)&&(B->y==C->y))||((B->x==D->x)&&(B->y==D->y)))
+ else if (((B->x == C->x) && (B->y == C->y)) || ((B->x == D->x) && (B->y == D->y)))
{
theP.x = B->x;
theP.y = B->y;
@@ -1911,23 +2006,21 @@
}
else
{
- theP.x = A->x+r*(B->x-A->x);
- theP.y = A->y+r*(B->y-A->y);
+ theP.x = A->x + r * (B->x - A->x);
+ theP.y = A->y + r * (B->y - A->y);
}
- dl->distance=0.0;
- dl->p1=theP;
- dl->p2=theP;
+ dl->distance = 0.0;
+ dl->p1 = theP;
+ dl->p2 = theP;
}
return LW_TRUE;
}
}
-
/*------------------------------------------------------------------------------------------------------------
End of Brute force functions
--------------------------------------------------------------------------------------------------------------*/
-
/*------------------------------------------------------------------------------------------------------------
New faster distance calculations
--------------------------------------------------------------------------------------------------------------*/
@@ -1940,22 +2033,22 @@
chosen selection of the points not all of them
*/
int
-lw_dist2d_fast_ptarray_ptarray(POINTARRAY *l1, POINTARRAY *l2,DISTPTS *dl, GBOX *box1, GBOX *box2)
+lw_dist2d_fast_ptarray_ptarray(POINTARRAY *l1, POINTARRAY *l2, DISTPTS *dl, GBOX *box1, GBOX *box2)
{
/*here we define two lists to hold our calculated "z"-values and the order number in the geometry*/
double k, thevalue;
- float deltaX, deltaY, c1m, c2m;
- POINT2D c1, c2;
+ float deltaX, deltaY, c1m, c2m;
+ POINT2D c1, c2;
const POINT2D *theP;
- float min1X, max1X, max1Y, min1Y,min2X, max2X, max2Y, min2Y;
+ float min1X, max1X, max1Y, min1Y, min2X, max2X, max2Y, min2Y;
int t;
int n1 = l1->npoints;
int n2 = l2->npoints;
LISTSTRUCT *list1, *list2;
- list1 = (LISTSTRUCT*)lwalloc(sizeof(LISTSTRUCT)*n1);
- list2 = (LISTSTRUCT*)lwalloc(sizeof(LISTSTRUCT)*n2);
+ list1 = (LISTSTRUCT *)lwalloc(sizeof(LISTSTRUCT) * n1);
+ list2 = (LISTSTRUCT *)lwalloc(sizeof(LISTSTRUCT) * n2);
LWDEBUG(2, "lw_dist2d_fast_ptarray_ptarray is called");
@@ -1968,64 +2061,62 @@
max2Y = box2->ymax;
min2Y = box2->ymin;
/*we want the center of the bboxes, and calculate the slope between the centerpoints*/
- c1.x = min1X + (max1X-min1X)/2;
- c1.y = min1Y + (max1Y-min1Y)/2;
- c2.x = min2X + (max2X-min2X)/2;
- c2.y = min2Y + (max2Y-min2Y)/2;
+ c1.x = min1X + (max1X - min1X) / 2;
+ c1.y = min1Y + (max1Y - min1Y) / 2;
+ c2.x = min2X + (max2X - min2X) / 2;
+ c2.y = min2Y + (max2Y - min2Y) / 2;
- deltaX=(c2.x-c1.x);
- deltaY=(c2.y-c1.y);
+ deltaX = (c2.x - c1.x);
+ deltaY = (c2.y - c1.y);
-
/*Here we calculate where the line perpendicular to the center-center line crosses the axes for each vertex
- if the center-center line is vertical the perpendicular line will be horizontal and we find it's crossing the Y-axes with z = y-kx */
- if ((deltaX*deltaX)<(deltaY*deltaY)) /*North or South*/
+ if the center-center line is vertical the perpendicular line will be horizontal and we find it's crossing the
+ Y-axes with z = y-kx */
+ if ((deltaX * deltaX) < (deltaY * deltaY)) /*North or South*/
{
- k = -deltaX/deltaY;
- for (t=0; t<n1; t++) /*for each segment in L1 */
+ k = -deltaX / deltaY;
+ for (t = 0; t < n1; t++) /*for each segment in L1 */
{
theP = getPoint2d_cp(l1, t);
thevalue = theP->y - (k * theP->x);
- list1[t].themeasure=thevalue;
- list1[t].pnr=t;
-
+ list1[t].themeasure = thevalue;
+ list1[t].pnr = t;
}
- for (t=0; t<n2; t++) /*for each segment in L2*/
+ for (t = 0; t < n2; t++) /*for each segment in L2*/
{
theP = getPoint2d_cp(l2, t);
thevalue = theP->y - (k * theP->x);
- list2[t].themeasure=thevalue;
- list2[t].pnr=t;
-
+ list2[t].themeasure = thevalue;
+ list2[t].pnr = t;
}
- c1m = c1.y-(k*c1.x);
- c2m = c2.y-(k*c2.x);
+ c1m = c1.y - (k * c1.x);
+ c2m = c2.y - (k * c2.x);
}
-
- /*if the center-center line is horizontal the perpendicular line will be vertical. To eliminate problems with dividing by zero we are here mirroring the coordinate-system
- and we find it's crossing the X-axes with z = x-(1/k)y */
- else /*West or East*/
+ /*if the center-center line is horizontal the perpendicular line will be vertical. To eliminate problems with
+ dividing by zero we are here mirroring the coordinate-system and we find it's crossing the X-axes with z =
+ x-(1/k)y */
+ else /*West or East*/
{
- k = -deltaY/deltaX;
- for (t=0; t<n1; t++) /*for each segment in L1 */
+ k = -deltaY / deltaX;
+ for (t = 0; t < n1; t++) /*for each segment in L1 */
{
theP = getPoint2d_cp(l1, t);
thevalue = theP->x - (k * theP->y);
- list1[t].themeasure=thevalue;
- list1[t].pnr=t;
+ list1[t].themeasure = thevalue;
+ list1[t].pnr = t;
/* lwnotice("l1 %d, measure=%f",t,thevalue ); */
}
- for (t=0; t<n2; t++) /*for each segment in L2*/
+ for (t = 0; t < n2; t++) /*for each segment in L2*/
{
theP = getPoint2d_cp(l2, t);
thevalue = theP->x - (k * theP->y);
- list2[t].themeasure=thevalue;
- list2[t].pnr=t;
+ list2[t].themeasure = thevalue;
+ list2[t].pnr = t;
/* lwnotice("l2 %d, measure=%f",t,thevalue ); */
}
- c1m = c1.x-(k*c1.y);
- c2m = c2.x-(k*c2.y);
+ c1m = c1.x - (k * c1.y);
+ c2m = c2.x - (k * c2.y);
}
/*we sort our lists by the calculated values*/
@@ -2034,7 +2125,7 @@
if (c1m < c2m)
{
- if (!lw_dist2d_pre_seg_seg(l1,l2,list1,list2,k,dl))
+ if (!lw_dist2d_pre_seg_seg(l1, l2, list1, list2, k, dl))
{
lwfree(list1);
lwfree(list2);
@@ -2043,8 +2134,8 @@
}
else
{
- dl->twisted= ((dl->twisted) * (-1));
- if (!lw_dist2d_pre_seg_seg(l2,l1,list2,list1,k,dl))
+ dl->twisted = ((dl->twisted) * (-1));
+ if (!lw_dist2d_pre_seg_seg(l2, l1, list2, list1, k, dl))
{
lwfree(list1);
lwfree(list2);
@@ -2059,22 +2150,19 @@
int
struct_cmp_by_measure(const void *a, const void *b)
{
- LISTSTRUCT *ia = (LISTSTRUCT*)a;
- LISTSTRUCT *ib = (LISTSTRUCT*)b;
- return
- (ia->themeasure > ib->themeasure) ? 1 : ((ia->themeasure < ib->themeasure) ? -1 : 0);
+ LISTSTRUCT *ia = (LISTSTRUCT *)a;
+ LISTSTRUCT *ib = (LISTSTRUCT *)b;
+ return (ia->themeasure > ib->themeasure) ? 1 : ((ia->themeasure < ib->themeasure) ? -1 : 0);
}
-/**
- preparation before lw_dist2d_seg_seg.
-*/
+/** preparation before lw_dist2d_seg_seg. */
int
-lw_dist2d_pre_seg_seg(POINTARRAY *l1, POINTARRAY *l2,LISTSTRUCT *list1, LISTSTRUCT *list2,double k, DISTPTS *dl)
+lw_dist2d_pre_seg_seg(POINTARRAY *l1, POINTARRAY *l2, LISTSTRUCT *list1, LISTSTRUCT *list2, double k, DISTPTS *dl)
{
const POINT2D *p1, *p2, *p3, *p4, *p01, *p02;
- int pnr1,pnr2,pnr3,pnr4, n1, n2, i, u, r, twist;
+ int pnr1, pnr2, pnr3, pnr4, n1, n2, i, u, r, twist;
double maxmeasure;
- n1= l1->npoints;
+ n1 = l1->npoints;
n2 = l2->npoints;
LWDEBUG(2, "lw_dist2d_pre_seg_seg is called");
@@ -2082,66 +2170,85 @@
p1 = getPoint2d_cp(l1, list1[0].pnr);
p3 = getPoint2d_cp(l2, list2[0].pnr);
lw_dist2d_pt_pt(p1, p3, dl);
- maxmeasure = sqrt(dl->distance*dl->distance + (dl->distance*dl->distance*k*k));
+ maxmeasure = sqrt(dl->distance * dl->distance + (dl->distance * dl->distance * k * k));
twist = dl->twisted; /*to keep the incoming order between iterations*/
- for (i =(n1-1); i>=0; --i)
+ for (i = (n1 - 1); i >= 0; --i)
{
- /*we break this iteration when we have checked every
- point closer to our perpendicular "checkline" than
- our shortest found distance*/
- if (((list2[0].themeasure-list1[i].themeasure)) > maxmeasure) break;
- for (r=-1; r<=1; r +=2) /*because we are not iterating in the original point order we have to check the segment before and after every point*/
+ /*we break this iteration when we have checked every point closer to our perpendicular "checkline" than
+ * our shortest found distance*/
+ if (((list2[0].themeasure - list1[i].themeasure)) > maxmeasure)
+ break;
+ /*because we are not iterating in the original point order we have to check the segment before and after
+ * every point*/
+ for (r = -1; r <= 1; r += 2)
{
pnr1 = list1[i].pnr;
p1 = getPoint2d_cp(l1, pnr1);
- if (pnr1+r<0)
+ if (pnr1 + r < 0)
{
- p01 = getPoint2d_cp(l1, (n1-1));
- if (( p1->x == p01->x) && (p1->y == p01->y)) pnr2 = (n1-1);
- else pnr2 = pnr1; /* if it is a line and the last and first point is not the same we avoid the edge between start and end this way*/
+ p01 = getPoint2d_cp(l1, (n1 - 1));
+ if ((p1->x == p01->x) && (p1->y == p01->y))
+ pnr2 = (n1 - 1);
+ else
+ pnr2 = pnr1; /* if it is a line and the last and first point is not the same we
+ avoid the edge between start and end this way*/
}
- else if (pnr1+r>(n1-1))
+ else if (pnr1 + r > (n1 - 1))
{
p01 = getPoint2d_cp(l1, 0);
- if (( p1->x == p01->x) && (p1->y == p01->y)) pnr2 = 0;
- else pnr2 = pnr1; /* if it is a line and the last and first point is not the same we avoid the edge between start and end this way*/
+ if ((p1->x == p01->x) && (p1->y == p01->y))
+ pnr2 = 0;
+ else
+ pnr2 = pnr1; /* if it is a line and the last and first point is not the same we
+ avoid the edge between start and end this way*/
}
- else pnr2 = pnr1+r;
+ else
+ pnr2 = pnr1 + r;
-
p2 = getPoint2d_cp(l1, pnr2);
- for (u=0; u<n2; ++u)
+ for (u = 0; u < n2; ++u)
{
- if (((list2[u].themeasure-list1[i].themeasure)) >= maxmeasure) break;
+ if (((list2[u].themeasure - list1[i].themeasure)) >= maxmeasure)
+ break;
pnr3 = list2[u].pnr;
p3 = getPoint2d_cp(l2, pnr3);
- if (pnr3==0)
+ if (pnr3 == 0)
{
- p02 = getPoint2d_cp(l2, (n2-1));
- if (( p3->x == p02->x) && (p3->y == p02->y)) pnr4 = (n2-1);
- else pnr4 = pnr3; /* if it is a line and the last and first point is not the same we avoid the edge between start and end this way*/
+ p02 = getPoint2d_cp(l2, (n2 - 1));
+ if ((p3->x == p02->x) && (p3->y == p02->y))
+ pnr4 = (n2 - 1);
+ else
+ pnr4 = pnr3; /* if it is a line and the last and first point is not the
+ same we avoid the edge between start and end this way*/
}
- else pnr4 = pnr3-1;
+ else
+ pnr4 = pnr3 - 1;
p4 = getPoint2d_cp(l2, pnr4);
- dl->twisted=twist;
- if (!lw_dist2d_selected_seg_seg(p1, p2, p3, p4, dl)) return LW_FALSE;
+ dl->twisted = twist;
+ if (!lw_dist2d_selected_seg_seg(p1, p2, p3, p4, dl))
+ return LW_FALSE;
- if (pnr3>=(n2-1))
+ if (pnr3 >= (n2 - 1))
{
p02 = getPoint2d_cp(l2, 0);
- if (( p3->x == p02->x) && (p3->y == p02->y)) pnr4 = 0;
- else pnr4 = pnr3; /* if it is a line and the last and first point is not the same we avoid the edge between start and end this way*/
+ if ((p3->x == p02->x) && (p3->y == p02->y))
+ pnr4 = 0;
+ else
+ pnr4 = pnr3; /* if it is a line and the last and first point is not the
+ same we avoid the edge between start and end this way*/
}
- else pnr4 = pnr3+1;
+ else
+ pnr4 = pnr3 + 1;
p4 = getPoint2d_cp(l2, pnr4);
- dl->twisted=twist; /*we reset the "twist" for each iteration*/
- if (!lw_dist2d_selected_seg_seg(p1, p2, p3, p4, dl)) return LW_FALSE;
-
- maxmeasure = sqrt(dl->distance*dl->distance + (dl->distance*dl->distance*k*k));/*here we "translate" the found mindistance so it can be compared to our "z"-values*/
+ dl->twisted = twist; /*we reset the "twist" for each iteration*/
+ if (!lw_dist2d_selected_seg_seg(p1, p2, p3, p4, dl))
+ return LW_FALSE;
+ /*here we "translate" the found mindistance so it can be compared to our "z"-values*/
+ maxmeasure = sqrt(dl->distance * dl->distance + (dl->distance * dl->distance * k * k));
}
}
}
@@ -2149,7 +2256,6 @@
return LW_TRUE;
}
-
/**
This is the same function as lw_dist2d_seg_seg but
without any calculations to determine intersection since we
@@ -2158,31 +2264,27 @@
int
lw_dist2d_selected_seg_seg(const POINT2D *A, const POINT2D *B, const POINT2D *C, const POINT2D *D, DISTPTS *dl)
{
- LWDEBUGF(2, "lw_dist2d_selected_seg_seg [%g,%g]->[%g,%g] by [%g,%g]->[%g,%g]",
- A->x,A->y,B->x,B->y, C->x,C->y, D->x, D->y);
-
/*A and B are the same point */
- if ( ( A->x == B->x) && (A->y == B->y) )
+ if ((A->x == B->x) && (A->y == B->y))
{
- return lw_dist2d_pt_seg(A,C,D,dl);
+ return lw_dist2d_pt_seg(A, C, D, dl);
}
/*U and V are the same point */
- if ( ( C->x == D->x) && (C->y == D->y) )
+ if ((C->x == D->x) && (C->y == D->y))
{
- dl->twisted= ((dl->twisted) * (-1));
- return lw_dist2d_pt_seg(D,A,B,dl);
+ dl->twisted *= -1;
+ return lw_dist2d_pt_seg(D, A, B, dl);
}
- if ((lw_dist2d_pt_seg(A,C,D,dl)) && (lw_dist2d_pt_seg(B,C,D,dl)))
+ if ((lw_dist2d_pt_seg(A, C, D, dl)) && (lw_dist2d_pt_seg(B, C, D, dl)))
{
- dl->twisted= ((dl->twisted) * (-1)); /*here we change the order of inputted geometrys and that we notice by changing sign on dl->twisted*/
- return ((lw_dist2d_pt_seg(C,A,B,dl)) && (lw_dist2d_pt_seg(D,A,B,dl))); /*if all is successful we return true*/
+ /* change the order of inputted geometries and that we notice by changing sign on dl->twisted */
+ dl->twisted *= -1;
+ return ((lw_dist2d_pt_seg(C, A, B, dl)) && (lw_dist2d_pt_seg(D, A, B, dl)));
}
else
- {
return LW_FALSE; /* if any of the calls to lw_dist2d_pt_seg goes wrong we return false*/
- }
}
/*------------------------------------------------------------------------------------------------------------
@@ -2189,7 +2291,6 @@
End of New faster distance calculations
--------------------------------------------------------------------------------------------------------------*/
-
/*------------------------------------------------------------------------------------------------------------
Functions in common for Brute force and new calculation
--------------------------------------------------------------------------------------------------------------*/
@@ -2205,12 +2306,11 @@
lw_dist2d_pt_seg(const POINT2D *p, const POINT2D *A, const POINT2D *B, DISTPTS *dl)
{
POINT2D c;
- double r;
+ double r;
/*if start==end, then use pt distance */
- if ( ( A->x == B->x) && (A->y == B->y) )
- {
- return lw_dist2d_pt_pt(p,A,dl);
- }
+ if ((A->x == B->x) && (A->y == B->y))
+ return lw_dist2d_pt_pt(p, A, dl);
+
/*
* otherwise, we use comp.graphics.algorithms
* Frequently Asked Questions method
@@ -2226,35 +2326,26 @@
* 0<r<1 P is interior to AB
*/
- r = ( (p->x-A->x) * (B->x-A->x) + (p->y-A->y) * (B->y-A->y) )/( (B->x-A->x)*(B->x-A->x) +(B->y-A->y)*(B->y-A->y) );
+ r = ((p->x - A->x) * (B->x - A->x) + (p->y - A->y) * (B->y - A->y)) /
+ ((B->x - A->x) * (B->x - A->x) + (B->y - A->y) * (B->y - A->y));
/*This is for finding the maxdistance.
- the maxdistance have to be between two vertexes,
- compared to mindistance which can be between
- two vertexes vertex.*/
+ the maxdistance have to be between two vertexes, compared to mindistance which can be between two vertexes.*/
if (dl->mode == DIST_MAX)
{
- if (r>=0.5)
- {
- return lw_dist2d_pt_pt(p,A,dl);
- }
- if (r<0.5)
- {
- return lw_dist2d_pt_pt(p,B,dl);
- }
+ if (r >= 0.5)
+ return lw_dist2d_pt_pt(p, A, dl);
+ else /* (r < 0.5) */
+ return lw_dist2d_pt_pt(p, B, dl);
}
- if (r<0) /*If p projected on the line is outside point A*/
- {
- return lw_dist2d_pt_pt(p,A,dl);
- }
- if (r>=1) /*If p projected on the line is outside point B or on point B*/
- {
- return lw_dist2d_pt_pt(p,B,dl);
- }
+ if (r < 0) /*If p projected on the line is outside point A*/
+ return lw_dist2d_pt_pt(p, A, dl);
+ if (r >= 1) /*If p projected on the line is outside point B or on point B*/
+ return lw_dist2d_pt_pt(p, B, dl);
/*If the point p is on the segment this is a more robust way to find out that*/
- if (( ((A->y-p->y)*(B->x-A->x)==(A->x-p->x)*(B->y-A->y) ) ) && (dl->mode == DIST_MIN))
+ if ((((A->y - p->y) * (B->x - A->x) == (A->x - p->x) * (B->y - A->y))) && (dl->mode == DIST_MIN))
{
dl->distance = 0.0;
dl->p1 = *p;
@@ -2263,32 +2354,29 @@
/*If the projection of point p on the segment is between A and B
then we find that "point on segment" and send it to lw_dist2d_pt_pt*/
- c.x=A->x + r * (B->x-A->x);
- c.y=A->y + r * (B->y-A->y);
+ c.x = A->x + r * (B->x - A->x);
+ c.y = A->y + r * (B->y - A->y);
- return lw_dist2d_pt_pt(p,&c,dl);
+ return lw_dist2d_pt_pt(p, &c, dl);
}
-
-/**
-
-Compares incoming points and
-stores the points closest to each other
-or most far away from each other
-depending on dl->mode (max or min)
-*/
+/** Compares incoming points and stores the points closest to each other or most far away from each other depending on
+ * dl->mode (max or min) */
int
lw_dist2d_pt_pt(const POINT2D *thep1, const POINT2D *thep2, DISTPTS *dl)
{
double hside = thep2->x - thep1->x;
double vside = thep2->y - thep1->y;
- double dist = sqrt ( hside*hside + vside*vside );
+ double dist = sqrt(hside * hside + vside * vside);
- if (((dl->distance - dist)*(dl->mode))>0) /*multiplication with mode to handle mindistance (mode=1) and maxdistance (mode = (-1)*/
+ /*multiplication with mode to handle mindistance (mode=1) and maxdistance (mode = (-1)*/
+ if (((dl->distance - dist) * (dl->mode)) > 0)
{
dl->distance = dist;
- if (dl->twisted>0) /*To get the points in right order. twisted is updated between 1 and (-1) every time the order is changed earlier in the chain*/
+ /* To get the points in right order. twisted is updated between 1 and (-1) every time the order is
+ * changed earlier in the chain*/
+ if (dl->twisted > 0)
{
dl->p1 = *thep1;
dl->p2 = *thep2;
@@ -2302,9 +2390,6 @@
return LW_TRUE;
}
-
-
-
/*------------------------------------------------------------------------------------------------------------
End of Functions in common for Brute force and new calculation
--------------------------------------------------------------------------------------------------------------*/
@@ -2321,11 +2406,11 @@
double
distance2d_pt_seg(const POINT2D *p, const POINT2D *A, const POINT2D *B)
{
- double r,s;
+ double r, s;
/*if start==end, then use pt distance */
- if ( ( A->x == B->x) && (A->y == B->y) )
- return distance2d_pt_pt(p,A);
+ if ((A->x == B->x) && (A->y == B->y))
+ return distance2d_pt_pt(p, A);
/*
* otherwise, we use comp.graphics.algorithms
@@ -2332,8 +2417,8 @@
* Frequently Asked Questions method
*
* (1) AC dot AB
- * r = ---------
- * ||AB||^2
+ * r = ---------
+ * ||AB||^2
* r has the following meaning:
* r=0 P = A
* r=1 P = B
@@ -2342,12 +2427,14 @@
* 0<r<1 P is interior to AB
*/
- r = ( (p->x-A->x) * (B->x-A->x) + (p->y-A->y) * (B->y-A->y) )/( (B->x-A->x)*(B->x-A->x) +(B->y-A->y)*(B->y-A->y) );
+ r = ((p->x - A->x) * (B->x - A->x) + (p->y - A->y) * (B->y - A->y)) /
+ ((B->x - A->x) * (B->x - A->x) + (B->y - A->y) * (B->y - A->y));
- if (r<0) return distance2d_pt_pt(p,A);
- if (r>1) return distance2d_pt_pt(p,B);
+ if (r < 0)
+ return distance2d_pt_pt(p, A);
+ if (r > 1)
+ return distance2d_pt_pt(p, B);
-
/*
* (2)
* (Ay-Cy)(Bx-Ax)-(Ax-Cx)(By-Ay)
@@ -2358,12 +2445,10 @@
*
*/
- s = ( (A->y-p->y)*(B->x-A->x)- (A->x-p->x)*(B->y-A->y) ) /
- ( (B->x-A->x)*(B->x-A->x) +(B->y-A->y)*(B->y-A->y) );
+ s = ((A->y - p->y) * (B->x - A->x) - (A->x - p->x) * (B->y - A->y)) /
+ ((B->x - A->x) * (B->x - A->x) + (B->y - A->y) * (B->y - A->y));
- return FP_ABS(s) * sqrt(
- (B->x-A->x)*(B->x-A->x) + (B->y-A->y)*(B->y-A->y)
- );
+ return FP_ABS(s) * sqrt((B->x - A->x) * (B->x - A->x) + (B->y - A->y) * (B->y - A->y));
}
/* return distance squared, useful to avoid sqrt calculations */
@@ -2370,17 +2455,19 @@
double
distance2d_sqr_pt_seg(const POINT2D *p, const POINT2D *A, const POINT2D *B)
{
- double r,s;
+ double r, s;
- if ( ( A->x == B->x) && (A->y == B->y) )
- return distance2d_sqr_pt_pt(p,A);
+ if ((A->x == B->x) && (A->y == B->y))
+ return distance2d_sqr_pt_pt(p, A);
- r = ( (p->x-A->x) * (B->x-A->x) + (p->y-A->y) * (B->y-A->y) )/( (B->x-A->x)*(B->x-A->x) +(B->y-A->y)*(B->y-A->y) );
+ r = ((p->x - A->x) * (B->x - A->x) + (p->y - A->y) * (B->y - A->y)) /
+ ((B->x - A->x) * (B->x - A->x) + (B->y - A->y) * (B->y - A->y));
- if (r<0) return distance2d_sqr_pt_pt(p,A);
- if (r>1) return distance2d_sqr_pt_pt(p,B);
+ if (r < 0)
+ return distance2d_sqr_pt_pt(p, A);
+ if (r > 1)
+ return distance2d_sqr_pt_pt(p, B);
-
/*
* (2)
* (Ay-Cy)(Bx-Ax)-(Ax-Cx)(By-Ay)
@@ -2391,14 +2478,12 @@
*
*/
- s = ( (A->y-p->y)*(B->x-A->x)- (A->x-p->x)*(B->y-A->y) ) /
- ( (B->x-A->x)*(B->x-A->x) +(B->y-A->y)*(B->y-A->y) );
+ s = ((A->y - p->y) * (B->x - A->x) - (A->x - p->x) * (B->y - A->y)) /
+ ((B->x - A->x) * (B->x - A->x) + (B->y - A->y) * (B->y - A->y));
- return s * s * ( (B->x-A->x)*(B->x-A->x) + (B->y-A->y)*(B->y-A->y) );
+ return s * s * ((B->x - A->x) * (B->x - A->x) + (B->y - A->y) * (B->y - A->y));
}
-
-
/**
* Compute the azimuth of segment AB in radians.
* Return 0 on exception (same point), 1 otherwise.
@@ -2406,50 +2491,8 @@
int
azimuth_pt_pt(const POINT2D *A, const POINT2D *B, double *d)
{
- if ( A->x == B->x )
- {
- if ( A->y < B->y ) *d=0.0;
- else if ( A->y > B->y ) *d=M_PI;
- else return 0;
- return 1;
- }
-
- if ( A->y == B->y )
- {
- if ( A->x < B->x ) *d=M_PI/2;
- else if ( A->x > B->x ) *d=M_PI+(M_PI/2);
- else return 0;
- return 1;
- }
-
- if ( A->x < B->x )
- {
- if ( A->y < B->y )
- {
- *d=atan(fabs(A->x - B->x) / fabs(A->y - B->y) );
- }
- else /* ( A->y > B->y ) - equality case handled above */
- {
- *d=atan(fabs(A->y - B->y) / fabs(A->x - B->x) )
- + (M_PI/2);
- }
- }
-
- else /* ( A->x > B->x ) - equality case handled above */
- {
- if ( A->y > B->y )
- {
- *d=atan(fabs(A->x - B->x) / fabs(A->y - B->y) )
- + M_PI;
- }
- else /* ( A->y < B->y ) - equality case handled above */
- {
- *d=atan(fabs(A->y - B->y) / fabs(A->x - B->x) )
- + (M_PI+(M_PI/2));
- }
- }
-
- return 1;
+ if (A->x == B->x && A->y == B->y)
+ return LW_FALSE;
+ *d = fmod(2 * M_PI + M_PI / 2 - atan2(B->y - A->y, B->x - A->x), 2 * M_PI);
+ return LW_TRUE;
}
-
-
Modified: trunk/liblwgeom/measures.h
===================================================================
--- trunk/liblwgeom/measures.h 2019-08-23 06:37:52 UTC (rev 17772)
+++ trunk/liblwgeom/measures.h 2019-08-23 09:56:37 UTC (rev 17773)
@@ -83,13 +83,19 @@
int lw_dist2d_ptarray_poly(POINTARRAY *pa, LWPOLY *poly, DISTPTS *dl);
int lw_dist2d_point_point(LWPOINT *point1, LWPOINT *point2, DISTPTS *dl);
int lw_dist2d_point_line(LWPOINT *point, LWLINE *line, DISTPTS *dl);
+int lw_dist2d_point_tri(LWPOINT *point, LWTRIANGLE *tri, DISTPTS *dl);
int lw_dist2d_point_circstring(LWPOINT *point, LWCIRCSTRING *circ, DISTPTS *dl);
int lw_dist2d_point_poly(LWPOINT *point, LWPOLY *poly, DISTPTS *dl);
int lw_dist2d_point_curvepoly(LWPOINT *point, LWCURVEPOLY *poly, DISTPTS *dl);
int lw_dist2d_line_line(LWLINE *line1, LWLINE *line2, DISTPTS *dl);
+int lw_dist2d_line_tri(LWLINE *line, LWTRIANGLE *tri, DISTPTS *dl);
int lw_dist2d_line_circstring(LWLINE *line1, LWCIRCSTRING *line2, DISTPTS *dl);
int lw_dist2d_line_poly(LWLINE *line, LWPOLY *poly, DISTPTS *dl);
int lw_dist2d_line_curvepoly(LWLINE *line, LWCURVEPOLY *poly, DISTPTS *dl);
+int lw_dist2d_tri_tri(LWTRIANGLE *tri1, LWTRIANGLE *tri2, DISTPTS *dl);
+int lw_dist2d_tri_circstring(LWTRIANGLE *tri, LWCIRCSTRING *line, DISTPTS *dl);
+int lw_dist2d_tri_poly(LWTRIANGLE *tri, LWPOLY *poly, DISTPTS *dl);
+int lw_dist2d_tri_curvepoly(LWTRIANGLE *tri, LWCURVEPOLY *poly, DISTPTS *dl);
int lw_dist2d_circstring_circstring(LWCIRCSTRING *line1, LWCIRCSTRING *line2, DISTPTS *dl);
int lw_dist2d_circstring_poly(LWCIRCSTRING *circ, LWPOLY *poly, DISTPTS *dl);
int lw_dist2d_circstring_curvepoly(LWCIRCSTRING *circ, LWCURVEPOLY *poly, DISTPTS *dl);
@@ -127,4 +133,4 @@
LWGEOM *lw_dist2d_distancepoint(const LWGEOM *lw1, const LWGEOM *lw2, int32_t srid, int mode);
LWGEOM *lw_dist2d_distanceline(const LWGEOM *lw1, const LWGEOM *lw2, int32_t srid, int mode);
-#endif /* !defined _MEASURES_H */
\ No newline at end of file
+#endif /* !defined _MEASURES_H */
Modified: trunk/regress/core/measures.sql
===================================================================
--- trunk/regress/core/measures.sql 2019-08-23 06:37:52 UTC (rev 17772)
+++ trunk/regress/core/measures.sql 2019-08-23 09:56:37 UTC (rev 17773)
@@ -282,3 +282,11 @@
-- cast to text and back as a way of getting rid of solid flag
select '#4278.3', ST_3DIntersects('BOX3D(0 0 0, 1 1 1)'::box3d::geometry::text::geometry, 'BOX3D(-1 -1 -1, 2 2 2)'::box3d::geometry::text::geometry);
select '#4278.4', ST_3DDistance('BOX3D(0 0 0, 1 1 1)'::box3d::geometry::text::geometry, 'BOX3D(-1 -1 -1, 2 2 2)'::box3d::geometry::text::geometry);
+
+
+SELECT '#4328.1', ST_Intersects('TIN(((0 0,1 0,0 1,0 0)))'::geometry, 'POINT(.1 .1)'::geometry), ST_3DIntersects('TIN(((0 0,1 0,0 1,0 0)))'::geometry, 'POINT(.1 .1)'::geometry);
+SELECT '#4328.2', ST_Intersects('TIN(((0 0,1 0,0 1,0 0)))'::geometry, 'LINESTRING(.1 .1, .2 .2)'::geometry), ST_3DIntersects('TIN(((0 0,1 0,0 1,0 0)))'::geometry, 'LINESTRING(.1 .1, .2 .2)'::geometry);
+SELECT '#4328.3', ST_Intersects('TIN(((0 0,1 0,0 1,0 0)))'::geometry, 'TRIANGLE((.1 .1, .2 .2, .2 .1, .1 .1))'::geometry), ST_3DIntersects('TIN(((0 0,1 0,0 1,0 0)))'::geometry, 'TRIANGLE((.1 .1, .2 .2, .2 .1, .1 .1))'::geometry);
+SELECT '#4328.4', ST_Intersects('TIN(((0 0,1 0,0 1,0 0)))'::geometry, 'POLYGON((.1 .1, .2 .2, .2 .1, .1 .1))'::geometry), ST_3DIntersects('TIN(((0 0,1 0,0 1,0 0)))'::geometry, 'POLYGON((.1 .1, .2 .2, .2 .1, .1 .1))'::geometry);
+SELECT '#4328.5', ST_Intersects('TIN(((0 0,3 0,0 3,0 0)))'::geometry, 'CIRCULARSTRING(1.1 1.1, 1.2 1.2, 1.2 1.1)'::geometry), ST_3DIntersects('TIN(((0 0,3 0,0 3,0 0)))'::geometry, 'CIRCULARSTRING(1.1 1.1, 1.2 1.2, 1.2 1.1)'::geometry);
+SELECT '#4328.6', ST_Intersects('TIN(((0 0,3 0,0 3,0 0)))'::geometry, 'CURVEPOLYGON(CIRCULARSTRING(1.1 1.1, 1.2 1.2, 1.2 1.1, 1.2 1.2, 1.1 1.1))'::geometry), ST_3DIntersects('TIN(((0 0,3 0,0 3,0 0)))'::geometry, 'CURVEPOLYGON(CIRCULARSTRING(1.1 1.1, 1.2 1.2, 1.2 1.1, 1.2 1.2, 1.1 1.1))'::geometry);
Modified: trunk/regress/core/measures_expected
===================================================================
--- trunk/regress/core/measures_expected 2019-08-23 06:37:52 UTC (rev 17772)
+++ trunk/regress/core/measures_expected 2019-08-23 09:56:37 UTC (rev 17773)
@@ -67,3 +67,15 @@
#4278.2|0
#4278.3|f
#4278.4|1
+NOTICE: One or both of the geometries is missing z-value. The unknown z-value will be regarded as "any value"
+#4328.1|t|t
+NOTICE: One or both of the geometries is missing z-value. The unknown z-value will be regarded as "any value"
+#4328.2|t|t
+NOTICE: One or both of the geometries is missing z-value. The unknown z-value will be regarded as "any value"
+#4328.3|t|t
+NOTICE: One or both of the geometries is missing z-value. The unknown z-value will be regarded as "any value"
+#4328.4|t|t
+NOTICE: One or both of the geometries is missing z-value. The unknown z-value will be regarded as "any value"
+#4328.5|t|t
+NOTICE: One or both of the geometries is missing z-value. The unknown z-value will be regarded as "any value"
+#4328.6|t|t
More information about the postgis-tickets
mailing list