Skip to content

Commit

Permalink
Fixes for python 3 compatibility (#8)
Browse files Browse the repository at this point in the history 
Add use_hydrstress, gswp3w535 options
Fix to load software environment for multi-site cases

Co-authored-by: Daniel Ricciuto <ac.ricciuto@chrlogin1.lcrc.anl.gov>
  • Loading branch information
@dmricciuto
dmricciuto and Daniel Ricciuto authored Mar 30, 2022
1 parent b914160 commit 99767ea
Show file tree
Hide file tree
Showing 9 changed files with 130 additions and 77 deletions.
3 changes: 2 additions & 1 deletion case_copy.py
View file
Original file line number Diff line number Diff line change
Expand Up @@ -85,7 +85,8 @@
s_out = s
elif ('drv_in' in f and 'restart_n' in s and int(options.nyears) > 0):
s_out = ' restart_n = '+str(options.nyears)+'\n'
elif ('lnd_in' in f and 'finidat =' in s and int(options.finyr) > 0):
elif ('lnd_in' in f and 'finidat =' in s and not ("finidat = ' '" in s) and int(options.finyr) > 0):
print(s)
year_orig = str((s.split('.')[-2:-1])[0])[0:4]
year_new = str(10000+int(options.finyr))[1:]
s_out = s.replace('.r.'+year_orig, '.r.'+year_new)
Expand Down
2 changes: 2 additions & 0 deletions ensemble_copy.py
View file
Original file line number Diff line number Diff line change
Expand Up @@ -242,6 +242,8 @@
print('Creting netcdf variable for '+p)
param = nffun.getvar(myfile,'flnr')
param[:] = parm_values[pnum]
elif (p == 'psi50'):
param[:,parm_indices[pnum]] = parm_values[pnum]
elif (parm_indices[pnum] > 0):
param[parm_indices[pnum]] = parm_values[pnum]
elif (parm_indices[pnum] == 0):
Expand Down
8 changes: 4 additions & 4 deletions makepointdata.py
View file
Original file line number Diff line number Diff line change
Expand Up @@ -223,7 +223,7 @@
elif ('f19' in options.res):
longxy = (numpy.cumsum(numpy.ones([145]))-1)*2.5-1.25
latixy_centers = (numpy.cumsum(numpy.ones([96]))-1)*(180.0/95) - 90.0
latixy = numpy.zeros([97], numpy.float)
latixy = numpy.zeros([97], float)
longxy[0] = 0
latixy[0] = -90
latixy[96] = 90
Expand All @@ -232,7 +232,7 @@
elif ('f09' in options.res):
longxy = (numpy.cumsum(numpy.ones([289]))-1)*1.25-0.625
latixy_centers = (numpy.cumsum(numpy.ones([192]))-1)*(180.0/191) - 90.0
latixy = numpy.zeros([193], numpy.float)
latixy = numpy.zeros([193], float)
longxy[0] = 0
latixy[0] = -90
latixy[192] = 90
Expand Down Expand Up @@ -563,7 +563,7 @@
npft_crop = 10

#read file for site-specific PFT information
mypft_frac = numpy.zeros([npft+npft_crop], numpy.float)
mypft_frac = numpy.zeros([npft+npft_crop], float)
mypct_sand = 0.0
mypct_clay = 0.0

Expand Down Expand Up @@ -684,7 +684,7 @@
'C3 non-arctic grass', 'C4 grass', 'Crop','xxx','xxx']
if options.marsh and n==1: # Set tidal channel column in marsh mode to zero PFT area
print('Setting PFT area in tidal column to zero')
mypft_frac = numpy.zeros([npft+npft_crop], numpy.float)
mypft_frac = numpy.zeros([npft+npft_crop], float)
mypft_frac[0]=100.0
if (options.mypft >= 0 and not (options.marsh and n==1)):
print('Setting PFT '+str(options.mypft)+'('+pft_names[int(options.mypft)]+') to 100%')
Expand Down
24 changes: 20 additions & 4 deletions manage_ensemble.py
View file
Original file line number Diff line number Diff line change
Expand Up @@ -207,7 +207,10 @@ def postproc(myvars, myyear_start, myyear_end, myday_start, myday_end, myavg, \
else: #Regular parameter file
mydata = nffun.getvar(pfname,p)
if (int(ppfts[pnum]) > 0):
parms[pnum] = mydata[int(ppfts[pnum])]
if (p == 'psi50'):
parms[pnum] = mydata[0,int(ppfts[pnum])]
else:
parms[pnum] = mydata[int(ppfts[pnum])]
elif(int(ppfts[pnum]) <= 0):
try:
parms[pnum] = mydata[0]
Expand Down Expand Up @@ -240,7 +243,7 @@ def postproc(myvars, myyear_start, myyear_end, myday_start, myday_end, myavg, \
myobs=[]
myobs_err=[]
mytreatment=[]
time.sleep(rank)
time.sleep(rank*0.2)
postproc_input = open(options.postproc_file,'r')
data_cols = 0
for s in postproc_input:
Expand Down Expand Up @@ -361,13 +364,26 @@ def postproc(myvars, myyear_start, myyear_end, myday_start, myday_end, myavg, \
obs_out.close()
myoutput = open(UQ_output+'/data/pnames.txt', 'w')
eden_header=''
pnum=0
for p in pnames:
myoutput.write(p+'\n')
if ((pnum == 0 and pnames[pnum+1] == p) or p == pnames[pnum-1]):
myoutput.write(p+'_'+str(mypft[pnum])+'\n')
else:
myoutput.write(p+'\n')
pnum=pnum+1
eden_header=eden_header+p+','
myoutput.close()
myoutput = open(UQ_output+'/data/outnames.txt', 'w')
vlast=''
for v in myvars:
myoutput.write(v+'\n')
if v != vlast:
vcount=0
vlast=v
if (myvars.count(v) > 1):
myoutput.write(v+'_'+str(vcount)+'\n')
vcount=vcount+1
else:
myoutput.write(v+'\n')
eden_header=eden_header+v+','
myoutput.close()
myoutput = open(UQ_output+'/data/param_range.txt', 'w')
Expand Down
15 changes: 7 additions & 8 deletions model_surrogate.py
View file
Original file line number Diff line number Diff line change
Expand Up @@ -12,13 +12,13 @@ def __init__(self,case=''):
self.nparms = self.ptrain.shape[1]
self.nobs = self.ytrain.shape[1]
self.ntrain=self.ptrain.shape[0]
self.yrange = np.zeros([2,self.nobs],np.float)
self.yrange = np.zeros([2,self.nobs],float)
for i in range(0,self.nobs):
self.yrange[0,i] = min(self.ytrain[:,i])
self.yrange[1,i] = max(self.ytrain[:,i])

self.pmin = np.zeros([self.nparms], np.float)
self.pmax= np.zeros([self.nparms], np.float)
self.pmin = np.zeros([self.nparms], float)
self.pmax= np.zeros([self.nparms], float)
for i in range(0,self.nparms):
self.pmin[i] = min(self.ptrain[:,i])
self.pmax[i] = max(self.ptrain[:,i])
Expand All @@ -29,12 +29,11 @@ def __init__(self,case=''):
self.parm_names.append(s[:-1])
pnamefile.close()
print(self.parm_names, self.nobs, self.nparms)
self.pdef = np.zeros([self.nparms], np.float)
self.pdef = np.zeros([self.nparms], float)
for i in range(0,self.nparms):
self.pdef[i] = (self.pmin[i]+self.pmax[i])/2

self.obs = np.zeros([self.nobs],np.float)
self.obs_err = np.zeros([self.nobs],np.float)
self.obs = np.zeros([self.nobs],float)
self.obs_err = np.zeros([self.nobs],float)
self.obs_name = []
obsfile = open(UQdir+'/data/obs.dat')
i=0
Expand All @@ -53,7 +52,7 @@ def __init__(self,case=''):
self.qoi_good = np.loadtxt(UQdir+'/NN_surrogate/qoi_good.txt').astype(int)

def run(self,parms):
parms_nn = np.zeros([1,self.nparms],np.float)
parms_nn = np.zeros([1,self.nparms],float)
for p in range(0,self.nparms):
parms_nn[0,p] = (parms[p]-self.pmin[p])/(self.pmax[p]-self.pmin[p])
self.output = np.zeros([self.nobs])
Expand Down
68 changes: 34 additions & 34 deletions plotcase.py
View file
Original file line number Diff line number Diff line change
Expand Up @@ -217,7 +217,7 @@ def getvar(fname, varname, npf, index, scale_factor):
#sys.exit(0)

obs = options.myobs
myobsdir = '/lustre/or-hydra/cades-ccsi/scratch/dmricciuto/fluxnet'
myobsdir = '/home/ac.ricciuto/fluxnet'

#get list of variables from varfile
myvars=[]
Expand Down Expand Up @@ -257,10 +257,10 @@ def getvar(fname, varname, npf, index, scale_factor):

#dirs=[]
nvar = len(myvars)
x_toplot = numpy.zeros([ncases, 2000000], numpy.float)
data_toplot = numpy.zeros([ncases, nvar, 2000000], numpy.float)
obs_toplot = numpy.zeros([ncases, nvar, 2000000], numpy.float)+numpy.NaN
err_toplot = numpy.zeros([ncases, nvar, 2000000], numpy.float)+numpy.NaN
x_toplot = numpy.zeros([ncases, 2000000], float)
data_toplot = numpy.zeros([ncases, nvar, 2000000], float)
obs_toplot = numpy.zeros([ncases, nvar, 2000000], float)+numpy.NaN
err_toplot = numpy.zeros([ncases, nvar, 2000000], float)+numpy.NaN
snum = numpy.zeros([ncases], numpy.int)

for c in range(0,ncases):
Expand Down Expand Up @@ -344,10 +344,10 @@ def getvar(fname, varname, npf, index, scale_factor):
nypf = max(1, nhtot/8760)

#initialize data arrays
mydata = numpy.zeros([nvar,2000000], numpy.float)
myobs = numpy.zeros([nvar,2000000], numpy.float)+numpy.NaN
myerr = numpy.zeros([nvar,2000000], numpy.float)+numpy.NaN
x = numpy.zeros([2000000], numpy.float)
mydata = numpy.zeros([nvar,2000000], float)
myobs = numpy.zeros([nvar,2000000], float)+numpy.NaN
myerr = numpy.zeros([nvar,2000000], float)+numpy.NaN
x = numpy.zeros([2000000], float)
nsteps=0

if (c == 0):
Expand All @@ -359,7 +359,7 @@ def getvar(fname, varname, npf, index, scale_factor):
if (obs):
myobsfiles = os.listdir(myobsdir+'/'+mytstep+'/')
for f in myobsfiles:
if mysites[c] in f and '.csv' in f:
if mysites[0] in f and '.csv' in f:
myobsfile = myobsdir+'/'+mytstep+'/'+f
avpd_obs = 1
if (mytstep == 'halfhourly' and '_HH_' in myobsfile):
Expand Down Expand Up @@ -391,40 +391,40 @@ def getvar(fname, varname, npf, index, scale_factor):
if (thisob == 0):
thisob = (int(myvals[0][0:4])-ystart)*npy*avpd_obs
if (thisob % avpd_obs == 0):
myobs[v,thisob/avpd_obs] = 0.0
myerr[v,thisob/avpd_obs] = 0.0
myobs[v,int(thisob/avpd_obs)] = 0.0
myerr[v,int(thisob/avpd_obs)] = 0.0
for h in header:
if (h.strip() == 'NEE_CUT_REF' and 'NEE' in myvars[v]):
myobs[v,thisob/avpd_obs] = myobs[v,thisob/avpd_obs] + float(myvals[thiscol])/avpd_obs
elif (h.strip () == 'NEE_CUT_REF_JOINTUNC' and \
'NEE' in myvars[v]):
myerr[v,thisob/avpd_obs] = myerr[v,thisob/avpd_obs] +float(myvals[thiscol])/avpd_obs
elif (h.strip() == 'GPP_NT_CUT_REF' and 'GPP' in myvars[v]):
myobs[v,thisob/avpd_obs] = myobs[v,thisob/avpd_obs] +float(myvals[thiscol])/avpd_obs
elif (h.strip() == 'GPP_NT_CUT_SE' and 'GPP' in myvars[v]):
myerr[v,thisob/avpd_obs] = myerr[v,thisob/avpd_obs] +float(myvals[thiscol])/avpd_obs
elif (h.strip() == 'GPP_NT_CUT_REF' and ('GPP' in myvars[v] or 'FPSN' in myvars[v])):
myobs[v,int(thisob/avpd_obs)] = myobs[v,int(thisob/avpd_obs)] +float(myvals[thiscol])/avpd_obs
elif (h.strip() == 'GPP_NT_CUT_SE' and ('GPP' in myvars[v] or 'FPSN' in myvars[v])):
myerr[v,int(thisob/avpd_obs)] = myerr[v,int(thisob/avpd_obs)] +float(myvals[thiscol])/avpd_obs
elif (h.strip() == 'RECO_NT_CUT_REF' and 'ER' in myvars[v]):
myobs[v,thisob/avpd_obs] = myobs[v,thisob/avpd_obs] +float(myvals[thiscol])/avpd_obs
elif (h.strip() == 'RECO_NET_CUT_SE' and 'ER' in myvars[v]):
myerr[v,thisob/avpd_obs] = myerr[v,thisob/avpd_obs] +float(myvals[thiscol])/avpd_obs
elif (h.strip() == 'LE_F_MDS' and 'EFLX_LH_TOT' in myvars[v]):
myobs[v,thisob/avpd_obs] = myobs[v,thisob/avpd_obs] +float(myvals[thiscol])/avpd_obs
myobs[v,int(thisob/avpd_obs)] = myobs[v,int(thisob/avpd_obs)] +float(myvals[thiscol])/avpd_obs
elif (h.strip() == 'LE_RANDUNC' and 'EFLX_LH_TOT' in myvars[v]):
myerr[v,thisob/avpd_obs] = myerr[v,thisob/avpd_obs] +float(myvals[thiscol])/avpd_obs
myerr[v,int(thisob/avpd_obs)] = myerr[v,int(thisob/avpd_obs)] +float(myvals[thiscol])/avpd_obs
elif (h.strip() == 'H_F_MDS' and 'FSH' in myvars[v]):
myobs[v,thisob/avpd_obs] = myobs[v,thisob/avpd_obs] +float(myvals[thiscol])/avpd_obs
elif (h.strip() == 'H_RANDUNC' and 'FSH' in myvars[v]):
myerr[v,thisob/avpd_obs] = myerr[v,thisob/avpd_obs] +float(myvals[thiscol])/avpd_obs
elif (h.strip() == 'TA_F_MDS' and 'TBOT' in myvars[v]):
myobs[v,thisob/avpd_obs] = myobs[v,thisob/avpd_obs] +float(myvals[thiscol])/avpd_obs
myobs[v,int(thisob/avpd_obs)] = myobs[v,thisob/avpd_obs] +float(myvals[thiscol])/avpd_obs
elif (h.strip() == 'SWIN_F_MDS' and 'FSDS' in myvars[v]):
myobs[v,thisob/avpd_obs] = myobs[v,thisob/avpd_obs] +float(myvals[thiscol])/avpd_obs
elif (h.strip() == 'WS_F' and 'WIND' in myvars[v]):
myobs[v,thisob/avpd_obs] = myobs[v,thisob/avpd_obs] +float(myvals[thiscol])/avpd_obs
elif (h.strip() == 'P_F' and 'RAIN' in myvars[v]):
myobs[v,thisob/avpd_obs] = myobs[v,thisob/avpd_obs] +float(myvals[thiscol])/avpd_obs
if myobs[v,thisob/avpd_obs] < -4000:
myobs[v,thisob/avpd_obs] = numpy.NaN
if myobs[v,int(thisob/avpd_obs)] < -4000:
myobs[v,int(thisob/avpd_obs)] = numpy.NaN
thiscol=thiscol+1
thisob=thisob+1
thisrow = thisrow+1
Expand Down Expand Up @@ -618,10 +618,10 @@ def getvar(fname, varname, npf, index, scale_factor):
if (avtype == 'diurnal'):
snum[c]=24
for v in range(0,nvar):
mysum = numpy.zeros(snum[c], numpy.float)
mysum_obs = numpy.zeros(snum[c], numpy.float)
myct = numpy.zeros(snum[c],numpy.float)
myct_obs = numpy.zeros(snum[c],numpy.float)
mysum = numpy.zeros(snum[c], float)
mysum_obs = numpy.zeros(snum[c], float)
myct = numpy.zeros(snum[c],float)
myct_obs = numpy.zeros(snum[c],float)
for y in range(0,(yend_all-ystart+1)):
for d in range (int(options.dstart),int(options.dend)):
for s in range(0,snum[c]):
Expand All @@ -646,8 +646,8 @@ def getvar(fname, varname, npf, index, scale_factor):
if (avtype == 'seasonal'):
for v in range(0,nvar):
snum[c] = 12
mysum=numpy.zeros(snum[c], numpy.float)
mysum_obs = numpy.zeros(snum[c], numpy.float)
mysum=numpy.zeros(snum[c], float)
mysum_obs = numpy.zeros(snum[c], float)
mycount_obs = numpy.zeros(snum[c], numpy.int)
for y in range(0,(yend_all-ystart+1)):
for s in range(0,snum[c]):
Expand Down Expand Up @@ -678,9 +678,9 @@ def getvar(fname, varname, npf, index, scale_factor):
else:
analysis_type=mytstep

rmse = numpy.zeros([len(myvars),ncases],numpy.float)
bias = numpy.zeros([len(myvars),ncases],numpy.float)
corr = numpy.zeros([len(myvars),ncases],numpy.float)
rmse = numpy.zeros([len(myvars),ncases],float)
bias = numpy.zeros([len(myvars),ncases],float)
corr = numpy.zeros([len(myvars),ncases],float)

options.nperpage=int(options.nperpage)
if (options.nperpage > 1):
Expand All @@ -695,7 +695,7 @@ def getvar(fname, varname, npf, index, scale_factor):
if (v % options.nperpage == 0):
fig = plt.figure(figsize=(11,8.5))
thisfig = v % options.nperpage+1
fignum = v/options.nperpage
fignum = int(v/options.nperpage)
thiscol = (thisfig -1) % ncol
thisrow = (thisfig -1) / ncol
ax = plt.subplot(nrow*100+ncol*10+thisfig)
Expand Down Expand Up @@ -779,7 +779,7 @@ def getvar(fname, varname, npf, index, scale_factor):
linestyle=styles[c], linewidth=3)
if (options.myobs and c == 0):
ax.errorbar(x_toplot[c, 0:snum[c]], obs_toplot[c,v,0:snum[c]], yerr = err_toplot[c,v,0:snum[c]], \
color=colors[c], fmt='o')
color='k', fmt='o')

if (options.noplot == False):
if (thisrow+1 == nrow):
Expand All @@ -793,8 +793,8 @@ def getvar(fname, varname, npf, index, scale_factor):
plt.ylabel(myvars[v]+' ('+var_units[v]+')')
box = ax.get_position()
ax.set_position([box.x0, box.y0, box.width * 0.8, box.height])
#if (v % options.nperpage == options.nperpage-1):
# ax.legend(loc='center left', bbox_to_anchor=(1, 0.5),prop={'size': 10})
if (v % options.nperpage == options.nperpage-1):
ax.legend(loc='center left', bbox_to_anchor=(1, 0.5),prop={'size': 10})
if (v % options.nperpage == 0):
ax.legend(loc='upper left', bbox_to_anchor=(0, 1.3), prop={'size': 8}, ncol=2)
plt.title(var_long_names[v]) #+' at '+mysites[0])
Expand Down
Loading

0 comments on commit 99767ea

Please sign in to comment.