3DGS投影变换矩阵全解析:从三维高斯到屏幕椭圆的数学推导

发布时间:2026/9/20 13:40:11
3DGS投影变换矩阵全解析:从三维高斯到屏幕椭圆的数学推导 第一次把3DGS整个渲染管线读通的时候我踩了一个特别蠢的坑我以为只需要把每个高斯中心当成普通点云用一个MVP矩阵投到屏幕上再叠上一个固定大小的圆斑当模糊效果就行。结果跑出来的图全是边缘发亮的空心圈和奇怪的条纹怎么看怎么不对。后来才意识到3DGS里的“高斯”不是一个点而是一个三维概率分布形状上是一个带方向、带缩放的三维椭球你要做的是把这个椭球“投影”到屏幕上得到一个二维高斯斑块然后再去参与alpha混合。这个从三维椭球到二维像素斑块的变换就是标题里那个投影变换矩阵真正要解决的事情。这篇文章我会按照自己的推导习惯把3DGS的投影链路完整走一遍。适合两类人看一类是正在复现3DGS、对着diff-gaussian-rasterization源码挠头的初学者另一类是想把3DGS改成自己的渲染器、或者换成非针孔相机模型全景、鱼眼的进阶玩家。我尽量把每个矩阵为什么长这样、每一步为什么非得这么做讲清楚而不是让读者背下来完事。1. 射影变换为什么不能直接套在高斯上问题就出在“分布”两个字1.1 高斯不是一个点从概率分布视角看渲染3DGS的场景表示说穿了就是一堆三维高斯分布。每个高斯由三样东西定义中心位置 μ一个三维坐标、协方差矩阵 Σ决定椭球的朝向和伸缩、以及一个不透明度 α 和一组颜色系数通常是球谐系数。单个三维高斯的概率密度函数可以写成G(x) exp(−0.5 * (x − μ)^T · Σ⁻¹ · (x − μ))注意这个“椭球”只是它的等值面形状真正参与渲染时高斯是对周围一片空间都有贡献的连续分布不是只有几何表面。渲染一张图本质上是把这一堆连续分布“画”到每个像素上看每个像素被多少高斯覆盖、每个高斯贡献多少颜色和不透明度。这里就出现了一个根上的矛盾点云的投影是把三维点 p(x,y,z) 经过旋转平移、透视除法、内参缩放变成一个二维像素坐标 (u,v)这是一条确定性的路径。但高斯的投影不能只投“中心点”因为高斯是一个有形状的分布。投影之后分布的形状也必须跟着变——中心投到哪、椭球投成什么样、长轴方向朝向哪里、在屏幕上覆盖多大范围这些全都要算。如果你只投中心点然后硬套一个固定半径的圆斑那等于抛弃了协方差矩阵里存的所有方向信息渲染结果自然就是那副空心圈的样子。1.2 非线性的射影变换为什么不能一步到位用 MVP 矩阵常规渲染里MVP矩阵把三维点从世界坐标变换到裁剪坐标再经过透视除法得到归一化设备坐标。这个流程对“点”是完全够用的因为点的坐标变换是线性的至少到裁剪坐标之前是线性的。但高斯分布不一样。高斯分布经过一个线性变换 y A·x 之后均值和协方差有现成的封闭形式μ_y A·μ_xΣ_y A·Σ_x·A^T注意协方差变换是 AΣA^T 而不是 AΣ这是线性代数里外积期望的必然结果。麻烦在于针孔相机的透视投影是 u fx·x/z, v fy·y/z这是一个非线性函数。非线性变换下高斯分布的形状不一定还是高斯。理论上有无数种可能但我们希望在“单个椭球局部范围内”用一阶泰勒展开把这个非线性变换近似成线性变换也就是用雅可比矩阵 J 代替 A于是Σ_2d ≈ J·Σ_cam·J^T这个近似叫作“局部线性化”。单个高斯椭球在相机坐标系里的空间范围其实很小所以一阶近似的误差在视觉上通常可以忽略。3DGS原文也是这么做的先做刚体变换世界到相机再做一个“近似线性”的射影变换最终在屏幕空间得到一个二维高斯协方差矩阵。整个推导的难点就不在“知道要用 JΣJ^T”而在“J 的每一项到底怎么求”以及“世界坐标、相机坐标、NDC、像素坐标这几套坐标系的矩阵到底怎么对位”。下面先把坐标系和符号钉死再动手推雅可比。2. 坐标系统一与符号约定把 T_w2c、K、雅可比的位置摆正2.1 四个坐标空间与变换关系3DGS的投影链路掰开揉碎就是四个坐标空间世界坐标系World、相机坐标系Camera、归一化坐标NDC以及像素坐标系Pixel。坐标空间输入变换输出World高斯中心 μ_w、协方差 Σ_w视图矩阵 T_w2c [R | t]CameraCameraμ_cam、Σ_cam透视除法u_ndc x/z, v_ndc y/zNDC归一化平面NDCu_ndc, v_ndc内参矩阵 Ku_pix fx·u_ndc cxPixelPixelu_pix, v_pix配合 cov2d 参与光栅化采样屏幕图像这里有一点很多人第一次看源码会懵3DGS原版并没有像传统渲染一样构造一个完整的投影矩阵而是把链路拆成了“刚体变换 雅可比近似 focal缩放”三段。为什么因为高斯的协方差矩阵是3×3的对称阵你要算的是它在屏幕空间被“挤成”什么样的2×2协方差直接去构造一个4×4的投影矩阵反而绕。更关键的是透视除法对协方差的变换不是线性的必须用雅可比来处理所以与其拼一个大矩阵不如把每个环节的矩阵拆开每一步都保持物理意义清晰。2.2 视图矩阵世界到相机旋转参与协方差变换平移只影响中心视图矩阵 T_w2c 是4×4的齐次矩阵包含旋转 R 和平移 t。高斯中心从世界坐标变到相机坐标很简单μ_cam R·μ_w t但协方差矩阵的变换里平移部分不参与。因为协方差描述的是分布相对中心的离散程度整体平移不影响形状。所以Σ_cam R·Σ_w·R^T很多人在这一步开始犯迷糊R 是3×3Σ_w 是3×3为什么不是 R·Σ_w回到前面那个 AΣA^T 推导因为协方差是“外积的期望”必须右侧也乘一个转置。如果你之前没有推导过这个公式建议自己拿一个一维分布算一遍如果所有点都乘以2方差会变成原来的4倍2²这就是两侧各乘一个 R 的直觉含义。2.3 内参矩阵相机到像素焦距负责缩放光心只做平移针孔相机模型里三维点在相机系下的坐标 (x,y,z) 投影到像素平面是u fx·x/z cx v fy·y/z cy其中 fx, fy 是焦距以像素为单位cx, cy 是光心坐标。在3DGS里畸变通常被忽略或者在做COLMAP标定时就先去畸变了所以内参矩阵 K 就是一个经典的上三角矩阵K [[fx, 0, cx], [0, fy, cy], [0, 0, 1]]注意cx 和 cy 在协方差的变换里会直接消失。为什么因为它们只是对坐标做了一个平移平移不改变分布的方差形状。所以在算 2D 协方差时你完全可以把 K 简化为只保留对角线上的 fx 和 fycx、cy 等最后算像素坐标中心时再加上去就行。把坐标系和符号搞整齐后面手推雅可比就不会乱套了。3. 手撕投影雅可比从透视除法开始逐项求偏导3.1 雅可比矩阵的几何直觉非线性函数的局部线性化如果你对“雅可比矩阵”这个名字有点发怵先忘掉教科书定义用一个生活类比理解。想象你站在一个山坡上脚下的地形是起伏不平的但你周围一米范围内基本都是斜平面。这个“斜平面”的坡度就是雅可比矩阵。你在这个小区域里走任何方向高度变化都可以用坡度乘以水平位移量来估算误差很小。同理透视投影在整个三维空间里是弯曲的x/z 是双曲型函数但在单个高斯椭球的局部范围内我们可以把它近似成一个线性变换这个线性变换的系数矩阵就是雅可比。所以在3DGS里我们要做的事是在相机坐标系里给定一个高斯中心把这个中心附近的射影变换做一阶展开得到一个2×3的雅可比矩阵 J然后用 J·Σ_cam·J^T 得到屏幕空间的2×2协方差矩阵。3.2 针孔相机投影的雅可比矩阵逐项推导设高斯中心在相机坐标系的位置是 (x, y, z)那么它投影到像素坐标的函数是u fx·x / z v fy·y / z这里我暂时忽略 cx, cy因为它们对形状没影响。接下来对四个可能的偏导方向逐一求导。先看 u 对 x、y、z 的偏导∂u/∂x fx / z∂u/∂y 0∂u/∂z − fx·x / z²再看 v 对 x、y、z 的偏导∂v/∂x 0∂v/∂y fy / z∂v/∂z − fy·y / z²把这些偏导按行排成矩阵就得到 2×3 的雅可比矩阵J [[fx/z, 0, -fx·x/z²], [ 0, fy/z, -fy·y/z²]]这个矩阵就是“从相机坐标系到屏幕像素坐标系”的投影雅可比。注意它的每一项都依赖高斯中心的实际位置 (x,y,z)离相机越近z 越小J 的数值越大说明近处的高斯在屏幕上的形变更剧烈离光轴越远x/y 越大∂u/∂z 这一项越大说明远离图像中心的高斯受透视“挤压”更明显。3.3 一个数值小例子让抽象代数落地光给公式还是有点悬我拿一个最简单的例子验证一下。假设世界系里有一个各向同性的球高斯Σ_w I单位阵R 也是单位阵所以 Σ_cam I。高斯中心在相机坐标 (2, 1, 5)焦距 fx fy 500。代入雅可比公式J [[500/5, 0, -500·2/25], [ 0, 500/5, -500·1/25]] [[100, 0, -40], [ 0, 100, -20]]然后Σ_2d J·I·J^T J·J^T [[100²40², 40·20], [40·20, 100²20²]] [[11600, 800], [ 800, 10400]]这告诉我们两件事。第一即使原来是一个各向同性的球高斯投影到屏幕上也是一个有方向的椭圆长轴和短轴都不是原来的轴向。第二这个 2D 协方差矩阵的对角线是 11600 和 10400标准差大约是 sqrt(11600)≈108 像素。3DGS光栅化时一般取 3σ 作为有效半径所以这个高斯在屏幕上会覆盖大约 300×300 像素的范围。一个球高斯都能占到这么大面积这就是为什么3DGS的每个splat都不是“一个点”而是一个实打实的椭球脚印。3.4 为什么原版代码里的雅可比是3×3如果你翻过 diff-gaussian-rasterization 的源码会发现它在 CUDA kernel 里构造的雅可比矩阵是3×3的J [[1/z, 0, -x/z²], [0, 1/z, -y/z²], [0, 0, 0 ]]这里跟我的推导有两个差异一是没乘 fx/fy二是多了个第三行。先解释第三行。3DGS在处理协方差变换时为了矩阵运算方便把2×3的雅可比补成3×3第三行全部为0。这个第三行对应的是齐次坐标的 w 分量透视除法之后 w 恒等于1所以 w 对 x,y,z 的偏导全是0。算完 Σ_2d J·Σ_cam·J^T 之后只需要取结果矩阵的左上2×2子矩阵就行第三行第三列不影响前两行的结果。至于焦距原版代码故意不把它乘进雅可比里。它先把相机坐标下的高斯用雅可比近似投影到归一化坐标未乘 fx/fy得到一个“归一化平面上的2×2协方差”之后在光栅化时通过 CUDA 的视野参数 FovX、FovY 或 focal 转成像素坐标。这样做的好处是雅可比和相机内参解耦——如果你想换相机只需要改后面的 focal 参数雅可比本身不用动。而我们自己写简化版时直接乘进 fx/fy 更直观数学上完全等价。4. 一次跑通从3D协方差到屏幕椭圆的完整 PyTorch 实现4.1 数据表示用四元数缩放构造正定协方差在代码里我们不会直接优化一个3×3的协方差矩阵因为那是自找麻烦。高斯协方差必须是半正定矩阵直接回归9个元素很容易在训练中途跑出一个不正定的矩阵渲染时出现“负方差”之类的鬼畜现象。3DGS的解法是参数化用四元数 q 表示旋转用三通道向量 s 表示缩放然后构造旋转矩阵 R_gauss注意这是高斯自身的旋转不是相机视图矩阵再构造缩放对角阵 S diag(exp(s))最终3D世界系协方差为Σ_w R_gauss·S·S^T·R_gauss^T R_gauss·S²·R_gauss^T其中对 s 做 exp 是为了保证缩放为正。这样构造出来的 Σ_w 天然是对称半正定的怎么优化都不会越界。这个技巧你在任何一本“高斯过程”或者“协方差参数化”的资料里都能看到但在3DGS里尤其重要因为每个高斯都要经历两次“旋转相似变换”世界到相机、相机到屏幕中间任何一步产生不正定矩阵画面就崩了。4.2 一次前向的完整 PyTorch 片段下面这段代码浓缩了前面所有推导。输入是一批高斯的世界系参数中心、四元数、缩放、视图矩阵和内参矩阵输出是每个高斯在相机系下的中心和屏幕空间的2×2协方差矩阵。import torch def quat_to_rotmat(q): 四元数 [w, x, y, z] - 旋转矩阵 [N, 3, 3] N q.shape[0] w, x, y, z q[:, 0], q[:, 1], q[:, 2], q[:, 3] R torch.zeros(N, 3, 3, dtypeq.dtype, deviceq.device) R[:, 0, 0] 1 - 2*(y*y z*z) R[:, 0, 1] 2*(x*y - w*z) R[:, 0, 2] 2*(x*z w*y) R[:, 1, 0] 2*(x*y w*z) R[:, 1, 1] 1 - 2*(x*x z*z) R[:, 1, 2] 2*(y*z - w*x) R[:, 2, 0] 2*(x*z - w*y) R[:, 2, 1] 2*(y*z w*x) R[:, 2, 2] 1 - 2*(x*x y*y) return R def build_cov2d(means, quats, scales, viewmatrix, K): # means: [N, 3] 世界坐标 # quats: [N, 4] 四元数 # scales: [N, 3] log缩放 # viewmatrix: [4, 4] 世界-相机 # K: [3, 3] 相机内参 N means.shape[0] R_w2c viewmatrix[:3, :3] t_w2c viewmatrix[:3, 3] # 1. 均值转到相机系 means_cam means R_w2c.T t_w2c # [N, 3] # 2. 构造世界系协方差 rot_gauss quat_to_rotmat(quats) # [N, 3, 3] scale_vec torch.exp(scales) # [N, 3] S torch.diag_embed(scale_vec) # [N, 3, 3] cov_world rot_gauss S.pow(2) rot_gauss.transpose(-1, -2) # 3. 转到相机系只乘旋转 cov_cam R_w2c[None] cov_world R_w2c[None].T # [N, 3, 3] # 4. 构造像素尺度的雅可比 J: [N, 2, 3] x, y, z means_cam[:, 0], means_cam[:, 1], means_cam[:, 2] fx, fy K[0, 0], K[1, 1] J torch.zeros(N, 2, 3, dtypemeans.dtype, devicemeans.device) J[:, 0, 0] fx / z J[:, 0, 2] -fx * x / (z * z) J[:, 1, 1] fy / z J[:, 1, 2] -fy * y / (z * z) # 5. 屏幕空间协方差 cov2d J cov_cam J.transpose(-1, -2) # [N, 2, 2] # 6. 数值保护防止cov2d奇异 eye2 torch.eye(2, dtypemeans.dtype, devicemeans.device) cov2d cov2d 0.1 * eye2 return means_cam, cov2d这段代码里有一个我强烈建议你保留的细节最后给 cov2d 加了一个 0.1×I。原因是当 z 很大或者高斯很扁时cov2d 可能接近奇异矩阵后面光栅化时求逆会爆炸加一个小的正则项让数值稳定。3DGS原版代码里没有这个操作但实际复现时如果遇到画面闪烁或梯度异常可以先从这个角度排查。4.3 和原版 diff-gaussian-rasterization 的差异说明原版代码为了速度是在 CUDA kernel 里逐高斯算 cov2d而且用一个技巧把 J 和 R_w2c 合并成一个2×3的矩阵 M J·R_w2c这样直接从世界系协方差跳到屏幕系协方差少一次矩阵乘法Σ_2d J·(R_w2c·Σ_w·R_w2c^T)·J^T (J·R_w2c)·Σ_w·(J·R_w2c)^T但这种合并写法可读性差和上面分步写的数学结果完全一样。我建议初学者用分步写等读懂了再去看源码里的合并写法。另外真正的光栅化器里还做了大量和投影无关的事情把3D高斯投影成二维椭圆之后要计算椭圆覆盖的 bounding box、在该区域内逐像素采样、计算不透明度、做深度排序和 alpha blending。这些环节和投影矩阵是解耦的所以我在这里先不展开第五部分再说。5. 光栅化落地时绕不开的矩阵细节与验证方法5.1 从2D高斯到屏幕上的椭圆如何确定覆盖范围拿到每个高斯的 cov2d 之后光栅化的第一步是判断“这个高斯在图像上到底占了哪些像素”。理论上二维高斯在整个平面上都有贡献但3σ之外的部分影响已经很小实践中会直接裁剪掉。对 cov2d 做特征值分解cov2d V·diag(λ1, λ2)·V^T两个特征值 λ1, λ2 的平方根就是椭圆两个轴向的标准差V 的列向量是轴向。于是椭圆长半轴近似为 3·sqrt(max(λ1, λ2))短半轴为 3·sqrt(min(λ1, λ2))。光栅化时就围绕高斯中心按这个椭圆的范围框出一个矩形区域只在这个区域里逐像素计算贡献。这种“按需计算”是3DGS实时性的关键——一幅1080p图几百万像素但每个高斯的 footprint 通常只有几十到几千像素省掉无关计算才能跑到实时帧率。这里有一个我在实现时踩过的细节坑特征值可能包含微小的负数导致 sqrt 出错。原因还是数值精度。所以在分解之前最好对 cov2d 做一次 clamp确保对角线非负或者直接执行 cov2d (cov2d cov2d.T) / 2 强制对称再加一点正则。5.2 深度排序与 alpha blending 的光栅化流程每个2D高斯在像素 p 处的贡献通常写成一个带不透明度的二维高斯衰减contribution α_i · exp(−0.5 · Δp^T · cov2d_i⁻¹ · Δp)其中 Δp 是像素坐标与高斯中心投影坐标的差。然后所有覆盖该像素的高斯按深度从近到远排序做 alpha 合成C Σ_i (T_i · contribution_i · c_i)其中 T_i 是前面所有高斯的透射率乘积即 T_i Π_{ji}(1 − contribution_j)。这个排序用的深度就是高斯中心在相机坐标系下的 z 值。看到这里你应该明白投影矩阵的推导在高斯中心排序时又用了一次——means_cam 的 z 坐标直接来自第4节的代码它和 cov2d 是同一个变换链路的产物。所以如果你把视图矩阵的 R 或 t 写错不仅椭圆形状错排序也错最后画面会同时出现“模糊错误”和“遮挡错误”非常难排查。5.3 梯度回传与 gradcheck验证投影矩阵实现的第一道防线如果你只是用现成的 diff-gaussian-rasterization梯度那部分已经被CUDA代码封装好了。但如果你想自己写渲染器或者在原版基础上改投影方式我强烈建议你写一个单元测试用 torch.autograd.gradcheck 验证你的解析梯度。gradcheck 的原理很简单它对输入加一个极小的扰动用数值差分估计梯度再和你代码里自动求导得到的梯度对比。只要投影链路里的矩阵变换有错误这里立刻暴露。一个小小的经验不要直接用整个 cov2d 做 gradcheck那个输出维度太高最好挑选一两个标量比如 cov2d[0,0,0]作为输出。我在最开始写自己的光栅化器时就靠 gradcheck 抓到了两个转置错误一个在视图矩阵乘 R_w2c 时方向反了另一个在 cov2d 的 J 矩阵里把 x/z² 的正负号写错。这些错误看一眼公式都对真的跑起来才露馅。6. 复现3DGS时的常见翻车点与参数心得6.1 环境与编译Ubuntu 20.04 下的常见问题很多人在环境这一关就卡住了。原版代码的依赖比较老在 Ubuntu 20.04 上经常会遇到几个老问题子模块没拉全必须用git clone --recursive否则 submodules 目录是空的编译直接失败。CUDA 版本不匹配原版环境用的是 PyTorch 1.13 CUDA 11.7如果你用更新的 PyTorch 2.x编译 diff-gaussian-rasterization 时可能报一堆稀奇古怪的错误。优先按 environment.yml 的版本装。GCC 版本太高新版 GCC 对旧 CUDA 代码兼容性不好如果编译报错可以试试用 GCC 9。我自己复现时的经验是别在宿主机上硬刚直接用官方 environment.yml 建一个 conda 环境然后按 README 顺序安装git clone https://github.com/graphdeco-inria/gaussian-splatting --recursive cd gaussian-splatting conda env create --file environment.yml conda activate gaussian_splatting pip install submodules/diff-gaussian-rasterization pip install submodules/simple-knn如果最后还是报编译错误优先去看diff-gaussian-rasterization的 README注意它要求的 PyTorch 和 CUDA 组合。很多坑其实不是代码问题是环境组合问题。6.2 显存与显卡4060 到底能不能跑3DGS这个问题在社区里被问过很多次。直接说结论4060 8GB 跑小场景完全够用但别碰大户外场景。我实测下来一个几百帧的室内场景1024×1024分辨率8GB显存可以正常训练但如果换成一两千帧的大场景或者分辨率拉到 1600×1600训练中爆显存几乎是一定的。主要的显存消耗不在前向而在反向传播的梯度存储——每个高斯在每个可见像素上的梯度都会被记录量大得惊人。如果显存不够能做几个事情降训练分辨率比如从原始分辨率降到1/2减少迭代次数3DGS不是非要跑满 30000 次迭代很多场景 15000 次已经能出效果关闭一些额外功能比如球谐阶数降到0也能省一点训练完成后导出模型推理阶段用同一个文件跑实时渲染显存占用会小很多。总之4060 是能入门3DGS的好卡但别把它当成能通吃所有数据集的卡。当作学习验证平台完全够。6.3 训练效果不行的常见原因先把矛头对准投影链路如果你跑出来的训练效果很差先别急着调学习率先确认投影链路本身没错。我排过很多此类问题最常见的几个原因相机内参对不上COLMAP 输出的 fx/fy 是以像素为单位的但如果是自己构造的数据集焦距可能被归一化到 [0,1]直接喂进去所有投影全错。判断方法就是拿特征点坐标反投影一下看投回图像的位置和原始帧能不能对上。高斯初始化太差3DGS 依赖 SfM 点云初始化。点太少了或者场景里大量区域没有初始点训练很久都长不出高斯。densification 被关掉3DGS 的自适应控制是核心每迭代100次进行分裂/克隆。如果设置不对高斯数量不增长场景细节死活出不来。视图矩阵方向写反如果你是自己从 COLMAP 读的位姿注意 COLMAP 的相机坐标系和 3DGS 代码默认的坐标系可能不同R 和 t 需要转换。这也是最常见的“图全花了”的原因。这里面最烦人的是内参问题因为外观看不出来。我给自己定的排错流程是先用原版数据集比如官方提供的小场景跑通再换自己的数据。如果原版数据正常、自己数据崩了那一定是数据预处理的问题跟3DGS代码无关。6.4 一个建议用暴力采样验证投影矩阵最后分享一个我每次修改投影相关代码都会做的“土办法”。我自己在实现阶段经常会写一个独立的测试函数在相机前方一定范围内均匀撒几千个三维点模拟某个高斯的内部样本对这些点做同样的世界到相机变换、雅可比近似、投影到像素平面然后统计这些样本在像素空间的协方差矩阵。这个统计出来的协方差应当和直接调用build_cov2d算出来的 cov2d 非常接近。如果两者差很多说明雅可比矩阵或视图矩阵有错。这个方法比对着论文推公式管用得多它能让你把“数学推导的笔误”和“代码实现的错误”区分开。我印象最深的一次是排查一个很难发现的 bug我在视图矩阵里用的 R 是从 COLMAP 位姿直接拆出来的但 COLMAP 的旋转表示是 world-to-camera而我代码里另一处用到了 camera-to-world一个转置之差让所有高斯的椭球方向全部错乱。用暴力采样一验误差立刻大得离谱才定位到问题。3DGS 的整个投影链路说穿了就是一堆矩阵乘法和一阶泰勒展开最怕的不是公式难而是各种约定不一致。把坐标系、前后乘顺序、转置都钉死后面写任何变体都会顺很多。