self-hosted/postgis
self-hosted/postgis
1
0
Fork 0
mirror of https://git.osgeo.org/gitea/postgis/postgis synced 2024-10-23 08:32:42 +00:00
Code Issues Packages Projects Releases Wiki Activity

Enhanced speed _ST_DWithin(g,g,d) that returns as soon as g and g are within d of each other, rather than using distance naively. Change ST_DWithin to use enhanced op. (Issue 20)

Browse source 
git-svn-id: http://svn.osgeo.org/postgis/trunk@2796 b70326c6-7e19-0410-871a-916f4a2858ee
This commit is contained in:
Paul Ramsey 2008-05-28 23:03:11 +00:00
parent 698d86d72f
commit 55a1ecae41
7 changed files with 89 additions and 29 deletions
Show stats Download patch file Download diff file Expand all files Collapse all files

4
lwgeom/box2d.c
View file

@ -37,12 +37,12 @@ box2d_same(BOX2DFLOAT4 *box1, BOX2DFLOAT4 *box2)
(box1->xmin==box2->xmin) &&
(box1->ymax==box2->ymax) &&
(box1->ymin==box2->ymin));
#if 0 // changed to exact equality
#if 0
return(FPeq(box1->xmax, box2->xmax) &&
FPeq(box1->xmin, box2->xmin) &&
FPeq(box1->ymax, box2->ymax) &&
FPeq(box1->ymin, box2->ymin));
#endif // 0
#endif
}
BOX2DFLOAT4 *

1
lwgeom/liblwgeom.h
View file

@ -946,6 +946,7 @@ extern double distance2d_poly_poly(LWPOLY *poly1, LWPOLY *poly2);
extern double distance2d_line_poly(LWLINE *line, LWPOLY *poly);
extern int azimuth_pt_pt(POINT2D *p1, POINT2D *p2, double *ret);
extern double lwgeom_mindistance2d_recursive(uchar *lw1, uchar *lw2);
extern double lwgeom_mindistance2d_recursive_tolerance(uchar *lw1, uchar *lw2, double tolerance);
extern int lwgeom_pt_inside_circle(POINT2D *p, double cx, double cy, double rad);
extern int32 lwgeom_npoints(uchar *serialized);
extern char ptarray_isccw(const POINTARRAY *pa);

26
lwgeom/lwgeom_functions_analytic.c
View file

@ -794,7 +794,7 @@ Datum LWGEOM_snaptogrid(PG_FUNCTION_ARGS)
PG_LWGEOM *out_geom = NULL;
LWGEOM *out_lwgeom;
gridspec grid;
//BOX3D box3d;
/* BOX3D box3d; */
if ( PG_ARGISNULL(0) ) PG_RETURN_NULL();
datum = PG_GETARG_DATUM(0);
@ -859,7 +859,7 @@ Datum LWGEOM_snaptogrid(PG_FUNCTION_ARGS)
out_lwgeom->bbox = box3d_to_box2df(&box3d);
}
#endif // 0
#endif /* 0 */
#if VERBOSE
elog(NOTICE, "SnapToGrid made a %s", lwgeom_typename(TYPE_GETTYPE(out_lwgeom->type)));
@ -880,7 +880,7 @@ Datum LWGEOM_snaptogrid_pointoff(PG_FUNCTION_ARGS)
PG_LWGEOM *out_geom = NULL;
LWGEOM *out_lwgeom;
gridspec grid;
//BOX3D box3d;
/* BOX3D box3d; */
POINT4D offsetpoint;
if ( PG_ARGISNULL(0) ) PG_RETURN_NULL();
@ -965,7 +965,7 @@ Datum LWGEOM_snaptogrid_pointoff(PG_FUNCTION_ARGS)
out_lwgeom->bbox = box3d_to_box2df(&box3d);
}
#endif // 0
#endif /* 0 */
#if VERBOSE
elog(NOTICE, "SnapToGrid made a %s", lwgeom_typename(TYPE_GETTYPE(out_lwgeom->type)));
@ -1018,7 +1018,7 @@ Datum LWGEOM_line_substring(PG_FUNCTION_ARGS)
opa = ptarray_substring(ipa, from, to);
if ( opa->npoints == 1 ) // Point returned
if ( opa->npoints == 1 ) /* Point returned */
olwgeom = (LWGEOM *)lwpoint_construct(iline->SRID, NULL, opa);
else
olwgeom = (LWGEOM *)lwline_construct(iline->SRID, NULL, opa);
@ -1366,8 +1366,8 @@ int point_in_polygon_deprecated(LWPOLY *polygon, LWPOINT *point)
/* assume bbox short-circuit has already been attempted */
ring = polygon->rings[0];
//root = createTree(ring);
//if(point_in_ring(root, &pt) != 1)
/* root = createTree(ring); */
/* if(point_in_ring(root, &pt) != 1) */
if(point_in_ring_deprecated(polygon->rings[0], &pt) != 1)
{
#ifdef PGIS_DEBUG
@ -1379,8 +1379,8 @@ int point_in_polygon_deprecated(LWPOLY *polygon, LWPOINT *point)
for(i=1; i<polygon->nrings; i++)
{
ring = polygon->rings[i];
//root = createTree(ring);
//if(point_in_ring(root, &pt) != -1)
/* root = createTree(ring); */
/* if(point_in_ring(root, &pt) != -1) */
if(point_in_ring_deprecated(polygon->rings[i], &pt) != -1)
{
#ifdef PGIS_DEBUG
@ -1447,8 +1447,8 @@ int point_outside_polygon_deprecated(LWPOLY *polygon, LWPOINT *point)
/* assume bbox short-circuit has already been attempted */
ring = polygon->rings[0];
//root = createTree(ring);
//if(point_in_ring(root, &pt) == -1)
/* root = createTree(ring); */
/* if(point_in_ring(root, &pt) == -1) */
if(point_in_ring_deprecated(ring, &pt) == -1)
{
#ifdef PGIS_DEBUG
@ -1460,8 +1460,8 @@ int point_outside_polygon_deprecated(LWPOLY *polygon, LWPOINT *point)
for(i=1; i<polygon->nrings; i++)
{
ring = polygon->rings[i];
//root = createTree(ring);
//if(point_in_ring(root, &pt) == 1)
/* root = createTree(ring); */
/* if(point_in_ring(root, &pt) == 1) */
if(point_in_ring_deprecated(ring, &pt) == 1)
{
#ifdef PGIS_DEBUG

45
lwgeom/lwgeom_functions_basic.c
View file

@ -34,6 +34,7 @@ Datum LWGEOM_summary(PG_FUNCTION_ARGS);
Datum LWGEOM_npoints(PG_FUNCTION_ARGS);
Datum LWGEOM_nrings(PG_FUNCTION_ARGS);
Datum LWGEOM_area_polygon(PG_FUNCTION_ARGS);
Datum LWGEOM_dwithin(PG_FUNCTION_ARGS);
Datum postgis_uses_stats(PG_FUNCTION_ARGS);
Datum postgis_autocache_bbox(PG_FUNCTION_ARGS);
Datum postgis_scripts_released(PG_FUNCTION_ARGS);
@ -1629,6 +1630,50 @@ Datum LWGEOM_mindistance2d(PG_FUNCTION_ARGS)
PG_RETURN_FLOAT8(mindist);
}
/* Minimum 2d distance between objects in geom1 and geom2. */
PG_FUNCTION_INFO_V1(LWGEOM_dwithin);
Datum LWGEOM_dwithin(PG_FUNCTION_ARGS)
{
PG_LWGEOM *geom1;
PG_LWGEOM *geom2;
double mindist, tolerance;
#ifdef PROFILE
profstart(PROF_QRUN);
#endif
geom1 = (PG_LWGEOM *)PG_DETOAST_DATUM(PG_GETARG_DATUM(0));
geom2 = (PG_LWGEOM *)PG_DETOAST_DATUM(PG_GETARG_DATUM(1));
tolerance = PG_GETARG_FLOAT8(2);
if( tolerance < 0 ) {
elog(ERROR,"Tolerance cannot be less than zero\n");
PG_RETURN_NULL();
}
if (pglwgeom_getSRID(geom1) != pglwgeom_getSRID(geom2))
{
elog(ERROR,"Operation on two GEOMETRIES with different SRIDs\n");
PG_RETURN_NULL();
}
mindist = lwgeom_mindistance2d_recursive_tolerance(
SERIALIZED_FORM(geom1),
SERIALIZED_FORM(geom2),
tolerance
);
#ifdef PROFILE
profstop(PROF_QRUN);
profreport("dist",geom1, geom2, NULL);
#endif
PG_FREE_IF_COPY(geom1, 0);
PG_FREE_IF_COPY(geom2, 1);
PG_RETURN_BOOL(tolerance >= mindist);
}
/*
* Maximum 2d distance between linestrings.
* Returns NULL if given geoms are not linestrings.

4
lwgeom/lwgeom_svg.c
View file

@ -280,7 +280,7 @@ print_svg_path_abs(char *result, POINTARRAY *pa, int precision, int polygonRing)
{
getPoint2d_p(pa, u, &pt);
//close PATH with 'Z' for polygon rings if last point equals first point
/* close PATH with 'Z' for polygon rings if last point equals first point */
if(u > 0 && u == (pa->npoints - 1) && polygonRing)
{
POINT2D firstPoint;
@ -333,7 +333,7 @@ print_svg_path_rel(char *result, POINTARRAY *pa, int precision, int polygonRing)
if(u == (pa->npoints - 1) && polygonRing)
{
//close PATH with 'z' if last point equals first point
/* close PATH with 'z' if last point equals first point */
POINT2D firstPoint;
getPoint2d_p(pa, 0, &firstPoint);
if(pt.x == firstPoint.x && pt.y == firstPoint.y)

8
lwgeom/lwpostgis.sql.in.c
View file

@ -4430,10 +4430,16 @@ CREATEFUNCTION ST_Touches(geometry,geometry)
AS 'SELECT $1 && $2 AND _ST_Touches($1,$2)'
LANGUAGE 'SQL' _IMMUTABLE; -- WITH (iscachable);
-- Availability: 1.3.4
CREATEFUNCTION _ST_DWithin(geometry,geometry,float8)
RETURNS boolean
AS '@MODULE_FILENAME@', 'LWGEOM_dwithin'
LANGUAGE 'C' _IMMUTABLE_STRICT; -- WITH (isstrict,iscachable);
-- Availability: 1.2.2
CREATEFUNCTION ST_DWithin(geometry, geometry, float8)
RETURNS boolean
AS 'SELECT $1 && ST_Expand($2,$3) AND $2 && ST_Expand($1,$3) AND ST_Distance($1, $2) < $3'
AS 'SELECT $1 && ST_Expand($2,$3) AND $2 && ST_Expand($1,$3) AND _ST_DWithin($1, $2, $3)'
LANGUAGE 'SQL' _IMMUTABLE; -- WITH (iscachable);
-- Deprecation in 1.2.3

30
lwgeom/measures.c
View file

@ -47,7 +47,7 @@ int pt_in_ring_2d(POINT2D *p, POINTARRAY *ring)
#ifdef PGIS_DEBUG
lwnotice("pt_in_ring_2d called with point: %g %g", p->x, p->y);
//printPA(ring);
/* printPA(ring); */
#endif
/* loop through all edges of the polygon */
@ -91,10 +91,10 @@ double distance2d_pt_pt(POINT2D *p1, POINT2D *p2)
return sqrt ( hside*hside + vside*vside );
// the above is more readable
//return sqrt(
// (p2->x-p1->x) * (p2->x-p1->x) + (p2->y-p1->y) * (p2->y-p1->y)
// );
/* the above is more readable
return sqrt(
(p2->x-p1->x) * (p2->x-p1->x) + (p2->y-p1->y) * (p2->y-p1->y)
); */
}
/*distance2d from p to line A->B */
@ -510,7 +510,7 @@ double distance2d_poly_poly(LWPOLY *poly1, LWPOLY *poly2)
poly2->rings[j]);
if ( d <= 0 ) return 0.0;
// mindist is -1 when not yet set
/* mindist is -1 when not yet set */
if (mindist > -1) mindist = LW_MIN(mindist, d);
else mindist = d;
#ifdef PGIS_DEBUG
@ -676,6 +676,12 @@ lwnotice("in lwgeom_polygon_perimeter (%d rings)", poly->nrings);
double
lwgeom_mindistance2d_recursive(uchar *lw1, uchar *lw2)
{
return lwgeom_mindistance2d_recursive_tolerance( lw1, lw2, 0.0 );
}
double
lwgeom_mindistance2d_recursive_tolerance(uchar *lw1, uchar *lw2, double tolerance)
{
LWGEOM_INSPECTED *in1, *in2;
int i, j;
@ -688,13 +694,13 @@ lwgeom_mindistance2d_recursive(uchar *lw1, uchar *lw2)
{
uchar *g1 = lwgeom_getsubgeometry_inspected(in1, i);
int t1 = lwgeom_getType(g1[0]);
double dist=0;
double dist=tolerance;
/* it's a multitype... recurse */
if ( t1 >= 4 )
{
dist = lwgeom_mindistance2d_recursive(g1, lw2);
if ( dist == 0 ) return 0.0; /* can't be closer */
dist = lwgeom_mindistance2d_recursive_tolerance(g1, lw2, tolerance);
if ( dist <= tolerance ) return tolerance; /* can't be closer */
if ( mindist == -1 ) mindist = dist;
else mindist = LW_MIN(dist, mindist);
continue;
@ -779,7 +785,7 @@ lwgeom_mindistance2d_recursive(uchar *lw1, uchar *lw2)
}
else /* it's a multitype... recurse */
{
dist = lwgeom_mindistance2d_recursive(g1, g2);
dist = lwgeom_mindistance2d_recursive_tolerance(g1, g2, tolerance);
}
if (mindist == -1 ) mindist = dist;
@ -791,7 +797,7 @@ lwgeom_mindistance2d_recursive(uchar *lw1, uchar *lw2)
#endif
if (mindist <= 0.0) return 0.0; /* can't be closer */
if (mindist <= tolerance) return tolerance; /* can't be closer */
}
@ -802,6 +808,8 @@ lwgeom_mindistance2d_recursive(uchar *lw1, uchar *lw2)
return mindist;
}
int
lwgeom_pt_inside_circle(POINT2D *p, double cx, double cy, double rad)
{