#!/usr/bin/env python # Wenchang Yang (wenchang@princeton.edu) # Thu Jun 4 11:50:59 EDT 2020 if __name__ == '__main__': from misc.timer import Timer tt = Timer(f'start {__file__}') import sys, os.path, os, glob import xarray as xr, numpy as np, pandas as pd import matplotlib.pyplot as plt #more imports import xaddon maindir = os.path.dirname(os.path.dirname(os.path.abspath(__file__))) if maindir not in sys.path: sys.path.append(maindir) from IBTrACS.data_ntc_cycle import get_cycle as cycle_ntc from ERA5.data_seedxp_map_noLP_basinsum_cycle import get_cycle as cycle_seedxp # if __name__ == '__main__': tt.check('end import') # #start from here basin = 'NA' years = slice('1980', '2018') life = 12 def scatterplot(x=None, y=None, data=None, ax=None, xlabel=None, ylabel=None, title=None, alpha=0.5, tag=None, **kws): if ax is None: ax = plt.gca() if tag is None: tag = '' ax.scatter(x=x, y=y, data=data, alpha=0.5, s=20, **kws) ax.set_xlim(0, None) ax.set_ylim(0, None) if xlabel is not None: ax.set_xlabel(xlabel) if ylabel is not None: ax.set_ylabel(ylabel) if title is not None: ax.set_title(title, loc='left') # linear regression xx, yy = data[x], data[y] rg = yy.linregress.on(xx) ax.plot(xx, rg.predicted) # text if rg.intercept.item()>=0: s = f'{tag}\ny={rg.slope.item():.1f}x $+$ {rg.intercept.item():.2f}\nr$^2$={rg.r.item()**2:.2f}' else: s = f'{tag}\ny={rg.slope.item():.1f}x $-$ {-rg.intercept.item():.2f}\nr$^2$={rg.r.item()**2:.2f}' ax.text(0.02, 1-0.02, s, ha='left', va='top', transform=ax.transAxes) #print(rg) def get_scatter_data(basin='NA', years=None): if years is None: years = slice('1980', '2018') ntc = cycle_ntc(basin=basin, years=years)['mclim'] seedxp = cycle_seedxp(basin=basin, years=years)['mclim'] ds = xr.Dataset(dict( ntc=ntc, seedxp=seedxp )) return ds if __name__ == '__main__': #from wyconfig import * #my plot settings ds = get_scatter_data(basin=basin, years=years) fig, ax = plt.subplots(figsize=(5,4)) figname = __file__.replace('.py', f'_{tt.today()}.png') scatterplot(ax=ax, xlabel='ERA5 ', ylabel='IBTrACS N_TC', x='seedxp', y='ntc', data=ds) plt.savefig(figname) print('[saved]:', figname) tt.check(f'**Done**') plt.show()