import pandas as pd import xarray as xr import numpy as np import glob def get_ens_ACCESS(exp_name,var_name,yrstart=2040,yrend=2060,concat=True): ''' read all 10 realizations of ACCESS return concatenated DataArray ''' var_concat = [] for ireal in range(1,11): dir_name = '/scratch/usr/shkhwwey/CDRMIP/ACCESS-ESM1-5/' + exp_name + '/' + var_name + '/real/' file_list = glob.glob(dir_name + '/*' + 'r' + str(ireal) + 'i*.nc') file_list.sort() print(file_list) ds = xr.open_mfdataset(file_list) yrini = 2015 istart = (yrstart-yrini) * 12 iend = (yrend-yrini)*12 var = ds[var_name][istart:iend+12] var_concat.append(var) if concat is True: return xr.concat(var_concat,dim='ireal') if concat is False: return var_concat def ens_values(df_list,var_name): # df_list: a list of DataFrame yrstart = 2015 diff_realizaions = [] # (nreal,nyear) for ireal in range(len(df_list)): diff_realizaions.append([]) for iyear in range(yrstart,2100): for ireal in range(len(df_list)): diff_realizaions[ireal].append(df_list[ireal][var_name].loc[iyear]) diff_realizaions = np.array(diff_realizaions) ens_mean = np.average(diff_realizaions,axis=0) ens_min = np.min(diff_realizaions,axis=0) ens_max = np.max(diff_realizaions,axis=0) ens_stddev = np.std(diff_realizaions,axis=0) yrindex = np.array(range(yrstart,2100)) df_ens = pd.DataFrame(data=ens_max - ens_min, index=yrindex, columns=['ens_range']) df_ens['ens_mean'] = ens_mean df_ens['ens_min'] = ens_min df_ens['ens_max'] = ens_max return df_ens def cOcean_from_fgco2(df,yrstart): cOcean = [] for year in range(yrstart,2100): cOcean.append(df['fgco2'].loc[yrstart:year].sum()) cOcean = np.array(cOcean) cOcean = pd.DataFrame(data=cOcean, index=range(yrstart,2100), columns=['cOcean']) cOcean.index.name = 'year' return cOcean class Figure1(): def __init__(self, var_name, exp): self.var_name = var_name self.exp = exp self.expnamelist = [] self.dflist = [] def read_FOCI(self,debug=False): self.model = 'FOCI' var_name = self.var_name model, exp = self.model, self.exp diff = True self.expnamelist = [] self.expnamelist.append(['CC104','CC105','CC106']) # 0 REF self.expnamelist.append(['HW002','HW005','HW006']) # 1 OAE self.expnamelist.append(['HW003','HW007','HW008']) # 2 AFF self.expnamelist.append(['HW013','HW012','HW011']) # 3 SynTraREF self.expnamelist.append(['HW022','HW023','HW024']) # 4 SynTraAR self.expnamelist.append(['HW028','HW032','HW033']) # 5 SynTraOAE self.expnamelist.append(['HW029','HW030','HW031']) # 6 SynTraBOTH self.expnamelist.append(['HW034','HW035','HW036']) # 7 SynTraHalfOAE self.expnamelist.append(['HW037','HW038','HW039']) # 8 SynTraHalfBOTH self.expnamelist.append(['HW040','HW041','HW042']) # 9 SynTraHalfAR if exp=='OAE': save_name = var_name + '_OAE' self.yrstart = 2015 if exp=='AFF': save_name = var_name + '_AFF' self.yrstart = 2015 if exp=='BOTH': save_name = var_name + '_BOTH' self.dflist = [] for iexpset in range(len(self.expnamelist)): self.dflist.append([]) for expname in self.expnamelist[iexpset]: if var_name=='cLand': try: self.dflist[iexpset].append( pd.read_csv('/home/shkhwwey/OAE/csv/' + expname + '_combined.csv',index_col='year')) except: try: self.dflist[iexpset].append( pd.read_csv('/home/shkhwwey/OAE/csv/' + expname + '_cLand.csv',index_col='year')) except: if debug is True: print(f'Error reading {expname} {var_name}!') else: try: self.dflist[iexpset].append( pd.read_csv('/home/shkhwwey/OAE/csv/' + expname + '_' + var_name + '.csv',index_col='year')) except: if debug is True: print(f'Error reading {expname} {var_name}!') def read_MPIESM(self,debug=False): self.model = 'MPI-ESM' var_name = self.var_name model, exp = self.model, self.exp diff = True self.expnamelist = [] if exp=='SynTra': self.expnamelist.append(['ymo1001','ymo1002','ymo1003']) # 0 REF self.expnamelist.append(['ymo1021','ymo1022','ymo1045']) # 1 AR self.expnamelist.append(['ymo1089','ymo1090','ymo1091']) # 2 OAE self.expnamelist.append(['ymo1092','ymo1093','ymo1094']) # 3 BOTH if exp=='OAE': # new experiments self.expnamelist.append(['essp585_0003','essp585_0004','essp585_0005']) self.expnamelist.append(['esm-ssp585-ocn-alk_EXP1','esm-ssp585-ocn-alk_EXP2','esm-ssp585-ocn-alk_EXP3']) # old experiments with wrong OAE intensity (too high in high latitudes) #self.expnamelist.append(['essp585_0001','essp585_0002','essp585_0003']) #self.expnamelist.append(['esm-ssp585-ocn-alk-r1','esm-ssp585-ocn-alk-r2','esm-ssp585-ocn-alk-r3']) self.expnamelist.append([]) self.yrstart = 2020 if exp=='AFF': self.expnamelist.append(['essp585_0004','essp585_0005','esm-ssp585']) self.expnamelist.append([]) self.expnamelist.append(['ymo_esm_ssp585_ssp126Lu_2','ymo_esm_ssp585_ssp126Lu_3','esm-ssp585-ssp126Lu']) self.yrstart = 2015 if exp=='BOTH': # new experiments self.expnamelist.append(['essp585_0003','essp585_0004','essp585_0005']) self.expnamelist.append(['esm-ssp585-ocn-alk_EXP1','esm-ssp585-ocn-alk_EXP2','esm-ssp585-ocn-alk_EXP3']) # old experiments #self.expnamelist.append(['essp585_0001','essp585_0002','essp585_0003']) #self.expnamelist.append(['esm-ssp585-ocn-alk-r1','esm-ssp585-ocn-alk-r2','esm-ssp585-ocn-alk-r3']) self.expnamelist.append(['essp585_0004','essp585_0005','esm-ssp585']) self.expnamelist.append(['ymo_esm_ssp585_ssp126Lu_2','ymo_esm_ssp585_ssp126Lu_3','esm-ssp585-ssp126Lu']) save_name = var_name + '_BOTH' if exp!='SynTra': base_dir = '/home/shkhwwey/OAE/csv/MPI-ESM1-2-LR/' if exp=='SynTra': base_dir = '/home/shkhwwey/OAE/SynTra_analysis/csv/MPIESM_csv/' self.dflist = [] for iexpset in range(len(self.expnamelist)): self.dflist.append([]) for expname in self.expnamelist[iexpset]: if var_name=='cLand': # self.dflist[iexpset].append( pd.read_csv(base_dir + 'MPI-ESM1-2-LR_' + expname + '_combined.csv',index_col='year')) self.dflist[iexpset].append( pd.read_csv(base_dir + 'MPI-ESM1-2-LR_' + expname + '_' + var_name + '.csv',index_col='year')) else: try: self.dflist[iexpset].append( pd.read_csv(base_dir + 'MPI-ESM1-2-LR_' + expname + '_' + var_name + '.csv',index_col='year')) except: if debug is True: print(f'Error reading {expname} {var_name}!') def read_NorESM_ensemble(self,debug=False): base_dir = '/home/shkhwwey/OAE/CDRMIP_csv/' self.model = 'NorESM2-LM' var_name = self.var_name self.expnamelist = ['esm-ssp585','esm-ssp585-ocn-alk'] for iexpset in range(len(self.expnamelist)): self.dflist.append([]) expname = self.expnamelist[iexpset] if expname=='esm-ssp585': real_list = [1] for ireal in real_list: try: self.dflist[iexpset].append(pd.read_csv(base_dir + 'NorESM2-LM_' + expname + '_' + var_name + '.csv',index_col='year')) except: if debug is True: print(f'Error reading {expname} {var_name} r{ireal}!') elif expname=='esm-ssp585-ocn-alk': real_list = [1,2,3] for ireal in real_list: try: self.dflist[iexpset].append(pd.read_csv(base_dir + 'NorESM2-LM_' + expname + '_' + var_name + '_r' + str(ireal) + '.csv',index_col='year')) except: if debug is True: print(f'Error reading {expname} {var_name} r{ireal}!') def read_ACCESS_ensemble(self,debug=False): base_dir = '/home/shkhwwey/OAE/CDRMIP_csv/' self.model = 'ACCESS-ESM1-5' var_name = self.var_name self.expnamelist = ['esm-ssp585','esm-ssp585-ssp126Lu'] for iexpset in range(len(self.expnamelist)): self.dflist.append([]) expname = self.expnamelist[iexpset] for ireal in range(1,11): try: self.dflist[iexpset].append(pd.read_csv(base_dir + 'ACCESS-ESM1-5_' + expname + '_' + var_name + '_r' + str(ireal) + '.csv',index_col='year')) except: if debug is True: print(f'Error reading {expname} {var_name} r{ireal}!') def read_model(self,model,debug=False): base_dir = '/home/shkhwwey/OAE/' self.model = model exp = self.exp var_name = self.var_name diff = True if exp=='OAE': save_name = var_name + '_OAE' self.expnamelist = [['esm-ssp585'],['esm-ssp585-ocn-alk'],['esm-ssp585-ssp126Lu']] #self.yrstart = 2015 if exp=='AFF': save_name = var_name + '_AFF' self.expnamelist = [['esm-ssp585'],['esm-ssp585-ocn-alk'],['esm-ssp585-ssp126Lu']] #self.yrstart = 2015 if exp=='BOTH': save_name = var_name + '_BOTH' self.dflist = [] for iexpset in range(len(self.expnamelist)): self.dflist.append([]) for expname in self.expnamelist[iexpset]: try: self.dflist[iexpset].append(pd.read_csv(base_dir + 'CDRMIP_csv/' + model + '_' + expname + '_' + var_name + '.csv',index_col='year')) except: if debug is True: print(f' Error reading {expname} {var_name} for {model}!') def read_NorESM(self): base_dir = '/home/shkhwwey/OAE/' self.model = 'NorESM' model, exp = self.model, self.exp var_name = self.var_name diff = True if exp=='OAE': save_name = var_name + '_OAE' self.expnamelist = [['esm-ssp585'],['esm-ssp585-ocn-alk'],['esm-ssp585-ssp126Lu']] self.yrstart = 2015 if exp=='AFF': save_name = var_name + '_AFF' self.expnamelist = [['esm-ssp585'],['esm-ssp585-ocn-alk'],['esm-ssp585-ssp126Lu']] self.yrstart = 2015 if exp=='BOTH': save_name = var_name + '_BOTH' self.dflist = [] for iexpset in range(len(self.expnamelist)): self.dflist.append([]) for expname in self.expnamelist[iexpset]: try: self.dflist[iexpset].append( pd.read_csv(base_dir + 'CDRMIP_csv/NorESM2-LM_' + expname + '_' + var_name + '.csv',index_col='year')) except: print(f'Error reading {expname} {var_name}!') def read_UKESM(self): self.model = 'UKESM' model, exp = self.model, self.exp var_name = self.var_name diff = True if exp=='OAE': save_name = var_name + '_OAE' self.expnamelist = [['esm-ssp585'],['esm-ssp585-ocn-alk'],['esm-ssp585-ssp126Lu']] self.yrstart = 2020 if exp=='AFF': save_name = var_name + '_AFF' self.expnamelist = [['esm-ssp585'],[],['esm-ssp585-ssp126Lu']] self.yrstart = 2015 if exp=='BOTH': save_name = var_name + '_BOTH' self.dflist = [] for iexpset in range(len(self.expnamelist)): self.dflist.append([]) for expname in self.expnamelist[iexpset]: if var_name in ['cLand','co2_concentration']: self.dflist[iexpset].append( pd.read_csv('csv/UKESM_' + expname + '_combined.csv',index_col='year')) elif var_name == 'tree': self.dflist[iexpset].append( pd.read_csv('csv/UKESM1-0-LL/UKESM1-0-LL_' + expname + '_' + var_name + '.csv',index_col='year')) else: try: self.dflist[iexpset].append( pd.read_csv('csv/UKESM_' + expname + '_' + var_name + '.csv',index_col='year')) except: print(f'Error reading {expname} {var_name}!') # if exp=='OAE': # ax.plot(self.dflist[1][0].index,self.dflist[1][0][var_name].loc[:]-self.dflist[0][0][var_name].loc[self.dflist[1][0].index],'-',linewidth=3,color='C2',label=model) # if exp=='AFF': # ax.plot(self.dflist[2][0].index,self.dflist[2][0][var_name].loc[:]-self.dflist[0][0][var_name].loc[self.dflist[2][0].index],'-',linewidth=3,color='C2',label=model) def plot_individual_run(self): model = self.model expnamelist = self.expnamelist dflist = self.dflist var_name = self.var_name diff = True exp = self.exp print('plot_individual_run()',model,var_name) for iexpset in range(len(expnamelist)): # how many exp sets? if diff and exp == 'OAE': if not (iexpset == 1): continue if diff and exp == 'AFF': if not (iexpset == 2): continue for iexp in range(len(expnamelist[0])): # how many realizations? # if not diff: # ax.plot(dflist[iexp].index,dflist[iexp][var_name],'--',color=collist[iexp],linewidth=1) ### if diff: if var_name!='cLand': if model=='FOCI': ax.plot(dflist[iexpset][iexp].index,dflist[iexpset][iexp][var_name].loc[:]-dflist[0][iexp][var_name].loc[dflist[iexpset][iexp].index],'--',linewidth=1,color='C0') if model=='MPI-ESM': ax.plot(dflist[iexpset][iexp].index,dflist[iexpset][iexp][var_name].loc[:]-dflist[0][iexp][var_name].loc[dflist[iexpset][iexp].index],'--',linewidth=1,color='C1') # print(f'{expnamelist[iexp]} {(dflist[iexp][var_name].loc[2090:2099]-dflist[iexp%3][var_name].loc[2090:2099]).mean():.3f} ++++++ individual ens member') # for itime in range (2,10): # temp1 = np.average(dflist[iexp][var_name].loc[2000+itime*10:2000+itime*10+10]) # temp2 = np.average(dflist[iexp%3][var_name].loc[2000+itime*10:2000+itime*10+10]) # temp1 = np.average(dflist[iexp][var_name].loc[2051:2060]) # temp2 = np.average(dflist[iexp%3][var_name].loc[2051:2060]) # print(f'2051-2060, EXP{iexp}, EXP: {temp1:.2f}, CTR: {temp2:.2f}, Diff: {(temp1-temp2):.2f}') # temp1 = np.average(dflist[iexp][var_name].loc[2071:2080]) # temp2 = np.average(dflist[iexp%3][var_name].loc[2071:2080]) # print(f'2071-2080, EXP{iexp}, EXP: {temp1:.2f}, CTR: {temp2:.2f}, Diff: {(temp1-temp2):.2f}') # print('---') if var_name =='cLand': ax.plot(dflist[iexp].index,dflist[iexp]['cLand'].loc[:]-dflist[iexp%3]['cLand'].loc[dflist[iexp].index],patlist[2],color=collist[iexp//3*3]) cOceanCTR = [] for iyear in dflist[0].index: cOceanCTR.append(sum(dflist[iexp%3]['fgco2'].loc[2015:iyear])) cOceanCTR = np.array(cOceanCTR) if dflist[iexp].index[0]==2015: cOcean = [] for iyear in dflist[iexp].index: cOcean.append(sum(dflist[iexp]['fgco2'].loc[2015:iyear])) cOcean = np.array(cOcean) ax.plot(dflist[iexp].index,cOcean-cOceanCTR,patlist[1],color=collist[iexp//3*3]) if dflist[iexp].index[0]==2070: cOcean = [] for iyear in dflist[0].index: if iyear<=2069: cOcean.append(sum(dflist[iexp-6]['fgco2'].loc[2015:iyear])) if iyear>=2070: cOcean.append(cOcean[-1] + dflist[iexp]['fgco2'].loc[iyear]) cOcean = np.array(cOcean) ax.plot(dflist[iexp].index,cOcean[55:85]-cOceanCTR[55:85],patlist[1],color=collist[iexp//3*3]) def plot_uncertainty(self): '''plot uncertainty range as shading area''' model = self.model expnamelist = self.expnamelist dflist = self.dflist var_name = self.var_name diff = True yrstart = self.yrstart exp = self.exp if (diff)and (exp=='BOTH') and (var_name=='cLand'): ### cLand of ocn_alk ens_min = [] ens_max = [] ens_mean = [] for iyear in range(yrstart,2100): diff1 = dflist[3]['cLand'].loc[iyear] - dflist[0]['cLand'].loc[iyear] diff2 = dflist[4]['cLand'].loc[iyear] - dflist[1]['cLand'].loc[iyear] diff3 = dflist[5]['cLand'].loc[iyear] - dflist[2]['cLand'].loc[iyear] ens_mean.append((diff1 + diff2 + diff3) / 3) ens_min.append( min(diff1,diff2,diff3) ) ens_max.append( max(diff1,diff2,diff3) ) ens_min = np.array(ens_min) ens_max = np.array(ens_max) ens_mean = np.array(ens_mean) ax.fill_between(dflist[0].index, ens_min, ens_max ,alpha=0.3, facecolor=collist[3]) print(f'OAE cLand ens_mean 2060-2015: {ens_mean[45]-ens_mean[0]:.2f}') ### cLand of A/R ens_min = [] ens_max = [] ens_mean = [] for iyear in range(yrstart,2100): diff1 = dflist[6]['cLand'].loc[iyear] - dflist[0]['cLand'].loc[iyear] diff2 = dflist[7]['cLand'].loc[iyear] - dflist[1]['cLand'].loc[iyear] diff3 = dflist[8]['cLand'].loc[iyear] - dflist[2]['cLand'].loc[iyear] ens_mean.append((diff1 + diff2 + diff3) / 3) ens_min.append( min(diff1,diff2,diff3) ) ens_max.append( max(diff1,diff2,diff3) ) ens_min = np.array(ens_min) ens_max = np.array(ens_max) ens_mean = np.array(ens_mean) ax.fill_between(dflist[0].index, ens_min, ens_max ,alpha=0.3, facecolor=collist[6]) print(f'A/R cLand ens_mean 2060-2015: {ens_mean[45]-ens_mean[0]:.2f}') ### cOcean of ocn_alk and A/R ens_min = [] ens_max = [] ens_mean = [] ens_min2 = [] ens_max2 = [] ens_mean2 = [] ens_min3 = [] ens_max3 = [] ens_mean3 = [] for iyear in dflist[0].index: cOceanCTR1 = (sum(dflist[0]['fgco2'].loc[2015:iyear])) cOceanCTR2 = (sum(dflist[1]['fgco2'].loc[2015:iyear])) cOceanCTR3 = (sum(dflist[2]['fgco2'].loc[2015:iyear])) cOcean1 = (sum(dflist[3]['fgco2'].loc[2015:iyear])) cOcean2 = (sum(dflist[4]['fgco2'].loc[2015:iyear])) cOcean3 = (sum(dflist[5]['fgco2'].loc[2015:iyear])) diff1 = cOcean1 - cOceanCTR1 diff2 = cOcean2 - cOceanCTR2 diff3 = cOcean3 - cOceanCTR3 ens_mean.append((diff1 + diff2 + diff3) / 3) ens_min.append( min(diff1,diff2,diff3) ) ens_max.append( max(diff1,diff2,diff3) ) ### ocn_alk_stop if OAE_STOP: if iyear>=2070: cOcean_stop1 = ( (sum(dflist[3]['fgco2'].loc[2015:2069])) + sum(dflist[9]['fgco2'].loc[2070:iyear] )) cOcean_stop2 = ( (sum(dflist[4]['fgco2'].loc[2015:2069])) + sum(dflist[10]['fgco2'].loc[2070:iyear])) cOcean_stop3 = ( (sum(dflist[5]['fgco2'].loc[2015:2069])) + sum(dflist[11]['fgco2'].loc[2070:iyear])) diff1 = cOcean_stop1 - cOceanCTR1 diff2 = cOcean_stop2 - cOceanCTR2 diff3 = cOcean_stop3 - cOceanCTR3 ens_mean3.append((diff1 + diff2 + diff3) / 3) ens_min3.append( min(diff1,diff2,diff3) ) ens_max3.append( max(diff1,diff2,diff3) ) ### A/R cOcean1 = (sum(dflist[6]['fgco2'].loc[2015:iyear])) cOcean2 = (sum(dflist[7]['fgco2'].loc[2015:iyear])) cOcean3 = (sum(dflist[8]['fgco2'].loc[2015:iyear])) diff1 = cOcean1 - cOceanCTR1 diff2 = cOcean2 - cOceanCTR2 diff3 = cOcean3 - cOceanCTR3 ens_mean2.append((diff1 + diff2 + diff3) / 3) ens_min2.append( min(diff1,diff2,diff3) ) ens_max2.append( max(diff1,diff2,diff3) ) ens_min = np.array(ens_min) ens_max = np.array(ens_max) ens_mean = np.array(ens_mean) ens_min2 = np.array(ens_min2) ens_max2 = np.array(ens_max2) ens_mean2 = np.array(ens_mean2) if OAE_STOP: ens_min3 = np.array(ens_min3) ens_max3 = np.array(ens_max3) ens_mean3 = np.array(ens_mean3) ax.fill_between(dflist[9].index, ens_min3, ens_max3 ,alpha=0.3, facecolor=collist[9],label='OAE-stop') ax.fill_between(dflist[0].index, ens_min, ens_max ,alpha=0.3, facecolor=collist[3],label='OAE') ax.fill_between(dflist[0].index, ens_min2, ens_max2 ,alpha=0.3, facecolor=collist[6],label='A/R') print(f'OAE cOcean ens_mean 2060-2015: {ens_mean[45]-ens_mean[0]:.2f}') print(f'A/R cOcean ens_mean 2060-2015: {ens_mean2[45]-ens_mean2[0]:.2f}') ##### end if (diff)and (exp=='BOTH') and (var_name=='cLand') if (diff) and (var_name!='cLand'): for iexpset in range(1,len(expnamelist[0])): if exp=='AFF' and iexpset==1: continue if exp=='OAE' and iexpset==2: continue # yrindex = [] diff_realizaions = [] # (nreal,nyear) for ireal in range(len(expnamelist[0])): diff_realizaions.append([]) for iyear in range(yrstart,2100): for ireal in range(len(expnamelist[0])): diff_realizaions[ireal].append(dflist[iexpset][ireal][var_name].loc[iyear] - dflist[0][ireal][var_name].loc[iyear]) # yrindex.append(iyear) diff_realizaions = np.array(diff_realizaions) ens_mean = np.average(diff_realizaions,axis=0) ens_min = np.min(diff_realizaions,axis=0) ens_max = np.max(diff_realizaions,axis=0) ens_stddev = np.std(diff_realizaions,axis=0) yrindex = np.array(range(yrstart,2100)) print(ens_mean.shape,ens_min.shape) ###### # mean of ensemble range from 2050 to 2059 df_ens = pd.DataFrame(data=ens_max - ens_min, index=yrindex, columns=['ens_range']) df_ens['ens_mean'] = ens_mean mean_range = df_ens['ens_range'].loc[2050:2059].mean() print('mean of range from 2050 to 2059') print(f'OAE {mean_range:.3f}') ###### # ax.plot(range(yrstart,2100),ens_mean,'-',color=collist[3],linewidth=3) # ax.fill_between(range(yrstart,2100), ens_min, ens_max ,alpha=0.3, facecolor=collist[3],label='OAE') if model=='MPI-ESM': ax.plot(range(yrstart,2100),ens_mean,'-',color='C1',linewidth=3,label=model) ax.fill_between(range(yrstart,2100), ens_min, ens_max ,alpha=0.3, facecolor='C1') if model=='FOCI': ax.plot(range(yrstart,2100),ens_mean,'-',color='C0',linewidth=3,label=model) ax.fill_between(range(yrstart,2100), ens_min, ens_max ,alpha=0.3, facecolor='C0') print(f'iexpset: {iexpset}') print(f'ens_mean 2060-2015: {ens_mean[45]-ens_mean[0]:.3f}') temp = df_ens['ens_mean'].loc[2050:2060].mean() print(f'ens_mean from 2050 to 2060 df {temp:.6f}') temp = df_ens['ens_mean'].loc[2090:2099].mean() print(f'ens_mean from 2090 to 2099 df {temp:.6f}') print(f'ens_mean from 2015 to 2060 {np.average(ens_mean[:46]):.3f}') if np.average(ens_mean[:46])<0: print('ensmean < 0',np.count_nonzero(ens_mean[:46]<0)) print('ens max < 0',np.count_nonzero(ens_max[:46]<0)) print('ensmean+std < 0',np.count_nonzero(ens_mean[:46]+ens_stddev[:46]<0)) else: # np.average(ens_mean[:46]>=0 print('ensmean > 0',np.count_nonzero(ens_mean[:46]>0)) print('ens min > 0',np.count_nonzero(ens_min[:46]>0)) print('ensmean-std > 0',np.count_nonzero(ens_mean[:46]-ens_stddev[:46]>0)) #def get_ens(expset,stream,yrstart=2090,yrend=2100,model='FOCI',exps='SynTra'): # foo = [] # for ireal in [0,1,2]: # if model=='FOCI': # foo.append(xr.open_mfdataset(get_file_list(expset,stream,yrstart=yrstart,yrend=yrend,real=[ireal],exps=exps))) # if model=='MPI-ESM': # foo.append(xr.open_mfdataset(get_file_list_mpiesm(expset,stream,yrstart=yrstart,yrend=yrend,real=[ireal]))) # return xr.concat(foo,dim='ireal') #def get_file_list(expset,stream,yrstart=2090,yrend=2100,real=[0],exps='SynTra'): # dirlist = read_dirlist(expset,exps=exps) # for ireal in real: ## print(f'{ireal}/2') # echam_all_years = [] # jsbach_all_years = [] # veg_all_years = [] # dirname = dirlist[ireal] # ptrc_T_all_years = [] # grid_T_all_years = [] # diad_T_all_years = [] # for year in range(yrstart,yrend): # # ptrc_T # foo = glob.glob(dirname + '../outdata/nemo/ym/*' + str(year)+'0101*_ptrc_T*.nc')[0] # ptrc_T_all_years.append(foo) # # grid_T # foo = glob.glob(dirname + '../outdata/nemo/ym/*' + str(year)+'0101*_grid_T*.nc')[0] # grid_T_all_years.append(foo) # # diad_T # foo = glob.glob(dirname + '../outdata/nemo/ym/*' + str(year)+'0101*_diad_T*.nc')[0] # diad_T_all_years.append(foo) # # echam # foo = glob.glob(dirname + 'echam6_echam_' + str(year) + '.nc')[0] # echam_all_years.append(foo) # # veg # foo = glob.glob(dirname + 'veg_' + str(year) + '.nc')[0] # veg_all_years.append(foo) # # jsbach # foo = glob.glob(dirname + 'jsbach_' + str(year) + '.nc')[0] # jsbach_all_years.append(foo) # if stream=='ptrc_T': # return ptrc_T_all_years # if stream=='diad_T': # return diad_T_all_years # if stream=='grid_T': # return grid_T_all_years # if stream=='echam': # return echam_all_years # if stream=='veg': # return veg_all_years # if stream=='jsbach': # return jsbach_all_years #def read_dirlist(expset,exps='SynTra'): # dirlist = [] # if exps=='SynTra': # if expset=='REF': # dirlist.append('/scratch/usr/shkhwwey/models/foci2.0/experiments/FOCI20.0-HW013_esm_ssp370_constLu/nc/') # dirlist.append('/scratch/usr/shkhwwey/models/foci2.0/experiments/FOCI20.0-HW012_esm_ssp370_constLu/nc/') # dirlist.append('/scratch/usr/shkhwwey/models/foci2.0/experiments/FOCI20.0-HW011_esm_ssp370_constLu/nc/') # if expset=='AFF': # dirlist.append('/scratch/usr/shkhwwey/models/foci2.0_maps_wh/experiments/FOCI20.0-HW022_esm_ssp370_AR/nc/') # dirlist.append('/scratch/usr/shkhwwey/models/foci2.0_maps_wh/experiments/FOCI20.0-HW023_esm_ssp370_AR/nc/') # dirlist.append('/scratch/usr/shkhwwey/models/foci2.0_maps_wh/experiments/FOCI20.0-HW024_esm_ssp370_AR/nc/') # if expset=='OAE': # dirlist.append('/scratch/usr/shkhwwey/models/foci2.0_maps_wh/experiments/FOCI20.0-HW028_esm_ssp370_OAE/nc/') # dirlist.append('/scratch/usr/shkhwwey/models/foci2.0_maps_wh/experiments/FOCI20.0-HW032_esm_ssp370_OAE/nc/') # dirlist.append('/scratch/usr/shkhwwey/models/foci2.0_maps_wh/experiments/FOCI20.0-HW033_esm_ssp370_OAE/nc/') # if expset=='BOTH': # dirlist.append('/scratch/usr/shkhwwey/models/foci2.0_maps_wh/experiments/FOCI20.0-HW029_esm_ssp370_AR_and_OAE/nc/') # dirlist.append('/scratch/usr/shkhwwey/models/foci2.0_maps_wh/experiments/FOCI20.0-HW030_esm_ssp370_AR_and_OAE/nc/') # dirlist.append('/scratch/usr/shkhwwey/models/foci2.0_maps_wh/experiments/FOCI20.0-HW031_esm_ssp370_AR_and_OAE/nc/') # return dirlist