r.watershed

Dan Epstein writes about a problem with r.watershed on a DEC 5100:

I am running GRASS4.0 on a DEC 5100. The watershed of concern is very
large - the elev map is roughly 1700 x 750 cells. In running
r.watershed for basin and stream delineation, I dont have enough memory
to run up front. As a result, the program resorts to swapping to disk.
I managed to free up enough disk space to allow it to do this.
However, after approx 3 days of running, the program is only 11% done
with Section 1b- Determining Offmap flow

SECTION 1b (of 5): Determining Offmap Flow. Percent Complete:
  11%

At this rate, the program might terminate in 2-3 weeks.

I am running using a min. outside sub-basin area of 150 km2. The map
is at roughly 90 meter resolution.

Any suggestions on how to speed this up ? Or, if I rescale the
elevation data to 270 meters ( for example), any ideas on how much
information will be lost ? At what resolution will the delination
algorithm not adequately describe the topography ? Will a 600 X 250
map run significantly faster ? Does anyone have experince with scale
and size in r.watershed ?

Dan Epstein

I had the same problem with r.watershed on an RS/6000 with a slightly
smaller mapset. The problem did not disappear when I reran it on the
mainframe with about 750M internal memory. After some
puzzling I found that there is a fault in memory allocation resulting
from am macro expansion in "ramseg.h", in init_vars.c, if I remember well. I
rewrote that macro as a function (don't forget to #define MAIN in main.c)
and all went well and in negligable time (5-10 minutes).
I enclose this function with the original macro commented out.

#define RAMSEG int
#define RAMSEGBITS 4
#define DOUBLEBITS 8 /* 2 * ramsegbits */
#define SEGLENLESS 15 /* 2 ^ ramsegbits - 1 */
/*
#define SEG_INDEX(s,r,c) (int) \
   ((((r) >> RAMSEGBITS) * (s) + (((c) >> RAMSEGBITS)) << DOUBLEBITS) \
    + (((r) & SEGLENLESS) << RAMSEGBITS) + ((c) & SEGLENLESS))
*/

#ifdef MAIN
SEG_INDEX(s,r,c)
int s,r,c;
{ return ((((r) >> RAMSEGBITS) * (s) + (((c) >> RAMSEGBITS)) << DOUBLEBITS)
    + (((r) & SEGLENLESS) << RAMSEGBITS) + ((c) & SEGLENLESS));
}
#endif

Michael Shapiro blamed the IBM prepocessor some time ago on this list.
This could be very well the case, (I don't understand what happens here
at all), but I am not able to test this out further, as this work
on watershed was just a quick hack to get some students working on a particular
problem. Since then we haven't worked with watershed analysis and the
dataset isn't available any more. Perhaps Dan Epstein's dataset
and a watershed script could be made available at moon for testing on
different architectures?

Jan Hartmann
Faculty of Environmental Sciences
University of Amsterdam