[GRASS-user] Read/Write numpy array

[Glynn Clements]

Johannes Radinger wrote:

> I try to do a mapcalculation with scipy (do to limited functionality of 
> usual mapcalc functions.
> 
> Therefore I want to read a map into numpy-array, do the calculation
> on it and write it back to a GRASS-map. I tried to follow the
> instructions in the wiki
> (http://grass.osgeo.org/wiki/GRASS_and_Python#Interfacing_with_NumPy)
> but I got stucked when it comes to the writing-process.
> Thats my code so far:
> 
>     x = garray.array()
>     x.read("distance_raster")
>     Density = garray.array()
>     Density = p * stats.norm.pdf(x,loc=m, scale=s1*K) + (1-p) * stats.norm.pdf(x, loc=m, scale=s2*K)
>     Density.write("Density")
> 
> The stats.norm.pdf is from the stats-scipy module. I get following
> error when I try to process the script:
> 
>    Density.write("Density")
> AttributeError: 'numpy.ndarray' object has no attribute
> 'write'
> 
> I don't know excactly what how I have to use the write command.

The write() method is part of garray.array, which is a subclass of
numpy.array (specifically, of numpy.memmap).

This line:

>     Density = garray.array()

creates an instance of garray.array and associates the name "Density"
with it.

This line:

     Density = p * stats.norm.pdf(x,loc=m, scale=s1*K) + (1-p) * stats.norm.pdf(x, loc=m, scale=s2*K)

creates a normal numpy.array and associates the name "Density" with
it, replacing the original association. Consequently, Density.write()
fails because Density now refers to a numpy.array and not a
garray.array().

Rather than making the name Density refer to the resulting
numpy.array, you have to assign the resulting numpy.array *into* the
existing garray.array to which the name Density already refers, i.e.:

     Density[...] = p * stats.norm.pdf(x,loc=m, scale=s1*K) + (1-p) * stats.norm.pdf(x, loc=m, scale=s2*K)

Note that the "..." (ellipsis) is literal. The above is equivalent to:

     Density[:,:] = p * stats.norm.pdf(x,loc=m, scale=s1*K) + (1-p) * stats.norm.pdf(x, loc=m, scale=s2*K)

> I am working on 6.5 and I read there were some problems with
> numpy-arrays. Are the solved now?

I believe that the issue with grass.script.array requiring the 7.0
version of r.{in,out}.bin has been solved by backporting the changes
into 6.5.

-- 
Glynn Clements <glynn at gclements.plus.com>