[GRASS5] Re: [GRASSLIST:1551] Re: Use of neighbourhood operator in r.mapcalc

[Hamish]

> > r.mapcalc buffmap = if(oldmap, oldmap, \
> >  if(oldmap[0,-1], oldmap[0,-1] \
> >   if(oldmap[0,1], oldmap[0,1] \
> >    if(oldmap[-1,0], oldmap[-1,0] \
> >     if(oldmap[1,0], oldmap[1,0])))))
>
> Aside: this operation seems to be sufficiently popular that we should
> consider adding an option to either r.buffer or r.grow to support it. 
> Using r.mapcalc quickly becomes impractical if you want a buffer
> larger than one cell.

[enter r.grow2]

Great, now I can get rid of my crappy script.

My use:
I've used this with a series of "region codes" at various site locations
which 'grow' outwards in a loop until the entire non-null region is
filled.

General method:
s.to.rast taking category number for raster value.
r.mapcalc base map so area of interest is 0 and no-go is NULL.
r.patch sites,basemap
r.grow_endless_loop
^C, check, r.grow_endless_loop [or {while find(==0)} might work..]

This allows me to 'grow' around corners; s.surf.rst+r.reclass has no
concept of impenetrable barriers, eg data from two parallel lakes
pollute each other.


One moderately important caveat:
For correct results I found I had to alternate the sign of the [-1,0]'s
(search order) after each iteration in order for it to make it around
certain corners and/or fill in some pixels which only have data on one
adjoining side. I cycled the search order but maybe randomly is more
appropriate.



thanks, this will be useful.
Hamish