目前所有的代码、数据、高斯点云模型、训练日志均已开源,欢迎大家来使用或是提交issue。
NeRF类重建方法由于需要密集地采集射线,并在射线上采样多个3D点,计算后通过再通过volume rendering渲染得到2D projection。这一过程非常消耗时间。
当前正火的3D Gaussian Splatting (3DGS)因为其高度平行化的渲染方法——Rasterization而有着比NeRF更快的渲染速度。然而3DGS是针对自然光成像设计的。
直接将3DGS用于X光成像会遇到两个问题:
(1)首先,如图2所示,自然光成像主要依靠于光线在物体表面的反射,这使得从不同角度看,物体表面的颜色会有差异。为了拟合这一各向异性的特点,3DGS采用球谐函数(Spherical Harmonics,SH)来模拟自然光的分布。然而在X光成像中,X光穿透物体并衰减,然后落在探测器上成像。X光无论从各个角度穿透同一物质点,其衰减都是一样的。直接使用SH很难拟合X光成像的这一各向同性的特点。
(2)其次,3DGS的初始化需要通过计算Structure-from-Motion(SfM)算法来得到各个视角的相机内外参数以及一个稀疏点云作为起始。这个算法十分耗时,增加了患者和医生的等待时间。
图2 3DGS 自然光成像与 X-Gaussian 进行 X 光成像对比
本文针对上述这些问题展开研究,做出了以下四点贡献:
针对X光新视角合成任务,本文提出首个基于3D Gaussian Splatting的技术框架——X-Gaussian;
设计了一个全新的辐射高斯点云模型(Radiative Gaussian Point Cloud Model),基于该模型,又设计了一个可微的辐射光栅化渲染方法(Differentiable Radiative Rasterization);
针对高斯点云模型,提出了一种初始化方法——Angle-pose Cuboid Uniform Initialization(ACUI),这种初始化方法能够通过X光扫描仪的设备参数和旋转角直接计算出相机内外参数和初始稀疏点云,这使得新方法免于计算 SfM,从而大幅提升训练速度。
X-Gaussian在性能上超过当前最好NeRF方法6.5dB的情况下,推理速度还达到了73倍。同时在传统算法上也验证了,通过新方法合成的新视角X光片能够提升CT重建的图像质量。
在圆形扫描轨迹锥形X光束扫描(circular cone-beam X-ray scanning)场景下研究三维重建问题。空间坐标系的变换关系如图3所示。被扫描物体的中心O为世界坐标系的原点。
扫描仪的中心S为相机坐标系的中心。探测器D的左上角为图像坐标系的原点。整个空间坐标系的变换遵循OpenCV三维视觉的标准。
图3 空间坐标系转换关系示意图
图4 X-Gaussian 算法框架流程图
算法的流程图如图4所示,首先通过图4(a)中的Angle-pose Cuboid Uniform Initialization(ACUI)来计算出X光源(Source)在对应旋转角𝜙下的相机内外参矩阵并计算出初始稀疏点云。然后,针对
X
光各向同性的成像特点设计了辐射高斯点云模型(Radiative Gaussian Point Cloud Model),如图4(b)所示。
针对这一点云模型,团队设计了一个可微的辐射光栅化(Differentiable Radiative Rasterization,DRR)渲染方法,用于三维高斯点云的泼溅渲染,如图4(c)所示。本节先介绍辐射高斯点云模型,然后是可微的辐射光栅化,最后介绍ACUI初始化方法。
辐射高斯点云模型
本小节首先回顾一下3DGS的基本知识。3DGS将一个物体或场景用𝑁𝑝个高斯点云表示如下:
其中的𝐺𝑖表示第𝑖个高斯点云,𝜇𝑖,Σ𝑖,𝛼𝑖分别表示高斯点云的中心位置,协方差,和不透明度。
协方差控制高斯点云椭球的三轴大小,即控制点云的形状。3DGS对每一个高斯点云采用球谐函数来拟合其颜色如下:
其中,𝑐表示颜色,𝑑=(𝜃,𝜙)表示观测视角,𝑘𝑙𝑚表示球谐函数系数,𝑌𝑙𝑚表示球谐函数,将球面上的点映射成一个实数值。然而,如前面的分析,球谐函数并不适合用来模拟各向同性的X光成像。
为此,团队设计了一个辐射强度响应函数(Radiation Intensity Response Function,RIRF)来替代球谐函数。
具体而言,让每一个高斯点云学一个特征向量𝑓其固有的辐射属性,如辐射密度等。然后该点云的辐射强度𝑖便可以被表示为:
其中𝜆表示一组常数。
因此,辐射高斯点云模型𝐺𝑥可以被表示为:
其中𝑓𝑖为可学习参数,表示分配给第𝑖个高斯点云的特征向量。
可微的辐射光栅化方法
基于提出的这个高斯点云模型,团队还设计了一个可微的辐射光栅化方法(Differentiable Radiative Rasterization,DRR)。
整DRR的过程𝐹𝐷𝑅𝑅总结如下:
其中𝐼表示被渲染的图像,𝑀𝑒𝑥𝑡和𝑀𝑖𝑛𝑡分别表示内外参矩阵。接着介绍𝐹𝐷𝑅𝑅的细节。
首先,计算第𝑖个高斯分布上的3D点𝑥的概率如下:
接着,将3D高斯点云从世界坐标系中投影到相机坐标系,进而再投影到图像坐标系上:
其中的𝑡𝑖=(𝑡𝑥,𝑡𝑦,𝑡𝑧)表示相机坐标,𝑢𝑖表示图像坐标。三维的协方差矩阵𝛴𝑖也被对应地投影到相机坐标系上:
其中𝐽𝑖是投影变换(projective transformation)的仿射近似的雅克比矩阵。
𝑊𝑖是viewing transformation。
其中的𝐿𝑆𝐷表示X光扫描仪中X光源(source)和探测器(detector)之间的距离,𝜙表示X光源的旋转角。
然后在图像坐标系下的二维协方差矩阵是直接取𝛴𝑖′的前两行前两列。
将
2D projection分割成互补重叠的
titles。
每一
个三维高斯点云都按照其对应投影所落在的位置分配到对应的
tiles
上。
这些
3D
高斯点云按照与二维探测器平面的距离进行排序。
那么,在2D projection上像素点𝑝上的辐射强度便是混合𝑁个与𝑝重叠的排好序的3D点得到的,如下公式所示:
其中的𝑥𝑗表示落在像素𝑝上的X射线与高斯点云之间的交点,𝑖𝑗表示𝑥𝑗的辐射强度。
模型训练的监督函数是一范数损失与SSIM损失之间的加权和:
其中的𝛾是加权稀疏,可调的超参。
角度位姿立方体均匀初始化
常规的3DGS使用SfM算法来计算每一个视角的相机内外参数以及初始的稀疏点云。
SfM算法的原理是检测不同视角投影之间的特征匹配点。对于X光片这种低对比度的图像来说,SfM的检测精度会降低。
同时运行SfM非常耗时,对几十张图像计算SfM可能需要耗费几个小时。这大大延长了病人和医生的等待时间。
为此,团队设计了角度位姿立方体均匀初始化(Angle-pose Cuboid Uniform Initialization,ACUI)算法。
ACUI首先直接使用X光扫描仪的参数来计算相机的内外参矩阵:
其中的𝑀𝑒𝑥𝑡表示外参矩阵,𝐿𝑆𝑂表示X光源与物体之间的距离。𝑀𝑖𝑛𝑡表示相机内参,𝑊,𝐻表示渲染图像的宽度和高度。