高斯泼溅
高斯泼溅原理与实现
概览
高斯泼溅是基于 Splatting 和机器学习的三维重建方法。根据多视角图像，从 SfM 生成的点云和相机参数矩阵出发，将点膨胀成带有位置、形状、颜色和不透明度参数的高斯椭球，并以此为基本单位来表示场景。为了实现对三维场景的高质量表达，高斯泼溅通过深度学习的方式优化高斯椭球的参数。具体来说，是将用高斯椭球表示的三维场景按照给定的视角投影到二维平面得到渲染图像，渲染图像与真值图像作 loss，反向传播到各个可学习参数上完成优化。整个流程归纳为下图，本文也会围绕此图，记录其中的原理和实现细节。(注，本文提到的高斯泼溅，均指在 SIGGRAPH 2023 上发表的原版 3DGS，以下统称高斯泼溅)
SfM
从多视角图像生成三维点云的方法主要分为 SfM + MVS 的几何方法和基于 VGGT 架构的 transformer 方法，高斯泼溅用到的是几何方法中 SfM 生成的稀疏点云和相机参数。SfM(Structure from Motion，运动恢复结构)是一种从多视角图像中恢复三维结构和相机运动的计算机视觉技术，通过特征点匹配和对极约束的纯几何方法，估计每张图像对应的相机位姿和各个特征点在三维中的位置，并进行全局的多视图几何优化，最终得到稀疏点云和相机的内外参。
SfM并不是高斯泼溅中研究的内容，高斯泼溅中直接使用COLMAP这个SfM工具，调用COLMAP的命令完成特征提取、特征匹配、稀疏重建和图像去畸变四个步骤，获得场景的稀疏点云和相机内外参。对COLMAP工具提供的读取函数进行重构，并根据后续算法需要，将读取出的相关数据加工组织如下：
- Intrinsics：两轴焦距、两轴主点偏移、图像的宽度和高度
- Extrinsics：旋转四元数、平移向量
- Viewpoints：两轴向的视场角、世界坐标系到相机坐标系的转换矩阵（numpy_array和torch_tensor）、投影矩阵（torch_tensor）、变换矩阵（torch_tensor）、相机位置（torch_tensor）
- ImageInfos：图像文件名、图像宽度和高度、图像本身（PIL_image和torch_tensor）
- SceneInfo：稀疏点云的位置和颜色、所有相机位置的中心到最远相机的距离
值得注意的是，从 COLMAP 中得到的 qvec 和 tvec 是 W2C 的，即 qvec 表示世界系的基底在相机系中的姿态，tvec 表示世界系的原点在相机系中的位置。以坐标变换来理解，设世界系为旧基，世界系下坐标为\(\mathbf{x}=(x_1, x_2, \cdots, x_n)^\top\)；相机系为新基，相机系下坐标为\(\mathbf{y}=(y_1, y_2, \cdots, y_n)^\top\)；过渡矩阵\(\mathbf{P}\)为用世界系表示相机系的变换矩阵，那么 qvec 对应的旋转矩阵和 tvec 对应的平移向量构成的变换矩阵是用相机系表示世界系的变换矩阵，即\(\mathbf{P}^{-1}\)。根据坐标变换公式\(\mathbf{y} = \mathbf{P}^{-1} \mathbf{x}\)，因此\(\mathbf{P}^{-1}\)记为 W2C 矩阵，即世界坐标系转换到相机坐标系的变换矩阵。另外，qvec已经做了归一化。
相机模型
相机模型是物理相机成型过程的数学模型，描述了三维坐标点到二维成像平面的几何映射关系；将这个映射过程中涉及的参数提取出来代表相机在成像中的作用，就是相机模型。相机模型中的参数分为内参和外参，不同种类相机的内参数量和形式是不同的，最简单的针孔相机内参包含两个轴向的焦距和主点偏移，其它相机模型如鱼眼相机、全向相机的内参还有各种畸变参数等；而外参则具有统一的形式，即相机的位置和姿态。高斯泼溅中使用的相机模型为针孔相机，故本文也以针孔相机为前提进行下面的介绍。为了使用线性代数的工具完成映射过程中的坐标变换，将参数组织成矩阵的形式，即内参矩阵和外参矩阵。 $$ \begin{aligned} 外参矩阵&: \mathbf{T}= \begin{bmatrix} \mathbf{R} & \mathbf{t}
\mathbf{0} & 1 \end{bmatrix}= \begin{bmatrix} \mathbf{I} & \mathbf{t}
\mathbf{0} & 1 \end{bmatrix} \times \begin{bmatrix} \mathbf{R} & \mathbf{0}
\mathbf{0} & 1 \end{bmatrix}= \begin{bmatrix} 1 & 0 & 0 & t_1
0 & 1 & 0 & t_2
0 & 0 & 1 & t_3
0 & 0 & 0 & 1 \end{bmatrix} \times \begin{bmatrix} r_{1,1} & r_{1,2} & r_{1,3} & 0
r_{2,1} & r_{2,2} & r_{2,3} & 0
r_{3,1} & r_{3,2} & r_{3,3} & 0
0 & 0 & 0 & 1 \end{bmatrix}\
内参矩阵&: \mathbf{K}= \begin{bmatrix} f_x&0&c_x
0&f_y&c_y
0&0&1 \end{bmatrix} \end{aligned} $$
整个成像的过程可以拆解成一系列坐标系之间的转化：世界坐标系–>相机坐标系–>归一化设备坐标系–>图像坐标系。以齐次坐标\([x, y, z, w]^\top\)表示三维空间中的点，外参矩阵用于将点从世界坐标系变换到相机坐标系；然后通过透视除法(即将齐次坐标除以深度分量\(z\))转化到归一化设备坐标系；最后再通过内参矩阵变换到图像坐标系。
高斯椭球
对于 SfM 获得的点云，用三维的高斯函数将它们膨胀成带有位置、形状、颜色和不透明度的椭球。这些参数都是可学习的，会在后续的流程中以深度学习的方法进行优化，且为了保证这些参数的取值符合自身的物理意义(主轴长度为正、旋转矩阵为正交矩阵、协方差矩阵为半正定对称矩阵、不透明度在0-1之间)，会将这些参数通过所谓的激活函数之后再使用。另外，这个所谓的膨胀，其实是用 KNN 取最近邻三个点的平均半径作为初始化高斯椭球的主轴长度，并以点的颜色作为椭球的初始颜色。
之所以选择高斯椭球作为场景的表达基元而不是立方体或球或其它形式，是因为高斯函数具有良好的数学性质，一是高斯函数对仿射变换封闭(渲染中需要用到仿射变换)，二是高斯函数可微，且三维高斯沿某一个轴积分成二维后仍为高斯函数(渲染实质上就是投影到二维平面)。
位置和形状
这个所谓的椭球形状，是由高斯函数的数学表达式决定的(高维高斯均假设各向异性，即不同方向梯度不同)
\[\begin{aligned} G(x) &= \frac{1}{\sqrt{2\pi\sigma^2}} \exp\left( -\frac{(x - \mu)^2}{2\sigma^2} \right)\\ G(\mathbf{x}) &= \frac{1}{2\pi \sqrt{|\mathbf{\Sigma}|}} \exp\left( -\frac{1}{2} (\mathbf{x} - \mathbf{\mu})^\top \mathbf{\Sigma}^{-1} (\mathbf{x} - \mathbf{\mu}) \right),\\ &\text{其中 } \mathbf{x} = \begin{bmatrix} x_1 \\ x_2 \end{bmatrix},\quad \mathbf{\mu} = \begin{bmatrix} \mu_1 \\ \mu_2 \end{bmatrix},\quad \mathbf{\Sigma} = \begin{bmatrix} \sigma_{11} & \sigma_{12} \\ \sigma_{21} & \sigma_{22} \end{bmatrix}\\ G(\mathbf{x}) &= \frac{1}{(2\pi)^{3/2} \sqrt{|\mathbf{\Sigma}|}} \exp\left( -\frac{1}{2} (\mathbf{x} - \mathbf{\mu})^\top \mathbf{\Sigma}^{-1} (\mathbf{x} - \mathbf{\mu}) \right),\\ &\text{其中 } \mathbf{x} = \begin{bmatrix} x_1 \\ x_2 \\ x_3 \end{bmatrix},\quad \mathbf{\mu} = \begin{bmatrix} \mu_1 \\ \mu_2 \\ \mu_3 \end{bmatrix},\quad \mathbf{\Sigma} = \begin{bmatrix} \sigma_{11} & \sigma_{12} & \sigma_{13} \\ \sigma_{21} & \sigma_{22} & \sigma_{23} \\ \sigma_{31} & \sigma_{32} & \sigma_{33} \end{bmatrix}\\ \mathbf{\mu}为均值&，决定分布的位置\\ \mathbf{\Sigma}为协方&差矩阵，为半正定的对称矩阵，对角线元素为方差，其余元素为协方差，决定分布的形状 \end{aligned}\]先以二维为例，令二维高斯函数为某一常数，提取各高次项系数并整理之后可以得到一个标准的椭圆方程；这一点从几何意义上理解更直观：想象一个平面与钟形曲面相交，截出的曲线就是椭圆。三维高斯函数同理，令其为常数并提取各高次项系数，整理之后就是一个标准的椭球面方程。在几何意义上为了可视化出来，可以想象用颜色而非欧式距离来代表函数值，这样，三维高斯函数就是一个实心椭球，从内到外就是颜色渐变的椭球壳。
高斯椭球的位置和形状分别由三维高斯分布中的均值和协方差矩阵确定，均值确定位置很好理解，下面解释协方差如何确定形状。
任意高斯分布可由标准高斯分布通过仿射变换得到，对高斯函数进行仿射变换满足如下规律
\[\begin{aligned} \mathbf{x} &\sim \mathcal{N}(\mathbf{\mu}, \mathbf{\Sigma})\\ \mathbf{y} &= \mathbf{A} \mathbf{x} + \mathbf{b}\\ \mathbf{y} &\sim \mathcal{N}(\mathbf{A} \mathbf{\mu} + \mathbf{b},\ \mathbf{A} \mathbf{\Sigma} \mathbf{A}^\top)\\ \end{aligned}\]由于仿射变换本身可以表达为旋转、缩放和平移(\(\mathbf{A} \mathbf{x} + \mathbf{b} = \mathbf{R}\mathbf{S} \mathbf{x} + \mathbf{b}\))的组合。标准正态分布的协方差矩阵\(\mathbf{\Sigma}=\mathbf{I}\)，因此变换之后的协方差矩阵\(\mathbf{\Sigma’}\)可以写成下式，即用旋转矩阵和缩放矩阵构成协方差矩阵。
\[\mathbf{\Sigma’} = \mathbf{A}\mathbf{I}\mathbf{A}^\top = \mathbf{R}\mathbf{S}\mathbf{S}^\top\mathbf{R}^\top\]之所以要这样分解，是为了在参数学习的过程中自然地满足协方差矩阵为半正定的对称矩阵这一约束(旋转矩阵为正交矩阵，满足转置等于逆的性质，而缩放矩阵为对角矩阵，两者合成为\(\mathbf{R}~ \text{diag}(\sigma_1^2\ \sigma_2^2\ \sigma_3^2) ~\mathbf{R}^\top\)，满足半正定对称矩阵的要求)。旋转矩阵控制了椭球绕三个坐标轴的旋转角度，缩放矩阵控制了椭球三个主轴的长度，因此它们构成的协方差矩阵决定了椭球的形状。另外，这一点如果从代数上理解，其实就是矩阵二次型对应的几何图像。
颜色
高斯椭球的颜色由球谐函数确定，高斯泼溅中以此实现不同视角下观察到不同颜色的效果。球谐函数是球坐标系下，定义在单位球面的一组具有正交完备性和旋转不变性的基函数。
之所以使用球坐标系函数，是为了将颜色与观测角度建立联系。而基函数是函数空间中的概念，基函数在其中的作用与向量空间中基向量的作用相同，即基函数可以通过线性组合张成函数空间中的任意函数（其实函数空间就是一种抽象的向量空间，同样满足向量空间中诸如对加法和数乘封闭等一系列性质）。球谐函数是球坐标系下的基函数，类似于笛卡尔坐标系下傅里叶级数中的\(\cos(nx)\)和\(\sin(nx)\)的概念，任何球面坐标系的函数都可以展开为无穷个球谐函数，换句话说，多个球谐函数的线性组合，可以近似表示任意的球坐标函数。
球谐函数本身，是球坐标系下的拉普拉斯方程在分离变量\(r\)后，角度部分通解的正交项。 \(\begin{aligned} 球坐标拉普拉斯方程:&\nabla^2 f = \frac{1}{r^2} \frac{\partial}{\partial r} \left( r^2 \frac{\partial f}{\partial r} \right) + \frac{1}{r^2 \sin\theta} \frac{\partial}{\partial \theta} \left( \sin\theta \frac{\partial f}{\partial \theta} \right) + \frac{1}{r^2 \sin^2\theta} \frac{\partial^2 f}{\partial \varphi^2}\\ 分离变量后角度部分的方程:& \frac{1}{\sin\theta} \frac{\partial}{\partial\theta} \left( \sin\theta \frac{\partial Y}{\partial\theta} \right) + \frac{1}{\sin^2\theta} \frac{\partial^2 Y}{\partial\varphi^2} + l(l + 1)Y = 0\quad (球函数方程)\\ \end{aligned}\) 球函数方程可以进一步分离变量，最终得到两个常微分方程，分别求解再回代可以得到球函数方程的解如下。其中关于方位角\(\theta\)的方程又称为连带勒让德方程，方程的解称为连带勒让德函数\(P_l^m (x)\)，全名为\(l\)次\(m\)阶连带勒让德函数(\(\text{degree} = l, \text{order} = m\))，在计算机中以递推的方式实现。 $$ \begin{aligned} 球函数方程:&Y(\theta, \varphi) = \Theta(\theta) \Phi(\varphi)
解得:&\Phi(\varphi) = \text{e}^{\text{i} m\varphi} ,\, m = 0, \pm 1, \cdots, \pm l &\quad(极角方程的解)
&\Theta(\theta) = P_l^m (\cos\theta) ,\, m = 0, \pm 1, \cdots, \pm l &\quad(方位角方程的解)
其中&P_l^m(x) = (-1)^m (1 - x^2)^{\frac{m}{2}} \frac{\text{d}^m}{\text{d}x^m} \left(P_l(x)\right)&\quad(递推公式)
&P_l(x) = \frac{1}{2^l \cdot l!} \frac{\text{d}^l}{\text{d}x^l} \left[(x^2 - 1)^l\right]&\quad(通项公式)\
球函数方程的解为:&Y(\theta, \varphi) = \sum_{l=0}^{\infty} \sum_{m=-l}^{l} P_l^m (\cos\theta) \text{e}^{\text{i} m\varphi}, \, m = 0, \pm 1, \pm 2, \cdots
\end{aligned} \(高斯泼溅中使用的是归一化的球谐函数的分量，归一化即在球函数方程解的基础上增加了归一化系数，分量是指在\)m>0\(时取\)\text{e}\(指数的实数部分，\)m<0\(时取\)\text{e}\(指数的虚数部分，最终得到高斯泼溅中球谐函数的表达式如下。\) \begin{aligned} Y_l^m (\theta, \varphi) &= K_l^m Y_{lm} (\theta, \varphi)
&=\begin{cases} \sqrt{2} K_l^m \cos(m\varphi) P_l^m(\cos\theta) & m > 0
\sqrt{2} K_l^m \sin(-m\varphi) P_l^{-m}(\cos\theta) & m < 0
K_l^0 P_l^0(\cos\theta) & m = 0 \end{cases}
其中K_l^m &= \sqrt{ \frac{2l + 1}{4\pi} \frac{(l - |m|)!}{(l + |m|)!} }
\end{aligned} $$ 解释完了球谐函数本身，再聊聊正交完备性和旋转不变性。正交完备性包括正交归一性和完备性两个性质。函数空间的正交归一性定义为不同函数在定义域上的内积为0，函数自身的内积为1。完备性是指函数展开为基函数表示时，当次数趋于无穷时，级数和平均收敛于原函数，这保证了展开结果会趋近于被展开的函数。而旋转不变性是指，当原函数发生旋转时，只需对生成原函数的线性组合系数(也称为广义傅里叶系数)进行相应的变换，就可以保证变换后的系数能等价还原出新函数，而无需对新函数进行展开并计算广义傅里叶系数(生成系数的计算开销很大)。这具体表现为，当光源发生旋转时，只要同步计算出变换后的广义傅里叶系数，就能保证画面的光照效果不会抖动跳变。
对于给定的球谐函数次数，球谐函数和对应的球谐系数的数量是确定的，用于表示颜色的原函数可以按如下式子展开(以高斯泼溅中支持的最高次数三次为例)，其中球谐函数线性组合的系数就称为球谐系数。由原函数生成线性组合系数的过程称为投影，而用系数和基函数组合成原函数这个逆过程称为重建。 \(\begin{aligned} f &= \sum_{l} \sum_{m=-l}^{l} \mathbf{c}_{l}^{m} y_{l}^{m}(\theta, \phi) \\ &= \mathbf{c}_{0}^{0} y_{0}^{0}+ \\ &\quad \mathbf{c}_{1}^{-1} y_{1}^{-1}+\mathbf{c}_{1}^{0} y_{1}^{0}+\mathbf{c}_{1}^{1} y_{1}^{1}+ \\ &\quad \mathbf{c}_{2}^{-2} y_{2}^{-2}+\mathbf{c}_{2}^{-1} y_{2}^{-1}+\mathbf{c}_{2}^{0} y_{2}^{0}+\mathbf{c}_{2}^{1} y_{2}^{1}+\mathbf{c}_{2}^{2} y_{2}^{2}+ \\ &\quad \mathbf{c}_{3}^{-3} y_{3}^{-3}+\mathbf{c}_{3}^{-2} y_{3}^{-2}+\mathbf{c}_{3}^{-1} y_{3}^{-1}+\mathbf{c}_{3}^{0} y_{3}^{0}+\mathbf{c}_{3}^{1} y_{3}^{1}+\mathbf{c}_{3}^{2} y_{3}^{2}+\mathbf{c}_{3}^{3} y_{3}^{3} \end{aligned}\) 与函数的傅里叶级数类似，球谐函数的次数\(l\)越高，实现的颜色表达就越细致。在三阶的情况下，按RGB三通道计算，完整的球谐系数存储为一个\(3\times 16\)的矩阵。
光栅化
在 3DGS 中，这个从三维场景到二维图像的渲染过程称为光栅化。总体来说，光栅化将每个 3D 高斯投影到 2D 屏幕上，并按深度顺序将其颜色和透明度累加到对应像素，生成最终图像。
渲染需要用到相机模型中的内参矩阵和外参矩阵，但转化的过程与物理相机成像并不完全相同，渲染过程使用视锥。视锥为一个四棱台，由从观察点出发连接图像四个顶点的四棱锥和近平面远平面相截形成。在图形学中，视锥由远近平面的深度znear、zfar、视场角fov和图像长宽比aspect_ratio确定。
光栅化将每个 3D 高斯投影到 2D 屏幕上，并按深度顺序将其颜色和透明度累加到对应像素，生成最终图像。3D 高斯椭球的投影具体来说可以细分为四个过程：观测变换、投影变换、视口变换、像素离散。
3D 高斯椭球的投影具体来说可以细分为四个过程：观测变换、投影变换、视口变换、像素离散。
观测变换是从世界坐标系转换到相机坐标系的过程，本质上是旋转、缩放和平移三种变换的组合。观测变换是从世界坐标系转换到相机坐标系的过程，本质上是仿射变换，通过乘相机外参矩阵实现。观测变换好比线性代数中基向量的变换，同一向量在不同基下的表达系数不同；同一物体在世界坐标系和相机坐标系中的表达也是不同的，观测变换就是两种坐标系之间的转换。
投影变换是从相机坐标系转换到图像平面的过程，先将视锥通过非线性变换化为正交空间，然后变换到归一化立的立方体空间内，最后正交投影忽略深度便转化到二维图像平面。投影变换是从相机坐标系转换到图像平面的过程，又分为透视投影和正交投影。正交投影即忽略深度，不存在近大远小的效应。透视投影可以看作是先将视锥通过非线性变换化为正交投影相同的情况。两种投影都会先将要观察的空间变换到归一化坐标系(NDC)，即\([-1, 1]^3\)范围内，然后正交投影忽略深度便转化到二维图像平面。透视投影符合人眼成像规律，应用更广泛，变换示意图和相应的变换矩阵如下
要观察的视锥按如下方法确定：SfM 中获得的内参矩阵中，fx 和 fy 分别为水平焦距和竖直焦距，与图像的宽度和高低一起可以算出 FOVx 和 FOVy，再根据指定的 zfar 和 znear，便可以确定要投影的区域和相应的投影矩阵。
视口变换是将投影得到的归一化平面还原到符合图像长宽的图像坐标系的过程，即将\([-1, 1]^2\)的矩形还原到\([0, w] \times [0, h]\)。
变换矩阵为
\[\mathbf{M}_{\text{viewport}} = \begin{bmatrix} \frac{w}{2} & 0 & 0 & \frac{w}{2} \\ 0 & \frac{h}{2} & 0 & \frac{h}{2} \\ 0 & 0 & 1 & 0 \\ 0 & 0 & 0 & 1 \end{bmatrix}\]
像素离散是将连续图像坐标系离散采样到像素坐标系的过程，根据像素中心点是否在图像中决定是否渲染该像素。
高斯球参数在变换中传递
设世界坐标系下的高斯核中心为\(\mathbf{t}_k= [t_{0}, t_{1}, t_{2}]^{\top}\)，协方差矩阵为\(\mathbf{V}''_k\)，高斯核为\(r''(k)=G_{\mathbf{V}''_k}(\mathbf{t} - \mathbf{t}_k)\)；
相机坐标系下的高斯核中心为\(\mathbf{u}_k= [u_{0}, u_{1}, u_{2}]^{\top}\)，协方差矩阵为\(\mathbf{V}'_k\)，高斯核为\(r'(k)=G_{\mathbf{V}'_k}(\mathbf{u} - \mathbf{u}_k)\)；
归一化坐标系的高斯核中心为\(\mathbf{x}_k= [x_{0}, x_{1}, x_{2}]^{\top}\)，协方差矩阵为\(\mathbf{V}_k\)，高斯核为\(r(k)=G_{\mathbf{V}_k}(\mathbf{x} - \mathbf{x}_k)\)；
观测变换为\(\mathbf{u}=\mathbf{W}\mathbf{t} + \mathbf{d}\)，则\(\mathbf{u}_k = \mathbf{W}\mathbf{t}_k + \mathbf{d}，\mathbf{V}'_k = \mathbf{W}\mathbf{V}''_k\mathbf{W}^\top\)；
投影变换为\(\mathbf{x}=m(\mathbf{u})\)，雅可比矩阵为\(\mathbf{J}\)，则\(\mathbf{x}_k = m(\mathbf{u}_k)，\mathbf{V}_k = \mathbf{J}\mathbf{V}'_k\mathbf{J}^\top = \mathbf{J}\mathbf{W}\mathbf{V}''_k\mathbf{W}^\top\mathbf{J}^\top\)。
由于投影变换为非线性变换，而希望高斯椭球只进行线性的仿射变换(非线性变换不能准确传播协方差矩阵)，对于均值这一个点可以直接应用非线性的投影变换，不存在形变的问题；而协方差则需要用雅可比矩阵在高斯核中心进行局部线性化，代替仿射变换矩阵\(\mathbf{W}\)。操作之后，均值在归一化坐标系中，后续再进行视口变换即可；而协方差在未缩放未平移的正交坐标系中，不需要作视口变换。
雅可比矩阵的具体计算方式如下：假设视锥中的点为\([x, y, z, 1]^\top\)(写成齐次形式)，可以求出变换到正交坐标系中的对应坐标，由此可求得雅可比矩阵
\[\begin{aligned} \begin{bmatrix} n & 0 & 0 & 0 \\ 0 & n & 0 & 0 \\ 0 & 0 & n + f & -nf \\ 0 & 0 & 1 & 0 \end{bmatrix} \begin{bmatrix} x \\ y \\ z \\ 1 \end{bmatrix} &= \begin{bmatrix} nx \\ ny \\ (n + f)z - nf \\ z \end{bmatrix} = \begin{bmatrix} \frac{nx}{z} \\ \frac{ny}{z} \\ (n + f)-\frac{nf}{z} \\ 1 \end{bmatrix} = \begin{bmatrix} f(x) \\ f(y) \\ f(z) \\ 1 \end{bmatrix}\\ \mathbf{J} &= \begin{bmatrix} \frac{\text{d}f(x)}{\text{d}x}&\frac{\text{d}f(x)}{\text{d}y}&\frac{\text{d}f(x)}{\text{d}z}\\ \frac{\text{d}f(y)}{\text{d}x}&\frac{\text{d}f(y)}{\text{d}y}&\frac{\text{d}f(y)}{\text{d}z}\\ \frac{\text{d}f(z)}{\text{d}x}&\frac{\text{d}f(z)}{\text{d}y}&\frac{\text{d}f(z)}{\text{d}z} \end{bmatrix} = \begin{bmatrix} \frac{n}{z}&0&-\frac{nx}{z^{2}}\\ 0&\frac{n}{z}&-\frac{ny}{z^{2}}\\ 0&0&\frac{n f}{z^{2}} \end{bmatrix} \end{aligned}\]图像合成
图像合成采用\(\alpha\)-blending 的方法，对于每一个区块，将高斯椭球按深度从近到远排序，由近到远遍历每个高斯椭球来填充像素的不透明度，直到该区块内所有像素的\(\alpha=1\)(完全不透明)为止；若混合完所有高斯椭球仍未达到\(\alpha=1\)则用背景颜色填充剩余不透明度。
渲染过程中将图像拆分成\(16\times16\)的区块(tile)，高斯椭球也随之被划分。GPU 为每个区块分配一个线程，区块之内共享内存。图像合成采用\(\alpha\)-blending 的方法，对于每一个区块，将高斯椭球按深度从近到远排序，由近到远遍历每个高斯椭球来填充像素的不透明度，直到该区块内所有像素的\(\alpha=1\)(完全不透明)为止；若混合完所有高斯椭球仍未达到\(\alpha=1\)则用背景颜色填充剩余不透明度。渲染过程中仍然执行以像素为单位的足迹渲染，但没有 NeRF 中从像素找粒子的过程。
前向渲染中，不透明度\(\alpha\)通过以下递推公式对图像产生作用，其中\(C_{acc}\)为之前累积混合得到的颜色，\(T_{acc}\)为之前累积的透射率(初始时\(T_0=1\))，\(C_{\text{final}}\)为当前这一层混合完成之后的颜色
\[\begin{aligned} \mathbf{C}_\text{final}&=\mathbf{C}_i\cdot\alpha_i+\mathbf{C}_\text{acc}\cdot(1 - \alpha_i)\\ T_\text{final}&=T_\text{acc}\cdot(1 - \alpha_i) \end{aligned}\]推导出等价的通项公式如下
\[\begin{aligned} \hat{C}(r)&=\sum_{i = 1}^{N} T_i\alpha_i c_i\\ T_i&= \prod_{j = 1}^{i - 1} (1 - \alpha_j) \end{aligned}\]
在公式上，\(\alpha\)-blending 的方法和体渲染是一致的，但实现却完全不同。3DGS 中的\(\alpha\)-blending 是将离散的几个高斯椭球进行投影和颜色混合；而 NeRF 中的体渲染是沿光线进行密集采样并求和(将沿光线的连续积分离散化)。
体渲染公式
连续域中的公式为一个第一类曲线积分(此处省略推导)，其中\(T(s)\)为在\(s\)之前，光线未被遮挡的概率；\(\sigma(s)\)为\(s\)处光线撞击粒子(被阻碍)的概率密度；\(C(s)\)为光线在\(s\)处造成的颜色；\(\hat{C}(s)\)为渲染得到的颜色
\[\begin{aligned} \hat{C}(s) &= \int_{0}^{+\infty} T(s) \sigma(s)C(s) \text{d}s \\ T(s) &= e^{-\int_{0}^{s} \sigma(t) dt} \end{aligned}\]离散域的公式如下(此处省略推导)，将光线\([0, s]\)划分为\(N\)个等距区间\([T_n\to T_{n+1}], n=0,1,\cdots,N\)，区间间隔长度\(\delta_n\)，区间足够小，认为区间内概率密度\(\sigma_n\)和颜色\(c_n\)固定
\[\begin{aligned} \hat{C}(r)&=\sum_{i = 1}^{N} T_i\alpha_i c_i\\ \alpha_i &=1 - \text{e}^{-\sigma_i\delta_i}\\ T_i&=\text{e}^{-\sum_{j = 1}^{i - 1}\sigma_j\delta_j}= \prod_{j = 1}^{i - 1} (1 - \alpha_j) \end{aligned}\]马氏距离
马氏距离是一种距离指标，可以应对高维线性分布的数据中各维度间非独立同分布的问题，修正了欧式距离中各维度尺度不一致且相关的问题，本质上是对PCA中的主成分进行标准化，将变量按主成分旋转（主成分就是特征向量的方向），使各维度之间相互独立，然后再标准化（缩放特征值）使它们同分布。马氏距离其实就是标准化之后的欧式距离 \(\begin{aligned} D_M(\mathbf{x}) &= \sqrt{(\mathbf{x} - \mathbf{\mu})^\top\Sigma^{-1}(\mathbf{x} - \mathbf{\mu})}\\ D_M(\mathbf{x}, \mathbf{y}) &= \sqrt{(\mathbf{x} - \mathbf{y})^\top\Sigma^{-1}(\mathbf{x} - \mathbf{y})} \end{aligned}\) EWA（椭圆加权平均）是一种高质量的采样和抗锯齿算法，适合处理椭圆形状分布的元素，其核心就是基于2d cov的马氏距离。
参数优化
通过以上步骤，可以获得给定视角下的渲染图像，与真值图像作 loss 之后，经过反向传播和梯度下降便可逐步调整高斯椭球的各个参数，使得渲染图像越发接近真值图像，从而完成对整个场景的重建。
3DGS 中，这个损失函数如下，其中\(L_1\)项代表逐像素绝对误差，而 D-SSIM 是一种衡量两张图像结构相似度的感知指标，考虑了亮度、对比度和结构信息，相应地，这项损失代表了图像整体结构上的误差。 \(L = (1-\lambda)L_1 + \lambda L_{D-SSIM}，原版高斯泼溅中取\lambda=0.2\)
补丁和工程 tricks
分图块以并行训练、自适应密度控制