[GRASS-user] r.stream.extract: map mis-match error

GRASS GRASS user
1

   On one project the r.stream.extract error, "Accumulation map does not
match elevation map!", was caused by a mask on the elevation map. With this
project the mask was removed, and no longer limits the DEM region.
anyway.

   I started over to ensure the elevation map was not masked. Because of the
large size, I used the NHD basin map's region to set the analytical region:

   g.region vect=basins align=dem_proj

Then changed resolution from 1 US Foot (0.305m) to 10m (30 US Feet).

   This is the same region and elevation map used in the process leading to
application of the r.stream.extract module, the accumulation map was
generated by r.watershed, and r.stream.extract has no problem running with
the unweighted accumulation map. Therefore, I assume the error arose in the
calculation of weight or applying that weight to the unweighted accumulation
map.

   If this assumption is correct, where do I look to find how to resolve the
problem?

Rich

2

On Thu, Apr 19, 2012 at 10:57 PM, Rich Shepard <rshepard@appl-ecosys.com> wrote:

On one project the r.stream.extract error, "Accumulation map does not
match elevation map!", was caused by a mask on the elevation map. With this
project the mask was removed, and no longer limits the DEM region.

It is perfectly ok to set the region to a subregion of the DEM and to
set a MASK. These settings should not be modified, however. As long as
r.watershed and r.stream.extract are used on the same region settings
and the same MASK, the map mismatch error should not happen. Just in
case, masking should be done by creating a MASK, not by modifying the
elevation values.

I started over to ensure the elevation map was not masked. Because of the
large size, I used the NHD basin map's region to set the analytical region:

   g\.region vect=basins align=dem\_proj

Then changed resolution from 1 US Foot (0.305m) to 10m (30 US Feet).

In your last posts you were using a resolution of 30m?

This is the same region and elevation map used in the process leading to
application of the r.stream.extract module, the accumulation map was
generated by r.watershed, and r.stream.extract has no problem running with
the unweighted accumulation map. Therefore, I assume the error arose in the
calculation of weight or applying that weight to the unweighted accumulation
map.

But in your last posts, applying the weight worked fine (map1.pdf,
map2.pdf, map3.pdf)? The weighed accumulation map seemed to produce
the best results. BTW, the vector output of r.stream.extract would
have been easier to see on the screenshots.

If this assumption is correct, where do I look to find how to resolve the
problem?

Try to set the region and MASK before the analysis, then don't modify
these settings during the analysis. Something like

g.region vect=basins align=dem_proj
g.region res=10m -a
r.resamp.stats in=dem_proj out=dem_10m method=average -w
v.to.rast in=basins out=MASK type=area use=val val=1
r.watershed elevation=dem_10m acc=dem_10m.acc basin=dem_10m.basin
threshold=150000
r.mapcalc "dem_10m.acc.w = dem_10m.acc * weight"
r.stream.extract elevation=dem_10m acc=dem_10m.acc.w
stream_rast=dem_10m.streams stream_vect=dem_10m_streams threshold=???

You may need to try different values for threshold with
r.stream.extract in order to get the desired detail.

HTH,

Markus M

3

On Fri, 20 Apr 2012, Markus Metz wrote:

It is perfectly ok to set the region to a subregion of the DEM and to set
a MASK. These settings should not be modified, however. As long as
r.watershed and r.stream.extract are used on the same region settings and
the same MASK, the map mismatch error should not happen. Just in case,
masking should be done by creating a MASK, not by modifying the elevation
values.

Markus,

   For the current there is no mask and the region is set to the same
parameters before the model is run. That's why I cannot find the source of
the error.

In your last posts you were using a resolution of 30m?

   No. The location units are US Feet so g.region reports them as 30. That's
approximately 10m and I use the latter as a descriptor. The original DEM has
a cell size of 1 x 1 feet (0.305 x 0.305 m) and is much too large to process
on a system with 4G RAM.

But in your last posts, applying the weight worked fine (map1.pdf,
map2.pdf, map3.pdf)? The weighed accumulation map seemed to produce the
best results. BTW, the vector output of r.stream.extract would have been
easier to see on the screenshots.

   This is true. The problem was the residual mask that caused two
comparatively large areas to not be defined as basins by r.watershed. That's
shen I removed that map and started over with the full-size DEM rather than
the one 'clipped' to the vector basin boundaries. Now I get too many small
basins but the former white areas with null elevations do contain identified
basins.

Try to set the region and MASK before the analysis, then don't modify
these settings during the analysis. Something like

g.region vect=basins align=dem_proj
g.region res=10m -a
r.resamp.stats in=dem_proj out=dem_10m method=average -w
v.to.rast in=basins out=MASK type=area use=val val=1
r.watershed elevation=dem_10m acc=dem_10m.acc basin=dem_10m.basin
threshold=150000
r.mapcalc "dem_10m.acc.w = dem_10m.acc * weight"
r.stream.extract elevation=dem_10m acc=dem_10m.acc.w
stream_rast=dem_10m.streams stream_vect=dem_10m_streams threshold=???

   This is different. I'll work on this today. The r.stream.extract,
r.stream.basins, r.stream.order, and r.basin modules are essential for the
two current big projects. So, I want to learn and better understand how to
use them to get meaningful results.

You may need to try different values for threshold with r.stream.extract
in order to get the desired detail.

   I've tried different thresholds for r.watershed but not -- yet -- for
r.stream.extract. I need to try the latter as well as looking at vector
output.

HTH,

   I'm sure it will.

Thank you very much,

Rich

4

On Fri, 20 Apr 2012, Markus Metz wrote:

g.region res=10m -a
r.resamp.stats in=dem_proj out=dem_10m method=average -w

   What is the reason to run r.resamp.stats if the region's resolution has
been changed to that cell size?

Thanks,

Rich

5

On Fri, Apr 20, 2012 at 6:22 PM, Rich Shepard <rshepard@appl-ecosys.com> wrote:

On Fri, 20 Apr 2012, Markus Metz wrote:

g.region res=10m -a

I meant
g.region res=30 -a

with 30 feet roughly equivalent to 10 meter

r.resamp.stats in=dem_proj out=dem_10m method=average -w

What is the reason to run r.resamp.stats if the region's resolution has
been changed to that cell size?

The resolution has been decreased, in this case several of the
original cells fall into one new cell. A reasonably accurate new value
can be obtained by using the average of all cells falling into a new
cell. If on the contrary resolution is changed from e.g. 30 to 10 feet
(increased), r.resamp.interp with bilinear or bicubic resampling would
be a good choice.

Markus M

6

On Fri, 20 Apr 2012, Markus Metz wrote:

The resolution has been decreased, in this case several of the original
cells fall into one new cell. A reasonably accurate new value can be
obtained by using the average of all cells falling into a new cell. If on
the contrary resolution is changed from e.g. 30 to 10 feet (increased),
r.resamp.interp with bilinear or bicubic resampling would be a good
choice.

Markus.

   Thank you for clarifying. I wasn't thinking about what happens to cells
when the resolution is made more coarse so I did not see the necessity of
running r.resamp.stats. Now I do.

Rich

7

On Fri, 20 Apr 2012, Markus Metz wrote:

Try to set the region and MASK before the analysis, then don't modify
these settings during the analysis. Something like

g.region vect=basins align=dem_proj
g.region res=10m -a
r.resamp.stats in=dem_proj out=dem_10m method=average -w
v.to.rast in=basins out=MASK type=area use=val val=1
r.watershed elevation=dem_10m acc=dem_10m.acc basin=dem_10m.basin
threshold=150000
r.mapcalc "dem_10m.acc.w = dem_10m.acc * weight"
r.stream.extract elevation=dem_10m acc=dem_10m.acc.w
stream_rast=dem_10m.streams stream_vect=dem_10m_streams threshold=???

   My only modifications to the above were to run r.param.scale after
r.watershed and r.mapcalc to calculate weight from the r.param.scale
tangential curve output maps before calling r.mapcalc to apply the weights.
Could these be problematic?

   The results are still

GRASS 6.5.svn :~/grassdata > r.stream.extract
elev=dem_10m acc=dem_10m.acc.w stream_rast=dem_10m.rstreams.w
stream_vect=dem_10m.vstreams.w thresh=1000 --o
Load input maps and get start points...
ERROR: Accumulation map does not match elevation map!

   I'm completely puzzled ... and a bit frustrated that the same error is
thrown so consistently.

Rich

8

On Fri, Apr 20, 2012 at 7:01 PM, Rich Shepard <rshepard@appl-ecosys.com> wrote:

On Fri, 20 Apr 2012, Markus Metz wrote:

Try to set the region and MASK before the analysis, then don't modify
these settings during the analysis. Something like

g.region vect=basins align=dem_proj
g.region res=10m -a
r.resamp.stats in=dem_proj out=dem_10m method=average -w
v.to.rast in=basins out=MASK type=area use=val val=1
r.watershed elevation=dem_10m acc=dem_10m.acc basin=dem_10m.basin
threshold=150000
r.mapcalc "dem_10m.acc.w = dem_10m.acc * weight"
r.stream.extract elevation=dem_10m acc=dem_10m.acc.w
stream_rast=dem_10m.streams stream_vect=dem_10m_streams threshold=???

My only modifications to the above were to run r.param.scale after
r.watershed and r.mapcalc to calculate weight from the r.param.scale
tangential curve output maps before calling r.mapcalc to apply the weights.
Could these be problematic?

I guess so because you said that r.stream.extract works with the
unweighed accumulation map. How exactly did you calculate the weight?
BTW, I checked the example script in the manual of r.stream.extract
and it works in 6.4.3 and 7.0.

Markus M

The results are still

GRASS 6.5.svn :~/grassdata > r.stream.extract
elev=dem_10m acc=dem_10m.acc.w stream_rast=dem_10m.rstreams.w
stream_vect=dem_10m.vstreams.w thresh=1000 --o
Load input maps and get start points...
ERROR: Accumulation map does not match elevation map!

I'm completely puzzled ... and a bit frustrated that the same error is
thrown so consistently.

Rich

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

9

On Sun, 22 Apr 2012, Markus Metz wrote:

I guess so because you said that r.stream.extract works with the
unweighed accumulation map. How exactly did you calculate the weight?

Markus,

r.param.scale in=dem_10m out=tancurv5 size=5 param=crosc
r.param.scale in=dem_10m out=tancurv7 size=7 param=crosc
r.param.scale in=dem_10m out=tancurv11 size=11 param=crosc
r.mapcalc "weight = if(tancurv5 < 0, -100 * tancurv5, if(tancurv7 < 0, -100
* tancurv7, if(tancurv11 < 0, -100 * tancurv11, 0.000001)))"

/* The r.mapcalc command is on a single line; wrapped here. */

r.mapcalc "dem_10m.acc.w = dem_10m.acc * weight"

BTW, I checked the example script in the manual of r.stream.extract
and it works in 6.4.3 and 7.0.

   With the Spearfish data set, correct? I could extract my data (dem_10m)
and send that for testing purposes. I'm using 6.5svn.

Rich

10

On Mon, Apr 23, 2012 at 12:17 AM, Rich Shepard <rshepard@appl-ecosys.com> wrote:

On Sun, 22 Apr 2012, Markus Metz wrote:

I guess so because you said that r.stream.extract works with the
unweighed accumulation map. How exactly did you calculate the weight?

Markus,

r.param.scale in=dem_10m out=tancurv5 size=5 param=crosc
r.param.scale in=dem_10m out=tancurv7 size=7 param=crosc
r.param.scale in=dem_10m out=tancurv11 size=11 param=crosc
r.mapcalc "weight = if(tancurv5 < 0, -100 * tancurv5, if(tancurv7 < 0, -100
* tancurv7, if(tancurv11 < 0, -100 * tancurv11, 0.000001)))"

/* The r.mapcalc command is on a single line; wrapped here. */

r.mapcalc "dem_10m.acc.w = dem_10m.acc * weight"

BTW, I checked the example script in the manual of r.stream.extract
and it works in 6.4.3 and 7.0.

With the Spearfish data set, correct?

Yes.

I could extract my data (dem_10m)
and send that for testing purposes.

That could shed some light on the mystery.

Markus M