Extracting Stream Networks

In [16]:

Copied!





import rioxarray as rxr
import geopandas as gpd
import matplotlib.pyplot as plt

from streamkit.terrain import flow_accumulation_workflow
from streamkit.extraction import heads_from_features
from streamkit.extraction import extract_channel_network
from streamkit.extraction import delineate_reaches
from streamkit.extraction import delineate_subbasins
from streamkit.network import vectorize_streams
import rioxarray as rxr
import geopandas as gpd
import matplotlib.pyplot as plt

from streamkit.terrain import flow_accumulation_workflow
from streamkit.extraction import heads_from_features
from streamkit.extraction import extract_channel_network
from streamkit.extraction import delineate_reaches
from streamkit.extraction import delineate_subbasins
from streamkit.network import vectorize_streams

In [7]:

Copied!

dem, flowlines = load_sample_data()
dem, flowlines = load_sample_data()

Loading sample data from cache at /Users/arthurkoehl/Library/Caches/streamkit...

Extracting the channel network from USGS 10m dem using NHD Flowlines as reference¶

In [14]:

Copied!

conditioned_dem, flow_dir, flow_acc = flow_accumulation_workflow(dem)
conditioned_dem, flow_dir, flow_acc = flow_accumulation_workflow(dem)

In [3]:

Copied!

channel_heads = heads_from_features(flowlines, dem)
channel_network = extract_channel_network(channel_heads, flow_dir)
channel_heads = heads_from_features(flowlines, dem)
channel_network = extract_channel_network(channel_heads, flow_dir)

In [4]:

Copied!

fig, ax = plt.subplots()
channel_network.plot(ax=ax)
fig, ax = plt.subplots()
channel_network.plot(ax=ax)

Out[4]:

<matplotlib.collections.QuadMesh at 0x139cb4200>

No description has been provided for this image

Computing Reaches and subbasins¶

In [15]:

Copied!

reaches = delineate_reaches(channel_network, dem, flow_dir, flow_acc, penalty=5, min_length=500, smooth_window=5, threshold_degrees=1)
reaches = delineate_reaches(channel_network, dem, flow_dir, flow_acc, penalty=5, min_length=500, smooth_window=5, threshold_degrees=1)

In [18]:

Copied!

# vectorize for the plot
vreaches = vectorize_streams(reaches, flow_dir, flow_acc)
vstreams = vectorize_streams(channel_network, flow_dir, flow_acc)
# vectorize for the plot
vreaches = vectorize_streams(reaches, flow_dir, flow_acc)
vstreams = vectorize_streams(channel_network, flow_dir, flow_acc)

In [21]:

Copied!





fig, axes = plt.subplots(1,2, figsize=(12,6))
vstreams['color_idx'] = range(len(vstreams))
vreaches['color_idx'] = range(len(vreaches))

dem.plot(ax=axes[0], add_colorbar=False, cmap='terrain')
vstreams.plot(ax=axes[0], column='color_idx', cmap='tab20', legend=False)
dem.plot(ax=axes[1], add_colorbar=False, cmap='terrain')
vreaches.plot(ax=axes[1], column='color_idx', cmap='tab20', legend=False)
fig, axes = plt.subplots(1,2, figsize=(12,6))
vstreams['color_idx'] = range(len(vstreams))
vreaches['color_idx'] = range(len(vreaches))

dem.plot(ax=axes[0], add_colorbar=False, cmap='terrain')
vstreams.plot(ax=axes[0], column='color_idx', cmap='tab20', legend=False)
dem.plot(ax=axes[1], add_colorbar=False, cmap='terrain')
vreaches.plot(ax=axes[1], column='color_idx', cmap='tab20', legend=False)

Out[21]:

<Axes: title={'center': 'band = 1, spatial_ref = 0'}, xlabel='x', ylabel='y'>

In [22]:

Copied!

subbasins = delineate_subbasins(reaches, flow_dir, flow_acc)
subbasins = delineate_subbasins(reaches, flow_dir, flow_acc)

In [26]:

Copied!





fig, ax = plt.subplots()
dem.plot(ax=ax, cmap='terrain',  add_colorbar=False,)
subbasins.plot(ax=ax, alpha=0.8, cmap='tab20', add_colorbar=False)
vreaches.plot(ax=ax, column='color_idx', cmap='tab20', legend=False)
fig, ax = plt.subplots()
dem.plot(ax=ax, cmap='terrain',  add_colorbar=False,)
subbasins.plot(ax=ax, alpha=0.8, cmap='tab20', add_colorbar=False)
vreaches.plot(ax=ax, column='color_idx', cmap='tab20', legend=False)

Out[26]:

<Axes: title={'center': 'spatial_ref = 0'}, xlabel='x coordinate of projection\n[metre]', ylabel='y coordinate of projection\n[metre]'>