Merge branch 'master' into runid-sort

.github/workflows/master.yml

@@ -20,7 +20,7 @@ jobs:
     if: always()
     strategy:
       matrix:
-        python-version: [3.7, 3.8, 3.9]
+        python-version: [3.8, 3.9, '3.10']
         pip-packages:
           - "setuptools pip pytest pytest-cov coverage codecov boutdata xarray numpy>=1.16.0"
       fail-fast: false
@@ -33,7 +33,7 @@ jobs:
         python-version: ${{ matrix.python-version }}
     - name: Install dependencies
       run: |
-        sudo apt-get install libhdf5-dev libnetcdf-dev
+        sudo apt-get update && sudo apt-get install libhdf5-dev libnetcdf-dev
         python -m pip install --upgrade pip
         pip install --upgrade ${{ matrix.pip-packages }}
         pip install -e .
@@ -53,7 +53,7 @@ jobs:
     if: always()
     strategy:
       matrix:
-        python-version: [3.7, 3.8]
+        python-version: [3.8]
         pip-packages:
           - "setuptools pip pytest pytest-cov coverage codecov boutdata==0.1.4 xarray==0.18.0 dask==2.10.0 numpy==1.18.0 natsort==5.5.0 matplotlib==3.1.1 animatplot==0.4.2 netcdf4==1.4.2 Pillow==6.1.0" # test with oldest supported version of packages. Note, using numpy==1.18.0 as a workaround because numpy==1.17.0 is not supported on Python-3.7, even though we should currently support numpy==1.17.0.
       fail-fast: false
@@ -66,7 +66,7 @@ jobs:
         python-version: ${{ matrix.python-version }}
     - name: Install dependencies
       run: |
-        sudo apt-get install libhdf5-dev libnetcdf-dev
+        sudo apt-get update && sudo apt-get install libhdf5-dev libnetcdf-dev
         python -m pip install --upgrade pip
         pip install --upgrade ${{ matrix.pip-packages }}
         pip install -e .

.github/workflows/pythonpackage.yml

@@ -33,7 +33,7 @@ jobs:
         python-version: ${{ matrix.python-version }}
     - name: Install dependencies
       run: |
-        sudo apt-get install libhdf5-dev libnetcdf-dev
+        sudo apt-get update && sudo apt-get install libhdf5-dev libnetcdf-dev
         python -m pip install --upgrade pip
         pip install --upgrade ${{ matrix.pip-packages }}
         pip install -e .
@@ -66,7 +66,7 @@ jobs:
         python-version: ${{ matrix.python-version }}
     - name: Install dependencies
       run: |
-        sudo apt-get install libhdf5-dev libnetcdf-dev
+        sudo apt-get update && sudo apt-get install libhdf5-dev libnetcdf-dev
         python -m pip install --upgrade pip
         pip install --upgrade ${{ matrix.pip-packages }}
         pip install -e .

.github/workflows/pythonpublish.yml

@@ -14,7 +14,7 @@ jobs:
     if: always()
     strategy:
       matrix:
-        python-version: [3.7, 3.8, 3.9]
+        python-version: [3.8, 3.9, '3.10']
         pip-packages:
           - "setuptools pip pytest pytest-cov coverage codecov boutdata xarray numpy>=1.16.0"
       fail-fast: true
@@ -27,7 +27,7 @@ jobs:
         python-version: ${{ matrix.python-version }}
     - name: Install dependencies
       run: |
-        sudo apt-get install libhdf5-dev libnetcdf-dev
+        sudo apt-get update && sudo apt-get install libhdf5-dev libnetcdf-dev
         python -m pip install --upgrade pip
         pip install --upgrade ${{ matrix.pip-packages }}
         pip install -e .
@@ -47,7 +47,7 @@ jobs:
     if: always()
     strategy:
       matrix:
-        python-version: [3.7, 3.8]
+        python-version: [3.8]
         pip-packages:
           - "setuptools pip pytest pytest-cov coverage codecov boutdata==0.1.4 xarray==0.18.0 dask==2.10.0 numpy==1.18.0 natsort==5.5.0 matplotlib==3.1.1 animatplot==0.4.2 netcdf4==1.4.2 Pillow==7.2.0" # test with oldest supported version of packages. Note, using numpy==1.18.0 as a workaround because numpy==1.17.0 is not supported on Python-3.7, even though we should currently support numpy==1.17.0.
       fail-fast: true
@@ -60,7 +60,7 @@ jobs:
         python-version: ${{ matrix.python-version }}
     - name: Install dependencies
       run: |
-        sudo apt-get install libhdf5-dev libnetcdf-dev
+        sudo apt-get update && sudo apt-get install libhdf5-dev libnetcdf-dev
         python -m pip install --upgrade pip
         pip install --upgrade ${{ matrix.pip-packages }}
         pip install -e .
@@ -123,7 +123,7 @@ jobs:
         python-version: '3.x'
     - name: Install dependencies
       run: |
-        sudo apt-get install libhdf5-dev libnetcdf-dev
+        sudo apt-get update && sudo apt-get install libhdf5-dev libnetcdf-dev
         python -m pip install --upgrade pip
         pip install --upgrade setuptools wheel twine
         pip install -e .

setup.cfg

@@ -28,7 +28,7 @@ setup_requires =
     setuptools_scm[toml]>=3.4
     setuptools_scm_git_archive
 install_requires =
-    xarray>=0.18.0,!=2022.9.0,!=2022.10.0
+    xarray>=0.18.0,!=2022.9.0,!=2022.10.0,!=2022.11.0,!=2022.12.0
     boutdata>=0.1.4
     dask[array]>=2.10.0
     gelidum>=0.5.3
@@ -50,6 +50,10 @@ calc =
     xrft
     xhistogram
 docs = sphinx >= 1.4
+3d_plot =
+    k3d >= 2.8.0
+    mayavi >= 4.7.2
+    wand
 
 [build_sphinx]
 project = $metadata.name

xbout/boutdataarray.py

@@ -385,6 +385,9 @@ def interpolate_parallel(
 
         return da
 
+    def add_cartesian_coordinates(self):
+        return _add_cartesian_coordinates(self.data)
+
     def add_cartesian_coordinates(self):
         """
         Add Cartesian (X,Y,Z) coordinates.
@@ -1077,3 +1080,20 @@ def plot_regions(self, ax=None, **kwargs):
         tokamak topology.
         """
         return plotfuncs.plot_regions(self.data, ax=ax, **kwargs)
+
+    def plot3d(self, ax=None, **kwargs):
+        """
+        Make a 3d plot
+
+        Warnings
+        --------
+
+        3d plotting functionality is still a bit of a work in progress. Bugs are likely, and
+        help developing is welcome!
+
+        Parameters
+        ----------
+
+        See plotfuncs.plot3d()
+        """
+        return plotfuncs.plot3d(self.data, **kwargs)

xbout/boutdataset.py

@@ -266,6 +266,9 @@ def find_with_dims(first_var, dims):
 
         return ds
 
+    def add_cartesian_coordinates(self):
+        return _add_cartesian_coordinates(self.data)
+
     def integrate_midpoints(self, variable, *, dims=None, cumulative_t=False):
         """
         Integrate using the midpoint rule for spatial dimensions, and trapezium rule for

xbout/load.py

@@ -35,6 +35,17 @@
     "MYPE",
 ]
 _BOUT_TIME_DEPENDENT_META_VARS = ["iteration", "hist_hi", "tt"]
+_BOUT_GEOMETRY_VARS = [
+    "ixseps1",
+    "ixseps2",
+    "jyseps1_1",
+    "jyseps2_1",
+    "jyseps1_2",
+    "jyseps2_2",
+    "nx",
+    "ny",
+    "ny_inner",
+]
 
 
 # This code should run whenever any function from this module is imported
@@ -350,6 +361,10 @@ def attrs_remove_section(obj, section):
                 pass
             else:
                 raise ValueError(msg)
+        for v in _BOUT_GEOMETRY_VARS:
+            if v not in ds.metadata and v in grid:
+                ds.metadata[v] = grid[v].values
+
     # Update coordinates to match particular geometry of grid
     ds = geometries.apply_geometry(ds, geometry, grid=grid)
 
@@ -365,6 +380,42 @@ def attrs_remove_section(obj, section):
     # BOUT++
     ds.bout.fine_interpolation_factor = 8
 
+    if ("dump" in input_type or "restart" in input_type) and ds.metadata[
+        "BOUT_VERSION"
+    ] < 4.0:
+        # Add workarounds for missing information or different conventions in data saved
+        # by BOUT++ v3.x.
+        for v in ds:
+            if ds.metadata["bout_zdim"] in ds[v].dims:
+                # All fields saved on aligned grid for BOUT-3
+                ds[v].attrs["direction_y"] = "Aligned"
+
+            added_location = False
+            if any(
+                d in ds[v].dims
+                for d in (
+                    ds.metadata["bout_xdim"],
+                    ds.metadata["bout_ydim"],
+                    ds.metadata["bout_zdim"],
+                )
+            ):
+                # zShift, etc. did not support staggered grids in BOUT++ v3 anyway, so
+                # just treat all variables as if they were at CELL_CENTRE
+                ds[v].attrs["cell_location"] = "CELL_CENTRE"
+                added_location = True
+            if added_location:
+                warn(
+                    "Detected data from BOUT++ v3.x. Treating all variables as being "
+                    "at `CELL_CENTRE`. Should be similar to what BOUT++ v3.x did, but "
+                    "if your code uses staggered grids, this may produce unexpected "
+                    "effects in some places."
+                )
+
+        if "nz" not in ds.metadata:
+            # `nz` used to be stored as `MZ` and `MZ` used to include an extra buffer
+            # point that was not used for data.
+            ds.metadata["nz"] = ds.metadata["MZ"] - 1
+
     if info == "terse":
         print("Read in dataset from {}".format(str(Path(datapath))))
     elif info:
@@ -600,17 +651,40 @@ def _auto_open_mfboutdataset(
 
         paths_grid, concat_dims = _arrange_for_concatenation(filepaths, nxpe, nype)
 
-        ds = xr.open_mfdataset(
-            paths_grid,
-            concat_dim=concat_dims,
-            combine="nested",
-            data_vars=data_vars,
-            preprocess=_preprocess,
-            engine=filetype,
-            chunks=chunks,
-            join="exact",
-            **kwargs,
-        )
+        try:
+            ds = xr.open_mfdataset(
+                paths_grid,
+                concat_dim=concat_dims,
+                combine="nested",
+                data_vars=data_vars,
+                preprocess=_preprocess,
+                engine=filetype,
+                chunks=chunks,
+                join="exact",
+                **kwargs,
+            )
+        except ValueError as e:
+            message_to_catch = (
+                "some variables in data_vars are not data variables on the first "
+                "dataset:"
+            )
+            if str(e)[: len(message_to_catch)] == message_to_catch:
+                # Open concatenating any variables that are different in
+                # different files as a work around to support opening older
+                # data.
+                ds = xr.open_mfdataset(
+                    paths_grid,
+                    concat_dim=concat_dims,
+                    combine="nested",
+                    data_vars="different",
+                    preprocess=_preprocess,
+                    engine=filetype,
+                    chunks=chunks,
+                    join="exact",
+                    **kwargs,
+                )
+            else:
+                raise
     else:
         # datapath was nested list of Datasets
 
@@ -744,8 +818,16 @@ def get_nonnegative_scalar(ds, key, default=1, info=True):
 
     # Check whether this is a single file squashed from the multiple output files of a
     # parallel run (i.e. NXPE*NYPE > 1 even though there is only a single file to read).
-    nx = ds["nx"].values
-    ny = ds["ny"].values
+    if "nx" in ds:
+        nx = ds["nx"].values
+    else:
+        # Workaround for older data files
+        nx = ds["MXSUB"].values * ds["NXPE"].values + 2 * ds["MXG"].values
+    if "ny" in ds:
+        ny = ds["ny"].values
+    else:
+        # Workaround for older data files
+        ny = ds["MYSUB"].values * ds["NYPE"].values
     nx_file = ds.dims["x"]
     ny_file = ds.dims["y"]
     is_squashed_doublenull = False
@@ -758,7 +840,10 @@ def get_nonnegative_scalar(ds, key, default=1, info=True):
                 mxg = 0
 
             # Check if there are two divertor targets present
-            if ds["jyseps1_2"] > ds["jyseps2_1"]:
+            # Note: if jyseps2_1 and jyseps1_2 are not in ds it probably
+            # indicates older data and likely the upper target boundary cells
+            # were not saved anyway, so continue as if they were not.
+            if "jyseps2_1" in ds and ds["jyseps1_2"] > ds["jyseps2_1"]:
                 upper_target_cells = myg
             else:
                 upper_target_cells = 0
@@ -771,7 +856,13 @@ def get_nonnegative_scalar(ds, key, default=1, info=True):
 
                 nxpe = 1
                 nype = 1
-                is_squashed_doublenull = (ds["jyseps2_1"] != ds["jyseps1_2"]).values
+                if "jyseps2_1" in ds:
+                    is_squashed_doublenull = (ds["jyseps2_1"] != ds["jyseps1_2"]).values
+                else:
+                    # For older data with no jyseps2_1 or jyseps1_2 in the
+                    # dataset, probably do not need to handle double null data
+                    # squashed with upper target points.
+                    is_squashed_doublenull = False
             elif ny_file == ny + 2 * myg:
                 # Older squashed file from double-null grid but containing only lower
                 # target boundary cells.

xbout/plotting/animate.py

@@ -600,7 +600,7 @@ def animate_line(
 
     # Check plot is the right orientation
     t_read, x_read = data.dims
-    if t_read is animate_over:
+    if t_read == animate_over:
         x = x_read
     else:
         data = data.transpose(animate_over, t_read, transpose_coords=True)