Repository: NCAS-CMS/cfdm

Issue:
Failure from `cfdm.write` when writing identical (auxiliary) coordinates to different data variables in different groups 
It should be possible to write to disk fields with identical dimension or auxiliary coordinates to different groups, but this causes a failure from `cfdm.write`, which complains that the netCDF data variable dimensions are not in the same group or a parent group as the data variable itself. 

It is good that `cfdm.write` complains about this, as it is a valid restriction, but it shouldn't have ended up in that situation in the first place. 

Many thanks to @ajelenak who noticed and diagnosed this problem, and provided a minimal reproducer, which I slightly adapted here to also  cover auxiliary coordinates: 

```python
import numpy as np
import cfdm

# Define dimension and auxiliary coordinates
station_name = cfdm.AuxiliaryCoordinate()
station_name.set_data(["station1", "station2"])
station_name.set_properties(
    {"long_name": "sensor id", "cf_role": "timeseries_id"}
)
station_name.nc_set_variable("station_name")

times = cfdm.DimensionCoordinate()
times.set_data(range(5))
times.set_properties(
    {"standard_name": "time", "units": "days since 1967-04-01"}
)
times.nc_set_variable("time")

file_content = list()
for mt in ("Dew Point", "Wind Speed"):
    grp_name = mt.replace(" ", "_")
    var_name = mt.replace(" ", "_")

    # Define domain axes
    stations = cfdm.DomainAxis(2)
    stations.nc_set_dimension("stations")
    stations.nc_set_dimension_groups([grp_name])

    time = cfdm.DomainAxis(5)
    time.nc_set_dimension("time")
    time.nc_set_dimension_groups([grp_name])

    # Assign new groups to dimension and auxiliary coordinates
    station_name = station_name.copy()
    station_name.nc_set_variable_groups([grp_name])

    times = times.copy()
    times.nc_set_variable_groups([grp_name])

    # Define the field in the same group as the coordinates
    f = cfdm.Field()

    station_axis = f.set_construct(stations)
    time_axis = f.set_construct(time)

    f.set_construct(times, axes=time_axis)
    f.set_construct(station_name, axes=station_axis)

    f.set_data(
        np.arange(10).reshape(5, 2), axes=[time_axis, station_axis]
    )

    f.nc_set_variable(var_name)
    f.nc_set_variable_groups([grp_name])

    file_content.append(f)

# This should not raise an exception, but does at v1.9.0.2
cfdm.write(file_content, "result.nc")
```
produces:
```python
Traceback (most recent call last):
  File "/home/david/alek.py", line 59, in <module>
    cfdm.write(file_content, "result.nc")
  File "/home/david/miniconda3.9/lib/python3.9/site-packages/cfdm/read_write/write.py", line 522, in write
    netcdf.write(
  File "/home/david/miniconda3.9/lib/python3.9/site-packages/cfdm/decorators.py", line 171, in verbose_override_wrapper
    return method_with_verbose_kwarg(*args, **kwargs)
  File "/home/david/miniconda3.9/lib/python3.9/site-packages/cfdm/read_write/netcdf/netcdfwrite.py", line 4809, in write
    self._file_io_iteration(
  File "/home/david/miniconda3.9/lib/python3.9/site-packages/cfdm/read_write/netcdf/netcdfwrite.py", line 5089, in _file_io_iteration
    self._write_field_or_domain(f)
  File "/home/david/miniconda3.9/lib/python3.9/site-packages/cfdm/read_write/netcdf/netcdfwrite.py", line 3972, in _write_field_or_domain
    self._write_netcdf_variable(
  File "/home/david/miniconda3.9/lib/python3.9/site-packages/cfdm/read_write/netcdf/netcdfwrite.py", line 2696, in _write_netcdf_variable
    raise ValueError(
ValueError: Can't create netCDF variable '/Wind_Speed/Wind_Speed' from <Field: ncvar%/Wind_Speed/Wind_Speed(time(5), cf_role=timeseries_id(2))> with dimension '/Dew_Point/time' that is not in the same group or a sub-group as the variable.
```



```python
>>> cfdm.environment(paths=False)
Platform: Linux-5.4.0-100-generic-x86_64-with-glibc2.31
HDF5 library: 1.12.1
netcdf library: 4.8.1
Python: 3.9.5
netCDF4: 1.5.8
numpy: 1.22.2
cfdm.core: 1.9.0.2
cftime: 1.5.2
netcdf_flattener: 1.2.0
cfdm: 1.9.0.2
```

Discussion:

Patch:
diff --git a/Changelog.rst b/Changelog.rst
index 01e43ed93..9c6d6fbfc 100644
--- a/Changelog.rst
+++ b/Changelog.rst
@@ -1,8 +1,11 @@
-Version 1.9.?.?
+Version 1.9.0.3
 ---------------
 
-**2022-0?-??**
+**2022-03-??**
 
+* Fixed bug that caused a failure from `cfdm.write` when writing
+  identical (auxiliary) coordinates to different data variables in
+  different groups (https://github.com/NCAS-CMS/cfdm/issues/177)
 * Fixed bug that caused `cf.Domain.__str__` to fail when a dimension
   coordinate construct does not have data
   (https://github.com/NCAS-CMS/cfdm/issues/174)
diff --git a/cfdm/core/__init__.py b/cfdm/core/__init__.py
index a06a41c34..d87c30e08 100644
--- a/cfdm/core/__init__.py
+++ b/cfdm/core/__init__.py
@@ -1,8 +1,8 @@
 """"""
 
-__date__ = "2022-01-31"
+__date__ = "2022-03-??"
 __cf_version__ = "1.9"
-__version__ = "1.9.0.2"
+__version__ = "1.9.0.3"
 
 from distutils.version import LooseVersion
 import platform
diff --git a/cfdm/read_write/netcdf/netcdfwrite.py b/cfdm/read_write/netcdf/netcdfwrite.py
index d255ce1e6..0fe26540a 100644
--- a/cfdm/read_write/netcdf/netcdfwrite.py
+++ b/cfdm/read_write/netcdf/netcdfwrite.py
@@ -78,7 +78,6 @@ class NetCDFWrite(IOWrite):
         return nc.createGroup(group_name)
 
     def _create_netcdf_variable_name(self, parent, default):
-        #                            force_use_existing=False):
         """Create an appropriate name for a netCDF variable.
 
         .. versionadded:: (cfdm) 1.7.0
@@ -97,12 +96,6 @@ class NetCDFWrite(IOWrite):
         """
         ncvar = self.implementation.nc_get_variable(parent, None)
 
-        #        if force_use_existing:
-        #            if ncvar is None:
-        #                raise ValueError()
-        #
-        #            return ncvar
-
         if ncvar is None:
             try:
                 ncvar = self.implementation.get_property(
@@ -150,12 +143,6 @@ class NetCDFWrite(IOWrite):
             if not role:
                 raise ValueError("Must supply role when providing dimsize")
 
-            #            if base in g['dimensions_with_role'].get(role, ()):
-            #                if base in ncdim_names and dimsize == g['ncdim_to_size'][base]:
-            #                    # Return the name of an existing netCDF dimension
-            #                    # with this name, this size, and matching the
-            #                    # given role.
-            #                    return base
             for ncdim in g["dimensions_with_role"].get(role, ()):
                 if g["ncdim_to_size"][ncdim] == dimsize:
                     # Return the name of an existing netCDF dimension
@@ -340,23 +327,7 @@ class NetCDFWrite(IOWrite):
 
         array = numpy.array(tuple(array.tobytes().decode("ascii")), dtype="S1")
 
-        #        else:
-        #            # dtype is 'U'
-        #            x = []
-        #            for s in array.flatten():
-        #                x.extend(tuple(s.ljust(N, '\x00')))
-        #
-        #            array = numpy.array(k, dtype='S1')
-
         array.resize(original_shape + (array.size // original_size,))
-        #        if masked:
-        #            array = numpy.ma.array(array, mask=mask, fill_value=fill_value)
-
-        #        if array.dtype.kind == 'U':
-        #            # Convert unicode to string
-        #            array = array.astype('S')
-        #
-        #        new = netCDF4.stringtochar(array, encoding='none')
 
         if masked:
             array = numpy.ma.masked_where(array == "", array)
@@ -365,17 +336,6 @@ class NetCDFWrite(IOWrite):
         if array.dtype.kind != "S":
             raise ValueError("Array must have string data type.")
 
-        #            new = numpy.ma.array(new, mask=mask, fill_value=fill_value)
-
-        #        new = numpy.ma.masked_all(shape + (strlen,), dtype='S1')
-        #
-        #        for index in numpy.ndindex(shape):
-        #            value = array[index]
-        #            if value is numpy.ma.masked:
-        #                new[index] = numpy.ma.masked
-        #            else:
-        #                new[index] = tuple(value.ljust(strlen, ' '))
-
         return array
 
     def _datatype(self, variable):
@@ -487,6 +447,10 @@ class NetCDFWrite(IOWrite):
             field: Field construct
 
             key: `str`
+                The construct identifier of the metadata construct.
+
+                *Parameter example:*
+                  ``'auxiliarycoordinate1'``
 
         :Returns:
 
@@ -678,7 +642,7 @@ class NetCDFWrite(IOWrite):
             coordinates: `list`
                This list may get updated in-place.
 
-               .. versionadded:: (cfdm) .8.7.0
+               .. versionadded:: (cfdm) 1.8.7.0
 
         :Returns:
 
@@ -687,14 +651,15 @@ class NetCDFWrite(IOWrite):
 
         """
         g = self.write_vars
-
         seen = g["seen"]
 
         data_axes = self.implementation.get_construct_data_axes(f, key)
         axis = data_axes[0]
 
+        already_in_file = self._already_in_file(coord)
+
         create = False
-        if not self._already_in_file(coord):
+        if not already_in_file:
             create = True
         elif seen[id(coord)]["ncdims"] != ():
             if seen[id(coord)]["ncvar"] != seen[id(coord)]["ncdims"][0]:
@@ -703,11 +668,20 @@ class NetCDFWrite(IOWrite):
                 # coordinate.
                 create = True
 
+        # If the dimension coordinate is already in the file but not
+        # in an approriate group then we have to create a new netCDF
+        # variable. This is to prevent a downstream error ocurring
+        # when the parent data variable tries to reference one of its
+        # netCDF dimensions that is not in the same group nor a parent
+        # group.
+        if already_in_file and not create:
+            ncvar = coord.nc_get_variable("")
+            groups = self._groups(seen[id(coord)]["ncvar"])
+            if not ncvar.startswith(groups):
+                create = True
+
         if create:
-            # ncvar = self._create_netcdf_variable_name(coord,
-            #                                           default='coordinate')
             ncvar = self._create_netcdf_variable_name(coord, default=None)
-
             if ncvar is None:
                 # No netCDF variable name has been set, so use the
                 # corresponding netCDF dimension name
@@ -726,38 +700,6 @@ class NetCDFWrite(IOWrite):
             # Create a new dimension
             unlimited = self.implementation.nc_is_unlimited_axis(f, axis)
 
-            #            if ncdim is None:
-            #                # A netCDF dimension name has NOT been specified, so
-            #                # put the dimension in the root group with the same
-            #                # name as the coordinate variable.
-            #                ncdim = self._remove_group_structure(ncvar)
-            #            elif ncdim in g['dimensions']:
-            #                # A netCDF dimension name has been specified, but
-            #                # matches one already in the file, so put the
-            #                # dimension in the root group with the same name as
-            #                # the coordinate variable.
-            #                ncdim = self._remove_group_structure(ncvar)
-            #            else:
-            #                ncdim = ncvar
-            #            if ncdim is not None:
-            #                # A netCDF dimension name HAS been specified, so make
-            #                # sure that the basename of the coordinate variable
-            #                # matches the basename of the dimension.
-            #                if g['group']:
-            #                    _, groups = self._remove_group_structure(
-            #                        ncvar, return_groups=True)
-            #                    ncvar = groups + self._remove_group_structure(ncdim)
-            #                else:
-            #                    ncvar = ncdim
-            #
-            #            if g['group']:
-            #
-            #                _, groups = self._remove_group_structure(
-            #                    ncvar, return_groups=True)
-            #               ncvar = groups + self._remove_group_structure(ncdim)
-            #            else:
-            #                ncvar = ncdim
-
             self._write_dimension(ncdim, f, axis, unlimited=unlimited)
 
             ncdimensions = self._netcdf_dimensions(f, key, coord)
@@ -998,11 +940,7 @@ class NetCDFWrite(IOWrite):
             geometry_dimension = g["key_to_ncdims"][key][0]
 
             geometry_id = (geometry_dimension, geometry_type)
-            gc.setdefault(
-                geometry_id,
-                {"geometry_type": geometry_type}
-                #                 'geometry_dimension': geometry_dimension,
-            )
+            gc.setdefault(geometry_id, {"geometry_type": geometry_type})
 
             # Nodes
             nodes_ncvar = g["seen"][id(nodes)]["ncvar"]
@@ -1156,6 +1094,13 @@ class NetCDFWrite(IOWrite):
 
             ncdims: `tuple`, optional
 
+                *Parameter example:*
+                  ``('x', 'y')``
+
+                *Parameter example:*
+                  Dimensions can be in a group other than the root
+                  group: ``('/Dew_Point/stations',)``
+
             ignore_type: `bool`, optional
 
         :Returns:
@@ -1385,17 +1330,6 @@ class NetCDFWrite(IOWrite):
             if not bounds_groups and coord_groups:
                 ncvar = coord_groups + ncvar
 
-            #            for ncdim in ncdimensions:
-            #                _, ncdim_groups = self._remove_group_structure(
-            #                    ncdim,
-            #                    return_groups=True)
-            #                if not bounds_groups.startswith(ncdim_groups):
-            #                    raise ValueError(
-            #                        "Can't find create a netCDF variable from {!r} "
-            #                        "with a dimension that is not in the same group or "
-            #                        "a sub-group as the variable: {}".format(bounds, ncdim)
-            #                    )
-
             # Note that, in a field, bounds always have equal units to
             # their parent coordinate
 
@@ -1892,9 +1826,6 @@ class NetCDFWrite(IOWrite):
         {}
 
         """
-        #        if self.implementation.get_data_ndim(bounds) < 3: # DCH
-        #            # No need for a part node count variable required
-        #            return {}
         if self.implementation.get_data_shape(bounds)[1] == 1:
             # No part node count variable required
             return {}
@@ -2163,7 +2094,11 @@ class NetCDFWrite(IOWrite):
         # The netCDF dimensions for the auxiliary coordinate variable
         ncdimensions = self._netcdf_dimensions(f, key, coord)
 
-        if self._already_in_file(coord, ncdimensions):
+        already_in_file = self._already_in_file(coord, ncdimensions)
+
+        create = False
+
+        if already_in_file:
             ncvar = g["seen"][id(coord)]["ncvar"]
 
             # Register that bounds as being the file, too. This is so
@@ -2179,6 +2114,9 @@ class NetCDFWrite(IOWrite):
                         "ncdims": None,
                     }
         else:
+            create = True
+
+        if create:
             if (
                 not self.implementation.get_properties(coord)
                 and self.implementation.get_data(coord, default=None) is None
@@ -2206,9 +2144,6 @@ class NetCDFWrite(IOWrite):
                         ncvar, ncdimensions, coord, extra=extra
                     )
 
-        #                g['key_to_ncvar'][key] = ncvar
-        #                g['key_to_ncdims'][key] = ncdimensions
-
         g["key_to_ncvar"][key] = ncvar
         g["key_to_ncdims"][key] = ncdimensions
 
@@ -2685,8 +2620,8 @@ class NetCDFWrite(IOWrite):
             )  # pragma: no cover
 
         # ------------------------------------------------------------
-        # Check that each dimension of the netCDF variable name is in
-        # the same group or else in a sub-group (CF>=1.8)
+        # Check that each dimension of the netCDF variable is in the
+        # same group or a parent group (CF>=1.8)
         # ------------------------------------------------------------
         if g["group"]:
             groups = self._groups(ncvar)
@@ -2695,8 +2630,8 @@ class NetCDFWrite(IOWrite):
                 if not groups.startswith(ncdim_groups):
                     raise ValueError(
                         f"Can't create netCDF variable {ncvar!r} from "
-                        f"{cfvar!r} with dimension {ncdim!r} that is not in "
-                        "the same group or a sub-group as the variable."
+                        f"{cfvar!r} with netCDF dimension {ncdim!r} that is "
+                        "not in the same group nor in a parent group."
                     )
 
         # ------------------------------------------------------------
@@ -2728,8 +2663,7 @@ class NetCDFWrite(IOWrite):
         kwargs = self._customize_createVariable(cfvar, kwargs)
 
         logger.info(
-            "        to netCDF variable: "
-            f"{ncvar}({ncvar, ', '.join(ncdimensions)})"
+            f"        to netCDF variable: {ncvar}({', '.join(ncdimensions)})"
         )  # pragma: no cover
 
         try:
@@ -3273,8 +3207,7 @@ class NetCDFWrite(IOWrite):
             ncdim = self.implementation.nc_get_dimension(
                 domain_axis, default=None
             )
-            #            if ncdim is not None:
-            #                ncdim = self._netcdf_name(ncdim)
+
             found_dimension_coordinate = False
             for key, dim_coord in dimension_coordinates.items():
                 if self.implementation.get_construct_data_axes(f, key) != (
@@ -3372,18 +3305,6 @@ class NetCDFWrite(IOWrite):
                         )
                         data_axes.append(axis)
 
-                #                spanning_constructs = self.implementation.get_constructs(
-                #                    f, axes=[axis])
-                #
-                #                if axis not in data_axes and spanning_constructs:
-                #                    # The data array doesn't span the domain axis but
-                #                    # an auxiliary coordinate, cell measure, domain
-                #                    # ancillary or field ancillary does, so expand the
-                #                    # data array to include it.
-                #                    f = self.implementation.field_insert_dimension(
-                #                            f, position=0, axis=axis)
-                #                    data_axes.append(axis)
-
                 # If the data array (now) spans this domain axis then
                 # create a netCDF dimension for it
                 if axis in data_axes:
@@ -3438,7 +3359,10 @@ class NetCDFWrite(IOWrite):
                                 if matched_construct:
                                     break
 
-                            if matched_construct:
+                            if (
+                                matched_construct
+                                and self._dimension_in_subgroup(f, ncdim1)
+                            ):
                                 use_existing_dimension = True
                                 break
 
@@ -3913,8 +3837,8 @@ class NetCDFWrite(IOWrite):
                         self.implementation.get_cell_method_qualifiers(cm)
                     ):
                         raise ValueError(
-                            "Can't write {!r}: Unknown cell method "
-                            "property: {!r}".format(org_f, cm.properties())
+                            f"Can't write {org_f!r}: Unknown cell method "
+                            f"property: {cm.properties()!r}"
                         )
 
                     axes = [
@@ -3931,7 +3855,7 @@ class NetCDFWrite(IOWrite):
                 cell_methods = " ".join(cell_methods_strings)
                 logger.info(
                     "    Writing cell_methods attribute to "
-                    "netCDF variable {}: {}".format(ncvar, cell_methods)
+                    f"netCDF variable {ncvar}: {cell_methods}"
                 )  # pragma: no cover
 
                 extra["cell_methods"] = cell_methods
@@ -3967,7 +3891,7 @@ class NetCDFWrite(IOWrite):
 
         # Note that for domain variables the ncdimensions parameter is
         # automatically changed to () within the
-        # _write_netcdf_variable method.  CF-1.9
+        # _write_netcdf_variable method. CF-1.9
         self._write_netcdf_variable(
             ncvar,
             ncdimensions,
@@ -5039,14 +4963,6 @@ class NetCDFWrite(IOWrite):
         g["filename"] = filename
         g["netcdf"] = self.file_open(filename, mode, fmt, fields)
 
-        #        # -----------------------------------------------------------
-        #        # Set the fill mode for a Dataset open for writing to
-        #        off. This # will prevent the data from being pre-filled with
-        #        fill values, # which may result in some performance
-        #        improvements.  #
-        #        -------------------------------------------------------------
-        #        g['netcdf'].set_fill_off()
-
         if not g["dry_run"]:
             # ------------------------------------------------------------
             # Write global properties to the file first. This is important
@@ -5146,3 +5062,39 @@ class NetCDFWrite(IOWrite):
             array = array.astype("int32", casting="same_kind")
 
         return array
+
+    def _dimension_in_subgroup(self, v, ncdim):
+        """Return True if the netCDF dimension is in a valid group.
+
+        Returns True if the dimension is in the same group, or a
+        parent group, as the group defined by the construct. Otherwise
+        return False.
+
+        .. versionadded:: (cfdm) 1.9.0.3
+
+        :Parameters:
+
+            v: Construct
+
+            ncdim: `str`
+                The netCDF dimension name.
+
+                *Parameter example:*
+                  ``'lat'``
+
+                *Parameter example:*
+                  ``'/group1/lat'``
+
+        :Returns:
+
+            `bool`
+                Whether or not the netCDF dimension is in a valid
+                group.
+
+        """
+        v_groups = self.implementation.nc_get_variable_groups(v)
+        v_groups = "/" + "/".join(v_groups)
+
+        ncdim_groups = self._groups(ncdim)
+
+        return v_groups.startswith(ncdim_groups)


Test patch:
diff --git a/cfdm/test/test_Domain.py b/cfdm/test/test_Domain.py
index 59b3ed1bc..bb6836c33 100644
--- a/cfdm/test/test_Domain.py
+++ b/cfdm/test/test_Domain.py
@@ -42,19 +42,12 @@ class DomainTest(unittest.TestCase):
         for title in (None, "title"):
             d.dump(display=False, _title=title)
 
-        # Test when any construct which can have data in fact has no data.
+        # Test when dimension coordinate has no data
         d = d.copy()
-        for identity in [
-            "time",  # a dimension coordinate
-            "latitude",  # an auxiliary coordinate
-            "measure:area",  # a cell measure
-            "surface_altitude",  # a domain ancillary
-        ]:
-            c = d.construct(identity)  # get relevant construct, type as above
-            c.del_data()
-            self.assertFalse(c.has_data())
-            str(d)
-            repr(d)
+        t = d.construct("time")
+        t.del_data()
+        self.assertFalse(t.has_data())
+        str(d)
 
     def test_Domain__init__(self):
         """Test the Domain constructor and source keyword."""
diff --git a/cfdm/test/test_Field.py b/cfdm/test/test_Field.py
index 56985ad63..e0d1f1b90 100644
--- a/cfdm/test/test_Field.py
+++ b/cfdm/test/test_Field.py
@@ -60,21 +60,6 @@ class FieldTest(unittest.TestCase):
         self.assertIsInstance(f.dump(display=False), str)
         self.assertEqual(f.construct_type, "field")
 
-        # Test when any construct which can have data in fact has no data.
-        f = f.copy()
-        for identity in [
-            "time",  # a dimension coordinate
-            "latitude",  # an auxiliary coordinate
-            "measure:area",  # a cell measure
-            "surface_altitude",  # a domain ancillary,
-            "air_temperature standard_error",  # a field ancillary
-        ]:
-            c = f.construct(identity)  # get relevant construct, type as above
-            c.del_data()
-            self.assertFalse(c.has_data())
-            str(f)
-            repr(f)
-
     def test_Field__init__(self):
         """Test the Field constructor and source keyword."""
         cfdm.Field(source="qwerty")
diff --git a/cfdm/test/test_groups.py b/cfdm/test/test_groups.py
index 86655f7d5..c8323bf96 100644
--- a/cfdm/test/test_groups.py
+++ b/cfdm/test/test_groups.py
@@ -8,10 +8,11 @@ import unittest
 faulthandler.enable()  # to debug seg faults and timeouts
 
 import netCDF4
+import numpy as np
 
 import cfdm
 
-n_tmpfiles = 9
+n_tmpfiles = 10
 tmpfiles = [
     tempfile.mkstemp("_test_groups.nc", dir=os.getcwd())[1]
     for i in range(n_tmpfiles)
@@ -26,6 +27,7 @@ tmpfiles = [
     grouped_file3,
     grouped_file4,
     grouped_file5,
+    grouped_file6,
 ) = tmpfiles
 
 
@@ -453,6 +455,65 @@ class GroupsTest(unittest.TestCase):
         h = h[0]
         self.assertTrue(f.equals(h, verbose=3))
 
+    def test_groups_identical_coordinates(self):
+        """Test for identical coordinates in different groups."""
+        # Define dimension and auxiliary coordinates
+        station_name = cfdm.AuxiliaryCoordinate()
+        station_name.set_data(["station1", "station2"])
+        station_name.set_properties(
+            {"long_name": "sensor id", "cf_role": "timeseries_id"}
+        )
+        station_name.nc_set_variable("station_name")
+
+        times = cfdm.DimensionCoordinate()
+        times.set_data(range(5))
+        times.set_properties(
+            {"standard_name": "time", "units": "days since 1967-04-01"}
+        )
+        times.nc_set_variable("time")
+
+        file_content = list()
+        for mt in ("Dew Point", "Wind Speed"):
+            grp_name = mt.replace(" ", "_")
+            var_name = mt.replace(" ", "_")
+
+            # Define domain axes
+            stations = cfdm.DomainAxis(2)
+            stations.nc_set_dimension("stations")
+            stations.nc_set_dimension_groups([grp_name])
+
+            time = cfdm.DomainAxis(5)
+            time.nc_set_dimension("time")
+            time.nc_set_dimension_groups([grp_name])
+
+            # Assign new groups to dimension and auxiliary coordinates
+            station_name = station_name.copy()
+            station_name.nc_set_variable_groups([grp_name])
+
+            times = times.copy()
+            times.nc_set_variable_groups([grp_name])
+
+            # Define the field in the same group as the coordinates
+            f = cfdm.Field()
+
+            station_axis = f.set_construct(stations)
+            time_axis = f.set_construct(time)
+
+            f.set_construct(times, axes=time_axis)
+            f.set_construct(station_name, axes=station_axis)
+
+            f.set_data(
+                np.arange(10).reshape(5, 2), axes=[time_axis, station_axis]
+            )
+
+            f.nc_set_variable(var_name)
+            f.nc_set_variable_groups([grp_name])
+
+            file_content.append(f)
+
+        # This should not raise an exception
+        cfdm.write(file_content, grouped_file6)
+
 
 if __name__ == "__main__":
     print("Run date:", datetime.datetime.now())