Add MFD flow length computation

README.md

@@ -294,6 +294,8 @@ In the GIS world, rasters are used for representing continuous phenomena (e.g. e
 | [Flow Direction (MFD)](xrspatial/flow_direction_mfd.py) | Partitions flow to all downslope neighbors with an adaptive exponent (Qin et al. 2007) | Qin et al. 2007 | ✅️ | ✅️ | ✅️ | ✅️ |
 | [Flow Accumulation (D8)](xrspatial/flow_accumulation.py) | Counts upstream cells draining through each cell in a D8 flow direction grid | Jenson & Domingue 1988 | ✅️ | ✅️ | ✅️ | ✅️ |
 | [Flow Accumulation (MFD)](xrspatial/flow_accumulation_mfd.py) | Accumulates upstream area through all MFD flow paths weighted by directional fractions | Qin et al. 2007 | ✅️ | ✅️ | ✅️ | 🔄 |
+| [Flow Length (D8)](xrspatial/flow_length.py) | Computes D8 flow path length from each cell to outlet (downstream) or from divide (upstream) | Standard (D8 tracing) | ✅️ | ✅️ | ✅️ | 🔄 |
+| [Flow Length (MFD)](xrspatial/flow_length_mfd.py) | Proportion-weighted flow path length using MFD fractions (downstream or upstream) | Qin et al. 2007 | ✅️ | ✅️ | ✅️ | 🔄 |
 | [Watershed](xrspatial/watershed.py) | Labels each cell with the pour point it drains to via D8 flow direction | Standard (D8 tracing) | ✅️ | ✅️ | ✅️ | ✅️ |
 | [Basins](xrspatial/watershed.py) | Delineates drainage basins by labeling each cell with its outlet ID | Standard (D8 tracing) | ✅️ | ✅️ | ✅️ | ✅️ |
 | [Stream Order](xrspatial/stream_order.py) | Assigns Strahler or Shreve stream order to cells in a drainage network | Strahler 1957, Shreve 1966 | ✅️ | ✅️ | ✅️ | ✅️ |

docs/source/reference/hydrology.rst

@@ -53,6 +53,13 @@ Flow Length
 
     xrspatial.flow_length.flow_length
 
+Flow Length (MFD)
+=================
+.. autosummary::
+    :toctree: _autosummary
+
+    xrspatial.flow_length_mfd.flow_length_mfd
+
 Flow Path
 =========
 .. autosummary::

examples/user_guide/32_Flow_Length_MFD.ipynb

@@ -0,0 +1,199 @@
+{
+ "cells": [
+  {
+   "cell_type": "markdown",
+   "id": "a1b2c3d4",
+   "metadata": {},
+   "source": [
+    "# Flow Length (MFD)\n",
+    "\n",
+    "Flow length measures the distance water travels along its flow path. The MFD (Multiple Flow Direction) variant computes proportion-weighted path lengths using MFD fraction grids, where each cell distributes flow to all downslope neighbors.\n",
+    "\n",
+    "Two modes are supported:\n",
+    "- **Downstream**: Expected (weighted-average) distance from each cell to its outlet\n",
+    "- **Upstream**: Longest flow path from any drainage divide to each cell"
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": null,
+   "id": "e5f6a7b8",
+   "metadata": {},
+   "outputs": [],
+   "source": [
+    "import numpy as np\n",
+    "import matplotlib.pyplot as plt\n",
+    "import xarray as xr\n",
+    "\n",
+    "from xrspatial import flow_direction_mfd, flow_length_mfd, flow_accumulation_mfd\n",
+    "from xrspatial.terrain import generate_terrain"
+   ]
+  },
+  {
+   "cell_type": "markdown",
+   "id": "c9d0e1f2",
+   "metadata": {},
+   "source": [
+    "## Generate terrain"
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": null,
+   "id": "a3b4c5d6",
+   "metadata": {},
+   "outputs": [],
+   "source": [
+    "W = 400\n",
+    "H = 400\n",
+    "terrain = generate_terrain(\n",
+    "    x_range=(-2, 2), y_range=(-2, 2), width=W, height=H, seed=42\n",
+    ")\n",
+    "\n",
+    "fig, ax = plt.subplots(figsize=(8, 6))\n",
+    "im = ax.imshow(terrain.values, cmap='terrain', origin='lower')\n",
+    "ax.set_title('Elevation')\n",
+    "plt.colorbar(im, ax=ax, label='meters')\n",
+    "plt.tight_layout()\n",
+    "plt.show()"
+   ]
+  },
+  {
+   "cell_type": "markdown",
+   "id": "e7f8a9b0",
+   "metadata": {},
+   "source": [
+    "## Compute MFD flow direction and accumulation"
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": null,
+   "id": "c1d2e3f4",
+   "metadata": {},
+   "outputs": [],
+   "source": [
+    "mfd_dir = flow_direction_mfd(terrain)\n",
+    "mfd_acc = flow_accumulation_mfd(mfd_dir)\n",
+    "\n",
+    "fig, ax = plt.subplots(figsize=(8, 6))\n",
+    "im = ax.imshow(np.log1p(mfd_acc.values), cmap='Blues', origin='lower')\n",
+    "ax.set_title('MFD flow accumulation (log scale)')\n",
+    "plt.colorbar(im, ax=ax, label='log(1 + accumulation)')\n",
+    "plt.tight_layout()\n",
+    "plt.show()"
+   ]
+  },
+  {
+   "cell_type": "markdown",
+   "id": "a5b6c7d8",
+   "metadata": {},
+   "source": [
+    "## Downstream flow length\n",
+    "\n",
+    "Downstream flow length gives the proportion-weighted average distance from each cell to the outlet it drains to. Cells near outlets have short distances; cells far upstream have long distances."
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": null,
+   "id": "e9f0a1b2",
+   "metadata": {},
+   "outputs": [],
+   "source": [
+    "downstream = flow_length_mfd(mfd_dir, direction='downstream')\n",
+    "\n",
+    "fig, ax = plt.subplots(figsize=(8, 6))\n",
+    "im = ax.imshow(downstream.values, cmap='viridis', origin='lower')\n",
+    "ax.set_title('MFD downstream flow length')\n",
+    "plt.colorbar(im, ax=ax, label='distance (coordinate units)')\n",
+    "plt.tight_layout()\n",
+    "plt.show()"
+   ]
+  },
+  {
+   "cell_type": "markdown",
+   "id": "c3d4e5f6",
+   "metadata": {},
+   "source": [
+    "## Upstream flow length\n",
+    "\n",
+    "Upstream flow length gives the longest flow path distance from any drainage divide to each cell. Cells on divides have zero length; cells at outlets accumulate the full watershed length."
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": null,
+   "id": "a7b8c9d0",
+   "metadata": {},
+   "outputs": [],
+   "source": [
+    "upstream = flow_length_mfd(mfd_dir, direction='upstream')\n",
+    "\n",
+    "fig, ax = plt.subplots(figsize=(8, 6))\n",
+    "im = ax.imshow(upstream.values, cmap='magma', origin='lower')\n",
+    "ax.set_title('MFD upstream flow length')\n",
+    "plt.colorbar(im, ax=ax, label='distance (coordinate units)')\n",
+    "plt.tight_layout()\n",
+    "plt.show()"
+   ]
+  },
+  {
+   "cell_type": "markdown",
+   "id": "e1f2a3b4",
+   "metadata": {},
+   "source": [
+    "## Comparison: D8 vs MFD flow length\n",
+    "\n",
+    "MFD distributes flow across multiple neighbors, producing smoother flow length fields than D8 which routes everything through the single steepest neighbor."
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": null,
+   "id": "c5d6e7f8",
+   "metadata": {},
+   "outputs": [],
+   "source": [
+    "from xrspatial import flow_direction, flow_length\n",
+    "\n",
+    "d8_dir = flow_direction(terrain)\n",
+    "d8_downstream = flow_length(d8_dir, direction='downstream')\n",
+    "\n",
+    "fig, axes = plt.subplots(1, 2, figsize=(14, 5))\n",
+    "\n",
+    "im0 = axes[0].imshow(d8_downstream.values, cmap='viridis', origin='lower')\n",
+    "axes[0].set_title('D8 downstream flow length')\n",
+    "plt.colorbar(im0, ax=axes[0], label='distance')\n",
+    "\n",
+    "im1 = axes[1].imshow(downstream.values, cmap='viridis', origin='lower')\n",
+    "axes[1].set_title('MFD downstream flow length')\n",
+    "plt.colorbar(im1, ax=axes[1], label='distance')\n",
+    "\n",
+    "plt.tight_layout()\n",
+    "plt.show()"
+   ]
+  }
+ ],
+ "metadata": {
+  "kernelspec": {
+   "display_name": "Python 3 (ipykernel)",
+   "language": "python",
+   "name": "python3"
+  },
+  "language_info": {
+   "codemirror_mode": {
+    "name": "ipython",
+    "version": 3
+   },
+   "file_extension": ".py",
+   "mimetype": "text/x-python",
+   "name": "python",
+   "nbconvert_exporter": "python",
+   "pygments_lexer": "ipython3",
+   "version": "3.14.2"
+  }
+ },
+ "nbformat": 4,
+ "nbformat_minor": 5
+}

xrspatial/__init__.py

@@ -45,6 +45,7 @@
 from xrspatial.flow_direction_dinf import flow_direction_dinf  # noqa
 from xrspatial.flow_direction_mfd import flow_direction_mfd  # noqa
 from xrspatial.flow_length import flow_length  # noqa
+from xrspatial.flow_length_mfd import flow_length_mfd  # noqa
 from xrspatial.flow_path import flow_path  # noqa
 from xrspatial.focal import mean  # noqa
 from xrspatial.glcm import glcm_texture  # noqa