[GRASS-user] Stream extraction from pseudo-elevation map using r.watershed

[Johannes Radinger]

Dear GRASS users,

I further elaborated on this issue of r.watershed and I found what is
actually causing this strange behaviour. It seems that r.watershed cannot
handle (pseudo)elevation values data are larger than >2100000. My procedure
is working value if the values are smaller, but I get an error message (see
previous email) if my elevation map contains e.g. 2200000.

Here small working example using the NC dataset to reproduce the behaviour
of getting stream segements (and stream order) from a raster river network
without having the original elevation map available.

##########################################
# Extract river network to get a working example
g.region raster=elev_ned_30m at PERMANENT
r.watershed -m --overwrite elevation=elev_ned_30m at PERMANENT threshold=250
drainage=drainage stream=stream
r.stream.basins direction=drainage coordinates="638985, 218135"
basins=selected_basin
r.mapcalc "selected_streams = if( selected_basin  && stream,1,null())"
--overwrite


# For the actual purpose only a raster stream network is available (no DEM)
# So we need to calculate a pseudoelevation map using r.cost
r.cost --overwrite input=selected_streams output=stream_accum
start_coordinates="638985, 218135"

# Create a buffer around stream (without outlet getting buffered). See
difference if the value is 2200000
r.grow.distance --overwrite input=stream_accum value=stream_raster_nearest
r.grow --overwrite input=stream_accum output=stream_accum_grow radius=2.01
old=1 new=2100000
r.mapcalc " stream_accum_grow_start = if( stream_raster_nearest  ==
0,null(), stream_accum_grow)" --overwrite
r.grow --overwrite input=stream_accum_grow_start
output=stream_raster_buffer old=2100000 new=2100000
r.patch --overwrite input=stream_accum,stream_raster_buffer
output=stream_accum_buffered

# Extract stream with buffer does  work only if buffer value <=2100000
r.watershed -m --overwrite elevation=stream_accum_buffered threshold=3
stream=stream_test2
#############################################

Of course I agree a value of 2200000 is pretty much extreme for a raster
elevation map. However, from the error message it is difficult to actually
know what is causing the problem. Is there any special reason for that
threshold? Maybe r.watershed can have a test before to check if this
threshold is not exceeded. Most probably this issue can be solved just by
dividing the elevation by a certain value (e.g. 1000) before using it with
r.watershed.

Cheers,
Johannes


On Thu, Jul 2, 2015 at 2:18 PM, Johannes Radinger <
johannesradinger at gmail.com> wrote:

> Hi,
>
> I tried something that was working at least several months ago, but
> somehow r.watershed has changed. What I am trying to do is to extract a
> stream network from a pseudo elevation map. This pseudo elevation map is
> created by getting the distance of each stream cell from to outlet in
> upstream direction (using r.cost). In addition, the streams are buffered (3
> cells wide) with cells of a value that is much greater than the actual
> stream value. This map has several features: 1) The cell values for the
> stream are certainly decreasing in downstream direction (=elevation is
> decreasing) 2) The stream is carved into a valley (important at river
> junctions where two streams meet). So this should actually meet the
> requirements for r.watershed (at least it was working some time ago) to
> extract the stream and get the drainage directions (e.g. for other tools
> like r.stream.order etc.).
>
> Now, when I try to run r.watershed (GRASS7.1) I get following warning
> several times:
> WARNING: MFD: cumulative proportion of flow distribution not 1.0 but
>          -0.000000
>
> Although I get a drainage map as output, no output is created for the
> extracted stream network. I am wondering what changed and why and if this
> is a bug in r.watershed or if I need to change settings differently.
>
> Attached you can find the pseudo elevation map as geotiff to try it
> yourself!
>
> The command I am using is:
> r.watershed -m --o  elevation="pseudo_elevation" drainage="drainage_test"
> stream="stream_test" threshold=3
>
> My system:
> GRASS version: 7.1.svn
>
> GRASS SVN revision: 65534
>
> Build date: 2015-07-02
>
> Build platform: i686-pc-linux-gnu
>
> GDAL: 1.10.0
>
> PROJ.4: 4.8.0
>
> GEOS: 3.4.2
>
> SQLite: 3.7.9
>
> Python: 2.7.3
>
> wxPython: 2.9.4.1
>
> Any ideas or suggestions?
>
> Best regards,
> Johannes
>
-------------- next part --------------
An HTML attachment was scrubbed...
URL: <http://lists.osgeo.org/pipermail/grass-user/attachments/20150713/65bed2ca/attachment.html>