|
@@ -2,22 +2,124 @@
|
|
|
|
|
|
@brief GRASS Python scripting module (2D and 3D raster with numpy)
|
|
|
|
|
|
-Functions to use GRASS rasters with NumPy.
|
|
|
+Functions to use GRASS 2D and 3D rasters with NumPy.
|
|
|
|
|
|
Usage:
|
|
|
|
|
|
@code
|
|
|
-from grass.script import array as garray
|
|
|
-
|
|
|
- map = 'elevation'
|
|
|
- x = garray.array()
|
|
|
- x.read(map)
|
|
|
- # calculate something on array
|
|
|
- x[...] = x / 50.
|
|
|
- x.write(map + ".new")
|
|
|
+
|
|
|
+>>> import grass.script as grass
|
|
|
+>>> from grass.script import array as garray
|
|
|
+>>>
|
|
|
+>>> # We create a temporary region that is only valid in this python session
|
|
|
+... grass.use_temp_region()
|
|
|
+>>> grass.run_command("g.region", n=80, e=120, t=60, s=0, w=0, b=0, res=20, res3=20)
|
|
|
+0
|
|
|
+>>>
|
|
|
+>>> # Lets create a raster map numpy array
|
|
|
+... # based at the current region settings
|
|
|
+... map2d_1 = garray.array()
|
|
|
+>>>
|
|
|
+>>> # Write some data
|
|
|
+... for y in range(map2d_1.shape[0]):
|
|
|
+... for x in range(map2d_1.shape[1]):
|
|
|
+... map2d_1[y][x] = y + x
|
|
|
+...
|
|
|
+>>> # Lets have a look at the array
|
|
|
+... print map2d_1
|
|
|
+[[ 0. 1. 2. 3. 4. 5.]
|
|
|
+ [ 1. 2. 3. 4. 5. 6.]
|
|
|
+ [ 2. 3. 4. 5. 6. 7.]
|
|
|
+ [ 3. 4. 5. 6. 7. 8.]]
|
|
|
+>>> # This will write the numpy array as GRASS raster map
|
|
|
+... # with name map2d_1
|
|
|
+... map2d_1.write(mapname="map2d_1", overwrite=True)
|
|
|
+ 100%
|
|
|
+0
|
|
|
+>>>
|
|
|
+>>> # We create a new array and read map2d_1 to modify it
|
|
|
+... map2d_2 = garray.array()
|
|
|
+>>> # Don't do map2d_2 = map2d_1 % 3
|
|
|
+... # because: this will overwrite the internal temporary filename
|
|
|
+... map2d_2.read("map2d_1")
|
|
|
+0
|
|
|
+>>> map2d_2 %= 3
|
|
|
+>>> # Show the result
|
|
|
+... print map2d_2
|
|
|
+[[ 0. 1. 2. 0. 1. 2.]
|
|
|
+ [ 1. 2. 0. 1. 2. 0.]
|
|
|
+ [ 2. 0. 1. 2. 0. 1.]
|
|
|
+ [ 0. 1. 2. 0. 1. 2.]]
|
|
|
+>>> # Write the result as new raster map with name map2d_2
|
|
|
+... map2d_2.write(mapname="map2d_2", overwrite=True)
|
|
|
+ 100%
|
|
|
+0
|
|
|
+>>>
|
|
|
+>>> # Here we create a 3D raster map numpy array
|
|
|
+... # based in the current region settings
|
|
|
+... map3d_1 = garray.array3d()
|
|
|
+>>>
|
|
|
+>>> # Write some data
|
|
|
+... # Note: the 3D array has map[depth][row][column] order
|
|
|
+... for z in range(map3d_1.shape[0]):
|
|
|
+... for y in range(map3d_1.shape[1]):
|
|
|
+... for x in range(map3d_1.shape[2]):
|
|
|
+... map3d_1[z][y][x] = z + y + x
|
|
|
+...
|
|
|
+>>> # Lets have a look at the 3D array
|
|
|
+... print map3d_1
|
|
|
+[[[ 0. 1. 2. 3. 4. 5.]
|
|
|
+ [ 1. 2. 3. 4. 5. 6.]
|
|
|
+ [ 2. 3. 4. 5. 6. 7.]
|
|
|
+ [ 3. 4. 5. 6. 7. 8.]]
|
|
|
+
|
|
|
+ [[ 1. 2. 3. 4. 5. 6.]
|
|
|
+ [ 2. 3. 4. 5. 6. 7.]
|
|
|
+ [ 3. 4. 5. 6. 7. 8.]
|
|
|
+ [ 4. 5. 6. 7. 8. 9.]]
|
|
|
+
|
|
|
+ [[ 2. 3. 4. 5. 6. 7.]
|
|
|
+ [ 3. 4. 5. 6. 7. 8.]
|
|
|
+ [ 4. 5. 6. 7. 8. 9.]
|
|
|
+ [ 5. 6. 7. 8. 9. 10.]]]
|
|
|
+>>> # This will write the numpy array as GRASS 3D raster map
|
|
|
+... # with name map3d_1
|
|
|
+... map3d_1.write(mapname="map3d_1", overwrite=True)
|
|
|
+Loading floating point data with 8 bytes ... (6x4x3)
|
|
|
+ 100%
|
|
|
+0
|
|
|
+>>> # We create a new 3D array and read map3d_1 to modify it
|
|
|
+... map3d_2 = garray.array3d()
|
|
|
+>>> # Don't do map3d_2 = map3d_1 % 3
|
|
|
+... # because: this will overwrite the internal temporary filename
|
|
|
+... map3d_2.read("map3d_1")
|
|
|
+0
|
|
|
+>>> map3d_2 %= 3
|
|
|
+>>> # Show the result
|
|
|
+... print map3d_2
|
|
|
+[[[ 0. 1. 2. 0. 1. 2.]
|
|
|
+ [ 1. 2. 0. 1. 2. 0.]
|
|
|
+ [ 2. 0. 1. 2. 0. 1.]
|
|
|
+ [ 0. 1. 2. 0. 1. 2.]]
|
|
|
+
|
|
|
+ [[ 1. 2. 0. 1. 2. 0.]
|
|
|
+ [ 2. 0. 1. 2. 0. 1.]
|
|
|
+ [ 0. 1. 2. 0. 1. 2.]
|
|
|
+ [ 1. 2. 0. 1. 2. 0.]]
|
|
|
+
|
|
|
+ [[ 2. 0. 1. 2. 0. 1.]
|
|
|
+ [ 0. 1. 2. 0. 1. 2.]
|
|
|
+ [ 1. 2. 0. 1. 2. 0.]
|
|
|
+ [ 2. 0. 1. 2. 0. 1.]]]
|
|
|
+>>> # Write the result as new 3D raster map with name map3d_2
|
|
|
+... map3d_2.write(mapname="map3d_2", overwrite=True)
|
|
|
+Loading floating point data with 8 bytes ... (6x4x3)
|
|
|
+ 100%
|
|
|
+0
|
|
|
+
|
|
|
@endcode
|
|
|
|
|
|
-(C) 2010-2011 by Glynn Clements and the GRASS Development Team
|
|
|
+(C) 2010-2012 by Glynn Clements and the GRASS Development Team
|
|
|
This program is free software under the GNU General Public
|
|
|
License (>=v2). Read the file COPYING that comes with GRASS
|
|
|
for details.
|
|
@@ -30,6 +132,7 @@ import numpy
|
|
|
|
|
|
import core as grass
|
|
|
|
|
|
+###############################################################################
|
|
|
|
|
|
class array(numpy.memmap):
|
|
|
def __new__(cls, dtype=numpy.double):
|
|
@@ -89,7 +192,8 @@ class array(numpy.memmap):
|
|
|
output=self.filename,
|
|
|
bytes=size,
|
|
|
null=null,
|
|
|
- quiet=True)
|
|
|
+ quiet=True,
|
|
|
+ overwrite=True)
|
|
|
|
|
|
def write(self, mapname, title=None, null=None, overwrite=None):
|
|
|
"""!Write array into raster map
|
|
@@ -139,6 +243,7 @@ class array(numpy.memmap):
|
|
|
rows=reg['rows'],
|
|
|
cols=reg['cols'])
|
|
|
|
|
|
+###############################################################################
|
|
|
|
|
|
class array3d(numpy.memmap):
|
|
|
def __new__(cls, dtype=numpy.double):
|
|
@@ -163,9 +268,11 @@ class array3d(numpy.memmap):
|
|
|
shape=shape)
|
|
|
|
|
|
self.filename = filename
|
|
|
+
|
|
|
return self
|
|
|
|
|
|
def _close(self):
|
|
|
+
|
|
|
numpy.memmap._close(self)
|
|
|
if isinstance(self, array):
|
|
|
grass.try_remove(self.filename)
|
|
@@ -199,7 +306,8 @@ class array3d(numpy.memmap):
|
|
|
output=self.filename,
|
|
|
bytes=size,
|
|
|
null=null,
|
|
|
- quiet=True)
|
|
|
+ quiet=True,
|
|
|
+ overwrite=True)
|
|
|
|
|
|
def write(self, mapname, null=None, overwrite=None):
|
|
|
"""!Write array into 3D raster map
|
|
@@ -213,6 +321,7 @@ class array3d(numpy.memmap):
|
|
|
"""
|
|
|
kind = self.dtype.kind
|
|
|
size = self.dtype.itemsize
|
|
|
+ flags = None
|
|
|
|
|
|
if kind == 'f':
|
|
|
if size != 4 and size != 8:
|
|
@@ -244,3 +353,4 @@ class array3d(numpy.memmap):
|
|
|
depths=reg['depths'],
|
|
|
rows=reg['rows3'],
|
|
|
cols=reg['cols3'])
|
|
|
+
|