3

如何找到使用 scipy 拟合的样条曲线的峰值曲率？（实际上，峰值二阶微分就足够了）

我使用我的 1d和向量计算了tck如下值：xsys

tck = splrep(xs, ys, s=0)

我知道我可以x根据我的任何选择评估二阶微分：

ddy = splev([x], tck, 2)

所以我可以遍历 的许多值x，计算曲率并取最大值。但我更愿意解释这些值tck以获得各个三次函数的系数，从而直接计算峰值曲率。但是，tck看起来相当不透明 - 我如何从中提取三次函数系数？

4

1 回答 1

0

只需der在函数上使用关键字参数splev：

ddy = splev(X, tck, der=2)

并且最好不要遍历 的许多值x，而是创建一个包含您要评估的每个值的 Nx1 数组X，以便返回一个值数组，而不是您必须放入一个序列中的单个值。

此外，非常建议您绘制结果作为调试它的一种方式。如果情节有意义，那么事情很可能会像您期望的那样工作（如果没有，它们肯定不会工作）。

编辑：如果使用 X 的插值仅给出一个近似值并且您想要 TRUE 最大值，则可以使用定义最大值（局部插值最大值及其邻居）的三个点的抛物线插值，考虑到样条是局部平滑的：

def parabolic_interpolation(p1, p2, p3):
    x1, y1 = p1
    x2, y2 = p2
    x3, y3 = p3

    denom = (x1-x2)*(x1-x3)*(x2-x3);
    a = (x3*(y2-y1)+x2*(y1-y3)+x1*(y3-y2))/denom
    b = (x3*x3*(y1-y2)+x2*x2*(y3-y1)+x1*x1*(y2-y3))/denom
    c = (x2*x3*(x2-x3)*y1+x3*x1*(x3-x1)*y2+x1*x2*(x1-x2)*y3)/denom

    xv = -b/(2*a)
    yv = c-b**2/(4*a)

    return (xv, yv)  # coordinates of the vertex

希望这可以帮助！

于 2012-11-27T18:58:55.653 回答