#857: add simple neighborhood functions to r.mapcalc
--------------------------+-------------------------------------------------
Reporter: dickeya | Owner: grass-dev@lists.osgeo.org
Type: enhancement | Status: new
Priority: normal | Milestone: 6.4.0
Component: Raster | Version: unspecified
Resolution: | Keywords: r.mapcalc
Platform: All | Cpu: All
--------------------------+-------------------------------------------------
Comment (by glynn):
Replying to [ticket:857 dickeya]:
> I find myself using a lot of neighbor functions in r.mapcalc, but find
the syntax to be tedious. Anything larger than a 3x3 neighborhood and
writing out "raster[-1,-1] + raster[-1,0] + raster[-1,-1] + ..." is quite
a task.
For a sum, use `r.mfilter[.fp]` or `r.neighbors ... weight=...`
> It's especially hard to do an average when nodata cells are present.
Convert nulls first in a separate pass.
> It would be great to have some built in neighborhood functions that
would handle the math, taking into account nodata cells, and working on a
rectangular or circular area. Something like:
>
{{{
neighborhood(raster, width, height, shape, function)
}}}
There area a couple of issues with this:
1. The functions which make up the bulk of r.mapcalc all operate on 1-D
arrays of size G_window_cols(). So we can't "leverage" the existing
functions for this purpose.
2. The above neighborhood() function cannot be implemented using the
existing function interface; it would need to be a new language construct,
recognised by the parser, the evaluator, and the I/O code.
If neighborhood() was a normal function, the above call would result in
its first argument being the current row of the specified raster, which is
no good for a neighborhood operation (and all of its other arguments would
also be row-sized arrays of CELL/FCELL/DCELL values; a numeric literal is
passed as a row-sized array filled with that value).
I suspect that it might be more fruitful to add new aggregates to
r.neighbors.
That still leaves a gap between what can be implemented with a combination
of r.mapcalc + r.neighbors + r.mapcalc and what could be implemented with
a more general language. E.g. there isn't any combination which will allow
you to calculate:
{{{
output[r,c] = sum<i,j>(input[r+i,c+j] <op> kernel[i,j])
}}}
for any <op> except multiplication.
Implementing that specific case in r.mapcalc (or r.neighbors) would be a
lot of work, but not entirely out of the question. My main concern is that
it would be the thin end of the wedge; either we draw a line in the sand,
or we end up with Perl.
IMHO, if you can't do what you want with r.mapcalc + r.neighbors, you'll
need a system with enough (virtual) memory to use Matlab/Octave/NumPy.
--
Ticket URL: <http://trac.osgeo.org/grass/ticket/857#comment:1>
GRASS GIS <http://grass.osgeo.org>