Scilab中的样条系数

Question

我想将分段三次样条曲线拟合到大量数据。我不希望使用B样条曲线，因为我希望这些样条曲线准确地通过数据点。我可以在Scilab中看到的唯一方法是使用splin和interp。

但是，我想确定结点之间的每个样条曲线的系数（因为我需要获取这些系数，并将其放入其他软件中）。所有splin给您的都是导数。有没有办法得到三次样条系数？还是有一种方法可以轻松地从一阶导数生成系数？

Answer 1

是，可以从您拥有的y值和splin返回的导数中获得系数。每个间隔[x（i），x（i + 1）]一个，您需要找到4个系数和4个等式：两端的值，两端的导数。最直接的方法就是告诉Scilab为每个子间隔解决这个4 x 4系统：花的时间不应该超过样条曲线本身的评估时间。下面的程序执行此操作。

x = [0,1,2,3,4,5]   // x values
y = [1,0,1,0,1,0]   // y values
d = splin(x,y)      
n = length(x)-1     // number of subintervals
cfs = zeros(4,n)    // matrix to store coefficients in
for i=1:n
    a = x(i)
    b = x(i+1)
    cfs(:,i) = [1,a,a^2,a^3; 1,b,b^2,b^3; 0,1,2*a,3*a^2; 0,1,2*b,3*b^2] \ [y(i);y(i+1);d(i);d(i+1)]
end

[前两个方程1,a,a^2,a^3; 1,b,b^2,b^3将多项式的值与y值相关；其他两个0,1,2*a,3*a^2; 0,1,2*b,3*b^2对于其导数也执行相同的操作。 （公式只是x的幂的导数。）

以上脚本的输出：

    1.     1.   - 8.6    13.     13.   
  - 3.4  - 3.4    11.  - 10.6  - 10.6  
    3.1    3.1  - 4.1    3.1     3.1   
  - 0.7  - 0.7    0.5  - 0.3   - 0.3   

每列具有四个系数：例如，样条的第一部分为1-3.4x+3.1x^2-0.7x^3。由于这是不打结的样条曲线，是splin命令的默认模式，因此第二部分与第一部分相同。并且倒数第二个与倒数第二个相同。

您可以通过绘制图块来检查其是否正常工作：

for i=1:n  
    t = linspace(x(i),x(i+1))
    plot(t,cfs(:,i)'*[ones(t); t; t.^2; t.^3])
end

---

话虽如此，以表格形式表示构成样条的多项式会更容易

 p(x) = y(i) + A*(x-x(i)) + B*(x-x(i))*(x-x(i+1)) + C*(x-x(i))^2*(x-x(i+1))

系数很容易在不求解线性系统的情况下一个接一个地找到：

• [A = (y(i+1)-y(i))/(x(i+1)-x(i))通过等于x（i + 1）处的值
• [B = (d(i)-A)/(x(i)-x(i+1))，通过在x（i）处等于导数]
• [C = (d(i+1)-A-B*(x(i+1)-x(i)))/(x(i+1)-x(i))^2，通过在x（i + 1）处等于导数

当然，这些系数应与上述适当的多项式一起使用。这是此替代版本

for i=1:n
    A = (y(i+1)-y(i))/(x(i+1)-x(i))
    B = (d(i)-A)/(x(i)-x(i+1))
    C = (d(i+1)-A-B*(x(i+1)-x(i)))/(x(i+1)-x(i))^2
    cfs(:,i) = [y(i);A;B;C]
end
// Again, plot for testing
for i=1:n
    t = linspace(x(i),x(i+1))
    plot(t,cfs(:,i)'*[ones(t); t-x(i); (t-x(i)).*(t-x(i+1)); ((t-x(i)).^2).*(t-x(i+1))])
end

Answer 2

我已经在scilab中编译了上述程序，输出与您的不一样，请告诉我我做错了什么。-> cfscfs =

1. 1. -20.333333 51.666667 51.666667

   -4.2222222 -4.2222222 27.777778 -44.222222 -44.222222

   4.3333333 4.3333333 -11.666667 12.333333 12.333333

   -1.1111111 -1.1111111 1.5555556 -1.1111111 -1.1111111