[GRASS-dev] [GRASS GIS] #1924: r.watershed - empty stream segment map

GRASS GIS GRASS GIS developer
1

#1924: r.watershed - empty stream segment map
-----------------------+----------------------------------------------------
Reporter: jradinger | Owner: grass-dev@…
     Type: defect | Status: new
Priority: normal | Milestone:
Component: Raster | Version: svn-trunk
Keywords: | Platform: Linux
      Cpu: x86-32 |
-----------------------+----------------------------------------------------
I am running GRASS 7 SVN (Rev55666, last revision of r.watershed) and face
following two problems with r.watershed:

1) When running with treshhold = 2 I don't get any errors (one warning,
see below), a flow direction map is produced correctly but the '''stream-
segments map is empty'''.
This was (and is still if nothing has changed since one week) in GRASS
6.5.
Here the debug output:

{{{
r.watershed --overwrite elevation=distance_raster_tmp_27369@Fidimo_test
threshold=2 drainage=flow_direction_test stream=segments_test
D1/1: Mode: All in RAM
SECTION 1a (of 5): Initiating Memory.
D1/1: G_find_raster2(): name=distance_raster_tmp_27369
mapset=Fidimo_test
D1/1: G_find_raster2(): name=distance_raster_tmp_27369
mapset=Fidimo_test
D1/1: G_find_raster2(): name=distance_raster_tmp_27369
mapset=Fidimo_test
D1/1: G_find_raster2(): name=distance_raster_tmp_27369
mapset=Fidimo_test
D1/1: G_find_raster2(): name=distance_raster_tmp_27369
mapset=Fidimo_test
D1/1: G_find_raster2(): name=distance_raster_tmp_27369
mapset=Fidimo_test
D1/1: G_find_raster2(): name=distance_raster_tmp_27369
mapset=Fidimo_test
D1/1: G_find_raster2(): name=distance_raster_tmp_27369
mapset=Fidimo_test
D1/1: G_find_raster(): name=MASK mapset=Fidimo_test
SECTION 1b (of 5): Determining Offmap Flow.
SECTION 2: A* Search.
SECTION 3a: Accumulating Surface Flow with MFD.
D1/1: MFD convergence factor set to 5.
D1/1: ns contour: 41.4214
D1/1: ew contour: 41.4214
D1/1: diag contour: 41.4214
SECTION 3b: Adjusting drainage directions.
SECTION 4: Watershed determination.
SECTION 5: Closing Maps.
D1/1: 0 basins created
WARNING: No basins were created. Verify threshold and region settings.
D1/1: Too many subbasins to reasonably check for color
brightness
D1/1: close segments_test compressed
D1/1: G_find_raster2(): name=segments_test
mapset=Fidimo_test
D1/1: G_find_raster2(): name=segments_test
mapset=Fidimo_test
D1/1: close flow_direction_test compressed
D1/1: G_find_raster2(): name=flow_direction_test
mapset=Fidimo_test
D1/1: G_find_raster2(): name=flow_direction_test
mapset=Fidimo_test
(Wed Apr 10 09:55:17 2013) Command finished (0 sec)
}}}

What does the warning mean? My region settings:

{{{
g.region -p
projection: 1 (UTM)
zone: 13
datum: nad27
ellipsoid: clark66
north: 4928000
south: 4914000
west: 590000
east: 609000
nsres: 100
ewres: 100
rows: 140
cols: 190
cells: 26600
}}}

2) When running with treshold = 0 (default), I get following Error:

{{{
r.watershed --overwrite elevation=distance_raster_tmp_27369@Fidimo_test
threshold=0 drainage=flow_direction_test stream=segments_test
D1/1: Mode: All in RAM
ERROR: USAGE for basin delineation:
/usr/local/grass-7.0.svn/etc/r.watershed/ram -4 elevation=elevation_map
threshold=swale_threshold [flow=overland_flow_map]
[drainage=drain_direction_map] [depression=depression_map]
[accumulation=accumulation_map] [basin=watershed_basin_map]
[stream=stream_segment_map]

USAGE for slope length determination:
/usr/local/grass-7.0.svn/etc/r.watershed/ram [-4] elevation=elevation_map
threshold=swale_threshold [drainage=drain_direction_map]
[depression=depression_map] [accumulation=accumulation_map]
[max_slope_length=max_slope_length] [blocking=overland_blocking_map]
[slope_steepness=slope_steepness_map] length_slope=length_slope_map
[disturbed_land=rill_erosion_map] [slope_deposition=slope_deposition value
or map]USAGE for ARMSED FILE creation:
/usr/local/grass-7.0.svn/etc/r.watershed/ram [-4] elevation=elevation_map
threshold=swale_threshold [flow=overland_flow_map]
[drainage=drain_direction_map] [depression=depression_map]
[accumulation=accumulation_map] [basin=watershed_basin_map]
[stream=stream_segment_map] [half_basin=half_basin_map]
ar=ARMSED_file_name

WARNING: Subprocess failed with exit code 1
}}}

I am especially interested in a solution of problem 1 to get also a
segment map as output of r.watershed.

--
Ticket URL: <http://trac.osgeo.org/grass/ticket/1924&gt;
GRASS GIS <http://grass.osgeo.org>

2

#1924: r.watershed - empty stream segment map
-----------------------+----------------------------------------------------
Reporter: jradinger | Owner: grass-dev@…
     Type: defect | Status: new
Priority: normal | Milestone:
Component: Raster | Version: svn-trunk
Keywords: | Platform: Linux
      Cpu: x86-32 |
-----------------------+----------------------------------------------------

Comment(by mmetz):

Does the region match the input raster? Can you post the complete output
of

{{{
r.info distance_raster_tmp_27369@Fidimo_test
}}}

?

The threshold must be > 0. With a threshold of 2, the stream network would
be very dense.

Markus M

--
Ticket URL: <http://trac.osgeo.org/grass/ticket/1924#comment:1&gt;
GRASS GIS <http://grass.osgeo.org>

3

#1924: r.watershed - empty stream segment map
-----------------------+----------------------------------------------------
Reporter: jradinger | Owner: grass-dev@…
     Type: defect | Status: new
Priority: normal | Milestone:
Component: Raster | Version: svn-trunk
Keywords: | Platform: Linux
      Cpu: x86-32 |
-----------------------+----------------------------------------------------

Comment(by jradinger):

I just checked an the region matches:

{{{
| Type of Map: raster Number of Categories: 0
|
  | Data Type: DCELL
|
  | Rows: 140
|
  | Columns: 190
|
  | Total Cells: 26600
|
  | Projection: UTM (zone 13)
|
  | N: 4928000 S: 4914000 Res: 100
|
  | E: 609000 W: 590000 Res: 100
|
  | Range of data: min = 0 max = 16433.8095116624
|
  |
|
  | Data Description:
|
  | generated by r.cost
}}}

Two more important things:
1) The elevation input map is accutally already a thinned river raster.
The raster
values represent a pseude elevation which was generated with r.cost from
the source.
Thus the value (elevation) is increasing upstream. r.watershed is just
used for getting
the flowdirection and segmenting the river raster. (if needed I can
provide the elevation-rcost
raster)

2) I tried r.watershed with treshold=1 and that is working perfectly

And maybe a warning/error message should be printed when threshold = 0,
instead of the (for me confusing error)

--
Ticket URL: <http://trac.osgeo.org/grass/ticket/1924#comment:2&gt;
GRASS GIS <http://grass.osgeo.org>

4

#1924: r.watershed - empty stream segment map
-----------------------+----------------------------------------------------
Reporter: jradinger | Owner: grass-dev@…
     Type: defect | Status: new
Priority: normal | Milestone:
Component: Raster | Version: svn-trunk
Keywords: | Platform: Linux
      Cpu: x86-32 |
-----------------------+----------------------------------------------------

Comment(by mmetz):

Replying to [comment:2 jradinger]:
>
> Two more important things:
> 1) The elevation input map is accutally already a thinned river raster.
The raster
> values represent a pseude elevation which was generated with r.cost
from the source.
> Thus the value (elevation) is increasing upstream. r.watershed is just
used for getting
> the flowdirection and segmenting the river raster. (if needed I can
provide the elevation-rcost
> raster)

r.watershed expects a surface as input. Using a thinned river raster as
input should give weird results: all non-NULL cells face at least one NULL
cell treated as unknown elevation, thus flow direction can not be
unambiguously determined. You can try to place a buffer around the rivers
and fill the buffer with pseudo-elevation values larger than the largest
cost value. For example, create a buffer using 300 meter (3 cells) as
distance, set a pseudo elevation value for the buffer zone, patch the
thinned river raster with the pseudo elevation, run r.watershed with
threshold=3. Works for me.

> And maybe a warning/error message should be printed when threshold = 0,
instead of the (for me confusing error)

Makes sense.

Markus M

--
Ticket URL: <http://trac.osgeo.org/grass/ticket/1924#comment:3&gt;
GRASS GIS <http://grass.osgeo.org>

5

Hi Markus,

Hi others,

I am coming back to the topic of running r.watershed on a rasterized river network.

As recommended I buffered now my river raster. This river raster is a thinned distance raster
with increasing values from the outlet to the branches (calculated via r.cost) and as

I said it is buffered now with a very large value (2 x max of the thinned distance raster).

So actually this should work as an input map to r.watershed. And it does more or
less perfectly: There’s just one issue concerning the outlet-stream segment.

Somehow this (in my case small segment, in the south-west corner) is missing in the stream-map output of r.watershed.
I am working on GRASS6.5. Attached there is a zipped GTIFF with the “pseudo-elevation” map

(this thinned and buffered distance raster) I am using as input to r.watershed.

The command I am using is:
r.watershed -f --overwrite elevation=buffer_river_raster stream=test_segments threshold=3

The threshold has been set to 3 as this is the number of raster cells that refer to the buffer (3 x resolution)

Can anyone reproduce that the outlet small river segment is missing in the produced “test_segments”??

What is causing that?

/Johannes

(attachments)

buffered_distance_raster.zip (99.5 KB)

···

On Wed, Apr 10, 2013 at 12:37 PM, GRASS GIS <trac@osgeo.org> wrote:

#1924: r.watershed - empty stream segment map
-----------------------±---------------------------------------------------
Reporter: jradinger | Owner: grass-dev@…
Type: defect | Status: new
Priority: normal | Milestone:
Component: Raster | Version: svn-trunk
Keywords: | Platform: Linux
Cpu: x86-32 |
-----------------------±---------------------------------------------------

Comment(by mmetz):

Replying to [comment:2 jradinger]:

Two more important things:

  1. The elevation input map is accutally already a thinned river raster.
    The raster
    values represent a pseude elevation which was generated with r.cost
    from the source.
    Thus the value (elevation) is increasing upstream. r.watershed is just
    used for getting
    the flowdirection and segmenting the river raster. (if needed I can
    provide the elevation-rcost
    raster)

r.watershed expects a surface as input. Using a thinned river raster as
input should give weird results: all non-NULL cells face at least one NULL
cell treated as unknown elevation, thus flow direction can not be
unambiguously determined. You can try to place a buffer around the rivers
and fill the buffer with pseudo-elevation values larger than the largest
cost value. For example, create a buffer using 300 meter (3 cells) as
distance, set a pseudo elevation value for the buffer zone, patch the
thinned river raster with the pseudo elevation, run r.watershed with
threshold=3. Works for me.

And maybe a warning/error message should be printed when threshold = 0,
instead of the (for me confusing error)

Makes sense.

Markus M

–
Ticket URL: <http://trac.osgeo.org/grass/ticket/1924#comment:3>
GRASS GIS <http://grass.osgeo.org>

6

On Mon, Jun 3, 2013 at 3:30 PM, Johannes Radinger
<johannesradinger@gmail.com> wrote:

Hi Markus,
Hi others,

I am coming back to the topic of running r.watershed on a rasterized river
network.

As recommended I buffered now my river raster. This river raster is a
thinned distance raster
with increasing values from the outlet to the branches (calculated via
r.cost) and as
I said it is buffered now with a very large value (2 x max of the thinned
distance raster).

So actually this should work as an input map to r.watershed. And it does
more or
less perfectly: There's just one issue concerning the outlet-stream segment.
Somehow this (in my case small segment, in the south-west corner) is missing
in the stream-map output of r.watershed.
I am working on GRASS6.5. Attached there is a zipped GTIFF with the
"pseudo-elevation" map
(this thinned and buffered distance raster) I am using as input to
r.watershed.

There is another error in the north-eastern part where a stream
segment extents into the buffer area. The reason for these artefacts
is that r.watershed does not know where the outlet is because all
river segments are inside the buffer. In this case, the outlet is the
starting point of r.cost which has in the river raster the value 0.
You could produce the desired output with

# buffer
r.grow.distance input=river_raster value=river_raster_nearest

# get value of nearest river cell
r.grow input=river_raster output=river_raster_grow radius=2.01 old=1 new=100000

# remove buffer for start point
r.mapcalc "river_raster_grow_start = if(river_raster_nearest == 0,
null(), river_raster_grow)"

# grow by one cell to make sure that the start point is the only cell
of the river bordering a NULL cell
r.grow input=river_raster_grow_start output=river_raster_buffer
radius=1.01 old=100000 new=100000

# path river raster with buffer
r.patch in=river_raster,river_raster_buffer out=buffered_distance_raster2

# get stream segments
r.watershed -s --overwrite elevation=buffered_distance_raster2
stream=test_segments threshold=3

HTH,

Markus

The command I am using is:
r.watershed -f --overwrite elevation=buffer_river_raster
stream=test_segments threshold=3
The threshold has been set to 3 as this is the number of raster cells that
refer to the buffer (3 x resolution)

Can anyone reproduce that the outlet small river segment is missing in the
produced "test_segments"??
What is causing that?

/Johannes

On Wed, Apr 10, 2013 at 12:37 PM, GRASS GIS <trac@osgeo.org> wrote:

#1924: r.watershed - empty stream segment map

-----------------------+----------------------------------------------------
Reporter: jradinger | Owner: grass-dev@…
     Type: defect | Status: new
Priority: normal | Milestone:
Component: Raster | Version: svn-trunk
Keywords: | Platform: Linux
      Cpu: x86-32 |

-----------------------+----------------------------------------------------

Comment(by mmetz):

Replying to [comment:2 jradinger]:
>
> Two more important things:
> 1) The elevation input map is accutally already a thinned river raster.
The raster
> values represent a pseude elevation which was generated with r.cost
from the source.
> Thus the value (elevation) is increasing upstream. r.watershed is just
used for getting
> the flowdirection and segmenting the river raster. (if needed I can
provide the elevation-rcost
> raster)

r.watershed expects a surface as input. Using a thinned river raster as
input should give weird results: all non-NULL cells face at least one
NULL
cell treated as unknown elevation, thus flow direction can not be
unambiguously determined. You can try to place a buffer around the rivers
and fill the buffer with pseudo-elevation values larger than the largest
cost value. For example, create a buffer using 300 meter (3 cells) as
distance, set a pseudo elevation value for the buffer zone, patch the
thinned river raster with the pseudo elevation, run r.watershed with
threshold=3. Works for me.

> And maybe a warning/error message should be printed when threshold = 0,
instead of the (for me confusing error)

Makes sense.

Markus M

--
Ticket URL: <http://trac.osgeo.org/grass/ticket/1924#comment:3&gt;
GRASS GIS <http://grass.osgeo.org>

_______________________________________________
grass-user mailing list
grass-user@lists.osgeo.org
http://lists.osgeo.org/mailman/listinfo/grass-user

7

Thank you Markus,

thats really a profound answer perfectly solving the problem. I’ve tested it, and
it is smoothly working in GRASS65. Just the variables for r.patch in python
are input and output instead of in and out…

cheers,
/johannes

···

On Mon, Jun 3, 2013 at 9:49 PM, Markus Metz <markus.metz.giswork@gmail.com> wrote:

On Mon, Jun 3, 2013 at 3:30 PM, Johannes Radinger
<johannesradinger@gmail.com> wrote:

Hi Markus,
Hi others,

I am coming back to the topic of running r.watershed on a rasterized river
network.

As recommended I buffered now my river raster. This river raster is a
thinned distance raster
with increasing values from the outlet to the branches (calculated via
r.cost) and as
I said it is buffered now with a very large value (2 x max of the thinned
distance raster).

So actually this should work as an input map to r.watershed. And it does
more or
less perfectly: There’s just one issue concerning the outlet-stream segment.
Somehow this (in my case small segment, in the south-west corner) is missing
in the stream-map output of r.watershed.
I am working on GRASS6.5. Attached there is a zipped GTIFF with the
“pseudo-elevation” map
(this thinned and buffered distance raster) I am using as input to
r.watershed.

There is another error in the north-eastern part where a stream
segment extents into the buffer area. The reason for these artefacts
is that r.watershed does not know where the outlet is because all
river segments are inside the buffer. In this case, the outlet is the
starting point of r.cost which has in the river raster the value 0.
You could produce the desired output with

buffer

r.grow.distance input=river_raster value=river_raster_nearest

get value of nearest river cell

r.grow input=river_raster output=river_raster_grow radius=2.01 old=1 new=100000

remove buffer for start point

r.mapcalc “river_raster_grow_start = if(river_raster_nearest == 0,
null(), river_raster_grow)”

grow by one cell to make sure that the start point is the only cell

of the river bordering a NULL cell
r.grow input=river_raster_grow_start output=river_raster_buffer
radius=1.01 old=100000 new=100000

path river raster with buffer

r.patch in=river_raster,river_raster_buffer out=buffered_distance_raster2

get stream segments

r.watershed -s --overwrite elevation=buffered_distance_raster2
stream=test_segments threshold=3

HTH,

Markus

The command I am using is:
r.watershed -f --overwrite elevation=buffer_river_raster
stream=test_segments threshold=3
The threshold has been set to 3 as this is the number of raster cells that
refer to the buffer (3 x resolution)

Can anyone reproduce that the outlet small river segment is missing in the
produced “test_segments”??
What is causing that?

/Johannes

On Wed, Apr 10, 2013 at 12:37 PM, GRASS GIS <trac@osgeo.org> wrote:

#1924: r.watershed - empty stream segment map

-----------------------±---------------------------------------------------
Reporter: jradinger | Owner: grass-dev@…
Type: defect | Status: new
Priority: normal | Milestone:
Component: Raster | Version: svn-trunk
Keywords: | Platform: Linux
Cpu: x86-32 |

-----------------------±---------------------------------------------------

Comment(by mmetz):

Replying to [comment:2 jradinger]:

Two more important things:

  1. The elevation input map is accutally already a thinned river raster.
    The raster
    values represent a pseude elevation which was generated with r.cost
    from the source.
    Thus the value (elevation) is increasing upstream. r.watershed is just
    used for getting
    the flowdirection and segmenting the river raster. (if needed I can
    provide the elevation-rcost
    raster)

r.watershed expects a surface as input. Using a thinned river raster as
input should give weird results: all non-NULL cells face at least one
NULL
cell treated as unknown elevation, thus flow direction can not be
unambiguously determined. You can try to place a buffer around the rivers
and fill the buffer with pseudo-elevation values larger than the largest
cost value. For example, create a buffer using 300 meter (3 cells) as
distance, set a pseudo elevation value for the buffer zone, patch the
thinned river raster with the pseudo elevation, run r.watershed with
threshold=3. Works for me.

And maybe a warning/error message should be printed when threshold = 0,
instead of the (for me confusing error)

Makes sense.

Markus M

–
Ticket URL: <http://trac.osgeo.org/grass/ticket/1924#comment:3>
GRASS GIS <http://grass.osgeo.org>

grass-user mailing list
grass-user@lists.osgeo.org
http://lists.osgeo.org/mailman/listinfo/grass-user

8

On Tue, Jun 4, 2013 at 11:00 AM, Johannes Radinger
<johannesradinger@gmail.com> wrote:

Thank you Markus,

thats really a profound answer perfectly solving the problem. I've tested
it, and
it is smoothly working in GRASS65. Just the variables for r.patch in python
are input and output instead of in and out....

Only shell supports option abbreviation. Sometimes shell is easier
(and faster to type).

The procedure I suggested will not work in all situations. Sometimes
it might be necessary to create a buffer, and then set the start point
(cumulative costs = 0) to NULL instead of modifying the buffer. This
will create segments leading up to the start point, but excluding the
start point.

Markus M

cheers,
/johannes

On Mon, Jun 3, 2013 at 9:49 PM, Markus Metz <markus.metz.giswork@gmail.com>
wrote:

On Mon, Jun 3, 2013 at 3:30 PM, Johannes Radinger
<johannesradinger@gmail.com> wrote:
> Hi Markus,
> Hi others,
>
> I am coming back to the topic of running r.watershed on a rasterized
> river
> network.
>
> As recommended I buffered now my river raster. This river raster is a
> thinned distance raster
> with increasing values from the outlet to the branches (calculated via
> r.cost) and as
> I said it is buffered now with a very large value (2 x max of the
> thinned
> distance raster).
>
> So actually this should work as an input map to r.watershed. And it does
> more or
> less perfectly: There's just one issue concerning the outlet-stream
> segment.
> Somehow this (in my case small segment, in the south-west corner) is
> missing
> in the stream-map output of r.watershed.
> I am working on GRASS6.5. Attached there is a zipped GTIFF with the
> "pseudo-elevation" map
> (this thinned and buffered distance raster) I am using as input to
> r.watershed.

There is another error in the north-eastern part where a stream
segment extents into the buffer area. The reason for these artefacts
is that r.watershed does not know where the outlet is because all
river segments are inside the buffer. In this case, the outlet is the
starting point of r.cost which has in the river raster the value 0.
You could produce the desired output with

# buffer
r.grow.distance input=river_raster value=river_raster_nearest

# get value of nearest river cell
r.grow input=river_raster output=river_raster_grow radius=2.01 old=1
new=100000

# remove buffer for start point
r.mapcalc "river_raster_grow_start = if(river_raster_nearest == 0,
null(), river_raster_grow)"

# grow by one cell to make sure that the start point is the only cell
of the river bordering a NULL cell
r.grow input=river_raster_grow_start output=river_raster_buffer
radius=1.01 old=100000 new=100000

# path river raster with buffer
r.patch in=river_raster,river_raster_buffer out=buffered_distance_raster2

# get stream segments
r.watershed -s --overwrite elevation=buffered_distance_raster2
stream=test_segments threshold=3

HTH,

Markus

>
> The command I am using is:
> r.watershed -f --overwrite elevation=buffer_river_raster
> stream=test_segments threshold=3
> The threshold has been set to 3 as this is the number of raster cells
> that
> refer to the buffer (3 x resolution)
>
> Can anyone reproduce that the outlet small river segment is missing in
> the
> produced "test_segments"??
> What is causing that?
>
> /Johannes
>
>
>
>
>
> On Wed, Apr 10, 2013 at 12:37 PM, GRASS GIS <trac@osgeo.org> wrote:
>>
>> #1924: r.watershed - empty stream segment map
>>
>>
>> -----------------------+----------------------------------------------------
>> Reporter: jradinger | Owner: grass-dev@…
>> Type: defect | Status: new
>> Priority: normal | Milestone:
>> Component: Raster | Version: svn-trunk
>> Keywords: | Platform: Linux
>> Cpu: x86-32 |
>>
>>
>> -----------------------+----------------------------------------------------
>>
>> Comment(by mmetz):
>>
>> Replying to [comment:2 jradinger]:
>> >
>> > Two more important things:
>> > 1) The elevation input map is accutally already a thinned river
>> raster.
>> The raster
>> > values represent a pseude elevation which was generated with r.cost
>> from the source.
>> > Thus the value (elevation) is increasing upstream. r.watershed is
>> just
>> used for getting
>> > the flowdirection and segmenting the river raster. (if needed I can
>> provide the elevation-rcost
>> > raster)
>>
>> r.watershed expects a surface as input. Using a thinned river raster
>> as
>> input should give weird results: all non-NULL cells face at least one
>> NULL
>> cell treated as unknown elevation, thus flow direction can not be
>> unambiguously determined. You can try to place a buffer around the
>> rivers
>> and fill the buffer with pseudo-elevation values larger than the
>> largest
>> cost value. For example, create a buffer using 300 meter (3 cells) as
>> distance, set a pseudo elevation value for the buffer zone, patch the
>> thinned river raster with the pseudo elevation, run r.watershed with
>> threshold=3. Works for me.
>>
>> > And maybe a warning/error message should be printed when threshold =
>> 0,
>> instead of the (for me confusing error)
>>
>> Makes sense.
>>
>> Markus M
>>
>> --
>> Ticket URL: <http://trac.osgeo.org/grass/ticket/1924#comment:3&gt;
>> GRASS GIS <http://grass.osgeo.org>
>>
>
>
> _______________________________________________
> grass-user mailing list
> grass-user@lists.osgeo.org
> http://lists.osgeo.org/mailman/listinfo/grass-user
>