[GIS] Raster Calculator Workaround in Arcpy for Suitabiliy Analysis

arcpypythonrasterraster-calculatorsuitability-analysis

This is my first attempt at a complex stand alone python code. I am self teaching the arcpy module, well at least attempting to.

I am having a great deal of difficulty getting my site suitability analysis code to run. It consistently crashes at the raster calculator step. I understand that the raster calculator tool is difficult to manage in arcpy, however I am equally unfamiliar with the method using map algebra.

The outputs/inputs are rasters. Everything is converted to rasters in prior steps. I get an error code: 000989 Python Syntax error: parsing error invalid syntax, referencing the line of code for raster calculator.

My code is included below. Can someone help me with a workaround for this issue?

# -*- coding: utf-8 -*-
# ---------------------------------------------------------------------------
# Megan E Brown
# ModelCode.py
# Created on: 2014-04-21 23:51:03.00000
#   (generated by ArcGIS/ModelBuilder)
# Usage: ModelCode <Scratch_Workspace> <Trails> <Trail_Proximity> <DEM> <Min_Elevation__m_> <Max_Slope____> <Final_raster> 
# ---------------------------------------------------------------------------

# Import arcpy module
import arcpy, arcview
#from arcpy.sa import *
#env.workspace = "C:\\temp"

# Check out any necessary licenses
arcpy.CheckOutExtension("3D")
arcpy.CheckOutExtension("spatial")

# Script arguments
Scratch_Workspace = arcpy.GetParameterAsText(0)
if Scratch_Workspace == '#' or not Scratch_Workspace:
    Scratch_Workspace = "C:\\temp" # provide a default value if unspecified

Trails = arcpy.GetParameterAsText(1)
if Trails == '#' or not Trails:
    Trails = "Z:\\GIS500\\FinalProject\\FPData\\WELP_trails.shp" # provide a default value if unspecified

Trail_Proximity = arcpy.GetParameterAsText(2)
if Trail_Proximity == '#' or not Trail_Proximity:
    Trail_Proximity = "200 Feet" # provide a default value if unspecified

DEM = arcpy.GetParameterAsText(3)
if DEM == '#' or not DEM:
    DEM = "Z:\\GIS500\\ned\\ned" # provide a default value if unspecified

Min_Elevation__m_ = arcpy.GetParameterAsText(4)
if Min_Elevation__m_ == '#' or not Min_Elevation__m_:
    Min_Elevation__m_ = "170" # provide a default value if unspecified

Max_Slope____ = arcpy.GetParameterAsText(5)
if Max_Slope____ == '#' or not Max_Slope____:
    Max_Slope____ = "10" # provide a default value if unspecified

Final_raster = arcpy.GetParameterAsText(6)
if Final_raster == '#' or not Final_raster:
    Final_raster = "Z:\\GIS500\\FinalProject\\FPData\\tempdata\\yay" # provide a default value if unspecified

# Local variables:
Percent_Slope = "C:\\temp\\Percent_Slope"            #slope output
Slope_LE = "C:\\temp\\Slope_LE"                      #"output of less than equal"
Calc = "C:\\temp\\Calc"                              #output from raster calculator, input to final step Equal to
Elevation_GE = "C:\\temp\\Elevation_GE"              #"greater than equal output"
Buff = "C:\\temp\\Buff"                              #"output of buffer"
buffer_raster = "C:\\temp\\buffer_raster"            #"output of feature to raster making buffer into polys"
Near_trails = "C:\\temp\\Near_trails"                #"output from not equal reclassing the buffer raster for raster calculation"
Cell_size_30m = 30
Value_0 = "0"
Value3 = "3"

# Process: Slope
arcpy.Slope_3d(DEM, Percent_Slope, "PERCENT_RISE", "9.0021749575298E-06")

# Process: Less Than Equal
arcpy.gp.LessThanEqual_sa(Percent_Slope, Max_Slope____, Slope_LE)

# Process: Greater Than Equal
arcpy.gp.GreaterThanEqual_sa(DEM, Min_Elevation__m_, Elevation_GE)

# Process: Buffer
arcpy.Buffer_analysis(Trails, Buff, Trail_Proximity, "FULL", "ROUND", "NONE", "")

# Process: Feature to Raster
arcpy.FeatureToRaster_conversion(Buff + ".shp", "FID", buffer_raster, Cell_size_30m)

# Process: Not Equal
arcpy.gp.NotEqual_sa(buffer_raster, Value_0, Near_trails)

# Process: Raster Calculator
arcpy.gp.RasterCalculator_sa(((Slope_LE) + (Elevation_GE) + (Near_trails)), Calc)

# Process: Equal To
arcpy.gp.EqualTo_sa(Calc, Value3, Final_raster)

Best Answer

Just looking around for an answer to this and I've found the following: On this answer to a previous post, Luke indicates that from the ArcGIS Help:

The Raster Calculator tool is intended for use in the ArcGIS Desktop application only as a GP tool dialog box or in ModelBuilder. It is not intended for use in scripting and is not available in the ArcPy Spatial Analyst module.

As an alternative to the Raster Calculator, you might need to investigate using the Math Toolset. Specifically, the Plus tool to do what you want.

I haven't tried this, but it looks like it might do what you want. Note, the Math Tools take two inputs, so you might need to create a temporary file. Something similar to this:

from arcpy.sa import *

PlusTemp = Plus(Slope_LE, Elevation_GE)
Calc = Plus(PlusTemp, Near_trails)
Related Question