Skip to content

Commit

Permalink
Add function to replicate "ccb" slice flag
Browse files Browse the repository at this point in the history
  • Loading branch information
tomvothecoder committed May 14, 2024
1 parent 23be59e commit 59466c4
Show file tree
Hide file tree
Showing 2 changed files with 99 additions and 0 deletions.
Original file line number Diff line number Diff line change
@@ -0,0 +1,64 @@
# %%
import xcdat as xc

# Parameters
# Time coordinates at the first day of the month (1850-01-01)
fn1 = "/global/cfs/cdirs/e3sm/e3sm_diags/test_model_data_for_acme_diags/time-series/E3SM_v1/FSNS_185001_201312.nc"
var = "FSNS"
start_time = "2011-01-15"
end_time = "2013-12-15"

slice_flag = "ccb"

# %%
ds_xcdat = xc.open_dataset(fn1, decode_times=True)

ds_xcdat_s = ds_xcdat.sel(time=slice(start_time, end_time))

# %%
# 35 2013-12-01 00:00:00
print(len(ds_xcdat_s.time), ds_xcdat_s.time.values[-1])

# [cftime.DatetimeNoLeap(2013, 11, 1, 0, 0, 0, 0, has_year_zero=True)
# cftime.DatetimeNoLeap(2013, 12, 1, 0, 0, 0, 0, has_year_zero=True)]
print(ds_xcdat_s.time_bnds.values[-1])


"""Possible solutions:
1. Use last time coordinate (2014, 1, 1) with subsetting
2. Use start value of the last bounds
- (2013, 12, 1), (2014, 1, 1)
"""

# 1. Get the last time coordinate as ds1
# 2. Combine with ds2 (subsetted dataset)

ds1 = ds_xcdat.isel(time=-1)

# %%
import xarray as xr

ds_final = xr.concat([ds_xcdat_s, ds1], dim="time")

# %%
import cdms2
import xcdat as xc

# Parameters
# Time coordinates at the first day of the month (1850-01-01)
fn1 = "/global/cfs/cdirs/e3sm/e3sm_diags/test_model_data_for_acme_diags/time-series/E3SM_v1/FSNS_185001_201312.nc"
var = "FSNS"
start_time = "2011-01-15"
end_time = "2013-12-15"

slice_flag = "ccb"


ds_cdat2 = cdms2.open(fn1)
var_s = ds_cdat2(var, time=(start_time, end_time, slice_flag))(squeeze=1)

# %%
import numpy as np

np.testing.assert_allclose(ds_final["FSNS"].values, var_s.data, rtol=0, atol=0)
# %%
35 changes: 35 additions & 0 deletions e3sm_diags/driver/utils/dataset_xr.py
Original file line number Diff line number Diff line change
Expand Up @@ -979,8 +979,43 @@ def _get_time_series_dataset_obj(self, var) -> xr.Dataset:
ds = xr.open_dataset(filepath, decode_times=True, use_cftime=True)
ds_subset = ds.sel(time=time_slice).squeeze()

if not self.is_sub_monthly:
ds_subset = self._extend_time_for_non_submonthly_data(ds_subset)

return ds_subset

def _extend_time_for_non_submonthly_data(self, ds: xr.Dataset) -> xr.Dataset:
"""Extend the time coordinates for non-sub-monthly data.
This function replicates the cdms2/cdutil "ccb" slice flag used for
subsetting. "ccb" only allows the right side to be closed. It will use
bounds to determine the right bound value to add the coordinate point.
For example, with "2013-12-01" and "ccb":
1. Actual end time: "2013-12-01"
2. Slice end time: "2013-12-15"
3. Bound values: ["2013-12-01", "2014-1-1"]
4. Use "2014-1-1" as last coordinate point.
More info can be found in the API docstrings here: https://github.com/CDAT/cdms/blob/3f8c7baa359f428628a666652ecf361764dc7b7a/Lib/axis.py#L392-L414
Parameters
----------
ds : xr.Dataset
The non-sub-monthly dataset.
Returns
-------
xr.Dataset
The non-sub-monthly dataset with an extra time coordinate point.12
"""
time_dim = xc.get_dim_keys(ds, axis="T")

ds_last_time = ds.isel({time_dim: -1})
ds_final = xr.concat([ds, ds_last_time], dim=time_dim)

return ds_final

def _get_timeseries_filepath(self, root_path: str, var_key: str) -> str:
"""Get the matching variable time series filepath.
Expand Down

0 comments on commit 59466c4

Please sign in to comment.