raphaeldussin / xoverturning

ocean overturning in xarray
19 stars 2 forks source link

xoverturning

MOM6 ocean overturning in xarray

QuickStart guide:

import xarray as xr
import matplotlib.pyplot as plt

ppdir = '/archive/Raphael.Dussin/FMS2019.01.03_devgfdl_20201120/CM4_piControl_c96_OM4p25_half_kdadd/gfdl.ncrc4-intel18-prod-openmp/pp/ocean_annual_z'
ds = xr.open_mfdataset([f"{ppdir}/ts/annual/10yr/ocean_annual_z.0021-0030.umo.nc",
                        f"{ppdir}/ts/annual/10yr/ocean_annual_z.0021-0030.vmo.nc",
                        f"{ppdir}/ts/annual/10yr/ocean_annual_z.0021-0030.uhml.nc",
                        f"{ppdir}/ts/annual/10yr/ocean_annual_z.0021-0030.vhml.nc"])

dsgrid = xr.open_dataset(f"{ppdir}/ocean_annual_z.static.nc")

# this is only useful for rotation
#hgrid = xr.open_dataset('/archive/gold/datasets/OM4_025/mosaic.v20170622.unpacked/ocean_hgrid.nc')
#ds['angle_dx'] = xr.DataArray(hgrid['angle_dx'].values[1::2,1::2], dims=('yh','xh'))
from xoverturning import calcmoc

# global MOC
moc = calcmoc(ds, dsgrid=dsgrid)
# atlantic
moc = calcmoc(ds, dsgrid=dsgrid, basin='atl-arc')
# remove thickness flux to mixed layer
moc = calcmoc(ds, dsgrid=dsgrid, basin='atl-arc', remove_hml=True)
# rotate velocities to True North
moc = calcmoc(ds, dsgrid=dsgrid, basin='atl-arc', rotate=True)
# mask the output (only for z)
moc = calcmoc(ds, dsgrid=dsgrid, basin='atl-arc', mask_output=True)
ppdir = '/archive/Raphael.Dussin/FMS2019.01.03_devgfdl_20201120/CM4_piControl_c96_OM4p25_half_kdadd/gfdl.ncrc4-intel18-prod-openmp/pp/ocean_annual_rho2'
ds = xr.open_mfdataset([f"{ppdir}/ts/annual/10yr/ocean_annual_rho2.0021-0030.umo.nc",
                        f"{ppdir}/ts/annual/10yr/ocean_annual_rho2.0021-0030.vmo.nc",
                        f"{ppdir}/ts/annual/10yr/ocean_annual_rho2.0021-0030.uhml.nc",
                        f"{ppdir}/ts/annual/10yr/ocean_annual_rho2.0021-0030.vhml.nc"])

dsgrid = xr.open_dataset(f"{ppdir}/ocean_annual_rho2.static.nc")

moc = calcmoc(ds, dsgrid=dsgrid, vertical='rho2')

Computing derived quantities:

maxmoc = moc.max(dim=['yq', 'z_i'])
maxmoc = moc.min(dim=['yq', 'z_i'])
maxmoc_265 = moc.sel(yq=slice(26.4,26.6)).max(dim=['z_i'])
maxsfn = moc.sel(yq=slice(20.0, 80.0), z_i=slice(500.0, 2500.0)).max(dim=['yq', 'z_i'])