Raster processing for climate change model











up vote
7
down vote

favorite












I am working on a project on how climate change may affect the development of spruce budworm larvae throughout each year (historical and projected) using different climate change models and emission scenarios.



The code that I have works, but would take a ridiculous amount of time to run: 38 years when extrapolated out. I'm using python 2.7.12 on a 64-bit windows machine with 24GB RAM and xeon processors. So, I don't think machine specs are an issue... it's how my script is written and organized.



I'm fairly new to Python, and have read through some articles on NumPy arrays and scipy packages, but do not have the knowledge/skills on how to implement them in this script to optimize the raster processing. Is numpy or scipy the way to go, and can those packages accomplish what I am trying to do?



Each time a new 4-h timestep raster is accessed, 27 conditional statements are executed, which probably is the reason the script runs so slow. Does anyone see a way to reduce the amount of conditional statements, but achieve the same goal?



Any help on how to improve the overall performance so that the script runs in a reasonable amount of time would be greatly appreciated!



import arcpy
from arcpy.sa import *
import math

arcpy.CheckOutExtension('spatial')
arcpy.env.overwriteOutput = True

ClimateModels = ["CanESM2", "CSIRO-Mk3-6-0", "HadGEM2-ES"]
ConcentrationPathway = ["rcp26", "rcp45", "rcp85"]

#100 random numbers between 0.4 and 2.5, representing the stochastic growth rates of individual spruce budworm larvae ....need to be the same each time it loops through GCM+RCP combination
BudwormList = [2.1, 2.3, 0.6, 2.1, 1.7, 1.1, 1.6, 2.4, 2.1, 2.2, 0.9, 1.6, 0.6, 2.0, 1.7, 0.6, 0.7, 1.8, 1.3, 2.4,
2.1, 0.8, 2.2, 1.6, 0.5, 1.0, 1.6, 2.4, 2.5, 2.5, 0.8, 0.5, 2.5, 0.9, 1.2, 0.7, 1.8, 1.8, 0.4, 1.9,
0.8, 1.3, 0.6, 2.2, 0.7, 1.6, 1.0, 0.7, 0.8, 2.4, 1.7, 0.5, 1.1, 2.2, 0.9, 1.6, 1.8, 0.8, 2.3, 0.9,
0.8, 1.0, 0.5, 1.6, 2.5, 2.2, 2.0, 1.3, 1.6, 1.2, 2.4, 0.9, 1.1, 0.8, 2.1, 2.2, 1.5, 2.4, 2.5, 2.2,
1.3, 0.4, 1.5, 1.3, 1.7, 1.9, 1.8, 0.4, 0.4, 0.5, 1.5, 0.5, 0.6, 0.7, 1.6, 0.9, 0.9, 2.5, 1.4, 1.3]
TimeStep = ['4h', '8h', '12h', '16h', '20h', '24h']

try:
for GCM in ClimateModels:
for RCP in ConcentrationPathway:
BudwormCount = 1
for budworm in BudwormList:
year = 1971
while year <= 2070:
if ((year%400 == 0) or ((year%4 == 0) and (year%100 != 0))):
days = 366
else:
days = 365
Dev_stage = Raster("C:\Users\Robert\Documents\PCIC_ClimateData\TEST\Dev_Stage.tif")
Dev_stage_female = Raster("C:\Users\Robert\Documents\PCIC_ClimateData\TEST\Dev_Stage_female.tif")
Dev_L2o = Raster("C:\Users\Robert\Documents\PCIC_ClimateData\TEST\Dev_L2o.tif")
Dev_L2 = Raster("C:\Users\Robert\Documents\PCIC_ClimateData\TEST\Dev_L2.tif")
Dev_L3 = Raster("C:\Users\Robert\Documents\PCIC_ClimateData\TEST\Dev_L3.tif")
Dev_L4 = Raster("C:\Users\Robert\Documents\PCIC_ClimateData\TEST\Dev_L4.tif")
Dev_L5 = Raster("C:\Users\Robert\Documents\PCIC_ClimateData\TEST\Dev_L5.tif")
Dev_L6_male = Raster("C:\Users\Robert\Documents\PCIC_ClimateData\TEST\Dev_L6_male.tif")
Dev_L6_female = Raster("C:\Users\Robert\Documents\PCIC_ClimateData\TEST\Dev_L6_female.tif")
Dev_pupa_male = Raster("C:\Users\Robert\Documents\PCIC_ClimateData\TEST\Dev_pupa_male.tif")
Dev_pupa_female = Raster("C:\Users\Robert\Documents\PCIC_ClimateData\TEST\Dev_pupa_female.tif")
L2_date = Raster("C:\Users\Robert\Documents\PCIC_ClimateData\TEST\Jul_L2.tif")
L3_date = Raster("C:\Users\Robert\Documents\PCIC_ClimateData\TEST\Jul_L3.tif")
L4_date = Raster("C:\Users\Robert\Documents\PCIC_ClimateData\TEST\Jul_L4.tif")
L5_date = Raster("C:\Users\Robert\Documents\PCIC_ClimateData\TEST\Jul_L5.tif")
L6_date = Raster("C:\Users\Robert\Documents\PCIC_ClimateData\TEST\Jul_L6.tif")
pupa_male_date = Raster("C:\Users\Robert\Documents\PCIC_ClimateData\TEST\Jul_pupa_male.tif")
pupa_female_date = Raster("C:\Users\Robert\Documents\PCIC_ClimateData\TEST\Jul_pupa_female.tif")
adult_male_date = Raster("C:\Users\Robert\Documents\PCIC_ClimateData\TEST\Jul_adult_male.tif")
adult_female_date = Raster("C:\Users\Robert\Documents\PCIC_ClimateData\TEST\Jul_adult_female.tif")

current_day = 1
while current_day <= days:
for time in TimeStep:
tempfile = Raster("C:\Users\Robert\Documents\PCIC_ClimateData\" + GCM + "\" + RCP + "\4hTempTimeStep\" + time + "_" + str(current_day) + "_" + str(year) + ".tif")
Dev_L2o = Con(((Dev_stage == 1) & ((tempfile > 2.5) & (tempfile < 35)) & (Dev_L2o <= 1)), Dev_L2o + (float(budworm))*((0.194*((1/(1+(math.e**(3.00-(5.84*((tempfile-2.5)/(35-2.5)))))))-(math.e**((((tempfile-2.5)/(35-2.5))-1)/0.034))))/6), Dev_L2o)
Dev_L2 = Con(((Dev_stage == 2) & ((tempfile > 4.4) & (tempfile < 38)) & (Dev_L2 <= 1)), Dev_L2 + (float(budworm))*((0.919*((1/(1+(math.e**(2.91-(5.32*((tempfile-4.4)/(38-4.4)))))))-(math.e**((((tempfile-4.4)/(38-4.4))-1)/0.061))))/6), Dev_L2)
Dev_L3 = Con(((Dev_stage == 3) & ((tempfile > 4.4) & (tempfile < 38)) & (Dev_L3 <= 1)), Dev_L3 + (float(budworm))*((0.438*((1/(1+(math.e**(3.06-(6.85*((tempfile-4.4)/(38-4.4)))))))-(math.e**((((tempfile-4.4)/(38-4.4))-1)/0.061))))/6), Dev_L3)
Dev_L4 = Con(((Dev_stage == 4) & ((tempfile > 4.4) & (tempfile < 38)) & (Dev_L4 <= 1)), Dev_L4 + (float(budworm))*((1.211*((1/(1+(math.e**(3.80-(7.55*((tempfile-4.4)/(38-4.4)))))))-(math.e**((((tempfile-4.4)/(38-4.4))-1)/0.148))))/6), Dev_L4)
Dev_L5 = Con(((Dev_stage == 5) & ((tempfile > 4.4) & (tempfile < 38)) & (Dev_L5 <= 1)), Dev_L5 + (float(budworm))*((0.269*((1/(1+(math.e**(3.02-(8.57*((tempfile-4.4)/(38-4.4)))))))-(math.e**((((tempfile-4.4)/(38-4.4))-1)/0.005))))/6), Dev_L5)
Dev_L6_male = Con(((Dev_stage == 6) & ((tempfile > 4.4) & (tempfile < 38)) & (Dev_L6_male <= 1)), Dev_L6_male + (float(budworm))*((0.288*((1/(1+(math.e**(2.67-(5.03*((tempfile-4.4)/(38-4.4)))))))-(math.e**((((tempfile-4.4)/(38-4.4))-1)/0.151))))/6), Dev_L6_male)
Dev_L6_female = Con(((Dev_stage_female == 1) & ((tempfile > 4.4) & (tempfile < 38)) & (Dev_L6_female <= 1)), Dev_L6_female + (float(budworm))*((0.317*((1/(1+(math.e**(3.06-(4.66*((tempfile-4.4)/(38-4.4)))))))-(math.e**((((tempfile-4.4)/(38-4.4))-1)/0.136))))/6), Dev_L6_female)
Dev_pupa_male = Con(((Dev_stage == 7) & ((tempfile > 4.4) & (tempfile < 35)) & (Dev_pupa_male <= 1)), Dev_pupa_male + ((0.259*((1/(1+(math.e**(2.75-(4.66*((tempfile-4.4)/(35-4.4)))))))-(math.e**((((tempfile-4.4)/(35-4.4))-1)/0.053))))/6), Dev_pupa_male)
Dev_pupa_female = Con(((Dev_stage_female == 2) & ((tempfile > 4.4) & (tempfile < 35)) & (Dev_pupa_female <= 1)), Dev_pupa_female + ((0.205*((1/(1+(math.e**(2.85-(6.28*((tempfile-4.4)/(35-4.4)))))))-(math.e**((((tempfile-4.4)/(35-4.4))-1)/0.044))))/6), Dev_pupa_female)

Dev_stage = Con(((Dev_stage == 1) & (Dev_L2o > 1)), Dev_stage + 1, Dev_stage)
Dev_stage = Con(((Dev_stage == 2) & (Dev_L2 > 1)), Dev_stage + 1, Dev_stage)
Dev_stage = Con(((Dev_stage == 3) & (Dev_L3 > 1)), Dev_stage + 1, Dev_stage)
Dev_stage = Con(((Dev_stage == 4) & (Dev_L4 > 1)), Dev_stage + 1, Dev_stage)
Dev_stage = Con(((Dev_stage == 5) & (Dev_L5 > 1)), Dev_stage + 1, Dev_stage)
Dev_stage = Con(((Dev_stage == 6) & (Dev_L6_male > 1)), Dev_stage + 1, Dev_stage)
Dev_stage = Con(((Dev_stage == 7) & (Dev_pupa_male > 1)), Dev_stage + 1, Dev_stage)
Dev_stage_female = Con(((Dev_stage_female == 1) & (Dev_L6_female > 1)), Dev_stage_female + 1, Dev_stage_female)
Dev_stage_female = Con(((Dev_stage_female == 2) & (Dev_pupa_female > 1)), Dev_stage_female + 1, Dev_stage_female)

L2_date = Con(((Dev_L2o > 1) & (L2_date == 0)), L2_date + current_day, L2_date)
L3_date = Con(((Dev_L2 > 1) & (L3_date == 0)), L3_date + current_day, L3_date)
L4_date = Con(((Dev_L3 > 1) & (L4_date == 0)), L4_date + current_day, L4_date)
L5_date = Con(((Dev_L4 > 1) & (L5_date == 0)), L5_date + current_day, L5_date)
L6_date = Con(((Dev_L5 > 1) & (L6_date == 0)), L6_date + current_day, L6_date)
pupa_male_date = Con(((Dev_L6_male > 1) & (pupa_male_date == 0)), pupa_male_date + current_day, pupa_male_date)
adult_male_date = Con(((Dev_pupa_male > 1) & (adult_male_date == 0)), adult_male_date + current_day, adult_male_date)
pupa_female_date = Con(((Dev_L6_female > 1) & (pupa_female_date == 0)), pupa_female_date + current_day, pupa_female_date)
adult_female_date = Con(((Dev_pupa_female > 1) & (adult_female_date == 0)), adult_female_date + current_day, adult_female_date)

current_day = current_day + 1
L2_date.save("C:\Users\Robert\Documents\PCIC_ClimateData\TEST\" + GCM + "_" + RCP + "_" + "Budworm_" + str(BudwormCount) + "_L2_" + str(year) + ".tif")
L3_date.save("C:\Users\Robert\Documents\PCIC_ClimateData\TEST\Budworm_" + str(BudwormCount) + "_L3_" + str(year) + ".tif")
L4_date.save("C:\Users\Robert\Documents\PCIC_ClimateData\TEST\Budworm_" + str(BudwormCount) + "_L4_" + str(year) + ".tif")
L5_date.save("C:\Users\Robert\Documents\PCIC_ClimateData\TEST\Budworm_" + str(BudwormCount) + "_L5_" + str(year) + ".tif")
L6_date.save("C:\Users\Robert\Documents\PCIC_ClimateData\TEST\Budworm_" + str(BudwormCount) + "_L6_" + str(year) + ".tif")
pupa_male_date.save("C:\Users\Robert\Documents\PCIC_ClimateData\TEST\Budworm_" + str(BudwormCount) + "_pupa_male_" + str(year) + ".tif")
adult_male_date.save("C:\Users\Robert\Documents\PCIC_ClimateData\TEST\Budworm_" + str(BudwormCount) + "_adult_male_" + str(year) + ".tif")
pupa_female_date.save("C:\Users\Robert\Documents\PCIC_ClimateData\TEST\Budworm_" + str(BudwormCount) + "_pupa_female_" + str(year) + ".tif")
adult_female_date.save("C:\Users\Robert\Documents\PCIC_ClimateData\TEST\Budworm_" + str(BudwormCount) + "_adult_female_" + str(year) + ".tif")

year = year + 1
BudwormCount = BudwormCount + 1
print GCM + " " + RCP + " has finished"
except Exception as e:
print "An error occurred"
print e

except arcpy.ExecuteError:
print "A geoprocessing error occurred"
print arcpy.GetMessages(2)

except IndexError:
print 'Error on line {}'.format(sys.exc_info()[-1].tb_lineno)

arcpy.CheckInExtension('spatial')









share|improve this question
























  • Cross-posted as gis.stackexchange.com/q/229266/115
    – PolyGeo
    Feb 21 '17 at 12:44















up vote
7
down vote

favorite












I am working on a project on how climate change may affect the development of spruce budworm larvae throughout each year (historical and projected) using different climate change models and emission scenarios.



The code that I have works, but would take a ridiculous amount of time to run: 38 years when extrapolated out. I'm using python 2.7.12 on a 64-bit windows machine with 24GB RAM and xeon processors. So, I don't think machine specs are an issue... it's how my script is written and organized.



I'm fairly new to Python, and have read through some articles on NumPy arrays and scipy packages, but do not have the knowledge/skills on how to implement them in this script to optimize the raster processing. Is numpy or scipy the way to go, and can those packages accomplish what I am trying to do?



Each time a new 4-h timestep raster is accessed, 27 conditional statements are executed, which probably is the reason the script runs so slow. Does anyone see a way to reduce the amount of conditional statements, but achieve the same goal?



Any help on how to improve the overall performance so that the script runs in a reasonable amount of time would be greatly appreciated!



import arcpy
from arcpy.sa import *
import math

arcpy.CheckOutExtension('spatial')
arcpy.env.overwriteOutput = True

ClimateModels = ["CanESM2", "CSIRO-Mk3-6-0", "HadGEM2-ES"]
ConcentrationPathway = ["rcp26", "rcp45", "rcp85"]

#100 random numbers between 0.4 and 2.5, representing the stochastic growth rates of individual spruce budworm larvae ....need to be the same each time it loops through GCM+RCP combination
BudwormList = [2.1, 2.3, 0.6, 2.1, 1.7, 1.1, 1.6, 2.4, 2.1, 2.2, 0.9, 1.6, 0.6, 2.0, 1.7, 0.6, 0.7, 1.8, 1.3, 2.4,
2.1, 0.8, 2.2, 1.6, 0.5, 1.0, 1.6, 2.4, 2.5, 2.5, 0.8, 0.5, 2.5, 0.9, 1.2, 0.7, 1.8, 1.8, 0.4, 1.9,
0.8, 1.3, 0.6, 2.2, 0.7, 1.6, 1.0, 0.7, 0.8, 2.4, 1.7, 0.5, 1.1, 2.2, 0.9, 1.6, 1.8, 0.8, 2.3, 0.9,
0.8, 1.0, 0.5, 1.6, 2.5, 2.2, 2.0, 1.3, 1.6, 1.2, 2.4, 0.9, 1.1, 0.8, 2.1, 2.2, 1.5, 2.4, 2.5, 2.2,
1.3, 0.4, 1.5, 1.3, 1.7, 1.9, 1.8, 0.4, 0.4, 0.5, 1.5, 0.5, 0.6, 0.7, 1.6, 0.9, 0.9, 2.5, 1.4, 1.3]
TimeStep = ['4h', '8h', '12h', '16h', '20h', '24h']

try:
for GCM in ClimateModels:
for RCP in ConcentrationPathway:
BudwormCount = 1
for budworm in BudwormList:
year = 1971
while year <= 2070:
if ((year%400 == 0) or ((year%4 == 0) and (year%100 != 0))):
days = 366
else:
days = 365
Dev_stage = Raster("C:\Users\Robert\Documents\PCIC_ClimateData\TEST\Dev_Stage.tif")
Dev_stage_female = Raster("C:\Users\Robert\Documents\PCIC_ClimateData\TEST\Dev_Stage_female.tif")
Dev_L2o = Raster("C:\Users\Robert\Documents\PCIC_ClimateData\TEST\Dev_L2o.tif")
Dev_L2 = Raster("C:\Users\Robert\Documents\PCIC_ClimateData\TEST\Dev_L2.tif")
Dev_L3 = Raster("C:\Users\Robert\Documents\PCIC_ClimateData\TEST\Dev_L3.tif")
Dev_L4 = Raster("C:\Users\Robert\Documents\PCIC_ClimateData\TEST\Dev_L4.tif")
Dev_L5 = Raster("C:\Users\Robert\Documents\PCIC_ClimateData\TEST\Dev_L5.tif")
Dev_L6_male = Raster("C:\Users\Robert\Documents\PCIC_ClimateData\TEST\Dev_L6_male.tif")
Dev_L6_female = Raster("C:\Users\Robert\Documents\PCIC_ClimateData\TEST\Dev_L6_female.tif")
Dev_pupa_male = Raster("C:\Users\Robert\Documents\PCIC_ClimateData\TEST\Dev_pupa_male.tif")
Dev_pupa_female = Raster("C:\Users\Robert\Documents\PCIC_ClimateData\TEST\Dev_pupa_female.tif")
L2_date = Raster("C:\Users\Robert\Documents\PCIC_ClimateData\TEST\Jul_L2.tif")
L3_date = Raster("C:\Users\Robert\Documents\PCIC_ClimateData\TEST\Jul_L3.tif")
L4_date = Raster("C:\Users\Robert\Documents\PCIC_ClimateData\TEST\Jul_L4.tif")
L5_date = Raster("C:\Users\Robert\Documents\PCIC_ClimateData\TEST\Jul_L5.tif")
L6_date = Raster("C:\Users\Robert\Documents\PCIC_ClimateData\TEST\Jul_L6.tif")
pupa_male_date = Raster("C:\Users\Robert\Documents\PCIC_ClimateData\TEST\Jul_pupa_male.tif")
pupa_female_date = Raster("C:\Users\Robert\Documents\PCIC_ClimateData\TEST\Jul_pupa_female.tif")
adult_male_date = Raster("C:\Users\Robert\Documents\PCIC_ClimateData\TEST\Jul_adult_male.tif")
adult_female_date = Raster("C:\Users\Robert\Documents\PCIC_ClimateData\TEST\Jul_adult_female.tif")

current_day = 1
while current_day <= days:
for time in TimeStep:
tempfile = Raster("C:\Users\Robert\Documents\PCIC_ClimateData\" + GCM + "\" + RCP + "\4hTempTimeStep\" + time + "_" + str(current_day) + "_" + str(year) + ".tif")
Dev_L2o = Con(((Dev_stage == 1) & ((tempfile > 2.5) & (tempfile < 35)) & (Dev_L2o <= 1)), Dev_L2o + (float(budworm))*((0.194*((1/(1+(math.e**(3.00-(5.84*((tempfile-2.5)/(35-2.5)))))))-(math.e**((((tempfile-2.5)/(35-2.5))-1)/0.034))))/6), Dev_L2o)
Dev_L2 = Con(((Dev_stage == 2) & ((tempfile > 4.4) & (tempfile < 38)) & (Dev_L2 <= 1)), Dev_L2 + (float(budworm))*((0.919*((1/(1+(math.e**(2.91-(5.32*((tempfile-4.4)/(38-4.4)))))))-(math.e**((((tempfile-4.4)/(38-4.4))-1)/0.061))))/6), Dev_L2)
Dev_L3 = Con(((Dev_stage == 3) & ((tempfile > 4.4) & (tempfile < 38)) & (Dev_L3 <= 1)), Dev_L3 + (float(budworm))*((0.438*((1/(1+(math.e**(3.06-(6.85*((tempfile-4.4)/(38-4.4)))))))-(math.e**((((tempfile-4.4)/(38-4.4))-1)/0.061))))/6), Dev_L3)
Dev_L4 = Con(((Dev_stage == 4) & ((tempfile > 4.4) & (tempfile < 38)) & (Dev_L4 <= 1)), Dev_L4 + (float(budworm))*((1.211*((1/(1+(math.e**(3.80-(7.55*((tempfile-4.4)/(38-4.4)))))))-(math.e**((((tempfile-4.4)/(38-4.4))-1)/0.148))))/6), Dev_L4)
Dev_L5 = Con(((Dev_stage == 5) & ((tempfile > 4.4) & (tempfile < 38)) & (Dev_L5 <= 1)), Dev_L5 + (float(budworm))*((0.269*((1/(1+(math.e**(3.02-(8.57*((tempfile-4.4)/(38-4.4)))))))-(math.e**((((tempfile-4.4)/(38-4.4))-1)/0.005))))/6), Dev_L5)
Dev_L6_male = Con(((Dev_stage == 6) & ((tempfile > 4.4) & (tempfile < 38)) & (Dev_L6_male <= 1)), Dev_L6_male + (float(budworm))*((0.288*((1/(1+(math.e**(2.67-(5.03*((tempfile-4.4)/(38-4.4)))))))-(math.e**((((tempfile-4.4)/(38-4.4))-1)/0.151))))/6), Dev_L6_male)
Dev_L6_female = Con(((Dev_stage_female == 1) & ((tempfile > 4.4) & (tempfile < 38)) & (Dev_L6_female <= 1)), Dev_L6_female + (float(budworm))*((0.317*((1/(1+(math.e**(3.06-(4.66*((tempfile-4.4)/(38-4.4)))))))-(math.e**((((tempfile-4.4)/(38-4.4))-1)/0.136))))/6), Dev_L6_female)
Dev_pupa_male = Con(((Dev_stage == 7) & ((tempfile > 4.4) & (tempfile < 35)) & (Dev_pupa_male <= 1)), Dev_pupa_male + ((0.259*((1/(1+(math.e**(2.75-(4.66*((tempfile-4.4)/(35-4.4)))))))-(math.e**((((tempfile-4.4)/(35-4.4))-1)/0.053))))/6), Dev_pupa_male)
Dev_pupa_female = Con(((Dev_stage_female == 2) & ((tempfile > 4.4) & (tempfile < 35)) & (Dev_pupa_female <= 1)), Dev_pupa_female + ((0.205*((1/(1+(math.e**(2.85-(6.28*((tempfile-4.4)/(35-4.4)))))))-(math.e**((((tempfile-4.4)/(35-4.4))-1)/0.044))))/6), Dev_pupa_female)

Dev_stage = Con(((Dev_stage == 1) & (Dev_L2o > 1)), Dev_stage + 1, Dev_stage)
Dev_stage = Con(((Dev_stage == 2) & (Dev_L2 > 1)), Dev_stage + 1, Dev_stage)
Dev_stage = Con(((Dev_stage == 3) & (Dev_L3 > 1)), Dev_stage + 1, Dev_stage)
Dev_stage = Con(((Dev_stage == 4) & (Dev_L4 > 1)), Dev_stage + 1, Dev_stage)
Dev_stage = Con(((Dev_stage == 5) & (Dev_L5 > 1)), Dev_stage + 1, Dev_stage)
Dev_stage = Con(((Dev_stage == 6) & (Dev_L6_male > 1)), Dev_stage + 1, Dev_stage)
Dev_stage = Con(((Dev_stage == 7) & (Dev_pupa_male > 1)), Dev_stage + 1, Dev_stage)
Dev_stage_female = Con(((Dev_stage_female == 1) & (Dev_L6_female > 1)), Dev_stage_female + 1, Dev_stage_female)
Dev_stage_female = Con(((Dev_stage_female == 2) & (Dev_pupa_female > 1)), Dev_stage_female + 1, Dev_stage_female)

L2_date = Con(((Dev_L2o > 1) & (L2_date == 0)), L2_date + current_day, L2_date)
L3_date = Con(((Dev_L2 > 1) & (L3_date == 0)), L3_date + current_day, L3_date)
L4_date = Con(((Dev_L3 > 1) & (L4_date == 0)), L4_date + current_day, L4_date)
L5_date = Con(((Dev_L4 > 1) & (L5_date == 0)), L5_date + current_day, L5_date)
L6_date = Con(((Dev_L5 > 1) & (L6_date == 0)), L6_date + current_day, L6_date)
pupa_male_date = Con(((Dev_L6_male > 1) & (pupa_male_date == 0)), pupa_male_date + current_day, pupa_male_date)
adult_male_date = Con(((Dev_pupa_male > 1) & (adult_male_date == 0)), adult_male_date + current_day, adult_male_date)
pupa_female_date = Con(((Dev_L6_female > 1) & (pupa_female_date == 0)), pupa_female_date + current_day, pupa_female_date)
adult_female_date = Con(((Dev_pupa_female > 1) & (adult_female_date == 0)), adult_female_date + current_day, adult_female_date)

current_day = current_day + 1
L2_date.save("C:\Users\Robert\Documents\PCIC_ClimateData\TEST\" + GCM + "_" + RCP + "_" + "Budworm_" + str(BudwormCount) + "_L2_" + str(year) + ".tif")
L3_date.save("C:\Users\Robert\Documents\PCIC_ClimateData\TEST\Budworm_" + str(BudwormCount) + "_L3_" + str(year) + ".tif")
L4_date.save("C:\Users\Robert\Documents\PCIC_ClimateData\TEST\Budworm_" + str(BudwormCount) + "_L4_" + str(year) + ".tif")
L5_date.save("C:\Users\Robert\Documents\PCIC_ClimateData\TEST\Budworm_" + str(BudwormCount) + "_L5_" + str(year) + ".tif")
L6_date.save("C:\Users\Robert\Documents\PCIC_ClimateData\TEST\Budworm_" + str(BudwormCount) + "_L6_" + str(year) + ".tif")
pupa_male_date.save("C:\Users\Robert\Documents\PCIC_ClimateData\TEST\Budworm_" + str(BudwormCount) + "_pupa_male_" + str(year) + ".tif")
adult_male_date.save("C:\Users\Robert\Documents\PCIC_ClimateData\TEST\Budworm_" + str(BudwormCount) + "_adult_male_" + str(year) + ".tif")
pupa_female_date.save("C:\Users\Robert\Documents\PCIC_ClimateData\TEST\Budworm_" + str(BudwormCount) + "_pupa_female_" + str(year) + ".tif")
adult_female_date.save("C:\Users\Robert\Documents\PCIC_ClimateData\TEST\Budworm_" + str(BudwormCount) + "_adult_female_" + str(year) + ".tif")

year = year + 1
BudwormCount = BudwormCount + 1
print GCM + " " + RCP + " has finished"
except Exception as e:
print "An error occurred"
print e

except arcpy.ExecuteError:
print "A geoprocessing error occurred"
print arcpy.GetMessages(2)

except IndexError:
print 'Error on line {}'.format(sys.exc_info()[-1].tb_lineno)

arcpy.CheckInExtension('spatial')









share|improve this question
























  • Cross-posted as gis.stackexchange.com/q/229266/115
    – PolyGeo
    Feb 21 '17 at 12:44













up vote
7
down vote

favorite









up vote
7
down vote

favorite











I am working on a project on how climate change may affect the development of spruce budworm larvae throughout each year (historical and projected) using different climate change models and emission scenarios.



The code that I have works, but would take a ridiculous amount of time to run: 38 years when extrapolated out. I'm using python 2.7.12 on a 64-bit windows machine with 24GB RAM and xeon processors. So, I don't think machine specs are an issue... it's how my script is written and organized.



I'm fairly new to Python, and have read through some articles on NumPy arrays and scipy packages, but do not have the knowledge/skills on how to implement them in this script to optimize the raster processing. Is numpy or scipy the way to go, and can those packages accomplish what I am trying to do?



Each time a new 4-h timestep raster is accessed, 27 conditional statements are executed, which probably is the reason the script runs so slow. Does anyone see a way to reduce the amount of conditional statements, but achieve the same goal?



Any help on how to improve the overall performance so that the script runs in a reasonable amount of time would be greatly appreciated!



import arcpy
from arcpy.sa import *
import math

arcpy.CheckOutExtension('spatial')
arcpy.env.overwriteOutput = True

ClimateModels = ["CanESM2", "CSIRO-Mk3-6-0", "HadGEM2-ES"]
ConcentrationPathway = ["rcp26", "rcp45", "rcp85"]

#100 random numbers between 0.4 and 2.5, representing the stochastic growth rates of individual spruce budworm larvae ....need to be the same each time it loops through GCM+RCP combination
BudwormList = [2.1, 2.3, 0.6, 2.1, 1.7, 1.1, 1.6, 2.4, 2.1, 2.2, 0.9, 1.6, 0.6, 2.0, 1.7, 0.6, 0.7, 1.8, 1.3, 2.4,
2.1, 0.8, 2.2, 1.6, 0.5, 1.0, 1.6, 2.4, 2.5, 2.5, 0.8, 0.5, 2.5, 0.9, 1.2, 0.7, 1.8, 1.8, 0.4, 1.9,
0.8, 1.3, 0.6, 2.2, 0.7, 1.6, 1.0, 0.7, 0.8, 2.4, 1.7, 0.5, 1.1, 2.2, 0.9, 1.6, 1.8, 0.8, 2.3, 0.9,
0.8, 1.0, 0.5, 1.6, 2.5, 2.2, 2.0, 1.3, 1.6, 1.2, 2.4, 0.9, 1.1, 0.8, 2.1, 2.2, 1.5, 2.4, 2.5, 2.2,
1.3, 0.4, 1.5, 1.3, 1.7, 1.9, 1.8, 0.4, 0.4, 0.5, 1.5, 0.5, 0.6, 0.7, 1.6, 0.9, 0.9, 2.5, 1.4, 1.3]
TimeStep = ['4h', '8h', '12h', '16h', '20h', '24h']

try:
for GCM in ClimateModels:
for RCP in ConcentrationPathway:
BudwormCount = 1
for budworm in BudwormList:
year = 1971
while year <= 2070:
if ((year%400 == 0) or ((year%4 == 0) and (year%100 != 0))):
days = 366
else:
days = 365
Dev_stage = Raster("C:\Users\Robert\Documents\PCIC_ClimateData\TEST\Dev_Stage.tif")
Dev_stage_female = Raster("C:\Users\Robert\Documents\PCIC_ClimateData\TEST\Dev_Stage_female.tif")
Dev_L2o = Raster("C:\Users\Robert\Documents\PCIC_ClimateData\TEST\Dev_L2o.tif")
Dev_L2 = Raster("C:\Users\Robert\Documents\PCIC_ClimateData\TEST\Dev_L2.tif")
Dev_L3 = Raster("C:\Users\Robert\Documents\PCIC_ClimateData\TEST\Dev_L3.tif")
Dev_L4 = Raster("C:\Users\Robert\Documents\PCIC_ClimateData\TEST\Dev_L4.tif")
Dev_L5 = Raster("C:\Users\Robert\Documents\PCIC_ClimateData\TEST\Dev_L5.tif")
Dev_L6_male = Raster("C:\Users\Robert\Documents\PCIC_ClimateData\TEST\Dev_L6_male.tif")
Dev_L6_female = Raster("C:\Users\Robert\Documents\PCIC_ClimateData\TEST\Dev_L6_female.tif")
Dev_pupa_male = Raster("C:\Users\Robert\Documents\PCIC_ClimateData\TEST\Dev_pupa_male.tif")
Dev_pupa_female = Raster("C:\Users\Robert\Documents\PCIC_ClimateData\TEST\Dev_pupa_female.tif")
L2_date = Raster("C:\Users\Robert\Documents\PCIC_ClimateData\TEST\Jul_L2.tif")
L3_date = Raster("C:\Users\Robert\Documents\PCIC_ClimateData\TEST\Jul_L3.tif")
L4_date = Raster("C:\Users\Robert\Documents\PCIC_ClimateData\TEST\Jul_L4.tif")
L5_date = Raster("C:\Users\Robert\Documents\PCIC_ClimateData\TEST\Jul_L5.tif")
L6_date = Raster("C:\Users\Robert\Documents\PCIC_ClimateData\TEST\Jul_L6.tif")
pupa_male_date = Raster("C:\Users\Robert\Documents\PCIC_ClimateData\TEST\Jul_pupa_male.tif")
pupa_female_date = Raster("C:\Users\Robert\Documents\PCIC_ClimateData\TEST\Jul_pupa_female.tif")
adult_male_date = Raster("C:\Users\Robert\Documents\PCIC_ClimateData\TEST\Jul_adult_male.tif")
adult_female_date = Raster("C:\Users\Robert\Documents\PCIC_ClimateData\TEST\Jul_adult_female.tif")

current_day = 1
while current_day <= days:
for time in TimeStep:
tempfile = Raster("C:\Users\Robert\Documents\PCIC_ClimateData\" + GCM + "\" + RCP + "\4hTempTimeStep\" + time + "_" + str(current_day) + "_" + str(year) + ".tif")
Dev_L2o = Con(((Dev_stage == 1) & ((tempfile > 2.5) & (tempfile < 35)) & (Dev_L2o <= 1)), Dev_L2o + (float(budworm))*((0.194*((1/(1+(math.e**(3.00-(5.84*((tempfile-2.5)/(35-2.5)))))))-(math.e**((((tempfile-2.5)/(35-2.5))-1)/0.034))))/6), Dev_L2o)
Dev_L2 = Con(((Dev_stage == 2) & ((tempfile > 4.4) & (tempfile < 38)) & (Dev_L2 <= 1)), Dev_L2 + (float(budworm))*((0.919*((1/(1+(math.e**(2.91-(5.32*((tempfile-4.4)/(38-4.4)))))))-(math.e**((((tempfile-4.4)/(38-4.4))-1)/0.061))))/6), Dev_L2)
Dev_L3 = Con(((Dev_stage == 3) & ((tempfile > 4.4) & (tempfile < 38)) & (Dev_L3 <= 1)), Dev_L3 + (float(budworm))*((0.438*((1/(1+(math.e**(3.06-(6.85*((tempfile-4.4)/(38-4.4)))))))-(math.e**((((tempfile-4.4)/(38-4.4))-1)/0.061))))/6), Dev_L3)
Dev_L4 = Con(((Dev_stage == 4) & ((tempfile > 4.4) & (tempfile < 38)) & (Dev_L4 <= 1)), Dev_L4 + (float(budworm))*((1.211*((1/(1+(math.e**(3.80-(7.55*((tempfile-4.4)/(38-4.4)))))))-(math.e**((((tempfile-4.4)/(38-4.4))-1)/0.148))))/6), Dev_L4)
Dev_L5 = Con(((Dev_stage == 5) & ((tempfile > 4.4) & (tempfile < 38)) & (Dev_L5 <= 1)), Dev_L5 + (float(budworm))*((0.269*((1/(1+(math.e**(3.02-(8.57*((tempfile-4.4)/(38-4.4)))))))-(math.e**((((tempfile-4.4)/(38-4.4))-1)/0.005))))/6), Dev_L5)
Dev_L6_male = Con(((Dev_stage == 6) & ((tempfile > 4.4) & (tempfile < 38)) & (Dev_L6_male <= 1)), Dev_L6_male + (float(budworm))*((0.288*((1/(1+(math.e**(2.67-(5.03*((tempfile-4.4)/(38-4.4)))))))-(math.e**((((tempfile-4.4)/(38-4.4))-1)/0.151))))/6), Dev_L6_male)
Dev_L6_female = Con(((Dev_stage_female == 1) & ((tempfile > 4.4) & (tempfile < 38)) & (Dev_L6_female <= 1)), Dev_L6_female + (float(budworm))*((0.317*((1/(1+(math.e**(3.06-(4.66*((tempfile-4.4)/(38-4.4)))))))-(math.e**((((tempfile-4.4)/(38-4.4))-1)/0.136))))/6), Dev_L6_female)
Dev_pupa_male = Con(((Dev_stage == 7) & ((tempfile > 4.4) & (tempfile < 35)) & (Dev_pupa_male <= 1)), Dev_pupa_male + ((0.259*((1/(1+(math.e**(2.75-(4.66*((tempfile-4.4)/(35-4.4)))))))-(math.e**((((tempfile-4.4)/(35-4.4))-1)/0.053))))/6), Dev_pupa_male)
Dev_pupa_female = Con(((Dev_stage_female == 2) & ((tempfile > 4.4) & (tempfile < 35)) & (Dev_pupa_female <= 1)), Dev_pupa_female + ((0.205*((1/(1+(math.e**(2.85-(6.28*((tempfile-4.4)/(35-4.4)))))))-(math.e**((((tempfile-4.4)/(35-4.4))-1)/0.044))))/6), Dev_pupa_female)

Dev_stage = Con(((Dev_stage == 1) & (Dev_L2o > 1)), Dev_stage + 1, Dev_stage)
Dev_stage = Con(((Dev_stage == 2) & (Dev_L2 > 1)), Dev_stage + 1, Dev_stage)
Dev_stage = Con(((Dev_stage == 3) & (Dev_L3 > 1)), Dev_stage + 1, Dev_stage)
Dev_stage = Con(((Dev_stage == 4) & (Dev_L4 > 1)), Dev_stage + 1, Dev_stage)
Dev_stage = Con(((Dev_stage == 5) & (Dev_L5 > 1)), Dev_stage + 1, Dev_stage)
Dev_stage = Con(((Dev_stage == 6) & (Dev_L6_male > 1)), Dev_stage + 1, Dev_stage)
Dev_stage = Con(((Dev_stage == 7) & (Dev_pupa_male > 1)), Dev_stage + 1, Dev_stage)
Dev_stage_female = Con(((Dev_stage_female == 1) & (Dev_L6_female > 1)), Dev_stage_female + 1, Dev_stage_female)
Dev_stage_female = Con(((Dev_stage_female == 2) & (Dev_pupa_female > 1)), Dev_stage_female + 1, Dev_stage_female)

L2_date = Con(((Dev_L2o > 1) & (L2_date == 0)), L2_date + current_day, L2_date)
L3_date = Con(((Dev_L2 > 1) & (L3_date == 0)), L3_date + current_day, L3_date)
L4_date = Con(((Dev_L3 > 1) & (L4_date == 0)), L4_date + current_day, L4_date)
L5_date = Con(((Dev_L4 > 1) & (L5_date == 0)), L5_date + current_day, L5_date)
L6_date = Con(((Dev_L5 > 1) & (L6_date == 0)), L6_date + current_day, L6_date)
pupa_male_date = Con(((Dev_L6_male > 1) & (pupa_male_date == 0)), pupa_male_date + current_day, pupa_male_date)
adult_male_date = Con(((Dev_pupa_male > 1) & (adult_male_date == 0)), adult_male_date + current_day, adult_male_date)
pupa_female_date = Con(((Dev_L6_female > 1) & (pupa_female_date == 0)), pupa_female_date + current_day, pupa_female_date)
adult_female_date = Con(((Dev_pupa_female > 1) & (adult_female_date == 0)), adult_female_date + current_day, adult_female_date)

current_day = current_day + 1
L2_date.save("C:\Users\Robert\Documents\PCIC_ClimateData\TEST\" + GCM + "_" + RCP + "_" + "Budworm_" + str(BudwormCount) + "_L2_" + str(year) + ".tif")
L3_date.save("C:\Users\Robert\Documents\PCIC_ClimateData\TEST\Budworm_" + str(BudwormCount) + "_L3_" + str(year) + ".tif")
L4_date.save("C:\Users\Robert\Documents\PCIC_ClimateData\TEST\Budworm_" + str(BudwormCount) + "_L4_" + str(year) + ".tif")
L5_date.save("C:\Users\Robert\Documents\PCIC_ClimateData\TEST\Budworm_" + str(BudwormCount) + "_L5_" + str(year) + ".tif")
L6_date.save("C:\Users\Robert\Documents\PCIC_ClimateData\TEST\Budworm_" + str(BudwormCount) + "_L6_" + str(year) + ".tif")
pupa_male_date.save("C:\Users\Robert\Documents\PCIC_ClimateData\TEST\Budworm_" + str(BudwormCount) + "_pupa_male_" + str(year) + ".tif")
adult_male_date.save("C:\Users\Robert\Documents\PCIC_ClimateData\TEST\Budworm_" + str(BudwormCount) + "_adult_male_" + str(year) + ".tif")
pupa_female_date.save("C:\Users\Robert\Documents\PCIC_ClimateData\TEST\Budworm_" + str(BudwormCount) + "_pupa_female_" + str(year) + ".tif")
adult_female_date.save("C:\Users\Robert\Documents\PCIC_ClimateData\TEST\Budworm_" + str(BudwormCount) + "_adult_female_" + str(year) + ".tif")

year = year + 1
BudwormCount = BudwormCount + 1
print GCM + " " + RCP + " has finished"
except Exception as e:
print "An error occurred"
print e

except arcpy.ExecuteError:
print "A geoprocessing error occurred"
print arcpy.GetMessages(2)

except IndexError:
print 'Error on line {}'.format(sys.exc_info()[-1].tb_lineno)

arcpy.CheckInExtension('spatial')









share|improve this question















I am working on a project on how climate change may affect the development of spruce budworm larvae throughout each year (historical and projected) using different climate change models and emission scenarios.



The code that I have works, but would take a ridiculous amount of time to run: 38 years when extrapolated out. I'm using python 2.7.12 on a 64-bit windows machine with 24GB RAM and xeon processors. So, I don't think machine specs are an issue... it's how my script is written and organized.



I'm fairly new to Python, and have read through some articles on NumPy arrays and scipy packages, but do not have the knowledge/skills on how to implement them in this script to optimize the raster processing. Is numpy or scipy the way to go, and can those packages accomplish what I am trying to do?



Each time a new 4-h timestep raster is accessed, 27 conditional statements are executed, which probably is the reason the script runs so slow. Does anyone see a way to reduce the amount of conditional statements, but achieve the same goal?



Any help on how to improve the overall performance so that the script runs in a reasonable amount of time would be greatly appreciated!



import arcpy
from arcpy.sa import *
import math

arcpy.CheckOutExtension('spatial')
arcpy.env.overwriteOutput = True

ClimateModels = ["CanESM2", "CSIRO-Mk3-6-0", "HadGEM2-ES"]
ConcentrationPathway = ["rcp26", "rcp45", "rcp85"]

#100 random numbers between 0.4 and 2.5, representing the stochastic growth rates of individual spruce budworm larvae ....need to be the same each time it loops through GCM+RCP combination
BudwormList = [2.1, 2.3, 0.6, 2.1, 1.7, 1.1, 1.6, 2.4, 2.1, 2.2, 0.9, 1.6, 0.6, 2.0, 1.7, 0.6, 0.7, 1.8, 1.3, 2.4,
2.1, 0.8, 2.2, 1.6, 0.5, 1.0, 1.6, 2.4, 2.5, 2.5, 0.8, 0.5, 2.5, 0.9, 1.2, 0.7, 1.8, 1.8, 0.4, 1.9,
0.8, 1.3, 0.6, 2.2, 0.7, 1.6, 1.0, 0.7, 0.8, 2.4, 1.7, 0.5, 1.1, 2.2, 0.9, 1.6, 1.8, 0.8, 2.3, 0.9,
0.8, 1.0, 0.5, 1.6, 2.5, 2.2, 2.0, 1.3, 1.6, 1.2, 2.4, 0.9, 1.1, 0.8, 2.1, 2.2, 1.5, 2.4, 2.5, 2.2,
1.3, 0.4, 1.5, 1.3, 1.7, 1.9, 1.8, 0.4, 0.4, 0.5, 1.5, 0.5, 0.6, 0.7, 1.6, 0.9, 0.9, 2.5, 1.4, 1.3]
TimeStep = ['4h', '8h', '12h', '16h', '20h', '24h']

try:
for GCM in ClimateModels:
for RCP in ConcentrationPathway:
BudwormCount = 1
for budworm in BudwormList:
year = 1971
while year <= 2070:
if ((year%400 == 0) or ((year%4 == 0) and (year%100 != 0))):
days = 366
else:
days = 365
Dev_stage = Raster("C:\Users\Robert\Documents\PCIC_ClimateData\TEST\Dev_Stage.tif")
Dev_stage_female = Raster("C:\Users\Robert\Documents\PCIC_ClimateData\TEST\Dev_Stage_female.tif")
Dev_L2o = Raster("C:\Users\Robert\Documents\PCIC_ClimateData\TEST\Dev_L2o.tif")
Dev_L2 = Raster("C:\Users\Robert\Documents\PCIC_ClimateData\TEST\Dev_L2.tif")
Dev_L3 = Raster("C:\Users\Robert\Documents\PCIC_ClimateData\TEST\Dev_L3.tif")
Dev_L4 = Raster("C:\Users\Robert\Documents\PCIC_ClimateData\TEST\Dev_L4.tif")
Dev_L5 = Raster("C:\Users\Robert\Documents\PCIC_ClimateData\TEST\Dev_L5.tif")
Dev_L6_male = Raster("C:\Users\Robert\Documents\PCIC_ClimateData\TEST\Dev_L6_male.tif")
Dev_L6_female = Raster("C:\Users\Robert\Documents\PCIC_ClimateData\TEST\Dev_L6_female.tif")
Dev_pupa_male = Raster("C:\Users\Robert\Documents\PCIC_ClimateData\TEST\Dev_pupa_male.tif")
Dev_pupa_female = Raster("C:\Users\Robert\Documents\PCIC_ClimateData\TEST\Dev_pupa_female.tif")
L2_date = Raster("C:\Users\Robert\Documents\PCIC_ClimateData\TEST\Jul_L2.tif")
L3_date = Raster("C:\Users\Robert\Documents\PCIC_ClimateData\TEST\Jul_L3.tif")
L4_date = Raster("C:\Users\Robert\Documents\PCIC_ClimateData\TEST\Jul_L4.tif")
L5_date = Raster("C:\Users\Robert\Documents\PCIC_ClimateData\TEST\Jul_L5.tif")
L6_date = Raster("C:\Users\Robert\Documents\PCIC_ClimateData\TEST\Jul_L6.tif")
pupa_male_date = Raster("C:\Users\Robert\Documents\PCIC_ClimateData\TEST\Jul_pupa_male.tif")
pupa_female_date = Raster("C:\Users\Robert\Documents\PCIC_ClimateData\TEST\Jul_pupa_female.tif")
adult_male_date = Raster("C:\Users\Robert\Documents\PCIC_ClimateData\TEST\Jul_adult_male.tif")
adult_female_date = Raster("C:\Users\Robert\Documents\PCIC_ClimateData\TEST\Jul_adult_female.tif")

current_day = 1
while current_day <= days:
for time in TimeStep:
tempfile = Raster("C:\Users\Robert\Documents\PCIC_ClimateData\" + GCM + "\" + RCP + "\4hTempTimeStep\" + time + "_" + str(current_day) + "_" + str(year) + ".tif")
Dev_L2o = Con(((Dev_stage == 1) & ((tempfile > 2.5) & (tempfile < 35)) & (Dev_L2o <= 1)), Dev_L2o + (float(budworm))*((0.194*((1/(1+(math.e**(3.00-(5.84*((tempfile-2.5)/(35-2.5)))))))-(math.e**((((tempfile-2.5)/(35-2.5))-1)/0.034))))/6), Dev_L2o)
Dev_L2 = Con(((Dev_stage == 2) & ((tempfile > 4.4) & (tempfile < 38)) & (Dev_L2 <= 1)), Dev_L2 + (float(budworm))*((0.919*((1/(1+(math.e**(2.91-(5.32*((tempfile-4.4)/(38-4.4)))))))-(math.e**((((tempfile-4.4)/(38-4.4))-1)/0.061))))/6), Dev_L2)
Dev_L3 = Con(((Dev_stage == 3) & ((tempfile > 4.4) & (tempfile < 38)) & (Dev_L3 <= 1)), Dev_L3 + (float(budworm))*((0.438*((1/(1+(math.e**(3.06-(6.85*((tempfile-4.4)/(38-4.4)))))))-(math.e**((((tempfile-4.4)/(38-4.4))-1)/0.061))))/6), Dev_L3)
Dev_L4 = Con(((Dev_stage == 4) & ((tempfile > 4.4) & (tempfile < 38)) & (Dev_L4 <= 1)), Dev_L4 + (float(budworm))*((1.211*((1/(1+(math.e**(3.80-(7.55*((tempfile-4.4)/(38-4.4)))))))-(math.e**((((tempfile-4.4)/(38-4.4))-1)/0.148))))/6), Dev_L4)
Dev_L5 = Con(((Dev_stage == 5) & ((tempfile > 4.4) & (tempfile < 38)) & (Dev_L5 <= 1)), Dev_L5 + (float(budworm))*((0.269*((1/(1+(math.e**(3.02-(8.57*((tempfile-4.4)/(38-4.4)))))))-(math.e**((((tempfile-4.4)/(38-4.4))-1)/0.005))))/6), Dev_L5)
Dev_L6_male = Con(((Dev_stage == 6) & ((tempfile > 4.4) & (tempfile < 38)) & (Dev_L6_male <= 1)), Dev_L6_male + (float(budworm))*((0.288*((1/(1+(math.e**(2.67-(5.03*((tempfile-4.4)/(38-4.4)))))))-(math.e**((((tempfile-4.4)/(38-4.4))-1)/0.151))))/6), Dev_L6_male)
Dev_L6_female = Con(((Dev_stage_female == 1) & ((tempfile > 4.4) & (tempfile < 38)) & (Dev_L6_female <= 1)), Dev_L6_female + (float(budworm))*((0.317*((1/(1+(math.e**(3.06-(4.66*((tempfile-4.4)/(38-4.4)))))))-(math.e**((((tempfile-4.4)/(38-4.4))-1)/0.136))))/6), Dev_L6_female)
Dev_pupa_male = Con(((Dev_stage == 7) & ((tempfile > 4.4) & (tempfile < 35)) & (Dev_pupa_male <= 1)), Dev_pupa_male + ((0.259*((1/(1+(math.e**(2.75-(4.66*((tempfile-4.4)/(35-4.4)))))))-(math.e**((((tempfile-4.4)/(35-4.4))-1)/0.053))))/6), Dev_pupa_male)
Dev_pupa_female = Con(((Dev_stage_female == 2) & ((tempfile > 4.4) & (tempfile < 35)) & (Dev_pupa_female <= 1)), Dev_pupa_female + ((0.205*((1/(1+(math.e**(2.85-(6.28*((tempfile-4.4)/(35-4.4)))))))-(math.e**((((tempfile-4.4)/(35-4.4))-1)/0.044))))/6), Dev_pupa_female)

Dev_stage = Con(((Dev_stage == 1) & (Dev_L2o > 1)), Dev_stage + 1, Dev_stage)
Dev_stage = Con(((Dev_stage == 2) & (Dev_L2 > 1)), Dev_stage + 1, Dev_stage)
Dev_stage = Con(((Dev_stage == 3) & (Dev_L3 > 1)), Dev_stage + 1, Dev_stage)
Dev_stage = Con(((Dev_stage == 4) & (Dev_L4 > 1)), Dev_stage + 1, Dev_stage)
Dev_stage = Con(((Dev_stage == 5) & (Dev_L5 > 1)), Dev_stage + 1, Dev_stage)
Dev_stage = Con(((Dev_stage == 6) & (Dev_L6_male > 1)), Dev_stage + 1, Dev_stage)
Dev_stage = Con(((Dev_stage == 7) & (Dev_pupa_male > 1)), Dev_stage + 1, Dev_stage)
Dev_stage_female = Con(((Dev_stage_female == 1) & (Dev_L6_female > 1)), Dev_stage_female + 1, Dev_stage_female)
Dev_stage_female = Con(((Dev_stage_female == 2) & (Dev_pupa_female > 1)), Dev_stage_female + 1, Dev_stage_female)

L2_date = Con(((Dev_L2o > 1) & (L2_date == 0)), L2_date + current_day, L2_date)
L3_date = Con(((Dev_L2 > 1) & (L3_date == 0)), L3_date + current_day, L3_date)
L4_date = Con(((Dev_L3 > 1) & (L4_date == 0)), L4_date + current_day, L4_date)
L5_date = Con(((Dev_L4 > 1) & (L5_date == 0)), L5_date + current_day, L5_date)
L6_date = Con(((Dev_L5 > 1) & (L6_date == 0)), L6_date + current_day, L6_date)
pupa_male_date = Con(((Dev_L6_male > 1) & (pupa_male_date == 0)), pupa_male_date + current_day, pupa_male_date)
adult_male_date = Con(((Dev_pupa_male > 1) & (adult_male_date == 0)), adult_male_date + current_day, adult_male_date)
pupa_female_date = Con(((Dev_L6_female > 1) & (pupa_female_date == 0)), pupa_female_date + current_day, pupa_female_date)
adult_female_date = Con(((Dev_pupa_female > 1) & (adult_female_date == 0)), adult_female_date + current_day, adult_female_date)

current_day = current_day + 1
L2_date.save("C:\Users\Robert\Documents\PCIC_ClimateData\TEST\" + GCM + "_" + RCP + "_" + "Budworm_" + str(BudwormCount) + "_L2_" + str(year) + ".tif")
L3_date.save("C:\Users\Robert\Documents\PCIC_ClimateData\TEST\Budworm_" + str(BudwormCount) + "_L3_" + str(year) + ".tif")
L4_date.save("C:\Users\Robert\Documents\PCIC_ClimateData\TEST\Budworm_" + str(BudwormCount) + "_L4_" + str(year) + ".tif")
L5_date.save("C:\Users\Robert\Documents\PCIC_ClimateData\TEST\Budworm_" + str(BudwormCount) + "_L5_" + str(year) + ".tif")
L6_date.save("C:\Users\Robert\Documents\PCIC_ClimateData\TEST\Budworm_" + str(BudwormCount) + "_L6_" + str(year) + ".tif")
pupa_male_date.save("C:\Users\Robert\Documents\PCIC_ClimateData\TEST\Budworm_" + str(BudwormCount) + "_pupa_male_" + str(year) + ".tif")
adult_male_date.save("C:\Users\Robert\Documents\PCIC_ClimateData\TEST\Budworm_" + str(BudwormCount) + "_adult_male_" + str(year) + ".tif")
pupa_female_date.save("C:\Users\Robert\Documents\PCIC_ClimateData\TEST\Budworm_" + str(BudwormCount) + "_pupa_female_" + str(year) + ".tif")
adult_female_date.save("C:\Users\Robert\Documents\PCIC_ClimateData\TEST\Budworm_" + str(BudwormCount) + "_adult_female_" + str(year) + ".tif")

year = year + 1
BudwormCount = BudwormCount + 1
print GCM + " " + RCP + " has finished"
except Exception as e:
print "An error occurred"
print e

except arcpy.ExecuteError:
print "A geoprocessing error occurred"
print arcpy.GetMessages(2)

except IndexError:
print 'Error on line {}'.format(sys.exc_info()[-1].tb_lineno)

arcpy.CheckInExtension('spatial')






python python-2.x geospatial scipy arcpy






share|improve this question















share|improve this question













share|improve this question




share|improve this question








edited Nov 16 at 13:01









rolfl

90.6k13190393




90.6k13190393










asked Feb 17 '17 at 16:30









bobby12345

363




363












  • Cross-posted as gis.stackexchange.com/q/229266/115
    – PolyGeo
    Feb 21 '17 at 12:44


















  • Cross-posted as gis.stackexchange.com/q/229266/115
    – PolyGeo
    Feb 21 '17 at 12:44
















Cross-posted as gis.stackexchange.com/q/229266/115
– PolyGeo
Feb 21 '17 at 12:44




Cross-posted as gis.stackexchange.com/q/229266/115
– PolyGeo
Feb 21 '17 at 12:44










1 Answer
1






active

oldest

votes

















up vote
3
down vote













First, fix the bug!



You're using sys.exc_info()[-1].tb_lineno but you're not importing the sys module anywhere, so do that first !



Now, I'll start my review with a couple of style guides. You can read more about them here.




  • use snake_case convention for your variable names


  • imports should be grouped in the following order:




    • standard library imports

    • related third party imports

    • local application/library specific imports




You should put a blank line between each group of imports. More, absolute imports are recommended, as they are usually more readable and tend to be better behaved (or at least give better error messages) if the import system is incorrectly configured




  • don't add redundant parentheses in conditional statements

  • use the print() function even if you're using Python 2.7.12

  • add a space before and after each operator

  • try to keep your lines no longer than 120 characters


  • use string formatting when you're concatenating strings: e.g:



    print("Skiwi hates stacks {} times more than I do!".format(100))



With all of the above in mind, we'll have the following code (the formatting is a bit ugly because you have really long formulas but heh..):



import math
import sys

import arcpy
from arcpy.sa import Con, Raster

arcpy.CheckOutExtension('spatial')
arcpy.env.overwriteOutput = True

climate_models = ["CanESM2", "CSIRO-Mk3-6-0", "HadGEM2-ES"]
concentration_pathway = ["rcp26", "rcp45", "rcp85"]

budworm_list = [2.1, 2.3, 0.6, 2.1, 1.7, 1.1, 1.6, 2.4, 2.1, 2.2, 0.9, 1.6, 0.6, 2.0, 1.7, 0.6, 0.7, 1.8, 1.3, 2.4,
2.1, 0.8, 2.2, 1.6, 0.5, 1.0, 1.6, 2.4, 2.5, 2.5, 0.8, 0.5, 2.5, 0.9, 1.2, 0.7, 1.8, 1.8, 0.4, 1.9,
0.8, 1.3, 0.6, 2.2, 0.7, 1.6, 1.0, 0.7, 0.8, 2.4, 1.7, 0.5, 1.1, 2.2, 0.9, 1.6, 1.8, 0.8, 2.3, 0.9,
0.8, 1.0, 0.5, 1.6, 2.5, 2.2, 2.0, 1.3, 1.6, 1.2, 2.4, 0.9, 1.1, 0.8, 2.1, 2.2, 1.5, 2.4, 2.5, 2.2,
1.3, 0.4, 1.5, 1.3, 1.7, 1.9, 1.8, 0.4, 0.4, 0.5, 1.5, 0.5, 0.6, 0.7, 1.6, 0.9, 0.9, 2.5, 1.4, 1.3]
time_step = ['4h', '8h', '12h', '16h', '20h', '24h']

try:
for gcm in climate_models:
for rcp in concentration_pathway:
budworm_count = 1
for budworm in budworm_list:
year = 1971
while year <= 2070:
if (year % 400 == 0) or ((year % 4 == 0) and (year % 100 != 0)):
days = 366
else:
days = 365
dev_stage = Raster("C:\Users\Robert\Documents\PCIC_ClimateData\TEST\Dev_Stage.tif")
dev_stage_female = Raster(
"C:\Users\Robert\Documents\PCIC_ClimateData\TEST\Dev_Stage_female.tif")
dev_l2o = Raster("C:\Users\Robert\Documents\PCIC_ClimateData\TEST\Dev_L2o.tif")
dev_l2 = Raster("C:\Users\Robert\Documents\PCIC_ClimateData\TEST\Dev_L2.tif")
dev_l3 = Raster("C:\Users\Robert\Documents\PCIC_ClimateData\TEST\Dev_L3.tif")
dev_l4 = Raster("C:\Users\Robert\Documents\PCIC_ClimateData\TEST\Dev_L4.tif")
dev_l5 = Raster("C:\Users\Robert\Documents\PCIC_ClimateData\TEST\Dev_L5.tif")
dev_l6_male = Raster("C:\Users\Robert\Documents\PCIC_ClimateData\TEST\Dev_L6_male.tif")
dev_l6_female = Raster("C:\Users\Robert\Documents\PCIC_ClimateData\TEST\Dev_L6_female.tif")
dev_pupa_male = Raster("C:\Users\Robert\Documents\PCIC_ClimateData\TEST\Dev_pupa_male.tif")
dev_pupa_female = Raster(
"C:\Users\Robert\Documents\PCIC_ClimateData\TEST\Dev_pupa_female.tif")
l2_date = Raster("C:\Users\Robert\Documents\PCIC_ClimateData\TEST\Jul_L2.tif")
l3_date = Raster("C:\Users\Robert\Documents\PCIC_ClimateData\TEST\Jul_L3.tif")
l4_date = Raster("C:\Users\Robert\Documents\PCIC_ClimateData\TEST\Jul_L4.tif")
l5_date = Raster("C:\Users\Robert\Documents\PCIC_ClimateData\TEST\Jul_L5.tif")
l6_date = Raster("C:\Users\Robert\Documents\PCIC_ClimateData\TEST\Jul_L6.tif")
pupa_male_date = Raster("C:\Users\Robert\Documents\PCIC_ClimateData\TEST\Jul_pupa_male.tif")
pupa_female_date = Raster(
"C:\Users\Robert\Documents\PCIC_ClimateData\TEST\Jul_pupa_female.tif")
adult_male_date = Raster("C:\Users\Robert\Documents\PCIC_ClimateData\TEST\Jul_adult_male.tif")
adult_female_date = Raster(
"C:\Users\Robert\Documents\PCIC_ClimateData\TEST\Jul_adult_female.tif")

current_day = 1
while current_day <= days:
for time in time_step:
tempfile = Raster(
"C:\Users\Robert\Documents\PCIC_ClimateData\{}\{}\4hTempTimeStep\{}_{}_{}.tif".format(
gcm, rcp, time, current_day, year))
dev_l2o = Con(((dev_stage == 1) & ((tempfile > 2.5) & (tempfile < 35)) & (dev_l2o <= 1)),
dev_l2o + (float(budworm)) * ((0.194 * (
(1 / (1 + (math.e ** (3.00 - (5.84 * ((tempfile - 2.5) / (35 - 2.5))))))) - (
math.e ** ((((tempfile - 2.5) / (35 - 2.5)) - 1) / 0.034)))) / 6), dev_l2o)
dev_l2 = Con(((dev_stage == 2) & ((tempfile > 4.4) & (tempfile < 38)) & (dev_l2 <= 1)),
dev_l2 + (float(budworm)) * ((0.919 * (
(1 / (1 + (math.e ** (2.91 - (5.32 * ((tempfile - 4.4) / (38 - 4.4))))))) - (
math.e ** ((((tempfile - 4.4) / (38 - 4.4)) - 1) / 0.061)))) / 6), dev_l2)
dev_l3 = Con(((dev_stage == 3) & ((tempfile > 4.4) & (tempfile < 38)) & (dev_l3 <= 1)),
dev_l3 + (float(budworm)) * ((0.438 * (
(1 / (1 + (math.e ** (3.06 - (6.85 * ((tempfile - 4.4) / (38 - 4.4))))))) - (
math.e ** ((((tempfile - 4.4) / (38 - 4.4)) - 1) / 0.061)))) / 6), dev_l3)
dev_l4 = Con(((dev_stage == 4) & ((tempfile > 4.4) & (tempfile < 38)) & (dev_l4 <= 1)),
dev_l4 + (float(budworm)) * ((1.211 * (
(1 / (1 + (math.e ** (3.80 - (7.55 * ((tempfile - 4.4) / (38 - 4.4))))))) - (
math.e ** ((((tempfile - 4.4) / (38 - 4.4)) - 1) / 0.148)))) / 6), dev_l4)
dev_l5 = Con(((dev_stage == 5) & ((tempfile > 4.4) & (tempfile < 38)) & (dev_l5 <= 1)),
dev_l5 + (float(budworm)) * ((0.269 * (
(1 / (1 + (math.e ** (3.02 - (8.57 * ((tempfile - 4.4) / (38 - 4.4))))))) - (
math.e ** ((((tempfile - 4.4) / (38 - 4.4)) - 1) / 0.005)))) / 6), dev_l5)
dev_l6_male = Con(
((dev_stage == 6) & ((tempfile > 4.4) & (tempfile < 38)) & (dev_l6_male <= 1)),
dev_l6_male + (float(budworm)) * ((0.288 * (
(1 / (1 + (math.e ** (2.67 - (5.03 * ((tempfile - 4.4) / (38 - 4.4))))))) - (
math.e ** ((((tempfile - 4.4) / (38 - 4.4)) - 1) / 0.151)))) / 6), dev_l6_male)
dev_l6_female = Con(
((dev_stage_female == 1) & ((tempfile > 4.4) & (tempfile < 38)) & (dev_l6_female <= 1)),
dev_l6_female + (float(budworm)) * ((0.317 * (
(1 / (1 + (math.e ** (3.06 - (4.66 * ((tempfile - 4.4) / (38 - 4.4))))))) - (
math.e ** ((((tempfile - 4.4) / (38 - 4.4)) - 1) / 0.136)))) / 6), dev_l6_female)
dev_pupa_male = Con(
((dev_stage == 7) & ((tempfile > 4.4) & (tempfile < 35)) & (dev_pupa_male <= 1)),
dev_pupa_male + ((0.259 * (
(1 / (1 + (math.e ** (2.75 - (4.66 * ((tempfile - 4.4) / (35 - 4.4))))))) - (
math.e ** ((((tempfile - 4.4) / (35 - 4.4)) - 1) / 0.053)))) / 6), dev_pupa_male)
dev_pupa_female = Con(((dev_stage_female == 2) & ((tempfile > 4.4) & (tempfile < 35)) & (
dev_pupa_female <= 1)), dev_pupa_female + ((0.205 * (
(1 / (1 + (math.e ** (2.85 - (6.28 * ((tempfile - 4.4) / (35 - 4.4))))))) - (
math.e ** ((((tempfile - 4.4) / (35 - 4.4)) - 1) / 0.044)))) / 6), dev_pupa_female)

dev_stage = Con(((dev_stage == 1) & (dev_l2o > 1)), dev_stage + 1, dev_stage)
dev_stage = Con(((dev_stage == 2) & (dev_l2 > 1)), dev_stage + 1, dev_stage)
dev_stage = Con(((dev_stage == 3) & (dev_l3 > 1)), dev_stage + 1, dev_stage)
dev_stage = Con(((dev_stage == 4) & (dev_l4 > 1)), dev_stage + 1, dev_stage)
dev_stage = Con(((dev_stage == 5) & (dev_l5 > 1)), dev_stage + 1, dev_stage)
dev_stage = Con(((dev_stage == 6) & (dev_l6_male > 1)), dev_stage + 1, dev_stage)
dev_stage = Con(((dev_stage == 7) & (dev_pupa_male > 1)), dev_stage + 1, dev_stage)
dev_stage_female = Con(((dev_stage_female == 1) & (dev_l6_female > 1)),
dev_stage_female + 1, dev_stage_female)
dev_stage_female = Con(((dev_stage_female == 2) & (dev_pupa_female > 1)),
dev_stage_female + 1, dev_stage_female)

l2_date = Con(((dev_l2o > 1) & (l2_date == 0)), l2_date + current_day, l2_date)
l3_date = Con(((dev_l2 > 1) & (l3_date == 0)), l3_date + current_day, l3_date)
l4_date = Con(((dev_l3 > 1) & (l4_date == 0)), l4_date + current_day, l4_date)
l5_date = Con(((dev_l4 > 1) & (l5_date == 0)), l5_date + current_day, l5_date)
l6_date = Con(((dev_l5 > 1) & (l6_date == 0)), l6_date + current_day, l6_date)
pupa_male_date = Con(((dev_l6_male > 1) & (pupa_male_date == 0)),
pupa_male_date + current_day, pupa_male_date)
adult_male_date = Con(((dev_pupa_male > 1) & (adult_male_date == 0)),
adult_male_date + current_day, adult_male_date)
pupa_female_date = Con(((dev_l6_female > 1) & (pupa_female_date == 0)),
pupa_female_date + current_day, pupa_female_date)
adult_female_date = Con(((dev_pupa_female > 1) & (adult_female_date == 0)),
adult_female_date + current_day, adult_female_date)

current_day = current_day + 1
l2_date.save(
"C:\Users\Robert\Documents\PCIC_ClimateData\TEST\{}_{}_Budworm_{}_L2_{}.tif".format(gcm,
rcp,
budworm_count,
year))
l3_date.save("C:\Users\Robert\Documents\PCIC_ClimateData\TEST\Budworm_{}_L3_{}.tif".format(
budworm_count, year))
l4_date.save("C:\Users\Robert\Documents\PCIC_ClimateData\TEST\Budworm_{}_L4_{}.tif".format(
budworm_count, year))
l5_date.save("C:\Users\Robert\Documents\PCIC_ClimateData\TEST\Budworm_{}_L5_{}.tif".format(
budworm_count, year))
l6_date.save("C:\Users\Robert\Documents\PCIC_ClimateData\TEST\Budworm_{}_L6_{}.tif".format(
budworm_count, year))
pupa_male_date.save(
"C:\Users\Robert\Documents\PCIC_ClimateData\TEST\Budworm_{}_pupa_male_{}.tif".format(
budworm_count, year))
adult_male_date.save(
"C:\Users\Robert\Documents\PCIC_ClimateData\TEST\Budworm_{}_adult_male_{}.tif".format(
budworm_count, year))
pupa_female_date.save(
"C:\Users\Robert\Documents\PCIC_ClimateData\TEST\Budworm_{}_pupa_female_{}.tif".format(
budworm_count, year))
adult_female_date.save(
"C:\Users\Robert\Documents\PCIC_ClimateData\TEST\Budworm_{}_adult_female_{}.tif".format(
budworm_count, year))

year += 1
budworm_count += 1
print("{} {} has finished".format(gcm, rcp))
except Exception as e:
print("An error occurred")
print(e)

except arcpy.ExecuteError:
print("A geoprocessing error occurred")
print(arcpy.GetMessages(2))

except IndexError:
print('Error on line {}'.format(sys.exc_info()[-1].tb_lineno))

arcpy.CheckInExtension('spatial')


Code related improvements





  • you said you want to generate a list of 100 random numbers (which should be the same while looping). That's quite easy to do if you're using the random module:



    def budworm_list(lower_limit, upper_limit, decimals):
    return [round(random.uniform(lower_limit, upper_limit), decimals) for _ in range(100)]


    Which can be latter called like this:



    budworm_list = generate_budworm_list(0.4, 2.5, 1)



  • use enumerate() to get rid of the usual counter. This Function can also take as a parameter the step which will let you start from 1.



    ...
    for budworm_count, budworm in enumerate(budworm_list, start=1):
    ...


  • don't initialize the year each time you're looping. Instead, declare it at the top of the program. More, that's a constant, which is usually declared in UPPER_CASE. (YEAR = 1971). Better yet, remove completely the while loop and add a for loop (see below).



  • there's a built-in which can tell you if an year is leap or not. You can find it in the calendar module:



    import calendar

    print calendar.isleap(1900)
    >> True


    So, let's make a little function which will return the number of days needed:



    def is_leap(year):
    return 366 if calendar.isleap(year) else 365


  • your program is a big main() where you've written all the logic. This makes the whole really hard to read! Instead, write small functions for each specific task.



  • your while loops can be re-written as for loops. This will save some variables allocation (extra incrementation):



    for year in range(1971, 2071):
    ...
    for current_day in range(1, days + 1):
    ...



I'd really like to review more of this, but unfortunately I can't test any of the code + you didn't add enough context (for me) to understand better what you'd like to achieve with each piece of code. For example, you create some variables in a for loop, and redefine them in another loop with the previous values from the previous loop... and so on.






share|improve this answer























  • Some more info. Simulating the development of 100 individual budworm throughout every year (each may have different growth rate, hence 100 random numbers). Calcs start Jan 1 and from then on if the temperature is between certain numbers (I have 4-hr timesteps through each day of each year), the development is incremented with the result of the equations. Once development is greater than 1, the budworm molts to the next larval stage, the julian date is recorded, and development is reset to 0. The next larval stage will have different parameters on how temperature affects the development.
    – bobby12345
    Feb 21 '17 at 15:34










  • Because temperature will vary between pixels, the development of a budworm will be spatially different, will molt to next larval stage at different times and will require the new parameters for calculating the development of next larval stage
    – bobby12345
    Feb 21 '17 at 15:37











Your Answer





StackExchange.ifUsing("editor", function () {
return StackExchange.using("mathjaxEditing", function () {
StackExchange.MarkdownEditor.creationCallbacks.add(function (editor, postfix) {
StackExchange.mathjaxEditing.prepareWmdForMathJax(editor, postfix, [["\$", "\$"]]);
});
});
}, "mathjax-editing");

StackExchange.ifUsing("editor", function () {
StackExchange.using("externalEditor", function () {
StackExchange.using("snippets", function () {
StackExchange.snippets.init();
});
});
}, "code-snippets");

StackExchange.ready(function() {
var channelOptions = {
tags: "".split(" "),
id: "196"
};
initTagRenderer("".split(" "), "".split(" "), channelOptions);

StackExchange.using("externalEditor", function() {
// Have to fire editor after snippets, if snippets enabled
if (StackExchange.settings.snippets.snippetsEnabled) {
StackExchange.using("snippets", function() {
createEditor();
});
}
else {
createEditor();
}
});

function createEditor() {
StackExchange.prepareEditor({
heartbeatType: 'answer',
convertImagesToLinks: false,
noModals: true,
showLowRepImageUploadWarning: true,
reputationToPostImages: null,
bindNavPrevention: true,
postfix: "",
imageUploader: {
brandingHtml: "Powered by u003ca class="icon-imgur-white" href="https://imgur.com/"u003eu003c/au003e",
contentPolicyHtml: "User contributions licensed under u003ca href="https://creativecommons.org/licenses/by-sa/3.0/"u003ecc by-sa 3.0 with attribution requiredu003c/au003e u003ca href="https://stackoverflow.com/legal/content-policy"u003e(content policy)u003c/au003e",
allowUrls: true
},
onDemand: true,
discardSelector: ".discard-answer"
,immediatelyShowMarkdownHelp:true
});


}
});














 

draft saved


draft discarded


















StackExchange.ready(
function () {
StackExchange.openid.initPostLogin('.new-post-login', 'https%3a%2f%2fcodereview.stackexchange.com%2fquestions%2f155610%2fraster-processing-for-climate-change-model%23new-answer', 'question_page');
}
);

Post as a guest















Required, but never shown

























1 Answer
1






active

oldest

votes








1 Answer
1






active

oldest

votes









active

oldest

votes






active

oldest

votes








up vote
3
down vote













First, fix the bug!



You're using sys.exc_info()[-1].tb_lineno but you're not importing the sys module anywhere, so do that first !



Now, I'll start my review with a couple of style guides. You can read more about them here.




  • use snake_case convention for your variable names


  • imports should be grouped in the following order:




    • standard library imports

    • related third party imports

    • local application/library specific imports




You should put a blank line between each group of imports. More, absolute imports are recommended, as they are usually more readable and tend to be better behaved (or at least give better error messages) if the import system is incorrectly configured




  • don't add redundant parentheses in conditional statements

  • use the print() function even if you're using Python 2.7.12

  • add a space before and after each operator

  • try to keep your lines no longer than 120 characters


  • use string formatting when you're concatenating strings: e.g:



    print("Skiwi hates stacks {} times more than I do!".format(100))



With all of the above in mind, we'll have the following code (the formatting is a bit ugly because you have really long formulas but heh..):



import math
import sys

import arcpy
from arcpy.sa import Con, Raster

arcpy.CheckOutExtension('spatial')
arcpy.env.overwriteOutput = True

climate_models = ["CanESM2", "CSIRO-Mk3-6-0", "HadGEM2-ES"]
concentration_pathway = ["rcp26", "rcp45", "rcp85"]

budworm_list = [2.1, 2.3, 0.6, 2.1, 1.7, 1.1, 1.6, 2.4, 2.1, 2.2, 0.9, 1.6, 0.6, 2.0, 1.7, 0.6, 0.7, 1.8, 1.3, 2.4,
2.1, 0.8, 2.2, 1.6, 0.5, 1.0, 1.6, 2.4, 2.5, 2.5, 0.8, 0.5, 2.5, 0.9, 1.2, 0.7, 1.8, 1.8, 0.4, 1.9,
0.8, 1.3, 0.6, 2.2, 0.7, 1.6, 1.0, 0.7, 0.8, 2.4, 1.7, 0.5, 1.1, 2.2, 0.9, 1.6, 1.8, 0.8, 2.3, 0.9,
0.8, 1.0, 0.5, 1.6, 2.5, 2.2, 2.0, 1.3, 1.6, 1.2, 2.4, 0.9, 1.1, 0.8, 2.1, 2.2, 1.5, 2.4, 2.5, 2.2,
1.3, 0.4, 1.5, 1.3, 1.7, 1.9, 1.8, 0.4, 0.4, 0.5, 1.5, 0.5, 0.6, 0.7, 1.6, 0.9, 0.9, 2.5, 1.4, 1.3]
time_step = ['4h', '8h', '12h', '16h', '20h', '24h']

try:
for gcm in climate_models:
for rcp in concentration_pathway:
budworm_count = 1
for budworm in budworm_list:
year = 1971
while year <= 2070:
if (year % 400 == 0) or ((year % 4 == 0) and (year % 100 != 0)):
days = 366
else:
days = 365
dev_stage = Raster("C:\Users\Robert\Documents\PCIC_ClimateData\TEST\Dev_Stage.tif")
dev_stage_female = Raster(
"C:\Users\Robert\Documents\PCIC_ClimateData\TEST\Dev_Stage_female.tif")
dev_l2o = Raster("C:\Users\Robert\Documents\PCIC_ClimateData\TEST\Dev_L2o.tif")
dev_l2 = Raster("C:\Users\Robert\Documents\PCIC_ClimateData\TEST\Dev_L2.tif")
dev_l3 = Raster("C:\Users\Robert\Documents\PCIC_ClimateData\TEST\Dev_L3.tif")
dev_l4 = Raster("C:\Users\Robert\Documents\PCIC_ClimateData\TEST\Dev_L4.tif")
dev_l5 = Raster("C:\Users\Robert\Documents\PCIC_ClimateData\TEST\Dev_L5.tif")
dev_l6_male = Raster("C:\Users\Robert\Documents\PCIC_ClimateData\TEST\Dev_L6_male.tif")
dev_l6_female = Raster("C:\Users\Robert\Documents\PCIC_ClimateData\TEST\Dev_L6_female.tif")
dev_pupa_male = Raster("C:\Users\Robert\Documents\PCIC_ClimateData\TEST\Dev_pupa_male.tif")
dev_pupa_female = Raster(
"C:\Users\Robert\Documents\PCIC_ClimateData\TEST\Dev_pupa_female.tif")
l2_date = Raster("C:\Users\Robert\Documents\PCIC_ClimateData\TEST\Jul_L2.tif")
l3_date = Raster("C:\Users\Robert\Documents\PCIC_ClimateData\TEST\Jul_L3.tif")
l4_date = Raster("C:\Users\Robert\Documents\PCIC_ClimateData\TEST\Jul_L4.tif")
l5_date = Raster("C:\Users\Robert\Documents\PCIC_ClimateData\TEST\Jul_L5.tif")
l6_date = Raster("C:\Users\Robert\Documents\PCIC_ClimateData\TEST\Jul_L6.tif")
pupa_male_date = Raster("C:\Users\Robert\Documents\PCIC_ClimateData\TEST\Jul_pupa_male.tif")
pupa_female_date = Raster(
"C:\Users\Robert\Documents\PCIC_ClimateData\TEST\Jul_pupa_female.tif")
adult_male_date = Raster("C:\Users\Robert\Documents\PCIC_ClimateData\TEST\Jul_adult_male.tif")
adult_female_date = Raster(
"C:\Users\Robert\Documents\PCIC_ClimateData\TEST\Jul_adult_female.tif")

current_day = 1
while current_day <= days:
for time in time_step:
tempfile = Raster(
"C:\Users\Robert\Documents\PCIC_ClimateData\{}\{}\4hTempTimeStep\{}_{}_{}.tif".format(
gcm, rcp, time, current_day, year))
dev_l2o = Con(((dev_stage == 1) & ((tempfile > 2.5) & (tempfile < 35)) & (dev_l2o <= 1)),
dev_l2o + (float(budworm)) * ((0.194 * (
(1 / (1 + (math.e ** (3.00 - (5.84 * ((tempfile - 2.5) / (35 - 2.5))))))) - (
math.e ** ((((tempfile - 2.5) / (35 - 2.5)) - 1) / 0.034)))) / 6), dev_l2o)
dev_l2 = Con(((dev_stage == 2) & ((tempfile > 4.4) & (tempfile < 38)) & (dev_l2 <= 1)),
dev_l2 + (float(budworm)) * ((0.919 * (
(1 / (1 + (math.e ** (2.91 - (5.32 * ((tempfile - 4.4) / (38 - 4.4))))))) - (
math.e ** ((((tempfile - 4.4) / (38 - 4.4)) - 1) / 0.061)))) / 6), dev_l2)
dev_l3 = Con(((dev_stage == 3) & ((tempfile > 4.4) & (tempfile < 38)) & (dev_l3 <= 1)),
dev_l3 + (float(budworm)) * ((0.438 * (
(1 / (1 + (math.e ** (3.06 - (6.85 * ((tempfile - 4.4) / (38 - 4.4))))))) - (
math.e ** ((((tempfile - 4.4) / (38 - 4.4)) - 1) / 0.061)))) / 6), dev_l3)
dev_l4 = Con(((dev_stage == 4) & ((tempfile > 4.4) & (tempfile < 38)) & (dev_l4 <= 1)),
dev_l4 + (float(budworm)) * ((1.211 * (
(1 / (1 + (math.e ** (3.80 - (7.55 * ((tempfile - 4.4) / (38 - 4.4))))))) - (
math.e ** ((((tempfile - 4.4) / (38 - 4.4)) - 1) / 0.148)))) / 6), dev_l4)
dev_l5 = Con(((dev_stage == 5) & ((tempfile > 4.4) & (tempfile < 38)) & (dev_l5 <= 1)),
dev_l5 + (float(budworm)) * ((0.269 * (
(1 / (1 + (math.e ** (3.02 - (8.57 * ((tempfile - 4.4) / (38 - 4.4))))))) - (
math.e ** ((((tempfile - 4.4) / (38 - 4.4)) - 1) / 0.005)))) / 6), dev_l5)
dev_l6_male = Con(
((dev_stage == 6) & ((tempfile > 4.4) & (tempfile < 38)) & (dev_l6_male <= 1)),
dev_l6_male + (float(budworm)) * ((0.288 * (
(1 / (1 + (math.e ** (2.67 - (5.03 * ((tempfile - 4.4) / (38 - 4.4))))))) - (
math.e ** ((((tempfile - 4.4) / (38 - 4.4)) - 1) / 0.151)))) / 6), dev_l6_male)
dev_l6_female = Con(
((dev_stage_female == 1) & ((tempfile > 4.4) & (tempfile < 38)) & (dev_l6_female <= 1)),
dev_l6_female + (float(budworm)) * ((0.317 * (
(1 / (1 + (math.e ** (3.06 - (4.66 * ((tempfile - 4.4) / (38 - 4.4))))))) - (
math.e ** ((((tempfile - 4.4) / (38 - 4.4)) - 1) / 0.136)))) / 6), dev_l6_female)
dev_pupa_male = Con(
((dev_stage == 7) & ((tempfile > 4.4) & (tempfile < 35)) & (dev_pupa_male <= 1)),
dev_pupa_male + ((0.259 * (
(1 / (1 + (math.e ** (2.75 - (4.66 * ((tempfile - 4.4) / (35 - 4.4))))))) - (
math.e ** ((((tempfile - 4.4) / (35 - 4.4)) - 1) / 0.053)))) / 6), dev_pupa_male)
dev_pupa_female = Con(((dev_stage_female == 2) & ((tempfile > 4.4) & (tempfile < 35)) & (
dev_pupa_female <= 1)), dev_pupa_female + ((0.205 * (
(1 / (1 + (math.e ** (2.85 - (6.28 * ((tempfile - 4.4) / (35 - 4.4))))))) - (
math.e ** ((((tempfile - 4.4) / (35 - 4.4)) - 1) / 0.044)))) / 6), dev_pupa_female)

dev_stage = Con(((dev_stage == 1) & (dev_l2o > 1)), dev_stage + 1, dev_stage)
dev_stage = Con(((dev_stage == 2) & (dev_l2 > 1)), dev_stage + 1, dev_stage)
dev_stage = Con(((dev_stage == 3) & (dev_l3 > 1)), dev_stage + 1, dev_stage)
dev_stage = Con(((dev_stage == 4) & (dev_l4 > 1)), dev_stage + 1, dev_stage)
dev_stage = Con(((dev_stage == 5) & (dev_l5 > 1)), dev_stage + 1, dev_stage)
dev_stage = Con(((dev_stage == 6) & (dev_l6_male > 1)), dev_stage + 1, dev_stage)
dev_stage = Con(((dev_stage == 7) & (dev_pupa_male > 1)), dev_stage + 1, dev_stage)
dev_stage_female = Con(((dev_stage_female == 1) & (dev_l6_female > 1)),
dev_stage_female + 1, dev_stage_female)
dev_stage_female = Con(((dev_stage_female == 2) & (dev_pupa_female > 1)),
dev_stage_female + 1, dev_stage_female)

l2_date = Con(((dev_l2o > 1) & (l2_date == 0)), l2_date + current_day, l2_date)
l3_date = Con(((dev_l2 > 1) & (l3_date == 0)), l3_date + current_day, l3_date)
l4_date = Con(((dev_l3 > 1) & (l4_date == 0)), l4_date + current_day, l4_date)
l5_date = Con(((dev_l4 > 1) & (l5_date == 0)), l5_date + current_day, l5_date)
l6_date = Con(((dev_l5 > 1) & (l6_date == 0)), l6_date + current_day, l6_date)
pupa_male_date = Con(((dev_l6_male > 1) & (pupa_male_date == 0)),
pupa_male_date + current_day, pupa_male_date)
adult_male_date = Con(((dev_pupa_male > 1) & (adult_male_date == 0)),
adult_male_date + current_day, adult_male_date)
pupa_female_date = Con(((dev_l6_female > 1) & (pupa_female_date == 0)),
pupa_female_date + current_day, pupa_female_date)
adult_female_date = Con(((dev_pupa_female > 1) & (adult_female_date == 0)),
adult_female_date + current_day, adult_female_date)

current_day = current_day + 1
l2_date.save(
"C:\Users\Robert\Documents\PCIC_ClimateData\TEST\{}_{}_Budworm_{}_L2_{}.tif".format(gcm,
rcp,
budworm_count,
year))
l3_date.save("C:\Users\Robert\Documents\PCIC_ClimateData\TEST\Budworm_{}_L3_{}.tif".format(
budworm_count, year))
l4_date.save("C:\Users\Robert\Documents\PCIC_ClimateData\TEST\Budworm_{}_L4_{}.tif".format(
budworm_count, year))
l5_date.save("C:\Users\Robert\Documents\PCIC_ClimateData\TEST\Budworm_{}_L5_{}.tif".format(
budworm_count, year))
l6_date.save("C:\Users\Robert\Documents\PCIC_ClimateData\TEST\Budworm_{}_L6_{}.tif".format(
budworm_count, year))
pupa_male_date.save(
"C:\Users\Robert\Documents\PCIC_ClimateData\TEST\Budworm_{}_pupa_male_{}.tif".format(
budworm_count, year))
adult_male_date.save(
"C:\Users\Robert\Documents\PCIC_ClimateData\TEST\Budworm_{}_adult_male_{}.tif".format(
budworm_count, year))
pupa_female_date.save(
"C:\Users\Robert\Documents\PCIC_ClimateData\TEST\Budworm_{}_pupa_female_{}.tif".format(
budworm_count, year))
adult_female_date.save(
"C:\Users\Robert\Documents\PCIC_ClimateData\TEST\Budworm_{}_adult_female_{}.tif".format(
budworm_count, year))

year += 1
budworm_count += 1
print("{} {} has finished".format(gcm, rcp))
except Exception as e:
print("An error occurred")
print(e)

except arcpy.ExecuteError:
print("A geoprocessing error occurred")
print(arcpy.GetMessages(2))

except IndexError:
print('Error on line {}'.format(sys.exc_info()[-1].tb_lineno))

arcpy.CheckInExtension('spatial')


Code related improvements





  • you said you want to generate a list of 100 random numbers (which should be the same while looping). That's quite easy to do if you're using the random module:



    def budworm_list(lower_limit, upper_limit, decimals):
    return [round(random.uniform(lower_limit, upper_limit), decimals) for _ in range(100)]


    Which can be latter called like this:



    budworm_list = generate_budworm_list(0.4, 2.5, 1)



  • use enumerate() to get rid of the usual counter. This Function can also take as a parameter the step which will let you start from 1.



    ...
    for budworm_count, budworm in enumerate(budworm_list, start=1):
    ...


  • don't initialize the year each time you're looping. Instead, declare it at the top of the program. More, that's a constant, which is usually declared in UPPER_CASE. (YEAR = 1971). Better yet, remove completely the while loop and add a for loop (see below).



  • there's a built-in which can tell you if an year is leap or not. You can find it in the calendar module:



    import calendar

    print calendar.isleap(1900)
    >> True


    So, let's make a little function which will return the number of days needed:



    def is_leap(year):
    return 366 if calendar.isleap(year) else 365


  • your program is a big main() where you've written all the logic. This makes the whole really hard to read! Instead, write small functions for each specific task.



  • your while loops can be re-written as for loops. This will save some variables allocation (extra incrementation):



    for year in range(1971, 2071):
    ...
    for current_day in range(1, days + 1):
    ...



I'd really like to review more of this, but unfortunately I can't test any of the code + you didn't add enough context (for me) to understand better what you'd like to achieve with each piece of code. For example, you create some variables in a for loop, and redefine them in another loop with the previous values from the previous loop... and so on.






share|improve this answer























  • Some more info. Simulating the development of 100 individual budworm throughout every year (each may have different growth rate, hence 100 random numbers). Calcs start Jan 1 and from then on if the temperature is between certain numbers (I have 4-hr timesteps through each day of each year), the development is incremented with the result of the equations. Once development is greater than 1, the budworm molts to the next larval stage, the julian date is recorded, and development is reset to 0. The next larval stage will have different parameters on how temperature affects the development.
    – bobby12345
    Feb 21 '17 at 15:34










  • Because temperature will vary between pixels, the development of a budworm will be spatially different, will molt to next larval stage at different times and will require the new parameters for calculating the development of next larval stage
    – bobby12345
    Feb 21 '17 at 15:37















up vote
3
down vote













First, fix the bug!



You're using sys.exc_info()[-1].tb_lineno but you're not importing the sys module anywhere, so do that first !



Now, I'll start my review with a couple of style guides. You can read more about them here.




  • use snake_case convention for your variable names


  • imports should be grouped in the following order:




    • standard library imports

    • related third party imports

    • local application/library specific imports




You should put a blank line between each group of imports. More, absolute imports are recommended, as they are usually more readable and tend to be better behaved (or at least give better error messages) if the import system is incorrectly configured




  • don't add redundant parentheses in conditional statements

  • use the print() function even if you're using Python 2.7.12

  • add a space before and after each operator

  • try to keep your lines no longer than 120 characters


  • use string formatting when you're concatenating strings: e.g:



    print("Skiwi hates stacks {} times more than I do!".format(100))



With all of the above in mind, we'll have the following code (the formatting is a bit ugly because you have really long formulas but heh..):



import math
import sys

import arcpy
from arcpy.sa import Con, Raster

arcpy.CheckOutExtension('spatial')
arcpy.env.overwriteOutput = True

climate_models = ["CanESM2", "CSIRO-Mk3-6-0", "HadGEM2-ES"]
concentration_pathway = ["rcp26", "rcp45", "rcp85"]

budworm_list = [2.1, 2.3, 0.6, 2.1, 1.7, 1.1, 1.6, 2.4, 2.1, 2.2, 0.9, 1.6, 0.6, 2.0, 1.7, 0.6, 0.7, 1.8, 1.3, 2.4,
2.1, 0.8, 2.2, 1.6, 0.5, 1.0, 1.6, 2.4, 2.5, 2.5, 0.8, 0.5, 2.5, 0.9, 1.2, 0.7, 1.8, 1.8, 0.4, 1.9,
0.8, 1.3, 0.6, 2.2, 0.7, 1.6, 1.0, 0.7, 0.8, 2.4, 1.7, 0.5, 1.1, 2.2, 0.9, 1.6, 1.8, 0.8, 2.3, 0.9,
0.8, 1.0, 0.5, 1.6, 2.5, 2.2, 2.0, 1.3, 1.6, 1.2, 2.4, 0.9, 1.1, 0.8, 2.1, 2.2, 1.5, 2.4, 2.5, 2.2,
1.3, 0.4, 1.5, 1.3, 1.7, 1.9, 1.8, 0.4, 0.4, 0.5, 1.5, 0.5, 0.6, 0.7, 1.6, 0.9, 0.9, 2.5, 1.4, 1.3]
time_step = ['4h', '8h', '12h', '16h', '20h', '24h']

try:
for gcm in climate_models:
for rcp in concentration_pathway:
budworm_count = 1
for budworm in budworm_list:
year = 1971
while year <= 2070:
if (year % 400 == 0) or ((year % 4 == 0) and (year % 100 != 0)):
days = 366
else:
days = 365
dev_stage = Raster("C:\Users\Robert\Documents\PCIC_ClimateData\TEST\Dev_Stage.tif")
dev_stage_female = Raster(
"C:\Users\Robert\Documents\PCIC_ClimateData\TEST\Dev_Stage_female.tif")
dev_l2o = Raster("C:\Users\Robert\Documents\PCIC_ClimateData\TEST\Dev_L2o.tif")
dev_l2 = Raster("C:\Users\Robert\Documents\PCIC_ClimateData\TEST\Dev_L2.tif")
dev_l3 = Raster("C:\Users\Robert\Documents\PCIC_ClimateData\TEST\Dev_L3.tif")
dev_l4 = Raster("C:\Users\Robert\Documents\PCIC_ClimateData\TEST\Dev_L4.tif")
dev_l5 = Raster("C:\Users\Robert\Documents\PCIC_ClimateData\TEST\Dev_L5.tif")
dev_l6_male = Raster("C:\Users\Robert\Documents\PCIC_ClimateData\TEST\Dev_L6_male.tif")
dev_l6_female = Raster("C:\Users\Robert\Documents\PCIC_ClimateData\TEST\Dev_L6_female.tif")
dev_pupa_male = Raster("C:\Users\Robert\Documents\PCIC_ClimateData\TEST\Dev_pupa_male.tif")
dev_pupa_female = Raster(
"C:\Users\Robert\Documents\PCIC_ClimateData\TEST\Dev_pupa_female.tif")
l2_date = Raster("C:\Users\Robert\Documents\PCIC_ClimateData\TEST\Jul_L2.tif")
l3_date = Raster("C:\Users\Robert\Documents\PCIC_ClimateData\TEST\Jul_L3.tif")
l4_date = Raster("C:\Users\Robert\Documents\PCIC_ClimateData\TEST\Jul_L4.tif")
l5_date = Raster("C:\Users\Robert\Documents\PCIC_ClimateData\TEST\Jul_L5.tif")
l6_date = Raster("C:\Users\Robert\Documents\PCIC_ClimateData\TEST\Jul_L6.tif")
pupa_male_date = Raster("C:\Users\Robert\Documents\PCIC_ClimateData\TEST\Jul_pupa_male.tif")
pupa_female_date = Raster(
"C:\Users\Robert\Documents\PCIC_ClimateData\TEST\Jul_pupa_female.tif")
adult_male_date = Raster("C:\Users\Robert\Documents\PCIC_ClimateData\TEST\Jul_adult_male.tif")
adult_female_date = Raster(
"C:\Users\Robert\Documents\PCIC_ClimateData\TEST\Jul_adult_female.tif")

current_day = 1
while current_day <= days:
for time in time_step:
tempfile = Raster(
"C:\Users\Robert\Documents\PCIC_ClimateData\{}\{}\4hTempTimeStep\{}_{}_{}.tif".format(
gcm, rcp, time, current_day, year))
dev_l2o = Con(((dev_stage == 1) & ((tempfile > 2.5) & (tempfile < 35)) & (dev_l2o <= 1)),
dev_l2o + (float(budworm)) * ((0.194 * (
(1 / (1 + (math.e ** (3.00 - (5.84 * ((tempfile - 2.5) / (35 - 2.5))))))) - (
math.e ** ((((tempfile - 2.5) / (35 - 2.5)) - 1) / 0.034)))) / 6), dev_l2o)
dev_l2 = Con(((dev_stage == 2) & ((tempfile > 4.4) & (tempfile < 38)) & (dev_l2 <= 1)),
dev_l2 + (float(budworm)) * ((0.919 * (
(1 / (1 + (math.e ** (2.91 - (5.32 * ((tempfile - 4.4) / (38 - 4.4))))))) - (
math.e ** ((((tempfile - 4.4) / (38 - 4.4)) - 1) / 0.061)))) / 6), dev_l2)
dev_l3 = Con(((dev_stage == 3) & ((tempfile > 4.4) & (tempfile < 38)) & (dev_l3 <= 1)),
dev_l3 + (float(budworm)) * ((0.438 * (
(1 / (1 + (math.e ** (3.06 - (6.85 * ((tempfile - 4.4) / (38 - 4.4))))))) - (
math.e ** ((((tempfile - 4.4) / (38 - 4.4)) - 1) / 0.061)))) / 6), dev_l3)
dev_l4 = Con(((dev_stage == 4) & ((tempfile > 4.4) & (tempfile < 38)) & (dev_l4 <= 1)),
dev_l4 + (float(budworm)) * ((1.211 * (
(1 / (1 + (math.e ** (3.80 - (7.55 * ((tempfile - 4.4) / (38 - 4.4))))))) - (
math.e ** ((((tempfile - 4.4) / (38 - 4.4)) - 1) / 0.148)))) / 6), dev_l4)
dev_l5 = Con(((dev_stage == 5) & ((tempfile > 4.4) & (tempfile < 38)) & (dev_l5 <= 1)),
dev_l5 + (float(budworm)) * ((0.269 * (
(1 / (1 + (math.e ** (3.02 - (8.57 * ((tempfile - 4.4) / (38 - 4.4))))))) - (
math.e ** ((((tempfile - 4.4) / (38 - 4.4)) - 1) / 0.005)))) / 6), dev_l5)
dev_l6_male = Con(
((dev_stage == 6) & ((tempfile > 4.4) & (tempfile < 38)) & (dev_l6_male <= 1)),
dev_l6_male + (float(budworm)) * ((0.288 * (
(1 / (1 + (math.e ** (2.67 - (5.03 * ((tempfile - 4.4) / (38 - 4.4))))))) - (
math.e ** ((((tempfile - 4.4) / (38 - 4.4)) - 1) / 0.151)))) / 6), dev_l6_male)
dev_l6_female = Con(
((dev_stage_female == 1) & ((tempfile > 4.4) & (tempfile < 38)) & (dev_l6_female <= 1)),
dev_l6_female + (float(budworm)) * ((0.317 * (
(1 / (1 + (math.e ** (3.06 - (4.66 * ((tempfile - 4.4) / (38 - 4.4))))))) - (
math.e ** ((((tempfile - 4.4) / (38 - 4.4)) - 1) / 0.136)))) / 6), dev_l6_female)
dev_pupa_male = Con(
((dev_stage == 7) & ((tempfile > 4.4) & (tempfile < 35)) & (dev_pupa_male <= 1)),
dev_pupa_male + ((0.259 * (
(1 / (1 + (math.e ** (2.75 - (4.66 * ((tempfile - 4.4) / (35 - 4.4))))))) - (
math.e ** ((((tempfile - 4.4) / (35 - 4.4)) - 1) / 0.053)))) / 6), dev_pupa_male)
dev_pupa_female = Con(((dev_stage_female == 2) & ((tempfile > 4.4) & (tempfile < 35)) & (
dev_pupa_female <= 1)), dev_pupa_female + ((0.205 * (
(1 / (1 + (math.e ** (2.85 - (6.28 * ((tempfile - 4.4) / (35 - 4.4))))))) - (
math.e ** ((((tempfile - 4.4) / (35 - 4.4)) - 1) / 0.044)))) / 6), dev_pupa_female)

dev_stage = Con(((dev_stage == 1) & (dev_l2o > 1)), dev_stage + 1, dev_stage)
dev_stage = Con(((dev_stage == 2) & (dev_l2 > 1)), dev_stage + 1, dev_stage)
dev_stage = Con(((dev_stage == 3) & (dev_l3 > 1)), dev_stage + 1, dev_stage)
dev_stage = Con(((dev_stage == 4) & (dev_l4 > 1)), dev_stage + 1, dev_stage)
dev_stage = Con(((dev_stage == 5) & (dev_l5 > 1)), dev_stage + 1, dev_stage)
dev_stage = Con(((dev_stage == 6) & (dev_l6_male > 1)), dev_stage + 1, dev_stage)
dev_stage = Con(((dev_stage == 7) & (dev_pupa_male > 1)), dev_stage + 1, dev_stage)
dev_stage_female = Con(((dev_stage_female == 1) & (dev_l6_female > 1)),
dev_stage_female + 1, dev_stage_female)
dev_stage_female = Con(((dev_stage_female == 2) & (dev_pupa_female > 1)),
dev_stage_female + 1, dev_stage_female)

l2_date = Con(((dev_l2o > 1) & (l2_date == 0)), l2_date + current_day, l2_date)
l3_date = Con(((dev_l2 > 1) & (l3_date == 0)), l3_date + current_day, l3_date)
l4_date = Con(((dev_l3 > 1) & (l4_date == 0)), l4_date + current_day, l4_date)
l5_date = Con(((dev_l4 > 1) & (l5_date == 0)), l5_date + current_day, l5_date)
l6_date = Con(((dev_l5 > 1) & (l6_date == 0)), l6_date + current_day, l6_date)
pupa_male_date = Con(((dev_l6_male > 1) & (pupa_male_date == 0)),
pupa_male_date + current_day, pupa_male_date)
adult_male_date = Con(((dev_pupa_male > 1) & (adult_male_date == 0)),
adult_male_date + current_day, adult_male_date)
pupa_female_date = Con(((dev_l6_female > 1) & (pupa_female_date == 0)),
pupa_female_date + current_day, pupa_female_date)
adult_female_date = Con(((dev_pupa_female > 1) & (adult_female_date == 0)),
adult_female_date + current_day, adult_female_date)

current_day = current_day + 1
l2_date.save(
"C:\Users\Robert\Documents\PCIC_ClimateData\TEST\{}_{}_Budworm_{}_L2_{}.tif".format(gcm,
rcp,
budworm_count,
year))
l3_date.save("C:\Users\Robert\Documents\PCIC_ClimateData\TEST\Budworm_{}_L3_{}.tif".format(
budworm_count, year))
l4_date.save("C:\Users\Robert\Documents\PCIC_ClimateData\TEST\Budworm_{}_L4_{}.tif".format(
budworm_count, year))
l5_date.save("C:\Users\Robert\Documents\PCIC_ClimateData\TEST\Budworm_{}_L5_{}.tif".format(
budworm_count, year))
l6_date.save("C:\Users\Robert\Documents\PCIC_ClimateData\TEST\Budworm_{}_L6_{}.tif".format(
budworm_count, year))
pupa_male_date.save(
"C:\Users\Robert\Documents\PCIC_ClimateData\TEST\Budworm_{}_pupa_male_{}.tif".format(
budworm_count, year))
adult_male_date.save(
"C:\Users\Robert\Documents\PCIC_ClimateData\TEST\Budworm_{}_adult_male_{}.tif".format(
budworm_count, year))
pupa_female_date.save(
"C:\Users\Robert\Documents\PCIC_ClimateData\TEST\Budworm_{}_pupa_female_{}.tif".format(
budworm_count, year))
adult_female_date.save(
"C:\Users\Robert\Documents\PCIC_ClimateData\TEST\Budworm_{}_adult_female_{}.tif".format(
budworm_count, year))

year += 1
budworm_count += 1
print("{} {} has finished".format(gcm, rcp))
except Exception as e:
print("An error occurred")
print(e)

except arcpy.ExecuteError:
print("A geoprocessing error occurred")
print(arcpy.GetMessages(2))

except IndexError:
print('Error on line {}'.format(sys.exc_info()[-1].tb_lineno))

arcpy.CheckInExtension('spatial')


Code related improvements





  • you said you want to generate a list of 100 random numbers (which should be the same while looping). That's quite easy to do if you're using the random module:



    def budworm_list(lower_limit, upper_limit, decimals):
    return [round(random.uniform(lower_limit, upper_limit), decimals) for _ in range(100)]


    Which can be latter called like this:



    budworm_list = generate_budworm_list(0.4, 2.5, 1)



  • use enumerate() to get rid of the usual counter. This Function can also take as a parameter the step which will let you start from 1.



    ...
    for budworm_count, budworm in enumerate(budworm_list, start=1):
    ...


  • don't initialize the year each time you're looping. Instead, declare it at the top of the program. More, that's a constant, which is usually declared in UPPER_CASE. (YEAR = 1971). Better yet, remove completely the while loop and add a for loop (see below).



  • there's a built-in which can tell you if an year is leap or not. You can find it in the calendar module:



    import calendar

    print calendar.isleap(1900)
    >> True


    So, let's make a little function which will return the number of days needed:



    def is_leap(year):
    return 366 if calendar.isleap(year) else 365


  • your program is a big main() where you've written all the logic. This makes the whole really hard to read! Instead, write small functions for each specific task.



  • your while loops can be re-written as for loops. This will save some variables allocation (extra incrementation):



    for year in range(1971, 2071):
    ...
    for current_day in range(1, days + 1):
    ...



I'd really like to review more of this, but unfortunately I can't test any of the code + you didn't add enough context (for me) to understand better what you'd like to achieve with each piece of code. For example, you create some variables in a for loop, and redefine them in another loop with the previous values from the previous loop... and so on.






share|improve this answer























  • Some more info. Simulating the development of 100 individual budworm throughout every year (each may have different growth rate, hence 100 random numbers). Calcs start Jan 1 and from then on if the temperature is between certain numbers (I have 4-hr timesteps through each day of each year), the development is incremented with the result of the equations. Once development is greater than 1, the budworm molts to the next larval stage, the julian date is recorded, and development is reset to 0. The next larval stage will have different parameters on how temperature affects the development.
    – bobby12345
    Feb 21 '17 at 15:34










  • Because temperature will vary between pixels, the development of a budworm will be spatially different, will molt to next larval stage at different times and will require the new parameters for calculating the development of next larval stage
    – bobby12345
    Feb 21 '17 at 15:37













up vote
3
down vote










up vote
3
down vote









First, fix the bug!



You're using sys.exc_info()[-1].tb_lineno but you're not importing the sys module anywhere, so do that first !



Now, I'll start my review with a couple of style guides. You can read more about them here.




  • use snake_case convention for your variable names


  • imports should be grouped in the following order:




    • standard library imports

    • related third party imports

    • local application/library specific imports




You should put a blank line between each group of imports. More, absolute imports are recommended, as they are usually more readable and tend to be better behaved (or at least give better error messages) if the import system is incorrectly configured




  • don't add redundant parentheses in conditional statements

  • use the print() function even if you're using Python 2.7.12

  • add a space before and after each operator

  • try to keep your lines no longer than 120 characters


  • use string formatting when you're concatenating strings: e.g:



    print("Skiwi hates stacks {} times more than I do!".format(100))



With all of the above in mind, we'll have the following code (the formatting is a bit ugly because you have really long formulas but heh..):



import math
import sys

import arcpy
from arcpy.sa import Con, Raster

arcpy.CheckOutExtension('spatial')
arcpy.env.overwriteOutput = True

climate_models = ["CanESM2", "CSIRO-Mk3-6-0", "HadGEM2-ES"]
concentration_pathway = ["rcp26", "rcp45", "rcp85"]

budworm_list = [2.1, 2.3, 0.6, 2.1, 1.7, 1.1, 1.6, 2.4, 2.1, 2.2, 0.9, 1.6, 0.6, 2.0, 1.7, 0.6, 0.7, 1.8, 1.3, 2.4,
2.1, 0.8, 2.2, 1.6, 0.5, 1.0, 1.6, 2.4, 2.5, 2.5, 0.8, 0.5, 2.5, 0.9, 1.2, 0.7, 1.8, 1.8, 0.4, 1.9,
0.8, 1.3, 0.6, 2.2, 0.7, 1.6, 1.0, 0.7, 0.8, 2.4, 1.7, 0.5, 1.1, 2.2, 0.9, 1.6, 1.8, 0.8, 2.3, 0.9,
0.8, 1.0, 0.5, 1.6, 2.5, 2.2, 2.0, 1.3, 1.6, 1.2, 2.4, 0.9, 1.1, 0.8, 2.1, 2.2, 1.5, 2.4, 2.5, 2.2,
1.3, 0.4, 1.5, 1.3, 1.7, 1.9, 1.8, 0.4, 0.4, 0.5, 1.5, 0.5, 0.6, 0.7, 1.6, 0.9, 0.9, 2.5, 1.4, 1.3]
time_step = ['4h', '8h', '12h', '16h', '20h', '24h']

try:
for gcm in climate_models:
for rcp in concentration_pathway:
budworm_count = 1
for budworm in budworm_list:
year = 1971
while year <= 2070:
if (year % 400 == 0) or ((year % 4 == 0) and (year % 100 != 0)):
days = 366
else:
days = 365
dev_stage = Raster("C:\Users\Robert\Documents\PCIC_ClimateData\TEST\Dev_Stage.tif")
dev_stage_female = Raster(
"C:\Users\Robert\Documents\PCIC_ClimateData\TEST\Dev_Stage_female.tif")
dev_l2o = Raster("C:\Users\Robert\Documents\PCIC_ClimateData\TEST\Dev_L2o.tif")
dev_l2 = Raster("C:\Users\Robert\Documents\PCIC_ClimateData\TEST\Dev_L2.tif")
dev_l3 = Raster("C:\Users\Robert\Documents\PCIC_ClimateData\TEST\Dev_L3.tif")
dev_l4 = Raster("C:\Users\Robert\Documents\PCIC_ClimateData\TEST\Dev_L4.tif")
dev_l5 = Raster("C:\Users\Robert\Documents\PCIC_ClimateData\TEST\Dev_L5.tif")
dev_l6_male = Raster("C:\Users\Robert\Documents\PCIC_ClimateData\TEST\Dev_L6_male.tif")
dev_l6_female = Raster("C:\Users\Robert\Documents\PCIC_ClimateData\TEST\Dev_L6_female.tif")
dev_pupa_male = Raster("C:\Users\Robert\Documents\PCIC_ClimateData\TEST\Dev_pupa_male.tif")
dev_pupa_female = Raster(
"C:\Users\Robert\Documents\PCIC_ClimateData\TEST\Dev_pupa_female.tif")
l2_date = Raster("C:\Users\Robert\Documents\PCIC_ClimateData\TEST\Jul_L2.tif")
l3_date = Raster("C:\Users\Robert\Documents\PCIC_ClimateData\TEST\Jul_L3.tif")
l4_date = Raster("C:\Users\Robert\Documents\PCIC_ClimateData\TEST\Jul_L4.tif")
l5_date = Raster("C:\Users\Robert\Documents\PCIC_ClimateData\TEST\Jul_L5.tif")
l6_date = Raster("C:\Users\Robert\Documents\PCIC_ClimateData\TEST\Jul_L6.tif")
pupa_male_date = Raster("C:\Users\Robert\Documents\PCIC_ClimateData\TEST\Jul_pupa_male.tif")
pupa_female_date = Raster(
"C:\Users\Robert\Documents\PCIC_ClimateData\TEST\Jul_pupa_female.tif")
adult_male_date = Raster("C:\Users\Robert\Documents\PCIC_ClimateData\TEST\Jul_adult_male.tif")
adult_female_date = Raster(
"C:\Users\Robert\Documents\PCIC_ClimateData\TEST\Jul_adult_female.tif")

current_day = 1
while current_day <= days:
for time in time_step:
tempfile = Raster(
"C:\Users\Robert\Documents\PCIC_ClimateData\{}\{}\4hTempTimeStep\{}_{}_{}.tif".format(
gcm, rcp, time, current_day, year))
dev_l2o = Con(((dev_stage == 1) & ((tempfile > 2.5) & (tempfile < 35)) & (dev_l2o <= 1)),
dev_l2o + (float(budworm)) * ((0.194 * (
(1 / (1 + (math.e ** (3.00 - (5.84 * ((tempfile - 2.5) / (35 - 2.5))))))) - (
math.e ** ((((tempfile - 2.5) / (35 - 2.5)) - 1) / 0.034)))) / 6), dev_l2o)
dev_l2 = Con(((dev_stage == 2) & ((tempfile > 4.4) & (tempfile < 38)) & (dev_l2 <= 1)),
dev_l2 + (float(budworm)) * ((0.919 * (
(1 / (1 + (math.e ** (2.91 - (5.32 * ((tempfile - 4.4) / (38 - 4.4))))))) - (
math.e ** ((((tempfile - 4.4) / (38 - 4.4)) - 1) / 0.061)))) / 6), dev_l2)
dev_l3 = Con(((dev_stage == 3) & ((tempfile > 4.4) & (tempfile < 38)) & (dev_l3 <= 1)),
dev_l3 + (float(budworm)) * ((0.438 * (
(1 / (1 + (math.e ** (3.06 - (6.85 * ((tempfile - 4.4) / (38 - 4.4))))))) - (
math.e ** ((((tempfile - 4.4) / (38 - 4.4)) - 1) / 0.061)))) / 6), dev_l3)
dev_l4 = Con(((dev_stage == 4) & ((tempfile > 4.4) & (tempfile < 38)) & (dev_l4 <= 1)),
dev_l4 + (float(budworm)) * ((1.211 * (
(1 / (1 + (math.e ** (3.80 - (7.55 * ((tempfile - 4.4) / (38 - 4.4))))))) - (
math.e ** ((((tempfile - 4.4) / (38 - 4.4)) - 1) / 0.148)))) / 6), dev_l4)
dev_l5 = Con(((dev_stage == 5) & ((tempfile > 4.4) & (tempfile < 38)) & (dev_l5 <= 1)),
dev_l5 + (float(budworm)) * ((0.269 * (
(1 / (1 + (math.e ** (3.02 - (8.57 * ((tempfile - 4.4) / (38 - 4.4))))))) - (
math.e ** ((((tempfile - 4.4) / (38 - 4.4)) - 1) / 0.005)))) / 6), dev_l5)
dev_l6_male = Con(
((dev_stage == 6) & ((tempfile > 4.4) & (tempfile < 38)) & (dev_l6_male <= 1)),
dev_l6_male + (float(budworm)) * ((0.288 * (
(1 / (1 + (math.e ** (2.67 - (5.03 * ((tempfile - 4.4) / (38 - 4.4))))))) - (
math.e ** ((((tempfile - 4.4) / (38 - 4.4)) - 1) / 0.151)))) / 6), dev_l6_male)
dev_l6_female = Con(
((dev_stage_female == 1) & ((tempfile > 4.4) & (tempfile < 38)) & (dev_l6_female <= 1)),
dev_l6_female + (float(budworm)) * ((0.317 * (
(1 / (1 + (math.e ** (3.06 - (4.66 * ((tempfile - 4.4) / (38 - 4.4))))))) - (
math.e ** ((((tempfile - 4.4) / (38 - 4.4)) - 1) / 0.136)))) / 6), dev_l6_female)
dev_pupa_male = Con(
((dev_stage == 7) & ((tempfile > 4.4) & (tempfile < 35)) & (dev_pupa_male <= 1)),
dev_pupa_male + ((0.259 * (
(1 / (1 + (math.e ** (2.75 - (4.66 * ((tempfile - 4.4) / (35 - 4.4))))))) - (
math.e ** ((((tempfile - 4.4) / (35 - 4.4)) - 1) / 0.053)))) / 6), dev_pupa_male)
dev_pupa_female = Con(((dev_stage_female == 2) & ((tempfile > 4.4) & (tempfile < 35)) & (
dev_pupa_female <= 1)), dev_pupa_female + ((0.205 * (
(1 / (1 + (math.e ** (2.85 - (6.28 * ((tempfile - 4.4) / (35 - 4.4))))))) - (
math.e ** ((((tempfile - 4.4) / (35 - 4.4)) - 1) / 0.044)))) / 6), dev_pupa_female)

dev_stage = Con(((dev_stage == 1) & (dev_l2o > 1)), dev_stage + 1, dev_stage)
dev_stage = Con(((dev_stage == 2) & (dev_l2 > 1)), dev_stage + 1, dev_stage)
dev_stage = Con(((dev_stage == 3) & (dev_l3 > 1)), dev_stage + 1, dev_stage)
dev_stage = Con(((dev_stage == 4) & (dev_l4 > 1)), dev_stage + 1, dev_stage)
dev_stage = Con(((dev_stage == 5) & (dev_l5 > 1)), dev_stage + 1, dev_stage)
dev_stage = Con(((dev_stage == 6) & (dev_l6_male > 1)), dev_stage + 1, dev_stage)
dev_stage = Con(((dev_stage == 7) & (dev_pupa_male > 1)), dev_stage + 1, dev_stage)
dev_stage_female = Con(((dev_stage_female == 1) & (dev_l6_female > 1)),
dev_stage_female + 1, dev_stage_female)
dev_stage_female = Con(((dev_stage_female == 2) & (dev_pupa_female > 1)),
dev_stage_female + 1, dev_stage_female)

l2_date = Con(((dev_l2o > 1) & (l2_date == 0)), l2_date + current_day, l2_date)
l3_date = Con(((dev_l2 > 1) & (l3_date == 0)), l3_date + current_day, l3_date)
l4_date = Con(((dev_l3 > 1) & (l4_date == 0)), l4_date + current_day, l4_date)
l5_date = Con(((dev_l4 > 1) & (l5_date == 0)), l5_date + current_day, l5_date)
l6_date = Con(((dev_l5 > 1) & (l6_date == 0)), l6_date + current_day, l6_date)
pupa_male_date = Con(((dev_l6_male > 1) & (pupa_male_date == 0)),
pupa_male_date + current_day, pupa_male_date)
adult_male_date = Con(((dev_pupa_male > 1) & (adult_male_date == 0)),
adult_male_date + current_day, adult_male_date)
pupa_female_date = Con(((dev_l6_female > 1) & (pupa_female_date == 0)),
pupa_female_date + current_day, pupa_female_date)
adult_female_date = Con(((dev_pupa_female > 1) & (adult_female_date == 0)),
adult_female_date + current_day, adult_female_date)

current_day = current_day + 1
l2_date.save(
"C:\Users\Robert\Documents\PCIC_ClimateData\TEST\{}_{}_Budworm_{}_L2_{}.tif".format(gcm,
rcp,
budworm_count,
year))
l3_date.save("C:\Users\Robert\Documents\PCIC_ClimateData\TEST\Budworm_{}_L3_{}.tif".format(
budworm_count, year))
l4_date.save("C:\Users\Robert\Documents\PCIC_ClimateData\TEST\Budworm_{}_L4_{}.tif".format(
budworm_count, year))
l5_date.save("C:\Users\Robert\Documents\PCIC_ClimateData\TEST\Budworm_{}_L5_{}.tif".format(
budworm_count, year))
l6_date.save("C:\Users\Robert\Documents\PCIC_ClimateData\TEST\Budworm_{}_L6_{}.tif".format(
budworm_count, year))
pupa_male_date.save(
"C:\Users\Robert\Documents\PCIC_ClimateData\TEST\Budworm_{}_pupa_male_{}.tif".format(
budworm_count, year))
adult_male_date.save(
"C:\Users\Robert\Documents\PCIC_ClimateData\TEST\Budworm_{}_adult_male_{}.tif".format(
budworm_count, year))
pupa_female_date.save(
"C:\Users\Robert\Documents\PCIC_ClimateData\TEST\Budworm_{}_pupa_female_{}.tif".format(
budworm_count, year))
adult_female_date.save(
"C:\Users\Robert\Documents\PCIC_ClimateData\TEST\Budworm_{}_adult_female_{}.tif".format(
budworm_count, year))

year += 1
budworm_count += 1
print("{} {} has finished".format(gcm, rcp))
except Exception as e:
print("An error occurred")
print(e)

except arcpy.ExecuteError:
print("A geoprocessing error occurred")
print(arcpy.GetMessages(2))

except IndexError:
print('Error on line {}'.format(sys.exc_info()[-1].tb_lineno))

arcpy.CheckInExtension('spatial')


Code related improvements





  • you said you want to generate a list of 100 random numbers (which should be the same while looping). That's quite easy to do if you're using the random module:



    def budworm_list(lower_limit, upper_limit, decimals):
    return [round(random.uniform(lower_limit, upper_limit), decimals) for _ in range(100)]


    Which can be latter called like this:



    budworm_list = generate_budworm_list(0.4, 2.5, 1)



  • use enumerate() to get rid of the usual counter. This Function can also take as a parameter the step which will let you start from 1.



    ...
    for budworm_count, budworm in enumerate(budworm_list, start=1):
    ...


  • don't initialize the year each time you're looping. Instead, declare it at the top of the program. More, that's a constant, which is usually declared in UPPER_CASE. (YEAR = 1971). Better yet, remove completely the while loop and add a for loop (see below).



  • there's a built-in which can tell you if an year is leap or not. You can find it in the calendar module:



    import calendar

    print calendar.isleap(1900)
    >> True


    So, let's make a little function which will return the number of days needed:



    def is_leap(year):
    return 366 if calendar.isleap(year) else 365


  • your program is a big main() where you've written all the logic. This makes the whole really hard to read! Instead, write small functions for each specific task.



  • your while loops can be re-written as for loops. This will save some variables allocation (extra incrementation):



    for year in range(1971, 2071):
    ...
    for current_day in range(1, days + 1):
    ...



I'd really like to review more of this, but unfortunately I can't test any of the code + you didn't add enough context (for me) to understand better what you'd like to achieve with each piece of code. For example, you create some variables in a for loop, and redefine them in another loop with the previous values from the previous loop... and so on.






share|improve this answer














First, fix the bug!



You're using sys.exc_info()[-1].tb_lineno but you're not importing the sys module anywhere, so do that first !



Now, I'll start my review with a couple of style guides. You can read more about them here.




  • use snake_case convention for your variable names


  • imports should be grouped in the following order:




    • standard library imports

    • related third party imports

    • local application/library specific imports




You should put a blank line between each group of imports. More, absolute imports are recommended, as they are usually more readable and tend to be better behaved (or at least give better error messages) if the import system is incorrectly configured




  • don't add redundant parentheses in conditional statements

  • use the print() function even if you're using Python 2.7.12

  • add a space before and after each operator

  • try to keep your lines no longer than 120 characters


  • use string formatting when you're concatenating strings: e.g:



    print("Skiwi hates stacks {} times more than I do!".format(100))



With all of the above in mind, we'll have the following code (the formatting is a bit ugly because you have really long formulas but heh..):



import math
import sys

import arcpy
from arcpy.sa import Con, Raster

arcpy.CheckOutExtension('spatial')
arcpy.env.overwriteOutput = True

climate_models = ["CanESM2", "CSIRO-Mk3-6-0", "HadGEM2-ES"]
concentration_pathway = ["rcp26", "rcp45", "rcp85"]

budworm_list = [2.1, 2.3, 0.6, 2.1, 1.7, 1.1, 1.6, 2.4, 2.1, 2.2, 0.9, 1.6, 0.6, 2.0, 1.7, 0.6, 0.7, 1.8, 1.3, 2.4,
2.1, 0.8, 2.2, 1.6, 0.5, 1.0, 1.6, 2.4, 2.5, 2.5, 0.8, 0.5, 2.5, 0.9, 1.2, 0.7, 1.8, 1.8, 0.4, 1.9,
0.8, 1.3, 0.6, 2.2, 0.7, 1.6, 1.0, 0.7, 0.8, 2.4, 1.7, 0.5, 1.1, 2.2, 0.9, 1.6, 1.8, 0.8, 2.3, 0.9,
0.8, 1.0, 0.5, 1.6, 2.5, 2.2, 2.0, 1.3, 1.6, 1.2, 2.4, 0.9, 1.1, 0.8, 2.1, 2.2, 1.5, 2.4, 2.5, 2.2,
1.3, 0.4, 1.5, 1.3, 1.7, 1.9, 1.8, 0.4, 0.4, 0.5, 1.5, 0.5, 0.6, 0.7, 1.6, 0.9, 0.9, 2.5, 1.4, 1.3]
time_step = ['4h', '8h', '12h', '16h', '20h', '24h']

try:
for gcm in climate_models:
for rcp in concentration_pathway:
budworm_count = 1
for budworm in budworm_list:
year = 1971
while year <= 2070:
if (year % 400 == 0) or ((year % 4 == 0) and (year % 100 != 0)):
days = 366
else:
days = 365
dev_stage = Raster("C:\Users\Robert\Documents\PCIC_ClimateData\TEST\Dev_Stage.tif")
dev_stage_female = Raster(
"C:\Users\Robert\Documents\PCIC_ClimateData\TEST\Dev_Stage_female.tif")
dev_l2o = Raster("C:\Users\Robert\Documents\PCIC_ClimateData\TEST\Dev_L2o.tif")
dev_l2 = Raster("C:\Users\Robert\Documents\PCIC_ClimateData\TEST\Dev_L2.tif")
dev_l3 = Raster("C:\Users\Robert\Documents\PCIC_ClimateData\TEST\Dev_L3.tif")
dev_l4 = Raster("C:\Users\Robert\Documents\PCIC_ClimateData\TEST\Dev_L4.tif")
dev_l5 = Raster("C:\Users\Robert\Documents\PCIC_ClimateData\TEST\Dev_L5.tif")
dev_l6_male = Raster("C:\Users\Robert\Documents\PCIC_ClimateData\TEST\Dev_L6_male.tif")
dev_l6_female = Raster("C:\Users\Robert\Documents\PCIC_ClimateData\TEST\Dev_L6_female.tif")
dev_pupa_male = Raster("C:\Users\Robert\Documents\PCIC_ClimateData\TEST\Dev_pupa_male.tif")
dev_pupa_female = Raster(
"C:\Users\Robert\Documents\PCIC_ClimateData\TEST\Dev_pupa_female.tif")
l2_date = Raster("C:\Users\Robert\Documents\PCIC_ClimateData\TEST\Jul_L2.tif")
l3_date = Raster("C:\Users\Robert\Documents\PCIC_ClimateData\TEST\Jul_L3.tif")
l4_date = Raster("C:\Users\Robert\Documents\PCIC_ClimateData\TEST\Jul_L4.tif")
l5_date = Raster("C:\Users\Robert\Documents\PCIC_ClimateData\TEST\Jul_L5.tif")
l6_date = Raster("C:\Users\Robert\Documents\PCIC_ClimateData\TEST\Jul_L6.tif")
pupa_male_date = Raster("C:\Users\Robert\Documents\PCIC_ClimateData\TEST\Jul_pupa_male.tif")
pupa_female_date = Raster(
"C:\Users\Robert\Documents\PCIC_ClimateData\TEST\Jul_pupa_female.tif")
adult_male_date = Raster("C:\Users\Robert\Documents\PCIC_ClimateData\TEST\Jul_adult_male.tif")
adult_female_date = Raster(
"C:\Users\Robert\Documents\PCIC_ClimateData\TEST\Jul_adult_female.tif")

current_day = 1
while current_day <= days:
for time in time_step:
tempfile = Raster(
"C:\Users\Robert\Documents\PCIC_ClimateData\{}\{}\4hTempTimeStep\{}_{}_{}.tif".format(
gcm, rcp, time, current_day, year))
dev_l2o = Con(((dev_stage == 1) & ((tempfile > 2.5) & (tempfile < 35)) & (dev_l2o <= 1)),
dev_l2o + (float(budworm)) * ((0.194 * (
(1 / (1 + (math.e ** (3.00 - (5.84 * ((tempfile - 2.5) / (35 - 2.5))))))) - (
math.e ** ((((tempfile - 2.5) / (35 - 2.5)) - 1) / 0.034)))) / 6), dev_l2o)
dev_l2 = Con(((dev_stage == 2) & ((tempfile > 4.4) & (tempfile < 38)) & (dev_l2 <= 1)),
dev_l2 + (float(budworm)) * ((0.919 * (
(1 / (1 + (math.e ** (2.91 - (5.32 * ((tempfile - 4.4) / (38 - 4.4))))))) - (
math.e ** ((((tempfile - 4.4) / (38 - 4.4)) - 1) / 0.061)))) / 6), dev_l2)
dev_l3 = Con(((dev_stage == 3) & ((tempfile > 4.4) & (tempfile < 38)) & (dev_l3 <= 1)),
dev_l3 + (float(budworm)) * ((0.438 * (
(1 / (1 + (math.e ** (3.06 - (6.85 * ((tempfile - 4.4) / (38 - 4.4))))))) - (
math.e ** ((((tempfile - 4.4) / (38 - 4.4)) - 1) / 0.061)))) / 6), dev_l3)
dev_l4 = Con(((dev_stage == 4) & ((tempfile > 4.4) & (tempfile < 38)) & (dev_l4 <= 1)),
dev_l4 + (float(budworm)) * ((1.211 * (
(1 / (1 + (math.e ** (3.80 - (7.55 * ((tempfile - 4.4) / (38 - 4.4))))))) - (
math.e ** ((((tempfile - 4.4) / (38 - 4.4)) - 1) / 0.148)))) / 6), dev_l4)
dev_l5 = Con(((dev_stage == 5) & ((tempfile > 4.4) & (tempfile < 38)) & (dev_l5 <= 1)),
dev_l5 + (float(budworm)) * ((0.269 * (
(1 / (1 + (math.e ** (3.02 - (8.57 * ((tempfile - 4.4) / (38 - 4.4))))))) - (
math.e ** ((((tempfile - 4.4) / (38 - 4.4)) - 1) / 0.005)))) / 6), dev_l5)
dev_l6_male = Con(
((dev_stage == 6) & ((tempfile > 4.4) & (tempfile < 38)) & (dev_l6_male <= 1)),
dev_l6_male + (float(budworm)) * ((0.288 * (
(1 / (1 + (math.e ** (2.67 - (5.03 * ((tempfile - 4.4) / (38 - 4.4))))))) - (
math.e ** ((((tempfile - 4.4) / (38 - 4.4)) - 1) / 0.151)))) / 6), dev_l6_male)
dev_l6_female = Con(
((dev_stage_female == 1) & ((tempfile > 4.4) & (tempfile < 38)) & (dev_l6_female <= 1)),
dev_l6_female + (float(budworm)) * ((0.317 * (
(1 / (1 + (math.e ** (3.06 - (4.66 * ((tempfile - 4.4) / (38 - 4.4))))))) - (
math.e ** ((((tempfile - 4.4) / (38 - 4.4)) - 1) / 0.136)))) / 6), dev_l6_female)
dev_pupa_male = Con(
((dev_stage == 7) & ((tempfile > 4.4) & (tempfile < 35)) & (dev_pupa_male <= 1)),
dev_pupa_male + ((0.259 * (
(1 / (1 + (math.e ** (2.75 - (4.66 * ((tempfile - 4.4) / (35 - 4.4))))))) - (
math.e ** ((((tempfile - 4.4) / (35 - 4.4)) - 1) / 0.053)))) / 6), dev_pupa_male)
dev_pupa_female = Con(((dev_stage_female == 2) & ((tempfile > 4.4) & (tempfile < 35)) & (
dev_pupa_female <= 1)), dev_pupa_female + ((0.205 * (
(1 / (1 + (math.e ** (2.85 - (6.28 * ((tempfile - 4.4) / (35 - 4.4))))))) - (
math.e ** ((((tempfile - 4.4) / (35 - 4.4)) - 1) / 0.044)))) / 6), dev_pupa_female)

dev_stage = Con(((dev_stage == 1) & (dev_l2o > 1)), dev_stage + 1, dev_stage)
dev_stage = Con(((dev_stage == 2) & (dev_l2 > 1)), dev_stage + 1, dev_stage)
dev_stage = Con(((dev_stage == 3) & (dev_l3 > 1)), dev_stage + 1, dev_stage)
dev_stage = Con(((dev_stage == 4) & (dev_l4 > 1)), dev_stage + 1, dev_stage)
dev_stage = Con(((dev_stage == 5) & (dev_l5 > 1)), dev_stage + 1, dev_stage)
dev_stage = Con(((dev_stage == 6) & (dev_l6_male > 1)), dev_stage + 1, dev_stage)
dev_stage = Con(((dev_stage == 7) & (dev_pupa_male > 1)), dev_stage + 1, dev_stage)
dev_stage_female = Con(((dev_stage_female == 1) & (dev_l6_female > 1)),
dev_stage_female + 1, dev_stage_female)
dev_stage_female = Con(((dev_stage_female == 2) & (dev_pupa_female > 1)),
dev_stage_female + 1, dev_stage_female)

l2_date = Con(((dev_l2o > 1) & (l2_date == 0)), l2_date + current_day, l2_date)
l3_date = Con(((dev_l2 > 1) & (l3_date == 0)), l3_date + current_day, l3_date)
l4_date = Con(((dev_l3 > 1) & (l4_date == 0)), l4_date + current_day, l4_date)
l5_date = Con(((dev_l4 > 1) & (l5_date == 0)), l5_date + current_day, l5_date)
l6_date = Con(((dev_l5 > 1) & (l6_date == 0)), l6_date + current_day, l6_date)
pupa_male_date = Con(((dev_l6_male > 1) & (pupa_male_date == 0)),
pupa_male_date + current_day, pupa_male_date)
adult_male_date = Con(((dev_pupa_male > 1) & (adult_male_date == 0)),
adult_male_date + current_day, adult_male_date)
pupa_female_date = Con(((dev_l6_female > 1) & (pupa_female_date == 0)),
pupa_female_date + current_day, pupa_female_date)
adult_female_date = Con(((dev_pupa_female > 1) & (adult_female_date == 0)),
adult_female_date + current_day, adult_female_date)

current_day = current_day + 1
l2_date.save(
"C:\Users\Robert\Documents\PCIC_ClimateData\TEST\{}_{}_Budworm_{}_L2_{}.tif".format(gcm,
rcp,
budworm_count,
year))
l3_date.save("C:\Users\Robert\Documents\PCIC_ClimateData\TEST\Budworm_{}_L3_{}.tif".format(
budworm_count, year))
l4_date.save("C:\Users\Robert\Documents\PCIC_ClimateData\TEST\Budworm_{}_L4_{}.tif".format(
budworm_count, year))
l5_date.save("C:\Users\Robert\Documents\PCIC_ClimateData\TEST\Budworm_{}_L5_{}.tif".format(
budworm_count, year))
l6_date.save("C:\Users\Robert\Documents\PCIC_ClimateData\TEST\Budworm_{}_L6_{}.tif".format(
budworm_count, year))
pupa_male_date.save(
"C:\Users\Robert\Documents\PCIC_ClimateData\TEST\Budworm_{}_pupa_male_{}.tif".format(
budworm_count, year))
adult_male_date.save(
"C:\Users\Robert\Documents\PCIC_ClimateData\TEST\Budworm_{}_adult_male_{}.tif".format(
budworm_count, year))
pupa_female_date.save(
"C:\Users\Robert\Documents\PCIC_ClimateData\TEST\Budworm_{}_pupa_female_{}.tif".format(
budworm_count, year))
adult_female_date.save(
"C:\Users\Robert\Documents\PCIC_ClimateData\TEST\Budworm_{}_adult_female_{}.tif".format(
budworm_count, year))

year += 1
budworm_count += 1
print("{} {} has finished".format(gcm, rcp))
except Exception as e:
print("An error occurred")
print(e)

except arcpy.ExecuteError:
print("A geoprocessing error occurred")
print(arcpy.GetMessages(2))

except IndexError:
print('Error on line {}'.format(sys.exc_info()[-1].tb_lineno))

arcpy.CheckInExtension('spatial')


Code related improvements





  • you said you want to generate a list of 100 random numbers (which should be the same while looping). That's quite easy to do if you're using the random module:



    def budworm_list(lower_limit, upper_limit, decimals):
    return [round(random.uniform(lower_limit, upper_limit), decimals) for _ in range(100)]


    Which can be latter called like this:



    budworm_list = generate_budworm_list(0.4, 2.5, 1)



  • use enumerate() to get rid of the usual counter. This Function can also take as a parameter the step which will let you start from 1.



    ...
    for budworm_count, budworm in enumerate(budworm_list, start=1):
    ...


  • don't initialize the year each time you're looping. Instead, declare it at the top of the program. More, that's a constant, which is usually declared in UPPER_CASE. (YEAR = 1971). Better yet, remove completely the while loop and add a for loop (see below).



  • there's a built-in which can tell you if an year is leap or not. You can find it in the calendar module:



    import calendar

    print calendar.isleap(1900)
    >> True


    So, let's make a little function which will return the number of days needed:



    def is_leap(year):
    return 366 if calendar.isleap(year) else 365


  • your program is a big main() where you've written all the logic. This makes the whole really hard to read! Instead, write small functions for each specific task.



  • your while loops can be re-written as for loops. This will save some variables allocation (extra incrementation):



    for year in range(1971, 2071):
    ...
    for current_day in range(1, days + 1):
    ...



I'd really like to review more of this, but unfortunately I can't test any of the code + you didn't add enough context (for me) to understand better what you'd like to achieve with each piece of code. For example, you create some variables in a for loop, and redefine them in another loop with the previous values from the previous loop... and so on.







share|improve this answer














share|improve this answer



share|improve this answer








edited Feb 21 '17 at 15:13

























answered Feb 21 '17 at 15:00









яүυк

7,07621854




7,07621854












  • Some more info. Simulating the development of 100 individual budworm throughout every year (each may have different growth rate, hence 100 random numbers). Calcs start Jan 1 and from then on if the temperature is between certain numbers (I have 4-hr timesteps through each day of each year), the development is incremented with the result of the equations. Once development is greater than 1, the budworm molts to the next larval stage, the julian date is recorded, and development is reset to 0. The next larval stage will have different parameters on how temperature affects the development.
    – bobby12345
    Feb 21 '17 at 15:34










  • Because temperature will vary between pixels, the development of a budworm will be spatially different, will molt to next larval stage at different times and will require the new parameters for calculating the development of next larval stage
    – bobby12345
    Feb 21 '17 at 15:37


















  • Some more info. Simulating the development of 100 individual budworm throughout every year (each may have different growth rate, hence 100 random numbers). Calcs start Jan 1 and from then on if the temperature is between certain numbers (I have 4-hr timesteps through each day of each year), the development is incremented with the result of the equations. Once development is greater than 1, the budworm molts to the next larval stage, the julian date is recorded, and development is reset to 0. The next larval stage will have different parameters on how temperature affects the development.
    – bobby12345
    Feb 21 '17 at 15:34










  • Because temperature will vary between pixels, the development of a budworm will be spatially different, will molt to next larval stage at different times and will require the new parameters for calculating the development of next larval stage
    – bobby12345
    Feb 21 '17 at 15:37
















Some more info. Simulating the development of 100 individual budworm throughout every year (each may have different growth rate, hence 100 random numbers). Calcs start Jan 1 and from then on if the temperature is between certain numbers (I have 4-hr timesteps through each day of each year), the development is incremented with the result of the equations. Once development is greater than 1, the budworm molts to the next larval stage, the julian date is recorded, and development is reset to 0. The next larval stage will have different parameters on how temperature affects the development.
– bobby12345
Feb 21 '17 at 15:34




Some more info. Simulating the development of 100 individual budworm throughout every year (each may have different growth rate, hence 100 random numbers). Calcs start Jan 1 and from then on if the temperature is between certain numbers (I have 4-hr timesteps through each day of each year), the development is incremented with the result of the equations. Once development is greater than 1, the budworm molts to the next larval stage, the julian date is recorded, and development is reset to 0. The next larval stage will have different parameters on how temperature affects the development.
– bobby12345
Feb 21 '17 at 15:34












Because temperature will vary between pixels, the development of a budworm will be spatially different, will molt to next larval stage at different times and will require the new parameters for calculating the development of next larval stage
– bobby12345
Feb 21 '17 at 15:37




Because temperature will vary between pixels, the development of a budworm will be spatially different, will molt to next larval stage at different times and will require the new parameters for calculating the development of next larval stage
– bobby12345
Feb 21 '17 at 15:37


















 

draft saved


draft discarded



















































 


draft saved


draft discarded














StackExchange.ready(
function () {
StackExchange.openid.initPostLogin('.new-post-login', 'https%3a%2f%2fcodereview.stackexchange.com%2fquestions%2f155610%2fraster-processing-for-climate-change-model%23new-answer', 'question_page');
}
);

Post as a guest















Required, but never shown





















































Required, but never shown














Required, but never shown












Required, but never shown







Required, but never shown

































Required, but never shown














Required, but never shown












Required, but never shown







Required, but never shown







Popular posts from this blog

Quarter-circle Tiles

build a pushdown automaton that recognizes the reverse language of a given pushdown automaton?

Mont Emei