Map algebra in Spatial Analyst

Available with Spatial Analyst license.

Map algebra allows you access to the Spatial Analyst tools, operators, functions, and classes through algebra. In its most basic form, an output raster is specified to the left of an equal sign (=), and the tools, operators, and their parameters are on the right. See the following example:

from arcpy.sa import *
outShade = Hillshade("inelevation.tif", 99, 33)

The above statement calculates a hillshade, determining the illumination with the sun being at an azimuth of 99 degrees and an altitude of 33 degrees, and creates a Raster object named outShade to store the results.

Map algebra can run simple statements, but the power of the language is realized when creating complex statements and models. As map algebra has been integrated in Python, all the functionality of Python and ArcPy and its extensions (modules, classes, functions, and properties) are available to you, the modeler.

As your needs grow, you can explore many of the facets of map algebra. The following quick tour will give you the essentials to get started.

Basics of running map algebra

There are three ways to use map algebra:

  • The Raster Calculator tool
  • The Python window
  • A Python integrated development environment (IDE)

Raster Calculator

The Raster Calculator tool runs map algebra expressions. The tool has a calculator interface from which most map algebra statements can be created by clicking buttons. Raster Calculator can be used as a stand-alone tool, but it can also be used in ModelBuilder. As a result, the tool allows map algebra to be integrated into ModelBuilder.

Raster Calculator user interface
The Raster Calculator user interface is shown.

In the above expression, three rasters are combined by multiplying the second and third rasters together and adding that result to the first. Note that operators follow a defined order of precedence.

The Raster Calculator tool is not intended to replace other Spatial Analyst tools. Continue to use the other tools for the appropriate calculations. For example, use the Slope tool to perform the slope calculations. The Raster Calculator tool is designed to run single-line algebraic statements.

Since Raster Calculator is a geoprocessing tool, like all tools, it can be integrated into ModelBuilder. See the following topics for more information:

Python window

The Python window allows you to use geoprocessing tools and Python functionality from within ArcGIS Pro. The Python functionality you run from this window can range from a single line to complex multiline blocks of code. The Python window also provides a place to access additional functionality using third-party Python modules and libraries.

To launch the Python window, on the Analysis tab, click the Python drop-down menu, and choose the Python window button Show Python window. Alternatively, on the View tab, in the Windows group, click the Python window.

Example of the Python window

In the above sequence of statements, the ArcPy site package, the geoprocessing environments, and the Spatial Analyst modules are imported; the workspace is set; the Slope tool is run; and the output is permanently saved. Once you press Enter at the end of a statement, that statement is immediately run.

Some features of the Python window include built-in line automatic completion, use of variables, and access to Python and ArcPy functionality.

Python integrated development environment

Even though there is no limit to the number of statements that can be entered in the Python window, it may be cumbersome to create more complex models. The Spatial Analyst modules' tools, operators, functions, and classes can also be accessed from your favorite integrated development environment. Start your preferred IDE and enter the desired statements.

In the following script, ArcPy, the geoprocessing environments, and the Spatial Analyst module are imported; variables are set; the extension is checked out; the Slope tool is run; and the output is saved.

# Name: Slope
# Description: Identifies the rate of maximum change
#               in z-value from each cell.
# Requirements: Spatial Analyst Extension


# Import system modules
import arcpy
from arcpy import env
from arcpy.sa import *

# Set environment settings
env.workspace = "C:/data"

# Set local variables
inRaster = "elevation"
outMeasurement = "DEGREE"
zFactor = 0.3043

# Check out the ArcGIS Spatial Analyst extension license
arcpy.CheckOutExtension("Spatial")

# Run Slope
outSlope = Slope(inRaster, outMeasurement, zFactor)

# Save the output
outSlope.save("C:/output/outslope02")

As is the case with the Python window, an IDE provides access to all available Python and ArcPy functionality.

Work with tools

All Spatial Analyst tools that output a raster are available in algebraic format. The dataset name can be used if it is in the Contents window or in the current workspace; otherwise, the full path must be entered.

# In the following statement, indem is either  
#   in the TOC or in the current workspace
outRas = Aspect("indem")

# In the following statement the full path is specified
outRas2 = Aspect("C:/Data/indem2.tif")

The output from one statement can be entered into a subsequent statement.

outRas = Select("inras.tif", "Value > 105")

# outRas is variable defined by the previous statement and is not quoted
outdist = EucDistance(outRas)

Work with operators

Map algebra supports a series of operators (for example, +, -, and *). These same operators also exist in Python but are modified for map algebra to handle Raster objects differently. For example, the following adds two numbers together into a variable:

# set outVar to 14 using the Python + operator
outVar = 5 + 9

To distinguish that the statement should work on rasters (that is, to use the Spatial Analyst operator), you must cast the dataset as a Raster. The following example uses the Spatial Analyst + operator to add two rasters together:

outRas = Raster("inras1.tif") + Raster("inras2.tif")

Operators can accept a mixture of rasters and numbers. For example, the following adds a constant value of 8 to all the cells in the input raster:

outRas = Raster("inras1.tif") + 8

Create complex expressions

Tools and operators can be strung together in a single statement. The following example runs several tools and operators in each expression:

outRas = Slope("indem" * 2) / 57
outdist = EucDistance(ExtractByAttributes("inras", "Value > 105"))

Parentheses can be used to control the order of processing. Consider the following two examples, which use the same operators but yield different results due to the use of parentheses:

outRas1 = (Raster("inras1") + Raster("inras2")) / Raster("inras3")

and

outRas2 = Raster("inras1") + Raster("inras2") / Raster("inras3")

In the first statement, inras1 is added to inras2, and the result is divided by inras3. Without the parentheses, as in the second statement, inras2 would be divided by inras3, and the result would be added to inras1.

When using multiple Boolean (~, &, ^, |) or Relational (<, <=, >, >=, ==, !=) operators consecutively in a single expression, parentheses should be used. For example, parentheses are required in the expression: (a>2) & (a<5). If parentheses are not used, the expression will result in an error: a>2 & a<5. The following expression will run because parentheses are used:

outRas = (Raster("a") > 2) & (Raster("a") < 5)
Dive-in:

Some expressions may not simply require parentheses but may instead have to be rewritten. For example, an expression in the form of a < b < c will not run, and adding parentheses will change the meaning of the expression. Therefore, to run successfully, this expression must be rewritten in the form of (a < b) & (b < c).

Use classes

Classes are used in map algebra tools for parameters that have multiple arguments. The use of classes for input parameters allows you to access the individual arguments of a parameter to query, alter, and add arguments. The following is an example showing the use of a class:

outRas = FocalStatistics("inRaster", NbrCircle(5, "CELL"), "SUM")

In the above statement, the sum is calculated for each cell within a five-cell circular neighborhood. NbrCircle is a class that creates a NbrCircle object.

The following is an example of a remap table class. Any number of values can be entered into a remap class.

outReclass = Reclassify("inRaster", "VALUE", RemapRange([[0, 1], [3, 10], [4, 8]]))

In the above statement, a class, RemapRange, is used to define the reclassification of the input values. The cells with a value of 0 in inRaster will be assigned to 1 in outReclass, and 3 will be assigned to 10, and 4 to 8.

Map algebra functions that output features, tables, or files

Only Spatial Analyst tools that produce a raster as output are implemented using the algebraic format. For Spatial Analyst tools that produce output that is not a raster (for example, features, tables, or text files), the output is specified as a parameter within the tool in the parentheses. Note the syntax in the following example, which creates contours as an output polyline feature dataset:

indem = "C:/Data/indem"
contourInterval = 100
Contour(indem, "C:/output/outcontours", contourInterval)

Suggestions for running map algebra statements

In all the map algebra examples shown below, the output is a Raster object. The Raster object points to a temporary raster dataset that, unless it is explicitly saved, will be removed when the ArcGIS session ends. To permanently save the temporary dataset, the save method is called on the Raster object (see the two examples below).

It is recommended that you set the appropriate analysis environments, in particular Current workspace, Extent, Cell size, Cell size projection method, Mask, and Snap, before implementing the map algebra tool or operator.

An example demonstrating the workspace environment is as follows:

import arcpy 
from arcpy import env 
from arcpy.sa import *

env.workspace = "C:/sapyexamples/data" 

outHillshade = Hillshade("elevation", 180, 75, "SHADOWS", 1)

outHillshade.save("outhillshd01")

In the above statement, the workspace is set, therefore, outhillshd01 will be saved in C:/sapyexamples/data.

It is advisable to set classes for any complex input to a map algebra tool to a variable and use the variable in the statement. In the statement below, a RemapRange class object is set to a variable, myRemapRange, and is used as input to the Reclassify tool.

import arcpy 
from arcpy import env
from arcpy.sa import *

env.workspace = "C:/sapyexamples/data" 
myRemapRange = RemapRange([[-3, 0, 0], [0, 1.75, 25], [1.75, 3.5, 50], 
                           [3.5, 5.25, 75], [5.25, 7, 100]]) 

outReclassRR = Reclassify("inreclass", "VALUE", myRemapRange)

outReclassRR.save("rclassremran")

Further reading

To gain a deeper understanding of ArcPy, explore these topics:

To obtain more information on geoprocessing in Python, the following may be helpful:

Related topics