[GRASS-SVN] r55649 - grass-addons/grass7/raster/r.stream.extract
svn_grass at osgeo.org
svn_grass at osgeo.org
Sat Apr 6 13:41:23 PDT 2013
Author: mmetz
Date: 2013-04-06 13:41:23 -0700 (Sat, 06 Apr 2013)
New Revision: 55649
Modified:
grass-addons/grass7/raster/r.stream.extract/close.c
grass-addons/grass7/raster/r.stream.extract/del_streams.c
grass-addons/grass7/raster/r.stream.extract/do_astar.c
grass-addons/grass7/raster/r.stream.extract/init_search.c
grass-addons/grass7/raster/r.stream.extract/load.c
grass-addons/grass7/raster/r.stream.extract/local_proto.h
grass-addons/grass7/raster/r.stream.extract/main.c
grass-addons/grass7/raster/r.stream.extract/streams.c
grass-addons/grass7/raster/r.stream.extract/thin.c
Log:
r.stream.extract: sync to r.watershed
Modified: grass-addons/grass7/raster/r.stream.extract/close.c
===================================================================
--- grass-addons/grass7/raster/r.stream.extract/close.c 2013-04-06 20:39:39 UTC (rev 55648)
+++ grass-addons/grass7/raster/r.stream.extract/close.c 2013-04-06 20:41:23 UTC (rev 55649)
@@ -6,9 +6,10 @@
int close_streamvect(char *stream_vect)
{
- int i, r, c, r_nbr, c_nbr, done;
+ int r, c, r_nbr, c_nbr, done;
+ GW_LARGE_INT i;
CELL stream_id, stream_nbr;
- char aspect;
+ ASP_FLAG af;
int next_node;
struct sstack
{
@@ -127,10 +128,10 @@
Vect_write_line(&Out, GV_POINT, Points, Cats);
- bseg_get(&asp, &aspect, r_nbr, c_nbr);
- while (aspect > 0) {
- r_nbr = r_nbr + asp_r[(int)aspect];
- c_nbr = c_nbr + asp_c[(int)aspect];
+ seg_get(&aspflag, (char *)&af, r_nbr, c_nbr);
+ while (af.asp > 0) {
+ r_nbr = r_nbr + asp_r[(int)af.asp];
+ c_nbr = c_nbr + asp_c[(int)af.asp];
cseg_get(&stream, &stream_nbr, r_nbr, c_nbr);
if (stream_nbr <= 0)
@@ -142,7 +143,7 @@
/* first point of parent stream */
break;
}
- bseg_get(&asp, &aspect, r_nbr, c_nbr);
+ seg_get(&aspflag, (char *)&af, r_nbr, c_nbr);
}
Vect_write_line(&Out, GV_LINE, Points, Cats);
@@ -153,7 +154,7 @@
}
G_percent(n_outlets, n_outlets, 1); /* finish it */
- G_message(_("Writing vector attribute table"));
+ G_message(_("Write vector attribute table"));
/* Prepeare strings for use in db_* calls */
db_init_string(&dbsql);
@@ -198,7 +199,7 @@
/* stream nodes */
for (i = 1; i <= n_stream_nodes; i++) {
- sprintf(buf, "insert into %s values ( %d, \'%s\', %d )",
+ sprintf(buf, "insert into %s values ( %lld, \'%s\', %d )",
Fi->table, i,
(stream_node[i].n_trib > 0 ? "intermediate" : "start"),
(stream_node[i].n_trib > 0));
@@ -236,16 +237,15 @@
int stream_fd, dir_fd, r, c, i;
CELL *cell_buf1, *cell_buf2;
struct History history;
- char flag_value;
CELL stream_id;
- char aspect;
+ ASP_FLAG af;
/* cheating... */
stream_fd = dir_fd = -1;
cell_buf1 = cell_buf2 = NULL;
G_message(_("Writing raster %s"),
- (stream_rast != NULL) + (dir_rast != NULL) > 1 ? _("maps") : _("map"));
+ (stream_rast != NULL) + (dir_rast != NULL) > 1 ? "maps" : "map");
/* write requested output rasters */
if (stream_rast) {
@@ -271,10 +271,9 @@
cell_buf1[c] = stream_id;
}
if (dir_rast) {
- bseg_get(&bitflags, &flag_value, r, c);
- if (!FLAG_GET(flag_value, NULLFLAG)) {
- bseg_get(&asp, &aspect, r, c);
- cell_buf2[c] = aspect;
+ seg_get(&aspflag, (char *)&af, r, c);
+ if (!FLAG_GET(af.flag, NULLFLAG)) {
+ cell_buf2[c] = af.asp;
}
}
@@ -294,11 +293,18 @@
Rast_write_history(stream_rast, &history);
}
if (dir_rast) {
+ struct Colors colors;
+
Rast_close(dir_fd);
G_free(cell_buf2);
+
Rast_short_history(dir_rast, "raster", &history);
Rast_command_history(&history);
Rast_write_history(dir_rast, &history);
+
+ Rast_init_colors(&colors);
+ Rast_make_aspect_colors(&colors, -8, 8);
+ Rast_write_colors(dir_rast, G_mapset(), &colors);
}
/* close stream vector */
Modified: grass-addons/grass7/raster/r.stream.extract/del_streams.c
===================================================================
--- grass-addons/grass7/raster/r.stream.extract/del_streams.c 2013-04-06 20:39:39 UTC (rev 55648)
+++ grass-addons/grass7/raster/r.stream.extract/del_streams.c 2013-04-06 20:41:23 UTC (rev 55649)
@@ -12,7 +12,7 @@
CELL curr_stream;
int asp_r[9] = { 0, -1, -1, -1, 0, 1, 1, 1, 0 };
int asp_c[9] = { 0, 1, 0, -1, -1, -1, 0, 1, 1 };
- char aspect;
+ ASP_FLAG af;
r = stream_node[stream_id].r;
c = stream_node[stream_id].c;
@@ -20,10 +20,10 @@
*next_stream_id = stream_id;
/* get next downstream point */
- bseg_get(&asp, &aspect, r, c);
- while (aspect > 0) {
- r_nbr = r + asp_r[(int)aspect];
- c_nbr = c + asp_c[(int)aspect];
+ seg_get(&aspflag, (char *)&af, r, c);
+ while (af.asp > 0) {
+ r_nbr = r + asp_r[(int)af.asp];
+ c_nbr = c + asp_c[(int)af.asp];
/* user-defined depression */
if (r_nbr == r && c_nbr == c)
@@ -40,7 +40,7 @@
slength++;
r = r_nbr;
c = c_nbr;
- bseg_get(&asp, &aspect, r, c);
+ seg_get(&aspflag, (char *)&af, r, c);
}
return slength;
@@ -54,7 +54,7 @@
CELL curr_stream;
int asp_r[9] = { 0, -1, -1, -1, 0, 1, 1, 1, 0 };
int asp_c[9] = { 0, 1, 0, -1, -1, -1, 0, 1, 1 };
- char aspect;
+ ASP_FLAG af;
r = stream_node[stream_id].r;
c = stream_node[stream_id].c;
@@ -65,10 +65,10 @@
curr_stream = stream_id;
/* get next downstream point */
- bseg_get(&asp, &aspect, r, c);
- while (aspect > 0) {
- r_nbr = r + asp_r[(int)aspect];
- c_nbr = c + asp_c[(int)aspect];
+ seg_get(&aspflag, (char *)&af, r, c);
+ while (af.asp > 0) {
+ r_nbr = r + asp_r[(int)af.asp];
+ c_nbr = c + asp_c[(int)af.asp];
/* user-defined depression */
if (r_nbr == r && c_nbr == c)
@@ -83,7 +83,7 @@
r = r_nbr;
c = c_nbr;
cseg_put(&stream, &new_stream, r, c);
- bseg_get(&asp, &aspect, r, c);
+ seg_get(&aspflag, (char *)&af, r, c);
}
if (curr_stream <= 0)
Modified: grass-addons/grass7/raster/r.stream.extract/do_astar.c
===================================================================
--- grass-addons/grass7/raster/r.stream.extract/do_astar.c 2013-04-06 20:39:39 UTC (rev 55648)
+++ grass-addons/grass7/raster/r.stream.extract/do_astar.c 2013-04-06 20:41:23 UTC (rev 55649)
@@ -8,8 +8,7 @@
#define GET_CHILD(p) (((p) << 3) - 6)
HEAP_PNT heap_drop(void);
-int sift_up(GW_LARGE_INT, HEAP_PNT);
-double get_slope(CELL, CELL, double);
+static double get_slope(CELL, CELL, double);
int do_astar(void)
{
@@ -19,8 +18,8 @@
int nextdc[8] = { 0, 0, -1, 1, 1, -1, 1, -1 };
CELL ele_val, ele_up, ele_nbr[8];
WAT_ALT wa;
- char asp_val;
- char flag_value, is_in_list, is_worked;
+ ASP_FLAG af;
+ char is_in_list, is_worked;
HEAP_PNT heap_p;
/* sides
* |7|1|4|
@@ -86,9 +85,9 @@
if (r_nbr < 0 || r_nbr >= nrows || c_nbr < 0 || c_nbr >= ncols)
continue;
- bseg_get(&bitflags, &flag_value, r_nbr, c_nbr);
- is_in_list = FLAG_GET(flag_value, INLISTFLAG);
- is_worked = FLAG_GET(flag_value, WORKEDFLAG);
+ seg_get(&aspflag, (char *)&af, r_nbr, c_nbr);
+ is_in_list = FLAG_GET(af.flag, INLISTFLAG);
+ is_worked = FLAG_GET(af.flag, WORKEDFLAG);
if (!is_worked) {
seg_get(&watalt, (char *)&wa, r_nbr, c_nbr);
ele_nbr[ct_dir] = wa.ele;
@@ -115,43 +114,41 @@
}
}
- if (is_in_list == 0 && skip_me == 0) {
- ele_up = ele_nbr[ct_dir];
- asp_val = drain[r_nbr - r + 1][c_nbr - c + 1];
- bseg_put(&asp, &asp_val, r_nbr, c_nbr);
- heap_add(r_nbr, c_nbr, ele_up);
- FLAG_SET(flag_value, INLISTFLAG);
- bseg_put(&bitflags, &flag_value, r_nbr, c_nbr);
- }
- else if (is_in_list && is_worked == 0) {
- if (FLAG_GET(flag_value, EDGEFLAG)) {
- /* neighbour is edge in list, not yet worked */
- bseg_get(&asp, &asp_val, r_nbr, c_nbr);
- if (asp_val < 0) {
- /* adjust flow direction for edge cell */
- asp_val = drain[r_nbr - r + 1][c_nbr - c + 1];
- bseg_put(&asp, &asp_val, r_nbr, c_nbr);
- }
+ if (!skip_me) {
+ if (is_in_list == 0) {
+ ele_up = ele_nbr[ct_dir];
+ af.asp = drain[r_nbr - r + 1][c_nbr - c + 1];
+ heap_add(r_nbr, c_nbr, ele_up);
+ FLAG_SET(af.flag, INLISTFLAG);
+ seg_put(&aspflag, (char *)&af, r_nbr, c_nbr);
}
- else if (FLAG_GET(flag_value, DEPRFLAG)) {
- G_debug(3, "real depression");
- /* neighbour is inside real depression, not yet worked */
- bseg_get(&asp, &asp_val, r_nbr, c_nbr);
- if (asp_val == 0 && ele_val <= ele_nbr[ct_dir]) {
- asp_val = drain[r_nbr - r + 1][c_nbr - c + 1];
- bseg_put(&asp, &asp_val, r_nbr, c_nbr);
- FLAG_UNSET(flag_value, DEPRFLAG);
- bseg_put(&bitflags, &flag_value, r_nbr, c_nbr);
+ else if (is_in_list && is_worked == 0) {
+ if (FLAG_GET(af.flag, EDGEFLAG)) {
+ /* neighbour is edge in list, not yet worked */
+ if (af.asp < 0) {
+ /* adjust flow direction for edge cell */
+ af.asp = drain[r_nbr - r + 1][c_nbr - c + 1];
+ seg_put(&aspflag, (char *)&af, r_nbr, c_nbr);
+ }
}
+ else if (FLAG_GET(af.flag, DEPRFLAG)) {
+ G_debug(3, "real depression");
+ /* neighbour is inside real depression, not yet worked */
+ if (af.asp == 0 && ele_val <= ele_nbr[ct_dir]) {
+ af.asp = drain[r_nbr - r + 1][c_nbr - c + 1];
+ FLAG_UNSET(af.flag, DEPRFLAG);
+ seg_put(&aspflag, (char *)&af, r_nbr, c_nbr);
+ }
+ }
}
}
} /* end neighbours */
/* add astar points to sorted list for flow accumulation and stream extraction */
first_cum--;
seg_put(&astar_pts, (char *)&heap_p.pnt, 0, first_cum);
- bseg_get(&bitflags, &flag_value, r, c);
- FLAG_SET(flag_value, WORKEDFLAG);
- bseg_put(&bitflags, &flag_value, r, c);
+ seg_get(&aspflag, (char *)&af, r, c);
+ FLAG_SET(af.flag, WORKEDFLAG);
+ seg_put(&aspflag, (char *)&af, r, c);
} /* end A* search */
G_percent(n_points, n_points, 1); /* finish it */
@@ -163,7 +160,7 @@
* compare function for heap
* returns 1 if point1 < point2 else 0
*/
-int heap_cmp(HEAP_PNT *a, HEAP_PNT *b)
+static int heap_cmp(HEAP_PNT *a, HEAP_PNT *b)
{
if (a->ele < b->ele)
return 1;
@@ -174,7 +171,7 @@
return 0;
}
-int sift_up(GW_LARGE_INT start, HEAP_PNT child_p)
+static int sift_up(GW_LARGE_INT start, HEAP_PNT child_p)
{
GW_LARGE_INT parent, child;
HEAP_PNT heap_p;
@@ -282,7 +279,7 @@
return root_p;
}
-double get_slope(CELL ele, CELL up_ele, double dist)
+static double get_slope(CELL ele, CELL up_ele, double dist)
{
if (ele >= up_ele)
return 0.0;
Modified: grass-addons/grass7/raster/r.stream.extract/init_search.c
===================================================================
--- grass-addons/grass7/raster/r.stream.extract/init_search.c 2013-04-06 20:39:39 UTC (rev 55648)
+++ grass-addons/grass7/raster/r.stream.extract/init_search.c 2013-04-06 20:41:23 UTC (rev 55649)
@@ -8,9 +8,9 @@
CELL *depr_buf, ele_value;
int nextdr[8] = { 1, -1, 0, 0, -1, 1, 1, -1 };
int nextdc[8] = { 0, 0, -1, 1, 1, -1, 1, -1 };
- char flag_value, flag_value_nbr, is_null;
+ char asp_value, is_null;
WAT_ALT wa;
- char asp_value;
+ ASP_FLAG af, af_nbr;
GW_LARGE_INT n_depr_cells = 0;
nxt_avail_pt = heap_size = 0;
@@ -31,8 +31,8 @@
for (c = 0; c < ncols; c++) {
- bseg_get(&bitflags, &flag_value, r, c);
- is_null = FLAG_GET(flag_value, NULLFLAG);
+ seg_get(&aspflag, (char *)&af, r, c);
+ is_null = FLAG_GET(af.flag, NULLFLAG);
if (is_null)
continue;
@@ -60,10 +60,10 @@
seg_get(&watalt, (char *)&wa, r, c);
ele_value = wa.ele;
heap_add(r, c, ele_value);
- FLAG_SET(flag_value, INLISTFLAG);
- FLAG_SET(flag_value, EDGEFLAG);
- bseg_put(&bitflags, &flag_value, r, c);
- bseg_put(&asp, &asp_value, r, c);
+ FLAG_SET(af.flag, INLISTFLAG);
+ FLAG_SET(af.flag, EDGEFLAG);
+ af.asp = asp_value;
+ seg_put(&aspflag, (char *)&af, r, c);
continue;
}
@@ -73,18 +73,18 @@
r_nbr = r + nextdr[ct_dir];
c_nbr = c + nextdc[ct_dir];
- bseg_get(&bitflags, &flag_value_nbr, r_nbr, c_nbr);
- is_null = FLAG_GET(flag_value_nbr, NULLFLAG);
+ seg_get(&aspflag, (char *)&af_nbr, r_nbr, c_nbr);
+ is_null = FLAG_GET(af_nbr.flag, NULLFLAG);
if (is_null) {
asp_value = -1 * drain[r - r_nbr + 1][c - c_nbr + 1];
seg_get(&watalt, (char *)&wa, r, c);
ele_value = wa.ele;
heap_add(r, c, ele_value);
- FLAG_SET(flag_value, INLISTFLAG);
- FLAG_SET(flag_value, EDGEFLAG);
- bseg_put(&bitflags, &flag_value, r, c);
- bseg_put(&asp, &asp_value, r, c);
+ FLAG_SET(af.flag, INLISTFLAG);
+ FLAG_SET(af.flag, EDGEFLAG);
+ af.asp = asp_value;
+ seg_put(&aspflag, (char *)&af, r, c);
break;
}
@@ -98,10 +98,10 @@
seg_get(&watalt, (char *)&wa, r, c);
ele_value = wa.ele;
heap_add(r, c, ele_value);
- FLAG_SET(flag_value, INLISTFLAG);
- FLAG_SET(flag_value, DEPRFLAG);
- bseg_put(&bitflags, &flag_value, r, c);
- bseg_put(&asp, &asp_value, r, c);
+ FLAG_SET(af.flag, INLISTFLAG);
+ FLAG_SET(af.flag, DEPRFLAG);
+ af.asp = asp_value;
+ seg_put(&aspflag, (char *)&af, r, c);
n_depr_cells++;
}
}
Modified: grass-addons/grass7/raster/r.stream.extract/load.c
===================================================================
--- grass-addons/grass7/raster/r.stream.extract/load.c 2013-04-06 20:39:39 UTC (rev 55648)
+++ grass-addons/grass7/raster/r.stream.extract/load.c 2013-04-06 20:41:23 UTC (rev 55649)
@@ -26,8 +26,9 @@
size_t ele_size, acc_size = 0;
int ele_map_type, acc_map_type = 0;
WAT_ALT *wabuf;
- char *flag_value_buf, *aspect;
+ ASP_FLAG *afbuf;
+
if (acc_fd < 0)
G_message(_("Loading elevation map..."));
else
@@ -35,7 +36,6 @@
n_search_points = n_points = 0;
- G_debug(1, "get buffers");
ele_map_type = Rast_get_map_type(ele_fd);
ele_size = Rast_cell_size(ele_map_type);
ele_buf = Rast_allocate_buf(ele_map_type);
@@ -60,9 +60,8 @@
ele_scale = 1000; /* should be enough to do the trick */
wabuf = G_malloc(ncols * sizeof(WAT_ALT));
- flag_value_buf = G_malloc(ncols * sizeof(char));
+ afbuf = G_malloc(ncols * sizeof(ASP_FLAG));
stream_id = G_malloc(ncols * sizeof(CELL));
- aspect = G_malloc(ncols * sizeof(char));
G_debug(1, "start loading %d rows, %d cols", nrows, ncols);
for (r = 0; r < nrows; r++) {
@@ -79,16 +78,16 @@
for (c = 0; c < ncols; c++) {
- flag_value_buf[c] = 0;
- aspect[c] = 0;
+ afbuf[c].flag = 0;
+ afbuf[c].asp = 0;
stream_id[c] = 0;
/* check for masked and NULL cells */
if (Rast_is_null_value(ptr, ele_map_type)) {
- FLAG_SET(flag_value_buf[c], NULLFLAG);
- FLAG_SET(flag_value_buf[c], INLISTFLAG);
- FLAG_SET(flag_value_buf[c], WORKEDFLAG);
- FLAG_SET(flag_value_buf[c], WORKED2FLAG);
+ FLAG_SET(afbuf[c].flag, NULLFLAG);
+ FLAG_SET(afbuf[c].flag, INLISTFLAG);
+ FLAG_SET(afbuf[c].flag, WORKEDFLAG);
+ FLAG_SET(afbuf[c].flag, WORKED2FLAG);
Rast_set_c_null_value(&ele_value, 1);
/* flow accumulation */
if (acc_fd >= 0) {
@@ -144,18 +143,16 @@
acc_ptr = G_incr_void_ptr(acc_ptr, acc_size);
}
seg_put_row(&watalt, (char *) wabuf, r);
- bseg_put_row(&asp, aspect, r);
+ seg_put_row(&aspflag, (char *) afbuf, r);
cseg_put_row(&stream, stream_id, r);
- bseg_put_row(&bitflags, flag_value_buf, r);
}
G_percent(nrows, nrows, 1); /* finish it */
Rast_close(ele_fd);
G_free(ele_buf);
G_free(wabuf);
- G_free(flag_value_buf);
+ G_free(afbuf);
G_free(stream_id);
- G_free(aspect);
if (acc_fd >= 0) {
Rast_close(acc_fd);
Modified: grass-addons/grass7/raster/r.stream.extract/local_proto.h
===================================================================
--- grass-addons/grass7/raster/r.stream.extract/local_proto.h 2013-04-06 20:39:39 UTC (rev 55648)
+++ grass-addons/grass7/raster/r.stream.extract/local_proto.h 2013-04-06 20:41:23 UTC (rev 55649)
@@ -29,6 +29,12 @@
DCELL wat;
};
+#define ASP_FLAG struct aspect_flag
+ASP_FLAG {
+ char asp;
+ char flag;
+};
+
struct snode
{
int r, c;
@@ -40,13 +46,14 @@
};
/* extern variables */
-extern struct snode *stream_node;
+
extern int nrows, ncols;
extern GW_LARGE_INT n_search_points, n_points, nxt_avail_pt;
extern GW_LARGE_INT heap_size;
-extern unsigned int n_stream_nodes, n_alloc_nodes;
+extern GW_LARGE_INT n_stream_nodes, n_alloc_nodes;
extern POINT *outlets;
-extern unsigned int n_outlets, n_alloc_outlets;
+extern struct snode *stream_node;
+extern GW_LARGE_INT n_outlets, n_alloc_outlets;
extern char drain[3][3];
extern char sides;
extern int c_fac;
@@ -55,9 +62,8 @@
extern SSEG search_heap;
extern SSEG astar_pts;
-extern BSEG bitflags;
-extern SSEG watalt;
-extern BSEG asp;
+extern SSEG watalt, aspflag;
+/* extern BSEG bitflags, asp; */
extern CSEG stream;
/* load.c */
@@ -72,7 +78,7 @@
/* streams.c */
int do_accum(double);
-int extract_streams(double, double, int, int);
+int extract_streams(double, double, int);
/* thin.c */
int thin_streams(void);
Modified: grass-addons/grass7/raster/r.stream.extract/main.c
===================================================================
--- grass-addons/grass7/raster/r.stream.extract/main.c 2013-04-06 20:39:39 UTC (rev 55648)
+++ grass-addons/grass7/raster/r.stream.extract/main.c 2013-04-06 20:41:23 UTC (rev 55649)
@@ -22,13 +22,13 @@
#include "local_proto.h"
/* global variables */
-struct snode *stream_node;
int nrows, ncols;
GW_LARGE_INT n_search_points, n_points, nxt_avail_pt;
GW_LARGE_INT heap_size;
-unsigned int n_stream_nodes, n_alloc_nodes;
+GW_LARGE_INT n_stream_nodes, n_alloc_nodes;
POINT *outlets;
-unsigned int n_outlets, n_alloc_outlets;
+struct snode *stream_node;
+GW_LARGE_INT n_outlets, n_alloc_outlets;
char drain[3][3] = { {7, 6, 5}, {8, 0, 4}, {1, 2, 3} };
char sides;
int c_fac;
@@ -37,9 +37,8 @@
SSEG search_heap;
SSEG astar_pts;
-BSEG bitflags;
-SSEG watalt;
-BSEG asp;
+SSEG watalt, aspflag;
+/* BSEG bitflags, asp; */
CSEG stream;
CELL *astar_order;
@@ -111,6 +110,7 @@
_("If accumulation is larger than d8cut, SFD is used instead of MFD."
" Applies only if no accumulation map is given.");
input.d8cut->required = NO;
+ input.d8cut->answer = "infinity";
input.d8cut->type = TYPE_DOUBLE;
input.mont_exp = G_define_option();
@@ -191,7 +191,7 @@
G_fatal_error(_("Threshold must be > 0 but is %f"), threshold);
/* d8cut */
- if (!input.d8cut->answer) {
+ if (strcmp(input.d8cut->answer, "infinity") == 0) {
d8cut = DBL_MAX;
}
else {
@@ -284,15 +284,14 @@
/* elevation + accumulation: * 2 */
memory_divisor = sizeof(WAT_ALT) * 2;
disk_space = sizeof(WAT_ALT);
- /* aspect: as is */
- memory_divisor += sizeof(char);
- disk_space += sizeof(char);
/* stream ids: / 2 */
memory_divisor += sizeof(int) / 2.;
disk_space += sizeof(int);
- /* flags: * 4 */
- memory_divisor += sizeof(char) * 4;
- disk_space += sizeof(char);
+
+ /* aspect and flags: * 2 */
+ memory_divisor += sizeof(ASP_FLAG) * 4;
+ disk_space += sizeof(ASP_FLAG);
+
/* astar_points: / 16 */
/* ideally only a few but large segments */
memory_divisor += sizeof(POINT) / 16.;
@@ -337,8 +336,14 @@
heap_mem += (num_open_segs * 2 - num_seg_total) * seg2kb *
sizeof(WAT_ALT) / 1024.;
cseg_open(&stream, seg_rows, seg_cols, num_open_segs / 2.);
+
+ seg_open(&aspflag, nrows, ncols, seg_rows, seg_cols, num_open_segs * 2,
+ sizeof(ASP_FLAG), 1);
+/*
bseg_open(&asp, seg_rows, seg_cols, num_open_segs);
bseg_open(&bitflags, seg_rows, seg_cols, num_open_segs * 4);
+*/
+
if (num_open_segs * 4 > num_seg_total)
heap_mem += (num_open_segs * 4 - num_seg_total) * seg2kb / 1024.;
@@ -411,8 +416,7 @@
}
/* extract streams */
- if (extract_streams
- (threshold, mont_exp, min_stream_length, acc_fd < 0) < 0)
+ if (extract_streams(threshold, mont_exp, acc_fd < 0) < 0)
G_fatal_error(_("Could not extract streams"));
seg_close(&astar_pts);
@@ -433,9 +437,8 @@
output.dir_rast->answer) < 0)
G_fatal_error(_("Could not write output maps"));
- bseg_close(&asp);
cseg_close(&stream);
- bseg_close(&bitflags);
+ seg_close(&aspflag);
exit(EXIT_SUCCESS);
}
Modified: grass-addons/grass7/raster/r.stream.extract/streams.c
===================================================================
--- grass-addons/grass7/raster/r.stream.extract/streams.c 2013-04-06 20:39:39 UTC (rev 55648)
+++ grass-addons/grass7/raster/r.stream.extract/streams.c 2013-04-06 20:41:23 UTC (rev 55649)
@@ -15,16 +15,15 @@
return result;
}
-int continue_stream(CELL stream_id, int r, int c, int r_max, int c_max,
+static int continue_stream(CELL stream_id, int r_max, int c_max,
int *stream_no)
{
- char aspect;
CELL curr_stream, stream_nbr, old_stream;
int r_nbr, c_nbr;
int asp_r[9] = { 0, -1, -1, -1, 0, 1, 1, 1, 0 };
int asp_c[9] = { 0, 1, 0, -1, -1, -1, 0, 1, 1 };
int stream_node_step = 1000;
- char this_flag_value;
+ ASP_FLAG af;
G_debug(3, "continue stream");
@@ -35,9 +34,9 @@
G_debug(3, "no confluence, just continue stream");
curr_stream = stream_id;
cseg_put(&stream, &curr_stream, r_max, c_max);
- bseg_get(&bitflags, &this_flag_value, r_max, c_max);
- FLAG_SET(this_flag_value, STREAMFLAG);
- bseg_put(&bitflags, &this_flag_value, r_max, c_max);
+ seg_get(&aspflag, (char *)&af, r_max, c_max);
+ FLAG_SET(af.flag, STREAMFLAG);
+ seg_put(&aspflag, (char *)&af, r_max, c_max);
return 0;
}
@@ -70,10 +69,9 @@
/* debug */
if (n_stream_nodes != *stream_no)
- G_warning(_("BUG: stream_no %d and n_stream_nodes %d out of sync"),
+ G_warning(_("BUG: stream_no %d and n_stream_nodes %lld out of sync"),
*stream_no, n_stream_nodes);
-
stream_node[*stream_no].n_alloc += 2;
new_size = stream_node[*stream_no].n_alloc * sizeof(int);
stream_node[*stream_no].trib =
@@ -93,17 +91,17 @@
old_stream = curr_stream;
curr_stream = *stream_no;
cseg_put(&stream, &curr_stream, r_nbr, c_nbr);
- bseg_get(&asp, &aspect, r_nbr, c_nbr);
+ seg_get(&aspflag, (char *)&af, r_nbr, c_nbr);
- while (aspect > 0) {
- r_nbr = r_nbr + asp_r[(int)aspect];
- c_nbr = c_nbr + asp_c[(int)aspect];
+ while (af.asp > 0) {
+ r_nbr = r_nbr + asp_r[(int)af.asp];
+ c_nbr = c_nbr + asp_c[(int)af.asp];
cseg_get(&stream, &stream_nbr, r_nbr, c_nbr);
if (stream_nbr != old_stream)
- aspect = -1;
+ af.asp = -1;
else {
cseg_put(&stream, &curr_stream, r_nbr, c_nbr);
- bseg_get(&asp, &aspect, r_nbr, c_nbr);
+ seg_get(&aspflag, (char *)&af, r_nbr, c_nbr);
}
}
}
@@ -147,17 +145,17 @@
double dx, dy;
int r_nbr, c_nbr, r_max, c_max, ct_dir, np_side;
int is_worked;
- char aspect;
double max_acc, prop;
int edge;
int asp_r[9] = { 0, -1, -1, -1, 0, 1, 1, 1, 0 };
int asp_c[9] = { 0, 1, 0, -1, -1, -1, 0, 1, 1 };
int nextdr[8] = { 1, -1, 0, 0, -1, 1, 1, -1 };
int nextdc[8] = { 0, 0, -1, 1, 1, -1, 1, -1 };
- unsigned int workedon, killer, count;
- char *flag_nbr, this_flag_value;
+ GW_LARGE_INT workedon, killer, count;
+ char *flag_nbr;
POINT astarpoint;
WAT_ALT wa;
+ ASP_FLAG af, af_nbr;
G_message(_("Calculating flow accumulation..."));
@@ -193,19 +191,18 @@
r = astarpoint.r;
c = astarpoint.c;
- bseg_get(&asp, &aspect, r, c);
+ seg_get(&aspflag, (char *)&af, r, c);
/* do not distribute flow along edges or out of real depressions */
- if (aspect <= 0) {
- bseg_get(&bitflags, &this_flag_value, r, c);
- FLAG_UNSET(this_flag_value, WORKEDFLAG);
- bseg_put(&bitflags, &this_flag_value, r, c);
+ if (af.asp <= 0) {
+ FLAG_UNSET(af.flag, WORKEDFLAG);
+ seg_put(&aspflag, (char *)&af, r, c);
continue;
}
- if (aspect) {
- dr = r + asp_r[abs((int)aspect)];
- dc = c + asp_c[abs((int)aspect)];
+ if (af.asp) {
+ dr = r + asp_r[abs((int)af.asp)];
+ dc = c + asp_c[abs((int)af.asp)];
}
r_max = dr;
@@ -216,9 +213,8 @@
/* WORKEDFLAG has been set during A* Search
* reversed meaning here: 0 = done, 1 = not yet done */
- bseg_get(&bitflags, &this_flag_value, r, c);
- FLAG_UNSET(this_flag_value, WORKEDFLAG);
- bseg_put(&bitflags, &this_flag_value, r, c);
+ FLAG_UNSET(af.flag, WORKEDFLAG);
+ seg_put(&aspflag, (char *)&af, r, c);
/***************************************/
/* get weights for flow distribution */
@@ -243,7 +239,8 @@
/* check that neighbour is within region */
if (r_nbr >= 0 && r_nbr < nrows && c_nbr >= 0 && c_nbr < ncols) {
- bseg_get(&bitflags, &flag_nbr[ct_dir], r_nbr, c_nbr);
+ seg_get(&aspflag, (char *)&af_nbr, r_nbr, c_nbr);
+ flag_nbr[ct_dir] = af_nbr.flag;
if ((edge = FLAG_GET(flag_nbr[ct_dir], NULLFLAG)))
break;
seg_get(&watalt, (char *)&wa, r_nbr, c_nbr);
@@ -367,7 +364,7 @@
/*
* extracts streams for threshold, accumulation is provided
*/
-int extract_streams(double threshold, double mont_exp, int min_length, int internal_acc)
+int extract_streams(double threshold, double mont_exp, int internal_acc)
{
int r, c, dr, dc;
CELL is_swale, ele_val, *ele_nbr;
@@ -378,7 +375,6 @@
double dx, dy;
int r_nbr, c_nbr, r_max, c_max, ct_dir, np_side, max_side;
int is_worked;
- char aspect;
double max_acc;
int edge, flat;
int asp_r[9] = { 0, -1, -1, -1, 0, 1, 1, 1, 0 };
@@ -394,9 +390,10 @@
GW_LARGE_INT workedon, killer, count;
int stream_no = 0, stream_node_step = 1000;
double slope, diag;
- char *flag_nbr, this_flag_value;
+ char *flag_nbr;
POINT astarpoint;
WAT_ALT wa;
+ ASP_FLAG af, af_nbr;
G_message(_("Extracting streams..."));
@@ -446,20 +443,18 @@
r = astarpoint.r;
c = astarpoint.c;
- bseg_get(&asp, &aspect, r, c);
-
- bseg_get(&bitflags, &this_flag_value, r, c);
+ seg_get(&aspflag, (char *)&af, r, c);
/* internal acc: SET, external acc: UNSET */
if (internal_acc)
- FLAG_SET(this_flag_value, WORKEDFLAG);
+ FLAG_SET(af.flag, WORKEDFLAG);
else
- FLAG_UNSET(this_flag_value, WORKEDFLAG);
- bseg_put(&bitflags, &this_flag_value, r, c);
+ FLAG_UNSET(af.flag, WORKEDFLAG);
+ seg_put(&aspflag, (char *)&af, r, c);
/* do not distribute flow along edges */
- if (aspect <= 0) {
+ if (af.asp <= 0) {
G_debug(3, "edge");
- is_swale = FLAG_GET(this_flag_value, STREAMFLAG);
+ is_swale = FLAG_GET(af.flag, STREAMFLAG);
if (is_swale) {
G_debug(2, "edge outlet");
/* add outlet point */
@@ -475,7 +470,7 @@
n_outlets++;
}
- if (aspect == 0) {
+ if (af.asp == 0) {
/* can only happen with real depressions */
if (!have_depressions)
G_fatal_error(_("Bug in stream extraction"));
@@ -484,9 +479,9 @@
continue;
}
- if (aspect) {
- dr = r + asp_r[abs((int)aspect)];
- dc = c + asp_c[abs((int)aspect)];
+ if (af.asp) {
+ dr = r + asp_r[abs((int)af.asp)];
+ dc = c + asp_c[abs((int)af.asp)];
}
else {
/* can only happen with real depressions,
@@ -529,7 +524,8 @@
if (dr == r_nbr && dc == c_nbr)
np_side = ct_dir;
- bseg_get(&bitflags, &flag_nbr[ct_dir], r_nbr, c_nbr);
+ seg_get(&aspflag, (char *)&af_nbr, r_nbr, c_nbr);
+ flag_nbr[ct_dir] = af_nbr.flag;
if ((edge = FLAG_GET(flag_nbr[ct_dir], NULLFLAG)))
break;
seg_get(&watalt, (char *)&wa, r_nbr, c_nbr);
@@ -583,7 +579,7 @@
break;
}
- is_swale = FLAG_GET(this_flag_value, STREAMFLAG);
+ is_swale = FLAG_GET(af.flag, STREAMFLAG);
/* do not continue streams along edges, these are artifacts */
if (edge) {
@@ -601,9 +597,9 @@
outlets[n_outlets].r = r;
outlets[n_outlets].c = c;
n_outlets++;
- if (aspect > 0) {
- aspect = -1 * drain[r - r_nbr + 1][c - c_nbr + 1];
- bseg_put(&asp, &aspect, r, c);
+ if (af.asp > 0) {
+ af.asp = -1 * drain[r - r_nbr + 1][c - c_nbr + 1];
+ seg_put(&aspflag, (char *)&af, r, c);
}
}
continue;
@@ -631,8 +627,8 @@
/* update aspect */
/* r_max == r && c_max == c should not happen */
if ((r_max != dr || c_max != dc) && (r_max != r || c_max != c)) {
- aspect = drain[r - r_max + 1][c - c_max + 1];
- bseg_put(&asp, &aspect, r, c);
+ af.asp = drain[r - r_max + 1][c - c_max + 1];
+ seg_put(&aspflag, (char *)&af, r, c);
}
/**********************/
@@ -660,8 +656,8 @@
G_debug(2, "start new stream");
is_swale = ++stream_no;
cseg_put(&stream, &is_swale, r, c);
- FLAG_SET(this_flag_value, STREAMFLAG);
- bseg_put(&bitflags, &this_flag_value, r, c);
+ FLAG_SET(af.flag, STREAMFLAG);
+ seg_put(&aspflag, (char *)&af, r, c);
/* add stream node */
if (stream_no >= n_alloc_nodes - 1) {
n_alloc_nodes += stream_node_step;
@@ -681,7 +677,7 @@
/* debug */
if (n_stream_nodes != stream_no)
- G_warning(_("BUG: stream_no %d and n_stream_nodes %d out of sync"),
+ G_warning(_("BUG: stream_no %d and n_stream_nodes %lld out of sync"),
stream_no, n_stream_nodes);
}
@@ -707,8 +703,7 @@
n_outlets++;
}
else {
- continue_stream(is_swale, r, c, r_max, c_max,
- &stream_no);
+ continue_stream(is_swale, r_max, c_max, &stream_no);
}
}
}
@@ -722,8 +717,8 @@
G_free(ele_nbr);
G_free(flag_nbr);
- G_debug(1, "%d outlets", n_outlets);
- G_debug(1, "%d nodes", n_stream_nodes);
+ G_debug(1, "%lld outlets", n_outlets);
+ G_debug(1, "%lld nodes", n_stream_nodes);
G_debug(1, "%d streams", stream_no);
return 1;
Modified: grass-addons/grass7/raster/r.stream.extract/thin.c
===================================================================
--- grass-addons/grass7/raster/r.stream.extract/thin.c 2013-04-06 20:39:39 UTC (rev 55648)
+++ grass-addons/grass7/raster/r.stream.extract/thin.c 2013-04-06 20:41:23 UTC (rev 55649)
@@ -10,28 +10,28 @@
CELL curr_stream, no_stream = 0;
int asp_r[9] = { 0, -1, -1, -1, 0, 1, 1, 1, 0 };
int asp_c[9] = { 0, 1, 0, -1, -1, -1, 0, 1, 1 };
- char flag_value, aspect;
+ ASP_FLAG af;
r = stream_node[stream_id].r;
c = stream_node[stream_id].c;
cseg_get(&stream, &curr_stream, r, c);
- bseg_get(&asp, &aspect, r, c);
- if (aspect > 0) {
+ seg_get(&aspflag, (char *)&af, r, c);
+ if (af.asp > 0) {
/* get downstream point */
- last_r = r + asp_r[(int)aspect];
- last_c = c + asp_c[(int)aspect];
+ last_r = r + asp_r[(int)af.asp];
+ last_c = c + asp_c[(int)af.asp];
cseg_get(&stream, &curr_stream, last_r, last_c);
if (curr_stream != stream_id)
return thinned;
/* get next downstream point */
- bseg_get(&asp, &aspect, last_r, last_c);
- while (aspect > 0) {
- r_nbr = last_r + asp_r[(int)aspect];
- c_nbr = last_c + asp_c[(int)aspect];
+ seg_get(&aspflag, (char *)&af, last_r, last_c);
+ while (af.asp > 0) {
+ r_nbr = last_r + asp_r[(int)af.asp];
+ c_nbr = last_c + asp_c[(int)af.asp];
if (r_nbr == last_r && c_nbr == last_c)
return thinned;
@@ -43,12 +43,12 @@
if (abs(r_nbr - r) < 2 && abs(c_nbr - c) < 2) {
/* eliminate last point */
cseg_put(&stream, &no_stream, last_r, last_c);
- bseg_get(&bitflags, &flag_value, last_r, last_c);
- FLAG_UNSET(flag_value, STREAMFLAG);
- bseg_put(&bitflags, &flag_value, last_r, last_c);
+ FLAG_UNSET(af.flag, STREAMFLAG);
+ seg_put(&aspflag, (char *)&af, last_r, last_c);
/* update start point */
- aspect = drain[r - r_nbr + 1][c - c_nbr + 1];
- bseg_put(&asp, &aspect, r, c);
+ seg_get(&aspflag, (char *)&af, r, c);
+ af.asp = drain[r - r_nbr + 1][c - c_nbr + 1];
+ seg_put(&aspflag, (char *)&af, r, c);
thinned = 1;
}
@@ -59,7 +59,7 @@
}
last_r = r_nbr;
last_c = c_nbr;
- bseg_get(&asp, &aspect, last_r, last_c);
+ seg_get(&aspflag, (char *)&af, last_r, last_c);
}
}
@@ -166,7 +166,7 @@
G_free(nodestack);
- G_verbose_message(_("%d of %d stream segments were thinned"), n_thinned, n_stream_nodes);
+ G_verbose_message(_("%d of %lld stream segments were thinned"), n_thinned, n_stream_nodes);
return 1;
}
More information about the grass-commit
mailing list