Updating rasters with r.mapcalc.

[Alun Evans]

I have an application in which a user digitizes a small block of cells
within a much larger raster layer.  I then take this block of cells,
which is represented in its own raster layer, and use it to update the
values in the larger raster.  The value I wish to write in the larger
raster could be 0 or non-zero.

I realise that I can do this using r.mapcalc with the command:

r.mapcalc base_map=if(block,block,base_map)

where base_map is the original, larger raster layer and block is the
raster layer holding the values to be used in updating.  In order for
this command to work properly I must have the current region set to
that of the base_map.  Given the size of this map (approx.  1500x1000
cells) this operation takes quite a few seconds.  Since it lies at the
centre of a loop that is performed very often I would really like to
reduce this time.

An obvious way of achieving this would be to set the region to
to just includes the values in the block layer (I can do this
using the zoom option to g.region).  mapcalc would then only
have to consider the values that were to be updated.  However,
given the current behaviour of r.mapcalc this does not work,
since the new base_map will inherit the smaller region and
all data outside this region is lost.

So, what I need is a modified version of r.mapcalc which first
checks if the layer being created already exists.  If so then
it simply updates the raster values specified by the current
region (which can be smaller than that of the raster's region)
and leaves all values outside this region unchanged.

If you have any ideas on how I might solve this problem then
I'd be very grateful.

Thanks,

Alun Evans Ph.D.
MacMillan Bloedel Limited
Email: evans at islandnet.com