[com1],[doc]: correct docu for raster inputs

avaframe/com1DFA/com1DFACfg.ini

@@ -238,7 +238,7 @@ meshCellSize = 5
 meshCellSizeThreshold = 0.001
 # clean DEMremeshed directory to ensure remeshing if chosen meshCellsize is different from rasters in Inputs/
 cleanRemeshedRasters = True
-# resize files read from Inputs/RASTERS to be resized to extent of DEM as resizeThreshold x meshCellSize
+# remesh if extent matches computational DEM within resizeThreshold x meshCellSize
 resizeThreshold = 3
 
 # Normal computation on rectangular grid

avaframe/com1DFA/deriveParameterSet.py

@@ -948,15 +948,8 @@ def checkExtentAndCellSize(cfg, inputFile, dem, fileType):
     cellSizeOld = inputField["header"]["cellsize"]
     demHeader = dem["header"]
 
-    # check if negative values or nan values
-    if np.any(np.isnan(inputField["rasterData"])):
-        message = "In %s file (%s) nan values found - this is not allowed" % (
-            fileType,
-            inputFile.name,
-        )
-        log.error(message)
-        raise AssertionError(message)
-    elif np.any(inputField["rasterData"] < 0):
+    # check if negative values in raster file
+    if np.any(inputField["rasterData"] < 0):
         message = "In %s file (%s) negative values found - this is not allowed" % (
             fileType,
             inputFile.name,
@@ -1025,6 +1018,36 @@ def checkExtentAndCellSize(cfg, inputFile, dem, fileType):
         returnStr = str(pathlib.Path("remeshedRasters", outFile.name))
         remeshedFlag = "Yes"
 
+    # check if no data values only where DEM also has no data values
+    # if remeshed - potentially nans at edges, if inputField has a different origin/extent
+    # for example if extent of DEM is larger -> nans in this region in remeshed input file
+    # first maks nans values of DEM as there, also inputField is allowed to have nans
+    nanDEMMasked = np.where(np.isnan(dem["rasterData"]), -9999, inputField["rasterData"])
+    # search for indices where nans come from remeshing use difference in extent prior to remeshing
+    # if diff is negative on Left side DEM is smaller
+    # if diff is negative on right side DEM is larger
+    nNeglectRowsMin = int(np.ceil(abs(min(0, diffY0 / dem["header"]["cellsize"]))))
+    nNeglectRowsMax = int(np.ceil(abs(max(0, diffY1 / dem["header"]["cellsize"]))))
+    # if diff is negative on lower side DEM is smaller
+    # if diff is negative on right side DEM is larger
+    nNeglectColsMin = int(np.ceil(abs(min(0, diffX0 / dem["header"]["cellsize"]))))
+    nNeglectColsMax = int(np.ceil(abs(max(0, diffX1 / dem["header"]["cellsize"]))))
+    maxRows = dem["header"]["nrows"] - nNeglectRowsMin
+    maxCols = dem["header"]["ncols"] - nNeglectColsMin
+
+    # add mask where nans come from remeshing
+    # change order because diff is with origin lower!!
+    nanDEMMaskedLimited = nanDEMMasked[nNeglectRowsMax:maxRows, nNeglectColsMax:maxCols]
+
+    # check if nan values in raster file data (excluding nans from remeshing and where also DEM has nans)
+    if np.any(np.isnan(nanDEMMaskedLimited)):
+        message = "In %s file (%s) nan values found inside DEM extent - this is not allowed" % (
+            fileType,
+            inputFile.name,
+        )
+        log.error(message)
+        raise AssertionError(message)
+
     return returnStr, outFile, remeshedFlag
 
 

avaframe/tests/test_deriveParameterSet.py

@@ -686,7 +686,7 @@ def test_checkExtentAndCellSize(tmp_path):
     assert newRaster2["header"]["xllcenter"] == 1.0
     assert newRaster2["header"]["yllcenter"] == 5.0
 
-    inputFile2 = inDirR / "inputFile1.asc"
+    inputFile2 = inDirR / "inputFile2.asc"
     headerInput2 = {
         "nrows": 5,
         "ncols": 5,
@@ -706,6 +706,102 @@ def test_checkExtentAndCellSize(tmp_path):
         assert dP.checkExtentAndCellSize(cfg, inputFile2, dem, "mu")
     assert "Lower left center coordinates of DEM and " in str(e.value)
 
+    demField = np.ones((4, 5))
+    dem = {
+        "header": {
+            "nrows": 4,
+            "ncols": 5,
+            "xllcenter": 1,
+            "yllcenter": 5,
+            "cellsize": 1,
+            "nodata_value": -9999,
+            "driver": "AAIGrid",
+        },
+        "rasterData": demField,
+    }
+
+    demField[0, :] = np.nan
+    dem["header"]["transform"] = IOf.transformFromASCHeader(dem["header"])
+    dem["header"]["crs"] = rasterio.crs.CRS()
+
+    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"
+
+    inputFile = inDirR / "inputFile4.asc"
+    inField = np.ones((4, 5))
+    inField[2, 2] = 10.0
+    inField[0, :] = np.nan
+    inField[1, 1] = 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)
+
+    inputFile = inDirR / "inputFile5.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, :] = 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 np.isnan(newRaster5["rasterData"][0, :]).all()
+    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)
+
+    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)
 
 # Produced by AI (test):
 

docs/moduleCom1DFA.rst

@@ -48,55 +48,66 @@ folder structure described below.
 
 
 In the directory ``Inputs``, the following files are required. Be aware that ALL inputs have to be provided in the same
-projection:
+projection. Allowed file types are either raster files, i.e. `ESRI grid format <https://desktop.arcgis.com/en/arcmap/10.3/manage-data/raster-and-images/esri-ascii-raster-format.htm>`_
+or GeoTIFF format, or shape files and specified for the respective input type below:
 
-* digital elevation model as raster file with either `ESRI grid format <https://desktop.arcgis.com/en/arcmap/10.3/manage-data/raster-and-images/esri-ascii-raster-format.htm>`_
-  or GeoTIFF format. The format of the DEM determines the format of the output files.
 
-* release area scenario as (multi-) polygon shapefile (in Inputs/REL; multiple features are possible)
-
-  - the release area polygon must not contain any "holes" or inner rings
-  - the release area name should not contain an underscore, if so '_AF' is added.
-  - recommended attributes are *name*, *thickness* (see :ref:`moduleCom1DFA:Release-, entrainment thickness settings`)
-    and *ci95* (see :ref:`moduleAna4Stats:probAna - Probability maps`)
-  - ALL features within one shapefile 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
+* **digital elevation model as raster file. The format of the DEM determines the format of the output files.**
 
+* release area scenario as (multi-) polygon shapefile OR raster file (in Inputs/REL; only shapefiles OR raster files)
+    - either polygon shapefile(s):
+        - the release area polygon must not contain any "holes" or inner rings
+        - multiple features are allowed
+        - recommended attributes are *name*, *thickness* (see :ref:`moduleCom1DFA:Release-, entrainment thickness settings`) and *ci95* (see :ref:`moduleAna4Stats:probAna - Probability maps`)
+    - 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)
+    - 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.
 
 and the following files are optional. Please note: in the standard configuration (i.e. ``simTypeList = available``) ,
 the *null* variant is always run! I.e. if a resistance and/or an entrainment file is given (as described below),
 at least two results are generated: the *null* variant and the variant with entrainment and/or resistance.
 
-* one entrainment area (multi-) polygon shapefile (in Inputs/ENT)
-
-  - marks the (multiple) areas where entrainment can occur.
-  - attribute *thickness* (see :ref:`moduleCom1DFA:Release-, entrainment thickness settings`)
-  - must not contain any "holes" or inner rings
-
-
-* one resistance area (multi-) polygon shapefile (in Inputs/RES)
-
-  - marks the (multiple) areas where resistance is considered
-  - resistance areas must not contain any "holes" or inner rings
-  - please consider the information about resistance below :ref:`moduleCom1DFA:Resistance setup`
-
-
-* one secondary release area (multi-) polygon shapefile (in Inputs/SECREL)
+* **one entrainment area as (multi-) polygon shapefile OR raster file (in Inputs/ENT)**
+    - marks the (multiple) areas where entrainment can occur.
+    - either polygon shapefile:
+        - attribute *thickness* (see :ref:`moduleCom1DFA:Release-, entrainment thickness settings`)
+        - must not contain any "holes" or inner rings
+    - 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)
+
+* **one resistance area as (multi-) polygon shapefile OR raster file (in Inputs/RES)**
+    - marks the (multiple) areas where resistance is considered
+    - please consider the information about resistance below :ref:`moduleCom1DFA:Resistance setup`
+    - either polygon shapefile:
+        - 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)
+        - if the cellsize does not match the requested meshCellSize, the file is
+          remeshed if within `resizeThreshold`
+
+* **one secondary release area (multi-) as polygon shapefile OR raster file (in Inputs/SECREL)**
 
   - can have multiple release areas, each as one feature
   - same setup as the release area scenario (see above)
   - features will release as soon as at least one particle enters its area
   - release area polygons must not contain any "holes" or inner rings
 
-* raster files for the Voellmy friction parameters :math:`\mu` and :math:`\xi` (in Inputs/RASTERS)
+* **raster files for the Voellmy friction parameters :math:`\mu` and :math:`\xi` (in Inputs/RASTERS)**
 
   - spatial field of :math:`\mu` and :math:`\xi` values with same extent as DEM
   - file names need to end with ``_mu.*`` and ``_xi.*``
   - only one file per parameter allowed
   - if ``meshCellSize`` is different from simulation ``meshCellSize`` fields will be remeshed
   - only used if ``frictionModel`` is set to ``spatialVoellmy``
 
-* one ``_cropshape.shp`` shape file (in Inputs/POLYGONS)
+* **one ``_cropshape.shp`` shape file (in Inputs/POLYGONS)**
 
   - provides a polygon located inside the DEM to define area for report plots of peak fields (bounds of polygon)
   - if not provided peak fields are shown for the extent where peak field values are nonzero
@@ -110,53 +121,51 @@ Release-, entrainment thickness settings
 .. Note::
     Thickness is unambiguous: it is measured normal to the slope.
 
-Release, entrainment and secondary release thickness can be specified in two different ways:
-
-1. Via **shape file**:
-
-  - add an attribute called `thickness` for each feature
-  - important: ALL features have to have a single thickness value, which can differ between features
-  - for entrainment area only: if the thickness value is missing, the thickness value is
-    taken from `entThIfMissingInShp` (default 0.3 m) in the configuration file. If multiple features are
-    in the entrainment file the thickness attribute has to be set either for ALL or NONE of the features.
-  - for backwards compatibility, the attribute 'd0' also works, but we suggest to use `thickness` in new projects
-  - set the flag `THICKNESSFromShp` (i.e. relThFromFile, entThFromFile,
-    secondaryRelthFromShp) to True in the configuration file (default is True)
-  - a parameter variation can be added with the `THICKNESSPercentVariation`
-    parameter in the configuration file in the form of
-    ``+-percentage$numberOfSteps``. Provided a `+` a positive variation will be
-    performed, if `-` is given, only a negative variation is performed. If no
-    sign is given: both directions will be used. Additionally, a variation can be
-    added with the `THICKNESSRangeVariation` parameter in the configuration file
-    in the form of ``+-range$numberOfSteps``. Provided a `+` a positive variation
-    will be performed, if `-` is given, only a negative variation is performed.
-    If no sign is given: both directions will be used. Furthermore, there is the
-    option to vary the thickness in a range of +- the 95% confidence interval
-    value, which is also read from the shape file (requires an attribute called
-    ci95). In order to use this variation, set the 'THICKESSRangeFromCiVariation'
-    to ``ci95$numberOfSteps``.
+Release, entrainment and secondary release thickness can be specified in two different ways, 1) directly from the provided
+input file (shape file or raster file) or 2) through the :py:mod:`com1DFA` configuration file:
+
+1. **Read from input file**:
+
+- Via **shape file**:
+
+    - add an attribute called `thickness` for each feature
+    - important: ALL features have to have a single thickness value, which can differ between features
+    - for entrainment area only: if the thickness value is missing, the thickness value is
+      taken from `entThIfMissingInShp` (default 0.3 m) in the configuration file. If multiple features are
+      in the entrainment file the thickness attribute has to be set either for ALL or NONE of the features.
+    - for backwards compatibility, the attribute 'd0' also works, but we suggest to use `thickness` in new projects
+    - set the flag `THICKNESSFromFile` (i.e. relThFromFile, entThFromFile,
+      secondaryRelThFromFile) to True in the configuration file (default is True)
+    - a parameter variation can be added with the `THICKNESSPercentVariation`
+      parameter in the configuration file in the form of
+      ``+-percentage$numberOfSteps``. Provided a `+` a positive variation will be
+      performed, if `-` is given, only a negative variation is performed. If no
+      sign is given: both directions will be used. Additionally, a variation can be
+      added with the `THICKNESSRangeVariation` parameter in the configuration file
+      in the form of ``+-range$numberOfSteps``. Provided a `+` a positive variation
+      will be performed, if `-` is given, only a negative variation is performed.
+      If no sign is given: both directions will be used. Furthermore, there is the
+      option to vary the thickness in a range of +- the 95% confidence interval
+      value, which is also read from the shape file (requires an attribute called
+      ci95). In order to use this variation, set the 'THICKESSRangeFromCiVariation'
+      to ``ci95$numberOfSteps``.
+
+- Via **raster file**:
+
+    - non-zero cell values are read as thickness
+    - set the flag `THICKNESSFromFile` (i.e. relThFromFile, entThFromFile,
+      secondaryRelThFromFile) to True in the configuration file (default is True)
+    - no thickness variation is available if input type is raster file and `THICKNESSFromFile` is *True*
 
 2. Via **configuration file (ini)**:
 
-  - set the flag 'THICKNESSFromShp' to False
-  - provide your desired thickness value in the respective THICKNESS parameter (i.e. relTh, entTh or secondaryRelth)
-  - in addition to the `THICKNESSPercentVariation` and `THICKNESSRangeVariation`
-    options (see option 1) and the standard variation options in
-    :ref:`configuration:Configuration`, you can also directly set e.g. `relTh =
-    1.$50$2`, ``referenceValue$+-percentage$numberOfSteps``, resulting in a
-    variation of relTh from 0.5 to 1.5m in two steps.
-
-Only available for release thickness:
-
-3. Via **release thickness file**:
-
-  - set the flag 'relThFromFile' to False
-  - set the flag 'relThFromFile' to True
-  - save a raster file with info on release thickness as raster file in
-    ``Inputs/RELTH`` the number of rows and columns must match the DEM raster
-    with desired meshCellSize (recommended)
-  - if the cellsize does not match the requested meshCellSize, the file is
-    remeshed.
+    - set the flag 'THICKNESSFromFile' to False
+    - provide your desired thickness value in the respective THICKNESS parameter (i.e. relTh, entTh or secondaryRelth)
+    - in addition to the `THICKNESSPercentVariation` and `THICKNESSRangeVariation`
+      options (see option 1) and the standard variation options in
+      :ref:`configuration:Configuration`, you can also directly set e.g. `relTh =
+      1.$50$2`, ``referenceValue$+-percentage$numberOfSteps``, resulting in a
+      variation of relTh from 0.5 to 1.5m in two steps.
 
 
 Friction parameters

pyproject.toml

@@ -32,7 +32,7 @@ dependencies = [
     "scipy",
     "cmcrameri",
     "seaborn",
-    "cython",
+    "cython<3.0", # cython 3.x compiles code that needs numpy > 2.0
     "pandas",
     "shapely",
     "configUpdater",