平滑算法:三次样条插值(Cubic Spline Interpolation)

https://blog.csdn.net/left_la/article/details/6347373

感谢强大的google翻译。

我从中认识到了航位推算dead reckoning，立方体样条Cubic Splines 算法。

我单独查找了 Cubic Splines ，里面的原理简单说明：

Cubic Splines 认为在 x 在[a, b]区间中，y对应是一条平滑的曲线，所以 y = f(x); 的一阶导函数和二阶导函数是平滑连续可导的。

拟定用三次方程，所以得出了一般的三次方程和一阶导数方程和二阶导数方程。

然后求各个分部的解。

这是三次样条的基本原理。

但文中最开始的链接中所得出的

x = At³ + Bt² + Ct + D
y = Et³ + Ft³ + Gt + H

t是percent（0~1）区间值，如果还有三维向量，我理解是同样的展开。

然后通过四个位置点来求出 A B C D … 各分部参数的值

A = x₃ – 3x₂ +3x₁ – x₀
B = 3x₂ – 6x₁ + 3x₀
C = 3x₁ – 3x₀
D = x₀

E = y₃ – 3y₂ +3y₁ – y₀
F = 3y₂ – 6y₁ + 3y₀
G = 3y₁ – 3y₀
H = y₀

…

相同分量展开。（如果有Z 分量的话）

学艺不精，无法从现有姿势推出这个分量求解过程。

实时运动游戏是通过预测其他玩家的位置来表现的，当服务器有新的输入的时候，本地玩家会发现其他玩家位置或状态发生一次跳变（瞬移）。

有两种思路，

一、预测未来

1. 通过当前位置和速度，通过预测未来精度（1s或者0.5s）推测出未来位置.
2. 得出公式参数，通过dt来平滑当前运动轨迹。

二、延迟渲染

1. 通过延迟渲染参数（延迟1s，0.5s来）来获得其他玩家的过去状态位置。
2. 得出公式参数，通过dt来平滑运动轨迹。

上述两种方案

1. 如果参数一致，速度不改，则运动轨迹跟预测一致，如果玩家输入多变，则永远不会是真实的位置。
2. 看到的玩家的过去位置，移动轨迹跟目标玩家运动轨迹基本保持一致。

https://gist.github.com/svdamani/1015c5c4b673c3297309#file-spline-c-L26

/** Numerical Analysis 9th ed - Burden, Faires (Ch. 3 Natural Cubic Spline, Pg. 149) */
#include <stdio.h>

int main() {
    /** Step 0 */
    int n, i, j;
    scanf("%d", &n);
    n--;
    float x[n + 1], a[n + 1], h[n], A[n], l[n + 1],
        u[n + 1], z[n + 1], c[n + 1], b[n], d[n];
    for (i = 0; i < n + 1; ++i) scanf("%f", &x[i]);
    for (i = 0; i < n + 1; ++i) scanf("%f", &a[i]);

    /** Step 1 */
    for (i = 0; i <= n - 1; ++i) h[i] = x[i + 1] - x[i];

    /** Step 2 */
    for (i = 1; i <= n - 1; ++i)
        A[i] = 3 * (a[i + 1] - a[i]) / h[i] - 3 * (a[i] - a[i - 1]) / h[i - 1];

    /** Step 3 */
    l[0] = 1;
    u[0] = 0;
    z[0] = 0;

    /** Step 4 */
    for (i = 1; i <= n - 1; ++i) {
        l[i] = 2 * (x[i + 1] - x[i - 1]) - h[i - 1] * u[i - 1];
        u[i] = h[i] / l[i];
        z[i] = (A[i] - h[i - 1] * z[i - 1]) / l[i];
    }

    /** Step 5 */
    l[n] = 1;
    z[n] = 0;
    c[n] = 0;

    /** Step 6 */
    for (j = n - 1; j >= 0; --j) {
        c[j] = z[j] - u[j] * c[j + 1];
        b[j] = (a[j + 1] - a[j]) / h[j] - h[j] * (c[j + 1] + 2 * c[j]) / 3;
        d[j] = (c[j + 1] - c[j]) / (3 * h[j]);
    }

    /** Step 7 */
    printf("%2s %8s %8s %8s %8s
", "i", "ai", "bi", "ci", "di");
    for (i = 0; i < n; ++i)
        printf("%2d %8.2f %8.2f %8.2f %8.2f
", i, a[i], b[i], c[i], d[i]);
    return 0;
}

这个上面根据 https://fac.ksu.edu.sa/sites/default/files/numerical_analysis_9th.pdf#page=167

实现了对应 x 求 y 的函数，这里x可以替换成 时间t，分别求 t 跟x 、y、z的abcd参数，最终求出s(t)函数。

INPUT

n; x0, x1, ... , xn;

a0 = f (x0), a1 = f (x1), ... , an = f (xn).

OUTPUT aj, bj, cj, dj for j = 0, 1, ... , n − 1.

(Note: S(x) = Sj(x) = aj + bj(x − xj) + cj(x − xj)2 + dj(x − xj)3 for xj ≤ x ≤ xj+1.).

最后通过x在哪个区间调用某个区间的 S(x) 函数。注意S(x)函数是一组函数，x多区间。

（或许上面两个文章介绍的其实不是一种算法 0.0）