射线与空间内三角形的相交检测算法（Möller-Trumbore）的推导与实践

背景介绍(学习算法之前需要先了解)

射线与空间内三角形的相交检测是游戏程序设计中一个常见的问题，最典型的应用就是拾取(Picking)，本文介绍一个最常见的方法，这个方法也是DirectX中采用的方法，该方法速度快，而且存储空间少。先讲述理论，然后文章末尾给出对应的代码实现与Unity中的显示。
简单而直观的方法是：先判断射线是否与三角形所在的平面相交，如果相交，再判断交点是否在三角形内。但这种方法效率并不高，因为多计算了三角形所在的平面。
Möller-Trumbore射线三角相交算法是一种快速计算射线与空间内三角形的交点的方法，通过向量与矩阵计算可以快速得出交点与重心坐标，而无需对包含三角形的平面方程进行预计算。它应用于计算机图形学中以实现涉及三角形网格的光线跟踪计算。算法名字是以发明者TomasMöller和Ben Trumbore的名字来命名的。

参数说明

设射线的方程Ray=O+td(O为起点，D为射线方向(方向向量)，t为权重)，一个点从起点O开始，沿着方向D移动任意长度，得到终点R，根据t值得不同，得到得R值也不同，所有这些不同的R值便构成了整条射线。
三角形三个顶点$P0$，$P1$，$P2$。$u$为$P_1$的权重，$v$为$P_2$的权重，而$1-u-v$是$P_0$的权重，可以理解为沿着边AC移动一段距离，然后再沿着边AB移动一段距离，最后求它们的和向量。至于移动多大距离，就是由参数u和v控制的，所表达的数学意义是三角形及其内部所有点的方程。

$$Ray:O+td{\quad}(t≥0)$$

$$(1-u-v)P0+uP1+vP2$$

推导过程

两个重要的定理

克莱姆法则

解线性方程组时

$$\left[ \begin{matrix} -D & E_1 & E_2\end{matrix} \right]\left[ \begin{matrix} t \\ u \\ v\end{matrix} \right]=S$$

D，E1，E2是含有三个参数的行列式，即可表达成

$$t=\frac{det{\left[ \begin{matrix} S & E_1 & E_2\end{matrix} \right]}}{det{\left[ \begin{matrix} -D & E_1 & E_2\end{matrix} \right]}}=\frac{\left[ \begin{matrix} S_x & E_{x1} & E_{x2} \\ S_y & E_{y1} & E_{y2} \\ S_{z} & E_{z1} & E_{z2}\end{matrix} \right]}{\left[ \begin{matrix} -D_x & E_{x1} & E_{x2} \\ -D_y & E_{y1} & E_{y2} \\ -D_{z} & E_{z1} & E_{z2}\end{matrix} \right]}$$

$$u=\frac{det{\left[ \begin{matrix} -D & S & E_2\end{matrix} \right]}}{det{\left[ \begin{matrix} -D & E_1 & E_2\end{matrix} \right]}}=\frac{\left[ \begin{matrix} -D_x & S_x & E_{x2} \\ -D_y & S_y & E_{y2} \\ -D_{z} & S_z & E_{z2}\end{matrix} \right]}{\left[ \begin{matrix} -D_x & E_{x1} & E_{x2} \\ -D_y & E_{y1} & E_{y2} \\ -D_{z} & E_{z1} & E_{z2}\end{matrix} \right]}$$

$$v=\frac{det{\left[ \begin{matrix} -D & E_1 & S\end{matrix} \right]}}{det{\left[ \begin{matrix} -D & E_1 & E_2\end{matrix} \right]}}=\frac{\left[ \begin{matrix} -D_x & E_{x1} & S_{x} \\ -D_y & E_{y1} & S_y \\ -D_{z} & E_{z1} & S_z\end{matrix} \right]}{\left[ \begin{matrix} -D_x & E_{x1} & E_{x2} \\ -D_y & E_{y1} & E_{y2} \\ -D_{z} & E_{z1} & E_{z2}\end{matrix} \right]}$$

向量混合积

$$a·(b×c)=b·(c×a)=c·(a×b)\\ a·(b×c)=-a·(c×b)\\ a·(b×c)=-b·(a×c)\\ a·(b×c)=-c·(b×a)$$

将方程联立

$$O+td=(1-u-v)P0+P1+P2{\quad}(u≥0,v≥0,u+v≤1)$$

化简

$$O-P0=u(P1-P0)+v(P2-P0)-td\\ S=uE_1+vE_2-td\\{\quad}(S=O-P0,E_1=P1-P0,E_2=P2-P0)$$

$$\left[ \begin{matrix} -D & E_1 & E_2\end{matrix} \right]\left[ \begin{matrix} t \\ u \\ v\end{matrix} \right]=S$$

克拉姆法则

$$t=\frac{det{\left[ \begin{matrix} S & E_1 & E_2\end{matrix} \right]}}{det{\left[ \begin{matrix} -D & E_1 & E_2\end{matrix} \right]}}$$

向量混合积
分母部分
运用向量混合积定理，D前的-号，被抵消了

$$det{\left[ \begin{matrix} -D & E_1 & E_2\end{matrix} \right]}=-D·(E_1×E_2)=E_1·(D×E_2)$$

令

$$S1=D×E_2$$

分子部分：

$$det{\left[ \begin{matrix} S & E_1 & E_2\end{matrix} \right]}=((S×E_1)·E_2)$$

令

$$S2=S×E_1$$

原式等于

$$det{\left[ \begin{matrix} S & E_1 & E_2\end{matrix} \right]}=S2·E_2$$

因此

$$t=\frac{S2·E_2}{E_1·S1}$$

同理可推得其他两个参数u，v

$$u=\frac{S1·S}{E_1·S1}$$

$$v=\frac{S2·D}{E_1·S1}$$

总结

我们在最后可以通过已知的数据，求出t、u、v三个参数，通过三个参数的范围限制判断是否相交，如果不满足则不相交，如果满足则相交

代码实现

    // Vector3 a b c triangle vertexs
    // orig is ray original point, dir is direction vector
    bool rayTriangleIntersect(Vector3 orig, Vector3 dir,
        Vector3 a, Vector3 b, Vector3 c, float t, float b1, float b2)
    {
        bool isIn = false;
        Vector3 E1 = b - a;
        Vector3 E2 = c - a;
        Vector3 S = orig - a;
        Vector3 S1 = Vector3.Cross(dir, E2);
        Vector3 S2 = Vector3.Cross(S, E1);

        // 共同系数
        float coeff = 1.0f / Vector3.Dot(S1, E1);
        t = coeff * Vector3.Dot(S2, E2);
        b1 = coeff * Vector3.Dot(S1, S);
        b2 = coeff * Vector3.Dot(S2, dir);

        Debug.Log($"t = {t}, b1 = {b1}, b2 = {b2}");

        if (t >= 0 && b1 >= 0 && b2 >= 0 && (1 - b1 - b2) >= 0)
        {
            isIn = true;
        }

        return isIn;
    }

Unity中的演示效果

当射线与三角形相交时，三角形的材质会变成红色，射线是采用LineRenderer的形式，三角形用了编辑材质顶点的脚本，过几天会分享出来~
项目地址:https://github.com/shadow-lr/RayTriangleIntersect

题外话

虽然射线和三角形的相交检测可以用来实现拾取(Picking)，但是大多数程序并不采用这个方法，原因是这个方法效率很低，我们可以设想，一个大型的3D游戏，某个模型的三角形数量很可能是百万级的，在此情况下，对模型上的每个三角形求交是一件极其耗费时间的事情。
所以一般可行的方法是，用包围球和包围盒（AABB、OBB、FDH）来代替，计算出能容纳模型的最小球体或者举行提，只要判断射线与包围球或者包围盒求交即可，只是精确度上有一定误差，但是足以满足多数程序的需要。