rivabar.visualization

rivabar.visualization

plot_im_and_lines(im, left_utm_x, right_utm_x, lower_utm_y, upper_utm_y, G_rook, G_primal, D_primal=None, plot_main_banklines=True, plot_lines=True, plot_image=True, smoothing=False, start_x=None, start_y=None, end_x=None, end_y=None, bankline_color='tab:blue', bankline_alpha=1.0, cmap='Grays', alpha=0.5, ax=None)

Plots an image with overlaid lines representing bank polygons and centerlines from graphs.

Parameters:

Name	Type	Description	Default
im	ndarray	The image to be plotted.	required
left_utm_x	float	The left boundary of the image in UTM coordinates.	required
right_utm_x	float	The right boundary of the image in UTM coordinates.	required
lower_utm_y	float	The lower boundary of the image in UTM coordinates.	required
upper_utm_y	float	The upper boundary of the image in UTM coordinates.	required
G_rook	Graph	A graph where nodes contain 'bank_polygon' attributes representing bank polygons.	required
G_primal	Graph	A graph where edges contain 'geometry' attributes representing edge geometries.	required
D_primal	DiGraph	A directed graph containing main channel bank coordinates in graph attributes.	None
plot_main_banklines	bool	If True, plot the main channel banklines (default is True).	True
plot_lines	bool	If True, plot the centerlines from G_primal (default is True).	True
plot_image	bool	If True, plot the image (default is True).	True
smoothing	bool	If True, apply smoothing to the lines (default is False).	False
start_x	float	The x-coordinate of the starting point for smoothing (default is None).	None
start_y	float	The y-coordinate of the starting point for smoothing (default is None).	None
end_x	float	The x-coordinate of the ending point for smoothing (default is None).	None
end_y	float	The y-coordinate of the ending point for smoothing (default is None).	None
bankline_color	str	The color of the banklines (default is 'tab:blue').	'tab:blue'
bankline_alpha	float	The transparency of the banklines (default is 1.0).	1.0
cmap	str	The colormap to use for the image (default is 'Grays').	'Grays'
alpha	float	The transparency of the image (default is 0.5).	0.5
ax	Axes	The axes on which to plot. If None, creates new figure and axes.	None

Returns:

Name	Type	Description
fig	Figure	The resulting figure.
ax	Axes	The axes object.

Source code in rivabar/visualization.py

def plot_im_and_lines(im, left_utm_x, right_utm_x, lower_utm_y, upper_utm_y, G_rook, 
                G_primal, D_primal=None, plot_main_banklines=True, plot_lines=True, plot_image=True, 
                smoothing=False, start_x=None, start_y=None, end_x=None, end_y=None,
                bankline_color='tab:blue', bankline_alpha=1.0,
                cmap='Grays', alpha=0.5, ax=None):
    """
    Plots an image with overlaid lines representing bank polygons and centerlines from graphs.

    Parameters
    ----------
    im : ndarray
        The image to be plotted.
    left_utm_x : float
        The left boundary of the image in UTM coordinates.
    right_utm_x : float
        The right boundary of the image in UTM coordinates.
    lower_utm_y : float
        The lower boundary of the image in UTM coordinates.
    upper_utm_y : float
        The upper boundary of the image in UTM coordinates.
    G_rook : networkx.Graph
        A graph where nodes contain 'bank_polygon' attributes representing bank polygons.
    G_primal : networkx.Graph
        A graph where edges contain 'geometry' attributes representing edge geometries.
    D_primal : networkx.DiGraph, optional
        A directed graph containing main channel bank coordinates in graph attributes.
    plot_main_banklines : bool, optional
        If True, plot the main channel banklines (default is True).
    plot_lines : bool, optional
        If True, plot the centerlines from G_primal (default is True).
    plot_image : bool, optional
        If True, plot the image (default is True).
    smoothing : bool, optional
        If True, apply smoothing to the lines (default is False).
    start_x : float, optional
        The x-coordinate of the starting point for smoothing (default is None).
    start_y : float, optional
        The y-coordinate of the starting point for smoothing (default is None).
    end_x : float, optional
        The x-coordinate of the ending point for smoothing (default is None).
    end_y : float, optional
        The y-coordinate of the ending point for smoothing (default is None).
    bankline_color : str, optional
        The color of the banklines (default is 'tab:blue').
    bankline_alpha : float, optional
        The transparency of the banklines (default is 1.0).
    cmap : str, optional
        The colormap to use for the image (default is 'Grays').
    alpha : float, optional
        The transparency of the image (default is 0.5).
    ax : matplotlib.axes.Axes, optional
        The axes on which to plot. If None, creates new figure and axes.

    Returns
    -------
    fig : matplotlib.figure.Figure
        The resulting figure.
    ax : matplotlib.axes.Axes
        The axes object.
    """
    if ax is None:
        fig, ax = plt.subplots()
    else:
        fig = ax.get_figure()
    if plot_image:
        if im is None:
            print('No image to plot!')
        else:
            plt.imshow(im, extent = [left_utm_x, right_utm_x, lower_utm_y, upper_utm_y], cmap=cmap, alpha=alpha)
    else:
        plt.axis('equal')

    # Plot main channel banklines from nodes 0 and 1
    for i in range(2):
        if type(G_rook.nodes()[i]['bank_polygon']) == Polygon:
            x = np.array(G_rook.nodes()[i]['bank_polygon'].exterior.xy[0])
            y = np.array(G_rook.nodes()[i]['bank_polygon'].exterior.xy[1])
            if smoothing and start_x is not None and end_x is not None:
                ind1 = find_closest_point(start_x, start_y, np.vstack((x, y)).T)
                ind2 = find_closest_point(end_x, end_y, np.vstack((x, y)).T)
                if ind1 < ind2:
                    x = x[ind1:ind2]
                    y = y[ind1:ind2]
                else:
                    x = x[ind2:ind1]
                    y = y[ind2:ind1]                
        else:
            x = np.array(G_rook.nodes()[i]['bank_polygon'].xy[0])
            y = np.array(G_rook.nodes()[i]['bank_polygon'].xy[1])
        if smoothing and len(x) > 31:
            x, y = smooth_line(x, y, spline_ds=100, spline_smoothing=10000, 
                              savgol_window=min(31, len(x)), savgol_poly_order=3)
        if plot_main_banklines:
            plt.plot(x, y, color=bankline_color, alpha=bankline_alpha)

    # Plot interior banklines (islands)
    for i in trange(2, len(G_rook.nodes)):
        x_interior = np.array(G_rook.nodes()[i]['bank_polygon'].exterior.xy[0])
        y_interior = np.array(G_rook.nodes()[i]['bank_polygon'].exterior.xy[1])
        if smoothing and len(x_interior) > 21:
            x_interior_sm = savgol_filter(x_interior, 21, 3)
            y_interior_sm = savgol_filter(y_interior, 21, 3)
            interior_lstr = LineString(zip(x_interior_sm, y_interior_sm))
            x_int, y_int = remove_endpoints(interior_lstr, remove_count=1)
            interior_lstr = LineString(zip(x_int, y_int)).simplify(3)
            x_interior = np.array(interior_lstr.xy[0])
            y_interior = np.array(interior_lstr.xy[1])
        plt.plot(x_interior, y_interior, '-', color=bankline_color, alpha=bankline_alpha)

    if plot_lines:
        for s,e,d in tqdm(G_primal.edges):
            x = G_primal[s][e][d]['geometry'].xy[0]
            y = G_primal[s][e][d]['geometry'].xy[1]
            plt.plot(x, y, 'k')
    return fig, ax

plot_colored_line(x, y, cmap='magma', linewidth=2, ax=None)

Plot a line with (x, y) coordinates colored according to a colormap.

Parameters:

Name	Type	Description	Default
x	(array - like, shape(n))	The x coordinates of the points.	required
y	(array - like, shape(n))	The y coordinates of the points.	required
cmap	str or Colormap	The colormap to use for coloring the line.	'viridis'
linewidth	float	The width of the line.	2
ax	Axes	The axes on which to plot the line. If None, a new figure and axes are created.	None

Returns:

Name	Type	Description
line	LineCollection	The LineCollection object representing the colored line.

Source code in rivabar/visualization.py

def plot_colored_line(x, y, cmap='magma', linewidth=2, ax=None):
    """
    Plot a line with (x, y) coordinates colored according to a colormap.

    Parameters
    ----------
    x : array-like, shape (n,)
        The x coordinates of the points.
    y : array-like, shape (n,)
        The y coordinates of the points.
    cmap : str or Colormap, optional, default: 'viridis'
        The colormap to use for coloring the line.
    linewidth : float, optional, default: 2
        The width of the line.
    ax : matplotlib.axes.Axes, optional
        The axes on which to plot the line. If None, a new figure and axes are created.

    Returns
    -------
    line : matplotlib.collections.LineCollection
        The LineCollection object representing the colored line.
    """

    # Create segments of the line
    points = np.array([x, y]).T.reshape(-1, 1, 2)
    segments = np.concatenate([points[:-1], points[1:]], axis=1)
    # Normalize the colormap
    norm = Normalize(0, len(segments))
    # Create the LineCollection
    lc = LineCollection(segments, cmap=plt.get_cmap(cmap), norm=norm, linewidths=linewidth, joinstyle='round', capstyle='round')
    lc.set_array(np.arange(len(segments)))  # Use the x-values for the color mapping
    # Plot the line using the LineCollection
    if not ax:
        fig, ax = plt.subplots()
    line = ax.add_collection(lc)
    return line

plot_directed_graph(D_primal, mndwi, left_utm_x, right_utm_x, lower_utm_y, upper_utm_y, edge_path=False, arrow_width=100, alpha=1.0)

Plots the 'D_primal' directed graph on a background image.

The nodes of the graph are plotted as red dots, and the edges are plotted as red lines. If a main path is provided, it is plotted in blue. Arrows are also drawn to indicate the direction of the edges.

Parameters:

Name	Type	Description	Default
D_primal	DiGraph	The directed graph to plot.	required
mndwi	ndarray	A 2D array representing the background image.	required
left_utm_x	float	The leftmost x-coordinate of the extent of the image.	required
right_utm_x	float	The rightmost x-coordinate of the extent of the image.	required
lower_utm_y	float	The lowest y-coordinate of the extent of the image.	required
upper_utm_y	float	The highest y-coordinate of the extent of the image.	required
edge_path	list	A list of edges defining the main path to plot. Defaults to False.	False
arrow_width	int	The width of the arrows indicating the direction of the edges. Defaults to 100.	100
alpha	float	The transparency of the edges and the main path. Defaults to 1.0.	1.0

Returns:

Name	Type	Description
fig	Figure	The figure containing the plot.
ax	Axes	The axes of the plot.

Source code in rivabar/visualization.py

def plot_directed_graph(D_primal, mndwi, left_utm_x, right_utm_x, lower_utm_y, upper_utm_y, edge_path=False, arrow_width=100, alpha=1.0):
    """
    Plots the 'D_primal' directed graph on a background image.

    The nodes of the graph are plotted as red dots, and the edges are plotted as red lines. If a main path is provided, it is plotted in blue. Arrows are also drawn to indicate the direction of the edges.

    Parameters
    ----------
    D_primal : networkx.classes.digraph.DiGraph
        The directed graph to plot.
    mndwi : numpy.ndarray
        A 2D array representing the background image.
    left_utm_x : float
        The leftmost x-coordinate of the extent of the image.
    right_utm_x : float
        The rightmost x-coordinate of the extent of the image.
    lower_utm_y : float
        The lowest y-coordinate of the extent of the image.
    upper_utm_y : float
        The highest y-coordinate of the extent of the image.
    edge_path : list, optional
        A list of edges defining the main path to plot. Defaults to False.
    arrow_width : int, optional
        The width of the arrows indicating the direction of the edges. Defaults to 100.
    alpha : float, optional
        The transparency of the edges and the main path. Defaults to 1.0.

    Returns
    -------
    fig : matplotlib.figure.Figure
        The figure containing the plot.
    ax : matplotlib.axes.Axes
        The axes of the plot.
    """

    fig, ax = plt.subplots()
    plt.imshow(mndwi, extent=[left_utm_x, right_utm_x, lower_utm_y, upper_utm_y], cmap='gray_r', alpha=0.5)
    # plot nodes:
    for node in D_primal.nodes:
        plt.plot(D_primal.nodes()[node]['geometry'].xy[0], D_primal.nodes()[node]['geometry'].xy[1], 'ro')
        plt.text(D_primal.nodes()[node]['geometry'].xy[0][0]+100, D_primal.nodes()[node]['geometry'].xy[1][0], str(node))
    # plot edges:    
    for s,e,d in D_primal.edges:
        if 'geometry' in D_primal[s][e][d].keys():
            x = D_primal[s][e][d]['geometry'].xy[0]
            y = D_primal[s][e][d]['geometry'].xy[1]
            plt.plot(x, y, 'r', alpha=alpha)
    # plot main path:
    if edge_path:
        for s,e,d in edge_path:
            x = D_primal[s][e][d]['geometry'].xy[0]
            y = D_primal[s][e][d]['geometry'].xy[1]
            plt.plot(x, y, 'b', alpha=alpha)
    # plot arrows:
    for s,e,d in D_primal.edges:
        x1 = D_primal.nodes()[s]['geometry'].xy[0][0]
        y1 = D_primal.nodes()[s]['geometry'].xy[1][0]
        x2 = D_primal.nodes()[e]['geometry'].xy[0][0]
        y2 = D_primal.nodes()[e]['geometry'].xy[1][0]
        plt.arrow(x1, y1, x2-x1, y2-y1, width=arrow_width, length_includes_head=True, linewidth=None)
    return fig, ax

read_and_plot_im(dirname, fname)

Reads a raster image from a file and plots it using matplotlib.

Parameters:

Name	Type	Description	Default
dirname	str	The directory name where the raster file is located.	required
fname	str	The filename of the raster file.	required

Returns:

Name	Type	Description
im	ndarray	The raster image data.
dataset	DatasetReader	The dataset object containing metadata and other information about the raster.
left_utm_x	float	The UTM x-coordinate of the left edge of the raster.
right_utm_x	float	The UTM x-coordinate of the right edge of the raster.
lower_utm_y	float	The UTM y-coordinate of the lower edge of the raster.
upper_utm_y	float	The UTM y-coordinate of the upper edge of the raster.

Source code in rivabar/visualization.py

def read_and_plot_im(dirname, fname):
    """
    Reads a raster image from a file and plots it using matplotlib.

    Parameters
    ----------
    dirname : str
        The directory name where the raster file is located.
    fname : str
        The filename of the raster file.

    Returns
    -------
    im : numpy.ndarray
        The raster image data.
    dataset : rasterio.io.DatasetReader
        The dataset object containing metadata and other information about the raster.
    left_utm_x : float
        The UTM x-coordinate of the left edge of the raster.
    right_utm_x : float
        The UTM x-coordinate of the right edge of the raster.
    lower_utm_y : float
        The UTM y-coordinate of the lower edge of the raster.
    upper_utm_y : float
        The UTM y-coordinate of the upper edge of the raster.
    """
    with rasterio.open(dirname+fname) as dataset:
        im = dataset.read(1)
        left_utm_x = dataset.transform[2]
        upper_utm_y = dataset.transform[5]
        delta_x = dataset.transform[0]
        delta_y = dataset.transform[4]
        nxpix = im.shape[1]
        nypix = im.shape[0]
        right_utm_x = left_utm_x + delta_x*nxpix
        lower_utm_y = upper_utm_y + delta_y*nypix  
    plt.figure()
    plt.imshow(im, extent = [left_utm_x, right_utm_x, lower_utm_y, upper_utm_y], cmap='gray')
    return im, dataset, left_utm_x, right_utm_x, lower_utm_y, upper_utm_y

plot_graph_mappings(D_primal_t1, G_rook_t1, D_primal_t2, G_rook_t2, mappings, ax=None, plot_rook_nodes=True, plot_primal_nodes=True, plot_primal_edges=True)

Plots the results of graph mapping between two time steps. Args: D_primal_t1, G_rook_t1: Graphs from the first time step. D_primal_t2, G_rook_t2: Graphs from the second time step. mappings (dict): The dictionary returned by map_graphs_over_time. ax (matplotlib.axes.Axes, optional): Matplotlib axes to plot on. If None, a new figure is created. plot_rook_nodes (bool): Whether to plot G_rook node mappings. plot_primal_nodes (bool): Whether to plot D_primal node mappings. plot_primal_edges (bool): Whether to plot D_primal edge mappings. Returns: matplotlib.axes.Axes: The axes object with the plot.

Source code in rivabar/visualization.py

def plot_graph_mappings(D_primal_t1, G_rook_t1, D_primal_t2, G_rook_t2, mappings, ax=None,
                        plot_rook_nodes=True, plot_primal_nodes=True, plot_primal_edges=True):
    """
    Plots the results of graph mapping between two time steps.
    Args:
        D_primal_t1, G_rook_t1: Graphs from the first time step.
        D_primal_t2, G_rook_t2: Graphs from the second time step.
        mappings (dict): The dictionary returned by map_graphs_over_time.
        ax (matplotlib.axes.Axes, optional): Matplotlib axes to plot on. If None, a new figure is created.
        plot_rook_nodes (bool): Whether to plot G_rook node mappings.
        plot_primal_nodes (bool): Whether to plot D_primal node mappings.
        plot_primal_edges (bool): Whether to plot D_primal edge mappings.
    Returns:
        matplotlib.axes.Axes: The axes object with the plot.
    """
    if ax is None:
        fig, ax = plt.subplots(figsize=(15, 15))

    # Plot graphs from t1
    for n, data in G_rook_t1.nodes(data=True):
        poly = data.get('cl_polygon')
        if poly and not poly.is_empty:
            x, y = poly.exterior.xy
            ax.plot(x, y, color='blue', alpha=0.5, label='t1 Polygons' if 't1 Polygons' not in [l.get_label() for l in ax.get_lines()] else "")
    for u, v, k, data in D_primal_t1.edges(keys=True, data=True):
        line = data['geometry']
        x, y = line.xy
        ax.plot(x, y, color='cyan', alpha=0.7, linewidth=2, label='t1 Centerline' if 't1 Centerline' not in [l.get_label() for l in ax.get_lines()] else "")

    # Plot graphs from t2
    for n, data in G_rook_t2.nodes(data=True):
        poly = data.get('cl_polygon')
        if poly and not poly.is_empty:
            x, y = poly.exterior.xy
            ax.plot(x, y, color='red', alpha=0.5, label='t2 Polygons' if 't2 Polygons' not in [l.get_label() for l in ax.get_lines()] else "")
    for u, v, k, data in D_primal_t2.edges(keys=True, data=True):
        line = data['geometry']
        x, y = line.xy
        ax.plot(x, y, color='magenta', alpha=0.7, linewidth=2, label='t2 Centerline' if 't2 Centerline' not in [l.get_label() for l in ax.get_lines()] else "")

    # Plot mappings
    # G_rook node mappings
    if plot_rook_nodes and 'rook_nodes' in mappings:
        for n1, n2 in mappings['rook_nodes'].items():
            if G_rook_t1.has_node(n1) and G_rook_t2.has_node(n2):
                p1 = G_rook_t1.nodes[n1]['cl_polygon'].centroid
                p2 = G_rook_t2.nodes[n2]['cl_polygon'].centroid
                ax.plot([p1.x, p2.x], [p1.y, p2.y], 'k--', alpha=0.6, linewidth=0.8, label='Rook Node Mapping' if 'Rook Node Mapping' not in [l.get_label() for l in ax.get_lines()] else "")

    # D_primal node mappings
    if plot_primal_nodes and 'primal_nodes' in mappings:
        for n1, n2 in mappings['primal_nodes'].items():
            if D_primal_t1.has_node(n1) and D_primal_t2.has_node(n2):
                p1 = D_primal_t1.nodes[n1]['geometry']
                p2 = D_primal_t2.nodes[n2]['geometry']
                ax.plot([p1.x, p2.x], [p1.y, p2.y], 'g-', alpha=0.7, linewidth=1, label='Primal Node Mapping' if 'Primal Node Mapping' not in [l.get_label() for l in ax.get_lines()] else "")

    # D_primal edge mappings
    if plot_primal_edges and 'primal_edges' in mappings:
        for e1, e2 in mappings['primal_edges'].items():
            if D_primal_t1.has_edge(*e1) and D_primal_t2.has_edge(*e2):
                line1 = D_primal_t1.edges[e1]['geometry']
                line2 = D_primal_t2.edges[e2]['geometry']
                c1 = line1.centroid
                c2 = line2.centroid
                ax.plot([c1.x, c2.x], [c1.y, c2.y], 'y:', alpha=0.8, linewidth=1.2, label='Primal Edge Mapping' if 'Primal Edge Mapping' not in [l.get_label() for l in ax.get_lines()] else "")

    ax.set_title("Graph Mappings Over Time")
    ax.legend()
    ax.set_aspect('equal', adjustable='box')
    return ax

plot_deviation_rose_diagram(deviations, classifications, weights=None, bins=12)

Plots a rose diagram (polar histogram) of node displacement deviations, separated for confluences and splits.

Args: deviations (dict): A dictionary of node IDs to deviation angles (in degrees). classifications (dict): A dictionary of node IDs to their classification. weights (dict, optional): A dictionary of node IDs to displacement lengths to weight the histogram. bins (int): The number of bins for the histogram.

Returns: matplotlib.figure.Figure, array of matplotlib.axes.Axes: The figure and axes objects.

Source code in rivabar/visualization.py

def plot_deviation_rose_diagram(deviations, classifications, weights=None, bins=12):
    """
    Plots a rose diagram (polar histogram) of node displacement deviations,
    separated for confluences and splits.

    Args:
        deviations (dict): A dictionary of node IDs to deviation angles (in degrees).
        classifications (dict): A dictionary of node IDs to their classification.
        weights (dict, optional): A dictionary of node IDs to displacement lengths to weight the histogram.
        bins (int): The number of bins for the histogram.

    Returns:
        matplotlib.figure.Figure, array of matplotlib.axes.Axes: The figure and axes objects.
    """
    confluence_deviations = [d for n, d in deviations.items() if classifications.get(n) == 'confluence' and d is not None]
    split_deviations = [d for n, d in deviations.items() if classifications.get(n) == 'split' and d is not None]

    confluence_weights = None
    if weights:
        confluence_weights = [weights.get(n) for n, d in deviations.items() if classifications.get(n) == 'confluence' and d is not None]

    split_weights = None
    if weights:
        split_weights = [weights.get(n) for n, d in deviations.items() if classifications.get(n) == 'split' and d is not None]

    fig, axes = plt.subplots(1, 2, figsize=(15, 7), subplot_kw={'projection': 'polar'})

    # Convert degrees to radians
    confluence_rad = np.radians(confluence_deviations)
    split_rad = np.radians(split_deviations)

    # Confluence plot
    if confluence_rad.size > 0:
        counts, bin_edges = np.histogram(confluence_rad, bins=np.linspace(0, 2 * np.pi, bins + 1), weights=confluence_weights)
        widths = np.diff(bin_edges)
        axes[0].bar(bin_edges[:-1], counts, width=widths, edgecolor='k', alpha=0.7, zorder=2)
    axes[0].set_title('Confluence Deviations')
    axes[0].set_theta_zero_location('N')
    axes[0].set_theta_direction(-1)
    axes[0].grid(True, zorder=0)

    # Split plot
    if split_rad.size > 0:
        counts, bin_edges = np.histogram(split_rad, bins=np.linspace(0, 2 * np.pi, bins + 1), weights=split_weights)
        widths = np.diff(bin_edges)
        axes[1].bar(bin_edges[:-1], counts, width=widths, edgecolor='k', alpha=0.7, color='C1', zorder=2)
    axes[1].set_title('Split Deviations')
    axes[1].set_theta_zero_location('N')
    axes[1].set_theta_direction(-1)
    axes[1].grid(True, zorder=0)

    fig.tight_layout()
    return fig, axes

plot_deviation_histogram(deviations, classifications, distances, bins=18)

Plots a histogram of node displacement deviations, separated for confluences and splits. Only includes nodes with a displacement distance greater than zero.

Args: deviations (dict): A dictionary of node IDs to deviation angles (in degrees). classifications (dict): A dictionary of node IDs to their classification. distances (dict): A dictionary of node IDs to displacement lengths. bins (int): The number of bins for the histogram (default is 18 for 10-degree bins).

Returns: matplotlib.figure.Figure, array of matplotlib.axes.Axes: The figure and axes objects.

Source code in rivabar/visualization.py

def plot_deviation_histogram(deviations, classifications, distances, bins=18):
    """
    Plots a histogram of node displacement deviations, separated for confluences and splits.
    Only includes nodes with a displacement distance greater than zero.

    Args:
        deviations (dict): A dictionary of node IDs to deviation angles (in degrees).
        classifications (dict): A dictionary of node IDs to their classification.
        distances (dict): A dictionary of node IDs to displacement lengths.
        bins (int): The number of bins for the histogram (default is 18 for 10-degree bins).

    Returns:
        matplotlib.figure.Figure, array of matplotlib.axes.Axes: The figure and axes objects.
    """
    confluence_deviations = [
        d for n, d in deviations.items()
        if d is not None and classifications.get(n) == 'confluence' and distances.get(n, 0) > 0
    ]
    split_deviations = [
        d for n, d in deviations.items()
        if d is not None and classifications.get(n) == 'split' and distances.get(n, 0) > 0
    ]

    fig, axes = plt.subplots(1, 2, figsize=(16, 6), sharey=True)

    # Confluence plot
    axes[0].hist(confluence_deviations, bins=bins, range=(0, 180), edgecolor='k')
    axes[0].set_xlim(0, 180)
    axes[0].set_title('Confluence Deviations')
    axes[0].set_xlabel('Deviation Angle (degrees)')
    axes[0].set_ylabel('Frequency')
    axes[0].grid(True, alpha=0.5)

    # Split plot
    axes[1].hist(split_deviations, bins=bins, range=(0, 180), edgecolor='k', color='C1')
    axes[1].set_xlim(0, 180)
    axes[1].set_title('Split Deviations')
    axes[1].set_xlabel('Deviation Angle (degrees)')
    axes[1].grid(True, alpha=0.5)

    fig.tight_layout()
    return fig, axes

plot_river_segments(rivers, common_confluences, segment_groups=None, ncols=4, figsize_per_panel=(4, 5), segment_cmap='tab10', rejected_color='0.75', confluence_color='k', confluence_marker='o', confluence_size=40, linewidth=2, show_mndwi=True, mndwi_cmap='Greys', mndwi_alpha=0.4, show_tributaries=True, tributary_color='0.5', tributary_linewidth=1.0)

Plot main-path segments for multiple rivers in a multi-panel figure.

Each panel shows one river's main path colored by segment (split at the common confluence points), with the confluence locations marked.

If segment_groups (output of :func:match_river_segments) is provided, only matched segments are colored; rejected segments are drawn in grey. Segment colors are keyed to the segment_index in each group so that the same physical reach has the same color across all panels.

Parameters:

Name	Type	Description	Default
rivers	list of River	Processed River objects.	required
common_confluences	list of dict	Output of :func:find_common_confluences. Each dict must have a 'utm_coords' key with (x, y) UTM coordinates.	required
segment_groups	list of dict	Output of :func:match_river_segments. When provided, colors are assigned by segment_index and only matched river/segment pairs are colored. Unmatched segments are drawn in rejected_color.	None
ncols	int	Number of columns in the panel grid (default 4).	4
figsize_per_panel	tuple	(width, height) in inches per panel (default (4, 5)).	(4, 5)
segment_cmap	str	Matplotlib colormap name for segment colours (default 'tab10').	'tab10'
rejected_color	str or tuple	Color for segments that were not matched (default '0.75', light grey).	'0.75'
confluence_color	str	Marker colour for confluence points (default 'k').	'k'
confluence_marker	str	Marker style for confluence points (default 'o').	'o'
confluence_size	float	Marker size for confluence points (default 40).	40
linewidth	float	Line width for centerline segments (default 2).	2
show_mndwi	bool	If True, show the MNDWI water mask as a background (default True).	True
mndwi_cmap	str	Colormap for the MNDWI background (default 'Greys').	'Greys'
mndwi_alpha	float	Alpha for the MNDWI background (default 0.4).	0.4
show_tributaries	bool	If True, plot tributary branches from each river for context (default True).	True
tributary_color	str or tuple	Color for tributary lines (default '0.5', medium grey).	'0.5'
tributary_linewidth	float	Line width for tributary branches (default 1.0).	1.0

Returns:

Type	Description
fig, axes : matplotlib Figure and array of Axes

Source code in rivabar/visualization.py

def plot_river_segments(rivers, common_confluences, segment_groups=None,
                        ncols=4, figsize_per_panel=(4, 5),
                        segment_cmap='tab10', rejected_color='0.75',
                        confluence_color='k',
                        confluence_marker='o', confluence_size=40, linewidth=2,
                        show_mndwi=True, mndwi_cmap='Greys', mndwi_alpha=0.4,
                        show_tributaries=True, tributary_color='0.5',
                        tributary_linewidth=1.0):
    """
    Plot main-path segments for multiple rivers in a multi-panel figure.

    Each panel shows one river's main path colored by segment (split at the
    common confluence points), with the confluence locations marked.

    If *segment_groups* (output of :func:`match_river_segments`) is provided,
    only matched segments are colored; rejected segments are drawn in grey.
    Segment colors are keyed to the ``segment_index`` in each group so that
    the same physical reach has the same color across all panels.

    Parameters
    ----------
    rivers : list of River
        Processed River objects.
    common_confluences : list of dict
        Output of :func:`find_common_confluences`. Each dict must have a
        ``'utm_coords'`` key with ``(x, y)`` UTM coordinates.
    segment_groups : list of dict, optional
        Output of :func:`match_river_segments`. When provided, colors are
        assigned by ``segment_index`` and only matched river/segment pairs
        are colored. Unmatched segments are drawn in *rejected_color*.
    ncols : int, optional
        Number of columns in the panel grid (default 4).
    figsize_per_panel : tuple, optional
        ``(width, height)`` in inches per panel (default ``(4, 5)``).
    segment_cmap : str, optional
        Matplotlib colormap name for segment colours (default ``'tab10'``).
    rejected_color : str or tuple, optional
        Color for segments that were not matched (default ``'0.75'``, light
        grey).
    confluence_color : str, optional
        Marker colour for confluence points (default ``'k'``).
    confluence_marker : str, optional
        Marker style for confluence points (default ``'o'``).
    confluence_size : float, optional
        Marker size for confluence points (default 40).
    linewidth : float, optional
        Line width for centerline segments (default 2).
    show_mndwi : bool, optional
        If True, show the MNDWI water mask as a background (default True).
    mndwi_cmap : str, optional
        Colormap for the MNDWI background (default ``'Greys'``).
    mndwi_alpha : float, optional
        Alpha for the MNDWI background (default 0.4).
    show_tributaries : bool, optional
        If True, plot tributary branches from each river for context
        (default True).
    tributary_color : str or tuple, optional
        Color for tributary lines (default ``'0.5'``, medium grey).
    tributary_linewidth : float, optional
        Line width for tributary branches (default 1.0).

    Returns
    -------
    fig, axes : matplotlib Figure and array of Axes
    """
    split_points = [c['utm_coords'] for c in common_confluences]

    # Build a lookup: river id -> {segment_index -> path} for matched segments
    matched_lookup = {}  # river id -> {seg_index: path}
    if segment_groups is not None:
        for g in segment_groups:
            seg_idx = g['segment_index']
            for river, path in zip(g['rivers'], g['paths']):
                rid = id(river)
                if rid not in matched_lookup:
                    matched_lookup[rid] = {}
                matched_lookup[rid][seg_idx] = path

    # Filter to valid rivers
    valid_rivers = [r for r in rivers if r._is_processed and r._processing_successful
                    and r.main_path is not None]
    if not valid_rivers:
        print("No valid processed rivers to plot.")
        return None, None

    nrows = int(np.ceil(len(valid_rivers) / ncols))
    fig_w = figsize_per_panel[0] * ncols
    fig_h = figsize_per_panel[1] * nrows
    fig, axes = plt.subplots(nrows, ncols, figsize=(fig_w, fig_h), squeeze=False)

    cmap = plt.get_cmap(segment_cmap)

    # Build a stable color mapping: segment_index -> color
    if segment_groups is not None:
        all_seg_indices = [g['segment_index'] for g in segment_groups]
        seg_color_map = {si: cmap(i % cmap.N) for i, si in enumerate(all_seg_indices)}

    for idx, river in enumerate(valid_rivers):
        row, col = divmod(idx, ncols)
        ax = axes[row][col]

        # Background
        if show_mndwi and river._mndwi is not None:
            ax.imshow(river._mndwi, extent=[river._left_utm_x, river._right_utm_x,
                                             river._lower_utm_y, river._upper_utm_y],
                      cmap=mndwi_cmap, alpha=mndwi_alpha, aspect='equal')

        rid = id(river)

        if segment_groups is not None and rid in matched_lookup:
            # Plot matched segments in color, everything else in grey
            # First draw the full path in grey as background
            for s, e, d in river.main_path:
                geom = river._D_primal[s][e][d]['geometry']
                xs, ys = geom.xy
                ax.plot(xs, ys, color=rejected_color, linewidth=linewidth)
            # Then overlay matched segments in their group color
            for seg_idx, path in matched_lookup[rid].items():
                color = seg_color_map[seg_idx]
                for s, e, d in path:
                    geom = river._D_primal[s][e][d]['geometry']
                    xs, ys = geom.xy
                    ax.plot(xs, ys, color=color, linewidth=linewidth + 0.5)
        else:
            # No segment_groups provided, or river not in any group:
            # fall back to coloring by split position
            try:
                segments, split_info = river.split_main_path_at_points(split_points)
            except Exception:
                segments = [river.main_path]

            if segment_groups is not None:
                # River didn't match — draw all grey
                for seg in segments:
                    for s, e, d in seg:
                        geom = river._D_primal[s][e][d]['geometry']
                        xs, ys = geom.xy
                        ax.plot(xs, ys, color=rejected_color, linewidth=linewidth)
            else:
                # No segment_groups at all — simple coloring by index
                for seg_idx, seg in enumerate(segments):
                    color = cmap(seg_idx % cmap.N)
                    for s, e, d in seg:
                        geom = river._D_primal[s][e][d]['geometry']
                        xs, ys = geom.xy
                        ax.plot(xs, ys, color=color, linewidth=linewidth)

        # Plot tributary branches for context
        if show_tributaries:
            for trib in river.tributary_confluences:
                if 'branch_utm_coords' in trib:
                    coords = trib['branch_utm_coords']
                    ax.plot(coords[:, 0], coords[:, 1], '-',
                            color=tributary_color, linewidth=tributary_linewidth)

        # Mark confluence points
        for c in common_confluences:
            cx, cy = c['utm_coords']
            ax.scatter(cx, cy, c=confluence_color, marker=confluence_marker,
                       s=confluence_size, zorder=5)

        # Label
        label = getattr(river, 'acquisition_date', None)
        if label is None:
            label = river.fname
        ax.set_title(str(label), fontsize=9)
        ax.set_aspect('equal')
        ax.tick_params(labelsize=7)

    # Compute global extent across all panels and apply uniform limits
    all_xlims = []
    all_ylims = []
    for idx in range(len(valid_rivers)):
        row, col = divmod(idx, ncols)
        ax = axes[row][col]
        all_xlims.append(ax.get_xlim())
        all_ylims.append(ax.get_ylim())
    if all_xlims:
        global_xlim = (min(xl[0] for xl in all_xlims), max(xl[1] for xl in all_xlims))
        global_ylim = (min(yl[0] for yl in all_ylims), max(yl[1] for yl in all_ylims))
        for idx in range(len(valid_rivers)):
            row, col = divmod(idx, ncols)
            axes[row][col].set_xlim(global_xlim)
            axes[row][col].set_ylim(global_ylim)

    # Hide unused panels
    for idx in range(len(valid_rivers), nrows * ncols):
        row, col = divmod(idx, ncols)
        axes[row][col].set_visible(False)

    fig.tight_layout()
    return fig, axes

plot_centerline_comparison(river, delta_s=100, smoothing_factor=1000000.0, path=None, pixel_size=None, ax=None, figsize=(14, 10))

Compare raw (D_primal) and smoothed centerlines against the banklines.

Plots three layers on top of the water mask: - banklines (blue) - raw centerline from D_primal geometry (black) - smoothed/resampled centerline from resample_and_smooth (red)

Parameters:

Name	Type	Description	Default
river	River	A processed River object.	required
delta_s	float	Resampling distance for smooth centerline (default 100).	100
smoothing_factor	float	Spline smoothing factor (default 1e6).	1000000.0
path	list of tuples	Edge path to use. If None, uses river.main_path.	None
pixel_size	float	Pixel size in metres. Inferred from river._dataset if None.	None
ax	Axes	Axes to plot on. Creates a new figure if None.	None
figsize	tuple	Figure size when creating a new figure (default (14, 10)).	(14, 10)

Returns:

Type	Description
(fig, ax)

Source code in rivabar/visualization.py

def plot_centerline_comparison(river, delta_s=100, smoothing_factor=1e6,
                                path=None, pixel_size=None, ax=None,
                                figsize=(14, 10)):
    """
    Compare raw (D_primal) and smoothed centerlines against the banklines.

    Plots three layers on top of the water mask:
    - banklines (blue)
    - raw centerline from D_primal geometry (black)
    - smoothed/resampled centerline from resample_and_smooth (red)

    Parameters
    ----------
    river : River
        A processed River object.
    delta_s : float, optional
        Resampling distance for smooth centerline (default 100).
    smoothing_factor : float, optional
        Spline smoothing factor (default 1e6).
    path : list of tuples, optional
        Edge path to use. If None, uses river.main_path.
    pixel_size : float, optional
        Pixel size in metres.  Inferred from river._dataset if None.
    ax : matplotlib.Axes, optional
        Axes to plot on.  Creates a new figure if None.
    figsize : tuple, optional
        Figure size when creating a new figure (default (14, 10)).

    Returns
    -------
    fig, ax
    """
    from .utils import resample_and_smooth

    if path is None:
        path = river.main_path
    if pixel_size is None:
        if river._dataset is not None:
            pixel_size = river._dataset.transform[0]
        else:
            pixel_size = 30.0

    # --- Extract raw centerline from D_primal edges ---
    x_raw, y_raw = [], []
    for s, e, d in path:
        xc = list(river._D_primal[s][e][d]['geometry'].xy[0])
        yc = list(river._D_primal[s][e][d]['geometry'].xy[1])
        x_raw += xc
        y_raw += yc
    x_raw, y_raw = np.array(x_raw), np.array(y_raw)

    # --- Smoothed centerline (same logic as get_width_and_curvature) ---
    x_smooth, y_smooth = resample_and_smooth(
        x_raw, y_raw, delta_s, smoothing_factor)

    # --- Banklines ---
    banks = river.main_channel_banks
    bank1 = banks['left_bank'] if banks is not None else None
    bank2 = banks['right_bank'] if banks is not None else None

    # --- Plot ---
    if ax is None:
        fig, ax = plt.subplots(figsize=figsize)
    else:
        fig = ax.get_figure()

    # Water mask background
    if river._mndwi is not None:
        ax.imshow(river._mndwi,
                  extent=[river._left_utm_x, river._right_utm_x,
                          river._lower_utm_y, river._upper_utm_y],
                  cmap='Grays', alpha=0.5)

    # Banklines
    if bank1 is not None:
        ax.plot(*bank1.xy, color='tab:blue', linewidth=1.2, label='banklines')
    if bank2 is not None:
        ax.plot(*bank2.xy, color='tab:blue', linewidth=1.2)

    # Raw centerline
    ax.plot(x_raw, y_raw, 'k-', linewidth=1.0, alpha=0.7,
            label='raw centerline (D_primal)')

    # Smoothed centerline
    ax.plot(x_smooth, y_smooth, 'r-', linewidth=1.5,
            label=f'smoothed (s={smoothing_factor:.0e}, Δs={delta_s})')

    ax.set_aspect('equal')
    ax.legend(loc='best', fontsize=9)
    ax.set_title(f'{river.fname}  –  raw vs smoothed centerline')

    return fig, ax

plot_pair(results, pair_idx, segment_group=None, km_interval=20, ax=None, figsize=(14, 10))

Plot banklines and smoothed centerlines for a pair from segment analysis results.

Parameters:

Name	Type	Description	Default
results	dict	Output of :func:~rivabar.utils.analyze_segment_group.	required
pair_idx	int	Index into results['pair_info'].	required
segment_group	dict	The segment group dict (from :func:~rivabar.temporal_analysis.match_river_segments). Used to look up river objects when results['rivers'] is not available (older results). If None, rivers are taken from results['rivers'].	None
km_interval	float	Spacing (in km) for distance markers along centerline 1 (default 20). Set to 0 or None to skip markers.	20
ax	Axes	Axes to plot on. Creates a new figure if None.	None
figsize	tuple	Figure size when creating a new figure (default (14, 10)).	(14, 10)

Returns:

Type	Description
(fig, ax)

Source code in rivabar/visualization.py

def plot_pair(results, pair_idx, segment_group=None, km_interval=20,
              ax=None, figsize=(14, 10)):
    """
    Plot banklines and smoothed centerlines for a pair from segment analysis results.

    Parameters
    ----------
    results : dict
        Output of :func:`~rivabar.utils.analyze_segment_group`.
    pair_idx : int
        Index into ``results['pair_info']``.
    segment_group : dict, optional
        The segment group dict (from :func:`~rivabar.temporal_analysis.match_river_segments`).
        Used to look up river objects when ``results['rivers']`` is not
        available (older results).  If *None*, rivers are taken from
        ``results['rivers']``.
    km_interval : float, optional
        Spacing (in km) for distance markers along centerline 1 (default 20).
        Set to 0 or *None* to skip markers.
    ax : matplotlib.Axes, optional
        Axes to plot on.  Creates a new figure if *None*.
    figsize : tuple, optional
        Figure size when creating a new figure (default (14, 10)).

    Returns
    -------
    fig, ax
    """
    pair_info = results['pair_info'][pair_idx]
    cc = results['centerline_coords'][pair_idx]

    # Resolve river objects ------------------------------------------------
    river1_idx = pair_info['river1_index']
    river2_idx = pair_info['river2_index']

    if 'rivers' in results:
        sg_rivers = results['rivers']
    elif segment_group is not None:
        sg_rivers = segment_group['rivers']
    else:
        raise ValueError(
            "Cannot resolve river objects: results has no 'rivers' key "
            "and no segment_group was provided.  Re-run "
            "analyze_segment_group to get updated results, or pass the "
            "segment_group explicitly.")

    river1 = sg_rivers[river1_idx]
    river2 = sg_rivers[river2_idx]

    # Create axes ----------------------------------------------------------
    if ax is None:
        fig, ax = plt.subplots(figsize=figsize)
    else:
        fig = ax.get_figure()

    # Channel polygons (filled) ---------------------------------------------
    from .polygon_processing import create_channel_nw_polygon
    for river, color in [(river1, 'tab:blue'), (river2, 'tab:red')]:
        try:
            ch_poly = create_channel_nw_polygon(
                river._G_rook, buffer=10, dataset=river._dataset)
            if isinstance(ch_poly, Polygon):
                # Plot exterior
                ax.fill(*ch_poly.exterior.xy, color=color, alpha=0.15)
                ax.plot(*ch_poly.exterior.xy, color=color, alpha=0.4,
                        linewidth=0.8)
                # Plot island holes
                for interior in ch_poly.interiors:
                    ax.fill(*interior.xy, color='white', alpha=1.0)
                    ax.plot(*interior.xy, color=color, alpha=0.4,
                            linewidth=0.8)
        except Exception:
            # Fallback: plot bank polygon outlines
            for node in range(len(river._G_rook.nodes())):
                poly = river._G_rook.nodes()[node]['bank_polygon']
                if isinstance(poly, Polygon):
                    ax.plot(*poly.exterior.xy, color=color, alpha=0.4,
                            linewidth=0.8)
                else:
                    ax.plot(*poly.xy, color=color, alpha=0.4,
                            linewidth=0.8)

    # Smoothed centerlines from the pair results ---------------------------
    d1 = str(pair_info['date1'])[:10]
    d2 = str(pair_info['date2'])[:10]

    ax.plot(cc['x1'], cc['y1'], color='tab:blue', linewidth=1.5,
            label=f'{d1} (river 1)')
    ax.plot(cc['x2'], cc['y2'], color='tab:red', linewidth=1.5,
            label=f'{d2} (river 2)')

    # Distance markers along centerline 1 ---------------------------------
    if km_interval:
        s1 = pair_info['s1']
        interval_m = km_interval * 1000
        ref_dists = np.arange(0, s1[-1], interval_m)
        for rd in ref_dists:
            idx = np.argmin(np.abs(s1 - rd))
            ax.plot(cc['x1'][idx], cc['y1'][idx], 'ko', markersize=7,
                    markerfacecolor='yellow', markeredgecolor='black',
                    zorder=5)
            ax.text(cc['x1'][idx], cc['y1'][idx], f'{rd/1000:.0f} km',
                    fontsize=9, ha='center', va='bottom',
                    bbox=dict(boxstyle='round,pad=0.2', facecolor='white',
                              alpha=0.8),
                    zorder=6)

    ax.set_aspect('equal')
    ax.legend(loc='best', fontsize=9)
    gap_yr = pair_info['time_gap_years']
    ax.set_title(f'Pair {pair_idx}: {d1} vs {d2}  '
                 f'({gap_yr:.1f} yr, segment {results.get("segment_index", "")})')

    return fig, ax

interactive_scatter(df, x, y, c=None, cmap=None, title=None, hover_cols=None)

Create an interactive scatter plot for selecting data points from a DataFrame.

Uses plotly for reliable interactive selection in Jupyter notebooks. Click on points to see their details in the hover tooltip. Use the box select or lasso select tools in the plotly toolbar to select multiple points. The selected row indices are returned.

Parameters:

Name	Type	Description	Default
df	DataFrame	The DataFrame containing the data (e.g. output of create_dataframe_from_results).	required
x	str	Column name for the x-axis.	required
y	str	Column name for the y-axis.	required
c	str or None	Column name used to colour the points. Works with both numeric and categorical columns.	None
cmap	dict or None	When c is categorical, a {category: color} mapping. Unmapped categories fall back to gray.	None
title	str or None	Plot title.	None
hover_cols	list of str or None	Extra column names to show in the hover tooltip. If None, shows date1, date2, and the colour column (if any).	None

Returns:

Name	Type	Description
selector	object	A selector object with a .selected attribute (list of int) that updates as you select points. Access selector.selected after making your selection.

Examples:

>>> sel = rb.interactive_scatter(
...     df, 'width_diff (m)', 'r_squared_curv_mr',
...     c='pair_class', title='Purus 2')
>>> # Use box/lasso select in the plotly toolbar, then:
>>> sel.selected  # [3, 17, 42]

Source code in rivabar/visualization.py

def interactive_scatter(df, x, y, c=None, cmap=None, title=None,
                        hover_cols=None):
    """
    Create an interactive scatter plot for selecting data points from a DataFrame.

    Uses plotly for reliable interactive selection in Jupyter notebooks.
    Click on points to see their details in the hover tooltip. Use the
    **box select** or **lasso select** tools in the plotly toolbar to
    select multiple points.  The selected row indices are returned.

    Parameters
    ----------
    df : pandas.DataFrame
        The DataFrame containing the data (e.g. output of
        ``create_dataframe_from_results``).
    x : str
        Column name for the x-axis.
    y : str
        Column name for the y-axis.
    c : str or None, optional
        Column name used to colour the points.  Works with both numeric
        and categorical columns.
    cmap : dict or None, optional
        When *c* is categorical, a ``{category: color}`` mapping.
        Unmapped categories fall back to gray.
    title : str or None, optional
        Plot title.
    hover_cols : list of str or None, optional
        Extra column names to show in the hover tooltip.  If *None*,
        shows ``date1``, ``date2``, and the colour column (if any).

    Returns
    -------
    selector : object
        A selector object with a ``.selected`` attribute (list of int)
        that updates as you select points.  Access
        ``selector.selected`` after making your selection.

    Examples
    --------
    >>> sel = rb.interactive_scatter(
    ...     df, 'width_diff (m)', 'r_squared_curv_mr',
    ...     c='pair_class', title='Purus 2')
    >>> # Use box/lasso select in the plotly toolbar, then:
    >>> sel.selected  # [3, 17, 42]
    """
    import plotly.graph_objects as go
    import plotly.express as px

    # Build hover text
    if hover_cols is None:
        hover_cols = []
        for col_name in ['date1', 'date2']:
            if col_name in df.columns:
                hover_cols.append(col_name)
        if c is not None and c in df.columns and c not in hover_cols:
            hover_cols.append(c)

    hover_template = f"<b>{x}</b>: %{{x:.3f}}<br><b>{y}</b>: %{{y:.3f}}"
    custom_data_cols = []
    for i, col_name in enumerate(hover_cols):
        if col_name in df.columns:
            custom_data_cols.append(col_name)
            hover_template += f"<br><b>{col_name}</b>: %{{customdata[{i}]}}"
    hover_template += "<br><b>row</b>: %{text}<extra></extra>"

    is_categorical = (c is not None and c in df.columns
                      and (df[c].dtype == object or hasattr(df[c], 'cat')))

    if is_categorical:
        categories = sorted(df[c].unique(), key=str)
        if cmap is None:
            colors_px = px.colors.qualitative.Plotly
            cmap = {cat: colors_px[i % len(colors_px)]
                    for i, cat in enumerate(categories)}

        fig = go.Figure()
        for cat in categories:
            mask = df[c] == cat
            sub = df[mask]
            fig.add_trace(go.Scatter(
                x=sub[x], y=sub[y],
                mode='markers',
                name=str(cat),
                marker=dict(size=8, color=cmap.get(cat, 'gray'),
                            line=dict(width=0.5, color='black')),
                text=sub.index.astype(str),
                customdata=sub[custom_data_cols].values if custom_data_cols else None,
                hovertemplate=hover_template,
                selectedpoints=[],
            ))
    else:
        marker_kw = dict(size=8, line=dict(width=0.5, color='black'))
        if c is not None and c in df.columns:
            marker_kw['color'] = df[c]
            marker_kw['colorscale'] = 'Viridis'
            marker_kw['colorbar'] = dict(title=c)

        fig = go.Figure(go.Scatter(
            x=df[x], y=df[y],
            mode='markers',
            marker=marker_kw,
            text=df.index.astype(str),
            customdata=df[custom_data_cols].values if custom_data_cols else None,
            hovertemplate=hover_template,
        ))

    fig.update_layout(
        title=title,
        xaxis_title=x,
        yaxis_title=y,
        dragmode='lasso',
        width=900, height=600,
    )

    widget = go.FigureWidget(fig)

    # Store row indices on each trace for retrieval
    _trace_row_indices = []
    for trace in widget.data:
        _trace_row_indices.append([int(t) for t in trace.text])

    class Selector:
        """Reads selected points from the FigureWidget traces."""

        def __init__(self, widget, trace_row_indices):
            self._widget = widget
            self._trace_row_indices = trace_row_indices

        @property
        def selected(self):
            """Return DataFrame row indices of currently selected points."""
            result = []
            for trace, row_indices in zip(
                    self._widget.data, self._trace_row_indices):
                sel_pts = trace.selectedpoints
                if sel_pts:
                    for idx in sel_pts:
                        result.append(row_indices[idx])
            return result

    selector = Selector(widget, _trace_row_indices)

    from IPython.display import display
    display(widget)

    print("Use lasso/box select in the toolbar to select points. "
          "Then read selector.selected to get the row indices.")
    return selector

plot_prediction_map(results, pair_idx, calibration, prediction=None, delta_s=50, pixel_size=None, ax=None, cmap='RdBu_r', vmax=None)

Map-view visualization of prediction residuals along the channel.

Colors the channel centerline by the difference between predicted and observed migration rate (red = overprediction, blue = underprediction). Also shows the observed and predicted centerlines.

Parameters:

Name	Type	Description	Default
results	dict	Results from analyze_river_pairs_filtered or analyze_segment_group.	required
pair_idx	int	Index of the pair to visualize.	required
calibration	dict	Output of calibrate_segment.	required
prediction	dict	Output of predict_forward for this pair's river1. If None, it is computed internally.	None
delta_s	float	Resampling interval (default 50).	50
pixel_size	float	Pixel size for width conversion.	None
ax	Axes	Axes to plot on. If None, a new figure is created.	None
cmap	str	Colormap for residuals (default 'RdBu_r').	'RdBu_r'
vmax	float	Symmetric color scale limit. If None, uses the 95th percentile of absolute residuals.	None

Returns:

Name	Type	Description
ax	Axes

Source code in rivabar/visualization.py

def plot_prediction_map(results, pair_idx, calibration, prediction=None,
                        delta_s=50, pixel_size=None, ax=None,
                        cmap='RdBu_r', vmax=None):
    """
    Map-view visualization of prediction residuals along the channel.

    Colors the channel centerline by the difference between predicted and
    observed migration rate (red = overprediction, blue = underprediction).
    Also shows the observed and predicted centerlines.

    Parameters
    ----------
    results : dict
        Results from ``analyze_river_pairs_filtered`` or
        ``analyze_segment_group``.
    pair_idx : int
        Index of the pair to visualize.
    calibration : dict
        Output of ``calibrate_segment``.
    prediction : dict, optional
        Output of ``predict_forward`` for this pair's river1. If *None*,
        it is computed internally.
    delta_s : float, optional
        Resampling interval (default 50).
    pixel_size : float, optional
        Pixel size for width conversion.
    ax : matplotlib.axes.Axes, optional
        Axes to plot on. If *None*, a new figure is created.
    cmap : str, optional
        Colormap for residuals (default 'RdBu_r').
    vmax : float, optional
        Symmetric color scale limit. If *None*, uses the 95th percentile
        of absolute residuals.

    Returns
    -------
    ax : matplotlib.axes.Axes
    """
    from .prediction import (nominal_migration_rate,
                             predicted_migration_rate)

    pair_info = results['pair_info'][pair_idx]
    curvature = results['curvatures'][pair_idx]
    distances = results['migration_distances'][pair_idx]
    s = results['along_channel_distances'][pair_idx]
    time_gap_years = pair_info['time_gap_years']
    W = np.mean(pair_info['width1'])
    observed_mr = distances / time_gap_years

    # centerline_coords x1/y1 are the same arrays as curvature/distances
    coords = results['centerline_coords'][pair_idx]
    x1, y1 = np.asarray(coords['x1']), np.asarray(coords['y1'])
    x2, y2 = np.asarray(coords['x2']), np.asarray(coords['y2'])

    # Compute predicted MR using calibrated parameters
    D = calibration['D']
    Cf = calibration['Cf_median']
    kl = calibration['kl_median']
    Omega = calibration.get('Omega', -1.0)
    Gamma = calibration.get('Gamma', 2.5)

    R0 = nominal_migration_rate(curvature, W, kl)
    R1 = predicted_migration_rate(R0, s, D, Cf, Omega, Gamma)

    # Residual = observed - predicted; positive means the model
    # underpredicted (actual migration was larger than predicted)
    residual = observed_mr - R1

    # Mask to stable segments only
    segment_results = pair_info.get('segment_results', [])
    if segment_results:
        mask = np.zeros(len(curvature), dtype=bool)
        for seg in segment_results:
            mask[seg['start_idx']:seg['end_idx'] + 1] = True
    else:
        mask = np.ones(len(curvature), dtype=bool)

    if ax is None:
        fig, ax = plt.subplots(figsize=(12, 12))

    if vmax is None:
        vmax = np.percentile(np.abs(residual[mask]), 95)

    # Plot residual colors on the predicted centerline if available,
    # otherwise fall back to x1/y1
    # Color the river1 centerline by residual — x1/y1, residual, and
    # mask are all on the same grid, no interpolation needed
    points = np.column_stack([x1, y1]).reshape(-1, 1, 2)
    segments = np.concatenate([points[:-1], points[1:]], axis=1)
    seg_colors = (residual[:-1] + residual[1:]) / 2

    # Set unstable segments to gray
    seg_mask = mask[:-1] & mask[1:]
    seg_colors_masked = np.where(seg_mask, seg_colors, 0.0)

    lc = LineCollection(segments, cmap=cmap,
                        norm=Normalize(vmin=-vmax, vmax=vmax),
                        linewidths=3)
    lc.set_array(seg_colors_masked)
    colors = lc.to_rgba(seg_colors_masked)
    colors[~seg_mask] = [0.7, 0.7, 0.7, 1.0]
    lc.set_colors(colors)
    ax.add_collection(lc)

    # Plot observed centerlines
    ax.plot(x2, y2, 'k-', lw=1, alpha=0.5,
            label=f'Observed ({pair_info["date2"].year})')

    # Plot predicted centerline if available
    if prediction is not None:
        ax.plot(prediction['predicted_x'], prediction['predicted_y'],
                'r--', lw=1, label='Predicted position')

    ax.set_aspect('equal')
    ax.autoscale_view()

    # Colorbar
    sm = plt.cm.ScalarMappable(cmap=cmap,
                                norm=Normalize(vmin=-vmax, vmax=vmax))
    sm.set_array([])
    cb = plt.colorbar(sm, ax=ax, shrink=0.6, pad=0.02)
    cb.set_label('Observed - Predicted MR (m/yr)\n'
                 'red = underpredicted, blue = overpredicted')

    ax.legend()
    ax.set_xlabel('Easting (m)')
    ax.set_ylabel('Northing (m)')

    return ax

select_aoi_interactive(image, extent, title=None)

Display an image and let the user select a rectangular area of interest.

Parameters:

Name	Type	Description	Default
image	ndarray	Image to display (e.g., a normalized false-color image).	required
extent	list	[left, right, bottom, top] coordinates of the image.	required
title	str	Title for the plot (e.g., the scene date).	None

Returns:

Type	Description
list	[left, right, bottom, top] AOI bounds in the image's coordinates.

Source code in rivabar/visualization.py

def select_aoi_interactive(image, extent, title=None):
    """
    Display an image and let the user select a rectangular area of interest.

    Parameters
    ----------
    image : numpy.ndarray
        Image to display (e.g., a normalized false-color image).
    extent : list
        [left, right, bottom, top] coordinates of the image.
    title : str, optional
        Title for the plot (e.g., the scene date).

    Returns
    -------
    list
        [left, right, bottom, top] AOI bounds in the image's coordinates.
    """
    from matplotlib.patches import Rectangle

    print("Interactive AOI selection: click two points to define a rectangle "
          "(two opposite corners)")

    fig, ax = plt.subplots(figsize=(12, 8))
    ax.imshow(image, extent=extent, origin='upper')
    ax.set_xlabel('UTM Easting (m)')
    ax.set_ylabel('UTM Northing (m)')
    ax.set_title((title + '\n' if title else '') +
                 'Click two points: two opposite corners of the AOI')
    ax.grid(True, alpha=0.3)
    plt.tight_layout()
    plt.show(block=False)

    points = plt.ginput(2, timeout=0)
    if len(points) != 2:
        plt.close(fig)
        raise ValueError("Need exactly 2 points to define the AOI")

    (x1, y1), (x2, y2) = points
    left, right = min(x1, x2), max(x1, x2)
    bottom, top = min(y1, y2), max(y1, y2)

    width, height = right - left, top - bottom
    rect = Rectangle((left, bottom), width, height,
                     linewidth=2, edgecolor='red', facecolor='none')
    ax.add_patch(rect)
    ax.set_title(f'Selected AOI: {width:.0f} m x {height:.0f} m')
    plt.draw()
    print(f"AOI selected: {width:.0f} m x {height:.0f} m, "
          f"bounds: ({left:.0f}, {bottom:.0f}) to ({right:.0f}, {top:.0f})")
    return [left, right, bottom, top]

make_river_evolution_frames(matched_rivers, image_stack, dates, aoi_bounds, output_dir='river_evolution_figures', select_dates=None, baseline_index=0, depo_cmap='plasma', erosion_cmap='plasma', alpha=1.0, dpi=300, scalebar=True, date_fontsize=14, frame_prefix='river_evolution')

Render animation frames showing river evolution through time.

Each frame has three panels: the (false-color) image, the image with cumulative deposition polygons, and the image with cumulative erosion polygons (colored by age via :func:create_and_plot_bars). Frames where deposition/erosion cannot be computed (the first frame, and any frame fewer than two time steps after baseline_index) show the plain channel polygon instead.

The three input lists must be aligned (same length, same order), as returned by :func:match_rivers_to_images and :func:prepare_image_stack.

Parameters:

Name	Type	Description	Default
matched_rivers	list of River	Processed rivers, in date order.	required
image_stack	list of numpy.ndarray or None	Cropped, normalized background images (None entries are skipped).	required
dates	list of datetime.datetime	Acquisition dates, used for frame labels and filenames.	required
aoi_bounds	list	[left, right, bottom, top] of the area of interest.	required
output_dir	str	Directory for the output PNGs (created if needed).	'river_evolution_figures'
select_dates	list of str or None	If given, only rivers whose acquisition_date ('YYYY-MM-DD') is in this list are rendered. None (default) renders all.	None
baseline_index	int	Index of the river used as time step 1 for the cumulative deposition/erosion maps (default 0).	0
depo_cmap	str	Colormaps for deposition and erosion polygons (default 'plasma').	'plasma'
erosion_cmap	str	Colormaps for deposition and erosion polygons (default 'plasma').	'plasma'
alpha	float	Transparency of the deposition/erosion polygons (default 1.0).	1.0
dpi	int	Resolution of the saved figures (default 300).	300
scalebar	bool	Whether to add scale bars (requires the optional matplotlib-scalebar package; silently skipped with a warning if not installed).	True
date_fontsize	int	Font size of the date label (default 14).	14
frame_prefix	str	Filename prefix for the frames (default 'river_evolution').	'river_evolution'

Returns:

Type	Description
list of str	Paths of the PNG frames that were written, in order.

Source code in rivabar/visualization.py

def make_river_evolution_frames(matched_rivers, image_stack, dates, aoi_bounds,
                                output_dir='river_evolution_figures',
                                select_dates=None, baseline_index=0,
                                depo_cmap='plasma', erosion_cmap='plasma',
                                alpha=1.0, dpi=300, scalebar=True,
                                date_fontsize=14, frame_prefix='river_evolution'):
    """
    Render animation frames showing river evolution through time.

    Each frame has three panels: the (false-color) image, the image with
    cumulative deposition polygons, and the image with cumulative erosion
    polygons (colored by age via :func:`create_and_plot_bars`). Frames where
    deposition/erosion cannot be computed (the first frame, and any frame
    fewer than two time steps after *baseline_index*) show the plain channel
    polygon instead.

    The three input lists must be aligned (same length, same order), as
    returned by :func:`match_rivers_to_images` and :func:`prepare_image_stack`.

    Parameters
    ----------
    matched_rivers : list of River
        Processed rivers, in date order.
    image_stack : list of numpy.ndarray or None
        Cropped, normalized background images (None entries are skipped).
    dates : list of datetime.datetime
        Acquisition dates, used for frame labels and filenames.
    aoi_bounds : list
        [left, right, bottom, top] of the area of interest.
    output_dir : str, optional
        Directory for the output PNGs (created if needed).
    select_dates : list of str or None, optional
        If given, only rivers whose ``acquisition_date`` ('YYYY-MM-DD') is in
        this list are rendered. None (default) renders all.
    baseline_index : int, optional
        Index of the river used as time step 1 for the cumulative
        deposition/erosion maps (default 0).
    depo_cmap, erosion_cmap : str, optional
        Colormaps for deposition and erosion polygons (default 'plasma').
    alpha : float, optional
        Transparency of the deposition/erosion polygons (default 1.0).
    dpi : int, optional
        Resolution of the saved figures (default 300).
    scalebar : bool, optional
        Whether to add scale bars (requires the optional matplotlib-scalebar
        package; silently skipped with a warning if not installed).
    date_fontsize : int, optional
        Font size of the date label (default 14).
    frame_prefix : str, optional
        Filename prefix for the frames (default 'river_evolution').

    Returns
    -------
    list of str
        Paths of the PNG frames that were written, in order.
    """
    import os
    import matplotlib.gridspec as gridspec
    from .temporal_analysis import create_and_plot_bars

    if not (len(matched_rivers) == len(image_stack) == len(dates)):
        raise ValueError(
            f"matched_rivers ({len(matched_rivers)}), image_stack "
            f"({len(image_stack)}), and dates ({len(dates)}) must have the "
            f"same length and order")

    ScaleBar = None
    if scalebar:
        try:
            from matplotlib_scalebar.scalebar import ScaleBar
        except ImportError:
            print("matplotlib-scalebar is not installed; skipping scale bars "
                  "(pip install matplotlib-scalebar)")

    os.makedirs(output_dir, exist_ok=True)
    written = []
    count = 0
    for i in range(len(matched_rivers)):
        if select_dates is not None and matched_rivers[i].acquisition_date not in select_dates:
            continue
        if image_stack[i] is None:
            print(f"Skipping {dates[i].strftime('%Y-%m-%d')}: no background image")
            continue

        fig_height = 1.6 * image_stack[i].shape[0] / 100
        fig_width = 5 * image_stack[i].shape[1] / 100
        fig = plt.figure(figsize=(fig_width, fig_height))
        gs = gridspec.GridSpec(1, 3, wspace=0.001, hspace=0.02)
        axes = [fig.add_subplot(gs[0, j]) for j in range(3)]
        for ax in axes:
            ax.imshow(image_stack[i], extent=aoi_bounds)
            ax.axis('off')

        axes[0].text(0.02, 0.015, dates[i].strftime('%Y-%m-%d'),
                     transform=axes[0].transAxes, fontsize=date_fontsize,
                     fontweight='bold', color='black',
                     bbox=dict(boxstyle='round,pad=0.3', facecolor='white', alpha=0.5))

        # Deposition/erosion needs at least two time steps after the baseline;
        # earlier frames (and failed computations) show the channel polygon
        plotted = False
        if i - baseline_index >= 2:
            chs, bars, erosions, aoi_dates, aoi_centerlines = create_and_plot_bars(
                matched_rivers, baseline_index, i,
                ax1=axes[1], ax2=axes[2],
                depo_cmap=depo_cmap, erosion_cmap=erosion_cmap,
                alpha=alpha, aoi=aoi_bounds, colorbar=False,
                color_scale_timestep=len(matched_rivers) - 1)
            plotted = chs is not None
        if not plotted:
            _plot_channel_polygon_on_axes(matched_rivers[i], aoi_bounds, axes[1:])

        if ScaleBar is not None:
            for ax in axes:
                ax.add_artist(ScaleBar(1, "m", length_fraction=0.25,
                                       location="lower right", border_pad=0.4,
                                       box_alpha=0.5, color='black'))

        output_filename = os.path.join(
            output_dir, f"{frame_prefix}_{count:03d}_{dates[i].strftime('%Y%m%d')}.png")
        fig.savefig(output_filename, dpi=dpi, bbox_inches='tight')
        plt.close(fig)
        written.append(output_filename)
        count += 1

    print(f"Wrote {len(written)} frames to {output_dir}")
    return written

assemble_movie(frames, output_path, fps=6, width=None, crf=18)

Assemble PNG frames into an MP4 movie using ffmpeg.

Parameters:

Name	Type	Description	Default
frames	str or list	Either a glob pattern (e.g. 'river_evolution_figures/river_evolution_*.png'), a directory containing the PNG frames, or an explicit list of frame paths (e.g., the return value of :func:make_river_evolution_frames, in which case their common directory and prefix are used).	required
output_path	str	Path of the output .mp4 file.	required
fps	int	Frames per second (default 6).	6
width	int	Output width in pixels; the height is scaled proportionally. If None (default), the full resolution is kept (trimmed to even dimensions, as required by the H.264 encoder).	None
crf	int	H.264 quality (lower = better/larger; default 18).	18

Returns:

Type	Description
str	The output path.

Source code in rivabar/visualization.py

def assemble_movie(frames, output_path, fps=6, width=None, crf=18):
    """
    Assemble PNG frames into an MP4 movie using ffmpeg.

    Parameters
    ----------
    frames : str or list
        Either a glob pattern (e.g. ``'river_evolution_figures/river_evolution_*.png'``),
        a directory containing the PNG frames, or an explicit list of frame
        paths (e.g., the return value of :func:`make_river_evolution_frames`,
        in which case their common directory and prefix are used).
    output_path : str
        Path of the output .mp4 file.
    fps : int, optional
        Frames per second (default 6).
    width : int, optional
        Output width in pixels; the height is scaled proportionally. If None
        (default), the full resolution is kept (trimmed to even dimensions,
        as required by the H.264 encoder).
    crf : int, optional
        H.264 quality (lower = better/larger; default 18).

    Returns
    -------
    str
        The output path.
    """
    import glob as glob_module
    import os
    import shutil
    import subprocess

    ffmpeg = shutil.which('ffmpeg')
    if ffmpeg is None:
        raise RuntimeError("ffmpeg not found on PATH; install it (e.g. "
                           "'brew install ffmpeg' or 'conda install ffmpeg') "
                           "to assemble movies")

    if isinstance(frames, (list, tuple)):
        if not frames:
            raise ValueError("No frames provided")
        dirname = os.path.dirname(frames[0])
        prefix = os.path.basename(frames[0]).rsplit('_', 2)[0]
        pattern = os.path.join(dirname, f'{prefix}_*.png')
    elif os.path.isdir(frames):
        pattern = os.path.join(frames, '*.png')
    else:
        pattern = frames
    matched = sorted(glob_module.glob(pattern))
    if not matched:
        raise ValueError(f"No frames match {pattern}")

    if width is not None:
        vf = f'scale={width}:-2'
    else:
        # H.264 with yuv420p requires even dimensions
        vf = 'crop=trunc(iw/2)*2:trunc(ih/2)*2'

    cmd = [ffmpeg, '-y', '-framerate', str(fps), '-pattern_type', 'glob',
           '-i', pattern, '-vf', vf, '-c:v', 'libx264', '-crf', str(crf),
           '-pix_fmt', 'yuv420p', '-movflags', '+faststart', output_path]
    result = subprocess.run(cmd, capture_output=True, text=True)
    if result.returncode != 0:
        raise RuntimeError(f"ffmpeg failed:\n{result.stderr[-2000:]}")
    print(f"Wrote {output_path} ({len(matched)} frames at {fps} fps, "
          f"{len(matched)/fps:.1f} s)")
    return output_path