#!/usr/bin/env python # Wenchang Yang (wenchang@princeton.edu) # Fri Feb 17 10:30:04 EST 2023 if __name__ == '__main__': import sys,os from misc.timer import Timer tt = Timer(f'[{os.getcwd()}] start ' + ' '.join(sys.argv)) import sys, os.path, os, glob, datetime import xarray as xr, numpy as np, pandas as pd, matplotlib.pyplot as plt #more imports from misc import get_kws_from_argv # if __name__ == '__main__': tt.check('end import') # #start from here #location: Nepal Lalitpur latitude longitude locname, lat0, lon0 = 'Lalitpur', 27.6588, 85.3247 dfname = '2m_temperature' dfname = get_kws_from_argv('dfname', dfname) daname = {'2m_temperature': 'tas', '2m_temperature_max': 'tasmax', '2m_temperature_min': 'tasmin', 'evaporation': 'evap'}[dfname] tstart = '1979-01' #param tstop = '2022-11' if daname=='tas' else '2022-10' #param trange = pd.date_range(tstart, tstop, freq='MS') year_month_vec = [(t.year, t.month) for t in trange] method = 'selnearest' #alternative is interp method = get_kws_from_argv('method', method) ofile = f'era5.{daname}.daily.{method}{locname}.{tstart}_{tstop}.nc' if os.path.exists(ofile): print('[exists]:', ofile) xr.open_dataarray(ofile).plot() sys.exit() das = [] for year,month in year_month_vec: print(f'{year}-{month:02} of {tstop}') #example: ifile = /tigress/wenchang/data/era5/analysis_wy/2m_temperature/daily/era5.2m_temperature.daily.2022-11.nc ifile = f'/tigress/wenchang/data/era5/analysis_wy/{dfname}/daily/era5.{dfname}.daily.{year}-{month:02d}.nc' print(ifile) da = xr.open_dataarray(ifile) units = da.attrs['units'] if method == 'selnearest': da = da.sel(latitude=lat0, longitude=lon0, method='nearest') elif method == 'interp': da = da.interp(latitude=lat0, longitude=lon0) da.load() das.append(da) print('concat...') da = xr.concat(das, dim='time') #convert units if dfname.startswith('2m_temperature'): da = da - 273.15 #K to C da.attrs['units'] = 'degC' else: da.attrs['units'] = units #save da.to_dataset(name=daname).to_netcdf(ofile) print('[saved]:', ofile) if __name__ == '__main__': from wyconfig import * #my plot settings da.plot() #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) tt.check(f'**Done**') print() if 'notshowfig' in sys.argv or 'n' in sys.argv: pass else: if 'plt' in globals(): plt.show()