[GRASS-user] Classify basins as "narrow"

[Ken Mankoff]

Hi List,

On 2020-09-02 at 04:27 -07, Ken Mankoff <mankoff at gmail.com> wrote...
> I'd like to detect "narrow" features in GRASS. The attached screenshot
> shows some basins (thick) and streams (thin) and some regions
> (hatched). These regions are spurious because the basin is narrow
> here. I'd like to estimate narrowness with an algorithm.
>
> I've looked into r.grow.distance r.distance and v.distance but haven't
> been able to imagine a solution yet. Can anyone on this list suggest
> something?

I've solved this, although it takes a lot of steps, and it only computes the width of the downstream region, not everywhere. The downstream/outflow region is what I'm interested in (but didn't specify this in my initial post) so this is OK with me. By this I mean that for a "Y" shaped catchment with three narrow regions and outflow at the bottom of the stem, the current algorithm will compute the width of one of the branches and the stem, but not the other branch.

I'm sharing this in case it helps someone else, or someone sees a way to do this more efficiently (perhaps all basins at once?). Here is the algorithm:

# elevation raster is input

# compute streams, outlets, and basins

r.stream.extract elevation=elevation threshold=11 memory=16384 direction=dir stream_raster=streams stream_vector=streams

r.mapcalc "outlets = if(dir < 0, 1, null())" # outlets

r.to.vect input=outlets output=outlets type=point

r.stream.basins -m direction=dir points=outlets basins=basins memory=16384 # basins

# compute distance from edge of each basin

r.to.vect -v input=basins output=basins type=area

v.to.lines input=basins output=bounds

v.to.rast -d input=bounds output=bounds use=val val=1

r.grow.distance input=bounds distance=edge_dist metric=euclidean

# Now the algorithm can only work on one basin at a time.
# Pick a basin with narrow outlet region manually, and mask to that basin.

basin_id=1174

r.mask raster=basins maskcats=${basin_id} --o

# Invert so center line is low, but make outlet lowest.
# Edge dist is 0 at edges, so this cost surface is a hole
# Make the outlet the lowest point so it pours from there.
r.mapcalc "cost_surf = if(isnull(outlets), -edge_dist, -max(edge_dist)-100)"

# Route a stream down the center line

r.stream.extract elevation=cost_surf threshold=11 memory=16384 direction=cost_dir stream_raster=cost_streams stream_vector=cost_streams

# Find the main stream (hack = 1)

r.stream.order -m stream_rast=cost_streams direction=cost_dir elevation=cost_surf hack=hack memory=16384

r.mapcalc "hack_1 = if(hack == 1, 1, null())"

# Grow center stream to edges

r.grow.distance -m input=hack_1 distance=center_dist

# DONE

Now the value of center_dist at the basin boundary is the distance to the center line, or 1/2 the width.