|
4 | 4 | "cell_type": "markdown", |
5 | 5 | "id": "title", |
6 | 6 | "metadata": {}, |
7 | | - "source": "# Xarray-spatial\n### User Guide: Hydrology tools\n-----\nThe Hydrology tools provide a complete workflow for analyzing drainage patterns and water flow across a terrain surface. Starting from a digital elevation model (DEM), you can compute flow directions, accumulate flow, delineate watersheds, and extract stream networks.\n\nThis guide walks through each tool in the order you would typically use them:\n\n[Sink Detection](#Sink-Detection): Identifies and labels depression cells in a flow direction grid.\n\n[Depression Filling](#Depression-Filling): Removes sinks from a DEM so water can flow continuously to the edges.\n\n[Flow Direction](#Flow-Direction): Computes D8 flow direction for each cell based on steepest descent.\n\n[Flow Accumulation](#Flow-Accumulation): Counts how many upstream cells drain through each cell.\n\n[Drainage Basins](#Drainage-Basins): Labels every cell with the ID of the outlet it drains to.\n\n[Stream Order](#Stream-Order): Assigns Strahler or Shreve stream order to cells in the drainage network.\n\n[Stream Link](#Stream-Link): Segments the stream network into individually labeled links.\n\n[Snap Pour Point](#Snap-Pour-Point): Moves user-placed pour points onto the nearest high-accumulation cell.\n\n[Watershed Delineation](#Watershed-Delineation): Labels every cell with the pour point it drains to.\n\n[Flow Path Tracing](#Flow-Path-Tracing): Traces downstream paths from selected start points to their outlets.\n\n-----------" |
| 7 | + "source": "# Xarray-spatial\n### User Guide: Hydrology tools\n-----\nThe Hydrology tools provide a complete workflow for analyzing drainage patterns and water flow across a terrain surface. Starting from a digital elevation model (DEM), you can compute flow directions, accumulate flow, delineate watersheds, and extract stream networks.\n\nThis guide walks through each tool in the order you would typically use them:\n\n[Sink Detection](#Sink-Detection): Identifies and labels depression cells in a flow direction grid.\n\n[Depression Filling](#Depression-Filling): Removes sinks from a DEM so water can flow continuously to the edges.\n\n[Flow Direction](#Flow-Direction): Computes D8 flow direction for each cell based on steepest descent.\n\n[Multiple Flow Direction (MFD)](#Multiple-Flow-Direction-(MFD)): Partitions flow to all downslope neighbors with an adaptive exponent.\n\n[Flow Accumulation](#Flow-Accumulation): Counts how many upstream cells drain through each cell.\n\n[Drainage Basins](#Drainage-Basins): Labels every cell with the ID of the outlet it drains to.\n\n[Stream Order](#Stream-Order): Assigns Strahler or Shreve stream order to cells in the drainage network.\n\n[Stream Link](#Stream-Link): Segments the stream network into individually labeled links.\n\n[Snap Pour Point](#Snap-Pour-Point): Moves user-placed pour points onto the nearest high-accumulation cell.\n\n[Watershed Delineation](#Watershed-Delineation): Labels every cell with the pour point it drains to.\n\n[Flow Path Tracing](#Flow-Path-Tracing): Traces downstream paths from selected start points to their outlets.\n\n-----------" |
8 | 8 | }, |
9 | 9 | { |
10 | 10 | "cell_type": "code", |
|
296 | 296 | } |
297 | 297 | ] |
298 | 298 | }, |
| 299 | + { |
| 300 | + "cell_type": "markdown", |
| 301 | + "id": "nbozw06cw1", |
| 302 | + "source": "## Multiple Flow Direction (MFD)\n\nD8 sends all flow to a single neighbor, which works well for channelized flow but not for dispersive hillslope flow. The `flow_direction_mfd` function partitions flow from each cell to **all** downslope neighbors. An adaptive exponent from Qin et al. (2007) adjusts per-cell: steep convergent terrain concentrates flow (closer to D8), while gentle slopes spread it out.\n\nThe output is a 3-D DataArray `(8, H, W)` where each band holds the fraction of flow directed to one of the 8 neighbors (E, SE, S, SW, W, NW, N, NE). Fractions sum to 1.0 at each cell.\n\nParameters:\n- **p**: flow-partition exponent. `None` (default) uses the adaptive exponent. A positive float sets a fixed exponent (e.g. `p=1.0` for Quinn et al. 1991).\n- **boundary**: same edge-handling options as `flow_direction`.", |
| 303 | + "metadata": {} |
| 304 | + }, |
| 305 | + { |
| 306 | + "cell_type": "code", |
| 307 | + "id": "48fq6hhjuwj", |
| 308 | + "source": "mfd_fracs = xrspatial.flow_direction_mfd(dem_filled)\nprint(f\"MFD output shape: {mfd_fracs.shape}\")\nprint(f\"Dims: {mfd_fracs.dims}\")\nprint(f\"Neighbors: {list(mfd_fracs.coords['neighbor'].values)}\")\n\n# Show the maximum fraction per cell as a measure of flow concentration.\n# Values near 1.0 = almost all flow to one neighbor (D8-like).\n# Values near 0.125 = flow spread nearly evenly across 8 neighbors.\nimport matplotlib.pyplot as plt\n\nmax_frac = mfd_fracs.max(dim='neighbor')\nn_receivers = (mfd_fracs > 0).sum(dim='neighbor').astype(float)\nn_receivers = n_receivers.where(~np.isnan(mfd_fracs.isel(neighbor=0)))\n\nfig, axes = plt.subplots(1, 2, figsize=(14, 5))\n\nmax_frac.plot(ax=axes[0], cmap='YlOrRd', vmin=0, vmax=1)\naxes[0].set_title('Max flow fraction (higher = more concentrated)')\n\nn_receivers.plot(ax=axes[1], cmap='viridis', vmin=0, vmax=8)\naxes[1].set_title('Number of receiving neighbors')\n\nplt.tight_layout()\nplt.show()", |
| 309 | + "metadata": {}, |
| 310 | + "execution_count": null, |
| 311 | + "outputs": [] |
| 312 | + }, |
| 313 | + { |
| 314 | + "cell_type": "markdown", |
| 315 | + "id": "8ua76wuod0q", |
| 316 | + "source": "### Comparing D8, D-infinity, and MFD\n\nThe three flow direction algorithms represent different trade-offs between simplicity and physical realism:\n\n- **D8**: all flow goes to the single steepest neighbor. Produces clean, deterministic channels but can't represent dispersive hillslope flow.\n- **D-infinity** (Tarboton 1997): flow direction is a continuous angle toward the steepest downslope facet. Splits flow between at most two neighbors. Better for smooth surfaces but still limited to two receivers.\n- **MFD** (Qin et al. 2007): partitions flow to all downslope neighbors. Most realistic for hillslope processes but produces 8 output bands instead of one.\n\nBelow we compare all three on the same DEM.", |
| 317 | + "metadata": {}, |
| 318 | + "execution_count": null, |
| 319 | + "outputs": [] |
| 320 | + }, |
| 321 | + { |
| 322 | + "cell_type": "code", |
| 323 | + "id": "tz72zbmimfk", |
| 324 | + "source": "# D-infinity flow direction (continuous angle in radians)\nflow_dir_dinf = xrspatial.flow_direction_dinf(dem_filled)\n\nfig, axes = plt.subplots(1, 3, figsize=(18, 5))\n\n# D8: color-code the 9 direction codes\nshade_d8 = flow_dir.copy()\nshade_d8.plot(ax=axes[0], cmap='hsv', add_colorbar=True)\naxes[0].set_title('D8 direction code\\n(one of 9 discrete values)')\n\n# Dinf: continuous angle 0-2pi\nflow_dir_dinf.plot(ax=axes[1], cmap='hsv', vmin=0, vmax=2*np.pi, add_colorbar=True)\naxes[1].set_title('D-infinity angle (radians)\\n(continuous 0 to 2pi)')\n\n# MFD: show max fraction as a proxy for how concentrated the flow is\nmax_frac.plot(ax=axes[2], cmap='YlOrRd', vmin=0, vmax=1, add_colorbar=True)\naxes[2].set_title('MFD max fraction\\n(1.0 = all flow to one neighbor)')\n\nfor ax in axes:\n ax.set_xlabel('')\n ax.set_ylabel('')\n\nplt.suptitle('Flow direction comparison: D8 vs D-infinity vs MFD')\nplt.tight_layout()\nplt.show()\n\n# Key differences in a small summary\nd8_codes = flow_dir.values[~np.isnan(flow_dir.values)]\nmfd_max = max_frac.values[~np.isnan(max_frac.values)]\nmfd_max_nonzero = mfd_max[mfd_max > 0]\n\nprint(f\"D8: {len(np.unique(d8_codes))} unique direction codes (always 1 receiver)\")\nprint(f\"Dinf: continuous angle, splits flow between at most 2 neighbors\")\nprint(f\"MFD: median max fraction = {np.median(mfd_max_nonzero):.3f}, \"\n f\"mean receivers = {(n_receivers.values[~np.isnan(n_receivers.values)]).mean():.1f}\")", |
| 325 | + "metadata": {}, |
| 326 | + "execution_count": null, |
| 327 | + "outputs": [] |
| 328 | + }, |
| 329 | + { |
| 330 | + "cell_type": "code", |
| 331 | + "id": "n9gfafn0u4b", |
| 332 | + "source": "# Compare MFD exponent modes\nmfd_p1 = xrspatial.flow_direction_mfd(dem_filled, p=1.0)\nmfd_p8 = xrspatial.flow_direction_mfd(dem_filled, p=8.0)\n\nfig, axes = plt.subplots(1, 3, figsize=(16, 5))\nfor ax, data, title in zip(\n axes,\n [mfd_p1.max(dim='neighbor'), max_frac, mfd_p8.max(dim='neighbor')],\n ['p=1.0 (dispersive)', 'Adaptive (Qin 2007)', 'p=8.0 (concentrated)'],\n):\n data.plot(ax=ax, cmap='YlOrRd', vmin=0, vmax=1)\n ax.set_title(title)\n ax.set_xlabel('')\n ax.set_ylabel('')\n\nplt.suptitle('MFD exponent comparison: lower p spreads flow, higher p concentrates it')\nplt.tight_layout()\nplt.show()", |
| 333 | + "metadata": {}, |
| 334 | + "execution_count": null, |
| 335 | + "outputs": [] |
| 336 | + }, |
299 | 337 | { |
300 | 338 | "cell_type": "markdown", |
301 | 339 | "id": "flowacc-md", |
|
423 | 461 | ")" |
424 | 462 | ] |
425 | 463 | }, |
| 464 | + { |
| 465 | + "cell_type": "markdown", |
| 466 | + "id": "httgn03vcco", |
| 467 | + "source": "## Flow Accumulation (MFD)\n\n`flow_accumulation_mfd` takes the 3-D fractional output from\n`flow_direction_mfd` and routes upstream contributing area through all\ndownslope paths at once. Where D8 accumulation produces sharp, single-pixel\ndrainage lines, MFD accumulation spreads flow across the landscape and\nproduces smoother contributing-area fields.", |
| 468 | + "metadata": {} |
| 469 | + }, |
| 470 | + { |
| 471 | + "cell_type": "code", |
| 472 | + "id": "btr592uiuve", |
| 473 | + "source": "flow_accum_mfd = xrspatial.flow_accumulation_mfd(mfd)\n\nfig, axes = plt.subplots(1, 2, figsize=(14, 5))\n\n# D8 accumulation (log scale)\nim0 = axes[0].imshow(np.log1p(flow_accum.values), cmap='Blues')\naxes[0].set_title('D8 flow accumulation (log)')\nfig.colorbar(im0, ax=axes[0], shrink=0.6)\n\n# MFD accumulation (log scale)\nim1 = axes[1].imshow(np.log1p(flow_accum_mfd.values), cmap='Blues')\naxes[1].set_title('MFD flow accumulation (log)')\nfig.colorbar(im1, ax=axes[1], shrink=0.6)\n\nplt.tight_layout()\nplt.show()", |
| 474 | + "metadata": {}, |
| 475 | + "execution_count": null, |
| 476 | + "outputs": [] |
| 477 | + }, |
426 | 478 | { |
427 | 479 | "cell_type": "markdown", |
428 | 480 | "id": "basin-md", |
|
911 | 963 | "cell_type": "markdown", |
912 | 964 | "id": "references", |
913 | 965 | "metadata": {}, |
914 | | - "source": [ |
915 | | - "## References\n", |
916 | | - "\n", |
917 | | - "- Jenson, S.K. and Domingue, J.O. (1988). Extracting Topographic Structure from Digital Elevation Data for Geographic Information System Analysis. *Photogrammetric Engineering and Remote Sensing*, 54(11), 1593-1600.\n", |
918 | | - "- Planchon, O. and Darboux, F. (2001). A fast, simple and versatile algorithm to fill the depressions of digital elevation models. *Catena*, 46(2-3), 159-176.\n", |
919 | | - "- Strahler, A.N. (1957). Quantitative analysis of watershed geomorphology. *Transactions of the American Geophysical Union*, 38(6), 913-920.\n", |
920 | | - "- Copernicus DEM on AWS: https://registry.opendata.aws/copernicus-dem/\n", |
921 | | - "- ESRI Hydrology Toolset: https://pro.arcgis.com/en/pro-app/tool-reference/spatial-analyst/an-overview-of-the-hydrology-tools.htm" |
922 | | - ] |
| 966 | + "source": "## References\n\n- Jenson, S.K. and Domingue, J.O. (1988). Extracting Topographic Structure from Digital Elevation Data for Geographic Information System Analysis. *Photogrammetric Engineering and Remote Sensing*, 54(11), 1593-1600.\n- Planchon, O. and Darboux, F. (2001). A fast, simple and versatile algorithm to fill the depressions of digital elevation models. *Catena*, 46(2-3), 159-176.\n- Quinn, P., Beven, K., Chevallier, P., and Planchon, O. (1991). The prediction of hillslope flow paths for distributed hydrological modelling using digital terrain models. *Hydrological Processes*, 5(1), 59-79.\n- Qin, C., Zhu, A.X., Pei, T., Li, B., Zhou, C., and Yang, L. (2007). An adaptive approach to selecting a flow-partition exponent for a multiple-flow-direction algorithm. *International Journal of Geographical Information Science*, 21(4), 443-458.\n- Strahler, A.N. (1957). Quantitative analysis of watershed geomorphology. *Transactions of the American Geophysical Union*, 38(6), 913-920.\n- Tarboton, D.G. (1997). A new method for the determination of flow directions and upslope areas in grid digital elevation models. *Water Resources Research*, 33(2), 309-319.\n- Copernicus DEM on AWS: https://registry.opendata.aws/copernicus-dem/\n- ESRI Hydrology Toolset: https://pro.arcgis.com/en/pro-app/tool-reference/spatial-analyst/an-overview-of-the-hydrology-tools.htm" |
923 | 967 | } |
924 | 968 | ], |
925 | 969 | "metadata": { |
|
0 commit comments