[GRASS-user] Read/Write numpy array

[Johannes Radinger]

Hello,

Am 06.04.2011 um 17:36 schrieb josef.pktd at gmail.com:

> On Wed, Apr 6, 2011 at 10:58 AM, Johannes Radinger <JRadinger at gmx.at> wrote:
>> Hello Glynn,
>> Hello Josef,
>> 
>> Thank you for your answer but I still get an error, although I changed the code
>> which is now:
>>    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 error says:
>>   Density[...] = p * stats.norm.pdf(x,loc=m, scale=s1*K) +
>> (1-p) * stats.norm.pdf(x, loc=m, scale=s2*K)
>>  File "/Library/Frameworks/Python.framework/Versions/2.6/li
>> b/python2.6/site-packages/scipy/stats/distributions.py",
>> line 1104, in pdf
>>    x = arr((x-loc)*1.0/scale)
>> TypeError: unsupported operand type(s) for /:
>> 'numpy.ndarray' and 'numpy.ndarray'
>> 
>> 
>> Propably it has something to do with the stats-norm function of scipy which is working
>> with numpy array (maybe not with the garray)..or what does that error mean?
> 
> I guess is that one of your variables doesn't have a type that is
> numpy compatible (or defines division).
> My guess is that the exception shows up too late in the code, I would
> think there should be an asarray earlier in the
> scipy.stats.distributions code. (If it's not then that might be a
> bug.)
> 
> you could try np.asarray on all your variables and check which one has
> not have been converted yet
O thank you for that, I didn't realize that the variables K and s weren't converted
into a float, they were still a string... so stupid from me that I didn't check that
before, I was only thinking about the input array(x)...It is now working!

anyway just some questions about the handling of NULL values. Like in that example: 
>>    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 distance_raster contains a lot of NULL-values. It seems that they are somehow transformed
into a real value (I don't now to which one exactly) because after the Scipy-calculation (Density[...] =)
and the rewrite process I a value of 0.00101776158703459 in every cell which was a NULL-cell before.

Two questions:

1) is it possible to exclude NULL values from calculation (to speed up the process)
2) What do I've to set that NULL cells remain NULL cells also after the processing?


/johannes



> 
> stats.norm.pdf(np.asarray(x), loc=np.asarray(m), scale=np.asarray(s1*K))
> 
> I don't know what methods garray.array() inherits or has, so I don't
> know any specific answer.
> 
> Cheers,
> 
> Josef
> 
>> 
>> 
>> PS: @ You once helped me with this scipy equation I now want to use within a python-GRASS-GIS script.
>> Transforming a raster-map into an array and apply the scipy-function on it and transform it back to a raster.
>> Maybe you have some ideas what is the problem...
>> 
>> 
>> Am 05.04.2011 um 19:43 schrieb 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>
>>> _______________________________________________
>>> grass-user mailing list
>>> grass-user at lists.osgeo.org
>>> http://lists.osgeo.org/mailman/listinfo/grass-user
>> 
>>