Merge pull request #1066 from xylar/add-ocean-region-names-coord

Add `regionNames` coordinate to NetCDF output with "standard" region names

Author: xylar

docs/developers_guide/api.rst

@@ -94,6 +94,9 @@ Ocean tasks
    WoaTransects
    WoceTransects
 
+Ocean subtasks
+--------------
+
 .. currentmodule:: mpas_analysis.ocean.compute_anomaly_subtask
 
 .. autosummary::
@@ -116,6 +119,19 @@ Ocean tasks
    PlotHovmollerSubtask
 
 
+Ocean utilities
+---------------
+
+.. currentmodule:: mpas_analysis.ocean.utility
+
+.. autosummary::
+   :toctree: generated/
+
+   add_standard_regions_and_subset
+   get_standard_region_names
+   compute_zmid
+
+
 Sea ice tasks
 -------------
 

docs/users_guide/config/regions.rst

@@ -15,10 +15,10 @@ within MPAS-Analysis using region mask files::
 
   # list of region names (needs to be in the same order as region indices in
   # time-series stats)
-  regions = ['arctic', 'equatorial', 'so', 'nino3', 'nino4', 'nino3.4', 'global']
-  # list of plot titles (needs to be in the same order as region indices in
-  # time-series stats)
-  plotTitles = ['Arctic', 'Equatorial (15S-15N)', 'Southern Ocean', 'Nino 3',
+  regionShortNames = ['arctic', 'equatorial', 'so', 'nino3', 'nino4', 'nino3.4',
+                      'global']
+  # list of full names (e.g. for plot titles) same order as regionShortNames
+  regionNames = ['Arctic', 'Equatorial (15S-15N)', 'Southern Ocean', 'Nino 3',
                 'Nino 4', 'Nino 3.4', 'Global Ocean']
 
 

docs/users_guide/tasks/timeSeriesOHCAnomaly.rst

@@ -21,8 +21,8 @@ The following configuration options are available for this task::
   ## options related to plotting time series of ocean heat content (OHC)
   ## anomalies from year 1
 
-  # list of regions to plot from the region list in [regions] below
-  regions = ['global']
+  # list of region shrot names to plot from the region list in [regions] above
+  regionShortNames = ['global']
 
   # approximate depths (m) separating plots of the upper, middle and lower ocean
   depths = [700, 2000]

docs/users_guide/tasks/timeSeriesSST.rst

@@ -20,8 +20,8 @@ The following configuration options are available for this task::
   [timeSeriesSST]
   ## options related to plotting time series of sea surface temperature (SST)
 
-  # list of regions to plot from the region list in [regions] below
-  regions = ['global']
+  # list of region shrot names to plot from the region list in [regions] above
+  regionShortNames = ['global']
 
   # Number of points over which to compute moving average (e.g., for monthly
   # output, movingAveragePoints=12 corresponds to a 12-month moving average

docs/users_guide/tasks/timeSeriesSalinityAnomaly.rst

@@ -20,8 +20,8 @@ The following configuration options are available for this task::
   [hovmollerSalinityAnomaly]
   ## options related to plotting time series of salinity vs. depth
 
-  # list of regions to plot from the region list in [regions] below
-  regions = ['global']
+  # list of region shrot names to plot from the region list in [regions] above
+  regionShortNames = ['global']
 
   # Number of points over which to compute moving average(e.g., for monthly
   # output, movingAveragePoints=12 corresponds to a 12-month moving average

docs/users_guide/tasks/timeSeriesTemperatureAnomaly.rst

@@ -20,8 +20,8 @@ The following configuration options are available for this task::
   [hovmollerTemperatureAnomaly]
   ## options related to plotting time series of potential temperature vs. depth
 
-  # list of regions to plot from the region list in [regions] below
-  regions = ['global']
+  # list of region shrot names to plot from the region list in [regions] above
+  regionShortNames = ['global']
 
   # Number of points over which to compute moving average(e.g., for monthly
   # output, movingAveragePoints=12 corresponds to a 12-month moving average

mpas_analysis/default.cfg

@@ -358,11 +358,11 @@ endYear = 20
 
 # list of region names (needs to be in the same order as region indices in
 # time-series stats)
-regions = ['arctic', 'equatorial', 'so', 'nino3', 'nino4', 'nino3.4', 'global']
-# list of plot titles (needs to be in the same order as region indices in
-# time-series stats)
-plotTitles = ['Arctic', 'Equatorial (15S-15N)', 'Southern Ocean', 'Nino 3',
-              'Nino 4', 'Nino 3.4', 'Global Ocean']
+regionShortNames = ['arctic', 'equatorial', 'so', 'nino3', 'nino4', 'nino3.4',
+                    'global']
+# list of full names (e.g. for plot titles) same order as regionShortNames
+regionNames = ['Arctic', 'Equatorial (15S-15N)', 'Southern Ocean', 'Nino 3',
+               'Nino 4', 'Nino 3.4', 'Global Ocean']
 
 
 [plot]
@@ -608,8 +608,8 @@ concentrationAltibergSH = Altiberg/Altiberg_1991-2017_20180308.nc
 ## options related to plotting time series of ocean heat content (OHC)
 ## anomalies from year 1
 
-# list of regions to plot from the region list in [regions] below
-regions = ['global']
+# list of region short names to plot from the region list in [regions] above
+regionShortNames = ['global']
 
 # approximate depths (m) separating plots of the upper, middle and lower ocean
 depths = [700, 2000]
@@ -672,8 +672,8 @@ movingAveragePoints = 12
 [hovmollerTemperatureAnomaly]
 ## options related to plotting time series of potential temperature vs. depth
 
-# list of regions to plot from the region list in [regions] below
-regions = ['global']
+# list of region short names to plot from the region list in [regions] above
+regionShortNames = ['global']
 
 # Number of points over which to compute moving average(e.g., for monthly
 # output, movingAveragePoints=12 corresponds to a 12-month moving average
@@ -706,8 +706,8 @@ contourLevelsResult = np.arange(-5.0, 5.51, 0.5)
 [hovmollerSalinityAnomaly]
 ## options related to plotting time series of salinity vs. depth
 
-# list of regions to plot from the region list in [regions] below
-regions = ['global']
+# list of region short names to plot from the region list in [regions] above
+regionShortNames = ['global']
 
 # Number of points over which to compute moving average(e.g., for monthly
 # output, movingAveragePoints=12 corresponds to a 12-month moving average
@@ -743,7 +743,7 @@ contourLevelsResult = np.arange(-2., 2.01, 0.2)
 ## options related to plotting time series of sea surface temperature (SST)
 
 # list of regions to plot from the region list in [regions] below
-regions = ['global']
+regionShortNames = ['global']
 
 # Number of points over which to compute moving average (e.g., for monthly
 # output, movingAveragePoints=12 corresponds to a 12-month moving average

mpas_analysis/ocean/compute_anomaly_subtask.py

@@ -28,6 +28,8 @@
 from mpas_analysis.shared.time_series import \
     compute_moving_avg_anomaly_from_start
 
+from mpas_analysis.ocean.utility import add_standard_regions_and_subset
+
 
 class ComputeAnomalySubtask(AnalysisTask):
     """
@@ -184,6 +186,9 @@ def run_task(self):
             movingAveragePoints=self.movingAveragePoints,
             alter_dataset=self.alter_dataset)
 
+        if 'nOceanRegions' in ds.dims or 'nOceanRegionsTmp' in ds.dims:
+            ds = add_standard_regions_and_subset(ds, config)
+
         outFileName = self.outFileName
         if not os.path.isabs(outFileName):
             baseDirectory = build_config_full_path(

mpas_analysis/ocean/index_nino34.py

@@ -38,6 +38,8 @@
 from mpas_analysis.shared import AnalysisTask
 from mpas_analysis.shared.html import write_image_xml
 
+from mpas_analysis.ocean.utility import add_standard_regions_and_subset
+
 
 class IndexNino34(AnalysisTask):
     """
@@ -129,7 +131,7 @@ def setup_and_check(self):
         regionToPlot = config.get('indexNino34', 'region')
 
         if regionToPlot not in ['nino3.4', 'nino3', 'nino4']:
-            raise ValueError('Unexpectes El Nino Index region {}'.format(
+            raise ValueError('Unexpected El Nino Index region {}'.format(
                 regionToPlot))
         ninoIndexNumber = regionToPlot[4:]
 
@@ -193,9 +195,7 @@ def run_task(self):
 
         # regionIndex should correspond to NINO34 in surface weighted Average
         # AM
-        regions = config.getexpression('regions', 'regions')
         regionToPlot = config.get('indexNino34', 'region')
-        regionIndex = regions.index(regionToPlot)
 
         # Load data:
         ds = open_mpas_dataset(fileName=self.inputFile,
@@ -204,6 +204,11 @@ def run_task(self):
                                startDate=startDate,
                                endDate=endDate)
 
+        ds = add_standard_regions_and_subset(ds, config,
+                                             regionShortNames=[regionToPlot])
+        # we want to collapse the nOceanRegions dimension
+        ds = ds.isel(nOceanRegions=0)
+
         # Observations have been processed to the nino34Index prior to reading
         dsObs = xr.open_dataset(dataPath, decode_cf=False, decode_times=False)
         # add the days between 0001-01-01 and the refDate so we have a new
@@ -215,7 +220,7 @@ def run_task(self):
         self.logger.info('  Compute El Nino {} Index...'.format(
             ninoIndexNumber))
         varName = self.variableList[0]
-        regionSST = ds[varName].isel(nOceanRegions=regionIndex)
+        regionSST = ds[varName]
         nino34Main = self._compute_nino34_index(regionSST, calendar)
 
         # Compute the observational index over the entire time range
@@ -262,7 +267,10 @@ def run_task(self):
                 calendar=calendar,
                 variableList=self.variableList)
 
-            regionSSTRef = dsRef[varName].isel(nOceanRegions=regionIndex)
+            dsRef = add_standard_regions_and_subset(
+                dsRef, self.controlConfig, regionShortNames=[regionToPlot])
+
+            regionSSTRef = dsRef[varName]
             nino34Ref = self._compute_nino34_index(regionSSTRef, calendar)
 
             nino34s = [nino34Subset, nino34Main[2:-3], nino34Ref[2:-3]]

mpas_analysis/ocean/plot_depth_integrated_time_series_subtask.py

@@ -31,6 +31,8 @@
 from mpas_analysis.shared.time_series import compute_moving_avg, \
     combine_time_series_with_ncrcat
 
+from mpas_analysis.ocean.utility import add_standard_regions_and_subset
+
 
 class PlotDepthIntegratedTimeSeriesSubtask(AnalysisTask):
     """
@@ -157,7 +159,8 @@ def __init__(self, parentTask, regionName, inFileName, outFileLabel,
 
         if subtaskName is None:
             suffix = regionName[0].upper() + regionName[1:]
-            subtaskName = 'plotDepthIntegratedTimeSeries{}'.format(suffix)
+            suffix = suffix.replace(' ', '_')
+            subtaskName = f'plotDepthIntegratedTimeSeries{suffix}'
 
         # first, call the constructor from the base class (AnalysisTask)
         super(PlotDepthIntegratedTimeSeriesSubtask, self).__init__(
@@ -205,7 +208,7 @@ def setup_and_check(self):
         if self.controlConfig is not None:
             # we need to know what file to read from the control run so
             # an absolute path won't work
-            assert(not os.path.isabs(self.inFileName))
+            assert not os.path.isabs(self.inFileName)
 
             baseDirectory = build_config_full_path(
                 self.controlConfig, 'output', 'timeSeriesSubdirectory')
@@ -217,7 +220,7 @@ def setup_and_check(self):
             'runs', 'preprocessedReferenceRunName')
         if preprocessedReferenceRunName != 'None':
 
-            assert(not os.path.isabs(self.inFileName))
+            assert not os.path.isabs(self.inFileName)
 
             baseDirectory = build_config_full_path(
                 config, 'output', 'timeSeriesSubdirectory')
@@ -263,11 +266,6 @@ def run_task(self):
 
         mainRunName = config.get('runs', 'mainRunName')
 
-        plotTitles = config.getexpression('regions', 'plotTitles')
-        allRegionNames = config.getexpression('regions', 'regions')
-        regionIndex = allRegionNames.index(self.regionName)
-        regionNameInTitle = plotTitles[regionIndex]
-
         startDate = config.get('timeSeries', 'startDate')
         endDate = config.get('timeSeries', 'endDate')
 
@@ -279,7 +277,17 @@ def run_task(self):
                                timeVariableNames=None,
                                startDate=startDate,
                                endDate=endDate)
-        ds = ds.isel(nOceanRegionsTmp=regionIndex)
+        if 'regionNames' in ds.coords:
+            # we added region names already
+            ds = ds.set_xindex('regionNames')
+            ds = ds.sel(regionNames=self.regionName)
+            regionNameInTitle = ds.regionNames.values
+        else:
+            # we need to add region names and select the right region by short
+            # name
+            ds = add_standard_regions_and_subset(
+                ds, config, regionShortNames=[self.regionName])
+            regionNameInTitle = ds.regionNames.values[0]
 
         depths = ds.depth.values
 
@@ -436,7 +444,17 @@ def run_task(self):
                                       timeVariableNames=None,
                                       startDate=controlStartDate,
                                       endDate=controlEndDate)
-            dsRef = dsRef.isel(nOceanRegionsTmp=regionIndex)
+
+            if 'regionNames' in dsRef.coords:
+                # we added region names already
+                dsRef = dsRef.set_xindex('regionNames')
+                dsRef = dsRef.sel(regionNames=self.regionName)
+            else:
+                # we need to add region names and select the right region by
+                # short name
+                dsRef = add_standard_regions_and_subset(
+                    dsRef,  self.controlConfig,
+                    regionShortNames=[self.regionName])
 
             color = config.get('timeSeries', 'controlColor')
 
@@ -468,8 +486,9 @@ def run_task(self):
 
         fig = timeseries_analysis_plot(
             config=config, dsvalues=timeSeries, calendar=calendar,
-            title=title, xlabel=xLabel, ylabel=yLabel, movingAveragePoints=None,
-            lineColors=lineColors, lineStyles=lineStyles, markers=lineMarkers,
+            title=title, xlabel=xLabel, ylabel=yLabel,
+            movingAveragePoints=None, lineColors=lineColors,
+            lineStyles=lineStyles, markers=lineMarkers,
             lineWidths=lineWidths, legendText=legendText, maxPoints=maxPoints,
             firstYearXTicks=firstYearXTicks, yearStrideXTicks=yearStrideXTicks)