Looping through equations in text file, Python script

Discussion created by dprosser on Feb 1, 2012
Latest reply on Feb 2, 2012 by dprosser
Hi all,

I am working on a script that will iterate through a number of geoprocessing tasks for each of 30+ animal species.
I am able to use os.listdir to list the inputs from different directories, then set up a for loop to have it use the species code (7 character code) in the file name to call the correct input files.
For my Single Output Map Algebra step, I would like to call the equation from a text file so that all 30+ equations are in a managable place, however I'm having trouble figuring out how to do this.

Attached is my current draft code; towards the bottom (see #Process 5) is the part where I'm getting stuck (see blue text). I'm attempting to open the text file, read the species code, and use the equation for the correct species for that iteration. As an example, here are the contents of the text file (just for 2 species and with truncated equations for testing)- not sure if I set this up right:

Example contents of the text file:
speciesCode equation
069swgo ( ild64p + ild41p + ild42p ) > (0)
074bhgo ( ild64p + ild41p ) > (0)

I'm using Arc10.0 SP3, Python 2.6, running script using PythonWin.
Finally, here is the output from the interactive window; the error message in red:

Printing shpRange
['r069swgo.shp', 'r074bhgo.shp']
printing valListShpBr
['v069swgo_br.shp', 'v074bhgo_br.shp']
Printing speciesCode: 069swgo
Printing currentValName: v069swgo_br
Printing shapfile: r069swgo.shp
Printing shapefileName: r069swgo
Running 50km buffer....
Clipping to validation layer....
Converting Vector Range to Raster (Polygon to Raster)...
Traceback (most recent call last):
File "C:\Python26\ArcGIS10.0\Lib\site-packages\pythonwin\pywin\framework\scriptutils.py", line 312, in RunScript
exec codeObject in __main__.__dict__
File "D:\GIS_Watf\Watf_Br.py", line 95, in <module>
print "Current Species Code from Text file is:" + " " + currentSppCode #just so I can see which one it's calling
NameError: name 'currentSppCode' is not defined
#DJP Updated 01 February 2012
#Revised for ArcGIS 10 on GIS desktop
#   slicing: general form is l[start:end:step]
#   to indent, use 4 spaces not a tab since tabs can be interpreted differently depending on the platform or tool being used
#   if, while, and for statements

# Import system modules
import sys, string, os, arcgisscripting

# Create the Geoprocessor object
gp = arcgisscripting.create()

# Set the necessary product code

# Check out any necessary licenses

#Python will overwrite exsisting outputfiles
gp.overwriteoutput = 1

# Load required toolboxes...(## commands below are for ArcGIS 9.*)
gp.AddToolbox("C:/ArcGIS/Desktop10.0/ArcToolbox/Toolboxes/Spatial Analyst Tools.tbx")
gp.AddToolbox("C:/ArcGIS/Desktop10.0/ArcToolbox/Toolboxes/Conversion Tools.tbx")
gp.AddToolbox("C:/ArcGIS/Desktop10.0/ArcToolbox/Toolboxes/Data Management Tools.tbx")
gp.AddToolbox("C:/ArcGIS/Desktop10.0/ArcToolbox/Toolboxes/Analysis Tools.tbx")
#gp.AddToolbox("C:/ArcGIS/ArcToolbox/Toolboxes/Spatial Analyst Tools.tbx")
#gp.AddToolbox("C:/ArcGIS/ArcToolbox/Toolboxes/Conversion Tools.tbx")
#gp.AddToolbox("C:/ArcGIS/ArcToolbox/Toolboxes/Data Management Tools.tbx")
#gp.AddToolbox("C:/ArcGIS/ArcToolbox/Toolboxes/Analysis Tools.tbx")

gp.workspace = "D:/GIS_Watf" #(**change back to C: for laptop)
outWorkspace = "D:/GIS_Watf/Output/"
inRange = "D:/GIS_Watf/Data/Ranges/" #path to input range maps

RangeList = os.listdir(inRange) #list what's in the Ranges folder

shpRange = [] #create empty variable
[shpRange.append(i) for i in RangeList if i[-4:] == '.shp'] #for each file in the RangeList list, append if it ends in ".shp"
print "Printing shpRange"
print shpRange

inValBr = "D:/GIS_Watf/Data/Valida/Breeding/" #create variable named inValBr that shows path to Breeding Validation data
inValWi = "D:/GIS_Watf/Data/Valida/Winter/" #create variable named inValBr that shows path to Breeding Validation data
valListFull = os.listdir(inValBr) #list all files in Valida/Breeding folder
valListShpBr = []
[valListShpBr.append(i) for i in valListFull if i[-4:] == '.shp'] #append to list only if they are shapefiles
print "printing valListShpBr"
print valListShpBr

for shapefile in shpRange: #starting a loop with new variable "shapefile" that loops thru shpRange
    shapefileName = shapefile[:-4] #new variable shapefileName that holds the shapefile name minus the .shp
    speciesCode = shapefileName[1:8] #new variable speciesCode that holds the 2nd thru 4th character from the file name (this is the species number)
    print "Printing speciesCode:" + " " + speciesCode #print message and the species number
    for valShapefile in valListShpBr: #starts new loop with variable valShapefile and loops thru the valListShp list
        if speciesCode in valShapefile: #chooses validation file having the speciesCode (from Range file)
            currentValFile = valShapefile
    currentValName = currentValFile[:-4] #currentValName is the validation shapefile (minus the ".shp")
    print "Printing currentValName:" + " " + currentValName
    print "Printing shapfile:" + " " + shapefile
    print "Printing shapefileName:" + " " + shapefileName

    # Process 1: BHGO Summer Range (Select)...
    selectOutput = outWorkspace + shapefileName + "_br.shp" #** change this to "_br.shp"
    selectInput = inRange + shapefile
    gp.Select_analysis(selectInput, selectOutput, "season = 'summer'")
    # Process 2: (Buffer) 50 km...
    print "Running 50km buffer...."
    bufferOutput = outWorkspace + shapefileName + "_buff50.shp"
    gp.Buffer_analysis(selectOutput, bufferOutput, "50 Kilometers", "FULL", "FLAT", "NONE", "")   
    # Process 3: Clip Winter Validation Points by Range...
    print "Clipping to validation layer...."
    valClipOutput = outWorkspace + currentValName + "_clip.shp"
    valInput = inValBr + currentValFile
    gp.Clip_analysis(valInput, bufferOutput, valClipOutput, "")

    # Process 4: Convert Vector Range to Raster  (Polygon to Raster)...'
    print "Converting Vector Range to Raster (Polygon to Raster)..."
    convertRangeOut = outWorkspace + shapefileName + "_br"
    gp.PolygonToRaster_conversion(bufferOutput, "season", convertRangeOut, "CELL_CENTER", "", "1000")

    # Process 5: Create Breeding Distribution Model for each species (Map Algebra)...
    # This process should find the correct species equation in a text file and perform the SOMA on it    SppEquaFile = "D:/GIS_Watf/Command/SpeciesEquations.txt"
    f = open(SppEquaFile)
    for speciesCode1 in f:
        if speciesCode1 in f:
            currentSppCode = speciesCode1 #Here I'm trying to get it to read the species code, move to the correct one, then use the equation that follows it; I realize I probably need another step here...    print "Current Species Code from Text file is:" + " " + currentSppCode #just so I can see which one it's calling

    v1_if_True = "1"
    v0_if_False = "0"
    SOMAoutput = outWorkspace + "a" + currentSppCode + "_br" #example output name is a069swgo_br
    gp.SingleOutputMapAlgebra_sa(ImNotSure, SOMAoutput) #here see 6 lines above where I'm getting stuck; ImNotSure will change to the correct input variable...

Thank you in advance for any help you can provide. As you might guess, I'm new to Python.