#!/usr/bin/env python # Wenchang Yang (wenchang@princeton.edu) # Thu Nov 28 10:06:59 PM EST 2024 if __name__ == '__main__': import sys,os try: from misc.timer import Timer tt = Timer(f'[{os.getcwd()}] start ' + ' '.join(sys.argv)) except: pass import sys, os.path, os, glob, datetime import xarray as xr, numpy as np, pandas as pd, matplotlib.pyplot as plt #more imports wython = '/tigress/wenchang/wython' if wython not in sys.path: sys.path.append(wython); print('added to python path:', wython) #from misc import get_kws_from_argv from misc.seasons import sel_season from xtc import wherebasin from misc.landmask import whereocean # if __name__ == '__main__': try: tt.check('end import') except: pass # #start from here ifile = 'Vshear_HIRAM_amipLMR2019SST0850ic_tigercpu_intelmpi_18_540PE_ens01_0850-1999.nc' da = xr.open_dataarray(ifile) timeslice_early = slice('1351', '1450') timeslice_late = slice('1451', '1550') season = 'ASO' basin = 'NA' with xr.set_options(keep_attrs=True): da_ref = da.sel(time=timeslice_early).pipe(sel_season, season).mean('time').mean('ens') #da_ref = da_ref.pipe(whereocean).pipe(wherebasin, basin).sel(lon=slice(250,360), lat=slice(0,60)) #da_ref = da_ref.pipe(wherebasin, basin).sel(lon=slice(250,360), lat=slice(0,60)) da_ref = da_ref.pipe(wherebasin, basin).sel(lon=slice(250,360), lat=slice(0,40)) daa = da.sel(time=timeslice_late).pipe(sel_season, season).mean('time').mean('ens') \ - da.sel(time=timeslice_early).pipe(sel_season, season).mean('time').mean('ens') #daa = daa.pipe(whereocean).pipe(wherebasin, basin).sel(lon=slice(250,360), lat=slice(0,60)) #daa = daa.pipe(wherebasin, basin).sel(lon=slice(250,360), lat=slice(0,60)) daa = daa.pipe(wherebasin, basin).sel(lon=slice(250,360), lat=slice(0,40)) daa = daa.assign_attrs(long_name=daa.attrs['long_name']+' change') if __name__ == '__main__': from wyconfig import * #my plot settings from geoplots import mapplot fig, ax = plt.subplots() daa.plot.contourf(levels=21, ax=ax) cs = da_ref.plot.contour(levels=11, ax=ax, colors='gray', linewidths=0.5) ax.clabel(cs) #mapplot(fill_continents=True) mapplot(coastlines_width=1) ax.set_title(f'HiRAM {season} wind shear change from {timeslice_early.start}-{timeslice_early.stop} to {timeslice_late.start}-{timeslice_late.stop}') ax.set_xlim(260,360) #savefig if 'savefig' in sys.argv or 's' in sys.argv: figname = __file__.replace('.py', f'.png') if 'overwritefig' in sys.argv or 'o' in sys.argv: wysavefig(figname, overwritefig=True) else: wysavefig(figname) try: tt.check(f'**Done**') except: pass print() if 'notshowfig' in sys.argv or 'n' in sys.argv: pass else: if 'plt' in globals(): plt.show()