[GRASS-SVN] r34821 - in grass/trunk/raster/r.watershed: front ram seg

svn_grass at osgeo.org svn_grass at osgeo.org
Wed Dec 10 10:13:06 EST 2008


Author: glynn
Date: 2008-12-10 10:13:06 -0500 (Wed, 10 Dec 2008)
New Revision: 34821

Modified:
   grass/trunk/raster/r.watershed/front/main.c
   grass/trunk/raster/r.watershed/front/r.watershed.html
   grass/trunk/raster/r.watershed/ram/do_cum.c
   grass/trunk/raster/r.watershed/ram/no_stream.c
   grass/trunk/raster/r.watershed/seg/close_maps.c
   grass/trunk/raster/r.watershed/seg/do_cum.c
   grass/trunk/raster/r.watershed/seg/no_stream.c
Log:
patch from Markus Metz


Modified: grass/trunk/raster/r.watershed/front/main.c
===================================================================
--- grass/trunk/raster/r.watershed/front/main.c	2008-12-10 05:35:27 UTC (rev 34820)
+++ grass/trunk/raster/r.watershed/front/main.c	2008-12-10 15:13:06 UTC (rev 34821)
@@ -368,16 +368,12 @@
 
     if (!flag_sfd->answer && opt16->answer) {
 	strcat(command, " conv=");
-	strcat(command, "\"");
 	strcat(command, opt16->answer);
-	strcat(command, "\"");
     }
 
     if (flag_seg->answer && opt17->answer) {
 	strcat(command, " mb=");
-	strcat(command, "\"");
 	strcat(command, opt17->answer);
-	strcat(command, "\"");
     }
 
     G_debug(1, "Mode: %s", flag_seg->answer ? "Segmented" : "All in RAM");

Modified: grass/trunk/raster/r.watershed/front/r.watershed.html
===================================================================
--- grass/trunk/raster/r.watershed/front/r.watershed.html	2008-12-10 05:35:27 UTC (rev 34820)
+++ grass/trunk/raster/r.watershed/front/r.watershed.html	2008-12-10 15:13:06 UTC (rev 34821)
@@ -224,9 +224,9 @@
 As a result, depressions and obstacles are overflown with a gracefull 
 flow convergence before the overflow. The convergence factor causes flow 
 accumulation to converge more strongly with higher values. The supported 
-range is 1 to 10, recommended is a convergence factor of 5. Although 
-basins and half-basins are calculated when using MFD, it is strongly 
-recommended to use SFD for basin and half-basin determination.
+range is 1 to 10, recommended is a convergence factor of 5. If many 
+small sliver basins are created with MFD, setting the convergence factor 
+to a higher value can reduce the amount of small sliver basins.
 <br>See example below with the North Carolina dataset for using MFD mode.
 
 <h4>In-memory mode and disk swap mode</h4>

Modified: grass/trunk/raster/r.watershed/ram/do_cum.c
===================================================================
--- grass/trunk/raster/r.watershed/ram/do_cum.c	2008-12-10 05:35:27 UTC (rev 34820)
+++ grass/trunk/raster/r.watershed/ram/do_cum.c	2008-12-10 15:13:06 UTC (rev 34821)
@@ -41,8 +41,9 @@
 		    valued = value - valued;
 	    }
 	    wat[SEG_INDEX(wat_seg, dr, dc)] = valued;
+	    valued = ABS(valued) + 0.5;
 	    is_swale = FLAG_GET(swale, r, c);
-	    if (is_swale || ABS(valued) >= threshold) {
+	    if (is_swale || ((int)valued) >= threshold) {
 		FLAG_SET(swale, dr, dc);
 	    }
 	    else {
@@ -83,7 +84,7 @@
 {
     int r, c, dr, dc;
     CELL is_swale;
-    DCELL value, valued, value_sfd;
+    DCELL value, valued;
     int killer, threshold, count;
 
     int mfd_cells, astar_not_set;
@@ -94,14 +95,15 @@
     double prop;
     int workedon;
 
-    G_message(_("SECTION 3: Accumulating Surface Flow with MFD."));
+    G_message(_("SECTION 3: Accumulating Surface Flow with MFD"));
+    G_debug(1, "MFD convergence factor set to %d.", c_fac);
 
     /* distances to neighbours */
     dist_to_nbr = (double *)G_malloc(sides * sizeof(double));
     weight = (double *)G_malloc(sides * sizeof(double));
 
     for (ct_dir = 0; ct_dir < sides; ct_dir++) {
-	/* get r, c (upr, upc) for neighbours */
+	/* get r, c (r_nbr, c_nbr) for neighbours */
 	r_nbr = nextdr[ct_dir];
 	c_nbr = nextdc[ct_dir];
 	/* account for rare cases when ns_res != ew_res */
@@ -111,8 +113,6 @@
 	    dist_to_nbr[ct_dir] = dx + dy;
 	else
 	    dist_to_nbr[ct_dir] = sqrt(dx * dx + dy * dy);
-
-
     }
 
     flag_clear_all(worked);
@@ -133,7 +133,7 @@
 	    c = astar_pts[killer].c;
 	    dc = astar_pts[killer].downc;
 	    value = wat[SEG_INDEX(wat_seg, r, c)];
-	    value_sfd = wat[SEG_INDEX(wat_seg, dr, dc)];
+	    valued = wat[SEG_INDEX(wat_seg, dr, dc)];
 
 	    /* disabled for MFD */
 	    /* if (ABS(value) > threshold)
@@ -150,7 +150,7 @@
 	    ele = alt[SEG_INDEX(alt_seg, r, c)];
 	    /* this loop is needed to get the sum of weights */
 	    for (ct_dir = 0; ct_dir < sides; ct_dir++) {
-		/* get r, c (upr, upc) for neighbours */
+		/* get r, c (r_nbr, c_nbr) for neighbours */
 		r_nbr = r + nextdr[ct_dir];
 		c_nbr = c + nextdc[ct_dir];
 		weight[ct_dir] = -1;
@@ -160,7 +160,6 @@
 		    in_val = FLAG_GET(worked, r_nbr, c_nbr);
 		    if (in_val == 0) {
 			ele_nbr = alt[SEG_INDEX(alt_seg, r_nbr, c_nbr)];
-			/* to consider neighbours with equal ele, change if() */
 			if (ele_nbr < ele) {
 			    weight[ct_dir] =
 				mfd_pow(((ele -
@@ -240,16 +239,13 @@
 		    }
 		}
 		if (ABS(prop - 1.0) > 5E-6f) {
-		    G_warning(_("Cumulative proportion not 1.0 but %f"),
+		    G_warning(_("MFD: cumulative proportion of flow distribution not 1.0 but %f"),
 			      prop);
 		}
-		valued =
-		    ABS(value_sfd) +
-		    ABS(value) * weight[np_side] / sum_weight;
+		valued = wat[SEG_INDEX(wat_seg, dr, dc)];
 	    }
 
 	    if (mfd_cells < 2) {
-		valued = value_sfd;
 		if (value > 0) {
 		    if (valued > 0)
 			valued += value;
@@ -270,7 +266,7 @@
 	    valued = ABS(valued) + 0.5;
 
 	    is_swale = FLAG_GET(swale, r, c);
-	    if (is_swale || (int)valued > threshold) {
+	    if (is_swale || ((int)valued) >= threshold) {
 		FLAG_SET(swale, dr, dc);
 	    }
 	    else {
@@ -282,8 +278,8 @@
     }
     G_percent(count, do_points, 1);	/* finish it */
     if (workedon)
-	G_message(_("A * path already processed: %d of %d"), workedon,
-		  do_points);
+	G_warning(_("MFD: A * path already processed when distributing flow: %d of %d"),
+		  workedon, do_points);
 
     G_free(astar_pts);
 

Modified: grass/trunk/raster/r.watershed/ram/no_stream.c
===================================================================
--- grass/trunk/raster/r.watershed/ram/no_stream.c	2008-12-10 05:35:27 UTC (rev 34820)
+++ grass/trunk/raster/r.watershed/ram/no_stream.c	2008-12-10 15:13:06 UTC (rev 34821)
@@ -6,7 +6,7 @@
     int r, rr, c, cc, uprow = 0, upcol = 0;
     double slope;
     CELL downdir, asp_value, hih_ele, new_ele, aspect;
-    DCELL value, max_drain;	/* flow acc is now DCELL */
+    DCELL dvalue, max_drain;	/* flow acc is now DCELL */
     SHORT updir, riteflag, leftflag, thisdir;
 
     while (1) {
@@ -16,14 +16,13 @@
 		if (r >= 0 && c >= 0 && r < nrows && c < ncols) {
 		    aspect = asp[SEG_INDEX(asp_seg, r, c)];
 		    if (aspect == drain[rr][cc]) {
-			value = wat[SEG_INDEX(wat_seg, r, c)];
-			/* here is the reason why MFD basins differ from SFD basins */
-			if (value < 0)
-			    value = -value;
-			if ((value - max_drain) > 5E-8f) {	/* floating point comparison problem workaround */
+			dvalue = wat[SEG_INDEX(wat_seg, r, c)];
+			if (dvalue < 0)
+			    dvalue = -dvalue;
+			if ((dvalue - max_drain) > 5E-8f) {	/* floating point comparison problem workaround */
 			    uprow = r;
 			    upcol = c;
-			    max_drain = value;
+			    max_drain = dvalue;
 			}
 		    }
 		}

Modified: grass/trunk/raster/r.watershed/seg/close_maps.c
===================================================================
--- grass/trunk/raster/r.watershed/seg/close_maps.c	2008-12-10 05:35:27 UTC (rev 34820)
+++ grass/trunk/raster/r.watershed/seg/close_maps.c	2008-12-10 15:13:06 UTC (rev 34821)
@@ -5,7 +5,7 @@
 {
     struct Colors colors, logcolors;
     int r, c;
-    CELL is_swale, value;
+    CELL is_swale;
     double dvalue;
     struct FPRange accRange;
     DCELL min, max;
@@ -100,7 +100,7 @@
     /* error in gislib.a
        G_free_colors(&colors);
      */
-    cseg_close(&wat);
+    dseg_close(&wat);
     if (ls_flag) {
 	dseg_write_cellfile(&l_s, ls_name);
 	dseg_close(&l_s);

Modified: grass/trunk/raster/r.watershed/seg/do_cum.c
===================================================================
--- grass/trunk/raster/r.watershed/seg/do_cum.c	2008-12-10 05:35:27 UTC (rev 34820)
+++ grass/trunk/raster/r.watershed/seg/do_cum.c	2008-12-10 15:13:06 UTC (rev 34821)
@@ -13,6 +13,7 @@
     int killer, threshold, count;
 
     G_message(_("SECTION 3: Accumulating Surface Flow with SFD."));
+    G_debug(1, "MFD convergence factor set to %d.", c_fac);
 
     count = 0;
     if (bas_thres <= 0)
@@ -87,7 +88,7 @@
 {
     int r, c, dr, dc;
     CELL is_swale;
-    DCELL value, valued, value_sfd;
+    DCELL value, valued;
     POINT point;
     int killer, threshold, count;
 
@@ -106,7 +107,7 @@
     weight = (double *)G_malloc(sides * sizeof(double));
 
     for (ct_dir = 0; ct_dir < sides; ct_dir++) {
-	/* get r, c (upr, upc) for neighbours */
+	/* get r, c (r_nbr, c_nbr) for neighbours */
 	r_nbr = nextdr[ct_dir];
 	c_nbr = nextdc[ct_dir];
 	/* account for rare cases when ns_res != ew_res */
@@ -116,8 +117,6 @@
 	    dist_to_nbr[ct_dir] = dx + dy;
 	else
 	    dist_to_nbr[ct_dir] = sqrt(dx * dx + dy * dy);
-
-
     }
 
     /* reset worked, takes time... */
@@ -145,7 +144,7 @@
 	    c = point.c;
 	    dc = point.downc;
 	    dseg_get(&wat, &value, r, c);
-	    dseg_get(&wat, &value_sfd, dr, dc);
+	    dseg_get(&wat, &valued, dr, dc);
 
 	    /* disabled for MFD */
 	    /* if (ABS(value) > threshold)
@@ -162,7 +161,7 @@
 	    cseg_get(&alt, &ele, r, c);
 	    /* this loop is needed to get the sum of weights */
 	    for (ct_dir = 0; ct_dir < sides; ct_dir++) {
-		/* get r, c (upr, upc) for neighbours */
+		/* get r, c (r_nbr, c_nbr) for neighbours */
 		r_nbr = r + nextdr[ct_dir];
 		c_nbr = c + nextdc[ct_dir];
 		weight[ct_dir] = -1;
@@ -172,7 +171,6 @@
 		    bseg_get(&worked, &in_val, r_nbr, c_nbr);
 		    if (in_val == 0) {
 			cseg_get(&alt, &ele_nbr, r_nbr, c_nbr);
-			/* to consider neighbours with equal ele, change if() */
 			if (ele_nbr < ele) {
 			    weight[ct_dir] =
 				mfd_pow(((ele -
@@ -252,16 +250,13 @@
 		    }
 		}
 		if (ABS(prop - 1.0) > 5E-6f) {
-		    G_warning(_("3 cumulative proportion not 1.0 but %f"),
+		    G_warning(_("MFD: cumulative proportion of flow distribution not 1.0 but %f"),
 			      prop);
 		}
-		valued =
-		    ABS(value_sfd) +
-		    ABS(value) * weight[np_side] / sum_weight;
+		dseg_get(&wat, &valued, dr, dc);
 	    }
 
 	    if (mfd_cells < 2) {
-		valued = value_sfd;
 		if (value > 0) {
 		    if (valued > 0)
 			valued += value;
@@ -282,7 +277,7 @@
 	    valued = ABS(valued) + 0.5;
 
 	    bseg_get(&swale, &is_swale, r, c);
-	    if (is_swale || (int)valued > threshold) {
+	    if (is_swale || (int)valued >= threshold) {
 		bseg_put(&swale, &one, dr, dc);
 	    }
 	    else {
@@ -294,8 +289,8 @@
     }
     G_percent(count, do_points, 1);	/* finish it */
     if (workedon)
-	G_message(_("A * path already processed: %d of %d"), workedon,
-		  do_points);
+	G_warning(_("MFD: A * path already processed when distributing flow: %d of %d"),
+		  workedon, do_points);
 
     seg_close(&astar_pts);
 

Modified: grass/trunk/raster/r.watershed/seg/no_stream.c
===================================================================
--- grass/trunk/raster/r.watershed/seg/no_stream.c	2008-12-10 05:35:27 UTC (rev 34820)
+++ grass/trunk/raster/r.watershed/seg/no_stream.c	2008-12-10 15:13:06 UTC (rev 34821)
@@ -6,8 +6,8 @@
 {
     int r, rr, c, cc, uprow = 0, upcol = 0;
     double slope;
-    CELL downdir, asp_value, hih_ele, new_ele, aspect;
-    DCELL value, max_drain;	/* flow acc is now DCELL */
+    CELL downdir, asp_value, hih_ele, new_ele, aspect, value;
+    DCELL dvalue, max_drain;	/* flow acc is now DCELL */
     SHORT updir, riteflag, leftflag, thisdir;
 
     while (1) {
@@ -18,14 +18,13 @@
 
 		    cseg_get(&asp, &aspect, r, c);
 		    if (aspect == drain[rr][cc]) {
-			dseg_get(&wat, &value, r, c);
-			/* here is the reason why MFD basins differ from SFD basins */
-			if (value < 0)
-			    value = -value;
-			if ((value - max_drain) > 5E-8f) {	/* floating point comparison problem workaround */
+			dseg_get(&wat, &dvalue, r, c);
+			if (dvalue < 0)
+			    dvalue = -dvalue;
+			if ((dvalue - max_drain) > 5E-8f) {	/* floating point comparison problem workaround */
 			    uprow = r;
 			    upcol = c;
-			    max_drain = value;
+			    max_drain = dvalue;
 			}
 		    }
 		}



More information about the grass-commit mailing list