Fix xarray dataset copies

conda_package/mpas_tools/mesh/conversion.py

@@ -246,7 +246,7 @@ def _masks_to_int(dsIn):
         'cullCell',
         'cellSeedMask',
     ]
-    dsOut = xr.Dataset(dsIn, attrs=dsIn.attrs)
+    dsOut = dsIn.copy()
     for var in var_list:
         if var in dsIn:
             dsOut[var] = dsIn[var].astype(np.int32)

conda_package/mpas_tools/viz/transect/horiz.py

@@ -894,7 +894,7 @@ def _sort_intersections(
     for (
         index,
         next_d,
-    ) in zip(sort_indices[1:], d_sorted[1:]):
+    ) in zip(sort_indices[1:], d_sorted[1:], strict=True):
         if next_d - d < epsilon:
             # this d value is effectively the same as the last, so we'll treat
             # it as the same
@@ -1001,7 +1001,7 @@ def _fix_periodic_tris(ds_tris, periodic_var, period):
         np.append(np.arange(0, n_triangles), pos_indices), neg_indices
     )
 
-    ds_new = xr.Dataset(ds_tris)
+    ds_new = ds_tris.copy()
     ds_new[periodic_var] = (('nTriangles', 'nNodes'), coord_node)
     ds_new = ds_new.isel(nTriangles=tri_indices)
     coord_node = ds_new[periodic_var].values

conda_package/pixi.toml

@@ -8,7 +8,7 @@ python = ">=3.10"
 cartopy = "*"
 cmocean = "*"
 dask = "*"
-geometric_features = ">=1.0.1,<2.0.0"
+geometric_features = ">=1.6.3,<2.0.0"
 h5py = "*"
 hdf5 = "*"
 inpoly = "*"

conda_package/recipe/recipe.yaml

@@ -46,7 +46,7 @@ requirements:
     - cartopy
     - cmocean
     - dask
-    - geometric_features >=1.0.1,<2.0.0
+    - geometric_features >=1.6.3,<2.0.0
     - h5py
     - inpoly
     - matplotlib-base >=3.9.0

conda_package/tests/test_conversion.py

@@ -4,7 +4,7 @@
 import numpy as np
 
 from mpas_tools.io import write_netcdf
-from mpas_tools.mesh.conversion import convert, cull, mask
+from mpas_tools.mesh.conversion import _masks_to_int, convert, cull, mask
 from mpas_tools.mesh.spherical import recompute_angle_edge
 
 from .util import get_test_data_file
@@ -47,6 +47,27 @@ def test_conversion_angle_edge():
     assert np.max(np.abs(angle_diff)) < 1.0e-10
 
 
+def test_masks_to_int_dataset_copy():
+    ds_in = xarray.Dataset(
+        data_vars={
+            'regionCellMasks': (
+                ('nCells', 'nRegions'),
+                np.array([[True, False], [False, True]]),
+            ),
+            'cullCell': (('nCells',), np.array([False, True])),
+            'xCell': (('nCells',), np.array([1.0, 2.0])),
+        },
+        attrs={'meshName': 'unit-test'},
+    )
+
+    ds_out = _masks_to_int(ds_in)
+
+    assert ds_out.regionCellMasks.dtype == np.int32
+    assert ds_out.cullCell.dtype == np.int32
+    assert ds_out.attrs == ds_in.attrs
+    assert np.array_equal(ds_out.xCell.values, ds_in.xCell.values)
+
+
 if __name__ == '__main__':
     test_conversion()
     test_conversion_angle_edge()

conda_package/tests/test_viz_transects.py

@@ -1,4 +1,5 @@
 #!/usr/bin/env python
+import numpy as np
 import xarray as xr
 
 from mpas_tools.cime.constants import constants
@@ -11,6 +12,7 @@
     find_spherical_transect_cells_and_weights,
     make_triangle_tree,
 )
+from mpas_tools.viz.transect.horiz import _fix_periodic_tris
 
 from .util import get_test_data_file
 
@@ -130,6 +132,24 @@ def test_find_planar_transect_cells_and_weights():
         assert var in ds_transect.data_vars
 
 
+def test_fix_periodic_tris_dataset_copy():
+    ds_tris = xr.Dataset(
+        data_vars={
+            'xNode': (
+                ('nTriangles', 'nNodes'),
+                np.array([[0.0, 8.0, -8.0]]),
+            ),
+        }
+    )
+
+    ds_new = _fix_periodic_tris(ds_tris, 'xNode', period=10.0)
+
+    assert ds_new.sizes['nTriangles'] == 3
+    assert np.allclose(ds_new.xNode.values[0, :], [0.0, -2.0, 2.0])
+    assert np.allclose(ds_new.xNode.values[1, :], [10.0, 8.0, 12.0])
+    assert np.allclose(ds_new.xNode.values[2, :], [-10.0, -12.0, -8.0])
+
+
 def _get_triangles():
     ds_mesh = xr.open_dataset(get_test_data_file('mesh.QU.1920km.151026.nc'))
     earth_radius = constants['SHR_CONST_REARTH']