Merge pull request #906 from xylar/remap_from_vertices

Add support for remapping from MPAS vertices

mpas_analysis/shared/climatology/climatology.py

@@ -35,7 +35,7 @@
 
 
 def get_remapper(config, sourceDescriptor, comparisonDescriptor,
-                 mappingFilePrefix, method, logger=None):
+                 mappingFilePrefix, method, logger=None, vertices=False):
     """
     Given config options and descriptions of the source and comparison grids,
     returns a ``pyremap.Remapper`` object that can be used to remap from source
@@ -49,10 +49,10 @@ def get_remapper(config, sourceDescriptor, comparisonDescriptor,
     config : mpas_tools.config.MpasConfigParser
         Contains configuration options
 
-    sourceDescriptor : ``MeshDescriptor`` subclass object
+    sourceDescriptor : pyremap.MeshDescriptor
         A description of the source mesh or grid
 
-    comparisonDescriptor : ``MeshDescriptor`` subclass object
+    comparisonDescriptor : pyremap.MeshDescriptor
         A description of the comparison grid
 
     mappingFilePrefix : str
@@ -61,12 +61,15 @@ def get_remapper(config, sourceDescriptor, comparisonDescriptor,
     method : {'bilinear', 'neareststod', 'conserve'}
         The method of interpolation used.
 
-    logger : ``logging.Logger``, optional
+    logger : logging.Logger, optional
         A logger to which ncclimo output should be redirected
 
+    vertices : bool, optional
+        Whether to remap from vertices, rather than cells
+
     Returns
     -------
-    remapper : ``pyremap.Remapper`` object
+    remapper : pyremap.Remapper
         A remapper that can be used to remap files or data sets from the source
         grid or mesh to the comparison grid.
     """
@@ -79,11 +82,14 @@ def get_remapper(config, sourceDescriptor, comparisonDescriptor,
     if not _matches_comparison(sourceDescriptor, comparisonDescriptor):
         # we need to remap because the grids don't match
 
-        mappingBaseName = '{}_{}_to_{}_{}.nc'.format(
-            mappingFilePrefix,
-            sourceDescriptor.meshName,
-            comparisonDescriptor.meshName,
-            method)
+        if vertices:
+            srcMeshName = f'{sourceDescriptor.meshName}_vertices'
+        else:
+            srcMeshName = sourceDescriptor.meshName
+        destMeshName = comparisonDescriptor.meshName
+
+        mappingBaseName = \
+            f'{mappingFilePrefix}_{srcMeshName}_to_{destMeshName}_{method}.nc'
 
         tryCustom = config.get('diagnostics', 'customDirectory') != 'none'
         if tryCustom:
@@ -92,17 +98,16 @@ def get_remapper(config, sourceDescriptor, comparisonDescriptor,
                 config, 'diagnostics', 'mappingSubdirectory',
                 baseDirectoryOption='customDirectory')
 
-            mappingFileName = '{}/{}'.format(mappingSubdirectory,
-                                             mappingBaseName)
+            mappingFileName = f'{mappingSubdirectory}/{mappingBaseName}'
+
         if not tryCustom or not os.path.exists(mappingFileName):
             # second see if mapping files are in the base directory
 
             mappingSubdirectory = build_config_full_path(
                 config, 'diagnostics', 'mappingSubdirectory',
                 baseDirectoryOption='base_path')
 
-            mappingFileName = '{}/{}'.format(mappingSubdirectory,
-                                             mappingBaseName)
+            mappingFileName = f'{mappingSubdirectory}/{mappingBaseName}'
 
         if not os.path.exists(mappingFileName):
             # we don't have a mapping file yet, so get ready to create one
@@ -111,8 +116,7 @@ def get_remapper(config, sourceDescriptor, comparisonDescriptor,
                 build_config_full_path(config, 'output',
                                        'mappingSubdirectory')
             make_directories(mappingSubdirectory)
-            mappingFileName = '{}/{}'.format(mappingSubdirectory,
-                                             mappingBaseName)
+            mappingFileName = f'{mappingSubdirectory}/{mappingBaseName}'
 
     remapper = Remapper(sourceDescriptor, comparisonDescriptor,
                         mappingFileName)

mpas_analysis/shared/climatology/remap_mpas_climatology_subtask.py

@@ -72,6 +72,8 @@ class RemapMpasClimatologySubtask(AnalysisTask):
     op : {'avg', 'min', 'max'}
          operator for monthly stats
 
+    vertices : bool
+        Whether to remap from vertices, rather than cells
     """
 
     # Authors
@@ -81,7 +83,7 @@ class RemapMpasClimatologySubtask(AnalysisTask):
     def __init__(self, mpasClimatologyTask, parentTask, climatologyName,
                  variableList, seasons, comparisonGridNames=None,
                  iselValues=None, subtaskName='remapMpasClimatology',
-                 useNcremap=None):
+                 useNcremap=None, vertices=False):
 
         """
         Construct the analysis task and adds it as a subtask of the
@@ -130,6 +132,9 @@ def __init__(self, mpasClimatologyTask, parentTask, climatologyName,
             if it is not explicitly given.  If a comparison grid other than
             ``latlon`` is given, ncremap is not supported so this flag is set
             to ``False``.
+
+        vertices : bool, optional
+            Whether to remap from vertices, rather than cells
         """
         # Authors
         # -------
@@ -174,6 +179,8 @@ def __init__(self, mpasClimatologyTask, parentTask, climatologyName,
         else:
             self.useNcremap = useNcremap
 
+        self.vertices = vertices
+
     def setup_and_check(self):
         """
         Perform steps to set up the analysis and check for errors in the setup.
@@ -405,17 +412,18 @@ def _setup_remappers(self):
             comparisonDescriptor = \
                 self.comparisonDescriptors[comparisonGridName]
             self.comparisonGridName = comparisonDescriptor.meshName
+            meshName = config.get('input', 'mpasMeshName')
             mpasDescriptor = MpasMeshDescriptor(
-                self.restartFileName, meshName=config.get('input',
-                                                          'mpasMeshName'))
+                self.restartFileName, meshName=meshName,
+                vertices=self.vertices)
             self.mpasMeshName = mpasDescriptor.meshName
 
             self.remappers[comparisonGridName] = get_remapper(
                 config=config, sourceDescriptor=mpasDescriptor,
                 comparisonDescriptor=comparisonDescriptor,
                 mappingFilePrefix=mappingFilePrefix,
                 method=config.get('climatology', 'mpasInterpolationMethod'),
-                logger=self.logger)
+                logger=self.logger, vertices=self.vertices)
 
     def _setup_file_names(self):
         """
@@ -523,10 +531,6 @@ def _mask_climatologies(self, season, dsMask):
             if len(iselValues.keys()) > 0:
                 climatology = climatology.isel(**iselValues)
 
-            # add valid mask as a variable, useful for remapping later
-            climatology['validMask'] = \
-                xr.DataArray(numpy.ones(climatology.sizes['nCells']),
-                             dims=['nCells'])
             # mask the data set
             for variableName in self.variableList:
                 climatology[variableName] = \
@@ -537,6 +541,15 @@ def _mask_climatologies(self, season, dsMask):
             climatology = self.customize_masked_climatology(climatology,
                                                             season)
 
+            if self.vertices:
+                dim = 'nVertices'
+            else:
+                dim = 'nCells'
+            # add valid mask as a variable, useful for remapping later
+            climatology['validMask'] = \
+                xr.DataArray(numpy.ones(climatology.sizes[dim]),
+                             dims=[dim])
+
             write_netcdf(climatology, maskedClimatologyFileName)
 
     def _remap(self, inFileName, outFileName, remapper, comparisonGridName,