Commit 538d7ff9 authored by Kennedy, Joseph H's avatar Kennedy, Joseph H
Browse files

Remove CESM2 table from energy extension

parent 2b0488b5
import os
# import numpy as np
import numpy.ma as ma
from netCDF4 import Dataset
def calculate_area_weighted_averages(cesm_path='/lustre/atlas1/cli115/world-shared/4ue/B_1PAL/'):
f_cesm_lnd_climo_jja = os.path.join(cesm_path, 'postproc/lnd/climos/B_1PAL_JJA_climo.nc')
f_cesm_lnd_climo_ann = os.path.join(cesm_path, 'postproc/lnd/climos/B_1PAL_ANN_climo.nc')
# f_cesm_lnd_climo_jja get FLDS (mean downwelling longwave rad),FIRA (mean net IR longwave rad, absorbed)
# FLDS(0,:,:), FIRA(0,:,:)
# Matrix gris_area(Greenland ice area) has values 0 and none
# Goal: print out area-weighted average of
# SW_d(268),
# SW_net(61),
# LW_d(235),
# LW_net(-46),
# R_net(15),
# SHF(7),
# LHF(-8),
# NOTE: these variables have size[192,288]
# FSDS: atmospheric incident solar radiation
# FSA: absorbed solar radiation
# FLDS: atmospheric longwave radiation
# FIRA: net infrared (longwave) radiation
# FIRE: emitted infrared (longwave) radiation
# FSH: sensible heat
# QSOIL: ground evaporation over ice landunits
# ---read input_files
ncid1 = Dataset(f_cesm_lnd_climo_jja, mode='r')
ncid3 = Dataset(f_cesm_lnd_climo_ann, mode='r')
fsds = ncid1.variables['FSDS'][0, :, :]
fsa = ncid1.variables['FSA'][0, :, :]
flds = ncid1.variables['FLDS'][0, :, :]
fira = ncid1.variables['FIRA'][0, :, :]
fsh = ncid1.variables['FSH'][0, :, :]
qsoil = ncid1.variables['QSOIL'][0, :, :]
rain = ncid3.variables['RAIN'][0, :, :]
snow = ncid3.variables['SNOW'][0, :, :]
qrunoff = ncid3.variables['QRUNOFF'][0, :, :]
# t2m = ncid11.variables['TREFHT'][0,:,:]
t2m = ncid1.variables['TSA'][0, :, :]
# CESM2
gris_mask = ncid1.variables['gris_mask'][:, :]
sw_d = fsds
sw_net = fsa
lw_d = flds
lw_net = -1.0 * fira
r_net = fsa - fira
shf = -1.0 * fsh
lhf = -1.0 * 2.835e6 * qsoil
# smb = QICE * 3.1567e7
smb = rain + snow - qrunoff - qsoil * 3.1567e7
t2m_c = t2m - 273.15
# ---grid cell area on Greenland. For other regions, the value is 9.999e35.
area = ncid1.variables['area'][:]
ncid1.close()
ncid3.close()
# Option I: use greenland mask to mask the area where the ice sheet fraction is zero
# Note: some non-greenland region has non-zero ice sheet fraction.
gris_mask = ma.masked_equal(gris_mask, 0)
area_maskice = ma.masked_array(area, mask=gris_mask.mask)
# Option II: use icesheet fraction mask to mask the area where the ice sheet fraction is zero
# ices_mask = ma.masked_equal(ices,0)
# area_maskice = ma.masked_array(area, mask=ices_mask.mask)
# Option III: first multiply area with ices, then use greenland mask to mask the area where the ice sheet
# fraction is zero
# NOTE: this option get the closest results to the Table 1 in the paper.
# areaices = np.multiply(area, ices)
# gris_mask = ma.masked_equal(gris_mask, 0)
# area_maskice = ma.masked_array(areaices, mask=gris_mask.mask)
aavgs = {
'Sw_d': (ma.average(sw_d, weights=area_maskice)),
'Sw_net': (ma.average(sw_net, weights=area_maskice)),
'Lw_d': (ma.average(lw_d, weights=area_maskice)),
'Lw_net': (ma.average(lw_net, weights=area_maskice)),
'Rnet': (ma.average(r_net, weights=area_maskice)),
'shf': (ma.average(shf, weights=area_maskice)),
'lhf': (ma.average(lhf, weights=area_maskice)),
'smb': (ma.average(smb, weights=area_maskice)),
't2c': (ma.average(t2m_c, weights=area_maskice)),
}
return aavgs
if __name__ == '__main__':
racmo_aavgs = calculate_area_weighted_averages()
for var, net in racmo_aavgs.items():
print('RACMO {}: {}'.format(var, net))
......@@ -51,7 +51,6 @@ with fn.temp_sys_path(os.path.dirname(__file__)):
import energy.cesm_racmo23_t2m_ann as t2m_ann
import energy.cesm_racmo23_t2m_djf as t2m_djf
import energy.cesm_racmo23_t2m_jja as t2m_jja
import energy.rad_aavg_cesm2 as cesm2_rad
import energy.rad_aavg_model as cesm_rad
import energy.rad_aavg_racmo as racmo_rad
import energy.timeseries_trefht as ts_trefht
......@@ -77,15 +76,12 @@ def run(name, config):
element_list = []
cesm_aavg = cesm_rad.calculate_area_weighted_averages()
cesm2_aavg = cesm2_rad.calculate_area_weighted_averages()
racmo_aavg = racmo_rad.calculate_area_weighted_averages()
cesm_aavg_title = """Area weighted averages of CESM's energy balance variables."""
cesm2_aavg_title = """Area weighted averages of CESM2's energy balance variables."""
racmo_aavg_title = """Area weighted averages of RACMO's energy balance variables."""
element_list.append(el.vtable(cesm_aavg_title, cesm_aavg.keys(), cesm_aavg))
element_list.append(el.vtable(cesm2_aavg_title, cesm2_aavg.keys(), cesm2_aavg))
element_list.append(el.vtable(racmo_aavg_title, racmo_aavg.keys(), racmo_aavg))
# PLOTS
......
Markdown is supported
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment