dtarb / TauDEM

Terrain Analysis Using Digital Elevation Models (TauDEM) software for hydrologic terrain analysis and channel network extraction.
http://hydrology.usu.edu/taudem
Other
222 stars 115 forks source link

Different DEM resolutions result in different stream networks for fixed threshold #230

Open jeremymatt opened 2 years ago

jeremymatt commented 2 years ago

TauDEM version: 5.3.7 using command-line tools

For computational efficiency purposes I need to downscale the resolution of the input DEM from a cell size of 1m per side (1sq-meter cell area) to 5m per side (25sq-meters cell area). When I do so, the TauDEM tools generate a less "detailed" stream network for the 5m cell size than for the 1m cell size for the same 2sq-mile flow accumulation threshold: tributaries that show up on the 1m map do not show up on the 5m map.

I've checked the results against the watersheds returned by the USGS Streamstats tool, and the stream network generated from the 1m DEM is correct.

I've done some testing and found that if I divide the threshold by the cell size before passing, the stream network is generated correctly regardless of the input DEM cell size. That is, if the intended threshold is 20,000 sq-meters and the cell size is 5m passing an adjusted threshold of 4,000 into the Threshold tool returns the correct stream network. Similarly if the cell size is 10m, passing an adjusted threshold of 2,000 into the Threshold tool returns the expected network.

My suspicion is that the area is not being calculated correctly; that the number of cells is being multiplied by the cell size (IE 5m) rather than the cell area of 25sq-meters. The work flow is as follows (is there a command-line flag I'm missing?)

For each DEM resolution (5m in this example): mpiexec -n 8 PitRemove -z "WIN_0501_5-m.tif" -fel "dem_filled.tif" mpiexec -n 8 DinfFlowDir -fel "dem_filled.tif" -ang "flowdir.tif" -slp "slope.tif" mpiexec -n 8 AreaDinf -ang "flowdir.tif" -sca "flow_accumulation.tif" -nc For each threshold value (2589988sq-meters in this example): adjusted_threshold = threshold/cell_size mpiexec -n 8 Threshold -ssa "flow_accumulation.tif" -src "1-th-sqmi_5-cs_2589988-th-meters_not-adj.tif" -thresh 2589988.0 mpiexec -n 8 Threshold -ssa "flow_accumulation.tif" -src "1-th-sqmi_5-cs_517998-th-meters_adj.tif" -thresh 517998.0

The attached PDF has some screenshots of the the varying outputs.
2021-07-12-TauDEM area issue.pdf

dtarb commented 2 years ago

The AreaDinf function returns its results as specific catchment area. It is the number of grid cells times cell area divided by unit contour length, here estimated as cell size. Let n be area in number of cells, and d cell size. The algorithm reports nd = nd^2/d. I believe that this explains your result.

See https://hydrology.usu.edu/taudem/taudem5/help53/DInfinityContributingArea.html where this is explained in a bit more detail.

Note that his approach is only used with Dinfinity, not D8 Contributing area.

jeremymatt commented 2 years ago

Thanks for the clarification! That makes sense now. I read the helpfile you linked above, but I glossed over the distinction between "area" and "specific catchment area". Assuming it isn't too hard to implement, allowing the user to choose either raw area or specific catchment area would be a nice addition (especially since this operation is not consistent with the D8 area function). Also, adding the information you gave above (that contour length is estimated as cell size and the Aspecific = nd = nd^2/d equation) to the AreaDinf documentation might help clarify the function operation for other inexperienced users.

In the meantime, I'm now confident that my workaround of adjusting the threshold to be a threshold of "specific catchment area" rather than an "area" threshold is valid for my application.