add postproc files to lex repo to prepare model files for analysis

postproc/README_scripts.md

+## Postprocessing of files of coupled or ice sheet only runs for use in LIVVkit
+
+This directory holds subdirectories with scripts that will postprocess model data for use in creating the plots and data used by LIVVkit. There are also scripts to processed observation data useful for any model analysis.
+
+These script environments take raw output from either a stand alone ice sheet model (e.g. cism-albany)
+or a coupled climate model (e.g. cesm1p2, e3sm (in beta testing))
+
+To enable or update the scripts to work on a given machine, alter the main scripts in the target model 
+subdirectory to load the necessary software.
+
+To create scripts for a model not yet included here, create a new subidirectory with README files and all 
+the necessary files to complete the postprocessing
+
+
+#### Some notes about nco. 
+The ncra and ncwa commands are mostly used for averaging, but in some cases they compute a sum. 
+Here are some alternative commands with ncra and ncwa using the -y op_typ flag that can be done within LIVVkit:
+* avg returns time average
+* min time min
+* max time max
+* ttl temporal sum
+* sdn temporal std deviation
+(see 3.35 in nco documentation for complete descrtion of op_types
+
+These features could also be used with the nces (ens evg feature)

postproc/cesm1p2/README_cesm1p2.md

+## Postprocessing of files of coupled CESM version CESM version 1.2 runs for use in LIVVkit
+(note: a test version for E3SM data 'main_e3sm.sh' is also here, but its not production yet)
+
+This directory holds subdirectories with scripts that will postprocess model data for use in creating the plots and data used by LIVVkit. 
+
+### What is needed to run postproc scripts:
+	*Create a working directory (with write permissions) to put all target files
+	*Refer to check list below to ensure all the files are accessable 
+	*Run on OLCFs rhea system, to ensure enough memory for postprocessing CISM-Albany data
+
+### Files needed to run coupled model postprocessing (right now a copy exists in the $WORLDWORK/cli115//4ue space)
+
+1. b.e10.BG20TRCN.f09_g16.002
+	* /lnd/hist
+	* /atm/hist
+	* /glc/hist
+2. surfdata_0.9x1.25_simyr1850_c110725.nc
+3. /obs_data
+	* /CLOUDSAT_*	for all months, DJF, JJA, and ANN
+	* /ISCCP_*    for all months, DJF, JJA, and ANN	
+4. /racmo23_GRN_monthly (all this data is already separated by variable )
+5. Greenland_5km_v1.1_SacksRev_c110629.nc 
+
+### Main script for analysis 
+
+#### Usage:
+	*Copy over main.sh to separate file (e.g. main_deal.sh) 
+	*Set parameters at the top of main_deal.sh to match your environment
+	*Below that, opt out, if desired, to postprocessing, remapping, and plots/tables
+
+```sh
+ bash main_deal.sh 
+```
+
+#### Postprocessing
+* loads computing environment for rhea, with hooks to add other machines
+* makes necessary mask files if they do not already exist
+* if 'postprocess_obs' is selected
+	* calls atm_obs.sh
+		* remaps ISCCP and CLOUDSAT climo data to model grid for selected variables for plots
+		* makes ice sheet area averaged climos from original climo data
+	* calls lnd_obs.sh 
+		* makes climos from original racmo climo data
+		* makes yearly time series from original racmo monthly files 
+		* makes ice sheet area average (but not ice sheet %) of yearly time series of racmo data
+* calls atm.sh
+	* makes climos and ice sheet area averages
+	* calls monthly_atm_ts.sh 
+		* make monthly and yearly averaged time series 
+		* make ice sheet area average of time series 
+* calls lnd.sh
+	* makes climos and ice area average of climos
+ 	* makes daily avergages over time record 
+	* submits rhea_tseries_lnd.pbs  
+		* makes time series of daily avg data
+		* makes ice sheet average of daily time series	
+	* makes annualized time series of daily avg data
+	* makes ice sheet average of annualized daily time series	
+	* calls tseries_lnd_h0.sh
+		* makes time series from monthly averaged data
+		* makes yearly averages of monthly averaged data
+		* makes ice sheet average of monthly averaged data 
+* calls glc.sh 
+	* select flag to compute velnorm (turn off if velocities do not exist)
+	* makes annual climos
+	* makes ice sheet average of climos
+	* makes time series of yearly averaged data
+	* area average time series (it is not masked right now)
+
+
+#### Make plots and tables
+* contour plots of 2m T MODEL, RACMO and difference for DJF and JJA
+* plots of ANN low, med, and high clouds for MODEL, ISCCP, and CLOUDSAT
+* plots of the annual cycle of low, med, and high clouds for MODEL, ISCCP, and CLOUDSAT
+* plot of time series of 2m T MODEL and RACMO
+* print ice sheet are and time averages of suite of variables for model and RACMO
+* time series of downscaled CISM SMB

postproc/cesm1p2/atm.sh

+#!/bin/bash
+
+# atmosphere
+
+if [ "$postprocess" = "true" ]; then
+
+# regridding step to regrid to lat lon, currently regrids all years even if not requested
+if [ "$version" == "e3sm" ]; then
+
+##ncwa -a -time $outpath/masks/mask_ne30np4.nc tmp.nc 
+##ncks -x -v time tmp.nc $outpath/masks/mask_ne30np4_notime.nc
+
+#for fl in $path/atm/hist/*h0*; do
+#echo "first add mask in original grid to monthly atm output in $fl"
+##ncks -A -C -v gris_mask $outpath/masks/mask_ne30np4_notime.nc $fl   
+#ncks -A -C -v gris_mask $mask_file $fl   
+#done
+
+h0atm_path=$outpath/atm/hist_regrid
+
+#mkdir -p "$h0atm_path"
+#ls $path/atm/hist/*h0* | ncremap -m $map_file -O $h0atm_path
+
+else
+h0atm_path=$path/atm/hist
+
+fi
+
+#run a check, is the requested data directory there?
+if test -d $h0atm_path; then
+        echo "Directory $h0atm_path exists"
+else
+        echo "Target directory $h0atm_path does not exist, you prob need to activate remapping"
+        exit
+fi
+
+echo "make atm climos"
+mkdir -p $outpath/atm/climos
+ncclimo -a sdd -c ${case} -s ${year_start} -e ${year_end} -m ${atm} -i ${h0atm_path} -o ${outpath}/atm/climos
+
+if test -a $outpath/atm/climos/${case}_ANN_climo.nc; then
+      	echo "File $outpath/atm/climos/${case}_ANN_climo.nc exists"
+else
+      	echo "File $outpath/atm/climos/${case}_ANN_climo.nc does not exist"
+	exit
+fi
+
+if [ "$version" == "cesm1" ]; then
+echo "add mask and area variables to monthly atm climos"
+ncks -A -v gris_mask,area $outpath/masks/$case.$lnd.gris_mask.nc $outpath/atm/climos/${case}_ANN_climo.nc
+ncks -A -v gris_mask,area $outpath/masks/$case.$lnd.gris_mask.nc $outpath/atm/climos/${case}_DJF_climo.nc
+ncks -A -v gris_mask,area $outpath/masks/$case.$lnd.gris_mask.nc $outpath/atm/climos/${case}_JJA_climo.nc
+ncks -A -v gris_mask,area $outpath/masks/$case.$lnd.gris_mask.nc $outpath/atm/climos/${case}_01_climo.nc
+ncks -A -v gris_mask,area $outpath/masks/$case.$lnd.gris_mask.nc $outpath/atm/climos/${case}_02_climo.nc
+ncks -A -v gris_mask,area $outpath/masks/$case.$lnd.gris_mask.nc $outpath/atm/climos/${case}_03_climo.nc
+ncks -A -v gris_mask,area $outpath/masks/$case.$lnd.gris_mask.nc $outpath/atm/climos/${case}_04_climo.nc
+ncks -A -v gris_mask,area $outpath/masks/$case.$lnd.gris_mask.nc $outpath/atm/climos/${case}_05_climo.nc
+ncks -A -v gris_mask,area $outpath/masks/$case.$lnd.gris_mask.nc $outpath/atm/climos/${case}_06_climo.nc
+ncks -A -v gris_mask,area $outpath/masks/$case.$lnd.gris_mask.nc $outpath/atm/climos/${case}_07_climo.nc
+ncks -A -v gris_mask,area $outpath/masks/$case.$lnd.gris_mask.nc $outpath/atm/climos/${case}_08_climo.nc
+ncks -A -v gris_mask,area $outpath/masks/$case.$lnd.gris_mask.nc $outpath/atm/climos/${case}_09_climo.nc
+ncks -A -v gris_mask,area $outpath/masks/$case.$lnd.gris_mask.nc $outpath/atm/climos/${case}_10_climo.nc
+ncks -A -v gris_mask,area $outpath/masks/$case.$lnd.gris_mask.nc $outpath/atm/climos/${case}_11_climo.nc
+ncks -A -v gris_mask,area $outpath/masks/$case.$lnd.gris_mask.nc $outpath/atm/climos/${case}_12_climo.nc
+fi
+
+echo "average over ice sheet region of atm climos "
+ncwa -O -a lat,lon -w area -B "gris_mask == 1" $outpath/atm/climos/${case}_ANN_climo.nc $outpath/atm/climos/${case}_ANN_aavg_climo.nc 
+ncwa -O -a lat,lon -w area -B "gris_mask == 1" $outpath/atm/climos/${case}_JJA_climo.nc $outpath/atm/climos/${case}_JJA_aavg_climo.nc 
+ncwa -O -a lat,lon -w area -B "gris_mask == 1" $outpath/atm/climos/${case}_DJF_climo.nc $outpath/atm/climos/${case}_DJF_aavg_climo.nc 
+ncwa -O -a lat,lon -w area -B "gris_mask == 1" $outpath/atm/climos/${case}_01_climo.nc $outpath/atm/climos/${case}_01_aavg_climo.nc 
+ncwa -O -a lat,lon -w area -B "gris_mask == 1" $outpath/atm/climos/${case}_02_climo.nc $outpath/atm/climos/${case}_02_aavg_climo.nc 
+ncwa -O -a lat,lon -w area -B "gris_mask == 1" $outpath/atm/climos/${case}_03_climo.nc $outpath/atm/climos/${case}_03_aavg_climo.nc 
+ncwa -O -a lat,lon -w area -B "gris_mask == 1" $outpath/atm/climos/${case}_04_climo.nc $outpath/atm/climos/${case}_04_aavg_climo.nc 
+ncwa -O -a lat,lon -w area -B "gris_mask == 1" $outpath/atm/climos/${case}_05_climo.nc $outpath/atm/climos/${case}_05_aavg_climo.nc 
+ncwa -O -a lat,lon -w area -B "gris_mask == 1" $outpath/atm/climos/${case}_06_climo.nc $outpath/atm/climos/${case}_06_aavg_climo.nc 
+ncwa -O -a lat,lon -w area -B "gris_mask == 1" $outpath/atm/climos/${case}_07_climo.nc $outpath/atm/climos/${case}_07_aavg_climo.nc 
+ncwa -O -a lat,lon -w area -B "gris_mask == 1" $outpath/atm/climos/${case}_08_climo.nc $outpath/atm/climos/${case}_08_aavg_climo.nc 
+ncwa -O -a lat,lon -w area -B "gris_mask == 1" $outpath/atm/climos/${case}_09_climo.nc $outpath/atm/climos/${case}_09_aavg_climo.nc 
+ncwa -O -a lat,lon -w area -B "gris_mask == 1" $outpath/atm/climos/${case}_10_climo.nc $outpath/atm/climos/${case}_10_aavg_climo.nc 
+ncwa -O -a lat,lon -w area -B "gris_mask == 1" $outpath/atm/climos/${case}_11_climo.nc $outpath/atm/climos/${case}_11_aavg_climo.nc 
+ncwa -O -a lat,lon -w area -B "gris_mask == 1" $outpath/atm/climos/${case}_12_climo.nc $outpath/atm/climos/${case}_12_aavg_climo.nc 
+
+echo "make yearly time series for atm h0 data"
+mkdir -p $outpath/atm/tseries
+
+for ivar in "${atm_vars[@]}"
+do
+
+. tseries_atm_h0.sh 
+
+done
+
+        echo "ATM POSTPROCESSING COMPLETE "
+        echo " "
+
+fi # to postprocess or not
+

postproc/cesm1p2/cesm2_rad_aavg.py

+##On my laptop I installed PyNGL and PyNIO to use Nio and Ngl module for reading *.nc file and for contour plot
+##On titan, I did not install PyNGL and PyNIO, so import Dataset from netCDF4 to read *.nc file
+##probably, I will the following command for contourplot in titan python "from mpl_toolkits.basemap import Basemap"
+
+
+#!/usr/bin/env python3
+
+import numpy as np
+import numpy.ma as ma
+import os
+from netCDF4 import Dataset
+
+
+# ---------------- Data sources ------------------------
+#CESM1
+#input_file1 = '/lustre/atlas1/cli115/world-shared/4ue/b.e10.BG20TRCN.f09_g16.002/postproc/lnd/climos/b.e10.BG20TRCN.f09_g16.002_JJA_climo.nc'
+#input_file11 = '/lustre/atlas1/cli115/world-shared/4ue/b.e10.BG20TRCN.f09_g16.002/postproc/atm/climos/b.e10.BG20TRCN.f09_g16.002_ANN_climo.nc'
+#input_file3 = '/lustre/atlas1/cli115/world-shared/4ue/b.e10.BG20TRCN.f09_g16.002/postproc/lnd/climos/b.e10.BG20TRCN.f09_g16.002_ANN_climo.nc'
+#CESM2
+input_file1 = '/lustre/atlas1/cli115/world-shared/4ue/B_1PAL/postproc/lnd/climos/B_1PAL_JJA_climo.nc'
+#input_file11 = '/lustre/atlas1/cli115/world-shared/4ue/B_1PAL/postproc/atm/climos/B_1PAL_ANN_climo.nc'
+input_file3 = '/lustre/atlas1/cli115/world-shared/4ue/B_1PAL/postproc/lnd/climos/B_1PAL_ANN_climo.nc'
+#for CESM1 and CESM2
+#input_file2 = '/lustre/atlas1/cli115/world-shared/4ue/surfdata_0.9x1.25_simyr1850_c110725.nc'
+
+# --------------------------------------------------------------
+# --------------------------------------------------------------
+
+# input_file1 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),
+# 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_file2, probably get area with ice sheet with nonzero values
+#ncid2 = Dataset(input_file2, mode='r')
+#ices = ncid2.variables['PCT_GLC_ICESHEET'][:] 
+#ncid2.close()
+
+# ---read input_files
+ncid1 = Dataset(input_file1, mode='r')
+#ncid11 = Dataset(input_file11, mode='r')
+ncid3 = Dataset(input_file3, 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,:,:] 
+
+#CESM1
+GrisMask = ncid1.variables['gris_mask'][0,:,:] 
+
+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
+T2mC = T2m - 273.15
+
+# ---grid cell area on Greenland. For other regions, the value is 9.999e35.
+area = ncid1.variables['area'][:] 
+
+ncid1.close()
+#ncid11.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(GrisMask,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
+## -- this option gives the closest results to the Table 1 in the paper.
+#areaices = np.multiply(area,ices)
+#Gris_mask = ma.masked_equal(GrisMask,0)
+#area_maskice = ma.masked_array(areaices, mask=Gris_mask.mask)
+
+
+swdout=(ma.average(SW_d,weights=area_maskice))
+swnetout=(ma.average(SW_net,weights=area_maskice))
+lwdout=(ma.average(LW_d,weights=area_maskice))
+lwnetout=(ma.average(LW_net,weights=area_maskice))
+rnetout=(ma.average(R_net,weights=area_maskice))
+shfout=(ma.average(SHF,weights=area_maskice))
+lhfout=(ma.average(LHF,weights=area_maskice))
+smbout=(ma.average(SMB,weights=area_maskice))
+t2cout=(ma.average(T2mC,weights=area_maskice))
+
+print "Sw_d model:", swdout
+print "Sw_net model:", swnetout
+print "Lw_d model:", lwdout
+print "Lw_net model:", lwnetout
+print "Rnet model:", rnetout
+print "SHF model:", shfout
+print "LHF model:", lhfout
+print "SMB model:", smbout
+#print "T2M C atm model:", t2cout
+print "T2M C land model:", t2cout
+
+

postproc/cesm1p2/e3sm_rad_aavg.py

+##On my laptop I installed PyNGL and PyNIO to use Nio and Ngl module for reading *.nc file and for contour plot
+##On titan, I did not install PyNGL and PyNIO, so import Dataset from netCDF4 to read *.nc file
+##probably, I will the following command for contourplot in titan python "from mpl_toolkits.basemap import Basemap"
+
+
+#!/usr/bin/env python3
+
+import numpy as np
+import numpy.ma as ma
+import os
+from netCDF4 import Dataset
+
+
+# ---------------- Data sources ------------------------
+#CESM1
+#input_file1 = '/lustre/atlas1/cli115/world-shared/4ue/b.e10.BG20TRCN.f09_g16.002/postproc/lnd/climos/b.e10.BG20TRCN.f09_g16.002_JJA_climo.nc'
+#input_file11 = '/lustre/atlas1/cli115/world-shared/4ue/b.e10.BG20TRCN.f09_g16.002/postproc/atm/climos/b.e10.BG20TRCN.f09_g16.002_ANN_climo.nc'
+#input_file3 = '/lustre/atlas1/cli115/world-shared/4ue/b.e10.BG20TRCN.f09_g16.002/postproc/lnd/climos/b.e10.BG20TRCN.f09_g16.002_ANN_climo.nc'
+#CESM2
+input_file1 = '/lustre/atlas1/cli115/world-shared/4ue/20171122.beta3rc10_1850.ne30_oECv3_ICG.edison/postproc/lnd/climos/20171122.beta3rc10_1850.ne30_oECv3_ICG.edison_JJA_climo.nc'
+#input_file11 = '/lustre/atlas1/cli115/world-shared/4ue/20171122.beta3rc10_1850.ne30_oECv3_ICG.edison/postproc/atm/climos/20171122.beta3rc10_1850.ne30_oECv3_ICG.edison_ANN_climo.nc'
+input_file3 = '/lustre/atlas1/cli115/world-shared/4ue/20171122.beta3rc10_1850.ne30_oECv3_ICG.edison/postproc/lnd/climos/20171122.beta3rc10_1850.ne30_oECv3_ICG.edison_ANN_climo.nc'
+#for CESM1 and CESM2
+#input_file2 = '/lustre/atlas1/cli115/world-shared/4ue/surfdata_0.9x1.25_simyr1850_c110725.nc'
+
+# --------------------------------------------------------------
+# --------------------------------------------------------------
+
+# input_file1 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),
+# 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_file2, probably get area with ice sheet with nonzero values
+#ncid2 = Dataset(input_file2, mode='r')
+#ices = ncid2.variables['PCT_GLC_ICESHEET'][:] 
+#ncid2.close()
+
+# ---read input_files
+ncid1 = Dataset(input_file1, mode='r')
+#ncid11 = Dataset(input_file11, mode='r')
+ncid3 = Dataset(input_file3, 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,:,:] 
+
+#E3SM
+GrisMask = 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
+T2mC = T2m - 273.15
+
+# ---grid cell area on Greenland. For other regions, the value is 9.999e35.
+area = ncid1.variables['area'][:] 
+
+ncid1.close()
+#ncid11.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(GrisMask,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
+## -- this option gives the closest results to the Table 1 in the paper.
+#areaices = np.multiply(area,ices)
+#Gris_mask = ma.masked_equal(GrisMask,0)
+#area_maskice = ma.masked_array(areaices, mask=Gris_mask.mask)
+
+
+swdout=(ma.average(SW_d,weights=area_maskice))
+swnetout=(ma.average(SW_net,weights=area_maskice))
+lwdout=(ma.average(LW_d,weights=area_maskice))
+lwnetout=(ma.average(LW_net,weights=area_maskice))
+rnetout=(ma.average(R_net,weights=area_maskice))
+shfout=(ma.average(SHF,weights=area_maskice))
+lhfout=(ma.average(LHF,weights=area_maskice))
+smbout=(ma.average(SMB,weights=area_maskice))
+t2cout=(ma.average(T2mC,weights=area_maskice))
+
+print "Sw_d model:", swdout
+print "Sw_net model:", swnetout
+print "Lw_d model:", lwdout
+print "Lw_net model:", lwnetout
+print "Rnet model:", rnetout
+print "SHF model:", shfout
+print "LHF model:", lhfout
+print "SMB model:", smbout
+#print "T2M C atm model:", t2cout
+print "T2M C land model:", t2cout
+
+

postproc/cesm1p2/glc.sh

+#!/bin/bash
+# ice sheet (CESM1 and CESM2 have the same naming conventions)
+
+if [ "$postprocess" = "true" ]; then
+calcvn=1 # create velnorm if uvel and vvel exist
+
+#run a check, is the requested data directory there?
+if test -d $path/glc/hist; then
+	echo "Directory $path/glc/hist exists"
+else
+	echo "Target directory $path/glc/hist does not exist"
+	exit 
+fi
+
+if [ $calcvn = 1 ]; then
+
+mkdir -p $outpath/glc/hist
+for (( year=${year_start}; year<=${year_end}; year++ ))
+do
+#echo "make velnorm from just uvel and vvel (memory constraints require separating uvel and vvel first)"
+#ncrcat -O -v time,level,uvel,vvel $path/glc/hist/${case}.$glc.h.${year}-01-01-00000.nc $outpath/glc/hist/${case}.$glc.h.${year}-01-01-00000.vels.nc
+#echo "created" $outpath/${case}.$glc.h.${year}-01-01-00000.vels.nc
+
+echo "creating $outpath/${case}.$glc.h.velnorm.${year}-01-01-00000.nc file with velnorm"
+ncap -O -s "velnorm=sqrt(uvel^2+vvel^2)" $path/glc/hist/${case}.$glc.h.${year}-01-01-00000.nc $outpath/glc/hist/${case}.$glc.h.velnorm.${year}-01-01-00000.nc
+done
+
+fi
+
+echo "make ice sheet climos"
+mkdir -p $outpath/glc/climos
+echo /ccs/home/zender/bin_rhea/ncclimo -C ann -c $case -s $year_start -e $year_end -m $glc -h h -i $path/glc/hist -o $outpath/glc/climos
+/ccs/home/zender/bin_rhea/ncclimo -C ann -c $case -s $year_start -e $year_end -m $glc -h h -i $path/glc/hist -o $outpath/glc/climos
+if test -a $outpath/glc/climos/${case}_ANN_${year_start}01_${year_end}12_climo.nc; then
+       	echo "File $outpath/glc/climos/${case}_ANN_${year_start}01_${year_end}12_climo.nc exists"
+else
+       	echo "File $outpath/glc/climos/${case}_ANN_${year_start}01_${year_end}12_climo.nc does not exist"
+	exit
+fi
+
+echo "add area variable to file for averaging/summing"
+ncks -A -v AREA $area_file $outpath/glc/climos/${case}_ANN_${year_start}01_${year_end}12_climo.nc
+
+echo "add masking variable to file for averaging/summing"
+ncks -A -v landcover $base/Greenland_5km_v1.1_SacksRev_c110629.nc $outpath/glc/climos/${case}_ANN_${year_start}01_${year_end}12_climo.nc
+
+echo "average annual ice sheet climos over ice sheet"
+ncwa -O -a x1,y1,lev,time -w AREA -B "landcover >= 1" $outpath/glc/climos/${case}_ANN_${year_start}01_${year_end}12_climo.nc $outpath/glc/climos/${case}_ANN_aavg.nc
+#
+if test -a $outpath/glc/climos/${case}_ANN_aavg.nc; then
+        echo "Global annual averages of $glc variables made"
+else
+        echo "Global annual averages of $glc variables failed"
+        exit
+fi
+
+
+echo "make time_series of annual averages"
+mkdir -p $outpath/glc/tseries
+
+for i in "${glc_vars[@]}"
+do
+
+echo "make time series of yearly averages"
+ncrcat -O -v $i $path/glc/hist/$case.$glc.h.*.nc $outpath/glc/tseries/$case.$glc.h.$i.nc
+#ncrcat -O -v time,$i $path/glc/hist/$case.$glc.h.*.nc $outpath/glc/tseries/$case.$glc.h.$i.nc
+#ncea -F -d time,first,last in.out out.nc
+echo "add area variable to file for averaging/summing"
+ncks -v AREA -A $area_file $outpath/glc/tseries/${case}.$glc.h.$i.nc
+# this line screws up time dimension
+#ncks -v landcover -A /lustre/atlas1/cli106/proj-shared/4ue/CESM-CISM/Greenland_5km_v1.1_SacksRev_c110629.nc $outpath/glc/tseries/${case}.$glc.h.$i.nc
+
+echo "average variable over whole ice sheet domain"
+ncwa -O -a x1,y1 -w AREA $outpath/glc/tseries/${case}.$glc.h.$i.nc $outpath/glc/tseries/${case}.$glc.h.${i}_aavg.nc
+#ncwa -O -a x1,y1 -w AREA -B "landcover >= 1" $outpath/glc/tseries/${case}.$glc.h.$i.nc $outpath/glc/tseries/${case}.$glc.h.${i}_aavg.nc
+
+
+if test -a $outpath/glc/tseries/$case.$glc.h.${i}_aavg.nc; then
+        echo "Timeseries of global annual averages of $glc for $i made"
+else
+        echo "Timeseries of global annual averages of $glc for $i failed"
+        exit
+fi
+
+
+done
+
+fi # to postprocess or not
+
+#if [ "$make_plots" = "true" ]; then
+
+#echo "make the plots"
+#mkdir -p "$outpath/glc/plots"
+
+#pwd
+#python plot_files/TimeSeries_acab.py
+
+#fi # to make plots or not

postproc/cesm1p2/lnd.sh

+#!/bin/bash
+
+# land
+
+if [ "$postprocess" = "true" ]; then
+
+if [ "$version" == "e3sm" ]; then
+#ncwa -a -time $outpath/masks/mask_ne30np4_lnd.nc $outpath/masks/mask_ne30np4_lnd_notime.nc 
+#for fl in $path/lnd/hist/*h0*; do
+#echo "first add mask in original grid to monthly lnd output in $fl"
+#ncks -A -C -v gris_mask $outpath/masks/mask_ne30np4_lnd_notime.nc $fl
+#done
+
+# regridding step to regrid to lat lon, currently regrids all years even if not requested
+#h0lnd_path=$outpath/lnd/hist_regrid
+#mkdir -p "$h0lnd_path"
+#ls $path/lnd/hist/*h0* | ncremap -m $map_file -O $h0lnd_path
+h0lnd_path=$path/archive/lnd/hist
+
+else 
+h0lnd_path=$path/lnd/hist
+fi
+
+echo "make land h0 climos"
+mkdir -p "$outpath/lnd/climos"
+
+ncclimo -c $case -s ${year_start} -e ${year_end} -m ${lnd} -i ${h0lnd_path} -o ${outpath}/lnd/climos
+if test -a $outpath/lnd/climos/${case}_ANN_climo.nc; then
+       	echo "File $outpath/lnd/climos/${case}_ANN_climo.nc exists"
+else
+       	echo "File $outpath/lnd/climos/${case}_ANN_climo.nc does not exist"
+	exit
+fi
+
+if [ "$version" == "cesm1" ]; then
+
+h0lnd_path=$outpath/lnd/climos
+
+echo "average land monthly climos over ice sheet region"
+ncwa -O -a lat,lon -w area -B "gris_mask == 1" $outpath/lnd/climos/${case}_ANN_climo.nc $outpath/lnd/climos/${case}_ANN_aavg_climo.nc
+ncwa -O -a lat,lon -w area -B "gris_mask == 1" $outpath/lnd/climos/${case}_JJA_climo.nc $outpath/lnd/climos/${case}_JJA_aavg_climo.nc
+ncwa -O -a lat,lon -w area -B "gris_mask == 1" $outpath/lnd/climos/${case}_DJF_climo.nc $outpath/lnd/climos/${case}_DJF_aavg_climo.nc
+
+mkdir -p "$outpath/lnd/davg"
+
+#echo "make annualized daily averages for land for all variables"
+#ls $path/lnd/hist/$case.$lnd.h1.*.nc |  /ccs/home/zender/bin_rhea/ncclimo --job_nbr=12 -c $case -C dly -s ${year_start} -e ${year_end} -o $outpath/lnd/davg -m $lnd 
+#/ccs/home/zender/bin_rhea/ncclimo --job_nbr=12 -c $case -C dly -s ${year_start} -e ${year_end} -o $outpath/lnd/davg -m $lnd -i $path/lnd/hist
+
+
+#echo "loop over each variable to create and process time series for h1 vars"
+#qsub rhea_tseries_lnd.sh 
+
+#for ivar in "${lnd_h1_vars[@]}"