310 lines
13 KiB
Markdown
Raw Permalink Blame History

This file contains ambiguous Unicode characters

This file contains Unicode characters that might be confused with other characters. If you think that this is intentional, you can safely ignore this warning. Use the Escape button to reveal them.

---
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
- picojsonJSON库
- rapidxmlXML库
- **nanoflann**是一个c++11标准库用于构建具有不同拓扑R2R3点云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。