Rationalising mesh dimension names

[bjlittle]

🐛 Bug Report

Given two valid CF-UGRID netCDF mesh files, the first of which has the following netCDF dimensions:

dimensions:
        time = 4 ;
        dim0 = 221184 ;
        num_node = 221186 ;
        num_vertices = 4 ;

The second of which has the following netCDF dimensions:

dimensions:
        time = 22 ;
        bnds = 2 ;
        dim0 = 221184 ;
        num_node = 221186 ;
        num_vertices = 4 ;

I easily can use iris to load these files into as single cube list and then save to a single combined netCDF file as follows:

import iris
from iris.experimental.ugrid import PARSE_UGRID_ON_LOAD

with PARSE_UGRID_ON_LOAD.context():
    cubes = iris.load(fnames)

iris.save(cubes, "combined.nc")

Inspecting the dimensions of combined.nc gives:

dimensions:
        num_node = 221186 ;
        dim0 = 221184 ;
        time = 22 ;
        dynamics_face_N_nodes = 4 ;
        bnds = 2 ;
        num_node_0 = 221186 ;
        dim0_0 = 221184 ;
        time_0 = 4 ;
        num_node_1 = 221186 ;
        dim0_1 = 221184 ;

😢

Consequently, attempting to then load combined.nc results in the following traceback:

Traceback (most recent call last):
  File "<snip>/script.py", line 36, in <module>
    combined = iris.load(fname) 
               ^^^^^^^^^^^^^^^^
  File "<snip>/iris/lib/iris/__init__.py", line 326, in load
    return _load_collection(uris, constraints, callback).merged().cubes()
           ^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^
  File "<snip>/iris/lib/iris/__init__.py", line 294, in _load_collection
    result = _CubeFilterCollection.from_cubes(cubes, constraints)
             ^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^
  File "<snip>/iris/lib/iris/cube.py", line 97, in from_cubes
    for cube in cubes:
  File "<snip>/iris/lib/iris/__init__.py", line 275, in _generate_cubes
    for cube in iris.io.load_files(part_names, callback, constraints):
  File "<snip>/iris/lib/iris/io/__init__.py", line 219, in load_files
    for cube in handling_format_spec.handler(fnames, callback, constraints):
  File "<snip>/iris/lib/iris/fileformats/netcdf/loader.py", line 634, in load_cubes
    mesh_coords, mesh_dim = _build_mesh_coords(mesh, cf_var)
                            ^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^
  File "<snip>/iris/lib/iris/experimental/ugrid/load.py", line 489, in _build_mesh_coords
    mesh_dim = cf_var.dimensions.index(mesh_dim_name)
               ^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^
ValueError: tuple.index(x): x not in tuple

This issue is a direct result of the dimension name mangling performed by iris.save.

In this particular case, the mesh_dim_name is 'dim0', and the cf_var.dimensions are ('time_0', 'dim0_0'). Hence the traceback from cf_var.dimensions.index.

As a quick temporary workaround, I can monkey-patch iris.experimental.ugrid.load._build_mesh_coords with this weak-contract:

def monkey_patch(mesh, cf_var):
    """this is a workaround"""
    element_dimensions = {
        "node": mesh.node_dimension,
        "edge": mesh.edge_dimension,
        "face": mesh.face_dimension,
    }
    mesh_dim_name = element_dimensions[cf_var.location]

    for i, dim_name in enumerate(cf_var.dimensions):
        if dim_name.startswith(mesh_dim_name):
            mesh_dim = i
            break
    else:
        raise ValueError("oops!")

    mesh_coords = mesh.to_MeshCoords(location=cf_var.location)
    return mesh_coords, mesh_dim

However, clearly we should address this issue in a more robust way.

Note that, I can supply data on request to recreate this issue.