diff --git a/docs/api.rst b/docs/api.rst index 93daa220a..1344745bb 100644 --- a/docs/api.rst +++ b/docs/api.rst @@ -780,16 +780,6 @@ Vector gis_utils.nearest_merge -PCRaster I/O -============ - -.. autosummary:: - :toctree: _generated - - gis_utils.write_map - gis_utils.write_clone - - .. _statistics: ===================================== diff --git a/docs/changelog.rst b/docs/changelog.rst index 00789d35f..e635900ae 100644 --- a/docs/changelog.rst +++ b/docs/changelog.rst @@ -31,7 +31,7 @@ Fixed Deprecated ---------- -- +- the dependencies ``pcraster`` and ``pygeos`` are no longer used and were removed. (PR #467) v0.8.0 (2023-07-18) diff --git a/hydromt/__init__.py b/hydromt/__init__.py index 6be7ef6b9..6dd041048 100644 --- a/hydromt/__init__.py +++ b/hydromt/__init__.py @@ -3,12 +3,6 @@ # version number without 'v' at start __version__ = "0.8.1.dev0" -# Set environment variables (this will be temporary) -# to use shapely 2.0 in favor of pygeos (if installed) -import os - -os.environ["USE_PYGEOS"] = "0" - # pkg_resource deprication warnings originate from dependencies # so silence them for now import warnings diff --git a/hydromt/_compat.py b/hydromt/_compat.py index a4db3e952..bc30b9848 100644 --- a/hydromt/_compat.py +++ b/hydromt/_compat.py @@ -3,10 +3,8 @@ __all__ = [] HAS_XUGRID = False -HAS_PCRASTER = False HAS_SHAPELY20 = False HAS_PYET = False -HAS_PYGEOS = False HAS_GCSFS = False HAS_S3FS = False HAS_OPENPYXL = False @@ -26,21 +24,6 @@ except ImportError: False -try: - import pygeos - - HAS_PYGEOS = True -except ImportError: - pass - -try: - import pcraster - - HAS_PCRASTER = True -except ImportError: - pass - - try: import xugrid diff --git a/hydromt/gis_utils.py b/hydromt/gis_utils.py index 7973b5067..d90ae387b 100644 --- a/hydromt/gis_utils.py +++ b/hydromt/gis_utils.py @@ -6,14 +6,12 @@ import logging import os import subprocess -from os.path import dirname, isfile, join +from os.path import dirname, join from typing import Optional, Tuple import geopandas as gpd import numpy as np -import rasterio import xarray as xr -from pyflwdir import core_conversion, core_d8, core_ldd from pyflwdir import gis_utils as gis from pyproj import CRS from rasterio.transform import Affine @@ -76,7 +74,6 @@ "jpg": "JPEG", "kro": "KRO", "lcp": "LCP", - "map": "PCRaster", "mbtiles": "MBTiles", "mpr/mpl": "ILWIS", "ntf": "NITF", @@ -498,121 +495,6 @@ def spread2d( return ds_out -## PCRASTER - - -def write_clone(tmpdir, gdal_transform, wkt_projection, shape): - """Write pcraster clone file to a tmpdir using gdal.""" - from osgeo import gdal - - gdal.AllRegister() - driver1 = gdal.GetDriverByName("GTiff") - driver2 = gdal.GetDriverByName("PCRaster") - fn = join(tmpdir, "clone.map") - # create temp tif file - fn_temp = join(tmpdir, "clone.tif") - TempDataset = driver1.Create(fn_temp, shape[1], shape[0], 1, gdal.GDT_Float32) - TempDataset.SetGeoTransform(gdal_transform) - if wkt_projection is not None: - TempDataset.SetProjection(wkt_projection) - # TODO set csr - # copy to pcraster format - driver2.CreateCopy(fn, TempDataset, 0) - # close and cleanup - TempDataset = None - return fn - - -def write_map( - data, - raster_path, - nodata, - transform, - crs=None, - clone_path=None, - pcr_vs="scalar", - **kwargs, -): - """Write pcraster map files using pcr.report functionality. - - A PCRaster clone map is written to a temporary directory if not provided. - For PCRaster types see https://www.gdal.org/frmt_various.html#PCRaster - - Parameters - ---------- - data : ndarray - Raster data - raster_path : str - Path to output map - nodata : int, float - no data value - transform : affine transform - Two dimensional affine transform for 2D linear mapping - clone_path : str, optional - Path to PCRaster clone map, by default None - pcr_vs : str, optional - pcraster type, by default "scalar" - **kwargs: - not used in this function, mainly here for compatability reasons. - crs: - The coordinate reference system of the data. - - - Raises - ------ - ImportError - pcraster package is required - ValueError - if invalid ldd - """ - if not _compat.HAS_PCRASTER: - raise ImportError("The pcraster package is required to write map files") - import tempfile - - import pcraster as pcr - - with tempfile.TemporaryDirectory() as tmpdir: - # deal with pcr clone map - if clone_path is None: - clone_path = write_clone( - tmpdir, - gdal_transform=transform.to_gdal(), - wkt_projection=None if crs is None else CRS.from_user_input(crs).wkt, - shape=data.shape, - ) - elif not isfile(clone_path): - raise IOError(f'clone_path: "{clone_path}" does not exist') - pcr.setclone(clone_path) - if nodata is None and pcr_vs != "ldd": - raise ValueError("nodata value required to write PCR map") - # write to pcrmap - if pcr_vs == "ldd": - # if d8 convert to ldd - data = data.astype(np.uint8) # force dtype - if core_d8.isvalid(data): - data = core_conversion.d8_to_ldd(data) - elif not core_ldd.isvalid(data): - raise ValueError("LDD data not understood") - mv = int(core_ldd._mv) - ldd = pcr.numpy2pcr(pcr.Ldd, data.astype(int), mv) - # make sure it is pcr sound - # NOTE this should not be necessary - pcrmap = pcr.lddrepair(ldd) - elif pcr_vs == "bool": - pcrmap = pcr.numpy2pcr(pcr.Boolean, data.astype(bool), np.bool_(nodata)) - elif pcr_vs == "scalar": - pcrmap = pcr.numpy2pcr(pcr.Scalar, data.astype(float), float(nodata)) - elif pcr_vs == "ordinal": - pcrmap = pcr.numpy2pcr(pcr.Ordinal, data.astype(int), int(nodata)) - elif pcr_vs == "nominal": - pcrmap = pcr.numpy2pcr(pcr.Nominal, data.astype(int), int(nodata)) - pcr.report(pcrmap, raster_path) - # set crs (pcrmap ignores this info from clone ??) - if crs is not None: - with rasterio.open(raster_path, "r+") as dst: - dst.crs = crs - - def create_vrt( vrt_path: str, file_list_path: str = None, diff --git a/hydromt/raster.py b/hydromt/raster.py index 6e4bfc057..308f36f03 100644 --- a/hydromt/raster.py +++ b/hydromt/raster.py @@ -12,7 +12,7 @@ import os import tempfile from itertools import product -from os.path import basename, dirname, isdir, join +from os.path import isdir, join from typing import Any, Optional, Union import dask @@ -2311,53 +2311,35 @@ def to_raster( count = da_out[dim0].size da_out = da_out.sortby(dim0) # write - if driver.lower() == "pcraster" and _compat.HAS_PCRASTER: - for i in range(count): - if dim0: - bname = basename(raster_path).split(".")[0] - bname = f"{bname[:8]:8s}".replace(" ", "0") - raster_path = join(dirname(raster_path), f"{bname}.{i+1:03d}") - data = da_out.isel({dim0: i}).load().squeeze().data + profile = dict( + driver=driver, + height=da_out.raster.height, + width=da_out.raster.width, + count=count, + dtype=str(da_out.dtype), + crs=da_out.raster.crs, + transform=da_out.raster.transform, + nodata=nodata, + **profile_kwargs, + ) + with rasterio.open(raster_path, "w", **profile) as dst: + if windowed: + window_iter = dst.block_windows(1) + else: + window_iter = [(None, None)] + for _, window in window_iter: + if window is not None: + row_slice, col_slice = window.toslices() + sel = {self.x_dim: col_slice, self.y_dim: row_slice} + data = da_out.isel(sel).load().values else: - data = da_out.load().data - gis_utils.write_map( - data, - raster_path, - crs=da_out.raster.crs, - transform=da_out.raster.transform, - nodata=nodata, - **profile_kwargs, - ) - else: - profile = dict( - driver=driver, - height=da_out.raster.height, - width=da_out.raster.width, - count=count, - dtype=str(da_out.dtype), - crs=da_out.raster.crs, - transform=da_out.raster.transform, - nodata=nodata, - **profile_kwargs, - ) - with rasterio.open(raster_path, "w", **profile) as dst: - if windowed: - window_iter = dst.block_windows(1) + data = da_out.load().values + if data.ndim == 2: + dst.write(data, 1, window=window) else: - window_iter = [(None, None)] - for _, window in window_iter: - if window is not None: - row_slice, col_slice = window.toslices() - sel = {self.x_dim: col_slice, self.y_dim: row_slice} - data = da_out.isel(sel).load().values - else: - data = da_out.load().values - if data.ndim == 2: - dst.write(data, 1, window=window) - else: - dst.write(data, window=window) - if tags is not None: - dst.update_tags(**tags) + dst.write(data, window=window) + if tags is not None: + dst.update_tags(**tags) def vectorize(self, connectivity=8): """Return geometry of grouped pixels with the same value in a DataArray object. diff --git a/pyproject.toml b/pyproject.toml index d46a9f6c2..e2a45be10 100644 --- a/pyproject.toml +++ b/pyproject.toml @@ -87,10 +87,6 @@ examples = [ "matplotlib-base", # plotting "notebook", # jupyter integration ] -deprecated = [ - "pcraster", # Old Wflow model was build on this - "pygeos", # Supersceded by shapely 2 -] full = ["hydromt[io,extra,dev,test,doc,examples]"] slim = ["hydromt[io,extra,examples]"] diff --git a/tests/conftest.py b/tests/conftest.py index 1cfe53bb6..a5b21c8f8 100644 --- a/tests/conftest.py +++ b/tests/conftest.py @@ -1,6 +1,3 @@ -import os - -os.environ["USE_PYGEOS"] = "0" import geopandas as gpd import numpy as np import pandas as pd diff --git a/tests/test_io.py b/tests/test_io.py index 00a1af4b2..9bc762f05 100644 --- a/tests/test_io.py +++ b/tests/test_io.py @@ -206,30 +206,3 @@ def test_rasterio_errors(tmpdir, rioda): da0.raster.to_raster(str(tmpdir.join("test2.tif")), count=3) with pytest.raises(ValueError, match="Extension unknown for driver"): da0.to_dataset().raster.to_mapstack(root=str(tmpdir), driver="unknown") - - -@pytest.mark.skipif(not _compat.HAS_PCRASTER, reason="PCRaster not installed.") -def test_io_pcr(tmpdir): - # test write ldd with clone - da = raster.full_from_transform( - [0.5, 0.0, 3.0, 0.0, -0.5, -9.0], (4, 6), nodata=247, dtype=np.uint8, name="ldd" - ) - fn_ldd = str(tmpdir.join("test_ldd.map")) - da.raster.to_raster(fn_ldd, driver="PCRaster", pcr_vs="ldd") - assert os.path.isfile(fn_ldd) - # test ordinal - da.raster.to_raster(fn_ldd, driver="PCRaster", pcr_vs="ordinal", clone_path=fn_ldd) - assert os.path.isfile(fn_ldd) - da.expand_dims("time").raster.to_raster( - tmpdir.join("testldd.map"), driver="PCRaster" - ) - assert os.path.isfile(tmpdir.join("testldd0.001")) - ds_in = hydromt.open_mfraster( - str(tmpdir.join("testldd*")), concat=True, mask_nodata=True - ) - assert "dim0" in ds_in.coords - # mapstack - prefix = "test_" - root = str(tmpdir) - assert np.all([np.isnan(ds_in[n].raster.nodata) for n in ds_in.raster.vars]) - ds_in.raster.to_mapstack(join(root, "pcr"), prefix=prefix, driver="PCRaster")