Merge pull request #1097 from xylar/fix-custom-thermal-forcing

Fix custom climatology maps

Author: xylar

mpas_analysis/default.cfg

@@ -2192,7 +2192,7 @@ availableVariables = {
          'mpas': ['timeMonthly_avg_activeTracers_temperature',
                   'timeMonthly_avg_activeTracers_salinity',
                   'timeMonthly_avg_density',
-                  'timeMonthly_avg_activeTracers_layerThickness']},
+                  'timeMonthly_avg_layerThickness']},
     'zonalVelocity':
         {'title': 'Zonal Velocity',
          'units': r'm s$$^{-1}$$',
@@ -2202,7 +2202,7 @@ availableVariables = {
          'units': r'm s$$^{-1}$$',
          'mpas': ['timeMonthly_avg_velocityMeridional']},
     'velocityMagnitude':
-        {'title': 'Zonal Velocity',
+        {'title': 'Velocity Magnitude',
          'units': r'm s$$^{-1}$$',
          'mpas': ['timeMonthly_avg_velocityZonal',
                   'timeMonthly_avg_velocityMeridional']},
@@ -2225,8 +2225,8 @@ availableVariables = {
          'has_depth': False},
     }
 
-# a list of fields top plot for each transect.  All supported fields are listed
-# below
+# a list of fields top plot for each depth slice.  All supported fields are
+# listed above
 variables = []
 
 

mpas_analysis/ocean/climatology_map_custom.py

@@ -10,6 +10,7 @@
 # https://raw.githubusercontent.com/MPAS-Dev/MPAS-Analysis/main/LICENSE
 
 import numpy as np
+import xarray as xr
 from mpas_tools.cime.constants import constants as cime_constants
 
 from mpas_analysis.ocean.remap_depth_slices_subtask import (
@@ -173,7 +174,7 @@ def __init__(self, config, mpasClimatologyTask, controlConfig=None):
                             subtaskName=subtaskName)
 
                         subtask.set_plot_info(
-                            outFileLabel=varName,
+                            outFileLabel=f'cust_{varName}',
                             fieldNameInTitle=title,
                             mpasFieldName=mpasVarName,
                             refFieldName=mpasVarName,
@@ -269,14 +270,19 @@ def customize_masked_climatology(self, climatology, season):
             the modified climatology data set
         """
 
-        # first, call the base class's version of this function so we extract
-        # the desired slices.
-        climatology = super().customize_masked_climatology(climatology,
-                                                           season)
+        # first, compute the derived variables, which may rely on having the
+        # full 3D variables available
+
         derivedVars = []
         self._add_vel_mag(climatology, derivedVars)
         self._add_thermal_forcing(climatology, derivedVars)
 
+        # then, call the superclass's version of this function so we extract
+        # the desired slices (but before renaming because it expects the
+        # original MPAS variable names)
+        climatology = super().customize_masked_climatology(climatology,
+                                                           season)
+        # finally, rename the variables and add metadata
         for varName, variable in self.variables.items():
             if varName not in derivedVars:
                 # rename variables from MPAS names to shorter names
@@ -295,10 +301,8 @@ def _add_vel_mag(self, climatology, derivedVars):
         """
         Add the velocity magnitude to the climatology if requested
         """
-        variables = self.variables
-
-        varName = 'verlocityMagnitude'
-        if varName not in variables:
+        varName = 'velocityMagnitude'
+        if varName not in self.variables:
             return
 
         derivedVars.append(varName)
@@ -311,10 +315,8 @@ def _add_thermal_forcing(self, climatology, derivedVars):
         """
         Add thermal forcing to the climatology if requested
         """
-        variables = self.variables
-
         varName = 'thermalForcing'
-        if varName not in variables:
+        if varName not in self.variables:
             return
 
         derivedVars.append(varName)
@@ -336,6 +338,12 @@ def _add_thermal_forcing(self, climatology, derivedVars):
         dp = cime_constants['SHR_CONST_G']*dens*thick
         press = dp.cumsum(dim='nVertLevels') - 0.5*dp
 
+        # add land ice pressure if available
+        ds_restart = xr.open_dataset(self.restartFileName)
+        ds_restart = ds_restart.isel(Time=0)
+        if 'landIcePressure' in ds_restart:
+            press += ds_restart.landIcePressure
+
         tempFreeze = c0 + cs*salin + cp*press + cps*press*salin
 
         climatology[varName] = temp - tempFreeze

mpas_analysis/ocean/remap_depth_slices_subtask.py

@@ -248,7 +248,9 @@ def customize_masked_climatology(self, climatology, season):
         interfaceIndices = self.dsSlice.interfaceIndices
         interfaceIndexMask = self.dsSlice.interfaceIndexMask
 
-        for variableName in self.variableList:
+        # iterate over all variables since some new ones may have been
+        # added by a subclass
+        for variableName in climatology.data_vars:
             if 'nVertLevels' in climatology[variableName].dims:
                 # mask only the values with the right vertical index
                 da = climatology[variableName].where(

suite/main.cfg

@@ -0,0 +1,19 @@
+[climatologyMapCustom]
+## options related to plotting climatology maps of any field at various depths
+## (if they include a depth dimension) without observatons for comparison
+
+# a list of fields top plot for each depth slice.  All supported fields are
+# listed above
+variables = [
+        'temperature',
+        'salinity',
+        'potentialDensity',
+        'thermalForcing',
+        'zonalVelocity',
+        'meridionalVelocity',
+        'velocityMagnitude',
+        'vertVelocity',
+        'vertDiff',
+        'vertVisc',
+        'mixedLayerDepth'
+    ]

suite/setup.py

@@ -146,6 +146,10 @@ def main():
             [config_from_job,
              os.path.join('..', '..', 'suite', f'{args.run}.cfg')])
 
+    if args.run.startswith('main_py'):
+        config_from_job = ' '.join(
+            [config_from_job, os.path.join('..', '..', 'suite', 'main.cfg')])
+
     if args.run not in ['main', 'ctrl']:
         try:
             os.makedirs(os.path.join(suite_path, args.run))