---
title: Untitled
date: 2024-01-01 18:57:57
excerpt:
tags:
rating: ⭐
---
# 前言
- https://github.com/graphdeco-inria/gaussian-splatting/tree/main/gaussian_renderer
基于Sibr渲染器制作的3D高斯查看器。
# 项目结构
- [x] gaussian
- render - sibr_gaussian
- apps - SIBR_gaussianViewer_app
- [x] diff-gaussian-rasterization(CUDA)
# render - sibr_gaussian
- picojson：JSON库
- rapidxml：XML库
- **nanoflann**：是一个c++11标准库，用于构建具有不同拓扑（R2，R3（点云），SO(2)和SO(3)（2D和3D旋转组））的KD树。
## GaussianSurfaceRenderer
>主要用于渲染椭圆体，估计是用于Debug用的。
### GaussianData
- GaussianData()：通过构造函数形参接受CPU端读取的高斯数据，再通过调用glCreateBuffers()、glNamedBufferStorage()创建GL缓存对象并且初始化，并使用GLuint进行记录（index）。
- render：给Shader绑定GL缓存，并且绘制数组实例。
### GaussianSurfaceRenderer
- GaussianSurfaceRenderer()：初始化相关变量。
- 初始化VS/Frag Shader。
- rayOrigin、MVP、alpha_limit、stage变量
- 创建idTexture、colorTexture贴图变量以及过滤器
- 创建fbo对象以及depthBuffer之后调用makeFBO()正式创建FBO
- 创建清屏Shader。
- makeFBX()：创建idTexture、colorTexture、depthBuffer FBO，用于将顶点数据传递到FragShader中。
- process()：整个渲染过程逻辑处理。
1. 清屏。
2. 判断如果分辨率与FBO大小不同，则重新创建FBO。
3. 获取绘制Buffer的Index，调用glDrawBuffers() 绘制colorTexture、idTexture。
4. 开启深度测试关闭Blend模式。
5. 给Shader绑定相关`_paramMVP`、`_paramCamPos`、`_paramLimit`、`_paramStage`变量，并且调用GaussianData.render()进行一次**不透明物体**的渲染。以小方盒的形式绘制点云数据。
6. 调用glDrawBuffers() 绘制colorTexture。
7. 关闭深度测试，开启透明Blend模式。
8. GaussianData.render()进行一次**透明物体**的渲染，融合模式**additive blendnig**。以小方盒的形式绘制点云数据。
9. 开启深度测试，关闭Blend模式。
10. 将结果显示在屏幕上？
## GaussianView
继承自sibr::ViewBase，用与调用渲染器以及显示结果。
### GaussianView
- GaussianView()：
- 初始化_pointbasedrenderer渲染器
- 初始化_copyRenderer渲染器
- 载入图片并且加入debug模式（应该sibr自带的那个多视角图片debug模式）
- 载入*.ply点云文件，函数为loadPly()。
- CUDA相关处理，应该是为了计算3D高斯结果所需的数据。
- 生成GaussianData指针变量gData。
- 初始化3D高斯渲染器对象_gaussianRenderer。
- 创建GL缓存对象imageBuffer。
- CUDA插值操作。
- 绑定3个geomBufferFunc、binningBufferFunc、imgBufferFunc仿函数，用来调整CUDA渲染时的缓存大小（创建或者回收内存空间）
- onRenderIBR()：View的渲染函数。
- Ellipsoids（椭圆体渲染）：使用_gaussianRenderer->process() 进行渲染。(OpenGL)
- Initial Points：`_pointbasedrenderer->process()`渲染点。
- Splats：使用CudaRasterizer::Rasterizer::forward()进行渲染。最后通过_copyRenderer->process()复制回imageBuffer缓存。
- onGUI()：GUI相关逻辑。
CUDA文件位于`SIBR_viewers\extlibs\CudaRasterizer\CudaRasterizer\cuda_rasterizer\rasterizer_impl.cu`以及`forward.cu`，这些为核心逻辑。
## Shader
可以理解为将点云渲染成一个个的椭圆体，每个椭圆体的颜色与点云数据中的颜色相关。
### VertexShader
1. 取得IndexID。
2. 使用IndexID从传入Shader的Buffer中获取的椭圆体中心、alpha、ellipsoidScale、q（四元数rotation），之后将rotation转成3x3矩阵 ellipsoidRotation。
3. 取得当前顶点Index并获得坐标。再乘以椭圆体旋转值并加上椭圆体中心坐标，取得最终的WorldPos（当前顶点的世界坐标）。
4. 使用IndexID从传入Shader的Buffer中取得**辐射照度?辐射强度?** 数据。
5. 将不符合要求的顶点堆到vec4(0,0,0,0)点。
6. 输出顶点数据到FragShader。
### FragShader
1. 计算摄像机=>当前顶点世界坐标的方向向量dir。
2. 调用closestEllipsoidIntersection()，计算与椭圆体的相交的坐标与相交点的法线。
1. 计算椭圆体空间的localRayOrigin与localRayDirection
2. 计算椭圆与直线相交的方程。
3. 计算摄像机朝向的椭圆体的外表面。如果是内表面最终颜色值 * 0.4。
4. 将相交的世界坐标乘以MVP矩阵，得到摄像机View坐标下的的世界坐标。
5. 计算深度缓存。
6. 计算Alpha。
7. 渲染`out_color = vec4(align * colorVert, a);` 也就是colorTexture
8. 渲染`out_id = boxID;`也就是idTexture
# CudaRasterizer
**本人没学过CUDA，以下仅仅是对代码的猜测。**
额外需要了解Tile渲染方式（具体可以看**Tiled-Based Deferred Rendering(TBDR)**) https://zhuanlan.zhihu.com/p/547943994
- 屏幕分成`16 * 16`的tile，每个tile进行单独计算。之后对每个像素进行计算。
- 取得对应tile中Start与End的位置，对已经排序完的高斯点进行计算，求微分。
- 计算当前像素的透明度T
- 2D协方差 => power => alpha。
- 每次循环都进行`float test_T = T * (1 - alpha)`，当test_T极小时（不透明）则停止循环。
- T = test_T。
- 计算当前像素的颜色，也就是计算各个方向接受的辐射照度。
- `for (int ch = 0; ch < CHANNELS; ch++)`
`C[ch] += features[collected_id[j] * CHANNELS + ch] * alpha * T;`
- 计算最终贡献值
- 如果当前像素在范围中则输出
- `final_T[pix_id]`最终透明度。
- `n_contrib[pix_id]`最终贡献值。
- `out_color[ch * H * W + pix_id]`最终颜色。`C[ch] + T * bg_color[ch]`
对屏幕分Tile
![[ScreenSpaceTile.jpg]]
以此减少需要遍历的点云数量。
![[TileRange.jpg|500]]
每个点云相当于空间中当前位置空间的辐射强度分布。
![[GS_radiation.jpg]]
一个像素的渲染会计算这个像素范围内所有的点云的辐射强度、透明度，最后求微分。下图两条横线内相当于一个像素的范围。
![[一个像素需要计算范围内所有电源的辐射强度.png|500]]
## rasterizer_impl.cu
- getHigherMsb()
- checkFrustum()：判断点云是否在视锥内，返回一个bool数组。
- duplicateWithKeys()
- identifyTileRanges()：确定每个Tile的工作起点与终点。
- markVisible()：标记高斯点云是否处于可视状态。
- GeometryState::fromChunk()：计算数据块的指针偏移，并且返回创建的GeometryState结构体对象。
- ImageState::fromChunk()：计算数据块的指针偏移，并且返回创建的ImageState结构体对象。
- BinningState::fromChunk()：计算数据块的指针偏移，并且返回创建的BinningState结构体对象。
- forward()：前向渲染可微分光栅化的高斯。具体见下文。
- backward()：生成优化所需的梯度数据，并传递到forward()。**该项目中目前未被调用**
相关数据结构体定义在rasterizer_impl.h中：
```c++
struct GeometryState
{
size_t scan_size;
float* depths;
char* scanning_space;
bool* clamped;
int* internal_radii;
float2* means2D;
float* cov3D;
float4* conic_opacity;
float* rgb;
uint32_t* point_offsets;
uint32_t* tiles_touched;
static GeometryState fromChunk(char*& chunk, size_t P);
};
struct ImageState
{
uint2* ranges;
uint32_t* n_contrib;
float* accum_alpha;
static ImageState fromChunk(char*& chunk, size_t N);
};
struct BinningState
{
size_t sorting_size;
uint64_t* point_list_keys_unsorted;
uint64_t* point_list_keys;
uint32_t* point_list_unsorted;
uint32_t* point_list;
char* list_sorting_space;
static BinningState fromChunk(char*& chunk, size_t P);
};
```
### forward()
1. 创建相关变量：GeometryState、ImageState、minn、maxx。
2. FORWARD::preprocess()
3. 计算所有tile的高斯点云总量。
4. 根据需要需要渲染的高斯点云总量来调整CUDA buffer大小。
5. 创建BinningState。
6. duplicateWithKeys()
7. getHigherMsb()
8. 对高斯点运行排序。
9. cudaMemset(imgState.ranges, 0, tile_grid.x * tile_grid.y * sizeof(uint2));
10. 调用identifyTileRanges()，确定每个Tile的工作起点与终点。
11. 取得点云颜色数组。
12. FORWARD::render()
## forward.cu
### preprocess()
在光栅化之前，对每个高斯进行初始化处理。
- 只处理在视锥中并且在盒子中的高斯。
- 使用投影矩阵对点云的点进行变换，并进行归一化，赋予给新变量p_proj。
- 计算协方差矩阵cov3D。
- 计算2D屏幕空间的协方差矩阵cov
- Invert covariance
- Compute extent in screen space (by finding eigenvalues of 2D covariance matrix). Use extent to compute a bounding rectangle of screen-space tiles that this Gaussian overlaps with. Quit if rectangle covers 0 tiles.
- 如果没有颜色数据则从球谐函数中计算辐射照度。
- 存储当前数据。
- `depths[idx]`
- `radii[idx]`
- `points_xy_image[idx]`
- `conic_opacity[idx]`
- `tiles_touched[idx]`
```c++
// Invert covariance (EWA algorithm)
float det = (cov.x * cov.z - cov.y * cov.y);
if (det == 0.0f)
return;
float det_inv = 1.f / det;
float3 conic = { cov.z * det_inv, -cov.y * det_inv, cov.x * det_inv };
// Compute extent in screen space (by finding eigenvalues of
// 2D covariance matrix). Use extent to compute a bounding rectangle
// of screen-space tiles that this Gaussian overlaps with. Quit if
// rectangle covers 0 tiles.
float mid = 0.5f * (cov.x + cov.z);
float lambda1 = mid + sqrt(max(0.1f, mid * mid - det));
float lambda2 = mid - sqrt(max(0.1f, mid * mid - det));
float my_radius = ceil(3.f * sqrt(max(lambda1, lambda2)));
float2 point_image = { ndc2Pix(p_proj.x, W), ndc2Pix(p_proj.y, H) };
uint2 rect_min, rect_max;
if (rects == nullptr) // More conservative
{
getRect(point_image, my_radius, rect_min, rect_max, grid);
}
else // Slightly more aggressive, might need a math cleanup
{
const int2 my_rect = { (int)ceil(3.f * sqrt(cov.x)), (int)ceil(3.f * sqrt(cov.z)) };
rects[idx] = my_rect;
getRect(point_image, my_rect, rect_min, rect_max, grid);
}
if ((rect_max.x - rect_min.x) * (rect_max.y - rect_min.y) == 0)
return;
```
### render()
对所有Tile进行并行计算。针对CUDA核心数量创建对应的Block以及对应数据。`int collected_id[BLOCK_SIZE]、float2 collected_xy[BLOCK_SIZE]、float4 collected_conic_opacity[BLOCK_SIZE]`。
递归所有的Block，计算透明度、Color以及贡献值（用于计算平均值）。
```c++
// Iterate over batches until all done or range is complete
for (int i = 0; i < rounds; i++, toDo -= BLOCK_SIZE)
{
// End if entire block votes that it is done rasterizing
int num_done = __syncthreads_count(done);
if (num_done == BLOCK_SIZE)
break;
// Collectively fetch per-Gaussian data from global to shared
int progress = i * BLOCK_SIZE + block.thread_rank();
if (range.x + progress < range.y)
{ int coll_id = point_list[range.x + progress];
collected_id[block.thread_rank()] = coll_id;
collected_xy[block.thread_rank()] = points_xy_image[coll_id];
collected_conic_opacity[block.thread_rank()] = conic_opacity[coll_id];
} block.sync();
// Iterate over current batch
for (int j = 0; !done && j < min(BLOCK_SIZE, toDo); j++)
{ // Keep track of current position in range
contributor++;
// Resample using conic matrix (cf. "Surface
// Splatting" by Zwicker et al., 2001)
float2 xy = collected_xy[j];
float2 d = { xy.x - pixf.x, xy.y - pixf.y };
float4 con_o = collected_conic_opacity[j];
float power = -0.5f * (con_o.x * d.x * d.x + con_o.z * d.y * d.y) - con_o.y * d.x * d.y;
if (power > 0.0f)
continue;
// Eq. (2) from 3D Gaussian splatting paper.
// Obtain alpha by multiplying with Gaussian opacity // and its exponential falloff from mean. // Avoid numerical instabilities (see paper appendix).float alpha = min(0.99f, con_o.w * exp(power));
if (alpha < 1.0f / 255.0f)
continue;
float test_T = T * (1 - alpha);
if (test_T < 0.0001f)
{ done = true;
continue;
}
// Eq. (3) from 3D Gaussian splatting paper.
for (int ch = 0; ch < CHANNELS; ch++)
C[ch] += features[collected_id[j] * CHANNELS + ch] * alpha * T;
T = test_T;
// Keep track of last range entry to update this
// pixel. last_contributor = contributor;
}}
```
```c++
// All threads that treat valid pixel write out their final
// rendering data to the frame and auxiliary buffers.
if (inside)
{
final_T[pix_id] = T;
n_contrib[pix_id] = last_contributor;
for (int ch = 0; ch < CHANNELS; ch++)
out_color[ch * H * W + pix_id] = C[ch] + T * bg_color[ch];
}
```
# apps - SIBR_gaussianViewer_app
调用`gaussianviewer/renderer/GaussianView.hpp`封装的App。