[com4FlowPy]: Calculate iteration with generations

avaframe/com4FlowPy/com4FlowPy.py

@@ -72,6 +72,7 @@ def com4FlowPyMain(cfgPath, cfgSetup):
 
     # Flag for use of old flux distribution version
     modelParameters["fluxDistOldVersionBool"] = cfgSetup.getboolean("fluxDistOldVersion")
+    modelParameters["calcGeneration"] = cfgSetup.getboolean("calcGeneration")
 
     # Tiling Parameters used for calculation of large model-domains
     tilingParameters = {}
@@ -262,10 +263,13 @@ def startLogging(modelParameters, forestParams, modelPaths, MPOptions):
     if modelParameters["fluxDistOldVersionBool"]:
         log.info("Calculation using old (BUGGY!!) version of flux distribution!")
         log.info("------------------------")
-    for param, value in MPOptions.items():
-        log.info(f"{'%s:'%param : <20}{value : <5}")
+    try:
+        for param, value in MPOptions.items():
+            log.info(f"{'%s:'%param : <20}{value : <5}")
         # log.info("{}:\t{}".format(param,value))
-    log.info("------------------------")
+        log.info("------------------------")
+    except:
+        log.info("------------------------")
     log.info(f"{'WorkDir:' : <12}{'%s'%modelPaths['workDir'] : <5}")
     log.info(f"{'ResultsDir:' : <12}{'%s'%modelPaths['resDir'] : <5}")
     # log.info("WorkDir: {}".format(modelPaths["workDir"]))

avaframe/com4FlowPy/com4FlowPyCfg.ini

@@ -154,6 +154,12 @@ skipForestCells = 1
 # be switched on by setting "fluxDistOldVersion = True".
 
 fluxDistOldVersion = False
+#++++++++++++ Calculate with generations
+# If calcGenerations = True, a different order of cells in a path are calculated.
+# The results can vary when computing with generations.
+# Additionally, the generation (iteration step) can be derived, which is required
+# to get thalweg information.
+calcGeneration = False
 
 #++++++++++++ Parameters for Tiling
 # tileSize: size of tiles in x and y direction in meters (if total size of) x

avaframe/com4FlowPy/flowClass.py

@@ -401,10 +401,10 @@ def calc_distribution(self):
         row_local, col_local = np.where(self.dist > threshold)
 
         return (
-            self.rowindex - 1 + row_local,
-            self.colindex - 1 + col_local,
-            self.dist[row_local, col_local],
-            self.z_delta_neighbour[row_local, col_local],
+            list(self.rowindex - 1 + row_local),
+            list(self.colindex - 1 + col_local),
+            list(self.dist[row_local, col_local]),
+            list(self.z_delta_neighbour[row_local, col_local]),
         )
 
     def forest_detrainment(self):

avaframe/data/avaArzlerAlm/Inputs/dem10m.prj

@@ -0,0 +1 @@
+PROJCS["MGI_Austria_GK_West",GEOGCS["GCS_MGI",DATUM["D_MGI",SPHEROID["Bessel_1841",6377397.155,299.1528128]],PRIMEM["Greenwich",0.0],UNIT["Degree",0.0174532925199433]],PROJECTION["Transverse_Mercator"],PARAMETER["False_Easting",0.0],PARAMETER["False_Northing",-5000000.0],PARAMETER["Central_Meridian",10.3333333333333],PARAMETER["Scale_Factor",1.0],PARAMETER["Latitude_Of_Origin",0.0],UNIT["Meter",1.0]]

avaframe/runStandardTestsCom4FlowPy.py

@@ -0,0 +1,130 @@
+"""
+    Run script for running the standard tests with com4FlowPy
+    in this test all the available tests tagged standardTest are performed
+"""
+
+# Load modules
+import time
+import pathlib
+import numpy as np
+import rasterio
+
+# Local imports
+from avaframe.com4FlowPy import com4FlowPy
+from avaframe.runCom4FlowPy import readFlowPyinputs
+from avaframe.ana1Tests import testUtilities as tU
+from avaframe.in3Utils import fileHandlerUtils as fU
+from avaframe.in3Utils import initializeProject as initProj
+from avaframe.in3Utils import cfgUtils
+from avaframe.in3Utils import logUtils
+
+
+def readRasters(path):
+    raster = rasterio.open(path)
+    output = raster.read(1)
+    return output
+
+
+def compareRasters(path, pathRef, rtol):
+    raster = readRasters(path)
+    rasterRef = readRasters(pathRef)
+    diff = rasterRef - raster
+    equal = np.array_equal(rasterRef, raster)
+    closeArray = np.isclose(raster, rasterRef, rtol=rtol, equal_nan=True)
+    closePercent = np.count_nonzero(closeArray[np.logical_or(raster > 0, rasterRef > 0)]) / rasterRef[np.logical_or(raster > 0, rasterRef > 0)].size
+    return diff, equal, closePercent
+
+
+# Which result types for comparison plots
+outputVariable = ['fpTravelAngle', 'zdelta']
+
+# log file name; leave empty to use default runLog.log
+logName = 'runStandardTestsCom4FlowPy'
+
+# Load settings from general configuration file
+cfgMain = cfgUtils.getGeneralConfig()
+
+# load all benchmark info as dictionaries from description files
+testDictList = tU.readAllBenchmarkDesDicts(info=False, inDir=pathlib.Path('..', 'benchmarksCom4FlowPy'))
+
+# filter benchmarks for tag standardTest
+#  filterType = 'TAGS'
+#  valuesList = ['resistance']
+filterType = 'TAGS'
+valuesList = ['standardTest']
+testList = tU.filterBenchmarks(testDictList, filterType, valuesList, condition='or')
+
+# Set directory for full standard test report
+outDir = pathlib.Path.cwd() / 'tests' / 'reportsCom4FlowPy'
+fU.makeADir(outDir)
+
+print(outDir)
+
+# Start writing markdown style report for standard tests
+reportFile = outDir / 'standardTestsReportCom4FlowPy.md'
+with open(reportFile, 'w') as pfile:
+
+    # Write header
+    pfile.write('# Standard Tests Report \n')
+    pfile.write('## Compare com4FlowPy simulations to benchmark results \n')
+
+log = logUtils.initiateLogger(outDir, logName)
+log.info('The following benchmark tests will be fetched ')
+for test in testList:
+    log.info('%s' % test['NAME'])
+
+# run Standard Tests sequentially
+
+for test in testList:
+
+    avaDir = test['AVADIR']
+    cfgMain['MAIN']['avalancheDir'] = avaDir
+    testName = test['NAME']
+
+    # Fetch benchmark test info
+    refDir = pathlib.Path('..', 'benchmarksCom4FlowPy', test['NAME'])
+
+    # Clean input directory(ies) of old work and output files
+    initProj.cleanSingleAvaDir(avaDir)
+
+    # Load input parameters from configuration file for standard tests
+    standardCfg = refDir / ('%s_com4FlowPyCfg.ini' % test['AVANAME'])
+    modName = 'com4FlowPy'
+    cfg = cfgUtils.getModuleConfig(com4FlowPy, fileOverride=standardCfg)
+    cfgGen = cfg["GENERAL"]
+    cfgGen["cpuCount"] = str(cfgUtils.getNumberOfProcesses(cfgMain, 9999))
+
+    compDir = pathlib.Path(avaDir, 'Outputs', modName, 'peakFiles')
+    avalancheDir = cfgMain["MAIN"]["avalancheDir"]
+    cfgPath = readFlowPyinputs(avalancheDir, cfg, log)
+    cfgPath["customDirs"] = False
+    cfgPath["resDir"] = compDir
+    fU.makeADir(cfgPath["resDir"])
+    cfgPath["tempDir"] = cfgPath["workDir"] / "temp"
+    fU.makeADir(cfgPath["tempDir"])
+    cfgPath["deleteTemp"] = "False"
+    cfgPath["uid"] = cfgUtils.cfgHash(cfg)
+    cfgPath["timeString"] = ""
+    cfgPath["outputFiles"] = "zDelta|travelLength|fpTravelAngle|flux"
+
+    # Set timing
+    startTime = time.time()
+    # call com4FlowPy run
+    com4FlowPy.com4FlowPyMain(cfgPath, cfgGen)
+    endTime = time.time()
+    timeNeeded = endTime - startTime
+    log.info(('Took %s seconds to calculate.' % (timeNeeded)))
+
+    rtol = 1e-01
+    for variable in outputVariable:
+        pathRasterRef = refDir / (f'com4_{cfgPath["uid"]}__{variable}.tif')
+        pathRaster = compDir / (f'com4_{cfgPath["uid"]}__{variable}.tif')
+        diff, eq, close = compareRasters(pathRaster, pathRasterRef, rtol)
+
+        if eq and np.sum(abs(diff[diff != 0])) == 0:
+            message = f'{testName}: for {variable}: rasters are equal \n'
+        else:
+            message = f'{testName}: for {variable}: rasters are *NOT* equal, but {np.round(close,4) * 100}% of the affected area is close (relative tolerance: {rtol}) \n'
+        log.info(message)
+        with open(reportFile, 'a') as pfile:
+            pfile.write(message)