[GRASS-SVN] r45838 - grass/branches/develbranch_6/lib/vector/Vlib
svn_grass at osgeo.org
svn_grass at osgeo.org
Mon Apr 4 04:31:30 EDT 2011
Author: mmetz
Date: 2011-04-04 01:31:30 -0700 (Mon, 04 Apr 2011)
New Revision: 45838
Modified:
grass/branches/develbranch_6/lib/vector/Vlib/buffer.c
grass/branches/develbranch_6/lib/vector/Vlib/buffer2.c
grass/branches/develbranch_6/lib/vector/Vlib/dgraph.c
Log:
improved cleaning, documentation
Modified: grass/branches/develbranch_6/lib/vector/Vlib/buffer.c
===================================================================
--- grass/branches/develbranch_6/lib/vector/Vlib/buffer.c 2011-04-03 12:54:46 UTC (rev 45837)
+++ grass/branches/develbranch_6/lib/vector/Vlib/buffer.c 2011-04-04 08:31:30 UTC (rev 45838)
@@ -15,9 +15,44 @@
\author Radim Blazek
\author Markus Metz
- \date 2001-2008
+ \date 2001-2011
*/
+/* BUGS
+ *
+ * buffering these shapes is unresolved for buffer distances which
+ * should create a closed loop in the inside
+ *
+ * ----------- ------------
+ * | \ /
+ * |
+ * |
+ * | / \
+ * ----------- ------------
+ *
+ *
+ * -----------
+ * | |
+ * |
+ * | |
+ * -----------
+ *
+ *
+ * for certain buffer distances, undefined behaviour for
+ *
+ * ---------
+ * | |
+ * |
+ * --------------
+ * |
+ * --------------
+ * |
+ * | |
+ * ---------
+ *
+ * MM April 2011
+ * */
+
#include <stdlib.h>
#include <math.h>
#include <grass/Vect.h>
@@ -25,8 +60,9 @@
#define LENGTH(DX, DY) ( sqrt( (DX*DX)+(DY*DY) ) )
#define PI M_PI
+#define D_MULT 0.99999999 /* distance multiplier for point_in_buf() */
-/* vector() calculates normalized vector form two points */
+/* vector() calculates normalized vector from two points */
static void vect(double x1, double y1, double x2, double y2, double *x,
double *y)
{
@@ -89,7 +125,7 @@
}
}
}
- G_debug(5, " no intersection");
+ G_debug(5, "find_cross() -> no intersection");
return 0;
}
@@ -101,14 +137,14 @@
** -1 found overlap
** 0 not found
*/
-static int find_cross2(struct line_pnts *Points, int s1, int s2, int s3,
+static int find_cross_reverse(struct line_pnts *Points, int s1, int s2, int s3,
int s4, int *s5, int *s6)
{
int i, j, jend, np, ret;
double *x, *y;
- G_debug(5,
- "find_cross2(): npoints = %d, s1 = %d, s2 = %d, s3 = %d, s4 = %d",
+ G_debug(5,
+ "find_cross_reverse(): npoints = %d, s1 = %d, s2 = %d, s3 = %d, s4 = %d",
Points->n_points, s1, s2, s3, s4);
x = Points->x;
@@ -139,11 +175,11 @@
}
}
}
- G_debug(5, " no intersection");
+ G_debug(5, "find_cross_reverse() -> no intersection");
return 0;
}
-/* find_cross3 find first crossing of segment s1 to s1 + 1 with all
+/* find_cross_from_start find first crossing of segment s1 to s1 + 1 with all
** segments from s1 + 1 to Points->n_points - 2
** ix is set to East crossing and iy to North crossing
** s2 is set to the intersecting segment
@@ -151,13 +187,14 @@
** returns: 1 found
** 0 not found
*/
-static int find_cross3(struct line_pnts *Points, int s1, int *s2, double *ix, double *iy)
+static int find_cross_from_start(struct line_pnts *Points, int s1, int *s2,
+ double *ix, double *iy)
{
int i, np, ret, found = 0;
double *x, *y, min_dist, new_dist, new_x, new_y, dx, dy;
- G_debug(5,
- "find_cross3(): npoints = %d, s1 = %d", Points->n_points, s1);
+ G_debug(5, "find_cross_from_start(): npoints = %d, s1 = %d",
+ Points->n_points, s1);
x = Points->x;
y = Points->y;
@@ -189,11 +226,12 @@
}
}
- G_debug(5, "find_cross3(): intersection %s", found ? "found" : "not found");
+ G_debug(5, "find_cross_from_start(): intersection %s",
+ found ? "found" : "not found");
return found;
}
-/* find_cross4 find first crossing of segment s1 to s1 + 1 with all
+/* find_cross_from_end find first crossing of segment s1 to s1 + 1 with all
** segments from 0 to s1 - 1
** ix is set to East crossing and iy to North crossing
** s2 is set to the intersecting segment
@@ -201,13 +239,14 @@
** returns: 1 found
** 0 not found
*/
-static int find_cross4(struct line_pnts *Points, int s1, int *s2, double *ix, double *iy)
+static int find_cross_from_end(struct line_pnts *Points, int s1, int *s2,
+ double *ix, double *iy)
{
int i, np, ret, found = 0;
double *x, *y, min_dist, new_dist, new_x, new_y, dx, dy;
- G_debug(5,
- "find_cross4(): npoints = %d, s1 = %d", Points->n_points, s1);
+ G_debug(5, "find_cross_from_end(): npoints = %d, s1 = %d",
+ Points->n_points, s1);
x = Points->x;
y = Points->y;
@@ -239,7 +278,8 @@
}
}
- G_debug(5, "find_cross4(): intersection %s", found ? "found" : "not found");
+ G_debug(5, "find_cross_from_end(): intersection %s",
+ found ? "found" : "not found");
return found;
}
@@ -254,7 +294,7 @@
double sd;
np = Points->n_points;
- d *= d;
+ /* d must be the squared distance */
for (i = 0; i < np - 1; i++) {
sd = dig_distance2_point_to_line(px, py, 0,
Points->x[i], Points->y[i], 0,
@@ -271,13 +311,14 @@
** - looking for loops and if loop doesn't contain any other loop
** and centroid of loop is in buffer removes this loop (repeated)
** - optionally removes all end points in buffer
+ ** - optionally closes output line if input line is looped
* parameters:
* Points - parallel line
* origPoints - original line
* d - offset
* rm_end - remove end points in buffer
- ** note1: on some lines (multiply selfcrossing; lines with end points
- ** in buffer of line other; some shapes of ends ) may create nosense
+ ** note1: on some lines (multiple selfcrossing; lines with end points
+ ** in buffer of other line; some shapes of ends ) may create nosense
** note2: this function is stupid and slow, somebody more clever
** than I am should write paralle_line + clean_parallel
** better; RB March 2000
@@ -292,6 +333,7 @@
int first = 0, current, last, lcount;
double *x, *y, px, py, ix, iy;
static struct line_pnts *sPoints = NULL;
+ double d2 = d * d * D_MULT;
G_debug(4, "clean_parallel(): npoints = %d, d = %f, rm_end = %d",
Points->n_points, d, rm_end);
@@ -303,46 +345,77 @@
if (np < 3)
return;
+
+ /* duplicate consecutive points cause problems for
+ * dig_test_for_intersection() and dig_find_intersection()
+ * -> output line must always be pruned */
if (rm_end) {
int inserted;
- /* remove points from start in buffer */
+ /* remove points from end in buffer */
j = inserted = 0;
- for (i = 0; i < Points->n_points - 1; i++) {
- px = (x[i] + x[i + 1]) / 2;
- py = (y[i] + y[i + 1]) / 2;
- if (point_in_buf(origPoints, x[i], y[i], d * 0.99999999)) {
+ for (i = Points->n_points - 1; i >= 1; i--) {
+ x = Points->x;
+ y = Points->y;
+ if (rm_end < 2) {
+ px = (x[i] + x[i - 1]) / 2;
+ py = (y[i] + y[i - 1]) / 2;
+ }
+ else {
+ /* this is for the last arc created by Vect_line_buffer() */
+ px = x[i - 1];
+ py = y[i - 1];
+ }
+ if (point_in_buf(origPoints, x[i], y[i], d2)) {
/* check for intersection of this segment */
- /* most not be the first or last point of this segment */
- current = i;
- if (find_cross3(Points, current, &sa, &ix, &iy) != 0) {
+ /* must not be the first or last point of this segment */
+ current = i - 1;
+ if (find_cross_from_end(Points, current, &sa, &ix, &iy) != 0) {
/* insert intersection point into this segment */
/* and adjust i and Points->n_points */
int k = 0;
-
- Vect_append_point(Points, Points->x[Points->n_points - 1],
- Points->y[Points->n_points - 1], 0);
- Vect_append_point(Points, Points->x[Points->n_points - 1],
- Points->y[Points->n_points - 1], 0);
-
- for (k = Points->n_points - 2; k > sa + 2; k--) {
- Points->x[k] = Points->x[k - 2];
- Points->y[k] = Points->y[k - 2];
+
+ if ((ix == x[sa] && iy == y[sa]) ||
+ (ix == x[sa + 1] && iy == y[sa + 1])) {
+
+ /* do not insert ix,iy into sa segment */
+ Vect_append_point(Points, Points->x[Points->n_points - 1],
+ Points->y[Points->n_points - 1], 0);
+
+ for (k = Points->n_points - 2; k > current + 1; k--) {
+ Points->x[k] = Points->x[k - 1];
+ Points->y[k] = Points->y[k - 1];
+ }
+ Points->x[current + 1] = ix;
+ Points->y[current + 1] = iy;
+ i += 2;
}
- Points->x[sa + 2] = ix;
- Points->y[sa + 2] = iy;
- for (k = sa + 1; k > current + 1; k--) {
- Points->x[k] = Points->x[k - 1];
- Points->y[k] = Points->y[k - 1];
+ else {
+ Vect_append_point(Points, Points->x[Points->n_points - 1],
+ Points->y[Points->n_points - 1], 0);
+ Vect_append_point(Points, Points->x[Points->n_points - 1],
+ Points->y[Points->n_points - 1], 0);
+
+ for (k = Points->n_points - 2; k > current + 2; k--) {
+ Points->x[k] = Points->x[k - 2];
+ Points->y[k] = Points->y[k - 2];
+ }
+ Points->x[current + 2] = ix;
+ Points->y[current + 2] = iy;
+ for (k = current + 1; k > sa + 1; k--) {
+ Points->x[k] = Points->x[k - 1];
+ Points->y[k] = Points->y[k - 1];
+ }
+ Points->x[sa + 1] = ix;
+ Points->y[sa + 1] = iy;
+ i += 3;
}
- Points->x[current + 1] = ix;
- Points->y[current + 1] = iy;
- i--;
+
inserted++;
}
- else if (point_in_buf(origPoints, px, py, d * 0.99999999)) {
+ else if (point_in_buf(origPoints, px, py, d2)) {
j++;
if (inserted)
close_loop = 0;
@@ -356,57 +429,64 @@
}
}
if (j > 0) {
- npn = 0;
- for (i = j; i < Points->n_points; i++) {
- x[npn] = x[i];
- y[npn] = y[i];
- npn++;
- }
- Points->n_points = npn;
+ Points->n_points -= j;
}
- /* remove points from end in buffer */
+
+ /* remove points from start in buffer */
j = inserted = 0;
- for (i = Points->n_points - 1; i >= 1; i--) {
- if (rm_end < 2) {
- px = (x[i] + x[i - 1]) / 2;
- py = (y[i] + y[i - 1]) / 2;
- }
- else {
- /* this is for the last arc created by Vect_line_buffer() */
- px = x[i - 1];
- py = y[i - 1];
- }
- if (point_in_buf(origPoints, x[i], y[i], d * 0.99999999)) {
+ for (i = 0; i < Points->n_points - 1; i++) {
+ x = Points->x;
+ y = Points->y;
+ px = (x[i] + x[i + 1]) / 2;
+ py = (y[i] + y[i + 1]) / 2;
+ if (point_in_buf(origPoints, x[i], y[i], d2)) {
/* check for intersection of this segment */
- /* most not be the first or last point of this segment */
- current = i - 1;
- if (find_cross4(Points, current, &sa, &ix, &iy) != 0) {
+ /* must not be the first or last point of this segment */
+ current = i;
+ if (find_cross_from_start(Points, current, &sa, &ix, &iy) != 0) {
/* insert intersection point into this segment */
/* and adjust i and Points->n_points */
int k = 0;
- Vect_append_point(Points, Points->x[Points->n_points - 1],
- Points->y[Points->n_points - 1], 0);
- Vect_append_point(Points, Points->x[Points->n_points - 1],
- Points->y[Points->n_points - 1], 0);
-
- for (k = Points->n_points - 2; k > current + 2; k--) {
- Points->x[k] = Points->x[k - 2];
- Points->y[k] = Points->y[k - 2];
+ if ((ix == x[sa] && iy == y[sa]) ||
+ (ix == x[sa + 1] && iy == y[sa + 1])) {
+
+ /* do not insert ix,iy into sa segment */
+ Vect_append_point(Points, Points->x[Points->n_points - 1],
+ Points->y[Points->n_points - 1], 0);
+
+ for (k = Points->n_points - 2; k > current + 1; k--) {
+ Points->x[k] = Points->x[k - 1];
+ Points->y[k] = Points->y[k - 1];
+ }
+ Points->x[current + 1] = ix;
+ Points->y[current + 1] = iy;
+ i--;
}
- Points->x[current + 2] = ix;
- Points->y[current + 2] = iy;
- for (k = current + 1; k > sa + 1; k--) {
- Points->x[k] = Points->x[k - 1];
- Points->y[k] = Points->y[k - 1];
+ else {
+ Vect_append_point(Points, Points->x[Points->n_points - 1],
+ Points->y[Points->n_points - 1], 0);
+ Vect_append_point(Points, Points->x[Points->n_points - 1],
+ Points->y[Points->n_points - 1], 0);
+
+ for (k = Points->n_points - 2; k > sa + 2; k--) {
+ Points->x[k] = Points->x[k - 2];
+ Points->y[k] = Points->y[k - 2];
+ }
+ Points->x[sa + 2] = ix;
+ Points->y[sa + 2] = iy;
+ for (k = sa + 1; k > current + 1; k--) {
+ Points->x[k] = Points->x[k - 1];
+ Points->y[k] = Points->y[k - 1];
+ }
+ Points->x[current + 1] = ix;
+ Points->y[current + 1] = iy;
+ i--;
}
- Points->x[sa + 1] = ix;
- Points->y[sa + 1] = iy;
- i += 3;
inserted++;
}
- else if (point_in_buf(origPoints, px, py, d * 0.99999999)) {
+ else if (point_in_buf(origPoints, px, py, d2)) {
j++;
if (inserted)
close_loop = 0;
@@ -419,8 +499,16 @@
break;
}
}
+ x = Points->x;
+ y = Points->y;
if (j > 0) {
- Points->n_points -= j;
+ npn = 0;
+ for (i = j; i < Points->n_points; i++) {
+ x[npn] = x[i];
+ y[npn] = y[i];
+ npn++;
+ }
+ Points->n_points = npn;
}
/* test for intersection of end segments */
@@ -460,7 +548,7 @@
npn = 1;
side = (int)(d / fabs(d));
-
+
/* remove loops */
while (np > 3 && first < np - 2) {
/* find first loop which doesn't contain any other loop */
@@ -486,7 +574,7 @@
} /* loop not found */
/* ensure sa is monotonically increasing, so npn doesn't reset low */
- /* disabled
+ /* disabled, otherwise nested loops are missed
if (sa > sa_max)
sa_max = sa;
if (sa < sa_max)
@@ -502,29 +590,47 @@
}
else {
double area_size;
+ int in_buffer = 0;
Vect_reset_line(sPoints);
dig_find_intersection(x[sa], y[sa], x[sa + 1], y[sa + 1], x[sb],
y[sb], x[sb + 1], y[sb + 1], &ix, &iy);
-
+
if (ix == x[0] && iy == y[0] && ix == x[np - 1] && iy == y[np - 1])
break;
+
Vect_append_point(sPoints, ix, iy, 0);
for (i = sa + 1; i < sb + 1; i++) { /* create loop polygon */
Vect_append_point(sPoints, x[i], y[i], 0);
}
+ /* close polygon */
+ Vect_append_point(sPoints, ix, iy, 0);
+ Vect_line_prune(sPoints);
- Vect_append_point(sPoints, ix, iy, 0);
+ /* is loop in buffer ? */
dig_find_area_poly(sPoints, &area_size);
+ in_buffer = area_size * side < 0;
+ if (!in_buffer && area_size) {
+ /* Vect_find_poly_centroid() may produce coords outside polygon */
+ Vect_find_poly_centroid(sPoints, &px, &py);
+ in_buffer = point_in_buf(origPoints, px, py, d2) != 0;
+ }
- /*Vect_find_poly_centroid(sPoints, &px, &py);
- if (point_in_buf(origPoints, px, py, d)) { */ /* is loop in buffer ? */
- if (area_size * side < 0) {
+ if (in_buffer) {
npn = sa + 1;
- x[npn] = ix;
- y[npn] = iy;
- j = sb + 1;
- npn++;
+ /* do not create zero-length segments */
+ if (ix != x[sa] && iy != y[sa]) {
+ x[npn] = ix;
+ y[npn] = iy;
+ npn++;
+ }
+ if (ix != x[sb + 1] && iy != y[sb + 1]) {
+ j = sb + 1;
+ }
+ else {
+ j = sb + 2;
+ }
+
/* lcount == 0 can not happen
if (lcount == 0) {
first = sa;
@@ -582,7 +688,8 @@
if (a > PI) {
G_debug(5, "search intersection");
/* there can be two intersections !!! */
- if (find_cross2(Points, 0, Points->n_points - 2, 1, Points->n_points - 1, &sa,
+ /* a > PI is probably handled by remove ends above */
+ if (find_cross_reverse(Points, 0, Points->n_points - 2, 1, Points->n_points - 1, &sa,
&sb) != 0) {
G_debug(5, "found intersection");
dig_find_intersection(x[sa], y[sa], x[sa + 1], y[sa + 1], x[sb],
@@ -748,6 +855,7 @@
int side, npoints;
static struct line_pnts *Points = NULL;
static struct line_pnts *PPoints = NULL;
+ double d2 = distance * distance * D_MULT;
distance = fabs(distance);
@@ -857,12 +965,12 @@
if (OutPoints->x[0] != OutPoints->x[OutPoints->n_points - 1] ||
OutPoints->y[0] != OutPoints->y[OutPoints->n_points - 1]) {
- if (point_in_buf(InPoints, OutPoints->x[0], OutPoints->y[0], distance * 0.99999999)) {
+ if (point_in_buf(InPoints, OutPoints->x[0], OutPoints->y[0], d2)) {
OutPoints->x[0] = OutPoints->x[OutPoints->n_points - 1];
OutPoints->y[0] = OutPoints->y[OutPoints->n_points - 1];
}
else if (point_in_buf(InPoints, OutPoints->x[OutPoints->n_points - 1],
- OutPoints->y[OutPoints->n_points - 1], distance * 0.99999999)) {
+ OutPoints->y[OutPoints->n_points - 1], d2)) {
OutPoints->x[OutPoints->n_points - 1] = OutPoints->x[0];
OutPoints->y[OutPoints->n_points - 1] = OutPoints->y[0];
}
Modified: grass/branches/develbranch_6/lib/vector/Vlib/buffer2.c
===================================================================
--- grass/branches/develbranch_6/lib/vector/Vlib/buffer2.c 2011-04-03 12:54:46 UTC (rev 45837)
+++ grass/branches/develbranch_6/lib/vector/Vlib/buffer2.c 2011-04-04 08:31:30 UTC (rev 45838)
@@ -485,7 +485,7 @@
/* close the output line */
Vect_append_point(nPoints, nPoints->x[0], nPoints->y[0], nPoints->z[0]);
- /* Vect_line_prune ( nPoints ); */
+ Vect_line_prune ( nPoints );
}
/*
@@ -533,7 +533,11 @@
eangle = atan2(vert->y - vert0->y, vert->x - vert0->x);
while (1) {
- Vect_append_point(nPoints, vert0->x, vert0->y, 0);
+ if (nPoints->n_points > 0 && (nPoints->x[nPoints->n_points - 1] != vert0->x ||
+ nPoints->y[nPoints->n_points - 1] != vert0->y))
+ Vect_append_point(nPoints, vert0->x, vert0->y, 0);
+ else
+ Vect_append_point(nPoints, vert0->x, vert0->y, 0);
G_debug(4, "ec: v0=%d, v=%d, eside=%d, edge->v1=%d, edge->v2=%d", v0,
v, eside, edge->v1, edge->v2);
G_debug(4, "ec: append point x=%.18f y=%.18f", vert0->x, vert0->y);
Modified: grass/branches/develbranch_6/lib/vector/Vlib/dgraph.c
===================================================================
--- grass/branches/develbranch_6/lib/vector/Vlib/dgraph.c 2011-04-03 12:54:46 UTC (rev 45837)
+++ grass/branches/develbranch_6/lib/vector/Vlib/dgraph.c 2011-04-04 08:31:30 UTC (rev 45838)
@@ -164,8 +164,8 @@
{
struct intersection_point *aa, *bb;
- aa = *((struct intersection_point **)a);
- bb = *((struct intersection_point **)b);
+ aa = *(struct intersection_point **)a;
+ bb = *(struct intersection_point **)b;
if (aa->x < bb->x)
return -1;
More information about the grass-commit
mailing list