CU_ASSERT_EQUAL(gh, rs);
}
-static void test_lwgeom_remote_repeated_points(void)
-{
- LWGEOM *g;
- char *ewkt;
-
- // g = lwgeom_from_wkt("MULTIPOINT(0 0, 10 0, 10 10, 10 10, 0 10, 0 10, 0 10, 0 0, 0 0, 0 0, 5 5, 0 0, 5 5)", LW_PARSER_CHECK_NONE);
- // lwgeom_remove_repeated_points_in_place(g, 1);
- // ewkt = lwgeom_to_ewkt(g);
- // CU_ASSERT_STRING_EQUAL(ewkt, "MULTIPOINT(0 0,10 0,10 10,0 10,5 5)");
- // printf("%s\n", ewkt);
- // lwgeom_free(g);
- // lwfree(ewkt);
-
- g = lwgeom_from_wkt("MULTIPOINT(0 0,10 0,10 10,10 10,0 10,0 10,0 10,0 0,0 0,0 0,5 5,5 5,5 8,8 8,8 8,8 8,8 5,8 5,5 5,5 5,5 5,5 5,5 5,50 50,50 50,50 50,50 60,50 60,50 60,60 60,60 50,60 50,50 50,55 55,55 58,58 58,58 55,58 55,55 55)", LW_PARSER_CHECK_NONE);
- lwgeom_remove_repeated_points_in_place(g, 1);
- ewkt = lwgeom_to_ewkt(g);
- CU_ASSERT_STRING_EQUAL(ewkt, "MULTIPOINT(0 0,10 0,10 10,0 10,5 5,5 8,8 8,8 5,50 50,50 60,60 60,60 50,55 55,55 58,58 58,58 55)");
- printf("%s\n", ewkt);
- lwgeom_free(g);
- lwfree(ewkt);
-}
-
static void test_lwgeom_simplify(void)
{
LWGEOM *l;
PG_ADD_TEST(suite,test_median_handles_3d_correctly);
PG_ADD_TEST(suite,test_median_robustness);
PG_ADD_TEST(suite,test_lwpoly_construct_circle);
- PG_ADD_TEST(suite,test_lwgeom_remote_repeated_points);
}
/*
* Repeated points
*/
+POINTARRAY *ptarray_remove_repeated_points_minpoints(const POINTARRAY *in, double tolerance, int minpoints);
POINTARRAY *ptarray_remove_repeated_points(const POINTARRAY *in, double tolerance);
+LWGEOM* lwmpoint_remove_repeated_points(const LWMPOINT *in, double tolerance);
LWGEOM* lwline_remove_repeated_points(const LWLINE *in, double tolerance);
-void ptarray_remove_repeated_points_in_place(POINTARRAY *pa, double tolerance, int min_points);
-void lwgeom_remove_repeated_points_in_place(LWGEOM *geom, double tolerance);
+LWGEOM* lwcollection_remove_repeated_points(const LWCOLLECTION *in, double tolerance);
+LWGEOM* lwpoly_remove_repeated_points(const LWPOLY *in, double tolerance);
/*
* Closure test
return outcol;
}
+LWGEOM*
+lwcollection_remove_repeated_points(const LWCOLLECTION *coll, double tolerance)
+{
+ uint32_t i;
+ LWGEOM **newgeoms;
+
+ newgeoms = lwalloc(sizeof(LWGEOM *)*coll->ngeoms);
+ for (i=0; i<coll->ngeoms; i++)
+ {
+ newgeoms[i] = lwgeom_remove_repeated_points(coll->geoms[i], tolerance);
+ }
+
+ return (LWGEOM*)lwcollection_construct(coll->type,
+ coll->srid, coll->bbox ? gbox_copy(coll->bbox) : NULL,
+ coll->ngeoms, newgeoms);
+}
+
LWCOLLECTION*
lwcollection_force_dims(const LWCOLLECTION *col, int hasz, int hasm)
extern LWGEOM* lwgeom_remove_repeated_points(const LWGEOM *in, double tolerance)
{
- LWGEOM *out = lwgeom_clone_deep(in);
- lwgeom_remove_repeated_points_in_place(out, tolerance);
- return out;
+ LWDEBUGF(4, "lwgeom_remove_repeated_points got type %s",
+ lwtype_name(in->type));
+
+ if(lwgeom_is_empty(in))
+ {
+ return lwgeom_clone_deep(in);
+ }
+
+ switch (in->type)
+ {
+ case MULTIPOINTTYPE:
+ return lwmpoint_remove_repeated_points((LWMPOINT*)in, tolerance);
+ break;
+ case LINETYPE:
+ return lwline_remove_repeated_points((LWLINE*)in, tolerance);
+
+ case MULTILINETYPE:
+ case COLLECTIONTYPE:
+ case MULTIPOLYGONTYPE:
+ case POLYHEDRALSURFACETYPE:
+ return lwcollection_remove_repeated_points((LWCOLLECTION *)in, tolerance);
+
+ case POLYGONTYPE:
+ return lwpoly_remove_repeated_points((LWPOLY *)in, tolerance);
+ break;
+
+ case POINTTYPE:
+ case TRIANGLETYPE:
+ case TINTYPE:
+ /* No point is repeated for a single point, or for Triangle or TIN */
+ return lwgeom_clone_deep(in);
+
+ case CIRCSTRINGTYPE:
+ case COMPOUNDTYPE:
+ case MULTICURVETYPE:
+ case CURVEPOLYTYPE:
+ case MULTISURFACETYPE:
+ /* Dunno how to handle these, will return untouched */
+ return lwgeom_clone_deep(in);
+
+ default:
+ lwnotice("%s: unsupported geometry type: %s",
+ __func__, lwtype_name(in->type));
+ return lwgeom_clone_deep(in);
+ break;
+ }
+ return 0;
}
void lwgeom_swap_ordinates(LWGEOM *in, LWORD o1, LWORD o2)
}
}
-
/**************************************************************/
-void
-lwgeom_remove_repeated_points_in_place(LWGEOM *geom, double tolerance)
-{
- switch (geom->type)
- {
- /* No-op! Cannot remote points */
- case POINTTYPE:
- return;
- case LINETYPE:
- {
- LWLINE *g = (LWLINE*)(geom);
- POINTARRAY *pa = g->points;
- ptarray_remove_repeated_points_in_place(pa, tolerance, 2);
- /* Invalid output */
- if (pa->npoints == 1 && pa->maxpoints > 1)
- {
- /* Use first point as last point */
- pa->npoints = 2;
- ptarray_copy_point(pa, 0, 1);
- }
- break;
- }
- case POLYGONTYPE:
- {
- int i, j = 0;
- LWPOLY *g = (LWPOLY*)(geom);
- for (i = 0; i < g->nrings; i++)
- {
- POINTARRAY *pa = g->rings[i];
- int minpoints = 4;
- /* Skip zero'ed out rings */
- if(!pa)
- continue;
- ptarray_remove_repeated_points_in_place(pa, tolerance, minpoints);
- /* Drop collapsed rings */
- if(pa->npoints < 4)
- {
- ptarray_free(pa);
- continue;
- }
- g->rings[j++] = pa;
- }
- /* Update ring count */
- g->nrings = j;
- break;
- }
- case MULTIPOINTTYPE:
- {
- static int out_stack_size = 32;
- double tolsq = tolerance*tolerance;
- int i, j, n = 0;
- LWMPOINT *mpt = (LWMPOINT *)(geom);
- LWPOINT **out;
- LWPOINT *out_stack[out_stack_size];
- int use_heap = (mpt->ngeoms > out_stack_size);
-
- /* No-op on empty */
- if (mpt->ngeoms == 0) return;
-
- /* We cannot write directly back to the multipoint */
- /* geoms array because we're reading out of it still */
- /* so we use a side array */
- if (use_heap)
- out = lwalloc(sizeof(LWMPOINT *) * mpt->ngeoms);
- else
- out = out_stack;
-
- /* Inefficient O(n^2) implementation */
- for (i = 0; i < mpt->ngeoms; i++)
- {
- int seen = 0;
- LWPOINT *p1 = mpt->geoms[i];
- const POINT2D *pt1 = getPoint2d_cp(p1->point, 0);
- for (j = 0; j < n; j++)
- {
- LWPOINT *p2 = out[j];
- const POINT2D *pt2 = getPoint2d_cp(p2->point, 0);
- if (distance2d_sqr_pt_pt(pt1, pt2) <= tolsq)
- {
- seen = 1;
- break;
- }
- }
- if (seen) continue;
- out[n++] = p1;
- }
-
- /* Copy remaining points back into the input */
- /* array */
- memcpy(mpt->geoms, out, sizeof(LWPOINT *) * n);
- mpt->ngeoms = n;
- if (use_heap) lwfree(out);
- return;
- }
-
- case CIRCSTRINGTYPE:
- /* Dunno how to handle these, will return untouched */
- return;
-
- /* Can process most multi* types as generic collection */
- case MULTILINETYPE:
- case MULTIPOLYGONTYPE:
- case COLLECTIONTYPE:
- /* Curve types we mostly ignore, but allow the linear */
- /* portions to be processed by recursing into them */
- case MULTICURVETYPE:
- case CURVEPOLYTYPE:
- case MULTISURFACETYPE:
- case COMPOUNDTYPE:
- {
- int i, j = 0;
- LWCOLLECTION *col = (LWCOLLECTION*)(geom);
- for (i = 0; i < col->ngeoms; i++)
- {
- LWGEOM *g = col->geoms[i];
- if (!g) continue;
- lwgeom_remove_repeated_points_in_place(g, tolerance);
- /* Drop zero'ed out geometries */
- if(lwgeom_is_empty(g))
- {
- lwgeom_free(g);
- continue;
- }
- col->geoms[j++] = g;
- }
- /* Update geometry count */
- col->ngeoms = j;
- break;
- }
- default:
- {
- lwerror("%s: unsupported geometry type: %s", __func__, lwtype_name(geom->type));
- break;
- }
- }
- return;
-}
-
-
-/**************************************************************/
-
void
lwgeom_simplify_in_place(LWGEOM *geom, double epsilon, int preserve_collapsed)
{
return;
}
+
+/**************************************************************/
+
LWGEOM* lwgeom_simplify(const LWGEOM *igeom, double dist, int preserve_collapsed)
{
LWGEOM *lwgeom_out = lwgeom_clone_deep(igeom);
return lwgeom_out;
}
-/**************************************************************/
-
double lwgeom_area(const LWGEOM *geom)
{
LWGEOM*
lwline_remove_repeated_points(const LWLINE *lwline, double tolerance)
{
- return lwgeom_remove_repeated_points((LWGEOM*)lwline, tolerance);
+ POINTARRAY* npts = ptarray_remove_repeated_points_minpoints(lwline->points, tolerance, 2);
+
+ LWDEBUGF(3, "%s: npts %p", __func__, npts);
+
+ return (LWGEOM*)lwline_construct(lwline->srid,
+ lwline->bbox ? gbox_copy(lwline->bbox) : 0,
+ npts);
}
int
lwfree(mpt);
}
+LWGEOM*
+lwmpoint_remove_repeated_points(const LWMPOINT *mpoint, double tolerance)
+{
+ uint32_t nnewgeoms;
+ uint32_t i, j;
+ LWGEOM **newgeoms;
+ LWGEOM *lwpt1, *lwpt2;
+
+ newgeoms = lwalloc(sizeof(LWGEOM *)*mpoint->ngeoms);
+ nnewgeoms = 0;
+ for (i=0; i<mpoint->ngeoms; ++i)
+ {
+ lwpt1 = (LWGEOM*)mpoint->geoms[i];
+ /* Brute force, may be optimized by building an index */
+ int seen=0;
+ for (j=0; j<nnewgeoms; ++j)
+ {
+ lwpt2 = (LWGEOM*)newgeoms[j];
+ if ( lwgeom_mindistance2d(lwpt1, lwpt2) <= tolerance )
+ {
+ seen=1;
+ break;
+ }
+ }
+ if ( seen ) continue;
+ newgeoms[nnewgeoms++] = lwgeom_clone_deep(lwpt1);
+ }
+
+ return (LWGEOM*)lwcollection_construct(mpoint->type,
+ mpoint->srid,
+ mpoint->bbox ? gbox_copy(mpoint->bbox) : NULL,
+ nnewgeoms, newgeoms);
+
+}
LWMPOINT*
lwmpoint_from_lwgeom(const LWGEOM *g)
return ret;
}
+LWGEOM*
+lwpoly_remove_repeated_points(const LWPOLY *poly, double tolerance)
+{
+ uint32_t i;
+ POINTARRAY **newrings;
+
+ newrings = lwalloc(sizeof(POINTARRAY *)*poly->nrings);
+ for (i=0; i<poly->nrings; i++)
+ {
+ newrings[i] = ptarray_remove_repeated_points_minpoints(poly->rings[i], tolerance, 4);
+ }
+
+ return (LWGEOM*)lwpoly_construct(poly->srid,
+ poly->bbox ? gbox_copy(poly->bbox) : NULL,
+ poly->nrings, newrings);
+
+}
+
+
LWPOLY*
lwpoly_force_dims(const LWPOLY *poly, int hasz, int hasm)
{
*
* Always returns a newly allocated object.
*/
-static POINTARRAY *
+POINTARRAY *
ptarray_remove_repeated_points_minpoints(const POINTARRAY *in, double tolerance, int minpoints)
{
- POINTARRAY *out = ptarray_clone_deep(in);
- ptarray_remove_repeated_points_in_place(out, tolerance, minpoints);
- return out;
-}
+ POINTARRAY *out;
+ size_t ptsize = ptarray_point_size(in);
+ int has_z = FLAGS_GET_Z(in->flags);
+ int has_m = FLAGS_GET_M(in->flags);
+ int ipn = 1, opn = 1;
+ const POINT2D *prev_point, *this_point;
+ uint8_t *p1, *p2;
+ double tolsq = tolerance * tolerance;
-POINTARRAY *
-ptarray_remove_repeated_points(const POINTARRAY *in, double tolerance)
-{
- return ptarray_remove_repeated_points_minpoints(in, tolerance, 2);
-}
+ LWDEBUGF(3, "%s called", __func__);
+ /* Single or zero point arrays can't have duplicates */
+ if ( in->npoints < 3 ) return ptarray_clone_deep(in);
-void
-ptarray_remove_repeated_points_in_place(POINTARRAY *pa, double tolerance, int min_points)
-{
- int i;
- double tolsq = tolerance * tolerance;
- const POINT2D *last = NULL;
- const POINT2D *pt;
- int n_points = pa->npoints;
- int n_points_out = 0;
- int pt_size = ptarray_point_size(pa);
- double dsq;
+ /* Condition minpoints */
+ if ( minpoints < 1 ) minpoints = 1;
- /* No-op on short inputs */
- if ( n_points <= 2 ) return;
+ /* Allocate enough output space for all points */
+ out = ptarray_construct(has_z, has_m, in->npoints);
- for (i = 0; i < n_points; i++)
- {
- int last_point = (i == n_points-1);
+ /* Keep the first point */
+ p1 = getPoint_internal(out, 0);
+ p2 = getPoint_internal(in, 0);
+ memcpy(p1, p2, ptsize);
- /* Look straight into the abyss */
- pt = getPoint2d_cp(pa, i);
+ /* Now fill up the actual points */
+ prev_point = getPoint2d_cp(in, 0);
+ LWDEBUGF(3, " first point copied, out points: %d", opn);
+ for ( ipn = 1; ipn < in->npoints; ipn++ )
+ {
+ this_point = getPoint2d_cp(in, ipn);
- /* Preserve first point always */
- if (last)
+ /*
+ * If number of points left > number of points we need
+ * then it's still OK to drop dupes
+ */
+ if ( in->npoints - ipn > minpoints - opn )
{
- /* Don't drop points if we are running short of points */
- if (n_points - i > min_points - n_points_out)
+ if (tolerance > 0.0)
{
- if (tolerance > 0.0)
- {
- /* Only drop points that are within our tolerance */
- dsq = distance2d_sqr_pt_pt(last, pt);
- /* Allow any point but the last one to be dropped */
- if (!last_point && dsq <= tolsq)
- {
- continue;
- }
- }
- else
- {
- /* At tolerance zero, only skip exact dupes */
- if (memcmp((char*)pt, (char*)last, pt_size) == 0)
- continue;
- }
+ /* within the removal tolerance? */
+ double dsq = distance2d_sqr_pt_pt(prev_point, this_point);
+ if (dsq <= tolsq)
+ /* then skip it */
+ continue;
+ }
+ else
+ {
+ /* exact duplicate? */
+ p1 = getPoint_internal(in, ipn-1);
+ p2 = getPoint_internal(in, ipn);
+ if (memcmp(p1, p2, ptsize) == 0)
+ /* then skip it */
+ continue;
}
}
+ /*
+ * The point is different (see above) from
+ * the previous point, so add it to output
+ */
+ p1 = getPoint_internal(out, opn++);
+ p2 = getPoint_internal(in, ipn);
+ memcpy(p1, p2, ptsize);
- /* Got to last point, and it's not very different from */
- /* the point that preceded it. We want to keep the last */
- /* point, not the second-to-last one, so we pull our write */
- /* index back one value */
- if (last_point && tolerance > 0.0 && dsq <= tolsq)
- {
- n_points_out--;
- }
-
- /* Compact all remaining values to front of array */
- ptarray_copy_point(pa, i, n_points_out++);
- last = pt;
+ prev_point = this_point;
}
- /* Adjust array length */
- pa->npoints = n_points_out;
- return;
+
+ /* Keep the last point */
+ p1 = getPoint_internal(out, opn-1); /* Last output point */
+ p2 = getPoint_internal(in, ipn-1); /* Last input point */
+ if ( memcmp(p1, p2, ptsize) != 0 )
+ memcpy(p1, p2, ptsize);
+
+ LWDEBUGF(3, " in:%d out:%d", out->npoints, opn);
+ out->npoints = opn;
+
+ return out;
}
+POINTARRAY *
+ptarray_remove_repeated_points(const POINTARRAY *in, double tolerance)
+{
+ return ptarray_remove_repeated_points_minpoints(in, tolerance, 2);
+}
/************************************************************************/
projects / postgis / commitdiff