[in2TopoHyd]: Module for computing hydraulic inital conditions

debrisframe/c1TIF/c1TIFCfg.ini

@@ -3,8 +3,11 @@
 
 
 [GENERAL]
-# so far no general parameters to be set as debris flow option is fully implemented in com1DFA
-# hence checkout the com1DFA_override section
+
+# if inputHydrograph is True, the in2TopoHyd module is computed before the c1TIF computation
+# Hydrograph input data need to be provided and the flag timeDependentRelease needs to be set to TRue
+inputHydrograph = False
+
 
 [com1DFA_com1DFA_override]
 
@@ -23,7 +26,7 @@ simTypeList = available
 modelType = dfa
 
 #+++++++++++++ Output++++++++++++
-# desired result Parameters (ppr, pft, pfv, pta, FT, FV, P, FM, Vx, Vy, Vz, TA, particles) - separated by |
+# desired result Parameters (ppr, pft, pfv, pta, FT, FV, P, FM, Vx, Vy, Vz, TA, dmDet, sfcChange, demAdapted, timeInfo, entrDepth, particles) - separated by |
 resType = pft|pfv|ppr
 
 #+++++++++Debris flow properties
@@ -38,6 +41,11 @@ rhoEnt = 2000
 # mass per particle (if MPPDIR is used) [kg]
 massPerPart = 280000.
 
+#+++++++++++++Mesh
+# remesh the input rasters or look for remeshed rasters
+# expected mesh size [m], default in AvaFrame: 5 m
+meshCellSize = 2
+
 #+++++++++++++general start conditions: time dependent release
 # if timeDependentRelease is True, provide the the timesteps, thickness and velocity
 # for a releases in a csv-table in the REL folder
@@ -78,6 +86,7 @@ xsivoellmy = 200
 #+++++++++++++Coulomb friction model
 mucoulomb = 0.47
 #+++++++++DEBRIS FLOW properties
+# only taken into account when using the O'Brien and Julien friction model
 # density of sediment [kg/m³]
 rhoSediment = 2650
 # mean sediment size [m]
@@ -141,7 +150,37 @@ entEroEnergy = 5000
 adaptSfcStopped = 0
 adaptSfcDetrainment = 0
 adaptSfcEntrainment = 0
+# deposition is only entrainable if it is an entrainment simulation and if adaptSfcStopped = 1 and adaptSfcEntrainment = 1
+# then mass that is stopped can be entrained (also if the entrainment area is not in this stopped area)
+entrainableDeposition = False
 # only adapt topography if changing height in at least one cell is > thresholdAdaptSfc [m]
 thresholdAdaptSfc = 0.1
 # use the adapted topography as background in the report plots
 adaptDemPlot = False
+
+# specify a particular release area scenario, provide name of shapefile with or without extension .shp (optional)
+releaseScenario =
+
+
+[in2TopoHyd_in2TopoHyd_override]
+# use default in2TopoHyd config as base configuration (True) and override following parameters
+# if False and local_in2TopoHydCfg is available use local
+defaultConfig = True
+
+#++++++++++++++++ Hydraulic Boundary Conditions
+# Elevation steps in vertical direction to calculate the flow area
+dElev = 0.1
+# method how the flow velocity as initial condition is calculated
+# currently only the method after Rickenmann (1999) is available
+velType = rickenmann
+# mean slope for the cross section in flow direction [m/m]
+# optional; if left blank, the slope is computed automatically
+slope =
+# normal distance to release line where two parallel cross sections are generated on each side of the release line
+# needed for computing the flow direction and the channel slope
+# optional; if left blank, the normal distance is computed automatically
+normalDist = 
+
+# export cell centers along the release line as points in a csv-file; for plausability check
+# True or False
+exportCrossSectionCells = True

debrisframe/data/debrisParabolaChanneled/Inputs/HYDR/hydrograph.csv

@@ -0,0 +1,92 @@
+timestep,discharge
+0,0.0
+1,2.83
+2,5.67
+3,8.5
+4,11.33
+5,14.17
+6,17.0
+7,19.83
+8,22.67
+9,25.5
+10,28.33
+11,31.17
+12,34.0
+13,36.83
+14,39.67
+15,42.5
+16,45.33
+17,48.17
+18,51.0
+19,53.83
+20,56.67
+21,59.5
+22,62.33
+23,65.17
+24,68.0
+25,70.83
+26,73.67
+27,76.5
+28,79.33
+29,82.17
+30,85.0
+31,83.58
+32,82.17
+33,80.75
+34,79.33
+35,77.92
+36,76.5
+37,75.08
+38,73.67
+39,72.25
+40,70.83
+41,69.42
+42,68.0
+43,66.58
+44,65.17
+45,63.75
+46,62.33
+47,60.92
+48,59.5
+49,58.08
+50,56.67
+51,55.25
+52,53.83
+53,52.42
+54,51.0
+55,49.58
+56,48.17
+57,46.75
+58,45.33
+59,43.92
+60,42.5
+61,41.08
+62,39.67
+63,38.25
+64,36.83
+65,35.42
+66,34.0
+67,32.58
+68,31.17
+69,29.75
+70,28.33
+71,26.92
+72,25.5
+73,24.08
+74,22.67
+75,21.25
+76,19.83
+77,18.42
+78,17.0
+79,15.58
+80,14.17
+81,12.75
+82,11.33
+83,9.92
+84,8.5
+85,7.08
+86,5.67
+87,4.25
+88,2.83
+89,1.42
+90,0.0

debrisframe/data/debrisParabolaChanneled/Inputs/LEVEE/crossSectionLevee.cpg

@@ -0,0 +1 @@
+UTF-8

debrisframe/data/debrisParabolaChanneled/Inputs/LEVEE/crossSectionLevee.prj

@@ -0,0 +1 @@
+PROJCS["ETRS_1989_Austria_Lambert",GEOGCS["GCS_ETRS_1989",DATUM["D_ETRS_1989",SPHEROID["GRS_1980",6378137.0,298.257222101]],PRIMEM["Greenwich",0.0],UNIT["Degree",0.0174532925199433]],PROJECTION["Lambert_Conformal_Conic"],PARAMETER["False_Easting",400000.0],PARAMETER["False_Northing",400000.0],PARAMETER["Central_Meridian",13.3333333333333],PARAMETER["Standard_Parallel_1",49.0],PARAMETER["Standard_Parallel_2",46.0],PARAMETER["Latitude_Of_Origin",47.5],UNIT["Meter",1.0]]

debrisframe/data/debrisParabolaChanneled/Inputs/XSECT/crossSection.cpg

@@ -0,0 +1 @@
+UTF-8

debrisframe/data/debrisParabolaChanneled/Inputs/XSECT/crossSection.prj

@@ -0,0 +1 @@
+PROJCS["ETRS_1989_Austria_Lambert",GEOGCS["GCS_ETRS_1989",DATUM["D_ETRS_1989",SPHEROID["GRS_1980",6378137.0,298.257222101]],PRIMEM["Greenwich",0.0],UNIT["Degree",0.0174532925199433]],PROJECTION["Lambert_Conformal_Conic"],PARAMETER["False_Easting",400000.0],PARAMETER["False_Northing",400000.0],PARAMETER["Central_Meridian",13.3333333333333],PARAMETER["Standard_Parallel_1",49.0],PARAMETER["Standard_Parallel_2",46.0],PARAMETER["Latitude_Of_Origin",47.5],UNIT["Meter",1.0]]

debrisframe/debrisframeCfg.ini

@@ -5,7 +5,7 @@
 [MAIN]
 # Path to avalanche directory
 #TODO: change avalancheDir to debrisDir
-avalancheDir = data/debrisGenericTopoTimeDep
+avalancheDir = data/debrisParabolaChanneled
 
 # number of CPU cores to use for the computation of com1DFA
 # possible values are:

debrisframe/in1Utils/fileUtils.py

@@ -0,0 +1,43 @@
+"""
+Directory and file handling helper functions
+"""
+
+import pathlib
+import shutil
+import logging
+import pandas as pd
+
+log = logging.getLogger("avaframe.debrisframe.in1Utils.fileUtils")
+
+
+def copyHydrToInput(debrisDir):
+    """
+    copy the output file from in2TopoHyd into Inputs/REL
+    to serve as input data for c1TIF
+
+    Parameters
+    -----------
+    debrisDir: str or pathlib path
+        path to debris-flow directory
+    """
+
+    debrisDir = pathlib.Path(debrisDir)
+
+    outputFile = debrisDir / "Outputs" / "in2TopoHyd" / "initCondHyd.csv"
+    inputsDir = debrisDir / "Inputs" / "REL"
+    inputsDir.mkdir(parents=True, exist_ok=True)
+
+    inputsHydFile = inputsDir / outputFile.name
+
+    if not inputsHydFile.exists():
+        shutil.copy2(outputFile, inputsHydFile)
+    else:
+        outputFileDF = pd.read_csv(outputFile)
+        inputsHydFileDF = pd.read_csv(inputsHydFile)
+
+        if not outputFileDF.equals(inputsHydFileDF):
+            log.info(f"{inputsHydFile} already exists and is inconsistent with {outputFile}.")
+            log.info(f"{inputsHydFile} is overwritten by {outputFile}!")
+            shutil.copy2(outputFile, inputsHydFile)
+        else:
+            log.info(f"{inputsHydFile} already exists and is consistent with {outputFile}")

debrisframe/in1Utils/plotUtils.py

@@ -0,0 +1,112 @@
+"""
+Plotting helper functions
+"""
+
+# load modules
+import matplotlib.pyplot as plt
+
+# local imports
+import avaframe.in3Utils.fileHandlerUtils as fU
+
+
+def plotCrossSection(crossSection, outputDir):
+    """
+    This function plots the terrain cross section
+
+    Parameters
+    -----------
+    crossSection: dict
+        dictionary containing
+        the elevation of the cross-section cells and levee points,
+        the path along (distance) along the cross-section cells
+    outputDir: str or Path
+        path to output directory of in2TopoHyd module
+
+    Returns
+    --------
+    Plot: .png
+        saves plot in debris-flow directory
+    """
+
+    # plot cross section
+    fig, ax = plt.subplots(figsize=(10, 8))
+
+    ax.plot(crossSection["s"], crossSection["elevation"])
+    ax.scatter(crossSection["sLevee"], crossSection["elevLevee"], color="red", label="Levee points")
+    ax.set_xlabel("distance [m]")
+    ax.set_ylabel("elevation [m]")
+    ax.grid(color="gray", linestyle="--", linewidth=0.5, alpha=0.6)
+    ax.set_title("Cross Section")
+
+    plt.legend()
+    plt.tight_layout()
+
+    # check if directory already exists
+    path = outputDir / "Plots"
+    fU.makeADir(path)
+
+    fig.savefig(path / "crossSection.png")
+
+
+def plotRatingCurve(ratingCurve, crossSection, outputDir):
+    """
+    This function plots the terrain cross section
+    and the corresponding rating curve
+
+    Parameters
+    -----------
+    ratingCurve: dict
+        dictionary containing
+        thickness values/surface elevations and corresponding flow area
+        min and max values of horizontal distance of flow area
+    crossSection: dict
+            dictionary containing
+            the elevation of the cross-section cells,
+            the path along (distance) along the cross-section cells
+    outputDir: str or Path
+            path to output directory of in2TopoHyd module
+
+    Returns
+    --------
+    Plot: .png
+        saves plot in debris-flow directory
+    """
+
+    thickness = ratingCurve["thickness"]
+    flowArea = ratingCurve["flowArea"]
+    surfElev = ratingCurve["surfElev"]
+    xmin = ratingCurve["xmin"]
+    xmax = ratingCurve["xmax"]
+
+    # plot cross section
+    fig, ax = plt.subplots(ncols=1, nrows=2, figsize=(10, 8))
+
+    ax[0].plot(crossSection["s"], crossSection["elevation"])
+    ax[0].scatter(crossSection["sLevee"], crossSection["elevLevee"], color="red", label="Levee points")
+    for i, (elev, minimum, maximum) in enumerate(zip(surfElev, xmin, xmax)):
+        elev = [elev] * len(minimum)
+        ax[0].hlines(
+            elev,
+            xmin=minimum,
+            xmax=maximum,
+            linestyles="--",
+            colors="black",
+            lw=0.5,
+            label="elevation increments" if i == 0 else "_nolegend_",
+        )
+    ax[0].set_xlabel("distance [m]"), ax[0].set_ylabel("elevation [m]")
+    ax[0].grid(color="gray", linestyle="--", linewidth=0.5, alpha=0.6)
+    ax[0].set_title("Cross Section")
+    ax[0].legend()
+
+    ax[1].plot(flowArea, thickness)
+    ax[1].set_xlabel("flow area [m²]")
+    ax[1].set_ylabel("flow thickness [m]")
+    ax[1].grid(color="gray", linestyle="--", linewidth=0.5, alpha=0.6)
+    ax[1].set_title("Rating Curve")
+
+    plt.tight_layout()
+
+    path = outputDir / "Plots"
+    fU.makeADir(path)
+    fig.savefig(path / "ratingCurve.png")