Commit 92562ed5 authored by Carl Schreck's avatar Carl Schreck

Mapping prism differences

parent ee938b7a
......@@ -15,7 +15,7 @@ import cjs # '~carl/lib/python'
cjs.tstamp('Here we go!')
# These are some parameters that could be useful to have up top
parser = argparse.ArgumentParser()
parser.add_argument('--var_name', default='tavg')
parser.add_argument('--var_name', default='pgt0')
parser.add_argument('--month', default=9)
args = parser.parse_args()
var_name = args.var_name
......@@ -29,9 +29,12 @@ ncg_var_to_prism_var = {
'tmin': 'tmin',
'tmax': 'tmax',
'tavg': 'tmean',
'precp': 'ppt',
'prcp': 'ppt',
}
ncg_var = var_name
if var_name == 'pgt0': # precip greater than zero
ncg_var = 'prcp'
else:
ncg_var = var_name
prism_var = ncg_var_to_prism_var[ncg_var]
for year in years:
......@@ -44,6 +47,9 @@ for year in years:
sub_prism = prism_ds.sel(lat=ncg_ds.lat, lon=ncg_ds.lon, method='nearest')
sub_prism['lat'] = ncg_ds.lat
sub_prism['lon'] = ncg_ds.lon
if var_name == 'pgt0': # precip greater than zero
ncg_ds[ncg_var] = xr.where(ncg_ds[ncg_var] == 0, 0, 1)
sub_prism[prism_var] = xr.where(sub_prism[prism_var] == 0, 0, 1)
if year == years[0]:
cjs.tstamp(year)
bias = ((ncg_ds[ncg_var] - sub_prism[prism_var]).mean(dim='time')
......@@ -59,12 +65,16 @@ for year in years:
cjs.tstamp('Writing')
bias.attrs['long_name'] = 'Bias ncg - prism'
bias.attrs['units'] = ncg_ds[ncg_var].attrs['units']
absdiff.attrs['long_name'] = 'Absolute difference ncg - prism'
absdiff.attrs['units'] = ncg_ds[ncg_var].attrs['units']
if var_name == 'pgt0':
bias.attrs['units'] = 'frequency'
absdiff.attrs['units'] = 'frequency'
else:
bias.attrs['units'] = ncg_ds[ncg_var].attrs['units']
absdiff.attrs['units'] = ncg_ds[ncg_var].attrs['units']
ds = xr.Dataset(data_vars={'bias': bias, 'absdiff': absdiff})
out_path = f'{os.environ["DATA_DIR"]}/nclimgrid/beta/compare/' \
f'{ncg_var}-prism-{month:02}.nc'
f'{var_name}-prism-{month:02}.nc'
ds.to_netcdf(out_path)
cjs.tstamp('Thank you, come again.')
......
""" Draw a map of comparison between nclimgrid and prism."""
__author__ = "Carl Schreck"
__email__ = "cjschrec@ncsu.edu"
__copyright__ = "Copyright 2020, North Carolina State University"
__license__ = "BSD-3.0"
import xarray as xr
import numpy as np
import os
import glob
from matplotlib import pyplot as plt
import cartopy.crs as ccrs
import cjs # '~carl/lib/python'
cjs.tstamp('Here we go!')
# These are some parameters that could be useful to have up top
var_name = 'tavg'
in_path = f'{os.environ["DATA_DIR"]}/nclimgrid/beta/compare/{var_name}*.nc'
var_name_to_contour_level = {
'tmin': 0.5,
'tmax': 0.5,
'tavg': 0.5,
'prcp': 0.5,
'pgt0': 0.5,
}
ds = xr.open_mfdataset(glob.glob(in_path), concat_dim='month', combine='nested')
ds = ds.mean(dim='month', keep_attrs=True)
ax_bias = plt.subplot(211,
projection=ccrs.AlbersEqualArea(
central_longitude=360-95,
standard_parallels=(40,60)))
plt.gca().outline_patch.set_visible(False)
ax_bias.set_extent([-118, -75, 23.5, 50.5])
shp_name = 'ne_10m_admin_1_states_provinces_lakes'
cjs.add_map(ax_bias, shp_name, edgecolor='none', facecolor='silver', zorder=1,
usa_only=True)
cjs.add_map(ax_bias, shp_name, edgecolor=(0.2, 0.2, 0.2), facecolor='none',
zorder=3, usa_only=True)
levels = list(np.arange(-5,6)*0.5)
ds['bias'].plot(ax=ax_bias, transform=ccrs.PlateCarree(),
cmap='RdBu_r', levels=levels, zorder=2, extend='both',
vmin=-np.inf, vmax=np.inf,
)
ax_bias.set_title('a) Bias', loc='left')
ax_bias.set_title(f'mean: {ds["bias"].mean().values:0.2f}', loc='right')
cjs.tstamp(f'mean bias: {ds["bias"].mean().values:0.2f}')
ax_absdiff = plt.subplot(212,
projection=ccrs.AlbersEqualArea(
central_longitude=360-95,
standard_parallels=(40,60)))
plt.gca().outline_patch.set_visible(False)
ax_absdiff.set_extent([-118, -75, 23.5, 50.5])
shp_name = 'ne_10m_admin_1_states_provinces_lakes'
cjs.add_map(ax_absdiff, shp_name, edgecolor='none', facecolor='silver',
zorder=1, usa_only=True)
cjs.add_map(ax_absdiff, shp_name, edgecolor=(0.2, 0.2, 0.2), facecolor='none',
zorder=3, usa_only=True)
levels = list(np.arange(1,8)*0.5)
ds['absdiff'].plot(ax=ax_absdiff, transform=ccrs.PlateCarree(),
levels=levels, zorder=2, extend='both',
vmin=-np.inf, vmax=np.inf,
)
ax_absdiff.set_title('b) Abs. diff.', loc='left')
ax_absdiff.set_title(f'mean: {ds["absdiff"].mean().values:0.2f}', loc='right')
cjs.tstamp(f'mean absdiff: {ds["absdiff"].mean().values:0.2f}')
plt.suptitle(var_name)
plt.savefig(f'figures/{var_name}.png', dpi=200)
plt.show()
cjs.tstamp('Thank you, come again.')
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