如何让Scipy LSQUnivariateSpline闭合样条首尾点满足切线拟合？

使用Scipy Interpolate模块的LSQUnivariateSpline方法拟合2D点生成B样条（用于输出DXF）时，整体效果良好，但在单元测试中发现：对于首尾点处于拐角的正方形这类闭合形状，拟合出的样条首尾点未像其他拐角那样进行切线拟合优化，需要找到让首尾点也纳入切线最佳拟合的方法，同时保留现有抗圆角特性。应用场景为激光扫描无序点云的2D切片处理，需先完成点排序、清理再拟合。

def b_spline_interpolation(np_xy_points: np.array) -> np.array:
    """Creates a b-spline based on the scipy module 'interpolate' and it's method 'LSQ Univariate Spline'."""

    # interpolate based on the 'LSQ Univariate Spline'
    ui = np.cumsum(np.r_[[0], np.linalg.norm(np.diff(np_xy_points, axis=0), axis=1)])
    knots = np.linspace(ui[0], ui[-1], 40)
    try:
        sx = interpolate.LSQUnivariateSpline(ui, np_xy_points[:, 0], knots[1:-1], k=3)
        sy = interpolate.LSQUnivariateSpline(ui, np_xy_points[:, 1], knots[1:-1], k=3)
    except ValueError as err:
        log.error(f"Failed to create the LSQ Univariate Spline, error: {err}")
        return np.array([])

    # sampling the resulting spline
    uu = np.linspace(ui[0], ui[-1], 150)
    xx = sx(uu)
    yy = sy(uu)

    # create the numpy 2D array
    np_xy_array = np.vstack((xx, yy)).T

    return np_xy_array

针对闭合形状首尾切线拟合的问题，有以下几种可行方案：

1. 给LSQUnivariateSpline添加边界切线约束 #

LSQUnivariateSpline默认采用自然样条边界条件（首尾二阶导数为0），这会导致首尾拐角处切线无法贴合原始点的趋势。可以手动计算首尾的切线方向，通过bc_type参数指定一阶导数约束：

def b_spline_interpolation(np_xy_points: np.array) -> np.array:
    ui = np.cumsum(np.r_[[0], np.linalg.norm(np.diff(np_xy_points, axis=0), axis=1)])
    knots = np.linspace(ui[0], ui[-1], 40)
    
    # 计算首尾切线（用相邻两点的斜率近似）
    # 首点切线：取前两个点的方向
    start_dir_x = np_xy_points[1,0] - np_xy_points[0,0]
    start_dir_y = np_xy_points[1,1] - np_xy_points[0,1]
    # 尾点切线：取后两个点的方向
    end_dir_x = np_xy_points[-1,0] - np_xy_points[-2,0]
    end_dir_y = np_xy_points[-1,1] - np_xy_points[-2,1]
    
    # 归一化切线（可选，确保导数幅度合理）
    start_norm = np.linalg.norm([start_dir_x, start_dir_y])
    end_norm = np.linalg.norm([end_dir_x, end_dir_y])
    start_tan_x = start_dir_x / start_norm if start_norm !=0 else 0
    start_tan_y = start_dir_y / start_norm if start_norm !=0 else 0
    end_tan_x = end_dir_x / end_norm if end_norm !=0 else 0
    end_tan_y = end_dir_y / end_norm if end_norm !=0 else 0
    
    try:
        # 给x、y分量分别设置一阶导数边界条件
        sx = interpolate.LSQUnivariateSpline(ui, np_xy_points[:, 0], knots[1:-1], k=3,
                                            bc_type=((1, start_tan_x), (1, end_tan_x)))
        sy = interpolate.LSQUnivariateSpline(ui, np_xy_points[:, 1], knots[1:-1], k=3,
                                            bc_type=((1, start_tan_y), (1, end_tan_y)))
    except ValueError as err:
        log.error(f"Failed to create the LSQ Univariate Spline, error: {err}")
        return np.array([])

    uu = np.linspace(ui[0], ui[-1], 150)
    xx = sx(uu)
    yy = sy(uu)
    np_xy_array = np.vstack((xx, yy)).T

    return np_xy_array

2. 扩展点云模拟闭合连续 #

对于闭合形状，将点云的前几个点复制到末尾（或末尾点复制到开头），让拟合算法在首尾衔接处有足够的点支撑切线拟合，最后采样时截断到原始点范围即可：

def b_spline_interpolation(np_xy_points: np.array) -> np.array:
    # 扩展点云：复制前3个点到末尾，模拟闭合连续
    extended_points = np.vstack([np_xy_points, np_xy_points[:3]])
    
    ui = np.cumsum(np.r_[[0], np.linalg.norm(np.diff(extended_points, axis=0), axis=1)])
    knots = np.linspace(ui[0], ui[-1], 40)
    try:
        sx = interpolate.LSQUnivariateSpline(ui, extended_points[:, 0], knots[1:-1], k=3)
        sy = interpolate.LSQUnivariateSpline(ui, extended_points[:, 1], knots[1:-1], k=3)
    except ValueError as err:
        log.error(f"Failed to create the LSQ Univariate Spline, error: {err}")
        return np.array([])

    # 采样时只取原始点对应的参数范围
    original_ui_max = np.cumsum(np.r_[[0], np.linalg.norm(np.diff(np_xy_points, axis=0), axis=1)])[-1]
    uu = np.linspace(ui[0], original_ui_max, 150)
    xx = sx(uu)
    yy = sy(uu)
    np_xy_array = np.vstack((xx, yy)).T

    return np_xy_array

3. 改用支持周期性的样条拟合工具 #

Scipy的splprep函数专门支持多变量样条拟合，且通过per=True参数可以开启周期性边界条件，自动处理首尾导数连续，更适合闭合曲线：

def closed_b_spline_interpolation(np_xy_points: np.array) -> np.array:
    """使用splprep处理闭合曲线样条拟合"""
    try:
        # per=True 开启周期性边界，s参数控制平滑度（0为过拟合，可根据需求调整）
        tck, u = interpolate.splprep([np_xy_points[:,0], np_xy_points[:,1]], s=0.5, per=True, k=3)
        uu = np.linspace(u.min(), u.max(), 150)
        xx, yy = interpolate.splev(uu, tck)
        np_xy_array = np.vstack((xx, yy)).T
        return np_xy_array
    except ValueError as err:
        log.error(f"Failed to create closed spline, error: {err}")
        return np.array([])

该方法无需手动处理边界条件，自动保证首尾切线连续，同时s参数可以灵活调整平滑度，保留抗圆角特性。

内容的提问来源于stack exchange，提问作者LVX