From 8b75640c7228cf3b97fc981a092e5fe731ebc0e8 Mon Sep 17 00:00:00 2001 From: Felix Oesterle <6945681+fso42@users.noreply.github.com> Date: Tue, 28 Apr 2026 09:14:06 +0200 Subject: [PATCH 1/2] fix(geoTrans): fill bilinear remesh edge NaNs with nearest fallback - Added additional test cases for resizing and interpolating raster data in `test_geoTrans`. - Adjusted assertions in `test_deriveParameterSet` to reflect new logic for handling `nan` values. - Enhanced bilinear interpolation in `geoTrans` to fill `nan` values using nearest-neighbor interpolation. --- avaframe/in3Utils/geoTrans.py | 10 +++++++++- avaframe/tests/test_deriveParameterSet.py | 10 ++++++---- avaframe/tests/test_geoTrans.py | 14 ++++++++++++++ 3 files changed, 29 insertions(+), 5 deletions(-) diff --git a/avaframe/in3Utils/geoTrans.py b/avaframe/in3Utils/geoTrans.py index eb415640e..fb3ca157a 100644 --- a/avaframe/in3Utils/geoTrans.py +++ b/avaframe/in3Utils/geoTrans.py @@ -201,7 +201,15 @@ def resizeData(raster, rasterRef): X, Y = np.meshgrid(xgrid, ygrid) Points = {"x": X, "y": Y} Points, _ = projectOnRaster(raster, Points, interp="bilinear") - raster["rasterData"] = Points["z"] + bilinearData = Points["z"] + + if np.isnan(bilinearData).any(): + PointsNearest = {"x": X, "y": Y} + PointsNearest, _ = projectOnRaster(raster, PointsNearest, interp="nearest") + nanMask = np.isnan(bilinearData) + bilinearData[nanMask] = PointsNearest["z"][nanMask] + + raster["rasterData"] = bilinearData return raster["rasterData"], rasterRef["rasterData"] diff --git a/avaframe/tests/test_deriveParameterSet.py b/avaframe/tests/test_deriveParameterSet.py index a6d6ab749..5870da628 100644 --- a/avaframe/tests/test_deriveParameterSet.py +++ b/avaframe/tests/test_deriveParameterSet.py @@ -781,7 +781,7 @@ def test_checkExtentAndCellSize(tmp_path): testFile5, outFile5, remeshedFlag5 = dP.checkExtentAndCellSize(cfg, inputFile, dem, "rel") newRaster5 = IOf.readRaster((inDir / testFile5)) assert remeshedFlag5 == "Yes" - assert np.isnan(newRaster5["rasterData"][0, :]).all() + assert not np.isnan(newRaster5["rasterData"][0, :]).any() assert np.isnan(newRaster5["rasterData"][-1, :]).all() inputFile = inDirR / "inputFile51.asc" @@ -802,9 +802,11 @@ def test_checkExtentAndCellSize(tmp_path): inField[0:2, :] = np.nan IOf.writeResultToRaster(headerInput, inField, inputFile.parent / inputFile.stem, flip=True) - with pytest.raises(AssertionError) as e: - assert dP.checkExtentAndCellSize(cfg, inputFile, dem, "rel") - assert "nan values found inside DEM extent - this is not allowed" in str(e.value) + testFile6, outFile6, remeshedFlag6 = dP.checkExtentAndCellSize(cfg, inputFile, dem, "rel") + newRaster6 = IOf.readRaster((inDir / testFile6)) + assert remeshedFlag6 == "Yes" + assert np.isnan(newRaster6["rasterData"][-1, :]).all() + # Produced by AI (test): diff --git a/avaframe/tests/test_geoTrans.py b/avaframe/tests/test_geoTrans.py index 6f20facad..95e9b578b 100644 --- a/avaframe/tests/test_geoTrans.py +++ b/avaframe/tests/test_geoTrans.py @@ -110,6 +110,20 @@ def test_resizeData(): testRes = np.allclose(np.nan_to_num(data1 - data2), zSol, atol=tol) assert testRes + raster3 = { + "header": {"ncols": 3, "nrows": 3, "xllcenter": 0.0, "yllcenter": 0.0, "cellsize": 10.0}, + "rasterData": np.array([[1.0, 2.0, 3.0], [4.0, 5.0, 6.0], [7.0, 8.0, 9.0]]), + } + raster4 = { + "header": {"ncols": 5, "nrows": 5, "xllcenter": 0.0, "yllcenter": 0.0, "cellsize": 5.0}, + "rasterData": np.zeros((5, 5)), + } + + data3, _ = geoTrans.resizeData(raster3, raster4) + + assert np.isnan(data3).sum() == 0 + assert np.isclose(data3[2, 2], 5.0) + # def test_findAngleProfile(capfd): def test_findAngleProfile(): From 903155eef7644e2befd9d14a33fd2645e9da9cff Mon Sep 17 00:00:00 2001 From: Felix Oesterle Date: Tue, 28 Apr 2026 15:27:44 +0200 Subject: [PATCH 2/2] fix(docs, tests, com1DFA): clarify handling of no-data values and enforce NaN validation - Updated documentation to explicitly state that no-data values are not allowed in configuration descriptions for shapefiles and raster files. - Added tests to validate that NaN values are not permitted in non-DEM input rasters, with appropriate error handling. - Adjusted `deriveParameterSet` logic to raise errors for NaN values in non-DEM input files while allowing them for DEM files. --- avaframe/com1DFA/deriveParameterSet.py | 27 +++++++--- avaframe/tests/test_deriveParameterSet.py | 61 ++++++++++------------- docs/moduleCom1DFA.rst | 6 +-- 3 files changed, 49 insertions(+), 45 deletions(-) diff --git a/avaframe/com1DFA/deriveParameterSet.py b/avaframe/com1DFA/deriveParameterSet.py index 8499d6d1a..c80bcc8ca 100644 --- a/avaframe/com1DFA/deriveParameterSet.py +++ b/avaframe/com1DFA/deriveParameterSet.py @@ -426,11 +426,16 @@ def checkThicknessSettings(cfg, thName, inputSimFiles): ) log.error(message) raise AssertionError(message) - if (cfg["GENERAL"].getboolean("timeDependentRelease") and inputSimFiles["entResInfo"]["relTh" + "FileType"] in [".asc", ".tif"]): + if cfg["GENERAL"].getboolean("timeDependentRelease") and inputSimFiles["entResInfo"][ + "relTh" + "FileType" + ] in [".asc", ".tif"]: message = "When release is time dependent, release file needs to be provided as shape file" log.error(message) raise FileNotFoundError(message) - if (cfg["GENERAL"].getboolean("timeDependentRelease") and cfg["GENERAL"].getboolean("relThFromFile") is False): + if ( + cfg["GENERAL"].getboolean("timeDependentRelease") + and cfg["GENERAL"].getboolean("relThFromFile") is False + ): message = "When release is time dependent, relThFromFile needs to be set to True" log.error(message) raise FileNotFoundError(message) @@ -890,11 +895,8 @@ def appendThicknessToCfg(cfg): for count, id in enumerate(idList): thNameId = thType + id if thNameId in cfg["GENERAL"].keys(): - log.info( - "Thickness value for %s already set in initial config file, \ - read from there not from shp file" - % thNameId - ) + log.info("Thickness value for %s already set in initial config file, \ + read from there not from shp file" % thNameId) else: cfgGen[thNameId] = str(float(thicknessList[count])) @@ -968,6 +970,17 @@ def checkExtentAndCellSize(cfg, inputFile, dem, fileType): """ inputField = IOf.readRaster(inputFile) + + # Check for NaN values before remeshing - not allowed in non-DEM input rasters + if fileType.upper() != "DEM": + if np.isnan(inputField["rasterData"]).any(): + message = "In %s file (%s) nan values found - this is not allowed" % ( + fileType, + inputFile.name, + ) + log.error(message) + raise AssertionError(message) + cellSizeOld = inputField["header"]["cellsize"] demHeader = dem["header"] diff --git a/avaframe/tests/test_deriveParameterSet.py b/avaframe/tests/test_deriveParameterSet.py index 5870da628..3c1315c88 100644 --- a/avaframe/tests/test_deriveParameterSet.py +++ b/avaframe/tests/test_deriveParameterSet.py @@ -740,15 +740,30 @@ def test_checkExtentAndCellSize(tmp_path): headerInput["transform"] = IOf.transformFromASCHeader(headerInput) headerInput["crs"] = rasterio.crs.CRS() + # Test that NaN values in non-DEM input files are not allowed + inputFile = inDirR / "inputFile3.asc" + headerInput = { + "nrows": 4, + "ncols": 5, + "xllcenter": 1, + "yllcenter": 5, + "cellsize": 1, + "nodata_value": -9999, + "driver": "AAIGrid", + } + headerInput["transform"] = IOf.transformFromASCHeader(headerInput) + headerInput["crs"] = rasterio.crs.CRS() + inField = np.ones((4, 5)) inField[2, 2] = 10.0 inField[0, :] = np.nan IOf.writeResultToRaster(headerInput, inField, inputFile.parent / inputFile.stem, flip=True) - testFile4, outFile4, remeshedFlag4 = dP.checkExtentAndCellSize(cfg, inputFile, dem, "rel") - - assert remeshedFlag4 == "No" + with pytest.raises(AssertionError) as e: + dP.checkExtentAndCellSize(cfg, inputFile, dem, "rel") + assert "nan values found - this is not allowed" in str(e.value) + # Test with interior NaN (not just on border) inputFile = inDirR / "inputFile4.asc" inField = np.ones((4, 5)) inField[2, 2] = 10.0 @@ -757,15 +772,16 @@ def test_checkExtentAndCellSize(tmp_path): IOf.writeResultToRaster(headerInput, inField, inputFile.parent / inputFile.stem, flip=True) with pytest.raises(AssertionError) as e: - assert dP.checkExtentAndCellSize(cfg, inputFile, dem, "rel") - assert "nan values found inside DEM extent - this is not allowed" in str(e.value) + dP.checkExtentAndCellSize(cfg, inputFile, dem, "rel") + assert "nan values found - this is not allowed" in str(e.value) + # Test that DEM files ARE allowed to have NaN values inputFile = inDirR / "inputFile5.asc" headerInput = { "nrows": 4, "ncols": 5, - "xllcenter": 1.3, - "yllcenter": 4.2, + "xllcenter": 1, + "yllcenter": 5, "cellsize": 1, "nodata_value": -9999, "driver": "AAIGrid", @@ -778,34 +794,9 @@ def test_checkExtentAndCellSize(tmp_path): inField[0, :] = np.nan IOf.writeResultToRaster(headerInput, inField, inputFile.parent / inputFile.stem, flip=True) - testFile5, outFile5, remeshedFlag5 = dP.checkExtentAndCellSize(cfg, inputFile, dem, "rel") - newRaster5 = IOf.readRaster((inDir / testFile5)) - assert remeshedFlag5 == "Yes" - assert not np.isnan(newRaster5["rasterData"][0, :]).any() - assert np.isnan(newRaster5["rasterData"][-1, :]).all() - - inputFile = inDirR / "inputFile51.asc" - headerInput = { - "nrows": 4, - "ncols": 5, - "xllcenter": 1.3, - "yllcenter": 4.2, - "cellsize": 1, - "nodata_value": -9999, - "driver": "AAIGrid", - } - headerInput["transform"] = IOf.transformFromASCHeader(headerInput) - headerInput["crs"] = rasterio.crs.CRS() - - inField = np.ones((4, 5)) - inField[2, 2] = 10.0 - inField[0:2, :] = np.nan - IOf.writeResultToRaster(headerInput, inField, inputFile.parent / inputFile.stem, flip=True) - - testFile6, outFile6, remeshedFlag6 = dP.checkExtentAndCellSize(cfg, inputFile, dem, "rel") - newRaster6 = IOf.readRaster((inDir / testFile6)) - assert remeshedFlag6 == "Yes" - assert np.isnan(newRaster6["rasterData"][-1, :]).all() + # fileType="DEM" should allow NaNs + testFile5, outFile5, remeshedFlag5 = dP.checkExtentAndCellSize(cfg, inputFile, dem, "DEM") + assert remeshedFlag5 == "No" # Produced by AI (test): diff --git a/docs/moduleCom1DFA.rst b/docs/moduleCom1DFA.rst index 509f16d47..81b12941c 100644 --- a/docs/moduleCom1DFA.rst +++ b/docs/moduleCom1DFA.rst @@ -63,7 +63,7 @@ or GeoTIFF format, or shape files and specified for the respective input type be - or raster file(s): - cells with non-zero values define the release area - cell value can be read as thickness (measured normal to the slope) (see :ref:`moduleCom1DFA:Release-, entrainment thickness settings`) - - negative values are not allowed (specified no-data values are not considered) + - negative values and no-data values are not allowed - ALL features within one shapefile or all cells with non-zero values in a raster file are released at the same time (and interact), this is what we refer to as *scenario* - if you want to simulate different scenarios with the same features, you have to copy them to separate shapefiles/raster files - the release area name should not contain an underscore, if so '_AF' is added. @@ -80,7 +80,7 @@ at least two results are generated: the *null* variant and the variant with entr - or raster file: - cells with non-zero values define where entrainment can occur - cell value can be read as thickness (measured normal to the slope) (see :ref:`moduleCom1DFA:Release-, entrainment thickness settings`) - - negative values are not allowed (specified no-data values are not considered) + - negative values and no-data values are not allowed * **one resistance area as (multi-) polygon shapefile OR raster file (in Inputs/RES)** - marks the (multiple) areas where resistance is considered @@ -89,7 +89,7 @@ at least two results are generated: the *null* variant and the variant with entr - resistance areas must not contain any "holes" or inner rings - or raster file: - cells with non-zero values define where entrainment can occur - - negative values are not allowed (specified no-data values are not considered) + - negative values and no-data values are not allowed - if the cellsize does not match the requested meshCellSize, the file is remeshed if within `resizeThreshold`