一种基于口腔CBCT重建点云的超分辨成像方法

技术领域

本发明涉及医学图像处理技术领域，具体涉及一种基于锥形束计算机断层扫描(Cone Beam Computred Tomography,CBCT)重建点云的超分辨成像方法。

背景技术

随着消费水平的升级，我国口腔医疗消费意识也快速觉醒。为保证临床牙科诊疗的安全性与高效性，医生在术前需要进行数字化辅助规划，对患者口腔内部状况进行初步评估，有效提高手术的预见性；患者也能因此获得更加舒适顺畅的就诊体验。

口腔数字化诊断已经成为必然趋势，其所依赖的计算机成像技术大致可归纳为三类：传统计算机断层扫描CT、锥形束计算机断层扫描(Cone Beam Computred Tomography,CBCT)和口内扫描仪。传统CT的投影是一维数据，重组得到的三维图像是由重建后的连续多个二维切片堆积而成，其金属伪影较重，影响医生对病灶区域的诊断。口内扫描仪虽然可以较好地实现快速获取口腔三维数字模型，但目前该设备在国内仍以外资品牌为主，价格昂贵，市场普及率和渗透率都比较低。CBCT是牙颌面部位的专属CT机，该项技术使用平板探测器和三维锥形X射线束围绕头部旋转扫描，得到多角度的序列二维投影切片，重建后可直接得到完整的三维牙体模型，方便医生全面了解患者牙齿及其周围的解剖结构。相比于传统CT，CBCT扫描仅用于面部照射，具有数据采集时间短、成像速度快、放射剂量小、体积小等优点，大大减少了由患者身体移动而导致的伪影，对牙和骨的空间分辨率较高。因此，CBCT是目前实用性最强、发展潜力最大的口腔影像设备。

CBCT技术已经在口腔正畸、口腔修复以及口腔种植等临床医学领域中得到广泛应用，但其口腔扫描的精度(150-200微米)相比三维牙体模型的分辨率(400微米左右)仍有一定差距。CBCT的空间分辨率最优可达到2.0lp/mm，在此基础上要获取高精度的牙齿完整形态结构仍有较大难度。因此，提高CBCT三维重建数据的分辨率对口腔医学的未来发展有着重要的临床意义与研究价值。

发明内容

针对上述CBCT技术现有问题，本发明提出了一种基于口腔CBCT重建点云的超分辨成像方法，基于CBCT数据构建稀疏分布的低分辨率(LR)点云数据，再通过LR点云预测密集分布的高分辨率(HR)点云数据，后者包含丰富细节信息和少量噪点的。本发明采用基于对抗残差图卷积深度学习网络(Adversarial Residual Graph Convolution Network,AR-GCN)的方法，考虑到LR和HR点云之间的相似性，将残差连接引入图卷积网络(GCN)，同时引用图对抗损失来减少过分依赖人类的先验知识的干扰，更易捕获HR点云如连续性等其他属性。

本发明的基于口腔CBCT重建点云的超分辨成像方法，包括以下步骤：

一、基于CBCT数据构建LR点云数据：

1)获取点云数据：

采用激光扫描仪提取牙冠三维模型，并将牙冠三维模型转换为高分辨率的激光扫描三维牙模型点云数据；采用锥形束计算机断层扫描仪获取CBCT数据，根据区域增长法从CBCT数据中提取出牙齿完整模型，并将牙齿完整模型转换为CBCT重建三维牙模型点云数据，并与激光扫描三维牙模型点云数据配准，得到三维牙体模型点云数据；

2)利用统计离群消除滤波器对三维牙体模型点云数据进行去噪处理，得到去噪后的三维牙体模型点云数据；

3)对去噪后的三维牙体模型点云数据进行下采样，在降低原始的点云数据量的同时，达到尽可能捕获牙体特征的目的，最终得到含有牙体局部精细细节和粗略全局形状结构的LR点云数据；

二、构建对抗残差图卷积深度学习网络AR-GCN，利用已有的训练集和测试集完成AR-GCN的训练：

1)搭建AR-GCN：

AR-GCN包括残差图卷积生成器G和图判别器D，残差图卷积生成器G包括特征网络、残差图卷积块和上池化块；图判别器D包括特征网络、池化块和残差图卷积块；

采用残差图卷积生成器G对点云数据进行上采样生成HR点云数据，通过图判别器D计算得到适应性图对抗损失函数L

2)构建损失函数并通过损失函数优化AR-GCN：

将倒角距离变式损失函数L

3)采用已有的训练集和测试集训练并评估AR-GCN网络：

从临床扫描得到的牙列CBCT数据中选取大部分作为训练集，占总数据的50％～90％，小部分作为测试集，占总数据的10％～50％，按照步骤一的方法将CBCT数据全部转化三维牙体模型点云数据，并构建LR点云数据，得到预处理后的训练集和测试集；然后将预处理后的训练集用于训练AR-GCN网络，网络训练完成后，由预处理后的测试集评估AR-GCN网络的泛化能力；

三、将步骤一中获得的LR点云数据输入到训练好的AR-GCN中，得到HR点云数据，

从而实现口腔三维点云的超分辨成像。

其中，在步骤一的2)中，利用统计离群消除滤波器对点云数据进行去噪处理，包括以下步骤：

a)采用k近邻算法计算点云数据中每个数据点到其最近k个邻域点的平均距离d，k为选取领域点个数，设置为1～10，；

b)计算平均距离d的期望值dm和标准差s；

c)根据期望值dm和标准差s，计算距离阈值dt的计算公式，

dt＝dm+η×s

其中，η为标准差参数，η取值为0.1～0.3；

d)将每一个特征点的平均距离d与距离阈值dt作比较，若d>dt，则滤除该点，否则保留。

在步骤一的3)中，分别对去噪后的三维牙体模型点云数据进行下采样，包括以下步骤：

(1)点云数据的搜索空间划分

确定三维牙体模型点云数据在空间中的大小，得到X、Y和Z坐标轴的最小值和最大值分别为Min_x、Max_x、Min_y、Max_y、Min_z和Max_z，构建牙齿点云的最大空间范围：[Min_x,Max_x]×[Min_y,Max_y]×[Min_z,Max_z]，对牙齿点云最大范围进行空间分割，得

到牙齿点云的最大包围盒L：

其中，β为调整牙齿最大包围盒的大小因子，取值为0.5～0.8，k为选取领域点的个数，设置为1～10，n

(2)点云数据的法向量估计

将位于最大包围盒L内的每个三维牙体模型点云数据的数据点的法向量作为该数据点的拟合切平面的法向量，计算每个数据点的拟合切平面，完成对每个三维牙体模型点云数据的数据点局部区域的平面拟合；

(3)点云数据的曲率估计

a)利用抛物面拟合算法分别对最大包围盒L内的完整模型表面点云数据进行曲率估计，完整模型表面点云数据的曲面主曲率与抛物面的主曲率相同，抛物面方程为：

Z＝bx

b)使抛物面方程的顶点x

AX＝Z

其中，A为抛物面内邻域点的x轴和y轴的坐标矩阵，X为抛物面方程中的系数向量，Z为抛物面内邻域点的z轴坐标向量，

c)利用奇异值分解(Singular Value Decomposition，SVD)方法求解公式AX＝Z得到系数b、c和d，从而得到抛物面方程；

(4)利用系数b、c和d计算三维牙体模型点云数据的平均曲率H：

H＝b+d

将计算三维牙体模型点云数据的平均曲率与设定的阈值进行对比，这一步骤中根据牙齿的几何形状及其凹凸性，将点云数据中牙冠和牙根分为：峰、谷和脊，其中峰对应牙体牙冠表面的牙尖特征，平均曲率H>0，即局部区域为凸状；谷对应牙体牙冠面上的沟，平均曲率H<0，即局部区域为凹状；脊对应牙体上各种脊线，根据邻点曲率确定凹凸状；对应的峰阈值设为30～50，过滤点云数据中峰的平均曲率低于峰阈值范围的点云数据；

对应的谷阈值区间设为-50～-30，过滤点云数据中谷的平均曲率大于谷阈值范围的点云数据。

在步骤一、3)的(2)中，数据点p的k个邻域点记为N(p)，ξ为数据点p的法向量，数据点的法向量ξ作为该数据点的拟合切平面的法向量，用最小二乘法得到拟合切平面Tp(p)，计算k个邻域点N(p)的协方差矩阵，根据N(p)的协方差矩阵计算数据点p的k个邻域点相应的法向量；对得到的拟合切平面的法向量进行方向统一，保持相邻平面的法向量方向一致，完成对每个三维牙体模型点云数据的数据点局部区域的平面拟合。

在步骤二的1)中，采用残差图卷积生成器G对点云数据进行上采样生成HR点云数据，包括以下步骤：

i.采用特征网络初步提取点云数据的特征：

对点云数据集x中每个大小为1×3的数据点p(p∈x)，点云数据集x位于欧几里得空间，获取数据点的k邻域P，其大小为k×3，然后经由带有最大池化的点卷积操作，将P中k个邻域点相对于数据点p的坐标

ii.采用残差图卷积块进一步提取点云数据的空间特征：

将欧几里得空间中的数据点p的k个邻域点记为N(p)，其坐标定义为x

其中每个残差层由修正线性单元(Rectified Linear Unit,ReLU)激活函数、图卷积G-conv和残差连接三部分组成：ReLu激活函数将AR-GCN深度神经网络中每个神经元的输入映射到输出端；图卷积G-conv利用第l层数据点p的特征值来计算第l+1层数据点的特征值，其中，l为G-conv图卷积的层数，取值范围为1～15，计算公式如下：

其中，f

iii.采用上池化块对点云数据进行上采样：

将上个步骤得到的点云数据x

在步骤二的1)中，通过图判别器D计算得到适应性图对抗损失函数L

i.将残差图卷积生成器G生成的HR点云数据

ii.通过池化块对输入点云数据进行下采样：

对于给定的输入的HR点云数据

其中，N[x

iii.通过残差图卷积块得到特征矩阵：

通过池化块得到的x

iv.将步骤一、2)中得到去噪后的三维牙体模型点云数据作为真实数据y输入到图判别器D中，重复步骤i～iii，得到降维后的特征矩阵D(y)；

v.由D(y)计算得到适应性图对抗损失函数L

其中L

vi.由适应性图对抗损失函数L

在步骤二的2)中，将倒角距离变式损失函数L

i.计算损失函数L

ii.计算联合损失函数Loss为：

其中，λ是平衡两个损失的参数；

iii.优化残差图卷积生成器G：

固定图判别器的参数，最小化Loss，图判别器将误差传到残差图卷积生成器，更新残差图卷积生成器的参数，残差图卷积生成器希望生成以假乱真的数据，即希望判别器向其分配接近1的概率，这样L

iv.优化图判别器：

对于真实数据y，图判别器尝试向其分配一个大于0.5的概率；对于残差图卷积生成器生成的“假数据”

进一步，还包括评估预测的HR点云数据

标准倒角距离(Chamfer Distance,CD)和陆地移动距离(Earth Mover’sDistance,EMD)作为评估预测的HR点云数据

采用标准倒角距离(Chamfer Distance,CD)、陆地移动距离(Earth Mover’sDistance,EMD)和豪斯多夫距离(Hausdorff Distanc,HD)衡量真实点云y和预测点云

采用标准化均匀系数(Normalized Uniformity Coefficient,NUC)度量预测点云

使用F-分数(F-score)作为评估指标：首先将点云超分辨看作一个分类问题，通过检测在特定阈值(0.01～0.05)内y或

本发明的优点：

本发明基于CBCT数据构建稀疏分布的低分辨率(LR)点云，再通过LR点云预测密集分布的高分辨率(HR)点云，后者包含丰富细节信息和少量噪点的；本发明采用基于对抗残差图卷积深度学习网络(Adversarial Residual Graph Convolution Network,AR-GCN)的方法，考虑到LR和HR点云之间的相似性，将残差连接引入图卷积网络(GCN)，同时引用图对抗损失来减少过分依赖人类的先验知识的干扰，更易捕获HR点云如连续性等其他属性。

附图说明

图1为本发明的基于口腔CBCT重建点云的超分辨成像方法的流程图；

图2为本发明的基于口腔CBCT重建点云的超分辨成像方法的下采样的流程图；

图3为根据本发明的基于口腔CBCT重建点云的超分辨成像方法得到的下采样前后的点云图；

图4为本发明的基于口腔CBCT重建点云的超分辨成像方法的AR-GCN深度学习网络结构图；

图5为本发明的基于口腔CBCT重建点云的超分辨成像方法的AR-GCN深度学习网络中的残差图卷积块结构图。

具体实施方式

下面结合附图，通过具体实施例，进一步阐述本发明。

如图1所示，本实施例的基于口腔CBCT重建点云的超分辨成像方法，包括以下步骤：一、基于CBCT数据构建LR点云数据：

1)获取不同来源的点云数据：

采用激光扫描仪提取牙冠三维模型，并将牙冠三维模型转换为高分辨率的激光扫描三维牙模型点云数据；采用锥形束计算机断层扫描仪获取CBCT数据，根据区域增长法从CBCT数据中提取出牙齿完整模型，并将牙齿完整模型转换为CBCT重建三维牙模型点云数据，并与激光扫描三维牙模型点云数据配准，得到三维牙体模型点云数据；

2)利用统计离群消除滤波器分别对三维牙体模型点云数据和牙冠模型点云数据进行去噪处理：

a)采用k近邻算法计算点云数据中每个数据点到其最近k个邻域点的平均距离d，k为选取领域点个数，设置为6；

b)计算平均距离d的期望值dm和标准差s；

c)根据期望值dm和标准差s，计算距离阈值dt的计算公式，

dt＝dm+η×s

其中，η为标准差参数，η为0.2；

d)将每一个特征点的平均距离d与距离阈值dt作比较，若d>dt，则滤除该点，否则保留，得到去噪后的三维牙体模型点云数据和牙冠模型点云数据；

3)分别对去噪后的三维牙体模型点云数据和牙冠模型点云数据进行下采样：

(1)点云数据的搜索空间划分

确定三维牙体模型点云数据和牙冠模型点云数据在空间中的大小，得到X、Y和Z坐标轴的最小值和最大值分别为Min_x、Max_x、Min_y、Max_y、Min_z和Max_z，构建牙齿点云的最大空间范围L：[Min_x,Max_x]×[Min_y,Max_y]×[Min_z,Max_z]，对最大空间范围进行空间分割得到最大包围盒L：

其中，β为调整牙齿最大包围盒的大小因子，k为选取领域点的个数，设置为1～10，n

(2)点云数据的法向量估计

三维扫描获取的初始采样点集只记录了各采样点的空间三维坐标，而不存在任何连接关系，法向量是三维点云数据具有的一个局部特性，在点云的许多处理中是必不可少的信息，求解法向量是处理点云数据的关键步骤。将位于最大包围盒L内的每个三维牙体模型点云数据和牙冠模型点云数据的数据点的法向量作为该数据点的拟合切平面的法向量，数据点p的k个邻域点记为N(p)，ξ为数据点p的法向量，用最小二乘法得到拟合切平面Tp(p)，根据N(p)协方差矩阵计算数据点p的k个邻域点相应的法向量；对得到的拟合切平面的法向量进行方向统一，保持相邻平面的法向量方向一致，完成对每个三维牙体模型点云数据和牙冠模型点云数据的数据点局部区域的平面拟合；

(3)点云数据的曲率估计

a)利用抛物面拟合算法分别对最大包围盒L内的完整模型表面点云数据和牙冠模型表面点云数据进行曲率估计，完整模型表面点云数据和牙冠模型表面点云数据的曲面主曲率与抛物面的主曲率相同，抛物面方程为：

Z＝bx

b)使抛物面方程的顶点p

AX＝Z

其中，

c)利用奇异值分解(Singular Value Decomposition，SVD)方法求解公式AX＝Z得到系数b，c，d，从而得到抛物面方程；

(4)利用b，c和d系数计算三维牙体模型点云数据和牙冠模型点云数据和平均曲率H：

H＝b+d

将计算三维牙体模型点云数据和牙冠模型点云数据的曲率与设定的阈值进行对比，这一步骤中根据牙齿的几何形状及其凹凸性，将牙体模型中牙冠、牙根分为：峰，谷，脊，其中峰对应牙体牙冠表面的牙尖特征，其高斯曲率H＞0(即点为椭圆点)，平均曲率K＞0(即局部区域为凸状)；谷对应牙体牙冠面上的沟，其高斯曲率H＞0(即点为椭圆点)，平均曲率K＜0(即局部区域为凹状)；脊对应牙体上各种脊线，根据邻点曲率确定凹凸状。对应的峰阈值设为50，过滤曲面曲率低于阈值范围的点集，保留曲率高于阈值范围的点集，对应的谷阈值设为-50，过滤曲面曲率大于阈值范围的点集，保留曲率小于阈值范围的点集，降低点云数据量的同时保留三维牙体的轮廓特征；

二、构建对抗残差图卷积深度学习网络AR-GCN，利用已有的训练集和测试集完成AR-GCN的训练：

1)搭建AR-GCN：

残差图卷积生成器G得到高分辨点云数据：

i.AR-GCN采用特征网络初步提取点云数据的特征：

对点云数据集x中每个大小为1×3的数据点p(p∈x)，点云数据集x位于欧几里得空间，获取数据点的k邻域P，其大小为k×3，然后经由带有最大池化的点卷积操作，将P中k个邻域点相对于数据点p的坐标

ii.采用残差图卷积块进一步提取点云数据的空间特征：

将欧几里得空间中的数据点p的k个邻域点记为N(p)，其坐标定义为x

其中每个残差层由ReLu激活函数、图卷积G-conv和残差连接三部分组成：ReLu激活函数将AR-GCN深度神经网络中每个神经元的输入映射到输出端，提高网络的非线性表达能力，防止梯度消失或梯度爆炸现象的发生；残差层中的图卷积G-conv是整个残差图卷积块的核心，它定义在图G＝(v,ε)上，并对数据点p和其邻点q进行不同处理，增强了残差图卷积生成器G的学习能力：具体来说，它可以利用第l层数据点p的特征值来计算第l+1层数据点AR-GCN网络的特征值，计算公式如下：

其中f

iii.采用上池化块对点云数据进行上采样：

将上个步骤得到的点云数据x

图判别器D负责区分高分辨率点云的真伪，包括以下步骤：

i.将残差图卷积生成器G生成的HR点云

ii.通过池化块对输入点云数据进行下采样：

对于给定的输入点云数据

其中N[x

iii.通过残差图卷积块得到特征矩阵：

通过池化块得到的x

iv.将步骤一、2)中得到去噪后的三维牙体模型点云数据作为真实数据y输入到图判别器D中，重复步骤i～iii，得到降维后的特征矩阵D(y)；

i.由D(y)计算得到适应性图对抗损失函数L

其中L

v.由适应性图对抗损失函数L

2)图判别器D构建损失函数并通过损失函数优化AR-GCN：

i.计算得到损失函数L

ii.联合损失函数为：

其中，λ是平衡两个损失的参数；

iii.优化残差图卷积生成器G：

固定图判别器的参数，最小化Loss，图判别器将误差传到残差图卷积生成器，更新残差图卷积生成器的参数，残差图卷积生成器希望生成以假乱真的数据，即希望判别器向其分配接近1的概率，这样L

iv.优化图判别器：

对于真实数据y，图判别器尝试向其分配一个接近1的概率；对于残差图卷积生成器生成的“假数据”

标准倒角距离(Chamfer Distance,CD)和陆地移动距离(Earth Mover’sDistance,EMD)作为评估预测点云

采用标准倒角距离(Chamfer Distance,CD)、陆地移动距离(Earth Mover’sDistance,EMD)和豪斯多夫距离(Hausdorff Distanc,HD)衡量真实点云y和预测点云

采用标准化均匀系数(Normalized Uniformity Coefficient,NUC)度量预测点云

使用F-score作为评估指标：首先将点云超分辨看作一个分类问题，通过检测在特定阈值(本实施例为0.01)内y或

3)采用已有的训练集和测试集训练并评估AR-GCN网络：

将临床扫描得到的牙列CBCT数据分成两部分，其中80％作为训练集，20％作为测试集，按照步骤一的方法将CBCT数据全部转化三维牙体模型点云数据，并构建LR点云数据，得到预处理有的训练集和测试集；然后将预处理后的训练集用于训练AR-GCN网络，网络训练完成后，由预处理后的测试集评估AR-GCN网络的泛化能力；

三、将步骤一中获得的LR点云数据输入到训练好的AR-GCN中，得到HR点云数据，从而实现口腔三维点云的超分辨成像。

最后需要注意的是，公布实施例的目的在于帮助进一步理解本发明，但是本领域的技术人员可以理解：在不脱离本发明及所附的权利要求的精神和范围内，各种替换和修改都是可能的。因此，本发明不应局限于实施例所公开的内容，本发明要求保护的范围以权利要求书界定的范围为准。