There is some pretty simple stuff you can do with numpy. Consider a slope with elevations ranging from 0 to 29 m and a cell size of 1m^2. What volume of earth would be removed if we removed all the material >=15 leaving behind a flat plateau behind.
In [15]: a = np.arange(5*6).reshape(5,6)
Out[16]: a
array([[ 0, 1, 2, 3, 4, 5],
[ 6, 7, 8, 9, 10, 11],
[12, 13, 14, 15, 16, 17],
[18, 19, 20, 21, 22, 23],
[24, 25, 26, 27, 28, 29]])
In [17]: a[a>15]
Out[17]: array([16, 17, 18, 19, 20, 21, 22, 23, 24, 25, 26, 27, 28, 29])
In [18]: np.sum(np.abs(15 - a[a>15]))
Out[18]: 105
So if you have a specific example, I am sure it can be done. You aren't limited to 1d or 2d arrays either. For instance, consider a tiny place (2x2) over 5 time periods. How do you calculate the per area averages over the time period, the averages per time and the overall average.
In [19]: t = np.arange(5*2*2).reshape(5,2,2)
In [20]: t
Out[20]:
array([[[ 0, 1],
[ 2, 3]],
[[ 4, 5],
[ 6, 7]],
[[ 8, 9],
[10, 11]],
[[12, 13],
[14, 15]],
[[16, 17],
[18, 19]]])
In [21]: np.average(t, axis=0)
Out[21]:
array([[ 8.000, 9.000],
[ 10.000, 11.000]])
In [22]: np.average(t, axis=1)
Out[22]:
array([[ 1.000, 2.000],
[ 5.000, 6.000],
[ 9.000, 10.000],
[ 13.000, 14.000],
[ 17.000, 18.000]])
In [23]: np.average(t)
Out[23]:
9.5
R and Numpy do share a lot... just a bit of difference in syntax.. there are a few translation tables out there converting between R, Numpy and Julia for instance. I have a large repository of 'stuff' on my blog if you are interested.
