Support for N axis dimensions mapped to N coordinates #343

tomvothecoder · 2022-09-12T21:57:07Z

Description

Related feature specification document: https://docs.google.com/document/d/1BnwNGzb-_pSUFzH5513x99Yv2IOanZlTtA7RvIF7hJA/edit#heading=h.disypbipz1ic

Closes [Bug]: AttributeError: ...multiple values for 'longitude' #183
Closes [Bug]: Issues opening up datasets with multiple coord vars for the same axis mapped to single dim/independent dims (IPSL ESM files) #285
Closes [Feature]: Use dimension coordinates from the input data variable in specific xCDAT APIs #314
Closes [Bug]: swap_lon_axis() renames data variables (non-bounds) #378

Overview:

xcdat APIs previously expected a 1:1 map between dimensions and coordinates. For example, the time dimension would be mapped to time coord var. However, datasets can have N dimensions mapped to N coordinate variables
Examples:
- E3SM datasets can have an N:N map, which might have a "Z" axis with the dimensions and coords "ilev" and "lev"
- IPSL datasets might have 1 time dimension mapped to multiple time coordinates ("coord vars" and "aux coord vars" in CF terms)
More information: https://docs.xarray.dev/en/stable/user-guide/data-structures.html#coordinates

Feature Updates:

Update bounds accessor methods to operate on the dataset or a variable within the dataset using the new var_key kwarg
Update add_missing_bounds() to loop over all axes and their coordinate vars and attempt to add bounds for each coordinate var if they do not exist
Update decoding function (decode_time()) to decode CF and non-CF compliant time coordiantes as cftime objects
- This function is now performed lazily, which means a significant decrease in performance runtime!
Update dataset longitude swapping function (swap_lon_axis()) to loop over longitude coords and their bounds to swap them
Update spatial and temporal accessor class methods to refer to the dimension coordinate variable on the data_var being operated on, rather than the parent dataset

Bug Fixes:

Fix add_bounds() not ignoring 0-dim singleton coords in addition to 1-dim singleton coords

Refactor:

Refactor open_dataset() and open_mfdataset() to remove legacy private functions for checking CF compliance before attempting to decode
- decode_time() now checks for CF compliance
Rename get_axis_coord() to get_dim_coords() and get_axis_dim() to get_dim_keys() for clarity
Update fixtures in fixtures.py to use cftime objects and add decode_times to generate_dataset()
Extract utility function _if_multidim_dask_array_then_load() and reuse in several places in the codebase, usually when manipulating bounds

Checklist

My code follows the style guidelines of this project
I have performed a self-review of my own code
My changes generate no new warnings
Any dependent changes have been merged and published in downstream modules

If applicable:

I have added tests that prove my fix is effective or that my feature works
New and existing unit tests pass with my changes (locally and CI/CD build)
I have commented my code, particularly in hard-to-understand areas
I have made corresponding changes to the documentation
I have noted that this is a breaking change for a major release (fix or feature that would cause existing functionality to not work as expected)

xcdat/axis.py

codecov-commenter · 2022-09-26T23:22:33Z

Codecov Report

Base: 100.00% // Head: 100.00% // No change to project coverage 👍

Coverage data is based on head (8575a09) compared to base (31c0ad1).
Patch coverage: 100.00% of modified lines in pull request are covered.

Additional details and impacted files

@@            Coverage Diff            @@
##              main      #343   +/-   ##
=========================================
  Coverage   100.00%   100.00%           
=========================================
  Files           13        13           
  Lines         1202      1242   +40     
=========================================
+ Hits          1202      1242   +40

Impacted Files	Coverage Δ
xcdat/axis.py	`100.00% <100.00%> (ø)`
xcdat/bounds.py	`100.00% <100.00%> (ø)`
xcdat/dataset.py	`100.00% <100.00%> (ø)`
xcdat/regridder/accessor.py	`100.00% <100.00%> (ø)`
xcdat/regridder/grid.py	`100.00% <100.00%> (ø)`
xcdat/spatial.py	`100.00% <100.00%> (ø)`
xcdat/temporal.py	`100.00% <100.00%> (ø)`
xcdat/utils.py	`100.00% <100.00%> (ø)`

Help us with your feedback. Take ten seconds to tell us how you rate us. Have a feature suggestion? Share it here.

☔ View full report at Codecov.
📢 Do you have feedback about the report comment? Let us know in this issue.

xcdat/dataset.py

tomvothecoder

My initial self-review of the PR.

tests/fixtures.py

tests/test_axis.py

xcdat/bounds.py

xcdat/temporal.py

xcdat/utils.py

tomvothecoder

Hi @pochedls, this PR is ready for review. We were spot-on with our estimate of about a month to get this ready.

Let me know if you would like me to tag others for review since this is a big PR, and I will also do a self-review. I marked areas of interest under this comment if you'd like to look over the implementation details.

Some helpful review resources:

Python script using the sample IPSL and E3SM datasets from the GH issues, sanity check: https://gist.github.com/tomvothecoder/7ae3b9ef34a4c9e8257eb71d9cb51866
Follow the order of Phases of the Algorithm Design in the feature spec document to streamline the review process.

Thanks!

xcdat/axis.py

xcdat/spatial.py

xcdat/temporal.py

xcdat/utils.py

pochedls · 2022-10-11T18:43:25Z

I am still reviewing, but I keep bumping into one major issue:

import xcdat as xc
fn = '/p/css03/esgf_publish/CMIP6/CMIP/CSIRO-ARCCSS/ACCESS-CM2/historical/r1i1p1f1/Amon/tas/gn/v20191108/tas_Amon_ACCESS-CM2_historical_r1i1p1f1_gn_185001-201412.nc`
ds = xc.open_dataset(fn)

2022-10-11 11:38:22,958 [WARNING]: dataset.py(decode_time:333) >> height does not have a 'calendar' attribute set, so it could not be decoded. Set the 'calendar' attribute (ds.height.attrs['calendar]) and try decoding them again

KeyError Traceback (most recent call last)
File .../bounds.py:271, in BoundsAccessor.add_bounds(self, axis, width)
270 try:
--> 271 bounds_key = ds[coord.name].attrs["bounds"]
272 ds[bounds_key]

KeyError: 'bounds'
...

File .../xcdat/bounds.py:352, in BoundsAccessor._create_bounds(self, axis, coord_var, width)
319 """Creates bounds for an axis using its coordinate points.
320
321 Parameters
(...)
349 .. [2] https://cf-xarray.readthedocs.io/en/latest/generated/xarray.Dataset.cf.add_bounds.html#
350 """
351 # Validate coordinate shape.
--> 352 if coord_var.shape[0] <= 1:
353 raise ValueError(
354 f"Cannot generate bounds for coordinate variable '{coord_var.name}'"
355 " which has a length <= 1."
356 )
358 # Retrieve coordinate dimension to calculate the diffs between points.

IndexError: tuple index out of range

I think this issue is related to a singleton height coordinate in this dataset.

If I open with xarray, things work okay (either with cftime=True or False).

docs/faqs.rst

pochedls · 2022-10-10T19:42:03Z

docs/index.rst

+  * Optional selection of single data variable to keep in the Dataset (bounds are also 
+    kept if they exist)


I'm happily surprised this is re-implemented – looking forward to checking this out.

This specific option has been available in open_dataset()/open_mfdataset() but not documented in the features list:

xcdat/xcdat/dataset.py

Lines 51 to 53 in 36a3539

data_var: Optional[str], optional

The key of the non-bounds data variable to keep in the Dataset,

alongside any existing bounds data variables, by default None.

I think you're referring to xcdat guessing the data variable to use in single variable datasets.

xcdat/dataset.py

xcdat/bounds.py

xcdat/axis.py

docs/faqs.rst

xcdat/axis.py

lee1043 · 2022-11-04T05:31:29Z

@pochedls great, thanks for confirming. I will work on this tomorrow or early next week.

pochedls · 2022-11-04T05:41:42Z

xcdat/axis.py

+        keys = keys + obj.cf.coordinates[cf_attrs["coordinate"]]
+    except KeyError:
+        pass
+


Suggested change

if CF_ATTR_MAP['T']['coordinate'] in obj.coords.keys():

keys = keys + [CF_ATTR_MAP[axis]['coordinate']]

This would fix this issue. I think this logic is reasonable and it seems like all unit tests still pass.

pochedls

@tomvothecoder - PR looks great. The lazy loading functionality you loaded goes above and beyond and is a really nice addition. Over all, I like how you cleaned up the open_*dataset and decode_time functionality.

I've verified that the current PR (with my one suggestion) passes all tests, resolves the linked issues/bugs, and spatial averaging and regridding look good.

I have one suggestion, which we can discuss if you're not a fan.

lee1043 · 2022-11-04T19:35:36Z

@tomvothecoder @pochedls - I've verified that current PR does not influence to the temporal and spatial average operations. My statistics (correlation, rmse, standard deviation, etc) that includes spatial and/or temporal average during its calculations were not changed by this PR.

xcdat/dataset.py

tomvothecoder · 2022-11-07T22:13:04Z

xcdat/axis.py

+# A dictionary that maps common variable names to coordinate variables. This
+# map is used as fall-back when coordinate variables don't have CF attributes
+# set for ``cf_xarray`` to interpret using `CF_ATTR_MAP`.
+VAR_NAME_MAP: Dict[CFAxisKey, List[str]] = {
+    "X": ["longitude", "lon"],
+    "Y": ["latitude", "lat"],
+    "T": ["time"],
+    "Z": coordinate_criteria["Z"]["standard_name"],
+}
+


This mapping table addresses the issue here #343 (comment).

I used cf_xarray's coordinate_criteria["Z"] because there are way too many to list ourselves. Let me know if there is a shorter static list we might want to use instead.

https://github.com/xarray-contrib/cf-xarray/blob/614b8ef448112ee51fce60f4db976b6050552dc5/cf_xarray/criteria.py#L40-L54

tomvothecoder

Final self-review code updates

xcdat/bounds.py

xcdat/dataset.py

xcdat/temporal.py

Overview: - `xcdat` APIs previously expected a 1:1 map between axis dimensions and coordinates. For example, the `time` axis dimension would be mapped to `time` coord var. However, datasets can have N dimensions mapped to N coordinate variables - For example, datasets that have N dims mapped to N coord vars include E3SM, which might have a "Z" axis with the dimensions and coords "ilev" and "lev" Feature Updates: - Update bounds accessor methods to operate on the dataset or a variable within the dataset using the new `var_key` kwarg - Update `add_missing_bounds()` to loop over all axes and their coordinate vars and attempt to add bounds for each coordinate var if they do not exist - Update decoding function (`decode_time()`) to decode CF and non-CF compliant time coordiantes as `cftime` objects - Update dataset longitude swapping function (`swap_lon_axis()`) to loop over longitude coords and their bounds to swap them - Update spatial and temporal accessor class methods to refer to the dimension coordinate variable on the `data_var` being operated on, rather than the parent dataset Bug Fixes: - Fix `add_bounds()` not ignoring 0-dim singleton coords in addition to 1-dim singleton coords Refactor: - Refactor `open_dataset()` and `open_mfdataset()` to remove unnecessary preprocessing functions - Rename `get_axis_coord()` to `get_dim_coords()` and `get_axis_dim()` to `get_dim_keys()` - Update fixtures in `fixtures.py` to use `cftime` objects and add `decode_times` to `generate_dataset()` - Extract utility function `_if_multidim_dask_array_then_load()` and reuse in several places in the codebase, usually when manipulating bounds

- This improves runtime performance, since `.values` performs type conversions to `np.array()` - Refactor `_get_cftime_coords()` to use `xr.coding.times.decode_cf_datetime()` for decoding CF compliant time units

- Update `_is_decoded()` to check `.encoding` attributes rather than the object types in the array for performance purposes, since accessing `.values` can involve type conversion

- Refactor imports and fix docstrings

- This is required for cases where time units are different between files in multi-file datasets, but the offsets are the same. More info in the `_preprocess()` docstring - Add raise ValueError for decoding datasets with time coords without CF attrs set

- Add `decode_time()` test for units in days

to interpret common var names - Update tests

tomvothecoder

Nit-picking updates

xcdat/dataset.py

xcdat/bounds.py

xcdat/dataset.py

tomvothecoder · 2022-11-08T00:26:10Z

Thank you @pochedls for your extensive review, and @lee1043 for jumping in to test the temporal APIs.

I performed a final review and now we are good to go. I will be merging this now!

pochedls · 2022-11-08T01:38:47Z

@tomvothecoder – this is an awesome PR – thanks for all of your work on it the past 2 months – it is a great advance for xcdat!

lee1043 · 2022-11-08T02:21:38Z

Excellent work @tomvothecoder and @pochedls, thank you for your hard work for this PR.

tomvothecoder added type: bug Inconsistencies or issues which will cause an issue or problem for users or implementors. type: enhancement New enhancement request Priority: High labels Sep 12, 2022

tomvothecoder self-assigned this Sep 12, 2022

tomvothecoder changed the title ~~Update get_axis_coord() for phase 1~~ Support datasets with axes represented by N dimensions mapped to N coordinates Sep 12, 2022

tomvothecoder force-pushed the bugfix/183-285-multiple-axis-vars branch 2 times, most recently from 963cec5 to f0b8715 Compare September 15, 2022 23:24

tomvothecoder changed the title ~~Support datasets with axes represented by N dimensions mapped to N coordinates~~ [Bug]: Support datasets with axes represented by N dimensions mapped to N coordinates Sep 16, 2022

tomvothecoder commented Sep 26, 2022

View reviewed changes

xcdat/axis.py Outdated Show resolved Hide resolved

tomvothecoder force-pushed the bugfix/183-285-multiple-axis-vars branch from 4214f91 to a592bca Compare September 27, 2022 16:26

This was referenced Sep 28, 2022

[Bug]: Issues opening up datasets with multiple coord vars for the same axis mapped to single dim/independent dims (IPSL ESM files) #285

Closed

[Bug]: AttributeError: ...multiple values for 'longitude' #183

Closed

tomvothecoder commented Sep 29, 2022

View reviewed changes

xcdat/dataset.py Show resolved Hide resolved

tomvothecoder commented Sep 30, 2022

View reviewed changes

tomvothecoder force-pushed the bugfix/183-285-multiple-axis-vars branch from c0ca5ab to 29efb0c Compare September 30, 2022 22:04

tomvothecoder marked this pull request as ready for review September 30, 2022 22:33

tomvothecoder commented Oct 6, 2022

View reviewed changes

pochedls reviewed Oct 11, 2022

View reviewed changes

tomvothecoder force-pushed the bugfix/183-285-multiple-axis-vars branch from 3c34fdb to c76f9fc Compare October 12, 2022 21:53

tomvothecoder commented Oct 12, 2022

View reviewed changes

docs/faqs.rst Outdated Show resolved Hide resolved

tomvothecoder commented Oct 12, 2022

View reviewed changes

xcdat/axis.py Outdated Show resolved Hide resolved

pochedls reviewed Nov 4, 2022

View reviewed changes

pochedls approved these changes Nov 4, 2022

View reviewed changes

tomvothecoder commented Nov 7, 2022

View reviewed changes

xcdat/dataset.py Outdated Show resolved Hide resolved

tomvothecoder force-pushed the bugfix/183-285-multiple-axis-vars branch from 2171067 to f622a8d Compare November 7, 2022 19:27

tomvothecoder commented Nov 7, 2022

View reviewed changes

tomvothecoder added 10 commits November 7, 2022 15:54

Refactor redundant calls to .values

23bcd2d

- This improves runtime performance, since `.values` performs type conversions to `np.array()` - Refactor `_get_cftime_coords()` to use `xr.coding.times.decode_cf_datetime()` for decoding CF compliant time units

Refactor decoding times to a lazy operation

2928e6a

- Update `_is_decoded()` to check `.encoding` attributes rather than the object types in the array for performance purposes, since accessing `.values` can involve type conversion

Add test for decoding with default calendar attr

c7f0e99

- Refactor imports and fix docstrings

Fix incorrect units arg passed to decode_cf_datetime()

1e57bc9

- Add `decode_time()` test for units in days

Update _get_cftime_coords() docstring

3a5d4fd

Update _get_all_coord_keys()

2920ffc

to interpret common var names - Update tests

Apply suggestions from code review

5431b0e

Remove unneeded portion of decode_time() docstring

6591295

tomvothecoder force-pushed the bugfix/183-285-multiple-axis-vars branch from 246f30b to 6591295 Compare November 7, 2022 23:54

tomvothecoder commented Nov 8, 2022

View reviewed changes

tomvothecoder added 3 commits November 7, 2022 16:01

Apply suggestions from code review

d5e5ccd

Update xcdat/bounds.py

ad587d4

Fix test warnings

8575a09

tomvothecoder changed the title ~~[Bug]: Support datasets with axes represented by N dimensions mapped to N coordinates~~ Support for N axis dimensions mapped to N coordinates Nov 8, 2022

tomvothecoder merged commit c6587b8 into main Nov 8, 2022

tomvothecoder deleted the bugfix/183-285-multiple-axis-vars branch November 8, 2022 00:29

Provide feedback

Saved searches

Use saved searches to filter your results more quickly

Support for N axis dimensions mapped to N coordinates #343

Support for N axis dimensions mapped to N coordinates #343

tomvothecoder commented Sep 12, 2022 •

edited

Loading

codecov-commenter commented Sep 26, 2022 •

edited

Loading

tomvothecoder left a comment

tomvothecoder left a comment •

edited

Loading

pochedls commented Oct 11, 2022 •

edited

Loading

pochedls Oct 10, 2022

tomvothecoder Oct 13, 2022

lee1043 commented Nov 4, 2022

pochedls Nov 4, 2022 •

edited

Loading

pochedls left a comment

lee1043 commented Nov 4, 2022

tomvothecoder Nov 7, 2022

tomvothecoder Nov 7, 2022 •

edited

Loading

tomvothecoder Nov 7, 2022

tomvothecoder left a comment

tomvothecoder left a comment

tomvothecoder commented Nov 8, 2022

pochedls commented Nov 8, 2022

lee1043 commented Nov 8, 2022

		* Optional selection of single data variable to keep in the Dataset (bounds are also
		kept if they exist)

	data_var: Optional[str], optional
	The key of the non-bounds data variable to keep in the Dataset,
	alongside any existing bounds data variables, by default None.



	if CF_ATTR_MAP['T']['coordinate'] in obj.coords.keys():
	keys = keys + [CF_ATTR_MAP[axis]['coordinate']]

Support for N axis dimensions mapped to N coordinates #343

Support for N axis dimensions mapped to N coordinates #343

Conversation

tomvothecoder commented Sep 12, 2022 • edited Loading

Description

Checklist

codecov-commenter commented Sep 26, 2022 • edited Loading

Codecov Report

tomvothecoder left a comment

Choose a reason for hiding this comment

tomvothecoder left a comment • edited Loading

Choose a reason for hiding this comment

pochedls commented Oct 11, 2022 • edited Loading

pochedls Oct 10, 2022

Choose a reason for hiding this comment

tomvothecoder Oct 13, 2022

Choose a reason for hiding this comment

lee1043 commented Nov 4, 2022

pochedls Nov 4, 2022 • edited Loading

Choose a reason for hiding this comment

pochedls left a comment

Choose a reason for hiding this comment

lee1043 commented Nov 4, 2022

tomvothecoder Nov 7, 2022

Choose a reason for hiding this comment

tomvothecoder Nov 7, 2022 • edited Loading

Choose a reason for hiding this comment

tomvothecoder Nov 7, 2022

Choose a reason for hiding this comment

tomvothecoder left a comment

Choose a reason for hiding this comment

tomvothecoder left a comment

Choose a reason for hiding this comment

tomvothecoder commented Nov 8, 2022

pochedls commented Nov 8, 2022

lee1043 commented Nov 8, 2022

tomvothecoder commented Sep 12, 2022 •

edited

Loading

codecov-commenter commented Sep 26, 2022 •

edited

Loading

tomvothecoder left a comment •

edited

Loading

pochedls commented Oct 11, 2022 •

edited

Loading

pochedls Nov 4, 2022 •

edited

Loading

tomvothecoder Nov 7, 2022 •

edited

Loading