公式中的矩阵操作如何应用于CocosCreator

前言

在使用CocosCreator开发的过程中，少不了进行一些数学运算。有一些是通用的图形学公式，去各种地方百度到的，但以一种数学论述的方式存在。猛然碰见，不是太好实现在程序之中。下载就以3次B样条曲线的公式为例来说明如何将矩阵操作“翻译成”CocosCreator中的JavaScript/TypeScript代码。

需求

这次我们发现需要画出3次B样条曲线，也就是输入n个点，来确定一条空间曲线路径，比如在游戏中需要飞行轨迹，弹道轨迹之类的需求。一般情况下，我们都去搜索引擎里面找相关的数学框架，其实就是数学结论。这次，我就找到了这些东西：
B样条：$P(u)=\sum _{i=0}^n{P}_i{B}_{i,k}(u)u\in [u_{k-1},u_{n+1}]$
$$\begin{gathered} B_{0,3}(t)=\frac{1}{6}(1-t{)}^3 \\ {B}_{1,3}(t)=\frac{1}{6}(3t^3-6t^2+4) \\ {B}_{2,3}(t)=\frac{1}{6}(-3t^3+3t^2+3t+1) \\ {B}_{3,3}(t)=\frac{1}{6}{t}^3 \end{gathered}$$
以及它的矩阵形式
$Q(t)=\sum _{i=0}^3{P}_i{B}_{i,3}(t)=\frac{1}{6}\left[\begin{array}{cccc}1& t& t^2& t^3\end{array}\right]\cdot \left[\begin{array}{cccc}1& 4& 1& 0\\ -3& 0& 3& 0\\ 3& -6& 3& 0\\ -1& 3& -3& 1\end{array}\right]\cdot \left[\begin{array}{c}{P}_i\\ P_{i+1}\\ P_{i+2}\\ P_{i+3}\end{array}\right]t\in [0,1]$

代数形式

其中代数形式不是本文想要解释的内容，这里随便实现了一下。其中P为cc.Vec3，N表示上面公式中的B函数组，循环就是连加，order是阶数，这里固定是3，addSelf是自身的每一个分量加上后面对应的分量，mul是向量每一个分量乘以某个值。

let P = cc.v3(0, 0, 0);
let N = [
    (t) => Math.pow(1 - t, 3) / 6,
    (t) => (3 * t * t * t - 6 * t * t + 4) / 6,
    (t) => (-3 * t * t * t + 3 * t * t + 3 * t + 1) / 6,
    (t) => t * t * t / 6
];
for (let i = 0; i < this.order; i++) {
    P.addSelf(this.pList[i + n].mul(N[i](t)));
}

矩阵形式

这里是文章的核心，直接用这里结论即可。
首先CocosCreator内置类型有相对完整的cc.Mat4，也有cc.Mat3，但是没有相关运算，不完整。我这里就都全部使用cc.Mat4作为矩阵的基本数据类型了。我发现的几点套路：

1. mat4当然是4阶矩阵，如果要算的矩阵行或列超过4阶，就不能用这个了。
2. 如果本身公式中的矩阵行列少于4个，那就空着，创建矩阵时自动保持单位矩阵中原来的值。一行都不够，那就补0。
   比如这个1x2矩阵$\left[\begin{array}{cccc}1& t& t^2& t^3\end{array}\right]$,表示为：

cc.mat4(1, t, t * t, t * t * t)

3. cc.mat4(…)创建的矩阵默认是个单位矩阵，其它矩阵与它相乘都不发生变化。
4. 点组成的列矩阵：$\left[\begin{array}{c}{P}_i\\ P_{i+1}\\ P_{i+2}\\ P_{i+3}\end{array}\right]$应该把每个点扩展为一个行比如：

cc.mat4(
     this.pList[n].x, this.pList[n].y, this.pList[n].z, 0,
     this.pList[n + 1].x, this.pList[n + 1].y, this.pList[n + 1].z, 0,
     this.pList[n + 2].x, this.pList[n + 2].y, this.pList[n + 2].z, 0,
     this.pList[n + 3].x, this.pList[n + 3].y, this.pList[n + 3].z, 0,
 );

5. 最后就是矩阵乘法，Mat4有个mul函数，参数是一个矩阵，实际意义就是原矩阵左乘参数中的矩阵，注意是左乘就对了，困惑了我好久。
6. 乘常量，Mat4里面的multiplyScalar。我发现还有个函数叫mulScalar，乘出来结果为null，不要用就好了。（我用的CocosCreator2.4.6）
7. 结果值，因为前面的矩阵补过0，最后的结果就只关注第一行，vec3的x,y,z就取结果矩阵的第一行前三个数即可。

公式的代码实现

$Q(t)=\sum _{i=0}^3{P}_i{B}_{i,3}(t)=\frac{1}{6}\left[\begin{array}{cccc}1& t& t^2& t^3\end{array}\right]\cdot \left[\begin{array}{cccc}1& 4& 1& 0\\ -3& 0& 3& 0\\ 3& -6& 3& 0\\ -1& 3& -3& 1\end{array}\right]\cdot \left[\begin{array}{c}{P}_i\\ P_{i+1}\\ P_{i+2}\\ P_{i+3}\end{array}\right]t\in [0,1]$

let entry1 = cc.mat4(1, t, t * t, t * t * t)
let entry2 = cc.mat4(1, 4, 1, 0, -3, 0, 3, 0, 3, -6, 3, 0, -1, 3, -3, 1);
let entry3 = cc.mat4(
    this.pList[n].x, this.pList[n].y, this.pList[n].z, 0,
    this.pList[n + 1].x, this.pList[n + 1].y, this.pList[n + 1].z, 0,
    this.pList[n + 2].x, this.pList[n + 2].y, this.pList[n + 2].z, 0,
    this.pList[n + 3].x, this.pList[n + 3].y, this.pList[n + 3].z, 0,
);
let result = entry3.mul(entry2).mul(entry1).multiplyScalar(1 / 6);
P = cc.v3(result.m[0], result.m[1], result.m[2]);

注意图中的n就是公式里面的i，由于需求给入的是多个连续点，n就在0~给定点数-阶数范围遍历取值。
由于mul表示左乘，那么就是用的 entry3.mul(entry2).mul(entry1)，而不是反过来，因为矩阵乘法不满足交换律。

观察

如果CocosCreator提供从 cc.vec3或者cc.vec4的数组直接创建cc.mat4就好了。代码里面最大的那一块就可以化简一些了。
至于算法复杂度，上面的公式法似乎已经把前两个矩阵乘完了。而矩阵法，第一个矩阵里面的$t^2,t^3$则已经算好存进了矩阵，且里面没有循环，顺序流程直接出结果。
当然，B样条曲线有专门的几何方法去生成。一般来说，计算机图形学中，几何方法总是比代数方法效率优化的多。B样条曲线的基础标准生成算方法就叫de Boor-cox。我最后把这个算法的代码贴在下面，这并不是写这篇文章的主要内容，可以感受一下简单明了的递归方法，而且还是任意阶次的。

class deBoorcox {

    n: number;
    pList: cc.Vec3[];
    k: number;
    u: number;

    constructor(pList: cc.Vec3[], k: number) {
        this.pList = pList;
        this.n = pList.length - 1;
        this.k = k;
        this.u = k - 1;
    }

    // u取值 [k-1, n+1)
    B(u, k, i) {
        if (k == 1) {
            if (i <= u && u < i + 1) return 1;
            else return 0;
        } else {
            let coef_0 = 0, coef_1 = 0;
            if (u - i == 0 && k - 1 == 0) {
            } else {
                coef_0 = (u - i) / (k - 1);
            }
            if (i + k - u == 0 && k - 1 == 0) {
            } else {
                coef_1 = (i + k - u) / (k - 1);
            }
            return coef_0 * this.B(u, k - 1, i) + coef_1 * this.B(u, k - 1, i + 1);
        }
    }

    solve(u) {
        let P = cc.v3(0, 0, 0);
        for (let i = 0; i <= this.n; i++) {
            P.addSelf(this.pList[i].mul(this.B(u, this.k, i)));
        }
        return P;
    }
}