Add support for transects plotted on the native MPAS-Ocean mesh

docs/api.rst

@@ -268,6 +268,7 @@ Plotting
    plot_polar_comparison
    plot_global_comparison
    plot_1D
+   plot_vertical_section_comparison
    plot_vertical_section
    colormap.setup_colormap
    ticks.plot_xtick_format

docs/config/transects.rst

@@ -8,46 +8,60 @@ the comparison grid for each transect::
 
     # The approximate horizontal resolution (in km) of each transect.  Latitude/
     # longitude between observation points will be subsampled at this interval.
-    # Use 'obs' to indicate no subsampling.
-    horizontalResolution = obs
+    # Use 'obs' to indicate no subsampling. Use 'mpas' to indicate plotting of
+    # model data on the native grid, in which case comparison with observations
+    # will take place on the observation grid.
+    #horizontalResolution = mpas
+    #horizontalResolution = obs
+    horizontalResolution = 5
 
     # The name of the vertical comparison grid.  Valid values are 'mpas' for the
     # MPAS vertical grid, 'obs' to use the locations of observations or
-    # any other name if the vertical grid is defined by 'verticalComparisonGrid'
-    # verticalComparisonGridName = obs
-    verticalComparisonGridName = uniform_0_to_4000m_at_10m
+    # any other name if the vertical grid is defined by 'verticalComparisonGrid'.
+    # If horizontalResolution is 'mpas', model data (both main and control) will be
+    # plotted on the MPAS vertical grid, regardless of the comparison grid.
     #verticalComparisonGridName = mpas
+    #verticalComparisonGridName = obs
+    verticalComparisonGridName = uniform_0_to_4000m_at_10m
 
     # The vertical comparison grid if 'verticalComparisonGridName' is not 'mpas' or
     # 'obs'.  This should be numpy array of (typically negative) elevations (in m).
     verticalComparisonGrid = numpy.linspace(0, -4000, 401)
 
-The ``horizontalResolution`` of all transects can be either ``obs`` or a number
-of kilometers.  If ``obs``, model data are sampled at latitute and longitude
-points corresponding to WOCE stations.  It a number of kilometers is given,
-linear interpolation between observation points is performed with approximately
-the requested resolution.  The distance between stations is always divided into
-an integer number of segments of equal length so the resolution may be slightly
-above or below ``horizontalResolution``.
+    # A range for the y axis (if any)
+    verticalBounds = []
 
+The ``horizontalResolution`` of all transects can be ``obs``, ``mpas`` or a
+number of kilometers.  If ``obs``, model data are sampled at latitude and
+longitude points corresponding to the observations.  If the horizontal grid
+is ``mpas``, then the native MPAS-Ocean mesh is used for both the horizontal and
+vertical grids.  If a number of kilometers is given, linear interpolation
+between observation points is performed with approximately the requested
+resolution.  The distance between observation points is always divided into an
+integer number of segments of equal length so the resolution may be slightly
+above or below ``horizontalResolution``.
 
 The vertical grid is determined by two parameters,
 ``verticalComparisonGridName`` and ``verticalComparisonGrid``.  If
-``verticalComparisonGridName = mpas``, the MPAS-Ocean vertical coordinate will
-be interpolated horitontally from grid cell centers to the latitude and
-longitude of each point along the transect, and the observations will be
-interpolated vertically to the resulting grid.  If
-``verticalComparisonGridName = obs``, the vertical grid of the observations
-is used instead.  If ``verticalComparisonGridName`` is anything else, it is
-taken to be the name of a user-defined vertical grid (best to make it
-descriptive and unique, e.g. ``uniform_0_to_4000m_at_10m``) and
+``verticalComparisonGridName = mpas``, but ``horizontalResoltuion`` is not
+``mpas``, the MPAS-Ocean vertical coordinate will be interpolated horizontally
+from grid cell centers to the latitude and longitude of each point along the
+transect, and the observations will be interpolated vertically to the resulting
+grid.  If ``verticalComparisonGridName = obs``, the vertical grid of the
+observations is used instead.  If ``verticalComparisonGridName`` is anything
+else, it is taken to be the name of a user-defined vertical grid (best to make
+it descriptive and unique, e.g. ``uniform_0_to_4000m_at_10m``) and
 ``verticalComparisonGrid`` should be assigned a valid array of positive-up
 depth values (in the form of a python list or numpy array), e.g.::
 
     verticalComparisonGrid = numpy.linspace(0, -4000, 401)
 
 produces points between 0 and -4000 m sampled every 10 m.
 
+``verticalBounds`` is a list of minimum and maximum limits for the vertical axis
+of the transect.  The default is an empty list, which means ``matplotlib``
+selects the axis limits to encompass the full range of the vertical grid.
+
 .. note::
 
   Some types of transects (e.g. those produce with geojson files) do not have

mpas_analysis/default.cfg

@@ -2332,7 +2332,10 @@ seasons =  ['ANN']
 
 # The approximate horizontal resolution (in km) of each transect.  Latitude/
 # longitude between observation points will be subsampled at this interval.
-# Use 'obs' to indicate no subsampling.
+# Use 'obs' to indicate no subsampling. Use 'mpas' to indicate plotting of
+# model data on the native grid, in which case comparison with observations
+# will take place on the observation grid.
+#horizontalResolution = mpas
 #horizontalResolution = obs
 horizontalResolution = 5
 
@@ -2345,15 +2348,22 @@ horizontalBounds = {'WOCE_A21': [],
 
 # The name of the vertical comparison grid.  Valid values are 'mpas' for the
 # MPAS vertical grid, 'obs' to use the locations of observations or
-# any other name if the vertical grid is defined by 'verticalComparisonGrid'
+# any other name if the vertical grid is defined by 'verticalComparisonGrid'.
+# If horizontalResolution is 'mpas', model data (both main and control) will be
+# plotted on the MPAS vertical grid, regardless of the comparison grid.
 #verticalComparisonGridName = mpas
 #verticalComparisonGridName = obs
 verticalComparisonGridName = uniform_0_to_4000m_at_10m
 
 # The vertical comparison grid if 'verticalComparisonGridName' is not 'mpas' or
 # 'obs'.  This should be numpy array of (typically negative) elevations (in m).
+# The first and last entries are used as axis bounds for 'mpas' and 'obs'
+# vertical comparison grids
 verticalComparisonGrid = numpy.linspace(0, -4000, 401)
 
+# A range for the y axis (if any)
+verticalBounds = []
+
 # The minimum weight of a destination cell after remapping. Any cell with
 # weights lower than this threshold will therefore be masked out.
 renormalizationThreshold = 0.01
@@ -2471,7 +2481,6 @@ normArgsDifference = {'vmin': -0.3, 'vmax': 0.3}
 contourLevelsDifference = []
 
 
-
 [geojsonTransects]
 ## options related to plotting model transects at points determined by a
 ## geojson file.  To generate your own geojson file, go to:
@@ -2524,7 +2533,9 @@ seasons =  ['ANN']
 
 # The approximate horizontal resolution (in km) of each transect.  Latitude/
 # longitude between observation points will be subsampled at this interval.
-# Use 'obs' to indicate no subsampling.
+# Use 'obs' to indicate no subsampling. Use 'mpas' to indicate plotting of
+# model data on the native grid.
+#horizontalResolution = mpas
 #horizontalResolution = obs
 horizontalResolution = 5
 
@@ -2536,8 +2547,13 @@ verticalComparisonGridName = uniform_0_to_4000m_at_10m
 
 # The vertical comparison grid if 'verticalComparisonGridName' is not 'mpas'.
 # This should be numpy array of (typically negative) elevations (in m).
+# The first and last entries are used as axis bounds for 'mpas' vertical
+# comparison grids
 verticalComparisonGrid = numpy.linspace(0, -4000, 401)
 
+# A range for the y axis (if any)
+verticalBounds = []
+
 # The minimum weight of a destination cell after remapping. Any cell with
 # weights lower than this threshold will therefore be masked out.
 renormalizationThreshold = 0.01
@@ -2713,7 +2729,10 @@ seasons =  ['ANN', 'JFM', 'JAS']
 
 # The approximate horizontal resolution (in km) of each transect.  Latitude/
 # longitude between observation points will be subsampled at this interval.
-# Use 'obs' to indicate no subsampling.
+# Use 'obs' to indicate no subsampling. Use 'mpas' to indicate plotting of
+# model data on the native grid, in which case comparison with observations
+# will take place on the observation grid.
+#horizontalResolution = mpas
 #horizontalResolution = obs
 horizontalResolution = 5
 
@@ -2726,8 +2745,13 @@ verticalComparisonGridName = uniform_10_to_1500m_at_10m
 
 # The vertical comparison grid if 'verticalComparisonGridName' is not 'mpas' or
 # 'obs'.  This should be numpy array of (typically negative) elevations (in m).
+# The first and last entries are used as axis bounds for 'mpas' and 'obs'
+# vertical comparison grids
 verticalComparisonGrid = numpy.linspace(-10, -1500, 150)
 
+# A range for the y axis (if any)
+verticalBounds = []
+
 # The minimum weight of a destination cell after remapping. Any cell with
 # weights lower than this threshold will therefore be masked out.
 renormalizationThreshold = 0.01

mpas_analysis/ocean/climatology_map_ohc_anomaly.py

@@ -349,7 +349,7 @@ def _compute_ohc(self, climatology):  # {{{
 
         nVertLevels = dsRestart.sizes['nVertLevels']
 
-        zMid = compute_zmid(dsRestart.bottomDepth, dsRestart.maxLevelCell,
+        zMid = compute_zmid(dsRestart.bottomDepth, dsRestart.maxLevelCell-1,
                             dsRestart.layerThickness)
 
         vertIndex = xarray.DataArray.from_dict(