Inconsistent results with BreachDepressionLeastCost from run to run

[jfbourdon]

Description

The current implementation of BreachDepressionLeastCost will not always give the same result even thought the same inputs are used. I found that if the DEM contains two pit cells with the same elevation, the order in which they will be process isn't always the same which can lead to different breaching solutions.

The code below reproduce the issue when used with this small DEM (dem.zip):

import os
from osgeo import gdal
import numpy as np
import subprocess

path_wbt = R"C:\Temp\whitebox_tools.exe"
outdir = "C:\Temp"
path_DEM = R"C:\Temp\dem.sdat"
max_iter = 20

def load_raster(path_raster):
    driver = gdal.GetDriverByName("SAGA")
    ds = gdal.Open(path_raster)
    band = ds.GetRasterBand(1)
    arr = band.ReadAsArray()
    return arr

# BreachDepressionsLeastCost command line
cmd_raw = " ".join([
    path_wbt,
    "--run=BreachDepressionsLeastCost",
    f"--dem={path_DEM}",
    "--min_dist",
    ])

# Run first iteration which will be used as a reference
ii = 0
outBDLC_ori = os.path.join(outdir, f"dem_BDLC_{ii}.sdat")
print(cmd_raw + f" --output={outBDLC_ori}")
subprocess.run(f"{cmd_raw} --output={outBDLC_ori}")
arr_ori = load_raster(outBDLC_ori)

# Rerun breaching until differences with the reference are found
while ii < max_iter:
    ii += 1
    outBDLC = os.path.join(outdir, f"dem_BDLC_{ii}.sdat")
    subprocess.run(f"{cmd_raw} --output={outBDLC}")
    arr_ii = load_raster(outBDLC)
    nb_diff = np.sum(arr_ori != arr_ii)
    print(f"{ii}/{max_iter} - Pixel differences: {nb_diff}")

    if nb_diff:
        print("Differences found!")
        break
    else:
        os.remove(outBDLC)
        os.remove(outBDLC[:-4] + "sgrd")

Metadata

• WhiteboxTools version: 1.3.0 and 2.4.0 (161cb56)
• Operating System version: Windows 10 Enterprise 22H2