[GRASS-SVN] r39080 - grass/branches/develbranch_6/vector/v.label.sa

svn_grass at osgeo.org svn_grass at osgeo.org
Mon Sep 7 17:21:25 EDT 2009


Author: wolf
Date: 2009-09-07 17:21:25 -0400 (Mon, 07 Sep 2009)
New Revision: 39080

Modified:
   grass/branches/develbranch_6/vector/v.label.sa/labels.c
Log:
Speedups by changing to a fast 2d segments_intersect function instead of the slower Vect_segment_intersect. Also eliminated extra calls to Vect_append_point in box_trans_rot, since we always have exactly 5 points.

Modified: grass/branches/develbranch_6/vector/v.label.sa/labels.c
===================================================================
--- grass/branches/develbranch_6/vector/v.label.sa/labels.c	2009-09-07 20:05:47 UTC (rev 39079)
+++ grass/branches/develbranch_6/vector/v.label.sa/labels.c	2009-09-07 21:21:25 UTC (rev 39080)
@@ -18,6 +18,9 @@
 			       label_point_t * p);
 static int box_overlap(BOUND_BOX * a, BOUND_BOX * b);
 static int box_overlap2(struct line_pnts *a, struct line_pnts *b);
+static int segments_intersect(double x1, double y1, double x2, double y2,
+			      double x3, double y3, double x4, double y4,
+			      double *x, double *y);
 
 /**
  * The font size in map units. A global variable because I'm lazy :P
@@ -375,10 +378,10 @@
 
     /* generate candidate location for each label based on feture type
      * see chapter 5 of MERL-TR-96-04 */
-	G_debug(3, "n_lables=%d", n_labels);
+	G_debug(3, "n_labels=%d", n_labels);
     fprintf(stderr, "Generating label candidates: ...");
     for (i = 0; i < n_labels; i++) {
-	G_percent(i, n_labels - 1, 1);
+	G_percent(i, n_labels, 1);
 	switch (labels[i].type) {
 	case GV_POINT:
 	    G_debug(3, "Point (%d): %s", i, labels[i].text);
@@ -397,6 +400,7 @@
 	    break;
 	}
     }
+    G_percent(n_labels, n_labels, 1);
     Vect_close(&Map);
     return;
 }
@@ -840,6 +844,8 @@
     return Points;
 }
 
+#define SET_POINT(P, i, X, Y) (P)->x[(i)]=(X); (P)->y[(i)]=(Y); (P)->z[(i)]=0
+
 struct line_pnts *box_trans_rot(BOUND_BOX * bb, label_point_t * p,
 				double angle)
 {
@@ -847,27 +853,35 @@
     double x0, y0, x1, y1, x2, y2;
 
     Points = Vect_new_line_struct();
+    Points->x = G_calloc(5, sizeof(double));
+    Points->y = G_calloc(5, sizeof(double));
+    Points->z = G_calloc(5, sizeof(double));
+    Points->n_points=5;
+    Points->alloc_points=5;
+    if((Points->x == NULL) || (Points->y == NULL)) {
+	return NULL;
+    }
 
     /* Lower Left, no rotation needed */
     x0 = p->x + bb->W;
     y0 = p->y + bb->S;
-    Vect_append_point(Points, x0, y0, 0);
+    SET_POINT(Points, 0, x0, y0);
     /* Lower Right */
     x1 = (bb->E - bb->W) * cos(angle);
     y1 = (bb->E - bb->W) * sin(angle);
-    Vect_append_point(Points, x0 + x1, y0 + y1, 0);
+    SET_POINT(Points, 1, x0 + x1, y0 + y1);
 
     /* Upper Right */
     x2 = (bb->N - bb->S) * sin(angle);
     y2 = (bb->N - bb->S) * cos(angle);
     /* First translate to LR, and then translate like UL */
-    Vect_append_point(Points, x0 + x1 - x2, y0 + y1 + y2, 0);
+    SET_POINT(Points, 2, x0 + x1 - x2, y0 + y1 + y2);
 
     /* Upper Left */
-    Vect_append_point(Points, x0 - x2, y0 + y2, 0);
+    SET_POINT(Points, 3, x0 - x2, y0 + y2);
 
     /* close polygon */
-    Vect_append_point(Points, x0, y0, 0);
+    SET_POINT(Points, 4, x0, y0);
 
     return Points;
 }
@@ -952,55 +966,49 @@
     for (i = 1; i < candidate->swathline->n_points; i++) {
 	int r;
 	double b, h;
-	double px1, py1, pz1, px2, py2, pz2;
+	double px, py, pz;
 
-	r = Vect_segment_intersection(x1, y1, 0, x2, y2, 0,
-				      candidate->swathline->x[i - 1],
-				      candidate->swathline->y[i - 1],
-				      0,
-				      candidate->swathline->x[i],
-				      candidate->swathline->y[i],
-				      0,
-				      &px1, &py1, &pz1, &px2, &py2, &pz2, 0);
+	r = segments_intersect(x1, y1, x2, y2,
+				candidate->swathline->x[i - 1],
+				candidate->swathline->y[i - 1],
+				candidate->swathline->x[i],
+				candidate->swathline->y[i],
+				&px, &py);
 	/* Now calculate the area between the swath and the line */
 	switch (r) {
 	case 0:		/* no intersection */
 	    dig_distance2_point_to_line(candidate->swathline->x[i],
 					candidate->swathline->y[i], 0,
 					x1, y1, 0, x2, y2, 0, 0,
-					&px1, &py1, &pz1, &h, NULL);
+					&px, &py, &pz, &h, NULL);
 	    h = (sqrt(pow(x1 - candidate->swathline->x[i - 1], 2.0) +
 		      pow(y1 - candidate->swathline->y[i - 1], 2.0)) +
 		 h) / 2.0;
-	    b = sqrt(pow(px1 - x1, 2) + pow(py1 - y1, 2));
+	    b = sqrt(pow(px - x1, 2) + pow(py - y1, 2));
 	    flatness += b * h;
-	    x1 = px1;
-	    y1 = py1;
+	    x1 = px;
+	    y1 = py;
 	    break;
 	case 1:
 	    h = sqrt(pow(x1 - candidate->swathline->x[i - 1], 2.0) +
 		     pow(y1 - candidate->swathline->y[i - 1], 2.0));
-	    b = sqrt(pow(px1 - x1, 2) + pow(py1 - y1, 2));
+	    b = sqrt(pow(px - x1, 2) + pow(py - y1, 2));
 	    flatness += b * h * 0.5;	/* the first triangle */
-	    x1 = px1;
-	    y1 = py1;
+	    x1 = px;
+	    y1 = py;
 	    dig_distance2_point_to_line(candidate->swathline->x[i],
 					candidate->swathline->y[i], 0,
 					x1, y1, 0, x2, y2, 0, 0,
-					&px1, &py1, &pz1, &h, NULL);
-	    b = sqrt(pow(px1 - x1, 2) + pow(py1 - y1, 2));
+					&px, &py, &pz, &h, NULL);
+	    b = sqrt(pow(px - x1, 2) + pow(py - y1, 2));
 	    flatness += b * h * 0.5;	/* the second triangle */
-	    x1 = px1;
-	    y1 = py1;
+	    x1 = px;
+	    y1 = py;
 	    break;
-	case 3:
-	    x1 = px2;
-	    y1 = py2;
+	case 2:
+	    x1 = px;
+	    y1 = py;
 	    break;
-	case 5:
-	    x1 = px2;
-	    y1 = py2;
-	    break;
 	default:
 	    G_fatal_error("Programming error!!\n");
 	    break;
@@ -1087,30 +1095,23 @@
 
 	    for (k = 1; k < trbb->n_points; k++) {
 		int r;
-		double x1, x2, y1, y2, z1, z2;
+		double x=0, y=0;
 
-		r = Vect_segment_intersection(trbb->x[k - 1], trbb->y[k - 1],
-					      0, trbb->x[k], trbb->y[k], 0,
-					      line->x[j - 1], line->y[j - 1],
-					      0, line->x[j], line->y[j], 0,
-					      &x1, &y1, &z1, &x2, &y2, &z2,
-					      0);
+		r = segments_intersect(trbb->x[k - 1], trbb->y[k - 1],
+				      trbb->x[k], trbb->y[k],
+				      line->x[j - 1], line->y[j - 1],
+				      line->x[j], line->y[j],
+				      &x, &y);
 		if (r > 0) {	/* intersection at one point */
 		    if (found == 0) {
 			found = 1;
-			v1.x = x1;
-			v1.y = y1;
+			v1.x = x;
+			v1.y = y;
 		    }
 		    else {
 			found++;
-			if (r > 1) {
-			    v2.x = x2;
-			    v2.y = y2;
-			}
-			else {
-			    v2.x = x1;
-			    v2.y = y1;
-			}
+			v2.x = x;
+			v2.y = y;
 		    }
 		}
 	    }
@@ -1260,7 +1261,7 @@
 		}
 	    }
 	}
-	G_debug(3, "i=%d n_lables=%d", i, n_labels);
+	G_debug(3, "i=%d n_labels=%d", i, n_labels);
 	G_percent(i, n_labels, 1);
     }
     G_percent(n_labels, n_labels, 1);
@@ -1304,18 +1305,77 @@
 	int j;
 
 	for (j = 0; j < (b->n_points - 1); j++) {
-	    double d[6];
+	    r += segments_intersect(a->x[i], a->y[i],
+				   a->x[i + 1], a->y[i + 1],
+				   b->x[j], b->y[j],
+				   b->x[j + 1], a->y[j + 1],
+				   NULL, NULL);
+	    if (r == 1) {
+		return 1;
+	    }
+	}
+    }
+    return 0;
+}
 
-	    r += Vect_segment_intersection(a->x[i], a->y[i], 0,
-					   a->x[i + 1], a->y[i + 1], 0,
-					   b->x[j], b->y[j], 0,
-					   b->x[j + 1], a->y[j + 1], 0,
-					   &d[0], &d[1], &d[2],
-					   &d[3], &d[4], &d[5], 0);
+/**
+ * This function checks if two line segments intersect. This function is based
+ * on the theory provided by Paul Burk.
+ * @sa http://ozviz.wasp.uwa.edu.au/~pbourke/geometry/lineline2d/
+ * @param x1 The X-coordinate of P1
+ * @param y1 The Y-coordinate of P1
+ * @param x2 The X-coordinate of P2
+ * @param y2 The Y-coordinate of P2
+ * @param x3 The X-coordinate of P3
+ * @param y3 The Y-coordinate of P3
+ * @param x4 The X-coordinate of P4
+ * @param y4 The Y-coordinate of P4
+ * @param x A pointer to a double which will be set to the X-coordinate of the
+            unique intersection point (Coincident segments: x4).
+ * @param y A pointer to a double which will be set to the Y-coordinate of the
+            unique intersection point (Coincident segments: y4).
+ * @return 0: No intersection. 1: Intersect. 2: The segments are coincident.
+ */
+static int segments_intersect(double x1, double y1, double x2, double y2,
+			      double x3, double y3, double x4, double y4,
+			      double *x, double *y)
+{
+    double u1numer = (x4-x3)*(y1-y3) - (y4-y3)*(x1-x3);
+    double u2numer = (x2-x1)*(y1-y3) - (y2-y1)*(x1-x3);
+    double denomin = (y4-y3)*(x2-x1) - (x4-x3)*(y2-y1);
+    double ua, ub;
+    int ret = 0;
+ /*
+    * If the denominator and numerator for the equations for ua and ub are 0 then the two lines are coincident.
+ */
+ if(((u1numer < GRASS_EPSILON) && (u1numer > -GRASS_EPSILON)) &&
+    ((u2numer < GRASS_EPSILON) && (u2numer > -GRASS_EPSILON)))
+    {
+	if((x != NULL) && (y != NULL)) {
+	    *x = x4;
+	    *y = y4;
 	}
+	ret = 2;
     }
-    if (r > 1)
-	return 1;
-    else
-	return 0;
+/*
+    * The denominators for the equations for ua and ub are the same.
+    * If the denominator for the equations for ua and ub is 0 then the two lines are parallel.
+ */
+    if((denomin < GRASS_EPSILON) && (denomin > -GRASS_EPSILON)) {
+	ret = 0; /* parallel but don't overlap */
+    }
+ /*
+    * The equations apply to lines, if the intersection of line segments is required then it is only necessary to test if ua and ub lie between 0 and 1. Whichever one lies within that range then the corresponding line segment contains the intersection point. If both lie within the range of 0 to 1 then the intersection point is within both line segments. 
+    */
+    ua = u1numer / denomin;
+    ub = u2numer / denomin;
+    if(((ua > 0.0) && (ua < 1.0)) && ((ub > 0.0) && (ub < 1.0))) {
+	ret = 1;
+	if((x != NULL) && (y != NULL)) {
+	    *x = x1 + ua * (x2 - x1);
+	    *y = y1 + ua * (y2 - y1);
+	}
+    }
+    return ret;
 }
+



More information about the grass-commit mailing list