refactor: migrate to pyramids-gis 0.32 and cleopatra 0.17 APIs

.gitignore

@@ -36,6 +36,9 @@ MANIFEST
 pip-log.txt
 pip-delete-this-directory.txt
 
+# GDAL statistics sidecar files generated when tests compute raster stats
+*.aux.xml
+
 # Unit test / coverage reports
 htmlcov/
 .tox/

examples/hydrological-model/HAND_HRUpy.py

@@ -2,6 +2,7 @@
 
 @author: Mostafa
 """
+
 # %links
 # from IPython import get_ipython  # to reset the variable explorer each time
 
@@ -18,40 +19,40 @@
 
 # %library
 import numpy as np
-from osgeo import gdal
+from pyramids.dataset import Dataset
 from scipy.stats import norm
 
 # functions
-from Hapi.rrm.parameters import Parameters as DP
+from hapi.rrm.parameters import Parameters as DP
 
-# from Hapi.catchment import GISCatchment as gc
+# from hapi.catchment import GISCatchment as gc
 
 # import matplotlib.pyplot as plt
 
 
 # %% to modify the basins raster
 # path="C:/Users/Mostafa/Desktop/delineation/Clipped/proj/basins.tif"
 # pathout="mask.tif"
-# basins=gdal.Open(path)
+# basins=Dataset.read_file(path)
 # gc.DeleteBasins(basins,pathout)
 # %%
 path = "C:/Users/Mostafa/Desktop/delineation/HRU/HAND/"
-DEM = gdal.Open(path + "DEM.tif")
+DEM = Dataset.read_file(path + "DEM.tif")
 # dem_A=DEM.ReadAsArray()
 # no_val=np.float32(DEM.GetRasterBand(1).GetNoDataValue())
 # rows=DEM.RasterYSize
 # cols=DEM.RasterXSize
-FD = gdal.Open(path + "FD.tif")
+FD = Dataset.read_file(path + "FD.tif")
 # fd_A=fd.ReadAsArray()
 # fd_val=[int(fd_A[i,j]) for i in range(rows) for j in range(cols) if fd_A[i,j] != no_val]
-FPL = gdal.Open(path + "fpl.tif")
+FPL = Dataset.read_file(path + "fpl.tif")
 # fpl_A=FPL.ReadAsArray()
-River = gdal.Open(path + "river.tif")
+River = Dataset.read_file(path + "river.tif")
 # river_A=River.ReadAsArray()
 # slope
-Slope = gdal.Open(path + "slope.tif")
-Slope_A = Slope.ReadAsArray()
-no_val_slope = np.float32(Slope.GetRasterBand(1).GetNoDataValue())
+Slope = Dataset.read_file(path + "slope.tif")
+Slope_A = Slope.read_array(band=0)
+no_val_slope = np.float32(Slope.no_data_value[0])
 # %% calculate HAND and DTND
 HAND, DTND = DP.hru_hand(DEM, FD, FPL, River)
 # %% calculate the cdf for

examples/hydrological-model/coello/calibration/coello-hrus-distributed-model-calibratation-muskingum.py

@@ -5,12 +5,12 @@
 import datetime as dt
 
 import numpy as np
-from osgeo import gdal
+from pyramids.dataset import Dataset
 from statista.descriptors import rmse
 
-from Hapi.calibration import Calibration
-from Hapi.rrm.hbv_bergestrom92 import HBVBergestrom92 as HBV
-from Hapi.rrm.parameters import Parameters as DP
+from hapi.calibration import Calibration
+from hapi.rrm.hbv_bergestrom92 import HBVBergestrom92 as HBV
+from hapi.rrm.parameters import Parameters as DP
 
 path = Comp + "/Coello/HAPI/Data/00inputs/"  # GIS/4000/
 SaveTo = Comp + "/Coello/Hapi/Model/results/"
@@ -61,7 +61,7 @@
 "01_rfcf","02_FC", "03_BETA", "04_ETF", "05_LP", "06_K0","07_K1", "08_K2",
 "09_UZL","10_PERC", + ["11_Kmuskingum", "12_Xmuskingum"]
 """
-raster = gdal.Open(path + "GIS/4000/subcatch_classes.tif")
+raster = Dataset.read_file(path + "GIS/4000/subcatch_classes.tif")
 no_parameters = 12
 no_lumped_par = 1
 lumped_par_pos = [7]

examples/hydrological-model/coello/calibration/coello-totally-distributed-model-calibratation-muskingum.py

@@ -1,12 +1,12 @@
 import datetime as dt
 
 import numpy as np
-from osgeo import gdal
+from pyramids.dataset import Dataset
 from statista.descriptors import rmse
 
-from Hapi.calibration import Calibration
-from Hapi.rrm.hbv_bergestrom92 import HBVBergestrom92 as HBV
-from Hapi.rrm.parameters import Parameters as DP
+from hapi.calibration import Calibration
+from hapi.rrm.hbv_bergestrom92 import HBVBergestrom92 as HBV
+from hapi.rrm.parameters import Parameters as DP
 
 # %% Paths
 Path = "examples/hydrological-model/data/distributed_model/"
@@ -54,7 +54,7 @@
 
 par_dist = SpatialVarFun(par,*SpatialVarArgs,kub=kub,klb=klb)
 """
-raster = gdal.Open(FlowAccPath)
+raster = Dataset.read_file(FlowAccPath)
 # -------------
 # for lumped catchment parameters
 no_parameters = 12

examples/hydrological-model/coello/prepare-Input-data/01-coello-prepare_inputdata.py

@@ -1,43 +1,47 @@
+"""Prepare the Coello meteorological and GIS input rasters.
+
+Make sure the working directory is set to the root of the Hapi repo:
+current_work_directory = Hapi/
 """
-Make sure the working directory is set to the examples folder in the Hapi repo"
-currunt_work_directory = Hapi/
-"""
+
+from __future__ import annotations
+
+from pathlib import Path
 
 import numpy as np
-from osgeo import gdal
-from pyramids.raster import Raster
+from pyramids.dataset import Dataset
 
-from Hapi.rrm.inputs import Inputs
+from hapi.inputs import Inputs
 
 metao_data_path = "examples/hydrological-model/data/meteo_data"
 gis_data_path = "examples/hydrological-model/data/gis_data"
 """
-this function prepare downloaded raster data to have the same align and
-nodatavalue from a GIS raster (DEM, flow accumulation, flow direction raster)
-and return a folder with the output rasters with a name “New_Rasters”
+prepare_inputs aligns and crops the downloaded raster data to have the same
+alignment and no-data value as a GIS raster (DEM, flow accumulation, or flow
+direction raster) and writes the prepared rasters to the given output folder.
 """
 
 dem_path = f"{gis_data_path}/acc4000.tif"
 
-In = Inputs(dem_path)
+inputs = Inputs(dem_path)
 outputpath = f"{metao_data_path}/meteodata_prepared/"
 
 # prec
 prec_in_path = f"{metao_data_path}/raw_data/prec/"
-In.prepareInputs(prec_in_path, outputpath + "prec0")
+inputs.prepare_inputs(prec_in_path, outputpath + "prec0")
 
 # evap
 evap_in_path = f"{metao_data_path}/raw_data/evap/"
-In.prepareInputs(evap_in_path, outputpath + "evap0")
+inputs.prepare_inputs(evap_in_path, outputpath + "evap0")
 # temp
 temp_in_path = f"{metao_data_path}/raw_data/temp/"
-In.prepareInputs(temp_in_path, outputpath + "temp0")
+inputs.prepare_inputs(temp_in_path, outputpath + "temp0")
 
 """
-in case you want to manipulate the value in all the rasters of one of the inputs
-for example evapotranspiration values in rasters downloaded from ECMWF are -ve
-and to change it to +ve in all rasters or if you want to operate any kind of function
-in all input rasters that are in the same folder FolderCalculator can do this task
+in case you want to manipulate the values in all the rasters of one of the inputs,
+for example evapotranspiration values in rasters downloaded from ECMWF are -ve,
+and to change them to +ve in all rasters (or apply any kind of function to all
+input rasters in the same folder), loop over the files and use `Dataset.apply`.
 
 "How can evaporation have both positive and negative values?
 Evaporation is normally negative due to the convention for fluxes.
@@ -46,50 +50,31 @@
 Link: https://confluence.ecmwf.int/pages/viewpage.action?pageId=111155327
 """
 evap_out_path = f"{metao_data_path}/meteodata_prepared/evap/"
-
-
-# define your function
-# this function is going to take the absolute value of the values in the raster
-# through MapAlgebra function then save the new raster to a given path with the same names
-def function(args):
-    # argument a list of two components
-    # first argument is the raster object [gdal object]
-    A = args[0]
-    # second argument is the path to save the resulted raster
-    path = args[1]
-    func = np.abs
-    # first function
-    B = Raster.mapAlgebra(A, func)
-    Raster.saveRaster(B, path)
-
-
-folder_path = evap_out_path
 new_folder_path = f"{metao_data_path}/meteodata_prepared/new_evap/"
-Raster.folderCalculator(folder_path, new_folder_path, function)
+
+Path(new_folder_path).mkdir(parents=True, exist_ok=True)
+for file in sorted(Path(evap_out_path).glob("*.tif")):
+    dataset = Dataset.read_file(str(file))
+    # apply operates on the domain cells only and keeps the no-data cells
+    dataset.apply(np.abs).to_file(f"{new_folder_path}/{file.name}")
 
 """
-in order to run the model all inputs have to have the same number of rows and columns
-for this purpose MatchRasterAlignment function was made to resample, change the coordinate
-system of the second raster and give it the same alignment like a source raster (DEM raster)
+in order to run the model all inputs have to have the same number of rows and columns.
+`Dataset.align` resamples a raster and copies the coordinate system, the number of
+rows/columns, and the cell size from a source raster (the DEM raster).
 """
 
 soil_path = f"{gis_data_path}/soil_classes.tif"
-DEM = gdal.Open(dem_path)
-dem_A = DEM.ReadAsArray()
-soil = gdal.Open(soil_path)
-soil_A = soil.ReadAsArray()
-
-# align
-aligned_soil = Raster.matchRasterAlignment(DEM, soil)
+dem = Dataset.read_file(dem_path)
+soil = Dataset.read_file(soil_path)
 
-# to check alignment of DEM raster compared to aligned_soil_A raster
-aligned_soil_A = aligned_soil.ReadAsArray()
+# align the soil raster to the DEM alignment
+aligned_soil = soil.align(dem)
 
-# nodatavalue is still different and some cells are no data value in the soil type raster but it is not in the dem raster
-# to match use Match MatchNoDataValue
-# match
-dst_Aligned_M = Raster.cropAlligned(DEM, aligned_soil)
-dst_Aligned_M_A = dst_Aligned_M.ReadAsArray()
+# the no-data cells are still different: some cells are no-data in the soil
+# raster but not in the DEM raster. Crop with the DEM as a raster mask so both
+# rasters share the same no-data cells.
+masked_soil = aligned_soil.crop(dem)
 
 # save the new raster
-Raster.saveRaster(dst_Aligned_M, f"{gis_data_path}/soil_classes_aligned.tif")
+masked_soil.to_file(f"{gis_data_path}/soil_classes_aligned.tif")

examples/hydrological-model/coello/prepare-Input-data/07-prepare-rasters.py

@@ -6,18 +6,19 @@
 currunt_work_directory = Hapi/Example
 """
 
-from osgeo import gdal
-from pyramids.raster import Raster
+from __future__ import annotations
+
+from pyramids.dataset import Dataset
 
 dem_path = "Data/GIS/Hapi_GIS_Data/acc4000.tif"
 SaveTo = "data/parameters/"
 # %%
-"craeate a raster typicall to the DEM and fill it with 1"
-K = 1
-src = gdal.Open(dem_path)
+# create a raster typical to the DEM and fill its domain cells with one value
+dem = Dataset.read_file(dem_path)
 
-Raster.rasterFill(src, K, SaveTo + "11_K_muskingum.tif")
+K = 1
+dem.fill(K, path=SaveTo + "11_K_muskingum.tif")
 # %%
 X = 0.2
 
-Raster.rasterFill(src, X, SaveTo + "12_X_muskingum.tif")
+dem.fill(X, path=SaveTo + "12_X_muskingum.tif")

examples/hydrological-model/coello/run/coello-distributed-model-run-maxbas.py

@@ -3,11 +3,11 @@
 import datetime as dt
 
 import pandas as pd
-from osgeo import gdal
+from pyramids.dataset import Dataset
 
-from Hapi.catchment import Catchment
-from Hapi.rrm.hbv_bergestrom92 import HBVBergestrom92 as HBV
-from Hapi.run import Run
+from hapi.catchment import Catchment
+from hapi.rrm.hbv_bergestrom92 import HBVBergestrom92 as HBV
+from hapi.run import Run
 
 # %% Paths
 Path = "examples/hydrological-model/data/distributed_model"
@@ -78,7 +78,7 @@
 Coello.plot_hydrograph(plotstart, plotend, gaugei)
 # %% store the result into rasters
 # create list of names
-src = gdal.Open(FlowAccPath)
+src = Dataset.read_file(FlowAccPath)
 s = dt.datetime(2012, 6, 14, 19, 00, 00)
 e = dt.datetime(2013, 12, 23, 00, 00, 00)
 index = pd.date_range(s, e, freq="1H")