diff --git a/.github/workflows/tests.yaml b/.github/workflows/tests.yaml index 8fd72a23..0cf55ed8 100644 --- a/.github/workflows/tests.yaml +++ b/.github/workflows/tests.yaml @@ -26,8 +26,8 @@ jobs: docker_tests: needs: linting runs-on: ubuntu-latest - name: Docker | python=${{ matrix.python-version }} | rasterio${{ matrix.rasterio-version }} | scipy ${{ matrix.run-with-scipy }} - container: osgeo/gdal:ubuntu-full-3.4.2 + name: Docker | GDAL=${{ matrix.gdal-version }} | python=${{ matrix.python-version }} | rasterio${{ matrix.rasterio-version }} | scipy ${{ matrix.run-with-scipy }} + container: osgeo/gdal:ubuntu-full-${{ matrix.gdal-version }} strategy: fail-fast: false matrix: @@ -35,23 +35,28 @@ jobs: rasterio-version: [''] xarray-version: [''] run-with-scipy: ['YES'] + gdal-version: ['3.5.3'] include: - python-version: '3.8' rasterio-version: '' xarray-version: '==0.17' run-with-scipy: 'YES' + gdal-version: '3.4.3' - python-version: '3.8' rasterio-version: '==1.1' xarray-version: '' run-with-scipy: 'YES' + gdal-version: '3.4.3' - python-version: '3.8' rasterio-version: '==1.2.1' xarray-version: '' run-with-scipy: 'YES' + gdal-version: '3.5.3' - python-version: '3.9' rasterio-version: '' xarray-version: '' run-with-scipy: 'NO' + gdal-version: '3.5.3' steps: - uses: actions/checkout@v3 diff --git a/docs/history.rst b/docs/history.rst index 1d24bdf4..07df3fb9 100644 --- a/docs/history.rst +++ b/docs/history.rst @@ -3,6 +3,7 @@ History Latest ------ +- BUG: Fix closing files manually (pull #607) 0.13.0 ------- diff --git a/rioxarray/_io.py b/rioxarray/_io.py index aff047a4..7a5b5756 100644 --- a/rioxarray/_io.py +++ b/rioxarray/_io.py @@ -13,6 +13,7 @@ import re import threading import warnings +from collections import defaultdict from typing import Any, Dict, Hashable, Iterable, List, Optional, Tuple, Union import numpy as np @@ -737,6 +738,30 @@ def _pop_global_netcdf_attrs_from_vars(dataset_to_clean: Dataset) -> Dataset: return dataset_to_clean +def _subdataset_groups_to_dataset( + dim_groups: Dict[Hashable, Dict[Hashable, DataArray]], global_tags: Dict +) -> Union[Dataset, List[Dataset]]: + if dim_groups: + dataset: Union[Dataset, List[Dataset]] = [] + for dim_group in dim_groups.values(): + dataset_group = _pop_global_netcdf_attrs_from_vars( + Dataset(dim_group, attrs=global_tags) + ) + + def _ds_close(): + # pylint: disable=cell-var-from-loop + for data_var in dim_group.values(): + data_var.close() + + dataset_group.set_close(_ds_close) + dataset.append(dataset_group) + if len(dataset) == 1: + dataset = dataset.pop() + else: + dataset = Dataset(attrs=global_tags) + return dataset + + def _load_subdatasets( riods: RasterioReader, group: Optional[Union[str, List[str], Tuple[str, ...]]], @@ -754,8 +779,7 @@ def _load_subdatasets( """ Load in rasterio subdatasets """ - global_tags = _parse_tags(riods.tags()) - dim_groups = {} + dim_groups: Dict[Hashable, Dict[Hashable, DataArray]] = defaultdict(dict) subdataset_filter = None if any((group, variable)): subdataset_filter = build_subdataset_filter(group, variable) @@ -777,23 +801,10 @@ def _load_subdatasets( decode_timedelta=decode_timedelta, **open_kwargs, ) - if shape not in dim_groups: - dim_groups[shape] = {rioda.name: rioda} - else: - dim_groups[shape][rioda.name] = rioda - - if len(dim_groups) > 1: - dataset: Union[Dataset, List[Dataset]] = [ - _pop_global_netcdf_attrs_from_vars(Dataset(dim_group, attrs=global_tags)) - for dim_group in dim_groups.values() - ] - elif not dim_groups: - dataset = Dataset(attrs=global_tags) - else: - dataset = _pop_global_netcdf_attrs_from_vars( - Dataset(list(dim_groups.values())[0], attrs=global_tags) - ) - return dataset + dim_groups[shape][rioda.name] = rioda + return _subdataset_groups_to_dataset( + dim_groups=dim_groups, global_tags=_parse_tags(riods.tags()) + ) def _load_bands_as_variables( @@ -836,7 +847,14 @@ def _load_bands_as_variables( .squeeze() # type: ignore .drop("band") # type: ignore ) - return Dataset(data_vars, attrs=global_tags) + dataset = Dataset(data_vars, attrs=global_tags) + + def _ds_close(): + for data_var in data_vars.values(): + data_var.close() + + dataset.set_close(_ds_close) + return dataset def _prepare_dask( @@ -1070,7 +1088,7 @@ def open_rasterio( captured_warnings = rio_warnings.copy() if band_as_variable: - return _load_bands_as_variables( + dataset_result = _load_bands_as_variables( riods=riods, parse_coordinates=parse_coordinates, chunks=chunks, @@ -1082,6 +1100,8 @@ def open_rasterio( decode_timedelta=decode_timedelta, **open_kwargs, ) + manager.close() + return dataset_result # raise the NotGeoreferencedWarning if applicable for rio_warning in captured_warnings: @@ -1092,7 +1112,7 @@ def open_rasterio( # open the subdatasets if they exist if riods.subdatasets: - return _load_subdatasets( + subdataset_result = _load_subdatasets( riods=riods, group=group, variable=variable, @@ -1106,6 +1126,8 @@ def open_rasterio( decode_timedelta=decode_timedelta, **open_kwargs, ) + manager.close() + return subdataset_result if vrt_params is not None: riods = WarpedVRT(riods, **vrt_params) @@ -1204,9 +1226,6 @@ def open_rasterio( if chunks is not None: result = _prepare_dask(result, riods, filename, chunks) - # Make the file closeable - result.set_close(manager.close) - result.rio._manager = manager # add file path to encoding result.encoding["source"] = riods.name result.encoding["rasterio_dtype"] = str(riods.dtypes[0]) @@ -1224,4 +1243,7 @@ def open_rasterio( for attr, value in result.attrs.items() if not attr.startswith(f"{result.name}#") } + # Make the file closeable + result.set_close(manager.close) + result.rio._manager = manager return result diff --git a/rioxarray/xarray_plugin.py b/rioxarray/xarray_plugin.py index 41639f59..01f6c4d1 100644 --- a/rioxarray/xarray_plugin.py +++ b/rioxarray/xarray_plugin.py @@ -70,8 +70,12 @@ def open_dataset( **open_kwargs, ) if isinstance(rds, xr.DataArray): - rds = rds.to_dataset() - if not isinstance(rds, xr.Dataset): + dataset = rds.to_dataset() + dataset.set_close(rds._close) + rds = dataset + if isinstance(rds, list): + for dataset in rds: + dataset.close() raise RioXarrayError( "Multiple resolution sets found. " "Use 'variable' or 'group' to filter." diff --git a/test/conftest.py b/test/conftest.py index e627af14..da6a49da 100644 --- a/test/conftest.py +++ b/test/conftest.py @@ -2,12 +2,15 @@ import os import pytest +import xarray from numpy.testing import assert_almost_equal, assert_array_equal from packaging import version import rioxarray from rioxarray.raster_array import UNWANTED_RIO_ATTRS +xarray.set_options(warn_for_unclosed_files=True) + TEST_DATA_DIR = os.path.join(os.path.dirname(__file__), "test_data") TEST_INPUT_DATA_DIR = os.path.join(TEST_DATA_DIR, "input") TEST_COMPARE_DATA_DIR = os.path.join(TEST_DATA_DIR, "compare") diff --git a/test/integration/test_integration__io.py b/test/integration/test_integration__io.py index f816d246..de9c1ae1 100644 --- a/test/integration/test_integration__io.py +++ b/test/integration/test_integration__io.py @@ -222,10 +222,11 @@ def test_open_multiple_resolution(): ) assert isinstance(rds_list, list) assert len(rds_list) == 2 - for rds in rds_list: - assert rds.attrs["SHORTNAME"] == "MOD09GA" assert rds_list[0].dims == {"y": 1200, "x": 1200, "band": 1} assert rds_list[1].dims == {"y": 2400, "x": 2400, "band": 1} + for rds in rds_list: + assert rds.attrs["SHORTNAME"] == "MOD09GA" + rds.close() def test_open_group_filter(open_rasterio): @@ -1060,17 +1061,17 @@ def test_rasterio_vrt_gcps__data_exists(): # NOTE: Eventually src_crs will not need to be provided # https://github.com/mapbox/rasterio/pull/2193 with rasterio.vrt.WarpedVRT(src, src_crs=crs) as vrt: - rds = rioxarray.open_rasterio(vrt) - assert rds.values.any() + with rioxarray.open_rasterio(vrt) as rds: + assert rds.values.any() @pytest.mark.parametrize("lock", [True, False]) def test_open_cog(lock): cog_file = os.path.join(TEST_INPUT_DATA_DIR, "cog.tif") - rdsm = rioxarray.open_rasterio(cog_file) - assert rdsm.shape == (1, 500, 500) - rdso = rioxarray.open_rasterio(cog_file, lock=lock, overview_level=0) - assert rdso.shape == (1, 250, 250) + with rioxarray.open_rasterio(cog_file) as rdsm: + assert rdsm.shape == (1, 500, 500) + with rioxarray.open_rasterio(cog_file, lock=lock, overview_level=0) as rdso: + assert rdso.shape == (1, 250, 250) def test_mask_and_scale(open_rasterio): @@ -1287,23 +1288,25 @@ def test_non_rectilinear(): def test_non_rectilinear__load_coords(open_rasterio): test_file = os.path.join(TEST_INPUT_DATA_DIR, "2d_test.tif") - xds = open_rasterio(test_file) - assert xds.rio.shape == (10, 10) - with rasterio.open(test_file) as rds: - assert rds.transform == xds.rio.transform() - for xi, yi in itertools.product(range(rds.width), range(rds.height)): - subset = xds.isel(x=xi, y=yi) - assert_almost_equal(rds.xy(yi, xi), (subset.xc.item(), subset.yc.item())) + with open_rasterio(test_file) as xds: + assert xds.rio.shape == (10, 10) + with rasterio.open(test_file) as rds: + assert rds.transform == xds.rio.transform() + for xi, yi in itertools.product(range(rds.width), range(rds.height)): + subset = xds.isel(x=xi, y=yi) + assert_almost_equal( + rds.xy(yi, xi), (subset.xc.item(), subset.yc.item()) + ) def test_non_rectilinear__skip_parse_coordinates(open_rasterio): test_file = os.path.join(TEST_INPUT_DATA_DIR, "2d_test.tif") - xds = open_rasterio(test_file, parse_coordinates=False) - assert "xc" not in xds.coords - assert "yc" not in xds.coords - assert xds.rio.shape == (10, 10) - with rasterio.open(test_file) as rds: - assert rds.transform == xds.rio.transform() + with open_rasterio(test_file, parse_coordinates=False) as xds: + assert "xc" not in xds.coords + assert "yc" not in xds.coords + assert xds.rio.shape == (10, 10) + with rasterio.open(test_file) as rds: + assert rds.transform == xds.rio.transform() def test_rotation_affine(): @@ -1331,14 +1334,14 @@ def test_rotation_affine(): @pytest.mark.parametrize("dtype", [None, "complex_int16"]) def test_cint16_dtype(dtype, tmp_path): test_file = os.path.join(TEST_INPUT_DATA_DIR, "cint16.tif") - xds = rioxarray.open_rasterio(test_file) - assert xds.rio.shape == (100, 100) - assert xds.dtype == "complex64" - assert xds.encoding["rasterio_dtype"] == "complex_int16" - - tmp_output = tmp_path / "tmp_cint16.tif" - with pytest.warns(NotGeoreferencedWarning): - xds.rio.to_raster(str(tmp_output), dtype=dtype) + with rioxarray.open_rasterio(test_file) as xds: + assert xds.rio.shape == (100, 100) + assert xds.dtype == "complex64" + assert xds.encoding["rasterio_dtype"] == "complex_int16" + + tmp_output = tmp_path / "tmp_cint16.tif" + with pytest.warns(NotGeoreferencedWarning): + xds.rio.to_raster(str(tmp_output), dtype=dtype) with rasterio.open(str(tmp_output)) as riofh: data = riofh.read() assert "complex_int16" in riofh.dtypes @@ -1351,20 +1354,20 @@ def test_cint16_dtype(dtype, tmp_path): ) def test_cint16_dtype_nodata(tmp_path): test_file = os.path.join(TEST_INPUT_DATA_DIR, "cint16.tif") - xds = rioxarray.open_rasterio(test_file) - assert xds.rio.nodata == 0 + with rioxarray.open_rasterio(test_file) as xds: + assert xds.rio.nodata == 0 - tmp_output = tmp_path / "tmp_cint16.tif" - with pytest.warns(NotGeoreferencedWarning): - xds.rio.to_raster(str(tmp_output), dtype="complex_int16") - with rasterio.open(str(tmp_output)) as riofh: - assert riofh.nodata == 0 + tmp_output = tmp_path / "tmp_cint16.tif" + with pytest.warns(NotGeoreferencedWarning): + xds.rio.to_raster(str(tmp_output), dtype="complex_int16") + with rasterio.open(str(tmp_output)) as riofh: + assert riofh.nodata == 0 - # Assign nodata=None - tmp_output = tmp_path / "tmp_cint16_nodata.tif" - xds.rio.write_nodata(None, inplace=True) - with pytest.warns(NotGeoreferencedWarning): - xds.rio.to_raster(str(tmp_output), dtype="complex_int16") + # Assign nodata=None + tmp_output = tmp_path / "tmp_cint16_nodata.tif" + xds.rio.write_nodata(None, inplace=True) + with pytest.warns(NotGeoreferencedWarning): + xds.rio.to_raster(str(tmp_output), dtype="complex_int16") with rasterio.open(str(tmp_output)) as riofh: assert riofh.nodata is None @@ -1373,15 +1376,15 @@ def test_cint16_dtype_nodata(tmp_path): @pytest.mark.parametrize("dtype", [None, "complex_int16"]) def test_cint16_dtype_masked(dtype, tmp_path): test_file = os.path.join(TEST_INPUT_DATA_DIR, "cint16.tif") - xds = rioxarray.open_rasterio(test_file, masked=True) - assert xds.rio.shape == (100, 100) - assert xds.dtype == "complex64" - assert xds.rio.encoded_nodata == 0 - assert np.isnan(xds.rio.nodata) - - tmp_output = tmp_path / "tmp_cint16.tif" - with pytest.warns(NotGeoreferencedWarning): - xds.rio.to_raster(str(tmp_output), dtype=dtype) + with rioxarray.open_rasterio(test_file, masked=True) as xds: + assert xds.rio.shape == (100, 100) + assert xds.dtype == "complex64" + assert xds.rio.encoded_nodata == 0 + assert np.isnan(xds.rio.nodata) + + tmp_output = tmp_path / "tmp_cint16.tif" + with pytest.warns(NotGeoreferencedWarning): + xds.rio.to_raster(str(tmp_output), dtype=dtype) with rasterio.open(str(tmp_output)) as riofh: data = riofh.read() assert "complex_int16" in riofh.dtypes @@ -1392,20 +1395,20 @@ def test_cint16_dtype_masked(dtype, tmp_path): @cint_skip def test_cint16_promote_dtype(tmp_path): test_file = os.path.join(TEST_INPUT_DATA_DIR, "cint16.tif") - xds = rioxarray.open_rasterio(test_file) + with rioxarray.open_rasterio(test_file) as xds: - tmp_output = tmp_path / "tmp_cfloat64.tif" - with pytest.warns(NotGeoreferencedWarning): - xds.rio.to_raster(str(tmp_output), dtype="complex64") - with rasterio.open(str(tmp_output)) as riofh: - data = riofh.read() - assert "complex64" in riofh.dtypes - assert riofh.nodata == 0 - assert data.dtype == "complex64" + tmp_output = tmp_path / "tmp_cfloat64.tif" + with pytest.warns(NotGeoreferencedWarning): + xds.rio.to_raster(str(tmp_output), dtype="complex64") + with rasterio.open(str(tmp_output)) as riofh: + data = riofh.read() + assert "complex64" in riofh.dtypes + assert riofh.nodata == 0 + assert data.dtype == "complex64" - tmp_output = tmp_path / "tmp_cfloat128.tif" - with pytest.warns(NotGeoreferencedWarning): - xds.rio.to_raster(str(tmp_output), dtype="complex128") + tmp_output = tmp_path / "tmp_cfloat128.tif" + with pytest.warns(NotGeoreferencedWarning): + xds.rio.to_raster(str(tmp_output), dtype="complex128") with rasterio.open(str(tmp_output)) as riofh: data = riofh.read() assert "complex128" in riofh.dtypes @@ -1469,13 +1472,15 @@ def test_read_file_handle_with_dask(): with open( os.path.join(TEST_COMPARE_DATA_DIR, "small_dem_3m_merged.tif"), "rb" ) as src: - rioxarray.open_rasterio(src, chunks=2048) + with rioxarray.open_rasterio(src, chunks=2048): + pass @cint_skip def test_read_cint16_with_dask(): test_file = os.path.join(TEST_INPUT_DATA_DIR, "cint16.tif") - rioxarray.open_rasterio(test_file, chunks=True) + with rioxarray.open_rasterio(test_file, chunks=True): + pass @python_vsi_skip diff --git a/test/integration/test_integration_merge.py b/test/integration/test_integration_merge.py index 1f0a225a..32bedf20 100644 --- a/test/integration/test_integration_merge.py +++ b/test/integration/test_integration_merge.py @@ -26,47 +26,47 @@ def test_merge_arrays(squeeze): arrays = [array.squeeze() for array in arrays] merged = merge_arrays(arrays) - if RASTERIO_GE_125: - assert_almost_equal( - merged.rio.bounds(), - (-7274009.6494863, 5003777.3385, -7227678.3778335, 5050108.6101528), - ) - assert merged.rio.shape == (200, 200) - assert_almost_equal(merged.sum(), 22865733) + if RASTERIO_GE_125: + assert_almost_equal( + merged.rio.bounds(), + (-7274009.6494863, 5003777.3385, -7227678.3778335, 5050108.6101528), + ) + assert merged.rio.shape == (200, 200) + assert_almost_equal(merged.sum(), 22865733) + + else: + assert_almost_equal( + merged.rio.bounds(), + ( + -7274009.649486291, + 5003545.682141737, + -7227446.721475236, + 5050108.61015275, + ), + ) + assert merged.rio.shape == (201, 201) + assert_almost_equal(merged.sum(), 11368261) - else: assert_almost_equal( - merged.rio.bounds(), + tuple(merged.rio.transform()), ( + 231.6563582639536, + 0.0, -7274009.649486291, - 5003545.682141737, - -7227446.721475236, + 0.0, + -231.65635826374404, 5050108.61015275, + 0.0, + 0.0, + 1.0, ), ) - assert merged.rio.shape == (201, 201) - assert_almost_equal(merged.sum(), 11368261) - - assert_almost_equal( - tuple(merged.rio.transform()), - ( - 231.6563582639536, - 0.0, - -7274009.649486291, - 0.0, - -231.65635826374404, - 5050108.61015275, - 0.0, - 0.0, - 1.0, - ), - ) - assert merged.rio._cached_transform() == merged.rio.transform() - assert merged.coords["band"].values == [1] - assert sorted(merged.coords) == ["band", "spatial_ref", "x", "y"] - assert merged.rio.crs == rds.rio.crs - assert merged.attrs == rds.attrs - assert merged.encoding["grid_mapping"] == "spatial_ref" + assert merged.rio._cached_transform() == merged.rio.transform() + assert merged.coords["band"].values == [1] + assert sorted(merged.coords) == ["band", "spatial_ref", "x", "y"] + assert merged.rio.crs == rds.rio.crs + assert merged.attrs == rds.attrs + assert merged.encoding["grid_mapping"] == "spatial_ref" @pytest.mark.parametrize("dataset", [True, False]) @@ -87,51 +87,51 @@ def test_merge__different_crs(dataset): else: merged = merge_arrays(arrays, crs=crs) - if dataset: - test_sum = merged[merged.rio.vars[0]].sum() - else: - test_sum = merged.sum() - if RASTERIO_GE_125: - assert_almost_equal( - merged.rio.bounds(), - (-7300984.0238134, 5003618.5908794, -7224054.1109682, 5050108.6101528), - ) - assert merged.rio.shape == (84, 139) - if RASTERIO_GE_13: - assert_almost_equal(test_sum, -131734881) + if dataset: + test_sum = merged[merged.rio.vars[0]].sum() else: - assert_almost_equal(test_sum, -128605446) - else: + test_sum = merged.sum() + if RASTERIO_GE_125: + assert_almost_equal( + merged.rio.bounds(), + (-7300984.0238134, 5003618.5908794, -7224054.1109682, 5050108.6101528), + ) + assert merged.rio.shape == (84, 139) + if RASTERIO_GE_13: + assert_almost_equal(test_sum, -131734881) + else: + assert_almost_equal(test_sum, -128605446) + else: + assert_almost_equal( + merged.rio.bounds(), + (-7300984.0238134, 5003618.5908794, -7223500.6583578, 5050108.6101528), + ) + assert merged.rio.shape == (84, 140) + assert_almost_equal(test_sum, -131013894) assert_almost_equal( - merged.rio.bounds(), - (-7300984.0238134, 5003618.5908794, -7223500.6583578, 5050108.6101528), + tuple(merged.rio.transform()), + ( + 553.4526103969893, + 0.0, + -7300984.023813409, + 0.0, + -553.4526103969796, + 5050108.610152751, + 0.0, + 0.0, + 1.0, + ), ) - assert merged.rio.shape == (84, 140) - assert_almost_equal(test_sum, -131013894) - assert_almost_equal( - tuple(merged.rio.transform()), - ( - 553.4526103969893, - 0.0, - -7300984.023813409, - 0.0, - -553.4526103969796, - 5050108.610152751, - 0.0, - 0.0, - 1.0, - ), - ) - assert merged.coords["band"].values == [1] - assert sorted(merged.coords) == ["band", "spatial_ref", "x", "y"] - assert merged.rio.crs == rds.rio.crs - if not dataset: - assert merged.attrs == { - "_FillValue": -28672, - "add_offset": 0.0, - "scale_factor": 1.0, - } - assert merged.encoding["grid_mapping"] == "spatial_ref" + assert merged.coords["band"].values == [1] + assert sorted(merged.coords) == ["band", "spatial_ref", "x", "y"] + assert merged.rio.crs == rds.rio.crs + if not dataset: + assert merged.attrs == { + "_FillValue": -28672, + "add_offset": 0.0, + "scale_factor": 1.0, + } + assert merged.encoding["grid_mapping"] == "spatial_ref" def test_merge_arrays__res(): @@ -148,32 +148,42 @@ def test_merge_arrays__res(): ] merged = merge_arrays(arrays, res=(300, 300)) - if RASTERIO_GE_125: - assert_almost_equal( - merged.rio.bounds(), - (-7274009.6494863, 5003608.6101528, -7227509.6494863, 5050108.6101528), - ) - assert merged.rio.shape == (155, 155) - else: + if RASTERIO_GE_125: + assert_almost_equal( + merged.rio.bounds(), + (-7274009.6494863, 5003608.6101528, -7227509.6494863, 5050108.6101528), + ) + assert merged.rio.shape == (155, 155) + else: + assert_almost_equal( + merged.rio.bounds(), + (-7274009.6494863, 5003308.6101528, -7227209.6494863, 5050108.6101528), + ) + assert merged.rio.shape == (156, 156) + assert_almost_equal( - merged.rio.bounds(), - (-7274009.6494863, 5003308.6101528, -7227209.6494863, 5050108.6101528), + tuple(merged.rio.transform()), + ( + 300.0, + 0.0, + -7274009.649486291, + 0.0, + -300.0, + 5050108.61015275, + 0.0, + 0.0, + 1.0, + ), ) - assert merged.rio.shape == (156, 156) - - assert_almost_equal( - tuple(merged.rio.transform()), - (300.0, 0.0, -7274009.649486291, 0.0, -300.0, 5050108.61015275, 0.0, 0.0, 1.0), - ) - assert merged.rio._cached_transform() == merged.rio.transform() - assert merged.coords["band"].values == [1] - assert sorted(merged.coords) == ["band", "spatial_ref", "x", "y"] - assert merged.rio.crs == rds.rio.crs - assert_almost_equal(merged.attrs.pop("_FillValue"), rds.attrs.pop("_FillValue")) - compare_attrs = dict(rds.attrs) - assert merged.attrs == compare_attrs - assert merged.encoding["grid_mapping"] == "spatial_ref" - assert_almost_equal(nansum(merged), 13760565) + assert merged.rio._cached_transform() == merged.rio.transform() + assert merged.coords["band"].values == [1] + assert sorted(merged.coords) == ["band", "spatial_ref", "x", "y"] + assert merged.rio.crs == rds.rio.crs + assert_almost_equal(merged.attrs.pop("_FillValue"), rds.attrs.pop("_FillValue")) + compare_attrs = dict(rds.attrs) + assert merged.attrs == compare_attrs + assert merged.encoding["grid_mapping"] == "spatial_ref" + assert_almost_equal(nansum(merged), 13760565) @pytest.mark.xfail(os.name == "nt", reason="On windows the merged data is different.") @@ -189,54 +199,59 @@ def test_merge_datasets(): rds.isel(x=slice(600, None), y=slice(600)), ] merged = merge_datasets(datasets) - data_vars = sorted(merged.data_vars) - assert data_vars == [ - "QC_500m_1", - "iobs_res_1", - "num_observations_500m", - "obscov_500m_1", - "sur_refl_b01_1", - "sur_refl_b02_1", - "sur_refl_b03_1", - "sur_refl_b04_1", - "sur_refl_b05_1", - "sur_refl_b06_1", - "sur_refl_b07_1", - ] - data_var = data_vars[0] - if RASTERIO_GE_125: - assert_almost_equal( - merged[data_var].rio.bounds(), - (-4447802.078667, -10007554.677, -3335851.559, -8895604.157333), - ) - assert merged.rio.shape == (2400, 2400) - assert_almost_equal(merged[data_var].sum(), 4539666606551516) - else: + data_vars = sorted(merged.data_vars) + assert data_vars == [ + "QC_500m_1", + "iobs_res_1", + "num_observations_500m", + "obscov_500m_1", + "sur_refl_b01_1", + "sur_refl_b02_1", + "sur_refl_b03_1", + "sur_refl_b04_1", + "sur_refl_b05_1", + "sur_refl_b06_1", + "sur_refl_b07_1", + ] + data_var = data_vars[0] + if RASTERIO_GE_125: + assert_almost_equal( + merged[data_var].rio.bounds(), + (-4447802.078667, -10007554.677, -3335851.559, -8895604.157333), + ) + assert merged.rio.shape == (2400, 2400) + assert_almost_equal(merged[data_var].sum(), 4539666606551516) + else: + assert_almost_equal( + merged[data_var].rio.bounds(), + ( + -4447802.078667, + -10008017.989716524, + -3335388.246283474, + -8895604.157333, + ), + ) + assert merged.rio.shape == (2401, 2401) + assert_almost_equal(merged[data_var].sum(), 4543446965182987) assert_almost_equal( - merged[data_var].rio.bounds(), - (-4447802.078667, -10008017.989716524, -3335388.246283474, -8895604.157333), + tuple(merged[data_var].rio.transform()), + ( + 463.3127165279158, + 0.0, + -4447802.078667, + 0.0, + -463.3127165279151, + -8895604.157333, + 0.0, + 0.0, + 1.0, + ), ) - assert merged.rio.shape == (2401, 2401) - assert_almost_equal(merged[data_var].sum(), 4543446965182987) - assert_almost_equal( - tuple(merged[data_var].rio.transform()), - ( - 463.3127165279158, - 0.0, - -4447802.078667, - 0.0, - -463.3127165279151, - -8895604.157333, - 0.0, - 0.0, - 1.0, - ), - ) - assert merged.coords["band"].values == [1] - assert sorted(merged.coords) == ["band", "spatial_ref", "x", "y"] - assert merged.rio.crs == rds.rio.crs - assert merged.attrs == rds.attrs - assert merged.encoding["grid_mapping"] == "spatial_ref" + assert merged.coords["band"].values == [1] + assert sorted(merged.coords) == ["band", "spatial_ref", "x", "y"] + assert merged.rio.crs == rds.rio.crs + assert merged.attrs == rds.attrs + assert merged.encoding["grid_mapping"] == "spatial_ref" @pytest.mark.xfail(os.name == "nt", reason="On windows the merged data is different.") @@ -252,33 +267,43 @@ def test_merge_datasets__res(): rds.isel(x=slice(1200, None), y=slice(1200)), ] merged = merge_datasets(datasets, res=1000) - data_vars = sorted(merged.data_vars) - assert data_vars == [ - "QC_500m_1", - "iobs_res_1", - "num_observations_500m", - "obscov_500m_1", - "sur_refl_b01_1", - "sur_refl_b02_1", - "sur_refl_b03_1", - "sur_refl_b04_1", - "sur_refl_b05_1", - "sur_refl_b06_1", - "sur_refl_b07_1", - ] - data_var = data_vars[0] - assert_almost_equal( - merged[data_var].rio.bounds(), - (-4447802.078667, -10007604.157333, -3335802.078667, -8895604.157333), - ) - assert_almost_equal( - tuple(merged[data_var].rio.transform()), - (1000.0, 0.0, -4447802.078667, 0.0, -1000.0, -8895604.157333, 0.0, 0.0, 1.0), - ) - assert merged.rio.shape == (1112, 1112) - assert merged.coords["band"].values == [1] - assert sorted(merged.coords) == ["band", "spatial_ref", "x", "y"] - assert merged.rio.crs == rds.rio.crs - assert merged.attrs == rds.attrs - assert merged.encoding["grid_mapping"] == "spatial_ref" - assert_almost_equal(merged[data_var].sum(), 974566547463955) + data_vars = sorted(merged.data_vars) + assert data_vars == [ + "QC_500m_1", + "iobs_res_1", + "num_observations_500m", + "obscov_500m_1", + "sur_refl_b01_1", + "sur_refl_b02_1", + "sur_refl_b03_1", + "sur_refl_b04_1", + "sur_refl_b05_1", + "sur_refl_b06_1", + "sur_refl_b07_1", + ] + data_var = data_vars[0] + assert_almost_equal( + merged[data_var].rio.bounds(), + (-4447802.078667, -10007604.157333, -3335802.078667, -8895604.157333), + ) + assert_almost_equal( + tuple(merged[data_var].rio.transform()), + ( + 1000.0, + 0.0, + -4447802.078667, + 0.0, + -1000.0, + -8895604.157333, + 0.0, + 0.0, + 1.0, + ), + ) + assert merged.rio.shape == (1112, 1112) + assert merged.coords["band"].values == [1] + assert sorted(merged.coords) == ["band", "spatial_ref", "x", "y"] + assert merged.rio.crs == rds.rio.crs + assert merged.attrs == rds.attrs + assert merged.encoding["grid_mapping"] == "spatial_ref" + assert_almost_equal(merged[data_var].sum(), 974566547463955) diff --git a/test/integration/test_integration_rioxarray.py b/test/integration/test_integration_rioxarray.py index 85cfa50e..983c0e8e 100644 --- a/test/integration/test_integration_rioxarray.py +++ b/test/integration/test_integration_rioxarray.py @@ -418,25 +418,25 @@ def test_clip_box__one_dimension_error(modis_clip): def test_clip_box__nodata_in_bounds(): - xds = rioxarray.open_rasterio( + with rioxarray.open_rasterio( os.path.join(TEST_INPUT_DATA_DIR, "clip_bbox__out_of_bounds.tif") - ) - with pytest.raises(NoDataInBounds): - xds.rio.clip_box( - 135.7821043660001123, - -17.1608065079999506, - 135.7849362810001139, - -17.1580839999999739, - crs="EPSG:4326", - ) + ) as xds: + with pytest.raises(NoDataInBounds): + xds.rio.clip_box( + 135.7821043660001123, + -17.1608065079999506, + 135.7849362810001139, + -17.1580839999999739, + crs="EPSG:4326", + ) - with pytest.raises(NoDataInBounds): - xds.rio.reproject("EPSG:4326").rio.clip_box( - 135.7821043660001123, - -17.1608065079999506, - 135.7849362810001139, - -17.1580839999999739, - ) + with pytest.raises(NoDataInBounds): + xds.rio.reproject("EPSG:4326").rio.clip_box( + 135.7821043660001123, + -17.1608065079999506, + 135.7849362810001139, + -17.1580839999999739, + ) def test_clip_box__reproject_bounds(modis_clip): @@ -1070,27 +1070,27 @@ def test_reproject__gcps_kwargs(tmp_path): ) as source: source.gcps = (src_gcps, crs) - rds = rioxarray.open_rasterio(tiffname) - rds.rio.write_crs(crs, inplace=True) - rds = rds.rio.reproject( - crs, - gcps=src_gcps, - ) - assert rds.rio.height == 923 - assert rds.rio.width == 1027 - assert rds.rio.crs == crs - assert rds.rio.transform().almost_equals( - Affine( - 216.8587081056465, - 0.0, - 115698.25, - 0.0, - -216.8587081056465, - 2818720.0, + with rioxarray.open_rasterio(tiffname) as rds: + rds.rio.write_crs(crs, inplace=True) + rds = rds.rio.reproject( + crs, + gcps=src_gcps, ) - ) - assert (rds.coords["x"].values > 11000).all() - assert (rds.coords["y"].values > 281700).all() + assert rds.rio.height == 923 + assert rds.rio.width == 1027 + assert rds.rio.crs == crs + assert rds.rio.transform().almost_equals( + Affine( + 216.8587081056465, + 0.0, + 115698.25, + 0.0, + -216.8587081056465, + 2818720.0, + ) + ) + assert (rds.coords["x"].values > 11000).all() + assert (rds.coords["y"].values > 281700).all() def test_reproject_match(modis_reproject_match): @@ -1803,11 +1803,11 @@ def test_to_raster__different_dtype(tmp_path, windowed): ), ): test_nd.rio.to_raster(tmpfile, dtype=numpy.uint8, windowed=windowed) - xds = rioxarray.open_rasterio(tmpfile) - assert str(xds.dtype) == "uint8" - assert xds.attrs["_FillValue"] == 255 - assert xds.rio.nodata == 255 - assert xds.squeeze().values[1, 1] == 255 + with rioxarray.open_rasterio(tmpfile) as xds: + assert str(xds.dtype) == "uint8" + assert xds.attrs["_FillValue"] == 255 + assert xds.rio.nodata == 255 + assert xds.squeeze().values[1, 1] == 255 def test_missing_spatial_dimensions(): @@ -2605,40 +2605,42 @@ def test_missing_crs_error_msg(): def test_missing_transform_bounds(): - xds = rioxarray.open_rasterio( + with rioxarray.open_rasterio( os.path.join(TEST_COMPARE_DATA_DIR, "small_dem_3m_merged.tif"), parse_coordinates=False, - ) - xds.coords["spatial_ref"].attrs.pop("GeoTransform") - with pytest.raises(DimensionMissingCoordinateError): - xds.rio.bounds() + ) as xds: + xds.coords["spatial_ref"].attrs.pop("GeoTransform") + with pytest.raises(DimensionMissingCoordinateError): + xds.rio.bounds() def test_missing_transform_resolution(): - xds = rioxarray.open_rasterio( + with rioxarray.open_rasterio( os.path.join(TEST_COMPARE_DATA_DIR, "small_dem_3m_merged.tif"), parse_coordinates=False, - ) - xds.coords["spatial_ref"].attrs.pop("GeoTransform") - with pytest.raises(DimensionMissingCoordinateError): - xds.rio.resolution() + ) as xds: + xds.coords["spatial_ref"].attrs.pop("GeoTransform") + with pytest.raises(DimensionMissingCoordinateError): + xds.rio.resolution() def test_shape_order(): - rds = rioxarray.open_rasterio(os.path.join(TEST_INPUT_DATA_DIR, "tmmx_20190121.nc")) - assert rds.air_temperature.rio.shape == (585, 1386) + with rioxarray.open_rasterio( + os.path.join(TEST_INPUT_DATA_DIR, "tmmx_20190121.nc") + ) as rds: + assert rds.air_temperature.rio.shape == (585, 1386) def test_write_transform__from_read(tmp_path): - xds = rioxarray.open_rasterio( + with rioxarray.open_rasterio( os.path.join(TEST_COMPARE_DATA_DIR, "small_dem_3m_merged.tif"), parse_coordinates=False, - ) - out_file = tmp_path / "test_geotransform.nc" - xds.to_netcdf(out_file) - xds2 = rioxarray.open_rasterio(out_file, parse_coordinates=False) - assert_almost_equal(tuple(xds2.rio.transform()), tuple(xds.rio.transform())) - assert xds.spatial_ref.GeoTransform == xds2.spatial_ref.GeoTransform + ) as xds: + out_file = tmp_path / "test_geotransform.nc" + xds.to_netcdf(out_file) + with rioxarray.open_rasterio(out_file, parse_coordinates=False) as xds2: + assert_almost_equal(tuple(xds2.rio.transform()), tuple(xds.rio.transform())) + assert xds.spatial_ref.GeoTransform == xds2.spatial_ref.GeoTransform def test_write_transform(): @@ -2829,18 +2831,18 @@ def test_grid_mapping_default(): def test_estimate_utm_crs(): - xds = rioxarray.open_rasterio( + with rioxarray.open_rasterio( os.path.join(TEST_INPUT_DATA_DIR, "cog.tif"), - ) - if PYPROJ_LT_3: - with pytest.raises(RuntimeError, match=r"pyproj 3\+ required"): - xds.rio.estimate_utm_crs() - else: - assert xds.rio.estimate_utm_crs().to_epsg() in (32618, 32655) - assert xds.rio.reproject("EPSG:4326").rio.estimate_utm_crs() == CRS.from_epsg( - 32618 - ) - assert xds.rio.estimate_utm_crs("WGS 72") in (32218, 32255) + ) as xds: + if PYPROJ_LT_3: + with pytest.raises(RuntimeError, match=r"pyproj 3\+ required"): + xds.rio.estimate_utm_crs() + else: + assert xds.rio.estimate_utm_crs().to_epsg() in (32618, 32655) + assert xds.rio.reproject( + "EPSG:4326" + ).rio.estimate_utm_crs() == CRS.from_epsg(32618) + assert xds.rio.estimate_utm_crs("WGS 72") in (32218, 32255) @pytest.mark.skipif(PYPROJ_LT_3, reason="pyproj 3+ required") @@ -2976,23 +2978,23 @@ def test_reproject__gcps_file(tmp_path): ) as source: source.gcps = (src_gcps, crs) - rds = rioxarray.open_rasterio(tiffname) - rds = rds.rio.reproject( - crs, - ) - assert rds.rio.height == 923 - assert rds.rio.width == 1027 - assert rds.rio.crs == crs - assert rds.rio.transform().almost_equals( - Affine( - 216.8587081056465, - 0.0, - 115698.25, - 0.0, - -216.8587081056465, - 2818720.0, + with rioxarray.open_rasterio(tiffname) as rds: + rds = rds.rio.reproject( + crs, + ) + assert rds.rio.height == 923 + assert rds.rio.width == 1027 + assert rds.rio.crs == crs + assert rds.rio.transform().almost_equals( + Affine( + 216.8587081056465, + 0.0, + 115698.25, + 0.0, + -216.8587081056465, + 2818720.0, + ) ) - ) def test_bounds__ordered__dataarray(): diff --git a/test/integration/test_integration_xarray_plugin.py b/test/integration/test_integration_xarray_plugin.py index 5d934a80..7cca368b 100644 --- a/test/integration/test_integration_xarray_plugin.py +++ b/test/integration/test_integration_xarray_plugin.py @@ -31,7 +31,7 @@ def test_xarray_open_dataset__drop_variables(): TEST_INPUT_DATA_DIR, "MOD09GA.A2008296.h14v17.006.2015181011753.hdf" ) - rds = xr.open_dataset( + with xr.open_dataset( input_file, engine="rasterio", group="MODIS_Grid_500m_2D", @@ -44,14 +44,13 @@ def test_xarray_open_dataset__drop_variables(): "sur_refl_b06_1", "sur_refl_b07_1", ], - ) - - assert sorted(rds.data_vars) == [ - "QC_500m_1", - "iobs_res_1", - "num_observations_500m", - "obscov_500m_1", - ] + ) as rds: + assert sorted(rds.data_vars) == [ + "QC_500m_1", + "iobs_res_1", + "num_observations_500m", + "obscov_500m_1", + ] def test_open_multiple_resolution():