def METHOD_NAME(
    gdf_stream: gpd.GeoDataFrame,
    da_rivmask: xr.DataArray,
    nmin=5,
) -> np.ndarray:
    """Return average river width along a segment based on a river mask raster.
    For each segment in gdf_stream the associated area is calculated from stream mask
    and divided by the segment length to obtain the average width.
    Parameters
    ----------
    gdf_stream : gpd.GeoDataFrame
        River segments
    da_rivmask : xr.DataArray
        Boolean river mask in projected grid.
    nmin : int, optional
        Minimum number of cells in rivmask to calculate the width, by default 5
    Returns
    -------
    rivwth : np.ndarray
        Average width per segment in gdf_stream
    """
    assert da_rivmask.raster.crs.is_projected
    gdf_stream = gdf_stream.copy()
    # get/check river length
    if "rivlen" not in gdf_stream.columns:
        gdf_stream["rivlen"] = gdf_stream.to_crs(da_rivmask.raster.crs).length
    # rasterize streams
    gdf_stream["segid"] = np.arange(1, gdf_stream.index.size + 1, dtype=np.int32)
    segid = da_rivmask.raster.rasterize(gdf_stream, "segid").astype(np.int32)
    segid.raster.set_nodata(0)
    segid.name = "segid"
    # remove islands to get total width of braided rivers
    da_mask = da_rivmask.copy()
    da_mask.data = ndimage.binary_fill_holes(da_mask.values)
    # find nearest stream segment for all river cells
    segid_spread = spread2d(segid, da_mask)
    # get average width based on da_rivmask area and segment length
    cellarea = abs(np.multiply(*da_rivmask.raster.res))
    seg_count = ndimage.sum(
        da_rivmask, segid_spread["segid"].values, gdf_stream["segid"].values
    )
    rivwth = seg_count * cellarea / gdf_stream["rivlen"]
    valid = np.logical_and(gdf_stream["rivlen"] > 0, seg_count > nmin)
    return np.where(valid, rivwth, -9999)