3d gaussian splatting公式推导

神气爱哥

已于 2024-11-12 09:44:21 修改

阅读量1.1k

点赞数 5

分类专栏：自动驾驶实战计算机图形学计算机视觉文章标签：自动驾驶计算机视觉算法

于 2024-09-22 12:08:28 首次发布

本文链接：https://blog.youkuaiyun.com/qingcaichongchong/article/details/142433156

版权

计算机视觉同时被 3 个专栏收录

16 篇文章

订阅专栏

自动驾驶实战

9 篇文章

订阅专栏

计算机图形学

6 篇文章

订阅专栏

1. 正向渲染公式推导

正向渲染是指从3D高斯球渲染到2D屏幕中某个像素值的颜色的过程。

1.1 3D高斯到2D高斯

3D高斯参数有高斯球的中心即均值 $\omega$ ，协方差 $\Sigma$ ，体密度 $\sigma$ ，颜色 $c$ ，2D高斯参数有高斯球的中心即均值 $\mu$ ，协方差 $\Omega$ 。

首先，对于高斯球上一个x，将其变换到相机坐标系下有： $p_{c}=T_{w2c}x$ ，然后再将其变换到2D高斯上一个u，有：

$u=\left\{\begin{matrix}p_c[0]/p_c[2] \\ p_c[1]/p_c[2] \\ \left \| p_c \right \| \end{matrix}\right.$

对于2D高斯的均值 $\mu$ 来说，则有： $p_{c}=T_{w2c}\omega$ ，

$\mu=\left\{\begin{matrix}p_c[0]/p_c[2] \\ p_c[1]/p_c[2] \\ \left \| p_c \right \| \end{matrix}\right.$

接下来是2D高斯协方差 $\Omega$ 的求法，三维高斯球 $\Sigma$ 是一个椭球体，它可以通过对半径为1的标准圆施加旋转R和缩放矩阵S生成，即 $\Sigma=RS\Lambda S^TR^T=RSS^TR^T$ , 这是因为标准正态分布的协方差 $\Lambda$ 就是一个半径为1的标准圆，椭球就是将这个圆进行旋转缩放变成的。 $\Lambda$ 是一个单位矩阵，所以可以省略。

由于3D转2D的变换是非线性变换，因此在3D协方差转2D协方差时需要做线性化，假设原先的3D高斯分布指数表示为 $(x-\omega)^T\Sigma^{-1}(x-\omega)$ ，2D高斯分布表示为 $(u-\mu)^T\Omega^{-1}(u-\mu)$ ，那么在 $\omega$ 处对x进行线性化，设 $u=Jx$ ， $J$ 就是2D坐标对3D坐标的雅可比矩阵，那么带入上式有：

$(Jx-\mu))^T\Omega^{-1}(Jx-\mu))=(x-J^{-1}\mu)J^T\Omega^{-1}{J}(x-J^{-1}\mu)$

其中 $J^{-1}\mu$ 是均值, 我们会用 $\omega$ 而不是 $J^{-1}\mu$ 。

对比发现： $\Sigma^{-1}=J^T\Omega^{-1}{J}$ ，即 $\Omega^{-1}=J^{-T}\Sigma^{-1}J^{-1}=(J\Sigma{J^T})^{-1}$ 。

于是， $\Omega=J\Sigma{J^T}$ , 对应论文公式（5），其中 $J=\frac{\partial u}{\partial x}$ ，为何少了一个旋转，是因为

$J=\frac{\partial u}{\partial x}=\frac{\partial u}{\partial p_c}\frac{\partial p_c}{\partial x}=\bigl(\begin{smallmatrix} \frac{1}{{p_c}[2]} & 0 &-\frac{pc[0]}{{{p_c}^2[2]}} \\ 0 & \frac{1}{p_c[2]} & -\frac{p_c[1]}{{p_c}^2[2]} \\ \frac{p_c[0]}{\left \| p_c \right \|} & \frac{p_c[1]}{\left \| p_c \right \|} & \frac{p_c[2]}{\left \| p_c \right \|}\end{smallmatrix}\bigr) * R_{w2c}$

其中， $R_{w2c}$ 是 $T_{w2c}$ 的旋转部分， $R_{w2c}$ 就是论文中所述的W。

由于2D高斯转换到2D高斯后只剩二维坐标了，所以通常取 $\Omega$ 的前两行两列作为2D高斯球的协方差矩阵。

1.2 2D高斯到像素颜色

gaussian splatting渲染实际是借鉴了nerf的积分方式，步骤是：

1. 图像中任选一个像素，遍历所有2D高斯球

2. 选择与该像素的距离不小于3倍方差的2D高斯球，则将其加入到候选集合

3. 按深度对候选集合排序，然后从近到远进行颜色的计算，计算公式如下：

颜色计算的离散公式推导：nerf中连续的积分渲染公式是：

$C(r)=\int_{0}^{t_i}T(t)\sigma(r(t))c(r(t), d)dt$

其中 $\sigma$ 是体密度，T(t)表示为被遮挡率，它的计算公式是：

$T(t)=exp(-\int_{0}^{t_i}\sigma(r(s))ds)$

那么转换为离散公式后有：

$T_i=T(t_i) =exp(-\sum_{j=0}^{j=i-1}\int_{t_j}^{t_{j+1}}\sigma(t)dt) =exp(-\sum_{j=0}^{j=i-1}\sigma_j\delta_j)$

其中， $\delta_j$ 代表j时刻的时间差，将其带入渲染公式：

$C(i)=\sum_{j=0}^{j=i-1}\int_{t_j}^{t_{j+1}}T(t)\sigma(t)c(t)dt \\=\sum_{j=0}^{j=i-1}\int_{t_j}^{t_{j+1}}exp(-\int_{0}^{t}\sigma(s)ds)\sigma_jc_jdt \\=\sum_{j=0}^{j={i-1}}\int_{t_j}^{t_{j+1}}exp(-\int_{0}^{t_j}\sigma(s)ds)exp(-\int_{t_j}^{t}\sigma(s)ds)\sigma_jc_jdt \\=\sum_{j=0}^{j=i-1}\int_{t_j}^{t_{j+1}}T_jexp(-\int_{t_j}^{t}\sigma(s)ds)\sigma_jc_jdt \\=\sum_{j=0}^{j=i-1}T_j\sigma_jc_j\int_{t_j}^{t_{j+1}}exp(-\sigma_j(t-t_j)dt \\=\sum_{j=0}^{j=i-1}T_j\sigma_jc_j(-\frac{1}{\sigma_j})(exp(-\sigma_j\delta_j)-1) \\=\sum_{j=0}^{j=i-1}T_jc_j(1-exp(-\sigma_j\delta_j))$

设透明度 $\alpha_i=(1-exp(-\sigma_i\delta_i))$

则被遮挡率 $T_i =\prod_{j=1}^{j=i-1}(1-\alpha_i)$

有 $C(i)=\sum_{j=0}^{j=i-1}c_j\alpha_j\prod_{j=1}^{j=i-1}(1-\alpha_j)$

而gaussian-splating的公式与nerf不一样的点在于，透明度 $\alpha_i$ 的计算是根据像素点处于2D高斯球的概率求得的：

$\alpha_i=\frac{1}{2\pi\sqrt{\left | \Omega \right |}}exp(-0.5(x-\mu)^T\Omega ^{-1}(x-\mu))*\sigma_i$

其中 $\mu$ 是3D高斯球溅射到屏幕上形成的2D高斯球的均值, x是屏幕坐标， $\Omega$ 是2D高斯球的协方差, $\sigma_i$ 是3D高斯球的体密度。

2. 反向传播公式推导

反向传播与正向渲染的方向是相反的，指从2D屏幕中某个像素值的颜色如何将梯度传导到3D高斯球的过程，从而更新3D高斯球的参数，需要更新的参数有：表示3D高斯球的协方差 $\Sigma$ 的参数R,S，体密度 $\sigma$ , 以及颜色c。

2.1 从像素到2D高斯

以任意一个像素i为例，其像素颜色可表示为：

$C(i)=\alpha_1c_1+\alpha_2(1-\alpha_1)c_2+...+\alpha_j\prod_{k=0}^{j-1}(1-\alpha_k)c_j+...$

那么，对于序号为j的高斯球来说，它的各项参数的导数的计算方法是，

对于颜色 $c_j$ 的导数显而易见：

$\frac{\partial C(i)}{\partial c_j}=\alpha_j\prod_{k=0}^{j-1}(1-\alpha_k)$

对于透明度 $\alpha_i$ 的导数可以求出来：

$\frac{\partial C(i)}{\partial \alpha_j}=\prod_{k=0}^{j-1}(1-\alpha_k)c_j-\frac{1}{1-\alpha_j}*\sum_{m=j+1}^{n}\prod_{k=0}^{m-1}(1-\alpha_k)\alpha_mc_m$

其中n是与该像素有交集的高斯球的个数，前面为什么除 $1-\alpha_j$ 是因为后面从k=0到m-1时也会乘以 $1-\alpha_j$ ,所以要把这个因素考虑在内。

根据链式求导法则，对于体密度 $\sigma_j$ 的导数也就可以求出来了:

$\frac{\partial C(i)}{\partial \sigma_j}=\frac{\partial C(i)}{\partial \alpha_j}*\frac{1}{2\pi\sqrt{\left | \Omega \right |}}exp(-0.5(x-\mu)^T\Omega ^{-1}(x-\mu))$

对于2D高斯协方差 $\Omega$ 的导数比较复杂，求解过程如下：

将 $\Omega$ 展开：

$\Omega=\begin{pmatrix} \Omega_{0,0} & \Omega_{0,1} \\ \Omega_{1,0} & \Omega_{1,1} \end{pmatrix}$

根据二维矩阵求逆公式有：

$\Omega^{-1}=\frac{\begin{pmatrix} \Omega_{1,1} & -\Omega_{0,1} \\ -\Omega_{1,0} & \Omega_{0,0} \end{pmatrix}}{\left | \Omega \right |}$

为了方便表示，设

$D=(x-\mu) = (D_{u}, D_{v})^T$ ， $\left | \Omega \right |=\Omega_{0, 0}\Omega_{1, 1}-\Omega_{0, 1}\Omega_{1, 0}$

$P_0=\frac{1}{2\pi\sqrt{\left | \Omega \right |}}$ ， $P_1=exp(-0.5D^T\Omega^{-1}D)$ ，

将 $\Omega^{-1}$ 按行列式展开得：

$D^T\Omega^{-1}D=\frac{\Omega_{1,1}D_u^2-\Omega_{0,1}D_uD_v-\Omega_{1,0}D_uD_v+\Omega_{0,0}D_v^2}{\left | \Omega \right |}\\=\frac{\Omega_{1,1}D_u^2-\Omega_{0,1}D_uD_v-\Omega_{1,0}D_uD_v+\Omega_{0,0}D_v^2}{\Omega_{0,0}\Omega_{1,1}-\Omega_{0,1}\Omega_{1,0}}=\frac{PP}{Det}$

那么，

$\frac{\partial C(i)}{\partial \Omega_{0, 0}}=\frac{\partial C(i)}{\partial \alpha_j} *\sigma_j * ( \frac{\partial C(i)}{\partial P_0}*\frac{\partial P_0}{\partial \Omega_{0, 0}} + \frac{\partial C(i)}{\partial P_1} * \frac{\partial P_1}{\partial \Omega_{0, 0}}) \\\\ = \frac{\partial C(i)}{\alpha_j} * \sigma_j * (-P_1*\frac{1}{4\pi\left | \Omega \right | ^{\frac{3}{2}}}*\Omega_{1,1}-0.5*P_0*exp(-\frac{D_u^T\Omega^{-1}D_u}{2})*(\frac{D_v^2*Det-PP*\Omega_{1, 1}}{Det^2})$

同理，

$\frac{\partial C(i)}{\partial \Omega_{0, 1}}=\frac{\partial C(i)}{\partial \alpha_j} *\sigma_j * (\frac{\partial C(i)}{\partial P_0}*\frac{\partial P_0}{\partial \Omega_{0, 1}} + \frac{\partial C(i)}{\partial P_1} * \frac{\partial P_1}{\partial \Omega_{0, 1}}) \\\\=\frac{\partial C(i)}{\alpha_j} *\sigma_j * (P_1*\frac{1}{4\pi\left | \Omega \right | ^{\frac{3}{2}}}*\Omega_{1,0}-0.5*P_0*exp(-\frac{D_u^T\Omega^{-1}D_u}{2})*(\frac{-D_uD_v*Det+PP*\Omega_{1, 0}}{Det^2})$

$\frac{\partial C(i)}{\partial \Omega_{1, 0}}=\frac{\partial C(i)}{\partial \alpha_j} *\sigma_j * (\frac{\partial C(i)}{\partial P_0}*\frac{\partial P_0}{\partial \Omega_{1, 0}} + \frac{\partial C(i)}{\partial P_1} * \frac{\partial P_1}{\partial \Omega_{1, 0}}) \\\\=\frac{\partial C(i)}{\alpha_j} *\sigma_j * (P_1*\frac{1}{4\pi\left | \Omega \right | ^{\frac{3}{2}}}*\Omega_{0,1}-0.5*P_0*exp(-\frac{D_u^T\Omega^{-1}D_u}{2})*(\frac{-D_uD_v*Det+PP*\Omega_{0, 1}}{Det^2})$

$\frac{\partial C(i)}{\partial \Omega_{1, 1}}=\frac{\partial C(i)}{\partial \alpha_j} *\sigma_j * (\frac{\partial C(i)}{\partial P_0}*\frac{\partial P_0}{\partial \Omega_{1, 1}} + \frac{\partial C(i)}{\partial P_1} * \frac{\partial P_1}{\partial \Omega_{1, 1}}) \\\\=\frac{\partial C(i)}{\alpha_j} *\sigma_j * (-P_1*\frac{1}{4\pi\left | \Omega \right | ^{\frac{3}{2}}}*\Omega_{0,0}-0.5*P_0*exp(-\frac{D_u^T\Omega^{-1}D_u}{2})*(\frac{D_u^2*Det-PP*\Omega_{0, 0}}{Det^2})$

2.2从2D高斯到3D高斯

从2D高斯到3D高斯的反向传播即将2D高斯中的协方差 $\Omega$ 的梯度如何传递给对表示3D高斯协方差的两个参数R,S的梯度，并对R,S进行更新的过程，即 $\frac{\partial C(i)}{\partial S}$ 和 $\frac{\partial C(i)}{\partial R}$ ，根据链式求导法则，其可以表示为

$\frac{\partial C(i)}{\partial S}=\frac{\partial C(i)}{\partial \Omega}\frac{\partial \Omega}{\partial \Sigma}\frac{\partial \Sigma}{\partial M}\frac{\partial M}{\partial S}$ 以及 $\frac{\partial C(i)}{\partial R}=\frac{\partial C(i)}{\partial \Omega}\frac{\partial \Omega}{\partial \Sigma}\frac{\partial \Sigma}{\partial M}\frac{\partial M}{\partial R}$ , 这里的M代表RS。

首先，求 $\frac{\partial \Omega}{\partial \Sigma}$ ，根据前面的推导有 $\Omega=J\Sigma{J^T}$ ，展开来：

$\begin{pmatrix} \Omega_{0, 0} & \Omega_{0, 1}\\ \Omega_{1, 0} & \Omega_{1, 1} \end{pmatrix} =\begin{pmatrix} J_{0,0} & J_{0,1} & J_{0,2} \\ J_{1,0} & J_{1,1} & J_{1,2} \end{pmatrix} \begin{pmatrix} \Sigma_{0, 0} & \Sigma_{0, 1} & \Sigma_{0, 2} \\ \Sigma_{1, 0} & \Sigma_{1, 1} & \Sigma_{1, 2} \\ \Sigma_{2, 0} & \Sigma_{2, 1} & \Sigma_{2, 2} \end{pmatrix} \begin{pmatrix} J_{0, 0} & J_{1, 0} \\ J_{0, 1} & J_{1, 1} \\ J_{0, 2} & J_{1, 2} \end{pmatrix}$

同理对于r从0到2，c从0到2，有：

$\frac{\partial C(i)}{\partial \Sigma_{r, c}}=\frac{\partial C(i)}{\partial \Omega_{0, 0}}\frac{\partial \Omega_{0, 0}}{\partial \Sigma_{r, c}}+\frac{\partial C(i)}{\partial \Omega_{0, 1}}\frac{\partial \Omega_{0, 1}}{\partial \Sigma_{r,c}}+\frac{\partial C(i)}{\partial \Omega_{1, 0}}\frac{\partial \Omega_{1, 0}}{\partial \Sigma_{r,c}}+\frac{\partial C(i)}{\partial \Omega_{1, 1}}\frac{\partial \Omega_{1, 1}}{\partial \Sigma_{r,c}}\\\\=\frac{\partial C(i)}{\partial \Omega_{0, 0}}J_{0, r}J_{0, c}+\frac{\partial C(i)}{\partial \Omega_{0, 1}}J_{0, r}J_{1, j}+\frac{\partial C(i)}{\partial \Omega_{1, 0}}J_{1, i}J_{0, c}+\frac{\partial C(i)}{\partial \Omega_{1, 1}}J_{1, i}J_{1, j}$

然后设 $M=RS$ ，那么 $\Sigma=MM^T$ ，将其按行列展开，

$\Sigma=MM^T=\begin{pmatrix} m_{0, 0}^2+m_{0,1}^2+m_{0,2}^2 & m_{0,0}m_{1,0}+m_{0,1}m_{1,1}+m_{0,2}m_{1, 2} & m_{0,0}m_{2,0}+m_{0,1}m_{2,1}+m_{0,2}m_{2,2}\\ m_{1,0}m_{0,0}+m_{1,1}m_{0,1}+m_{1,2}m_{0,2} &m_{1,0}^2+m_{1,1}^2+m_{1,2}^2 & m_{1,0}m_{2,0}+m_{1,1}m_{2,1}+m_{1,2}m_{2,2} \\ m_{2,0}m_{0,0}+m_{2,1}m_{0,1}+m_{2,2}m_{0,2} & m_{2, 0}m_{1,0}+m_{2, 1}m_{1, 1}+m_{2, 2}m_{1,2} & m_{2,0}^2+m_{2,1}^2+m_{2,2}^2 \end{pmatrix}$

那么对于 $\frac{\partial C(i)}{\partial M}$ ，可得：

$\frac{\partial C(i)}{\partial m_{r,c}}=\sum_{p=0}^{2}\sum_{q=0}^{2} \frac{\partial C(i)}{\partial \Sigma_{p, q}}*\frac{\Sigma_{p, q}}{\partial m_{r, c}} =\sum_{k=0}^{2}(\frac{\partial C(i)}{\partial \Sigma_{k, r}}+ \frac{\partial C(i)}{\partial \Sigma_{r, k}})*m_{k, c}$

最后，求 $\frac{\partial C(i)}{\partial S}$ 和 $\frac{\partial C(i)}{\partial R}$ ,

$M=\begin{pmatrix} m_{0, 0} & m_{0, 1} & m_{0, 2} \\ m_{1, 0} & m_{1, 1} & m_{2, 1} \\ m_{2, 0} & m_{2, 1} & m_{2, 2} \end{pmatrix}= \begin{pmatrix} R_{0,0} & R_{0,1} & R_{0,2}\\ R_{1,0} & R_{1,1} & R_{1,2}\\ R_{2,0} & R_{2,1} & R_{2,2} \end{pmatrix} * \begin{pmatrix} S_0 & & \\ & S_1 & \\ & & S_2 \end{pmatrix}$

对于 $\frac{\partial C(i)}{\partial S}$ ：

$\frac{\partial C(i)}{\partial S_j}|_{j=0,1,2}=\frac {\partial C(i)}{\partial m_{0,j}}*R_{0, j}+\frac {\partial C(i)}{\partial m_{1, j}}* R_{1, j} + \frac{\partial C(i)}{\partial m_{2, j}}*R_{2,j}$

对于 $\frac{\partial C(i)}{\partial R}$ ，先将R写为四元数的方式：

$R=2\begin{pmatrix} 0.5-(q_j^2+q_k^2) & (q_iq_j-q_rq_k) & (q_iq_k+q_rq_j) \\ (q_iq_j+q_rq_k) & 0.5-(q_i^2+q_k^2) & (q_jq_k-q_rq_i)\\ (q_iq_k-q_rq_j) & (q_jq_k+q_rq_i) & 0.5-(q_i^2+q_j^2) \end{pmatrix}$

那么有：

$D^r =\frac{\partial M}{\partial q_r} = 2\begin{pmatrix} 0 & -q_k & q_j \\ q_k & 0 & -q_i \\ -q_j & q_i & 0 \end{pmatrix} * S =2\begin{pmatrix} 0 & -q_ks_1 & q_js_2 \\ q_ks_0 & 0 & -q_is_2 \\ -q_js_0 & q_is_1 & 0 \end{pmatrix}$

那么，

$\frac{\partial C(i)}{\partial q_r}=\sum_{p=0}^{2}\sum_{q=0}^{2}\frac{\partial C(i)}{\partial m_{p,q}}*D^r_{pq}$

同理，对于 $q_i$ ,

$D^i=\frac{\partial M}{\partial q_i} = 2\begin{pmatrix} 0 & q_js_1 & q_ks_2 \\ q_js_0 & -2q_is_1 & -q_rs_2 \\ q_ks_0 & q_rs_1 & -2q_is_2 \end{pmatrix}$

$\frac{\partial C(i)}{\partial q_i}=\sum_{p=0}^{2}\sum_{q=0}^{2}\frac{\partial C(i)}{\partial m_{p,q}}*D^i_{p,q}$

对于 $q_j$ ,

$D^j=\frac{\partial M}{\partial q_j} = 2\begin{pmatrix} -2q_js_0 & q_is_1 & q_rs_2 \\ q_is_0 & 0 & q_ks_2 \\ -q_rs_0 & q_ks_1 & -2q_js_2 \end{pmatrix}$

$\frac{\partial C(i)}{\partial q_i}=\sum_{p=0}^{2}\sum_{q=0}^{2}\frac{\partial C(i)}{\partial m_{p,q}}*D^j_{p,q}$

对于 $q_k$ ,

$D^k=\frac{\partial M}{\partial q_k} = 2\begin{pmatrix} -2q_ks_0 & -q_rs_1 & q_is_2 \\ q_rs_0 & -2q_ks_1 & q_js_2 \\ q_is_0 & q_js_1 & 0 \end{pmatrix}$