ARTICLE DETAIL

资讯详情

深耕郑州网站建设与运营推广的一线实战洞察。

CT重建(三) | R²-Gaussian 知识点全解:当 3D Gaussian Splatting 遇上 CT 重建

CT重建(三) | R²-Gaussian 知识点全解:当 3D Gaussian Splatting 遇上 CT 重建 R²-Gaussian 知识点全解项目R²-Gaussian: Rectifying Radiative Gaussian Splatting for Tomographic ReconstructionNeurIPS 2024仓库r2_gaussian 论文arXiv 2405.20693本文档整理自代码通读 论文理解用于系统掌握该项目的原理、实现与工程要点。目录项目是什么整体流程总览背景一CT 重建与采集几何背景二3D Gaussian Splatting 怎么做核心思想把体积表示成 3D 高斯R²-Gaussian 的核心创新关键沿射线对密度做线积分公式对照3DGS vs R²-Gaussian初始化为什么用 FDK 而不是 SfM损失函数与正则化自适应密度控制相似点与差异点总览数据格式与坐标系代码地图工程与运行术语表1. 项目是什么R²-Gaussian是把3D Gaussian Splatting3DGS改造用于CT计算机断层成像重建的官方实现核心是用一组 3D 高斯椭球表示被测物体的密度体积并用可微的 X 射线渲染直接从投影数据优化这些高斯。关键点任务CT 重建 从多角度 X 射线投影正弦图反推 3D 密度分布是一个病态逆问题。方法把体积表示成 3D 高斯渲染 沿射线对密度做线积分符合 Beer–Lambert 物理。特点速度快、直接、可微擅长稀疏视角 / 有限角度场景。2. 整体流程总览投影数据 扫描仪几何 │ ▼ FDK 初始化 → 采样点云 → 3D 高斯表示 (xyz/scale/rot/density) │ ▼ 可微渲染 ┬─ 光栅化渲染 X 射线投影 └─ 体素化积分成 3D 体积供 TV 正则 / 3D 评估 │ ▼ 损失 L1 SSIM TV → 反向传播 │ ▼ 更新高斯参数 自适应密度控制 (clone / split / prune) │ └──────────► 回到渲染循环 30000 次3. 背景一CT 重建与采集几何CT 成像本质是求解逆问题已知多个角度下探测器测到的 X 射线线积分反推物体内部的密度分布。探测器每个像素的值满足Beer–Lambert 定律I I₀ · exp( −∫ ρ(l) dl ) ⟹ p −ln(I/I₀) ∫ ρ(l) dl即一个像素 一条射线穿过的整条路径上密度的累积。两种采集几何项平行束 (parallel beam)锥束 (cone beam)代码mode01射线形态互相平行准直/同步辐射从一个点光源锥形发散探测器线/面理想化面探测器存在放大与倾斜重建算法FBP 有简洁解析式FDK 近似解析代码分支computeCov2D中雅可比近似正交投影透视仿射近似额外加第三行t/l几何参数来自meta_data.jsonDSO源到物体距离、DSD源到探测器距离、探测器像素间距、旋转角度序列等。几何必须与数据匹配否则重建会模糊或错位。4. 背景二3D Gaussian Splatting 怎么做3DGSKerbl 等SIGGRAPH 2023的核心用一堆 3D 各向异性高斯椭球表示场景投影到屏幕、排序、alpha 混合渲染全程可微用图像反传优化高斯。① 场景表示每个高斯有位置μ、协方差Σ形状朝向、不透明度σ、视角相关颜色球谐 SH。② 协方差构造Σ R S Sᵀ RᵀS diag(sₓ, s_y, s_z)控制三个轴长短R用四元数表示朝向。好处天然半正定避免直接优化Σ的约束冲突各向异性可贴合薄壁、细边缘。③ 投影到 2DEWA splattingΣ J W Σ Wᵀ JᵀW是世界→相机旋转J是投影变换的仿射雅可比。对Σ特征分解得到屏幕椭圆主轴。④ 渲染遮挡式 alpha 合成按深度排序后C Σ_i c_i · α_i · Π_{ji}(1 − α_j), α_i σ_i exp(−½ dᵀΣ⁻¹d)前景遮挡背景。⑤ 优化L (1−λ)·L1 λ·D-SSIM。⑥ 自适应密度控制Clone复制小高斯补细节、Split分裂大高斯、Prune剪枝。高斯向需要细节处聚集。⑦ 初始化用 SfM 稀疏点云。为什么快NeRF 沿光线密采样过 MLP秒级3DGS 只在有高斯覆盖的 tile 上做解析投影/混合可达实时。5. 核心思想把体积表示成 3D 高斯体积不是体素网格而是一组3D 高斯椭球。每个高斯 4 个参数位置xyz、尺度scaling、旋转rotation、密度densityr2_gaussian/gaussian/gaussian_model.py。相比体素网格无网格分辨率限制、内存友好、可微、天然支持从少量投影恢复。本项目参数化的特殊点对照原版 3DGS参数3DGSR²-Gaussian幅值不透明度σsigmoid → (0,1)密度ρ₀softplus → [0,∞)尺度exp支持scale_bound相对体积大小的百分比约束颜色球谐 SH视角相关无单通道灰度均值/旋转四元数相同6. R²-Gaussian 的核心创新辐射高斯光栅化radiative splatting把 3DGS 的表面 alpha 合成改造成沿射线对密度的线积分物理匹配 Beer–Lambert。标题 “Rectifying” 即指修正原版 3DGS 不适用的渲染模型。解析射线积分因子 μ论文 Eq.7不是数值采样而是用闭式积分偏置因子由投影后的 3D 协方差解析算出。代码里即computeCov2D返回的第四分量cov.w。高斯体素化voxelization同一套高斯用第二个 CUDA 内核直接积分成 3D 体素体积用于体积域 TV 正则和 3D 指标评估。完整可训练框架FDK 初始化 密度阈值剪枝 尺度归一化约束 3D TV 正则 自适应密度控制适配 CT 稀疏视角/噪声/有限角度。速度与质量稀疏视角下比 NeRF 类SAX-NeRF 等快很多且质量更好支持平行束与锥束。7. 关键沿射线对密度做线积分物理直觉探测器不看见表面点而是记录穿过整条射线所有物质的累积吸收。没有遮挡前面的物质衰减一点后面的继续贡献——多个高斯对一条射线的贡献是线性叠加。单个高斯对一条射线的积分设高斯密度场ρ(x) ρ₀ exp(−½ (x−μ)ᵀ Σ⁻¹ (x−μ))射线x(t) o t·d则∫ ρ(otd) dt ρ₀ · √(2π) · σ∥ · exp(−½ · d⊥² / σ⊥²)沿射线积分后横向仍是 2D 高斯形状但幅值乘了一个与射线方向、协方差有关的解析系数。结论原版 3DGS 在像素中心的采样值 ≠ 该积分值差一个解析系数 μ。直觉又薄又长且顺着射线的高斯采样值小、积分却大必须显式修正。代码实现cuda_rasterizer/forward.cu的renderCUDA//! We compute mu to consider integration.constfloatalphacon_o.w*mu*exp(power);// 密度 × 积分因子 × 2D 高斯形状if(alpha0.00001f)continue;//! Simply add all alphas for X-ray imaging.for(intch0;chCHANNELS;ch)C[ch]alpha;// 所有高斯直接累加无透射率μ 的解析形式把投影后射线坐标系下的 3D 协方差记分量为a…f| a b c | cov | b d e | | c e f | diamond a·d − b² circ a·d·f 2·b·c·e − a·e² − f·b² − d·c² μ² 2π · circ / diamond ⟹ μ √(2π · circ / diamond)其中circ与 3D 协方差行列式相关diamond是它沿射线积掉一维后的 2D 行列式两者之比就是积分后保留的幅值系数。反向传播在backward.cu里用dL_dmu继续回传。8. 公式对照3DGS vs R²-Gaussian原版 3DGS表面型带遮挡α_i σ_i · exp( −½ (x − μ_i)ᵀ Σ⁻¹_i (x − μ_i) ) C(x) Σ_i c_i · α_i · Π_{ji} (1 − α_j)R²-Gaussian穿透型线积分C(x) Σ_i ρ₀,i · μ_i · exp( −½ (x − μ_i)ᵀ Σ⁻¹_i (x − μ_i) )差异一览维度3DGSR²-Gaussian遮挡有透射率Π(1−α)无直接求和幅值σ_iρ₀,i · μ_i多积分因子 μ合成加权遮挡和简单累加颜色每高斯 SH 颜色c_i无单通道通道RGBNUM_CHANNELS 1物理表面覆盖 / 视角合成密度线积分 / CT 投影9. 初始化为什么用 FDK 而不是 SfMSfM3DGS 用的处理普通照片特征匹配 三角测量输出稀疏表面点云 相机位姿。前提有纹理、有可匹配角点、像素对应单一可见物点。为什么 CT 用不了 SfMX 射线是穿透叠加不同深度的结构会叠到同一个像素没有单一物点对应关系无法特征匹配/三角测量。CT 的相机位姿角度、DSO/DSD、探测器几何是已知标定的不需要估计真正缺的是体积本身的初始估计而 SfM 给不了。FDK本项目用的FDK 是锥束 CT 的经典解析重建算法加权 → 斜坡滤波 → 反投影。流程initialize_pcd.pyTIGRE 跑 FDK 得到体积vol密度阈值density_thresh0.05筛出有效体素随机采样n_points50000个有效体素 → 高斯初始μ与ρ。效果高斯直接落在真实结构上好初始化对病态 CT 问题至关重要。本质区别维度SfMFDK输出稀疏表面点云稠密体积场依赖纹理 / 可见性已知扫描几何 投影精度点位置较准有伪影噪声但结构正确用途3DGS 初始骨架R²-Gaussian 初始骨架代码支持--recon_method fdk默认质量好与random调试/消融参数random_density_max。10. 损失函数与正则化train.py中的损失utils/loss_utils.pyL_total L1(render, gt) λ_dssim · (1 − SSIM) λ_tv · TV_3D(vol_pred)L1投影域像素级差异。SSIMλ_dssim 0.25结构相似性保护结构。3D TVλ_tv 0.05tv_vol_size 32随机采 32³ 小块算三维全变分抑制噪声与伪影。为什么需要CT 是病态逆问题稀疏视角、有限角度、噪声仅 2D 监督会过拟合必须加 3D 先验稳定重建。这是Rectifying提升质量的关键之一。11. 自适应密度控制在gaussian_model.py中每densification_interval步执行一次Clone复制梯度大且小的高斯 → 在欠重建区域复制增加细节。Split分裂梯度大且大的高斯 → 拆成两个更小的提升局部表达力。Prune剪枝密度过低density_min_threshold、屏幕尺寸过大、或出界的高斯删除。CT 场景的特殊改造阈值基于世界尺度归一化volume_to_world、densify_scale_threshold因为体积被归一到[-1,1]³。增加密度阈值剪枝表示密度场而非不透明度与bbox 约束高斯不能漂出感兴趣体积。新增max_num_gaussians上限。默认参数参数值iterations30000densify_from_iter500densify_until_iter15000densification_interval100densify_grad_threshold5e-5max_num_gaussians500000lambda_dssim0.25lambda_tv0.05tv_vol_size3212. 相似点与差异点总览继承自 3DGS 的相似点表示3D 各向异性高斯xyz scale rotation 幅值。协方差Σ R S Sᵀ Rᵀ四元数归一化。投影EWA splattingΣ J W Σ Wᵀ Jᵀconic 逆协方差exp(−½ dᵀΣ⁻¹d)形状。光栅化结构tile16×16 每线程一像素 共享内存批量取数 深度排序点表。训练可微渲染 梯度下降 各参数独立学习率与指数衰减。自适应密度控制clone / split / prune统计量xyz_gradient_accum、denom、max_radii2D。工程PyTorch CUDA 扩展、simple-knn、pickle/ply 存取、triansetup/capture/restore。代码同源文件头保留 Inria 版权auxiliary.h注明 “Modified from diff-gaussian-rasterization”。为 CT 改造的差异点渲染物理表面 alpha 合成 →密度线积分无遮挡、多一个 μ。参数语义不透明度 →密度softplus、可 1 尺度约束。去掉颜色/SH单通道。投影几何相机透视 →平行束 / 锥束保留 3D 协方差第三行以算 μ。损失增加3D TV 正则。初始化SfM →FDK 重建采样。数据图像/位姿 →投影 meta_data.json扫描仪配置体积归一到[-1,1]³。新增能力体素化内核3D 评估与正则。13. 数据格式与坐标系采用 NeRF 风格目录也支持 NAF / SAX-NeRF 的*.pickledata/synthetic_dataset/cone_ntrain_50_angle_360/0_chest_cone/ ├── proj_train/ # 训练投影 *.npy ├── proj_test/ # 评估投影 *.npy ├── init_*.npy # 初始化点云 (xyz density) ├── meta_data.json # 扫描仪配置 / 投影参数 └── vol_gt.npy # 真值体积meta_data.json关键字段nVoxel / sVoxel / dVoxel体积格子数与尺寸、nDetector / sDetector / dDetector探测器、DSD / DSO距离、offOrigin / offDetector偏移。坐标系归一化把感兴趣体积缩放到[-1,1]³scene_scale 2 / max(sVoxel)学习率与尺度阈值都随volume_to_world缩放。14. 代码地图r2_gaussian/ ├── arguments/__init__.py # ModelParams / OptimizationParams / PipelineParams所有超参 ├── gaussian/ │ ├── gaussian_model.py # 高斯模型参数、激活、clone/split/prune、存取 │ ├── initialize.py # 从点云npy/ply创建高斯 │ └── render_query.py # render() 光栅化渲染 / query() 体素化体积 ├── dataset/ │ ├── dataset_readers.py # 读 blender 风格与 NAF 格式几何归一化 │ └── cameras.py # Camera 定义与投影矩阵 ├── utils/ │ ├── loss_utils.py # l1 / ssim / tv_3d │ ├── gaussian_utils.py # 旋转/协方差/学习率调度 │ ├── ct_utils.py # TIGRE 几何与重建调用 │ ├── camera_utils.py # 相机/投影 │ ├── image_utils.py # metric_vol / metric_projPSNR/SSIM │ └── plot_utils.py # 可视化可交互显示体积/高斯 └── submodules/ ├── simple-knn/ # 最近邻距离初始化尺度来自 3DGS └── xray-gaussian-rasterization-voxelization/ ├── cuda_rasterizer/ # 投影渲染内核forward/backward ├── cuda_voxelizer/ # 体积体素化内核 └── third_party/glm顶层脚本initialize_pcd.pyFDK / random 初始化点云可选评估重建。train.py训练主循环渲染→损失→反传→自适应控制→保存/日志。test.py加载模型评估投影域与体积域指标。scripts/visualize_scene.py可视化数据场景。15. 工程与运行环境Ubuntu 20.04 RTX 3090 测试PyTorch 2.1.2 CUDA 11.8TIGRE v2.3。安装pip 路线gitclone https://github.com/Ruyi-Zha/r2_gaussian.git--recursiveconda create-nr2_gaussianpython3.9-yconda activate r2_gaussian pipinstalltorch2.1.2cu118torchvision0.16.2cu118 --extra-index-url https://download.pytorch.org/whl/cu118 pipinstall-rrequirements.txt pipinstall-er2_gaussian/submodules/simple-knn pipinstall-er2_gaussian/submodules/xray-gaussian-rasterization-voxelizationwgethttps://github.com/CERN/TIGRE/archive/refs/tags/v2.3.zipunzipv2.3.zip pipinstallTIGRE-2.3/Python --no-build-isolation运行python initialize_pcd.py-sdata_pathpython train.py-sdata_path-moutput_pathpython test.py-sdata_path-moutput_path注意事项必须结合--recursive克隆或git submodule update --init --recursive否则simple-knn、glm等子模块为空CUDA 扩展无法编译。需要 NVIDIA GPU CUDA 环境。数据与预训练模型见 README 中的 Google Drive 链接。16. 术语表术语含义CT计算机断层成像从投影反推内部密度正弦图 (sinogram)各角度投影组成的二维数据FDK锥束 CT 的经典解析重建算法FBP滤波反投影二维/平行束的解析重建DSO / DSD源到物体 / 源到探测器距离Beer–LambertX 射线强度按指数衰减的物理定律线积分沿射线对密度求积分等于投影值3DGS3D Gaussian Splatting原版高斯泼溅EWA splatting把 3D 高斯投影为 2D 椭圆的经典方法协方差 Σ决定高斯椭球的形状与朝向透射率 Talpha 合成中Π(1−α)代表遮挡μ (mu)本项目的解析射线积分因子Eq.7自适应密度控制clone / split / prune 动态增删高斯体素化 (voxelization)把高斯积分成 3D 体素体积TV 正则全变分平滑先验抑制噪声伪影SfM从运动恢复结构输出稀疏点云PSNR / SSIM常用图像/体积重建质量指标写在最后R²-Gaussian 的价值不在于把 3DGS 搬到 CT 上这一句话本身而在于它精准地指出了两者之间的物理不兼容点并给出了一个干净、解析、可微的修正方案。原版 3DGS 的渲染模型建立在表面 遮挡的假设之上光线碰到最近的表面就被挡住后面的东西不再贡献。而 X 射线成像的物理本质恰恰相反——它是穿透 累积的一条射线穿过多少物质就记录多少吸收前面的物质衰减一点后面的继续叠加。这两套逻辑在数学形式上只差一个透射率连乘和积分因子 μ但在物理上却是完全不同的世界。围绕这个核心洞察整个项目串起了一条完整的技术链表示用 3D 各向异性高斯代替体素网格摆脱分辨率限制渲染用解析线积分代替 alpha 合成把 μ 这个采样值 ≠ 积分值的偏差显式补上初始化用 FDK 代替 SfM因为 CT 没有纹理、没有单一可见物点却有已知几何和解析重建正则用 3D TV 压制病态逆问题的噪声与伪影控制把 clone / split / prune 的阈值重新锚定到归一化的世界尺度上评估用同一套高斯再做一次体素化得到体积域的监督与指标。如果只记一句话可以记这句R²-Gaussian 把 3DGS 的遮挡式加权和改写成了穿透式线积分和其余的一切改造——密度参数、μ 因子、FDK 初始化、TV 正则——都是为这个物理正确性服务的。
返回列表