I have an ArcPy code that I am using to find various characteristics about certain watersheds in Indiana. I start with a shapefile of all 800 watersheds and go through it row by row and perform a bunch of calculations on each one. Right now the way it is set up I have one large loop to go through all of these calculations and get all of the outputs into lists to get the output as an excel file. This process takes an extremely long time, and I am new to Python and do not know how to improve the code so that it will be faster.
#search cursor through each row of county subwatersheds file
rows = arcpy.SearchCursor(all_subwatersheds)
for row in rows:
subwatershed_number = row.getValue("HUC12")
subwatershed_list.append(subwatershed_number)
print(subwatershed_number, " is the subwatershed to be worked on")
#select subwatershed from shapefile with all of them
subwatershed_selection = arcpy.SelectLayerByAttribute_management(all_subwatersheds, "NEW_SELECTION",
"HUC12 = " + "'"+subwatershed_number+"'")
#copy selected subwatershed to it's own file
subwatershed = arcpy.CopyFeatures_management(subwatershed_selection, "subwatershed")
#clip dem to subwatershed area
dem_clip = arcpy.Clip_management(dem_Indiana, "#", "dem_clip.tif", subwatershed, "#" , "ClippingGeometry", "NO_MAINTAIN_EXTENT")
print("clipped dem to subwatershed area")
#calculate area of subwatershed
arcpy.AddField_management(subwatershed, "AREA", "DOUBLE")
area = arcpy.CalculateGeometryAttributes_management(subwatershed, "AREA AREA_GEODESIC", '',
"SQUARE_KILOMETERS",
None)
#print area
rows = arcpy.SearchCursor(area)
for row in rows:
area = row.getValue("AREA")
area_list.append(area)
print(area, " square kilometers")
#calculate perimeter of subwatershed
arcpy.AddField_management(subwatershed, "PERIMETER", "DOUBLE")
perimeter = arcpy.CalculateGeometryAttributes_management(subwatershed, "PERIMETER PERIMETER_LENGTH_GEODESIC",
"KILOMETERS")
print("perimeter calculated")
#print perimeter
rows = arcpy.SearchCursor(perimeter)
for row in rows:
perimeter = row.getValue("PERIMETER")
perimeter_list.append(perimeter)
print(perimeter, " kilometers")
#get average slope using zonal statistics from slope dem
avg_slope = arcpy.sa.ZonalStatistics(subwatershed, "FID", dem_slope, "MEAN")
#get slope value
avg_slope_result = arcpy.GetRasterProperties_management(avg_slope, "MAXIMUM")
avg_slope_value = avg_slope_result.getOutput(0)
avg_slope_list.append(avg_slope_value)
print(avg_slope_value, " = avg slope (%)")