Skip to content

to_xgcm_grid_dataset doesn't make z_w_top a Z axis  #97

@ALDepp

Description

@ALDepp

Hi there,

there is a problem with using to_xgcm_grid_dataset when opening a file only partially (e.g. with the function read_dat)

def read_dat(files, variables, pop=False, chunks={'time':-1, 'nlat':100, 'nlon':100, 'z_t':-1}):
    def preprocess(ds):
        return ds[variables].reset_coords(drop=True) # reset coords means they are reset as variables
    ds = xr.open_mfdataset(files, parallel=True, preprocess=preprocess,
                           chunks=chunks,
                           combine='by_coords')
    if pop==True:
        file0 = xr.open_dataset(files[0], chunks=chunks)
        ds.update(file0[['ULONG', 'ULAT', 'TLONG', 'TLAT']])
        file0.close()

    ds
    return ds

My dataset does contain z_w_top:

ds_sort.dims
Frozen(SortedKeysDict({'nlat': 305, 'nlon': 1301, 'time': 73, 'z_t': 62, 'z_w_top': 62, 'z_w_bot': 62}))

but when I convert to a xgcm compatible dataset and grid

metrics = {
    ("X",): ["DXU", "DXT"],  # X distances
    ("Y",): ["DYU", "DYT"],  # Y distances
    ("Z",): ["DZU", "DZT"],  # Z distances
    ("X", "Y"): ["UAREA", "TAREA"],
}

# here we get the xgcm compatible dataset
gridxgcm, dsxgcm = pop_tools.to_xgcm_grid_dataset(
    ds_sort,
    periodic=False,
    metrics=metrics,
    boundary={"X": "extend", "Y": "extend", "Z": "extend"},
)

I get this grid:

<xgcm.Grid>
Y Axis (not periodic, boundary='extend'):
  * center   nlat_t --> right
  * right    nlat_u --> center
X Axis (not periodic, boundary='extend'):
  * center   nlon_t --> right
  * right    nlon_u --> center
Z Axis (not periodic, boundary='extend'):
  * center   z_t --> right
  * right    z_w_bot --> center

While converting a dataset loaded with xr.open_dataset() and converting in the same way results in this grid:

<xgcm.Grid>
Y Axis (not periodic, boundary='extend'):
  * center   nlat_t --> right
  * right    nlat_u --> center
X Axis (not periodic, boundary='extend'):
  * center   nlon_t --> right
  * right    nlon_u --> center
Z Axis (not periodic, boundary='extend'):
  * center   z_t --> left
  * right    z_w_bot --> center
  * left     z_w_top --> center

Why is z_w_top not converted in the upper case?

Metadata

Metadata

Assignees

No one assigned

    Labels

    No labels
    No labels

    Type

    No type

    Projects

    Status

    No status

    Milestone

    No milestone

    Relationships

    None yet

    Development

    No branches or pull requests

    Issue actions