[GRASS-SVN] r47347 - grass/trunk/raster/r.watershed/ram
svn_grass at osgeo.org
svn_grass at osgeo.org
Tue Aug 2 07:54:07 EDT 2011
Author: mmetz
Date: 2011-08-02 04:54:07 -0700 (Tue, 02 Aug 2011)
New Revision: 47347
Modified:
grass/trunk/raster/r.watershed/ram/Gwater.h
grass/trunk/raster/r.watershed/ram/close_maps.c
grass/trunk/raster/r.watershed/ram/def_basin.c
grass/trunk/raster/r.watershed/ram/do_astar.c
grass/trunk/raster/r.watershed/ram/do_cum.c
grass/trunk/raster/r.watershed/ram/find_pour.c
grass/trunk/raster/r.watershed/ram/init_vars.c
grass/trunk/raster/r.watershed/ram/over_cells.c
grass/trunk/raster/r.watershed/ram/split_str.c
Log:
fix depression handling (ram)
Modified: grass/trunk/raster/r.watershed/ram/Gwater.h
===================================================================
--- grass/trunk/raster/r.watershed/ram/Gwater.h 2011-08-02 10:25:10 UTC (rev 47346)
+++ grass/trunk/raster/r.watershed/ram/Gwater.h 2011-08-02 11:54:07 UTC (rev 47347)
@@ -94,7 +94,7 @@
/* do_astar.c */
int do_astar(void);
-int add_pt(int, int, CELL, CELL);
+int add_pt(int, int, CELL);
int drop_pt(void);
int sift_up(int, CELL);
double get_slope(int, int, int, int, CELL, CELL);
Modified: grass/trunk/raster/r.watershed/ram/close_maps.c
===================================================================
--- grass/trunk/raster/r.watershed/ram/close_maps.c 2011-08-02 10:25:10 UTC (rev 47346)
+++ grass/trunk/raster/r.watershed/ram/close_maps.c 2011-08-02 11:54:07 UTC (rev 47347)
@@ -27,157 +27,128 @@
sum = sum_sqr = stddev = 0.0;
if (wat_flag) {
fd = Rast_open_new(wat_name, DCELL_TYPE);
- {
- if (abs_acc) {
- G_warning("Writing out only positive flow accumulation values.");
- G_warning("Cells with a likely underestimate for flow accumulation can no longer be identified.");
- for (r = 0; r < nrows; r++) {
- Rast_set_d_null_value(dbuf, ncols); /* reset row to all NULL */
- for (c = 0; c < ncols; c++) {
- dvalue = wat[SEG_INDEX(wat_seg, r, c)];
- if (Rast_is_d_null_value(&dvalue) == 0 && dvalue) {
- dvalue = ABS(dvalue);
- dbuf[c] = dvalue;
- sum += dvalue;
- sum_sqr += dvalue * dvalue;
- }
+ if (abs_acc) {
+ G_warning("Writing out only positive flow accumulation values.");
+ G_warning("Cells with a likely underestimate for flow accumulation can no longer be identified.");
+ for (r = 0; r < nrows; r++) {
+ Rast_set_d_null_value(dbuf, ncols); /* reset row to all NULL */
+ for (c = 0; c < ncols; c++) {
+ dvalue = wat[SEG_INDEX(wat_seg, r, c)];
+ if (Rast_is_d_null_value(&dvalue) == 0 && dvalue) {
+ dvalue = ABS(dvalue);
+ dbuf[c] = dvalue;
+ sum += dvalue;
+ sum_sqr += dvalue * dvalue;
}
- Rast_put_row(fd, dbuf, DCELL_TYPE);
}
+ Rast_put_row(fd, dbuf, DCELL_TYPE);
}
- else {
- for (r = 0; r < nrows; r++) {
- Rast_set_d_null_value(dbuf, ncols); /* reset row to all NULL */
- for (c = 0; c < ncols; c++) {
- dvalue = wat[SEG_INDEX(wat_seg, r, c)];
- if (Rast_is_d_null_value(&dvalue) == 0 && dvalue) {
- dbuf[c] = dvalue;
- dvalue = ABS(dvalue);
- sum += dvalue;
- sum_sqr += dvalue * dvalue;
- }
+ }
+ else {
+ for (r = 0; r < nrows; r++) {
+ Rast_set_d_null_value(dbuf, ncols); /* reset row to all NULL */
+ for (c = 0; c < ncols; c++) {
+ dvalue = wat[SEG_INDEX(wat_seg, r, c)];
+ if (Rast_is_d_null_value(&dvalue) == 0 && dvalue) {
+ dbuf[c] = dvalue;
+ dvalue = ABS(dvalue);
+ sum += dvalue;
+ sum_sqr += dvalue * dvalue;
}
- Rast_put_row(fd, dbuf, DCELL_TYPE);
}
+ Rast_put_row(fd, dbuf, DCELL_TYPE);
}
- Rast_close(fd);
+ }
+ Rast_close(fd);
- stddev =
- sqrt((sum_sqr - (sum + sum / do_points)) / (do_points - 1));
- G_debug(1, "stddev: %f", stddev);
+ stddev =
+ sqrt((sum_sqr - (sum + sum / do_points)) / (do_points - 1));
+ G_debug(1, "stddev: %f", stddev);
- /* set nice color rules: yellow, green, cyan, blue, black */
+ /* set nice color rules: yellow, green, cyan, blue, black */
- lstddev = log(stddev);
+ lstddev = log(stddev);
- Rast_read_fp_range(wat_name, this_mapset, &accRange);
- min = max = 0;
- Rast_get_fp_range_min_max(&accRange, &min, &max);
+ Rast_read_fp_range(wat_name, this_mapset, &accRange);
+ min = max = 0;
+ Rast_get_fp_range_min_max(&accRange, &min, &max);
- Rast_init_colors(&colors);
+ Rast_init_colors(&colors);
- if (min < 0) {
- if (min < (-stddev - 1)) {
- clr_min = min - 1;
- clr_max = -stddev - 1;
- Rast_add_d_color_rule(&clr_min, 0, 0, 0, &clr_max, 0,
- 0, 0, &colors);
- }
- clr_min = -stddev - 1.;
- clr_max = -1. * exp(lstddev * 0.75);
+ if (min < 0) {
+ if (min < (-stddev - 1)) {
+ clr_min = min - 1;
+ clr_max = -stddev - 1;
Rast_add_d_color_rule(&clr_min, 0, 0, 0, &clr_max, 0,
- 0, 255, &colors);
- clr_min = clr_max;
- clr_max = -1. * exp(lstddev * 0.5);
- Rast_add_d_color_rule(&clr_min, 0, 0, 255, &clr_max, 0,
- 255, 255, &colors);
- clr_min = clr_max;
- clr_max = -1. * exp(lstddev * 0.35);
- Rast_add_d_color_rule(&clr_min, 0, 255, 255, &clr_max, 0,
- 255, 0, &colors);
- clr_min = clr_max;
- clr_max = -1.;
- Rast_add_d_color_rule(&clr_min, 0, 255, 0, &clr_max, 255,
- 255, 0, &colors);
+ 0, 0, &colors);
}
- clr_min = -1.;
- clr_max = 1.;
- Rast_add_d_color_rule(&clr_min, 255, 255, 0, &clr_max, 255,
- 255, 0, &colors);
- clr_min = 1;
- clr_max = exp(lstddev * 0.35);
- Rast_add_d_color_rule(&clr_min, 255, 255, 0, &clr_max, 0,
- 255, 0, &colors);
+ clr_min = -stddev - 1.;
+ clr_max = -1. * exp(lstddev * 0.75);
+ Rast_add_d_color_rule(&clr_min, 0, 0, 0, &clr_max, 0,
+ 0, 255, &colors);
clr_min = clr_max;
- clr_max = exp(lstddev * 0.5);
- Rast_add_d_color_rule(&clr_min, 0, 255, 0, &clr_max, 0,
+ clr_max = -1. * exp(lstddev * 0.5);
+ Rast_add_d_color_rule(&clr_min, 0, 0, 255, &clr_max, 0,
255, 255, &colors);
clr_min = clr_max;
- clr_max = exp(lstddev * 0.75);
+ clr_max = -1. * exp(lstddev * 0.35);
Rast_add_d_color_rule(&clr_min, 0, 255, 255, &clr_max, 0,
- 0, 255, &colors);
+ 255, 0, &colors);
clr_min = clr_max;
- clr_max = stddev + 1.;
- Rast_add_d_color_rule(&clr_min, 0, 0, 255, &clr_max, 0, 0,
- 0, &colors);
+ clr_max = -1.;
+ Rast_add_d_color_rule(&clr_min, 0, 255, 0, &clr_max, 255,
+ 255, 0, &colors);
+ }
+ clr_min = -1.;
+ clr_max = 1.;
+ Rast_add_d_color_rule(&clr_min, 255, 255, 0, &clr_max, 255,
+ 255, 0, &colors);
+ clr_min = 1;
+ clr_max = exp(lstddev * 0.35);
+ Rast_add_d_color_rule(&clr_min, 255, 255, 0, &clr_max, 0,
+ 255, 0, &colors);
+ clr_min = clr_max;
+ clr_max = exp(lstddev * 0.5);
+ Rast_add_d_color_rule(&clr_min, 0, 255, 0, &clr_max, 0,
+ 255, 255, &colors);
+ clr_min = clr_max;
+ clr_max = exp(lstddev * 0.75);
+ Rast_add_d_color_rule(&clr_min, 0, 255, 255, &clr_max, 0,
+ 0, 255, &colors);
+ clr_min = clr_max;
+ clr_max = stddev + 1.;
+ Rast_add_d_color_rule(&clr_min, 0, 0, 255, &clr_max, 0, 0,
+ 0, &colors);
- if (max > 0 && max > stddev + 1) {
- clr_min = stddev + 1;
- clr_max = max + 1;
- Rast_add_d_color_rule(&clr_min, 0, 0, 0, &clr_max, 0, 0,
- 0, &colors);
- }
- Rast_write_colors(wat_name, this_mapset, &colors);
+ if (max > 0 && max > stddev + 1) {
+ clr_min = stddev + 1;
+ clr_max = max + 1;
+ Rast_add_d_color_rule(&clr_min, 0, 0, 0, &clr_max, 0, 0,
+ 0, &colors);
}
+ Rast_write_colors(wat_name, this_mapset, &colors);
}
/* TODO: elevation == NULL -> drainage direction == NULL (wat == 0 where ele == NULL) */
/* keep drainage direction == 0 for real depressions */
if (asp_flag) {
fd = Rast_open_c_new(asp_name);
- {
- for (r = 0; r < nrows; r++) {
- Rast_set_c_null_value(buf, ncols); /* reset row to all NULL */
- for (c = 0; c < ncols; c++) {
+ for (r = 0; r < nrows; r++) {
+ Rast_set_c_null_value(buf, ncols); /* reset row to all NULL */
+ for (c = 0; c < ncols; c++) {
+ dvalue = wat[SEG_INDEX(wat_seg, r, c)];
+ if (!Rast_is_d_null_value(&dvalue))
buf[c] = asp[SEG_INDEX(asp_seg, r, c)];
- }
- Rast_put_row(fd, buf, CELL_TYPE);
}
- Rast_close(fd);
+ Rast_put_row(fd, buf, CELL_TYPE);
}
+ Rast_close(fd);
Rast_init_colors(&colors);
Rast_make_aspect_colors(&colors, 0, 8);
Rast_write_colors(asp_name, this_mapset, &colors);
}
G_free(asp);
- /* visual output no longer needed */
- if (dis_flag) {
- fd = Rast_open_c_new(dis_name);
- {
- if (bas_thres <= 0)
- bas_thres = 60;
- for (r = 0; r < nrows; r++) {
- for (c = 0; c < ncols; c++) {
- buf[c] = wat[SEG_INDEX(wat_seg, r, c)];
- if (buf[c] < 0) {
- buf[c] = 0;
- }
- else {
- value = FLAG_GET(swale, r, c);
- if (value) {
- buf[c] = bas_thres;
- }
- }
- }
- Rast_put_row(fd, buf, CELL_TYPE);
- }
- Rast_close(fd);
- }
- Rast_init_colors(&colors);
- Rast_make_rainbow_colors(&colors, 1, 120);
- Rast_write_colors(dis_name, this_mapset, &colors);
- }
flag_destroy(swale);
/*
Rast_free_colors(&colors);
@@ -186,31 +157,27 @@
if (ls_flag) {
fd = Rast_open_new(ls_name, DCELL_TYPE);
- {
- for (r = 0; r < nrows; r++) {
- for (c = 0; c < ncols; c++) {
- dbuf[c] = l_s[SEG_INDEX(l_s_seg, r, c)];
- }
- Rast_put_row(fd, dbuf, DCELL_TYPE);
+ for (r = 0; r < nrows; r++) {
+ for (c = 0; c < ncols; c++) {
+ dbuf[c] = l_s[SEG_INDEX(l_s_seg, r, c)];
}
- Rast_close(fd);
+ Rast_put_row(fd, dbuf, DCELL_TYPE);
}
+ Rast_close(fd);
G_free(l_s);
}
if (sl_flag) {
fd = Rast_open_new(sl_name, DCELL_TYPE);
- {
- for (r = 0; r < nrows; r++) {
- for (c = 0; c < ncols; c++) {
- dbuf[c] = s_l[SEG_INDEX(s_l_seg, r, c)];
- if (dbuf[c] > max_length)
- dbuf[c] = max_length;
- }
- Rast_put_row(fd, dbuf, DCELL_TYPE);
+ for (r = 0; r < nrows; r++) {
+ for (c = 0; c < ncols; c++) {
+ dbuf[c] = s_l[SEG_INDEX(s_l_seg, r, c)];
+ if (dbuf[c] > max_length)
+ dbuf[c] = max_length;
}
- Rast_close(fd);
+ Rast_put_row(fd, dbuf, DCELL_TYPE);
}
+ Rast_close(fd);
}
if (sl_flag || ls_flag || sg_flag)
@@ -218,15 +185,13 @@
if (sg_flag) {
fd = Rast_open_new(sg_name, DCELL_TYPE);
- {
- for (r = 0; r < nrows; r++) {
- for (c = 0; c < ncols; c++) {
- dbuf[c] = s_g[SEG_INDEX(s_g_seg, r, c)];
- }
- Rast_put_row(fd, dbuf, DCELL_TYPE);
+ for (r = 0; r < nrows; r++) {
+ for (c = 0; c < ncols; c++) {
+ dbuf[c] = s_g[SEG_INDEX(s_g_seg, r, c)];
}
- Rast_close(fd);
+ Rast_put_row(fd, dbuf, DCELL_TYPE);
}
+ Rast_close(fd);
G_free(s_g);
}
Modified: grass/trunk/raster/r.watershed/ram/def_basin.c
===================================================================
--- grass/trunk/raster/r.watershed/ram/def_basin.c 2011-08-02 10:25:10 UTC (rev 47346)
+++ grass/trunk/raster/r.watershed/ram/def_basin.c 2011-08-02 11:54:07 UTC (rev 47347)
@@ -14,6 +14,8 @@
for (r = row - 1, rr = 0; rr < 3; r++, rr++) {
for (c = col - 1, cc = 0; cc < 3; c++, cc++) {
if (r >= 0 && c >= 0 && r < nrows && c < ncols) {
+ if (r == row && c == col)
+ continue;
value = asp[SEG_INDEX(asp_seg, r, c)];
if (value < -1)
value = -value;
@@ -44,6 +46,8 @@
for (r = row - 1, rr = 0; rr < 3; r++, rr++) {
for (c = col - 1, cc = 0; cc < 3; c++, cc++) {
if (r >= 0 && c >= 0 && r < nrows && c < ncols) {
+ if (r == row && c == col)
+ continue;
direction = asp[SEG_INDEX(asp_seg, r, c)];
if (direction == drain[rr][cc]) {
thisdir = updrain[rr][cc];
Modified: grass/trunk/raster/r.watershed/ram/do_astar.c
===================================================================
--- grass/trunk/raster/r.watershed/ram/do_astar.c 2011-08-02 10:25:10 UTC (rev 47346)
+++ grass/trunk/raster/r.watershed/ram/do_astar.c 2011-08-02 11:54:07 UTC (rev 47347)
@@ -153,7 +153,7 @@
/* add neighbour as new point if not in the list */
if (is_in_list == 0 && skip_diag == 0) {
- add_pt(upr, upc, alt_nbr[ct_dir], alt_val);
+ add_pt(upr, upc, alt_nbr[ct_dir]);
/* set flow direction */
asp[index_up] = drain[upr - r + 1][upc - c + 1];
}
@@ -165,6 +165,10 @@
if (wat[index_doer] > 0)
wat[index_doer] = -wat[index_doer];
}
+ /* neighbour is inside real depression, not yet worked */
+ else if (asp[index_up] == 0) {
+ asp[index_up] = drain[upr - r + 1][upc - c + 1];
+ }
}
} /* end if in region */
} /* end sides */
@@ -205,7 +209,7 @@
}
/* add point routine for min heap */
-int add_pt(int r, int c, CELL ele, CELL downe)
+int add_pt(int r, int c, CELL ele)
{
FLAG_SET(in_list, r, c);
Modified: grass/trunk/raster/r.watershed/ram/do_cum.c
===================================================================
--- grass/trunk/raster/r.watershed/ram/do_cum.c 2011-08-02 10:25:10 UTC (rev 47346)
+++ grass/trunk/raster/r.watershed/ram/do_cum.c 2011-08-02 11:54:07 UTC (rev 47347)
@@ -53,13 +53,6 @@
wat[down_index] = valued;
valued = ABS(valued) + 0.5;
is_swale = FLAG_GET(swale, r, c);
- /* update asp for depression */
- if (is_swale && pit_flag) {
- if (aspect > 0 && asp[down_index] == 0) {
- aspect = -aspect;
- asp[this_index] = aspect;
- }
- }
if (is_swale || ((int)valued) >= threshold) {
FLAG_SET(swale, dr, dc);
}
@@ -340,13 +333,6 @@
asp[this_index] = aspect;
}
is_swale = FLAG_GET(swale, r, c);
- /* update asp for depression */
- if (is_swale && pit_flag) {
- if (aspect > 0 && asp[SEG_INDEX(asp_seg, r_max, c_max)] == 0) {
- aspect = -aspect;
- asp[this_index] = aspect;
- }
- }
/* start new stream */
value = ABS(value) + 0.5;
if (!is_swale && (int)value >= threshold && stream_cells < 1 &&
@@ -380,14 +366,10 @@
double mfd_pow(double base, int exp)
{
- int x;
- double result;
+ int i;
+ double result = base;
- result = base;
- if (exp == 1)
- return result;
-
- for (x = 2; x <= exp; x++) {
+ for (i = 2; i <= exp; i++) {
result *= base;
}
return result;
Modified: grass/trunk/raster/r.watershed/ram/find_pour.c
===================================================================
--- grass/trunk/raster/r.watershed/ram/find_pour.c 2011-08-02 10:25:10 UTC (rev 47346)
+++ grass/trunk/raster/r.watershed/ram/find_pour.c 2011-08-02 11:54:07 UTC (rev 47347)
@@ -16,7 +16,7 @@
northing = window.north - (row + .5) * window.ns_res;
for (col = 0; col < ncols; col++) {
value = FLAG_GET(swale, row, col);
- if (value && asp[SEG_INDEX(asp_seg, row, col)] < 0) {
+ if (value && asp[SEG_INDEX(asp_seg, row, col)] <= 0) {
basin_num += 2;
if (arm_flag) {
easting = window.west + (col + .5) * window.ew_res;
Modified: grass/trunk/raster/r.watershed/ram/init_vars.c
===================================================================
--- grass/trunk/raster/r.watershed/ram/init_vars.c 2011-08-02 10:25:10 UTC (rev 47346)
+++ grass/trunk/raster/r.watershed/ram/init_vars.c 2011-08-02 11:54:07 UTC (rev 47347)
@@ -238,23 +238,7 @@
G_free(dbuf);
}
- /* depression: drainage direction will be set to zero later */
- if (pit_flag) {
- buf = Rast_allocate_c_buf();
- fd = Rast_open_old(pit_name, "");
- for (r = 0; r < nrows; r++) {
- Rast_get_c_row(fd, buf, r);
- for (c = 0; c < ncols; c++) {
- asp_value = buf[c];
- if (!Rast_is_c_null_value(&asp_value) && asp_value)
- asp[SEG_INDEX(asp_seg, r, c)] = 1;
- }
- }
- Rast_close(fd);
- G_free(buf);
- }
-
- /* this is also creating streams... */
+ /* overland blocking map; this is also creating streams... */
if (ob_flag) {
buf = Rast_allocate_c_buf();
fd = Rast_open_old(ob_name, "");
@@ -302,129 +286,32 @@
/* heap is empty */
heap_size = 0;
+ if (pit_flag) {
+ buf = Rast_allocate_c_buf();
+ fd = Rast_open_old(pit_name, "");
+ }
+ else
+ buf = NULL;
first_astar = first_cum = -1;
- if (MASK_flag) {
- for (r = 0; r < nrows; r++) {
- G_percent(r, nrows, 3);
- for (c = 0; c < ncols; c++) {
- seg_idx = SEG_INDEX(wat_seg, r, c);
- if (FLAG_GET(worked, r, c)) {
- wat[seg_idx] = 0;
- }
- else {
- if (er_flag)
- s_l[seg_idx] = half_res;
- asp_value = asp[seg_idx];
- if (r == 0 || c == 0 || r == nrows - 1 ||
- c == ncols - 1 || asp_value != 0) {
- wat_value = wat[seg_idx];
- if (wat_value > 0)
- wat[seg_idx] = -wat_value;
- /* set depression */
- if (asp_value) {
- asp_value = 0;
- wat[seg_idx] = ABS(wat_value);
- }
- else if (r == 0)
- asp_value = -2;
- else if (c == 0)
- asp_value = -4;
- else if (r == nrows - 1)
- asp_value = -6;
- else if (c == ncols - 1)
- asp_value = -8;
- asp[seg_idx] = asp_value;
- alt_value = alt[seg_idx];
- add_pt(r, c, alt_value, alt_value);
- }
- else if (FLAG_GET(worked, r - 1, c)) {
- alt_value = alt[seg_idx];
- add_pt(r, c, alt_value, alt_value);
- asp[seg_idx] = -2;
- wat_value = wat[seg_idx];
- if (wat_value > 0)
- wat[seg_idx] = -wat_value;
- }
- else if (FLAG_GET(worked, r + 1, c)) {
- alt_value = alt[seg_idx];
- add_pt(r, c, alt_value, alt_value);
- asp[seg_idx] = -6;
- wat_value = wat[seg_idx];
- if (wat_value > 0)
- wat[seg_idx] = -wat_value;
- }
- else if (FLAG_GET(worked, r, c - 1)) {
- alt_value = alt[seg_idx];
- add_pt(r, c, alt_value, alt_value);
- asp[seg_idx] = -4;
- wat_value = wat[seg_idx];
- if (wat_value > 0)
- wat[seg_idx] = -wat_value;
- }
- else if (FLAG_GET(worked, r, c + 1)) {
- alt_value = alt[seg_idx];
- add_pt(r, c, alt_value, alt_value);
- asp[seg_idx] = -8;
- wat_value = wat[seg_idx];
- if (wat_value > 0)
- wat[seg_idx] = -wat_value;
- }
- else if (sides == 8 && FLAG_GET(worked, r - 1, c - 1)) {
- alt_value = alt[seg_idx];
- add_pt(r, c, alt_value, alt_value);
- asp[seg_idx] = -3;
- wat_value = wat[seg_idx];
- if (wat_value > 0)
- wat[seg_idx] = -wat_value;
- }
- else if (sides == 8 && FLAG_GET(worked, r - 1, c + 1)) {
- alt_value = alt[seg_idx];
- add_pt(r, c, alt_value, alt_value);
- asp[seg_idx] = -1;
- wat_value = wat[seg_idx];
- if (wat_value > 0)
- wat[seg_idx] = -wat_value;
- }
- else if (sides == 8 && FLAG_GET(worked, r + 1, c - 1)) {
- alt_value = alt[seg_idx];
- add_pt(r, c, alt_value, alt_value);
- asp[seg_idx] = -5;
- wat_value = wat[seg_idx];
- if (wat_value > 0)
- wat[seg_idx] = -wat_value;
- }
- else if (sides == 8 && FLAG_GET(worked, r + 1, c + 1)) {
- alt_value = alt[seg_idx];
- add_pt(r, c, alt_value, alt_value);
- asp[seg_idx] = -7;
- wat_value = wat[seg_idx];
- if (wat_value > 0)
- wat[seg_idx] = -wat_value;
- }
- }
+ for (r = 0; r < nrows; r++) {
+ G_percent(r, nrows, 3);
+ if (pit_flag)
+ Rast_get_c_row(fd, buf, r);
+ for (c = 0; c < ncols; c++) {
+ seg_idx = SEG_INDEX(wat_seg, r, c);
+ if (FLAG_GET(worked, r, c)) {
+ wat[seg_idx] = 0;
}
- }
- }
- else {
- for (r = 0; r < nrows; r++) {
- G_percent(r, nrows, 3);
- for (c = 0; c < ncols; c++) {
- seg_idx = SEG_INDEX(wat_seg, r, c);
+ else {
+ asp_value = asp[seg_idx];
if (er_flag)
s_l[seg_idx] = half_res;
- asp_value = asp[seg_idx];
if (r == 0 || c == 0 || r == nrows - 1 ||
- c == ncols - 1 || asp_value != 0) {
+ c == ncols - 1) {
wat_value = wat[seg_idx];
- if (wat_value > 0) {
+ if (wat_value > 0)
wat[seg_idx] = -wat_value;
- }
- /* set depression */
- if (asp_value) {
- asp_value = 0;
- wat[seg_idx] = ABS(wat_value);
- }
- else if (r == 0)
+ if (r == 0)
asp_value = -2;
else if (c == 0)
asp_value = -4;
@@ -434,14 +321,91 @@
asp_value = -8;
asp[seg_idx] = asp_value;
alt_value = alt[seg_idx];
- add_pt(r, c, alt_value, alt_value);
+ add_pt(r, c, alt_value);
}
+ else if (FLAG_GET(worked, r - 1, c)) {
+ alt_value = alt[seg_idx];
+ add_pt(r, c, alt_value);
+ asp[seg_idx] = -2;
+ wat_value = wat[seg_idx];
+ if (wat_value > 0)
+ wat[seg_idx] = -wat_value;
+ }
+ else if (FLAG_GET(worked, r + 1, c)) {
+ alt_value = alt[seg_idx];
+ add_pt(r, c, alt_value);
+ asp[seg_idx] = -6;
+ wat_value = wat[seg_idx];
+ if (wat_value > 0)
+ wat[seg_idx] = -wat_value;
+ }
+ else if (FLAG_GET(worked, r, c - 1)) {
+ alt_value = alt[seg_idx];
+ add_pt(r, c, alt_value);
+ asp[seg_idx] = -4;
+ wat_value = wat[seg_idx];
+ if (wat_value > 0)
+ wat[seg_idx] = -wat_value;
+ }
+ else if (FLAG_GET(worked, r, c + 1)) {
+ alt_value = alt[seg_idx];
+ add_pt(r, c, alt_value);
+ asp[seg_idx] = -8;
+ wat_value = wat[seg_idx];
+ if (wat_value > 0)
+ wat[seg_idx] = -wat_value;
+ }
+ else if (sides == 8 && FLAG_GET(worked, r - 1, c - 1)) {
+ alt_value = alt[seg_idx];
+ add_pt(r, c, alt_value);
+ asp[seg_idx] = -3;
+ wat_value = wat[seg_idx];
+ if (wat_value > 0)
+ wat[seg_idx] = -wat_value;
+ }
+ else if (sides == 8 && FLAG_GET(worked, r - 1, c + 1)) {
+ alt_value = alt[seg_idx];
+ add_pt(r, c, alt_value);
+ asp[seg_idx] = -1;
+ wat_value = wat[seg_idx];
+ if (wat_value > 0)
+ wat[seg_idx] = -wat_value;
+ }
+ else if (sides == 8 && FLAG_GET(worked, r + 1, c - 1)) {
+ alt_value = alt[seg_idx];
+ add_pt(r, c, alt_value);
+ asp[seg_idx] = -5;
+ wat_value = wat[seg_idx];
+ if (wat_value > 0)
+ wat[seg_idx] = -wat_value;
+ }
+ else if (sides == 8 && FLAG_GET(worked, r + 1, c + 1)) {
+ alt_value = alt[seg_idx];
+ add_pt(r, c, alt_value);
+ asp[seg_idx] = -7;
+ wat_value = wat[seg_idx];
+ if (wat_value > 0)
+ wat[seg_idx] = -wat_value;
+ }
+
+ /* real depression ? */
+ if (pit_flag && asp[seg_idx] == 0) {
+ if (!Rast_is_c_null_value(&buf[c]) && buf[c] != 0) {
+ alt_value = alt[seg_idx];
+ add_pt(r, c, alt_value);
+ }
+ }
}
}
}
G_percent(r, nrows, 1); /* finish it */
+ if (pit_flag) {
+ Rast_close(fd);
+ G_free(buf);
+ }
+
return 0;
}
Modified: grass/trunk/raster/r.watershed/ram/over_cells.c
===================================================================
--- grass/trunk/raster/r.watershed/ram/over_cells.c 2011-08-02 10:25:10 UTC (rev 47346)
+++ grass/trunk/raster/r.watershed/ram/over_cells.c 2011-08-02 11:54:07 UTC (rev 47347)
@@ -13,6 +13,8 @@
for (r = row - 1, rr = 0; r <= row + 1; r++, rr++) {
for (c = col - 1, cc = 0; c <= col + 1; c++, cc++) {
if (r >= 0 && c >= 0 && r < nrows && c < ncols) {
+ if (r == row && c == col)
+ continue;
value = asp[SEG_INDEX(asp_seg, r, c)];
if (value == drain[rr][cc]) {
overland_cells(r, c, basin_num, haf_num, &new_ele);
@@ -58,6 +60,8 @@
for (r = next_r - 1, rr = 0; r <= next_r + 1; r++, rr++) {
for (c = next_c - 1, cc = 0; c <= next_c + 1; c++, cc++) {
if (r >= 0 && c >= 0 && r < nrows && c < ncols) {
+ if (r == row && c == col)
+ continue;
idx = SEG_INDEX(bas_seg, r, c);
if (asp[idx] == drain[rr][cc]) {
if (top >= ocs_alloced) {
Modified: grass/trunk/raster/r.watershed/ram/split_str.c
===================================================================
--- grass/trunk/raster/r.watershed/ram/split_str.c 2011-08-02 10:25:10 UTC (rev 47346)
+++ grass/trunk/raster/r.watershed/ram/split_str.c 2011-08-02 11:54:07 UTC (rev 47347)
@@ -22,6 +22,8 @@
for (r = row - 1, rr = 0; rr < 3; r++, rr++) {
for (c = col - 1, cc = 0; cc < 3; c++, cc++) {
if (r >= 0 && c >= 0 && r < nrows && c < ncols) {
+ if (r == row && c == col)
+ continue;
aspect = asp[SEG_INDEX(asp_seg, r, c)];
if (aspect == drain[rr][cc]) {
doit = 1;
More information about the grass-commit
mailing list