#!/usr/bin/env python # Wenchang Yang (wenchang@princeton.edu) # Mon Nov 18 11:56:55 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) import geoxarray from misc import get_kws_from_argv # if __name__ == '__main__': try: tt.check('end import') except: pass # #start from here #ifile = 'Data-added-Ingrid-Jun2025_subset.csv' ifile = 'Updated-data-added-Ingrid-Jun28-2025_subset.csv' #updated data df = pd.read_csv(ifile, index_col=0) lat = xr.DataArray(df['Latitude'], dims='location', coords=[df.location]) lon = xr.DataArray(df['Longitude'], dims='location', coords=[df.location]) lon = lon.where(lon>0, other=lon+360) #0-180 to 180 -> 0 to 360 to be consistent with ERA5 func = lambda da: da.interp(latitude=lat, longitude=lon) dsname = get_kws_from_argv('dsname', 'skin_temperature') if dsname == 'skin_temperature': daname = 'skt' elif dsname == '2m_temperature': daname = 't2m' elif dsname == '2m_temperature_max': daname = 'mx2t' elif dsname == '2m_temperature_min': daname = 'mn2t' else: daname = dsname expname = 'era5' dates = pd.date_range('1979-01', '2024-12', freq='MS') ofile = f'{expname}.{daname}.daily.Ingrid.{dates.year[0]:04d}-{dates.year[-1]:04d}.nc' #era5.t2m.daily.Colorado.1979-2024.nc print('ofile:', ofile) if True: #concat over the year dim das = [] for d in dates: year, month = d.year, d.month #ifile = f'/scratch/gpfs/GEOCLIM/wenchang/analysis/TC/{model}/{expname}/en{en:02d}/analysis_wy/TCIplev/{year:04d}0101.atmos_month.{dsname}.nc' #ifile = f'/tigress/wenchang/MODEL_OUT/{model}/{expname}/en{en:02d}/analysis_wy/TCIplev/{year:04d}0101.atmos_month.{dsname}.nc' #ifile = f'/projects/w/wenchang/data/era5/analysis_wy/TCI_v202201/era5.monthly.{year:04d}.{dsname}.nc' #/projects/w/wenchang/data/era5/analysis_wy/2m_temperature/daily/era5.2m_temperature.daily.2024-12.nc ifile = f'/projects/w/wenchang/data/era5/analysis_wy/{dsname}/daily/era5.{dsname}.daily.{year:04d}-{month:02d}.nc' print(f'{expname} year {year:04d}-{month:02d} {dsname}: {ifile}') ds = xr.open_dataset(ifile) ds = ds[[daname]] #only keep the data variable if 'expver' in ds.coords: ds = ds.drop('expver') #drop the coords that cause problems if 'valid_time' in ds.dims: ds = ds.rename(valid_time='time') print('[renamed]: valid_time -> time') if 'pressure_level' in ds.dims: levels_new = ds.pressure_level.values.astype(int) ds = ds.rename(pressure_level='level').assign_coords(level=levels_new) print('[converted]: pressure_level -> level') da = ds[daname].load() da0 = da.isel(longitude=slice(0,1)) lon0 = da0.longitude.values da0 = da0.assign_coords(longitude=lon0+360) da = xr.concat([da, da0], dim='longitude') if year == 1979 and month == 1: print(da.longitude.values) da = da.pipe(func) #da.isel(time=0).plot(); plt.show(); sys.exit() das.append(da) print('concat over time...') da = xr.concat(das, dim='time') # dss_ens.append(ds) #print('concat over ens...') #ds = xr.concat(dss_ens, dim=pd.Index(ens, name='ens')) #saveing print('saving...') ds = da.to_dataset(name=daname) ds.to_netcdf(ofile) print('[saved]:', ofile) if __name__ == '__main__': #from wyconfig import * #my plot settings #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()