理解并编码如何在 3D 高斯溅射中使用高斯
欢迎来到雲闪世界现在开始讨论高斯!这是每个人最喜欢的分布。如果您刚刚加入我们,我们已经在第1 部分中介绍了如何根据相机的位置获取 3D 点并将其转换为 2D 。在本文中,我们将讨论高斯分布的高斯部分。我们将使用GitHub中的 part_2.ipynb 。
我们在这里要做出的一点小改动是,我们将使用透视投影,它利用的内部矩阵与上一篇文章中所示的不同。然而,当将点投影到 2D 时,这两种方法是等效的,我发现第 1 部分中介绍的第一种方法更容易理解,但是我们改变了我们的方法,以便在 Python 中尽可能多地复制作者的代码。具体来说,我们的“内部”矩阵现在将由此处显示的 OpenGL 投影矩阵给出,乘法顺序现在将是 points @ external.transpose() @ internal。
对于那些想了解这个新的内部矩阵的人(否则请随意跳过本段)r 和 l 是右侧和左侧的裁剪平面,本质上是相对于照片宽度可以看到哪些点,t 和 b 是顶部和底部裁剪平面。N 是近裁剪平面(点将被投影到该平面),f 是远裁剪平面。有关更多信息,及完整代码可联系博主。这还会返回标准化设备坐标(介于 -1 和 1 之间)中的点,然后我们将其投影到像素坐标。撇开题外话,任务保持不变,即取 3D 中的点并投影到 2D 图像平面上。但是,在本教程的这一部分中,我们现在使用高斯而不是点。
def getIntinsicMatrix(
focal_x: torch.Tensor,
focal_y: torch.Tensor,
height: torch.Tensor,
width: torch.Tensor,
znear: torch.Tensor = torch.Tensor([100.0]),
zfar: torch.Tensor = torch.Tensor([0.001]),,
) -> torch.Tensor:
"""
Gets the internal perspective projection matrix
znear: near plane set by user
zfar: far plane set by user
fovX: field of view in x, calculated from the focal length
fovY: field of view in y, calculated from the focal length
"""
fovX = torch.Tensor([2 * math.atan(width / (2 * focal_x))])
fovY = torch.Tensor([2 * math.atan(height / (2 * focal_y))])
tanHalfFovY = math.tan((fovY / 2))
tanHalfFovX = math.tan((fovX / 2))
top = tanHalfFovY * znear
bottom = -top
right = tanHalfFovX * znear
left = -right
P = torch.zeros(4, 4)
z_sign = 1.0
P[0, 0] = 2.0 * znear / (right - left)
P[1, 1] = 2.0 * znear / (top - bottom)
P[0, 2] = (right + left) / (right - left)
P[1, 2] = (top + bottom) / (top - bottom)
P[3, 2] = z_sign
P[2, 2] = z_sign * zfar / (zfar - znear)
P[2, 3] = -(zfar * znear) / (zfar - znear)
return P
3D 高斯图由 x、y 和 z 坐标以及相关的协方差矩阵组成。正如作者所指出的:“一种显而易见的方法是直接优化协方差矩阵 Σ 以获得表示辐射场的 3D 高斯。然而,协方差矩阵只有当它们是半正定的时候才具有物理意义。对于我们所有参数的优化,我们使用梯度下降法,这种梯度下降法不容易被限制来产生这样的有效矩阵,而更新步骤和梯度很容易产生无效的协方差矩阵。”¹
因此,作者使用协方差矩阵的分解,该分解将始终产生正半定协方差矩阵。具体来说,他们使用 3 个“尺度”参数和 4 个四元数,将其转换为 3x3 旋转矩阵 (R)。然后,协方差矩阵由下式给出
请注意,在转换为旋转矩阵之前,必须对四元数向量进行归一化,才能获得有效的旋转矩阵。因此,在我们的实现中,高斯点由以下参数组成:坐标(3x1 向量)、四元数(4x1 向量)、比例(3x1 向量)以及与不透明度(splat 的透明度)相关的最终浮点值。现在我们需要做的就是优化这 11 个参数来获得我们的场景 — 很简单!
但实际上情况比这要复杂一些。如果你还记得高中数学,高斯在特定点的强度由以下公式给出:
然而,我们关心的是 2D 中 3D 高斯的强度,即在图像平面中。但你可能会说,我们知道如何将点投影到 2D!尽管如此,我们还没有讨论将协方差矩阵投影到 2D,所以如果我们还没有找到 2D 协方差矩阵,我们就不可能找到 2D 协方差矩阵的逆。
现在,这是最有趣的部分(取决于你如何看待它)。3D 高斯分布作者的一篇论文参考文献 EWA Splatting 精确地展示了如何将 3D 协方差矩阵投影到 2D。² 但是,这需要了解雅可比仿射变换矩阵,我们将在下面计算它。我发现代码在理解一个困难的概念时最有帮助,因此我在下面提供了一些代码,以举例说明如何从 3D 协方差矩阵转换为 2D。
def compute_2d_covariance(
points: torch.Tensor,
external_matrix: torch.Tensor,
covariance_3d: torch.Tensor,
tan_fovY: torch.Tensor,
tan_fovX: torch.Tensor,
focal_x: torch.Tensor,
focal_y: torch.Tensor,
) -> torch.Tensor:
"""
Compute the 2D covariance matrix for each gaussian
"""
points = torch.cat(
[points, torch.ones(points.shape[0], 1, device=points.device)], dim=1
)
points_transformed = (points @ external_matrix)[:, :3]
limx = 1.3 * tan_fovX
limy = 1.3 * tan_fovY
x = points_transformed[:, 0] / points_transformed[:, 2]
y = points_transformed[:, 1] / points_transformed[:, 2]
z = points_transformed[:, 2]
x = torch.clamp(x, -limx, limx) * z
y = torch.clamp(y, -limy, limy) * z
J = torch.zeros((points_transformed.shape[0], 3, 3), device=covariance_3d.device)
J[:, 0, 0] = focal_x / z
J[:, 0, 2] = -(focal_x * x) / (z**2)
J[:, 1, 1] = focal_y / z
J[:, 1, 2] = -(focal_y * y) / (z**2)
# transpose as originally set up for perspective projection
# so we now transform back
W = external_matrix[:3, :3].T
return (J @ W @ covariance_3d @ W.T @ J.transpose(1, 2))[:, :2, :2]
首先,tan_fovY 和 tan_fovX 是半个视场角的切线。我们使用这些值来限制我们的投影,防止任何狂野的屏幕外投影影响我们的渲染。我们可以从 3D 到 2D 的变换中推导出雅可比矩阵,就像我们在第 1 部分中介绍的初始正向变换一样,但我为您省去了麻烦,并在上面展示了预期的推导。最后,如果您还记得,我们在上面转置了旋转矩阵以适应术语的重新排列,因此我们在返回最终协方差计算之前在倒数第二行转置回来。正如 EWA 溅射论文所述,我们可以忽略第三行和第三列,因为我们只关心 2D 图像平面。您可能会想,为什么我们不能从一开始就这样做?好吧,协方差矩阵参数将根据您从哪个角度观察它而变化,因为在大多数情况下它不会是一个完美的球体!现在我们已经转换到正确的视点,协方差 z 轴信息就没有用了,可以丢弃。
假设我们有 2D 协方差矩阵,我们几乎能够计算出每个高斯对图像中任何随机像素的影响,我们只需要找到逆协方差矩阵。再次回想一下线性代数,要找到 2x2 矩阵的逆,你只需要找到行列式,然后对项进行一些重新排列。以下是一些代码,可帮助你完成该过程。
def compute_inverted_covariance(covariance_2d: torch.Tensor) -> torch.Tensor:
"""
Compute the inverse covariance matrix
For a 2x2 matrix
given as
[[a, b],
[c, d]]
the determinant is ad - bc
To get the inverse matrix reshuffle the terms like so
and multiply by 1/determinant
[[d, -b],
[-c, a]] * (1 / determinant)
"""
determinant = (
covariance_2d[:, 0, 0] * covariance_2d[:, 1, 1]
- covariance_2d[:, 0, 1] * covariance_2d[:, 1, 0]
)
determinant = torch.clamp(determinant, min=1e-3)
inverse_covariance = torch.zeros_like(covariance_2d)
inverse_covariance[:, 0, 0] = covariance_2d[:, 1, 1] / determinant
inverse_covariance[:, 1, 1] = covariance_2d[:, 0, 0] / determinant
inverse_covariance[:, 0, 1] = -covariance_2d[:, 0, 1] / determinant
inverse_covariance[:, 1, 0] = -covariance_2d[:, 1, 0] / determinant
return inverse_covariance
好了,现在我们可以计算图像中每个像素的像素强度了。但是,这样做非常慢而且没有必要。例如,除非协方差矩阵非常大,否则我们真的不需要浪费计算能力来弄清楚 (0,0) 处的 splat 如何影响 (1000, 1000) 处的像素。因此,作者选择计算他们所谓的每个 splat 的“半径”。如下面的代码所示,我们计算沿每个轴的特征值(请记住,特征值显示变化)。然后,我们取最大特征值的平方根以获得标准偏差测量值并将其乘以 3.0,这涵盖了 3 个标准偏差内的 99.7% 的分布。这个半径有助于我们找出 splat 接触的最小和最大 x 和 y 值。渲染时,我们只计算这些范围内像素的 splat 强度,从而节省了大量不必要的计算。很聪明,对吧?
def compute_extent_and_radius(covariance_2d: torch.Tensor):
mid = 0.5 * (covariance_2d[:, 0, 0] + covariance_2d[:, 1, 1])
det = covariance_2d[:, 0, 0] * covariance_2d[:, 1, 1] - covariance_2d[:, 0, 1] ** 2
intermediate_matrix = (mid * mid - det).view(-1, 1)
intermediate_matrix = torch.cat(
[intermediate_matrix, torch.ones_like(intermediate_matrix) * 0.1], dim=1
)
max_values = torch.max(intermediate_matrix, dim=1).values
lambda1 = mid + torch.sqrt(max_values)
lambda2 = mid - torch.sqrt(max_values)
# now we have the eigenvalues, we can calculate the max radius
max_radius = torch.ceil(3.0 * torch.sqrt(torch.max(lambda1, lambda2)))
return max_radius
上述所有步骤都为我们提供了预处理场景,然后可以在渲染步骤中使用。总结一下,我们现在有了 2D 中的点、与这些点相关的颜色、2D 中的协方差、2D 中的逆协方差、排序的深度顺序、每个 splat 的最小 x、最小 y、最大 x、最大 y 值以及相关的不透明度。有了所有这些组件,我们终于可以开始渲染图像了!
感谢关注雲闪世界。(亚马逊aws和谷歌GCP服务协助解决云计算及产业相关解决方案)
订阅频道(https://t.me/awsgoogvps_Host)
TG交流群(t.me/awsgoogvpsHost)