Refactor area_mean_time_series

[forsyth2]

Description

Refactor the area_mean_time_series set. Following the directions at https://github.com/E3SM-Project/e3sm_diags/wiki/CDAT-Migration-FY24-%E2%80%90-General-Guide#getting-started:

git fetch upstream
git checkout -b refactor/662-area_mean_time_series-set upstream/cdat-migration-fy24
mamba env create -f conda-env/dev.yml -n e3sm_diags_dev_662
conda activate e3sm_diags_dev_662
python -m pip install .

• Closes CDAT Migration Phase 2: Refactor area_mean_time_series set #662

TODO

5/9/24

• Refactor area_mean_time_series_driver.py
• Refactor area_mean_time_series_plot.py

5/13/24

• Fix time subsetting issue with string not conforming to ISO-8601 (e.g., "50-01-01" should be "0050-01-01")
• Fix QFLX not deriving correctly -- due to incorrect exception with ref dataset that does not exist
• Perform regression test to ensure results are within relative tolerance (.json and .nc)

Notebook

regression-test-json.ipynb

Viewer links:

• https://portal.nersc.gov/project/e3sm/cdat-migration-fy24/662-area-mean-time-series/viewer/area_mean_time_series/index.html
• https://portal.nersc.gov/project/e3sm/cdat-migration-fy24/main-area-mean-time-series/viewer/area_mean_time_series/index.html

Results

• 770 / 792 metrics are within diff tolerance
• 22 / 792 metrics are not within diff tolerance
• Debug large diffs, specifically these >10%
  • 3 FSNTOA 90S50S.json 0 14.17% 99.968042 114.128953
  • 12 SWCF 90S50S.json 0 10.06% -39.872953 -43.884584
  • 15 TREFHT 50N90N.json 0 27.07% -3.814642 -4.847081

Possible explanation for 22 / 792 metric diffs

1. The cdms2.open slice flag ("ccb") adds an additional time coordinate point when opening up time series datasets ([Refactor]: CDAT Migration: cdms2.open slice flag adds additional time coordinate which affects downstream time series operations and diffs with Xarray code #759).

Here's the logic from dataset.py:

if sub_monthly:
     start_time = "{}-01-01".format(start_year)
     end_time = "{}-01-01".format(str(int(end_year) + 1))
     slice_flag = "co"
 else:
     start_time = "{}-01-15".format(start_year)
     end_time = "{}-12-15".format(end_year)
     slice_flag = "ccb"

Result: Adding code to handle this does not change the results

2. The way masking is applied might be slightly different between xESMF and CDAT's ESMF?

5/24/24

• Add code to handle [Refactor]: CDAT Migration: cdms2.open slice flag adds additional time coordinate which affects downstream time series operations and diffs with Xarray code #759 (need to test in a separate PR)

5/29/24

• Figure out why 22/792 metrics are mismatching
  • Output yearly averages to .nc on both branches
  • Run regression testing notebook for .nc files
  • Debug the mismatching time coordinates for yearly averages (all 72 files mismatch) -- different months too
  • For some reason, cdutil.YEAR changes the month and day of the grouped time coordinates to July 2nd.

    # %%
    import cdms2
    import cdutil
    import xcdat as xc
    
    # %%
    # CDAT
    # 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_key = "FSNS"
    
    ds_cdat2 = cdms2.open(fn1)
    var = ds_cdat2(var_key)
    
    var_avg = cdutil.averager(var, axis="xy")
    var_avg_year = cdutil.YEAR(var_avg)
    
    print(var_avg_year.getTime().asComponentTime()[-3:-1])
    # [2011-7-2 12:0:0.0, 2012-7-2 12:0:0.0]
    
    # %%
    # xCDAT
    ds_xcdat = xc.open_dataset(fn1, decode_times=True)
    ds_xcdat_avg = ds_xcdat.spatial.average(var_key, axis=["X", "Y"])
    ds_xcdat_avg_year = ds_xcdat_avg.temporal.group_average(var_key, freq="year")
    
    var_avg_year_xc = ds_xcdat_avg_year[var_key]
    
    print(var_avg_year_xc.time.values[-3:-1])
    # [cftime.DatetimeNoLeap(2012, 1, 1, 0, 0, 0, 0, has_year_zero=True)
    #  cftime.DatetimeNoLeap(2013, 1, 1, 0, 0, 0, 0, has_year_zero=True)]

  • The main branch has a function called adjust_time_from_time_bounds() which centers time coordinates using time bounds for non-submonthly data (code block). The logic on cdat-migration-fy24 was incorrect and did not center time coordinates correctly, resulting in mismatching time coordinates and in some cases an extra time coordinate point after calculating annual averages. This commits fixes the logic to center time bounds on the dev branch.
  • Fix viewer plot description -- update to driver to get "long_name" attribute from the data variable
  • Fix qflxconvert_units() not maintaining attributes, including "long_name"

All NetCDF files and metrics json are within relative tolerance (1e-5). This PR is good to go.

Plots are also identical (example dev and main).

Next steps

In a separate PR, validate fix for #759 using ex1.py for lat_lon set and time series files.

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:

• New and existing unit tests pass with my changes (locally and CI/CD build)
• I have added tests that prove my fix is effective or that my feature works
• 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)

[forsyth2]

Using this script:

diags_set=area_mean_time_series

cd e3sm_diags
# https://anaconda.org/cdat/repo?sort=_name&sort_order=asc
for cdat_package in cd77 cdat cdat_compute_graph cdat_info cdat_jupyter_vcdat cdms2 cdp cdtime cdutil cibots compute_graph django-allow-cidr dv3d e3sm_nex esgf-compute-api esgf-compute-cert esgf-search ezget genutil image-compare jupyter-vcdat jupyterlab-git lats libcdms lmoments mesalib myproxyclient output_viewer selenium_testlib testsrunner thermo toolchain vcdat vcs vcs-js vcsaddons vtk vtk-cdat wk xmgrace
do
    git grep -n ${cdat_package} e3sm_diags/driver/${diags_set}_driver.py > /dev/null
    if [ $? == 0 ]; then
	echo ""
	echo "git grep -n ${cdat_package} e3sm_diags/driver/${diags_set}_driver.py"
	git grep -n ${cdat_package} e3sm_diags/driver/${diags_set}_driver.py
    fi
    git grep -n ${cdat_package} e3sm_diags/plot/cartopy/${diags_set}_plot.py > /dev/null
    if [ $? == 0 ]; then
	echo ""
	echo "git grep -n ${cdat_package} e3sm_diags/plot/cartopy/${diags_set}_plot.py"
	git grep -n ${cdat_package} e3sm_diags/plot/cartopy/${diags_set}_plot.py
    fi
done

I get:

git grep -n cdms2 e3sm_diags/driver/area_mean_time_series_driver.py
e3sm_diags/driver/area_mean_time_series_driver.py:8:import cdms2
e3sm_diags/driver/area_mean_time_series_driver.py:70:        with cdms2.open(mask_path) as f:

git grep -n cdutil e3sm_diags/driver/area_mean_time_series_driver.py
e3sm_diags/driver/area_mean_time_series_driver.py:9:import cdutil
e3sm_diags/driver/area_mean_time_series_driver.py:97:            test_domain = cdutil.averager(test_domain, axis="xy")
e3sm_diags/driver/area_mean_time_series_driver.py:98:            cdutil.setTimeBoundsMonthly(test_domain)
e3sm_diags/driver/area_mean_time_series_driver.py:99:            test_domain_year = cdutil.YEAR(test_domain)
e3sm_diags/driver/area_mean_time_series_driver.py:100:            # add back attributes since they got lost after applying cdutil.YEAR
e3sm_diags/driver/area_mean_time_series_driver.py:123:                    ref_domain = cdutil.averager(ref_domain, axis="xy")
e3sm_diags/driver/area_mean_time_series_driver.py:124:                    cdutil.setTimeBoundsMonthly(ref_domain)
e3sm_diags/driver/area_mean_time_series_driver.py:133:                    ref_domain_year = cdutil.YEAR(ref_domain)

[tomvothecoder]

I recently squashed the cdat-migration-fy24 branch to a few commits and just rebased this branch on top of it.

Can you make sure to checkout the latest version of his branch? I also recommend stashing and reapplying any changes you might have.

[forsyth2]

@tomvothecoder Do you know what parameters are best for testing area_mean_time_series? I tried python run_area_mean_time_series.py -d diags.cfg using all the variables listed on the defaults but I keep getting this error:

RuntimeError: No parameters we able to be extracted. Please check the parameters you defined.

[forsyth2]

Initial analysis of refactoring needs:

Driver

e3sm_diags/e3sm_diags/driver/area_mean_time_series_driver.py

Lines 67 to 72 in f6c4fdf

	mask_path = os.path.join(
	e3sm_diags.INSTALL_PATH, "acme_ne30_ocean_land_mask.nc"
	)
	with cdms2.open(mask_path) as f:
	land_frac = f("LANDFRAC")
	ocean_frac = f("OCNFRAC")

=> ds_land_sea_mask: xr.Dataset = test_ds._get_land_sea_mask(season) (as in https://github.com/E3SM-Project/e3sm_diags/pull/749/files#r1467834914)

e3sm_diags/e3sm_diags/driver/area_mean_time_series_driver.py

Lines 97 to 100 in f6c4fdf

	test_domain = cdutil.averager(test_domain, axis="xy")
	cdutil.setTimeBoundsMonthly(test_domain)
	test_domain_year = cdutil.YEAR(test_domain)
	# add back attributes since they got lost after applying cdutil.YEAR

e3sm_diags/e3sm_diags/driver/area_mean_time_series_driver.py

Lines 123 to 124 in f6c4fdf

	ref_domain = cdutil.averager(ref_domain, axis="xy")
	cdutil.setTimeBoundsMonthly(ref_domain)

e3sm_diags/e3sm_diags/driver/area_mean_time_series_driver.py

Line 133 in f6c4fdf

ref_domain_year = cdutil.YEAR(ref_domain)

=> All the cdutil references need to be changed out here. I'm not sure with what. I'm not really seeing any cdutil in the polar set (#749), and not much relevant in the lat_lon set (#677), to compare to.

Plot

e3sm_diags/e3sm_diags/plot/cartopy/area_mean_time_series_plot.py

Lines 26 to 36 in f6c4fdf

	panel = [
	(0.1, 0.70, 0.25, 0.225),
	(0.4, 0.70, 0.25, 0.225),
	(0.7, 0.70, 0.25, 0.225),
	(0.1, 0.38, 0.25, 0.225),
	(0.4, 0.38, 0.25, 0.225),
	(0.7, 0.38, 0.25, 0.225),
	(0.1, 0.06, 0.25, 0.225),
	(0.4, 0.06, 0.25, 0.225),
	(0.7, 0.06, 0.25, 0.225),
	]

=> Rename to PANEL_CFG

e3sm_diags/e3sm_diags/plot/cartopy/area_mean_time_series_plot.py

Line 40 in f6c4fdf

border = (-0.047, -0.06, 0.006, 0.03)

=> Rename to BORDER_PADDING

The two changes above will be passed into:

e3sm_diags/e3sm_diags/plot/utils.py

Line 59 in f6c4fdf

def _save_plot(

e3sm_diags/e3sm_diags/plot/cartopy/area_mean_time_series_plot.py

Lines 100 to 142 in f6c4fdf

	# Save the figure.
	output_file_name = var
	for f in parameter.output_format:
	f = f.lower().split(".")[-1]
	fnm = os.path.join(
	get_output_dir(parameter.current_set, parameter),
	output_file_name + "." + f,
	)
	plt.savefig(fnm)
	# Get the filename that the user has passed in and display that.
	fnm = os.path.join(
	get_output_dir(parameter.current_set, parameter),
	output_file_name + "." + f,
	)
	logger.info(f"Plot saved in: {fnm}")
	# Save individual subplots
	for f in parameter.output_format_subplot:
	fnm = os.path.join(
	get_output_dir(parameter.current_set, parameter),
	parameter.output_file,
	)
	page = fig.get_size_inches()
	i = 0
	for p in panel:
	# Extent of subplot
	subpage = np.array(p).reshape(2, 2)
	subpage[1, :] = subpage[0, :] + subpage[1, :]
	subpage = subpage + np.array(border).reshape(2, 2)
	subpage = list(((subpage) * page).flatten()) # type: ignore
	extent = matplotlib.transforms.Bbox.from_extents(*subpage)
	# Save subplot
	fname = fnm + ".%i." % (i) + f
	plt.savefig(fname, bbox_inches=extent)
	orig_fnm = os.path.join(
	get_output_dir(parameter.current_set, parameter),
	parameter.output_file,
	)
	fname = orig_fnm + ".%i." % (i) + f
	logger.info(f"Sub-plot saved in: {fname}")
	i += 1

=> Replace with a call to

e3sm_diags/e3sm_diags/plot/utils.py

Line 59 in f6c4fdf

def _save_plot(

[tomvothecoder]

This function can stored in a general utility module later if needed.