I will attempt to answer my own question - dun dun dun.
I used SAGA GIS to examine the differences in filled watersheds using their Planchon and Darboux (PD) based filling tool ( and their Wang and Liu (WL) based filling tool for 6 different watersheds. (Here I only show case two sets of results - they were similar across all 6 watersheds) I say "based", because there is always the question as to whether differences are due to the algorithm or the specific implementation of the algorithm.
Watershed DEMs were generated by clipping mosaiced NED 30 m data using USGS provided watershed shapefiles. For each base DEM, the two tools were run; there is only one option for each tool, the minimum enforced slope, which was set in both tools to 0.01.
After the watersheds were filled, I used the raster calculator to determine the differences in the resulting grids - these differences should only be due to the different behaviors of the two algorithms.
Images representing the differences or lack of differences (basically the calculated difference raster) are presented below. The formula used in calculating differences was: (((PD_Filled - WL_Filled) / PD_Filled) * 100) - give the percent difference on a cell by cell basis. Cells grey in color show now difference, with cells redder in color indicating the resulting PD elevation was greater, and cells greener in color indicating the resulting WL elevation was greater.
1st Watershed: Clear Watershed, Wyoming
Here's the legend for these images:
Differences only range from -0.0915% to +0.0910%. Differences seem to be focused around peaks and narrow stream channels, with the WL algorithm slightly higher in the channels and PD slightly higher around the localized peaks.
Clear Watershed, Wyoming, Zoom 1
Clear Watershed, Wyoming, Zoom 2
2nd Watershed: Winnipesaukee River, NH
Here's the legend for these images:
Winnipesaukee River, NH, Zoom 1
Differences only range from -0.323% to +0.315%. Differences seem to be focused around peaks and narrow stream channels, with (as before) the WL algorithm slightly higher in the channels and PD slightly higher around the localized peaks.
Sooooooo, thoughts? To me, the differences seem trivial are probably aren't likely to affect further calculations; anyone agree? I am checking by completing my workflow for these six watersheds.
Edit: More information. It seems that the WL algorithm leads to wider less distinct channels, causing high topographic index values (my final derivative data set). The image on the left below is the PD algorithm, the image on the right is the WL algorithm.
These images show the difference in topographic index at the same locations - wider wetter areas (more channel - redder, higher TI) in the WL pic on the right; narrower channels (less wet area - less red, narrower red area, lower TI in area) in the PD pic on the left.
Additionally, here is how PD handled (left) a depression and how WL handled it (right) - notice the raised orange (lower Topographic index) segment/line crossing through the depression in the WL filled output?
So the differences, however small, do seem to trickle through the additional analyses.
Here is my Python script if anyone is interested:
#! /usr/bin/env python
# ----------------------------------------------------------------------
# Create Fill Algorithm Comparison
# Author: T. Taggart
# ----------------------------------------------------------------------
import os, sys, subprocess, time
# function definitions
def runCommand_logged (cmd, logstd, logerr):
p = subprocess.call(cmd, stdout=logstd, stderr=logerr)
# ----------------------------------------------------------------------
# ----------------------------------------------------------------------
# environmental variables/paths
if (os.name == "posix"):
os.environ["PATH"] += os.pathsep + "/usr/local/bin"
else:
os.environ["PATH"] += os.pathsep + "C:\program files (x86)\SAGA-GIS"
# ----------------------------------------------------------------------
# ----------------------------------------------------------------------
# global variables
WORKDIR = "D:\TomTaggart\DepressionFillingTest\Ran_DEMs"
# This directory is the toplevel directoru (i.e. DEM_8)
INPUTDIR = "D:\TomTaggart\DepressionFillingTest\Ran_DEMs"
STDLOG = WORKDIR + os.sep + "processing.log"
ERRLOG = WORKDIR + os.sep + "processing.error.log"
# ----------------------------------------------------------------------
# ----------------------------------------------------------------------
# open logfiles (append in case files are already existing)
logstd = open(STDLOG, "a")
logerr = open(ERRLOG, "a")
# ----------------------------------------------------------------------
# ----------------------------------------------------------------------
# initialize
t0 = time.time()
# ----------------------------------------------------------------------
# ----------------------------------------------------------------------
# loop over files, import them and calculate TWI
# this for loops walks through and identifies all the folder, sub folders, and so on.....and all the files, in the directory
# location that is passed to it - in this case the INPUTDIR
for dirname, dirnames, filenames in os.walk(INPUTDIR):
# print path to all subdirectories first.
#for subdirname in dirnames:
#print os.path.join(dirname, subdirname)
# print path to all filenames.
for filename in filenames:
#print os.path.join(dirname, filename)
filename_front, fileext = os.path.splitext(filename)
#print filename
if filename_front == "w001001":
#if fileext == ".adf":
# Resetting the working directory to the current directory
os.chdir(dirname)
# Outputting the working directory
print "\n\nCurrently in Directory: " + os.getcwd()
# Creating new Outputs directory
os.mkdir("Outputs")
# Checks
#print dirname + os.sep + filename_front
#print dirname + os.sep + "Outputs" + os.sep + ".sgrd"
# IMPORTING Files
# --------------------------------------------------------------
cmd = ['saga_cmd', '-f=q', 'io_gdal', 'GDAL: Import Raster',
'-FILES', filename,
'-GRIDS', dirname + os.sep + "Outputs" + os.sep + filename_front + ".sgrd",
#'-SELECT', '1',
'-TRANSFORM',
'-INTERPOL', '1'
]
print "Beginning to Import Files"
try:
runCommand_logged(cmd, logstd, logerr)
except Exception, e:
logerr.write("Exception thrown while processing file: " + filename + "\n")
logerr.write("ERROR: %s\n" % e)
print "Finished importing Files"
# --------------------------------------------------------------
# Resetting the working directory to the ouputs directory
os.chdir(dirname + os.sep + "Outputs")
# Depression Filling - Wang & Liu
# --------------------------------------------------------------
cmd = ['saga_cmd', '-f=q', 'ta_preprocessor', 'Fill Sinks (Wang & Liu)',
'-ELEV', filename_front + ".sgrd",
'-FILLED', filename_front + "_WL_filled.sgrd", # output - NOT optional grid
'-FDIR', filename_front + "_WL_filled_Dir.sgrd", # output - NOT optional grid
'-WSHED', filename_front + "_WL_filled_Wshed.sgrd", # output - NOT optional grid
'-MINSLOPE', '0.0100000',
]
print "Beginning Depression Filling - Wang & Liu"
try:
runCommand_logged(cmd, logstd, logerr)
except Exception, e:
logerr.write("Exception thrown while processing file: " + filename + "\n")
logerr.write("ERROR: %s\n" % e)
print "Done Depression Filling - Wang & Liu"
# Depression Filling - Planchon & Darboux
# --------------------------------------------------------------
cmd = ['saga_cmd', '-f=q', 'ta_preprocessor', 'Fill Sinks (Planchon/Darboux, 2001)',
'-DEM', filename_front + ".sgrd",
'-RESULT', filename_front + "_PD_filled.sgrd", # output - NOT optional grid
'-MINSLOPE', '0.0100000',
]
print "Beginning Depression Filling - Planchon & Darboux"
try:
runCommand_logged(cmd, logstd, logerr)
except Exception, e:
logerr.write("Exception thrown while processing file: " + filename + "\n")
logerr.write("ERROR: %s\n" % e)
print "Done Depression Filling - Planchon & Darboux"
# Raster Calculator - DIff between Planchon & Darboux and Wang & Liu
# --------------------------------------------------------------
cmd = ['saga_cmd', '-f=q', 'grid_calculus', 'Grid Calculator',
'-GRIDS', filename_front + "_PD_filled.sgrd",
'-XGRIDS', filename_front + "_WL_filled.sgrd",
'-RESULT', filename_front + "_DepFillDiff.sgrd", # output - NOT optional grid
'-FORMULA', "(((g1-h1)/g1)*100)",
'-NAME', 'Calculation',
'-FNAME',
'-TYPE', '8',
]
print "Depression Filling - Diff Calc"
try:
runCommand_logged(cmd, logstd, logerr)
except Exception, e:
logerr.write("Exception thrown while processing file: " + filename + "\n")
logerr.write("ERROR: %s\n" % e)
print "Done Depression Filling - Diff Calc"
# ----------------------------------------------------------------------
# ----------------------------------------------------------------------
# finalize
logstd.write("\n\nProcessing finished in " + str(int(time.time() - t0)) + " seconds.\n")
logstd.close
logerr.close
# ----------------------------------------------------------------------