<jupyter_start><jupyter_text>Adding Back the Tides - CIOPSAn effort to make the daily files more accurate as they are currently lacking the tidal pumping that is so important to the flow of the Salish Sea<jupyter_code>import xarray as xr
from pathlib import Path
import numpy as np
import datetime as dt
float(np.nan)<jupyter_output><empty_output><jupyter_text>unlike what we did with SSC in this file make all add back tides days seperate to reduce space requirement just doing 2018 for nowplus 30 days on the end for running purposes<jupyter_code>startday = [dt.datetime(2017,12,31)+dt.timedelta(days=i) for i in range(int(399))]
print(startday[-1])
# for i in range(131,len(startday)-1):
#     print(i)
folders = [dt.datetime(2017,12,31)+dt.timedelta(days=7*(i+1)) for i in range(int(57))]
print(folders[-1])
folders = np.repeat(folders,7)
# start = dt.datetime(2015,11,22) #start day of your run

#dates for each run
# numdays = 2 #15 for all except the last run (7) 
# date_list = [start + dt.timedelta(days=x) for x in range(numdays)]
# folderday = [start+ dt.timedelta(days=7),start+ dt.timedelta(days=7*2)]#,start+ dt.timedelta(days=7*3)]
# folderday = np.repeat(folderday,7)
# folderday = folderday[:-6]
date_list_daily = [startday[0],startday[1]]
folderday_daily = [folders[0], folders[1]]

# for whatever reason the naming it raelly really weird and the daily and the hourly files names for the SAME DAY aren't\
# on the same day.. must take hourly files a day ahead woooooooooooo
date_list_hourly = [startday[1],startday[2]]
folderday_hourly = [folders[1], folders[2]]
"{:%Y%m%d}00/BC12_1h_grid_U_2D_{:%Y%m%d}_{:%Y%m%d}.nc".format(folderday_hourly[0], date_list_hourly[0], date_list_hourly[0])
path = Path("/ocean/mdunphy/CIOPSW-BC12/")

drop_vars = (
    "sbu", "tauuo", "time_counter_bounds", "time_instant_bounds", "uos", "time_instant",
)

files = [sorted(path.glob("{:%Y%m%d}00/BC12_1h_grid_U_2D_{:%Y%m%d}_{:%Y%m%d}.nc".format(folderday_hourly[i], date_list_hourly[i], date_list_hourly[i]))) for i in range(len(date_list_hourly))]

mydata = xr.open_mfdataset(files, drop_variables=drop_vars)
ut_h = mydata['ubar']
files
ut_h
ut_h[0,646,267].values
drop_vars = (
    "sbv", "tauvo", "time_counter_bounds", "time_instant_bounds", "vos", "time_instant",
)

files = [sorted(path.glob("{:%Y%m%d}00/BC12_1h_grid_V_2D_{:%Y%m%d}_{:%Y%m%d}.nc".format(folderday_hourly[i], date_list_hourly[i], date_list_hourly[i]))) for i in range(len(date_list_hourly))]

mydata = xr.open_mfdataset(files, drop_variables=drop_vars)
vt_h = mydata['vbar']
# 2D hourly outputs are alraedy in bartropic form so you dont need to do any conversion!
# Now get the required data from the daily files
drop_vars = (
    "depthu_bounds", "nav_lat", "nav_lon", 'time_counter_bounds', 'time_instant',
    'time_instant_bounds', 
)

files = [sorted(path.glob("{:%Y%m%d}00/BC12_1d_grid_U_{:%Y%m%d}_{:%Y%m%d}.nc".format(folderday_daily[i], date_list_daily[i], date_list_daily[i]))) for i in range(len(date_list_daily))]

mydata = xr.open_mfdataset(files, drop_variables=drop_vars)
u_d = mydata['uo']

drop_vars = (
    "depthv_bounds", "nav_lat", "nav_lon", 'time_counter_bounds', 'time_instant',
    'time_instant_bounds', 
)

files = [sorted(path.glob("{:%Y%m%d}00/BC12_1d_grid_V_{:%Y%m%d}_{:%Y%m%d}.nc".format(folderday_daily[i], date_list_daily[i], date_list_daily[i]))) for i in range(len(date_list_daily))]

mydata = xr.open_mfdataset(files, drop_variables=drop_vars)
v_d = mydata['vo']
u_d = u_d.where(u_d != 0)
u_d[0,:,267,646].values
u_d[0,:,240,646].values
u_d
# we WILL need to do some conversions here, so get e3t from the mesh_mask file
mydata = xr.open_dataset("/ocean/mdunphy/CIOPSW-BC12/grid/mesh_mask_Bathymetry_NEP36_714x1020_SRTM30v11_NOAA3sec_WCTSS_JdeFSalSea.nc")
e3t = mydata['e3t_0']
e3t
# convert e3t to e3u and to e3v
e3t_xshift = e3t.shift(x=-1,fill_value=0)
e3u = e3t_xshift+e3t
e3u = e3u*0.5
e3u = e3u.rename({'z': 'depthu'})

e3t_yshift = e3t.shift(y=-1,fill_value=0)
e3v = e3t_yshift+e3t
e3v = e3v*0.5
e3v = e3v.rename({'z': 'depthv'})
#calcuate bartropic component
ut_d = (u_d*e3u[0,:,:,:]).sum(dim='depthu')/e3u[0,:,:,:].sum(dim='depthu')
ut_d[0,646,267].values
#subtract from u to get baroclinic component
uc_d = u_d-ut_d #does this work even though their ut_d lacks the depth dimension?<jupyter_output><empty_output><jupyter_text>interpolate + resample uc_d to get it in an hourly format<jupyter_code>uc_d[0,:,646,267].values
offset = dt.timedelta(hours=1) # daily avg timestamp as next day time 0, offset by -23 hours to line up with hourly data
uc_h_interp = uc_d.resample(time_counter="1H", loffset=offset).interpolate("linear")
#instead of taking 12hours off each side like when done with SSC, this method takes off 24hours on the end!
# , loffset=offset
uc_h_interp
uc_h_interp[0,0,646,267].values
u_new = ut_h  + uc_h_interp
u_new = u_new.isel(time_counter = np.arange(0,24,1))
u_new = u_new.rename('vozocrtx')
u_new
u_new[0,646,267,0].values
encoding={
          "vozocrtx": {"zlib": True, "complevel": 4, "_FillValue": np.nan}
}
encoding
# encoding = {
#     "depthu": {
#         "zlib": True, "complevel": 4,
#         "chunksizes": {
#             "depthu": 75,
#         }
#     },
#     "vozocrtx": {
#         "zlib": True, "complevel": 4,
#         "chunksizes": {
#             "time_counter": 1, "y": 1020, "x": 714, "depthu": 75,
#         }
#     }
# }

# encoding = {
#     "depthu": {
#         "zlib": True, "complevel": 4,
#         "chunksizes": [75],
#     },
#     "time_counter": {
#         "zlib": True, "complevel": 4,
#         "chunksizes": [1],
#         "units": "hours since 1970-01-01 00:00:00"
#     },
#     "nav_lat": {
#         "zlib": True, "complevel": 4,
#         "chunksizes": {
#         "y": 1020, "x": 714,
#         }
#     },
#     "nav_lon": {
#         "zlib": True, "complevel": 4,
#         "chunksizes": {
#         "y": 1020, "x": 714,
#         }
#     },
#     "vozocrtx": {
#         "zlib": True, "complevel": 4,
#         "chunksizes": {
#         "time_counter": 1, "y": 1020, "x": 714, "depthu": 75,
#         }
#     }
# }
path = "/ocean/rbeutel/data/"
u_new.to_netcdf(str(path)+'u_new_{:%d%b%y}_fill.nc'.format(date_list_hourly[0],date_list_hourly[0]), encoding=encoding)
print('u_new_{:%d%b%y}_{:%d%b%y}.nc'.format(date_list_hourly[0],date_list_hourly[0]))
#calcuate bartropic component
vt_d = (v_d*e3v[0,:,:,:]).sum(dim='depthv')/e3v[0,:,:,:].sum(dim='depthv')
#subtract from v to get baroclinic component
vc_d = v_d-vt_d 
vc_h_interp = vc_d.resample(time_counter="1H", loffset=offset).interpolate("linear")
v_new = vt_h  + vc_h_interp
v_new = v_new.isel(time_counter = np.arange(0,24,1))
v_new = v_new.rename('vomecrty')
v_new
v_new.to_netcdf(str(path)+'v_new_{:%d%b%y}_{:%d%b%y}.nc'.format(date_list[0],date_list[-1]))<jupyter_output><empty_output>
{'fraction_non_alphanumeric': 0.1714997070884593, 'fraction_numerical': 0.042765084944346804, 'mean_word_length': 4.958987783595114, 'pattern_counts': {'":': 46, '<': 13, '<?xml version=': 0, '>': 13, 'https://': 0, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 1, 'word_list_counts': {'cursed_substrings.json': 1, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
{'chain_length': '3', 'dir': 'jupyter-structured-clean-dedup', 'id': '502156', 'license': 'Apache-2.0', 'n_tokens_mistral': 3048, 'n_tokens_neox': 2854, 'n_words': 554, 'path': 'ArtificialTides/BaroclinicBarotropic-BC12.ipynb', 'repo_name': 'UBC-MOAD/analysis-becca', 'size': '6828'}