Merge PR #377 (Bug fixes for stretched grid rst file generation/benchmarking)

This merge brings PR #377 (Bug fixes for stretched grid restart file
generation and benchmarking, by @lizziel) into the GCPy 1.7.0
development stream.

PR #377 adds the following updates:

1. Remove requirement to pass additional arguments to compare_single_level
   and compare_zonal_mean functions for the case of stretched grid datasets.
   The stretch parameters are now read from the ref and dev datasets passed
   to the functions. If not found the default values are used which correspond
   to no stretching. This enables using stretched grid files in benchmarking
   which previously did not work.

2. Call function to drop GC-Classic variables when regridding a GC-Classic
   restart file to GCHP. This fixes a bug where ilev appeared in GCHP
   restart files, inadvertently triggering MAPL to vertically flip all
   GCHP restart files upon read.

Signed-off-by: Bob Yantosca <yantosca@seas.harvard.edu>

Author: yantosca

CHANGELOG.md

@@ -11,6 +11,7 @@ The format is based on [Keep a Changelog](https://keepachangelog.com/en/1.0.0/),
 - Added "Lint" GitHub Action to check other actions for security issues
 - Added `gcpy_environment_py314.ym1` to specify the GCPy environment packages with Python 3.14 
 - Added GitHub action `build-gcpy-environment-py314.yml` to test building the GCPy environment with Python 3.14
+- Added call to drop GC-Classic variables when regridding a GC-Classic restart file to cubed-sphere
 
 ### Changed
 - Modified criteria for terminating read of log files in `benchmark_scrape_gcclassic_timers.py` to avoid being spoofed by  output that is attached by Intel VTune

gcpy/file_regrid.py

@@ -3,6 +3,7 @@
 (including stretched grids).
 """
 import argparse
+import logging
 import os
 import warnings
 import numpy as np
@@ -1294,6 +1295,7 @@ def drop_classic_vars(
     """
     with xr.set_options(keep_attrs=True):
         if towards_gchp:
+            logging.info("Dropping GC-Classic variables")
             dset = dset.drop_vars(
                 ["P0",
                  "hyam",

gcpy/plot/compare_single_level.py

@@ -65,8 +65,6 @@ def compare_single_level(
         second_ref=None,
         second_dev=None,
         spcdb_dir=os.path.dirname(__file__),
-        sg_ref_path='',
-        sg_dev_path='',
         ll_plot_func='imshow',
         **extra_plot_args
 ):
@@ -177,14 +175,6 @@ def compare_single_level(
         spcdb_dir: str
             Directory containing species_database.yml file.
             Default value: Path of GCPy code repository
-        sg_ref_path: str
-            Path to NetCDF file containing stretched-grid info
-            (in attributes) for the ref dataset
-            Default value: '' (will not be read in)
-        sg_dev_path: str
-            Path to NetCDF file containing stretched-grid info
-            (in attributes) for the dev dataset
-            Default value: '' (will not be read in)
         ll_plot_func: str
             Function to use for lat/lon single level plotting with
             possible values 'imshow' and 'pcolormesh'. imshow is much
@@ -249,40 +239,39 @@ def compare_single_level(
             quiet=True
         )
 
-    sg_ref_params = [1, 170, -90]
-    sg_dev_params = [1, 170, -90]
-    # Get stretched-grid info if passed
-    if sg_ref_path != '':
-        sg_ref_attrs = xr.open_dataset(sg_ref_path).attrs
-        if 'stretch_factor' in sg_ref_attrs:
-            sg_ref_params = [
-                sg_ref_attrs['stretch_factor'],
-                sg_ref_attrs['target_lon'],
-                sg_ref_attrs['target_lat']]
-        elif 'STRETCH_FACTOR' in sg_ref_attrs:
-            sg_ref_params = [
-                sg_ref_attrs['STRETCH_FACTOR'],
-                sg_ref_attrs['TARGET_LON'],
-                sg_ref_attrs['TARGET_LAT']]
-        else:
-            msg = "Stretched grid global parameters missing from ref file"
-            raise ValueError(msg)
-
-    if sg_dev_path != '':
-        sg_dev_attrs = xr.open_dataset(sg_dev_path).attrs
-        if 'stretch_factor' in sg_dev_attrs:
-            sg_dev_params = [
-                sg_dev_attrs['stretch_factor'],
-                sg_dev_attrs['target_lon'],
-                sg_dev_attrs['target_lat']]
-        elif 'STRETCH_FACTOR' in sg_dev_attrs:
-            sg_dev_params = [
-                sg_dev_attrs['STRETCH_FACTOR'],
-                sg_dev_attrs['TARGET_LON'],
-                sg_dev_attrs['TARGET_LAT']]
-        else:
-            msg = "Stretched grid global parameters missing from dev file"
-            raise ValueError(msg)
+    # Get stretched grid info, if any.
+    # Parameter order is stretch factor, target longitude, target latitude.
+    # Stretch factor 1 corresponds with no stretch.
+
+    # Ref stretch attributes
+    print(refdata.attrs)
+    if 'stretch_factor' in refdata.attrs:
+        sg_ref_params = [
+            refdata.attrs['stretch_factor'],
+            refdata.attrs['target_lon'],
+            refdata.attrs['target_lat']]
+    elif 'STRETCH_FACTOR' in refdata.attrs:
+        sg_ref_params = [
+            refdata.attrs['STRETCH_FACTOR'],
+            refdata.attrs['TARGET_LON'],
+            refdata.attrs['TARGET_LAT']]
+    else:
+        sg_ref_params = [1, 170, -90]
+
+    # Dev stretch attributes
+    print(devdata.attrs)
+    if 'stretch_factor' in devdata.attrs:
+        sg_dev_params = [
+            devdata.attrs['stretch_factor'],
+            devdata.attrs['target_lon'],
+            devdata.attrs['target_lat']]
+    elif 'STRETCH_FACTOR' in devdata.attrs:
+        sg_dev_params = [
+            devdata.attrs['STRETCH_FACTOR'],
+            devdata.attrs['TARGET_LON'],
+            devdata.attrs['TARGET_LAT']]
+    else:
+        sg_dev_params = [1, 170, -90]
 
     # Get grid info and regrid if necessary
     [refres, refgridtype, devres, devgridtype, cmpres, cmpgridtype, regridref,

gcpy/plot/compare_zonal_mean.py

@@ -65,8 +65,6 @@ def compare_zonal_mean(
         second_ref=None,
         second_dev=None,
         spcdb_dir=os.path.dirname(__file__),
-        sg_ref_path='',
-        sg_dev_path='',
         ref_vert_params=None,
         dev_vert_params=None,
         **extra_plot_args
@@ -179,14 +177,6 @@ def compare_zonal_mean(
         spcdb_dir: str
             Directory containing species_database.yml file.
             Default value: Path of GCPy code repository
-        sg_ref_path: str
-            Path to NetCDF file containing stretched-grid info
-            (in attributes) for the ref dataset
-            Default value: '' (will not be read in)
-        sg_dev_path: str
-            Path to NetCDF file containing stretched-grid info
-            (in attributes) for the dev dataset
-            Default value: '' (will not be read in)
         ref_vert_params: list(AP, BP) of list-like types
             Hybrid grid parameter A in hPa and B (unitless).
             Needed if ref grid is not 47 or 72 levels.
@@ -312,40 +302,39 @@ def compare_zonal_mean(
         fracrefdata = fracrefdata.isel(lev=ref_pmid_ind)
         fracdevdata = fracdevdata.isel(lev=dev_pmid_ind)
 
-    sg_ref_params = [1, 170, -90]
-    sg_dev_params = [1, 170, -90]
-    # Get stretched-grid info if passed
-    if sg_ref_path != '':
-        sg_ref_attrs = xr.open_dataset(sg_ref_path).attrs
-        if 'stretch_factor' in sg_ref_attrs:
-            sg_ref_params = [
-                sg_ref_attrs['stretch_factor'],
-                sg_ref_attrs['target_lon'],
-                sg_ref_attrs['target_lat']]
-        elif 'STRETCH_FACTOR' in sg_ref_attrs:
-            sg_ref_params = [
-                sg_ref_attrs['STRETCH_FACTOR'],
-                sg_ref_attrs['TARGET_LON'],
-                sg_ref_attrs['TARGET_LAT']]
-        else:
-            msg = "Stretched grid global parameters missing from ref file"
-            raise ValueError(msg)
-
-    if sg_dev_path != '':
-        sg_dev_attrs = xr.open_dataset(sg_dev_path).attrs
-        if 'stretch_factor' in sg_dev_attrs:
-            sg_dev_params = [
-                sg_dev_attrs['stretch_factor'],
-                sg_dev_attrs['target_lon'],
-                sg_dev_attrs['target_lat']]
-        elif 'STRETCH_FACTOR' in sg_dev_attrs:
-            sg_dev_params = [
-                sg_dev_attrs['STRETCH_FACTOR'],
-                sg_dev_attrs['TARGET_LON'],
-                sg_dev_attrs['TARGET_LAT']]
-        else:
-            msg = "Stretched grid global parameters missing from dev file"
-            raise ValueError(msg)
+    # Get stretched grid info, if any.
+    # Parameter order is stretch factor, target longitude, target latitude.
+    # Stretch factor 1 corresponds with no stretch.
+
+    # Ref stretch attributes
+    print(refdata.attrs)
+    if 'stretch_factor' in refdata.attrs:
+        sg_ref_params = [
+            refdata.attrs['stretch_factor'],
+            refdata.attrs['target_lon'],
+            refdata.attrs['target_lat']]
+    elif 'STRETCH_FACTOR' in refdata.attrs:
+        sg_ref_params = [
+            refdata.attrs['STRETCH_FACTOR'],
+            refdata.attrs['TARGET_LON'],
+            refdata.attrs['TARGET_LAT']]
+    else:
+        sg_ref_params = [1, 170, -90]
+
+    # Dev stretch attributes
+    print(devdata.attrs)
+    if 'stretch_factor' in devdata.attrs:
+        sg_dev_params = [
+            devdata.attrs['stretch_factor'],
+            devdata.attrs['target_lon'],
+            devdata.attrs['target_lat']]
+    elif 'STRETCH_FACTOR' in devdata.attrs:
+        sg_dev_params = [
+            devdata.attrs['STRETCH_FACTOR'],
+            devdata.attrs['TARGET_LON'],
+            devdata.attrs['TARGET_LAT']]
+    else:
+        sg_dev_params = [1, 170, -90]
 
     [refres, refgridtype, devres, devgridtype, cmpres, cmpgridtype,
      regridref, regriddev, regridany, refgrid, devgrid, cmpgrid,

gcpy/regrid_restart_file.py

@@ -36,7 +36,7 @@
 import sparselt.esmf
 import sparselt.xr
 import requests
-
+from gcpy.file_regrid import drop_classic_vars
 
 TEMP_FILES = []
 
@@ -515,6 +515,7 @@ def regrid_restart_file(
     is_conversion = input_is_gchp != output_is_gchp
     if is_conversion:
         to_gchp = output_is_gchp
+        dataset = drop_classic_vars(dataset, to_gchp)
         dataset = rename_variables(dataset, to_gchp)
         dataset = reverse_lev(dataset, to_gchp)