arcpy: For loop to sequence Processing of two items in list using ArcPy?

Question

I have multiple rasters, organized by years. I.e. raster y1 = from year 1 raster y2 = from year 2, etc...).

All of rasters contain values 0-11, representing forest (0-6) and disturbance (7-11). I would like to:

• Identify the forest on scene in y1 (by Reclassify tool)

• Use the forest from y1 to Extract by mask raster y2

• Identify the forest on scene in y2 (by Reclassify tool)

• Use the forest from y2 to Extract by mask raster y3

• Identify the forest on scene in y3 (by Reclassify tool)

• Use the forest from y3 to Extract by mask raster y4

etc for every year

However, I can't imagine how can I automate this process? I tried to run for loop and advance the raster index value by one to move to following raster in a list of rasters. Instead I have created unfinite loop, where the for loop always start with the first raster in a list... Or should I use another type of loop?

My approach:

# Import system modules
import arcpy
from arcpy.sa import *
# Set environment settings
wd = "C:/Users/Projects/2017_deforest"
env.workspace = wd + "/raw/SAO"
arcpy.env.overwriteOutput = True
# Set local variables
rasters = arcpy.ListRasters("sao*", "IMG")
reclassField = "Value"
remap = RemapRange([[0,0, "NODATA"],
 [1,6,1], # forest in that year
 [7,11, "NODATA"], # disturbance 
 [12,255, "NODATA"]])
# Check out the ArcGIS Spatial Analyst extension licence
arcpy.CheckOutExtension("Spatial")
# Execute Reclassify
for raster in rasters:
 # Identify forest in raster y1
 tempForest = Reclassify(raster, reclassField, remap, "NODATA")
 # Create list index to move to following raster in a list
 rasNumber = 1 
 inRaster2 = rasters[rasNumber]
 # Execute Extract by mask
 outName = wd + "/outputSAO/" + inRaster2
 outExtractByMask = ExtractByMask(inRaster2, tempForest)
 # advance the index by one to move to another element in a list
 rasNumber = rasNumber + 1 
 # Save the output
 outExtractByMask.save(outName)
# Check in the ArcGIS Spatial Analyst extension license
arcpy.CheckInExtension("Spatial")

Answer 1

The problem is that you re-define

rasNumber = 1 

in each iteration, so you will be accessing the same raster in each iteration. This would not result in infinite loop though, because the for loop will exhaust once all raster files are iterated (of course you would get wrong results as the same - the first - raster was used). You should move this line outside of the for loop.

rasNumber = 1 
for raster in rasters:

When you will run your modified script, you would get an IndexError on the line inRaster2 = rasters[rasNumber] when you are processing the last file. This is because there is no raster that comes after it, so you cannot get from the list. This can be solved either by wrapping the indexing operation into a try/except statement:

try:
 inRaster2 = rasters[rasNumber]
except:
 break

Alternatively, you could just check the item being accessed:

if rasNumber > len(rasters):
 break

The for loop will break when you would start processing the last raster (you cannot extract by mask the non-existing raster).

Apart from that, some tips on the code that will make it safer:

• Don't use string concatenation as in outName = wd + "/outputSAO/" + inRaster2 and env.workspace = wd + "/raw/SAO". Instead, use the os.path.join to safely join folder paths and file names.

• Don't import all objects from the arcpy.sa module using from arcpy.sa import * because this takes extra time every time you run the script.

• Alex Tereshenkov: No worries! I've updated my answer to explain on this. You may add print rasNumber statement in the loop to see what is the number so you get the logic right for rasNumber > len(rasters). It maybe rasNumber == len(rasters) depending on what number is used as a starter.
• maycca: Cool, thank you! please, when I run the modified code, I got an error at the end ' IndexError: list index out of range' and my inRaster2 = rasters[rasNumber]. How can I nicely handle this in my script, should I just use while loop until the last raster? Thank you for your tips, but can you please be little bit more specific? how could I correctly include your suggestions in my script?
• Alex Tereshenkov: Just change the inRaster2 = rasters[rasNumber] line to be if rasNumber > len(rasters): break as I have outlined in the answer. Does it help?
• maycca: yes, definitely !! thank you! I am just checking the documentation to correctly use the os.path.join and if I can I replace from arcpy.sa import * to make is less memory demanding :-D those are very helpful tips!
• Alex Tereshenkov: No worries. Excellent reading on paths - community.esri.com/blogs/dan_patterson/2016/08/14/…. You can just write arcpy.sa.Reclassify after you have done import arcpy. this will suffice

Answer 2

Thanks to @alextereshenkov answer I have put together my code. As it took me a while where exactly to put individual statements, here I attach the whole script for the record:

# Import system modules
import arcpy, os
# Set environment settings
wd = "C:/Users/Projects/2017_deforest"
arcpy.env.workspace = os.path.join(wd,"raw/SAO")
arcpy.env.overwriteOutput = True
# Set local variables
rasters = arcpy.ListRasters("sao*", "IMG")
reclassField = "Value"
remap = arcpy.sa.RemapRange([[0, 0, "NODATA"], # need to specify the Spatial extension 
 [1, 6, 2], # forest in that year
 [7, 11, "NODATA"], # disturbance 
 [12, 255, "NODATA"]])
# Check out the ArcGIS Spatial Analyst extension licence
arcpy.CheckOutExtension("Spatial")
# define the index of the first raster in the list
rasNumber = 0
for raster in rasters:
 try:
 # Identify forest in raster y1
 tempForest = arcpy.sa.Reclassify(raster, reclassField, remap, "NODATA")
 # Execute Extract by mask 
 # print "(1) rasNumber is " + str(rasNumber)
 # Define the position of the second raster in the list
 inRaster2 = rasters[rasNumber]
 print "inRaster2 is " + str(inRaster2)
 outExtractByMask = ExtractByMask(inRaster2, tempForest)
 # Advance raster position
 rasNumber = rasNumber + 1
 print " (2)rasNumber is " + str(rasNumber)
 # Save the output
 outName = os.path.join(wd, "outputSAO", inRaster2)
 outExtractByMask.save(outName)
 except:
 # Break into code after overapassing the number of rasters in a list
 break 
# Check in the ArcGIS Spatial Analyst extension license
arcpy.CheckInExtension("Spatial")
print "DONE"