Merge pull request #662 from xylar/fix-write-netcdf-without-time

Remove encoding unlimited_dims if Time not present

Author: xylar

conda_package/mpas_tools/io.py

@@ -135,6 +135,8 @@ def write_netcdf(
         # make sure the Time dimension is unlimited because MPAS has trouble
         # reading Time otherwise
         ds.encoding['unlimited_dims'] = {'Time'}
+    else:
+        ds.encoding['unlimited_dims'] = None
 
     # for performance, we have to handle this as a special case
     convert = format == 'NETCDF3_64BIT_DATA'

conda_package/mpas_tools/mesh/creation/sort_mesh.py

@@ -1,11 +1,13 @@
+import argparse
+import os
 
 import numpy as np
-import os
 import xarray
-import argparse
 from scipy.sparse import csr_matrix
 from scipy.sparse.csgraph import reverse_cuthill_mckee
 
+from mpas_tools.io import write_netcdf
+
 
 def sort_mesh(mesh):
     """
@@ -24,9 +26,9 @@ def sort_mesh(mesh):
     """
     # Authors: Darren Engwirda
 
-    ncells = mesh.sizes["nCells"]
-    nedges = mesh.sizes["nEdges"]
-    nduals = mesh.sizes["nVertices"]
+    ncells = mesh.sizes['nCells']
+    nedges = mesh.sizes['nEdges']
+    nduals = mesh.sizes['nVertices']
 
     cell_fwd = np.arange(0, ncells) + 1
     cell_rev = np.arange(0, ncells) + 1
@@ -42,23 +44,20 @@ def sort_mesh(mesh):
     cell_rev = np.zeros(ncells, dtype=np.int32)
     cell_rev[cell_fwd - 1] = np.arange(ncells) + 1
 
-    mesh["cellsOnCell"][:] = \
-        _sort_rev(mesh["cellsOnCell"], cell_rev)
-    mesh["cellsOnEdge"][:] = \
-        _sort_rev(mesh["cellsOnEdge"], cell_rev)
-    mesh["cellsOnVertex"][:] = \
-        _sort_rev(mesh["cellsOnVertex"], cell_rev)
+    mesh['cellsOnCell'][:] = _sort_rev(mesh['cellsOnCell'], cell_rev)
+    mesh['cellsOnEdge'][:] = _sort_rev(mesh['cellsOnEdge'], cell_rev)
+    mesh['cellsOnVertex'][:] = _sort_rev(mesh['cellsOnVertex'], cell_rev)
 
     for var in mesh.keys():
         dims = mesh.variables[var].dims
-        if ("nCells" in dims):
+        if 'nCells' in dims:
             mesh[var][:] = _sort_fwd(mesh[var], cell_fwd)
 
-    mesh["indexToCellID"][:] = np.arange(ncells) + 1
+    mesh['indexToCellID'][:] = np.arange(ncells) + 1
 
     # sort duals via pseudo-linear cell-wise ordering
 
-    dual_fwd = np.ravel(mesh["verticesOnCell"].values)
+    dual_fwd = np.ravel(mesh['verticesOnCell'].values)
     dual_fwd = dual_fwd[dual_fwd > 0]
 
     __, imap = np.unique(dual_fwd, return_index=True)
@@ -68,22 +67,20 @@ def sort_mesh(mesh):
     dual_rev = np.zeros(nduals, dtype=np.int32)
     dual_rev[dual_fwd - 1] = np.arange(nduals) + 1
 
-    mesh["verticesOnCell"][:] = \
-        _sort_rev(mesh["verticesOnCell"], dual_rev)
+    mesh['verticesOnCell'][:] = _sort_rev(mesh['verticesOnCell'], dual_rev)
 
-    mesh["verticesOnEdge"][:] = \
-        _sort_rev(mesh["verticesOnEdge"], dual_rev)
+    mesh['verticesOnEdge'][:] = _sort_rev(mesh['verticesOnEdge'], dual_rev)
 
     for var in mesh.keys():
         dims = mesh.variables[var].dims
-        if ("nVertices" in dims):
+        if 'nVertices' in dims:
             mesh[var][:] = _sort_fwd(mesh[var], dual_fwd)
 
-    mesh["indexToVertexID"][:] = np.arange(nduals) + 1
+    mesh['indexToVertexID'][:] = np.arange(nduals) + 1
 
     # sort edges via pseudo-linear cell-wise ordering
 
-    edge_fwd = np.ravel(mesh["edgesOnCell"].values)
+    edge_fwd = np.ravel(mesh['edgesOnCell'].values)
     edge_fwd = edge_fwd[edge_fwd > 0]
 
     __, imap = np.unique(edge_fwd, return_index=True)
@@ -93,60 +90,75 @@ def sort_mesh(mesh):
     edge_rev = np.zeros(nedges, dtype=np.int32)
     edge_rev[edge_fwd - 1] = np.arange(nedges) + 1
 
-    mesh["edgesOnCell"][:] = \
-        _sort_rev(mesh["edgesOnCell"], edge_rev)
+    mesh['edgesOnCell'][:] = _sort_rev(mesh['edgesOnCell'], edge_rev)
 
-    mesh["edgesOnEdge"][:] = \
-        _sort_rev(mesh["edgesOnEdge"], edge_rev)
+    mesh['edgesOnEdge'][:] = _sort_rev(mesh['edgesOnEdge'], edge_rev)
 
-    mesh["edgesOnVertex"][:] = \
-        _sort_rev(mesh["edgesOnVertex"], edge_rev)
+    mesh['edgesOnVertex'][:] = _sort_rev(mesh['edgesOnVertex'], edge_rev)
 
     for var in mesh.keys():
         dims = mesh.variables[var].dims
-        if ("nEdges" in dims):
+        if 'nEdges' in dims:
             mesh[var][:] = _sort_fwd(mesh[var], edge_fwd)
 
-    mesh["indexToEdgeID"][:] = np.arange(nedges) + 1
+    mesh['indexToEdgeID'][:] = np.arange(nedges) + 1
 
     return mesh
 
 
 def main():
     parser = argparse.ArgumentParser(
-        description=__doc__,
-        formatter_class=argparse.RawTextHelpFormatter)
+        description=__doc__, formatter_class=argparse.RawTextHelpFormatter
+    )
+
+    parser.add_argument(
+        '--mesh-file',
+        dest='mesh_file',
+        type=str,
+        required=True,
+        help='Path+name to unsorted mesh file.',
+    )
 
     parser.add_argument(
-        "--mesh-file", dest="mesh_file", type=str,
-        required=True, help="Path+name to unsorted mesh file.")
+        '--sort-file',
+        dest='sort_file',
+        type=str,
+        required=True,
+        help='Path+name to sorted output file.',
+    )
 
     parser.add_argument(
-        "--sort-file", dest="sort_file", type=str,
-        required=True, help="Path+name to sorted output file.")
+        '--format',
+        dest='format',
+        type=str,
+        required=False,
+        default='NETCDF4',
+        help='Output format for the sorted mesh file.',
+    )
 
     args = parser.parse_args()
 
     mesh = xarray.open_dataset(args.mesh_file)
 
     sort_mesh(mesh)
 
-    with open(os.path.join(os.path.dirname(
-              args.sort_file), "graph.info"), "w") as fptr:
-        cellsOnCell = mesh["cellsOnCell"].values
+    with open(
+        os.path.join(os.path.dirname(args.sort_file), 'graph.info'), 'w'
+    ) as fptr:
+        cellsOnCell = mesh['cellsOnCell'].values
 
-        ncells = mesh.sizes["nCells"]
+        ncells = mesh.sizes['nCells']
         nedges = np.count_nonzero(cellsOnCell) // 2
 
-        fptr.write(f"{ncells} {nedges}\n")
+        fptr.write(f'{ncells} {nedges}\n')
         for cell in range(ncells):
             data = cellsOnCell[cell, :]
             data = data[data > 0]
             for item in data:
-                fptr.write(f"{item} ")
-            fptr.write("\n")
+                fptr.write(f'{item} ')
+            fptr.write('\n')
 
-    mesh.to_netcdf(args.sort_file, format="NETCDF4")
+    write_netcdf(mesh, args.sort_file, format=args.format)
 
 
 def _sort_fwd(data, fwd):
@@ -212,11 +224,10 @@ def _cell_del2(mesh):
     ivec = np.array([], dtype=np.int32)
     jvec = np.array([], dtype=np.int32)
 
-    topolOnCell = mesh["nEdgesOnCell"].values
-    cellsOnCell = mesh["cellsOnCell"].values
+    topolOnCell = mesh['nEdgesOnCell'].values
+    cellsOnCell = mesh['cellsOnCell'].values
 
     for edge in range(np.max(topolOnCell)):
-
         # cell-to-cell pairs, if edges exist
         mask = topolOnCell > edge
         idx_self = np.argwhere(mask).ravel()
@@ -237,5 +248,5 @@ def _cell_del2(mesh):
     return csr_matrix((xvec, (ivec, jvec)))
 
 
-if (__name__ == "__main__"):\
+if __name__ == '__main__':
     main()