You can use the Raster Cell Iterator (RCI) to perform custom raster analysis. It allows you to visit individual raster cell locations, and provides iteration control to query and modify raster cell values, all in a Python environment.
Basic operations using Raster Cell Iterator
In this section, you will gain a conceptual understanding of the basic operations that can be performed using RCI. You will learn how to: create an implicit RCI, query and assign values to raster cells, create an empty raster with metadata, index a raster object, and handle NoData. Each section includes code examples that demonstrate how to use RCI from a Python window. It is recommended that you read this document from top-down to understand these code examples.
Implicit Raster Cell Iterator
An RCI can be invoked on a raster object implicitly to iterate through the row and column indices of a raster. The code example shown below creates a raster object called myRas from an existing raster dataset "myras". An implicit raster iterator is defined on the raster object that enumerates row and column index pairs across the raster in a loop. Within the iterator, row index i, column index j, and cell value myRas[i,j] are printed at each cell location. This is an example of querying a cell value at a given cell location using index notation.
from arcpy.sa import *
myRas = Raster("myras")
for i,j in myRas:
print(i, j, myRas[i,j])
Assign values to cells
Using RCI, you can assign values to raster cells and thereby modify the raster dataset. By default, raster objects created from an existing raster dataset have a readOnly property that is set to True. To modify raster datasets, the readOnly property should first be set to False. Once the raster cells have been modified, the save() method should be called on the raster object to persist these changes. For more information on the readOnly property and saving rasters, see the Raster object help topic.
myRas.readOnly = False
for i,j in myRas:
if myRas[i,j] == 0:
myRas[i,j] = 128
myRas.save()
Creating empty rasters
The getRasterInfo() method allows you to easily create an empty raster object and copy over metadata from an existing raster object. The getRasterInfo() method is called on an existing raster object, which returns a rasterInfo() object for this raster. The rasterInfo() object is then passed as an input parameter to the Raster() class to instantiate an empty raster. The empty raster object has the same metadata as the raster object from which it is created, but has NoData in all raster cells. Raster metadata consists of bandCount, cell size, extent, spatialReference, pixelType, and noDataValues, among other properties. When a raster object is created using the getRasterInfo() method, the readOnly property is set to False by default so that values can be assigned to this empty raster object.
outRas = Raster(myRas.getRasterInfo())
You can also create an empty rasterInfo class object and populate raster metadata using the FromJSONString() method on the rasterInfo() object. Subsequently, this object can be used to create an empty raster as shown below.
ri = arcpy.RasterInfo()
ri_as_jsonstring = '{"bandCount":1,' \
'"extent":{"xmin":0,"ymin":0,"xmax":10,"ymax":10,"spatialReference":{"wkid":26912}},' \
'"pixelSizeX":1,' \
'"pixelSizeY":1,' \
'"pixelType":"F32"}'
ri.fromJSONString(ri_as_jsonstring)
outRas = Raster(ri)
Indexing a raster object
Indexing a raster object allows read and write access to a specific cell in a raster dataset. Indexing can be used within or outside a Raster Cell Iterator.
For a single-band raster, a raster object has two indices, [i, j], which are the row and column indices, respectively. A minimum of two indices are required for querying and assigning values. All indices are zero-based, starting from zero. Row indices have a value range, defined by 0 to n(rows) - 1, where n(rows) is the number of rows in that raster. Similarly, column indices have a value range, defined by 0 to n(cols) - 1, where n(cols) is the number of columns in that raster.
Indexing starts from the upper left corner of a raster dataset, where the upper leftmost cell has an index of [0, 0], and the lower rightmost cell has an index of [n(rows) - 1, n(cols) - 1].
Handling NoData cells
When a cell with NoData is queried using index notation, within or outside the context of an RCI, the value returned in Python is NaN. For example, if it is known that myRas has a NoData value for the row and a column index of [2, 2], querying this cell value would return NaN. This can be verified using the math.isnan() method, which would return a Boolean value of True for this cell.
import math
math.isnan(myRas[2,2])
out: True
If a cell is assigned NaN in Python, the raster cell becomes NoData, regardless of the format or NoDataValues for that raster dataset. This is the correct way of turning a cell to NoData. For example, assigning NaN to myRas at the row and column index of [3, 3] would modify that cell to NoData.
myRas[3,3] = math.nan
It is recommended that you do not turn a cell to NoData by assigning the noDataValues of that raster dataset directly to a cell, but instead follow the method described above.
Advanced operations using Raster Cell Iterator
In this section, you will gain a conceptual understanding of advanced operations that can be performed using RCI. You will learn how to create an explicit Raster Cell Iterator to iterate over multiple rasters, and how to use the padding and skipNoData options. You will learn indexing of multiband rasters and use it within an RCI. You will learn how to handle different raster pixel types and how to deal with out-of-bound values.
Explicit Raster Cell Iterator
Another way of invoking an RCI is using an explicit RasterCellIterator class object. While the implicit iterator is suited for iterating over a single raster dataset, the explicit raster iterator is designed to iterate over multiple rasters and is optimized for this purpose. While iterating over multiple rasters, you can perform analysis with input rasters whose cell values are used in the computation and output rasters whose cell values are written to during the computation.
The following example demonstrates using explicit RCI to iterate over multiple input and output rasters. The input rasters myRas1 and myRas2 have the same number of rows and columns and other properties such as cell size, extent, and spatial reference. The output rasters outRas1 and outRas2 are created using the rasterInfo class object of the first input raster. At each raster cell, a focal operation is performed to compute mean values for a 3 by 3 cell neighborhood. for myRas1 and myRas2. At this cell location, outRas1 is assigned a cell value as the minimum of the means, and outRas2 is assigned a cell value as the maximum of the means.
myRas1 = Raster("myras1")
myRas2 = Raster("myras2")
outRas1 = Raster(myRas1.getRasterInfo())
outRas2 = Raster(myRas1.getRasterInfo())
with RasterCellIterator({'rasters':[myRas1, myRas2, outRas1, outRas2]}) as rci:
for i,j in rci:
meanMyRas1 = (myRas1[i-1,j-1] + myRas1[i-1, j] + myRas1[i-1, j+1] + \
myRas1[i,j-1] + myRas1[i, j] + myRas1[i, j+1] + \
myRas1[i+1,j-1] + myRas1[i+1, j] + myRas1[i+1, j+1]) / 9
meanMyRas2 = (myRas2[i-1,j-1] + myRas2[i-1, j] + myRas2[i-1, j+1] + \
myRas2[i,j-1] + myRas2[i, j] + myRas2[i, j+1] + \
myRas2[i+1,j-1] + myRas2[i+1, j] + myRas2[i+1, j+1]) / 9
outRas1[i,j] = min(meanMyRas1, meanMyRas2)
outRas2[i,j] = max(meanMyRas1, meanMyRas2)
outRas1.save()
outRas2.save()
When multiple rasters are used in the iterator, the first raster in the list defines the iterator's raster analysis environment. All other rasters in the list are expected to match in spatial reference, cell size, and extent with that of the first raster. If these properties of the second and subsequent rasters in this list don't match with the first raster, the analysis environment of the first raster is applied to them. It is recommended that you create all output rasters using the rasterInfo() object of the first raster so that they inherit its spatial reference, cell size, and extent.
Additionally, if you want to honor environment settings, they should be applied to the first raster in this list using the ApplyEnvironment function before using it in the iterator. This is demonstrated in the example below.
myRas1 = Raster("myras1")
myRas2 = Raster("myras2")
#Apply environment settings to myRas1
myRas1_env = arcpy.sa.ApplyEnvironment(myRas1)
outRas1_env = Raster(myRas1_env.getRasterInfo())
with RasterCellIterator({'rasters':[myRas1_env, myRas2, outRas1_env]}) as rci:
for i,j in rci:
meanMyRas1 = (myRas1_env[i-1,j-1] + myRas1_env[i-1, j] + myRas1_env[i-1, j+1] + \
myRas1_env[i,j-1] + myRas1_env[i, j] + myRas1_env[i, j+1] + \
myRas1_env[i+1,j-1] + myRas1_env[i+1, j] + myRas1_env[i+1, j+1]) / 9
meanMyRas2 = (myRas2[i-1,j-1] + myRas2[i-1, j] + myRas2[i-1, j+1] + \
myRas2[i,j-1] + myRas2[i, j] + myRas2[i, j+1] + \
myRas2[i+1,j-1] + myRas2[i+1, j] + myRas2[i+1, j+1]) / 9
outRas1_env[i,j] = min(meanMyRas1, meanMyRas2)
outRas1_env.save()
Padding option
The explicit iterator has an option to provide padding. Padding does not change the output in any way, but it improves performance if you are accessing neighborhood cell values while iterating through cell locations. For example, if a neighborhood operation has a kernel that accesses cell values up to two row cells or two column cells away from the cell location under iteration, provide a padding of 2.
with RasterCellIterator({'rasters':[myRas1, outRas1], 'padding': 2}) as rci_padded:
for i,j in rci_padded:
outRas1[i,j] = (myRas1[i-2,j-2] + myRas1[i-2, j] + myRas1[i-2, j+2] + \
myRas1[i,j-2] + myRas1[i, j] + myRas1[i, j+2] + \
myRas1[i+2,j-2] + myRas1[i+2, j] + myRas1[i+2, j+2]) / 9
outRas1.save()
Skip NoData option
The explicit iterator has the capability to skip NoData cells entirely, instead of handling them on a per-cell basis. Use the skipNoData key to provide a list of rasters that would determine which raster cells are skipped by the iterator. For any given cell, if any of the rasters in this list has a NoData value, that cell would be skipped. Since the state of NoData cells in these rasters can change over the course of an iteration, whether or not skipping is applied is decided by the state of the rasters at each cell, as it is being visited by the iterator. The primary benefit of this option is the improvement in performance when iterating over sparse raster datasets, since many NoData cells no longer need to be visited.
In the following example, a sparse raster streamRas has values only at cells representing river streams, and NoData everywhere else. A flow direction raster, flowDirRas, has values everywhere within the raster extent without any NoData cells. While using RCI to extract the flow direction raster values along the stream cells, you can use the skipNoData key such that only stream cells with values are visited by the iterator, minimizing the number of total cells visited.
streamRas = Raster("streams")
flowDirRas = Raster("flowdir")
outRas = Raster(flowDirRas.getRasterInfo())
with RasterCellIterator({'rasters':[streamRas, flowDirRas, outRas],'skipNoData':[flowDirRas, streamRas]}) as rci_skip:
for i,j in rci_skip:
outRas[i,j] = flowDirRas[i,j]
outRas.save()
Indexing of multiband rasters
For multiband rasters, you can provide two or three indices. The differences in behavior between them is discussed below.
If two indices are provided, these are interpreted as row and column indices, and a tuple of values is returned from each band for the [row, column] index specified. For example, a multiband raster "mbras" with three bands can be queried with two indices, say [2,2]. These are interpreted as row and column indices, and values from each of the three bands are returned as a tuple.
multibandRas = Raster("mbras")
multibandRas[2,2]
out: (4,5,6)
If three indices are specified, these are interpreted as [band, row, column] indices, and the cell value is returned for that band, row, and column cell location. Band index is zero-based and its value range for a multiband raster is determined by 0 to n(bands) - 1, where n(bands) is the total number of bands.
multibandRas[0,2,2]
out: 4
Handling out of bound values
All rasters have a pixel type and a corresponding bit depth. Refer to Bit depth capacity for raster dataset cells in the ArcMap online help for a complete list of bit depths and range of values that each cell can contain. When a value is assigned to a raster cell that exceeds the range of acceptable values for its bit depth, an out-of-bound situation occurs. In such situations, the cell is assigned NoData. For example, the pixel type of a raster object myRas1 is 'U8', which is an unsigned 8-bit integer. This raster can contain integer values in the range 0-255. If a value of 260 is assigned to a cell in myRas1, it is assigned NoData. Following this assignment, querying the value at this cell would return NaN.
myRas1 = Raster("myras1")
myRas1Info = myRas1.getRasterInfo()
pixelType_MyRas1 = myRas1Info.getPixelType()
print(pixelType_MyRas1)
out: 'U8'
myRas1[3,3] = 260
math.isnan(myRas1[3,3])
out: True