[发明专利]一种超声或CT医学影像三维重建方法

权利要求书

1.一种超声或CT医学影像三维重建方法，其特征在于，该方法输入一个超声或者CT影像序列，其影像分辨率为M×N，100≤M≤1500，100≤N≤1500，三维重建的过程具体包括以下步骤：

步骤1：构建数据集

(a)构建自然图像数据集D

选取一个自然图像网站，要求具有图像序列及对应的摄像机内部参数，从所述自然图像网站下载a个图像序列及序列对应的内部参数，1≤a≤20，对于每个图像序列，每相邻3帧图像记为图像b、图像c和图像d，将图像b和图像d按照颜色通道进行拼接，得到图像τ，由图像c与图像τ构成一个数据元素，图像c为自然目标图像，图像c的采样视点作为目标视点，图像b、图像c和图像d的内部参数均为e_t(t＝1，2，3，4)，其中e₁为水平焦距，e₂为垂直焦距，e₃及e₄是主点坐标的两个分量；如果同一图像序列中最后剩余图像不足3帧，则舍弃；利用所有序列构建数据集D，数据集D有f个元素，而且3000≤f≤20000；

(b)构建超声影像数据集E

采样g个超声影像序列，其中1≤g≤20，对于每个序列，每相邻3帧影像记为影像i、影像j和影像k，将影像i和影像k按照颜色通道进行拼接得到影像π，由影像j与影像π构成一个数据元素，影像j为超声目标影像，影像j的采样视点作为目标视点，如果同一影像序列中最后剩余影像不足3帧，则舍弃，利用所有序列构建数据集E，数据集E有F个元素，而且1000≤F≤20000；

(c)构建CT影像数据集G

采样h个CT影像序列，其中1≤h≤20，对于每个序列，每相邻3帧记为影像l、影像m和影像n，将影像l和影像n按照颜色通道进行拼接得到影像σ，由影像m与影像σ构成一个数据元素，影像m为CT目标影像，影像m的采样视点作为目标视点，如果同一影像序列中最后剩余影像不足3帧，则舍弃，利用所有序列构建数据集G，数据集G有ξ个元素，而且1000≤ξ≤20000；步骤2：构建神经网络

神经网络处理的图像或影像的分辨率均为p×o，p为宽度，o为高度，以像素为单位，100≤o≤2000，100≤p≤2000；

(1)网络A的结构

张量H作为输入，尺度为α×o×p×3，张量I作为输出，尺度为α×o×p×1，α为批次数量；

网络A由编码器和解码器组成，对于张量H，依次经过编码和解码处理后，获得输出张量I；

编码器由5个残差单元组成，第1至5个单元分别包括2，3，4，6，3个残差模块，每个残差模块进行3次卷积，卷积核的形状均为3×3，卷积核的个数分别是64，64，128，256，512，其中，第一个残差单元之后，包含一个最大池化层；

解码器由6个解码单元组成，每个解码单元均包括反卷积和卷积两步处理，反卷积和卷积处理的卷积核形状、个数都相同，第1至6个解码单元卷积核的形状均为3×3，卷积核的个数分别是512，256，128，64，32，16，编码器与解码器的网络层之间进行跨层连接，跨层连接的对应关系为：1与4、2与3、3与2、4与1；

(2)网络B的结构

张量J和张量K作为输入，尺度分别为α×o×p×3和α×o×p×6，张量L和张量O作为输出，尺度分别为α×2×6和α×4×1，α为批次数量；

网络B由模块P和模块Q构成，共有11层卷积单元，首先，将张量J和张量K按照最后一个通道进行拼接，获得尺度为α×o×p×9的张量，对于所述张量经过模块P和模块Q处理后，分别获得输出张量L和张量O；

模块Q与模块P共享前4层卷积单元，前4层卷积单元结构为：前两层单元中卷积核尺度分别为7×7、5×5，第3层到第4层的卷积核尺度均为3×3，1至4层的卷积核个数依次为16、32、64、128；

对于模块P，除了共享4层外，占据网络B的第5层到第7层卷积单元，卷积核尺度均为3×3，卷积核个数均为256，第7层的处理结果利用12个3×3的卷积核进行卷积处理后，将12个结果顺次排成2行，得到张量L的结果；

对于模块Q，除了共享网络B的1至4层外，还占据网络B的第8至11层卷积单元，网络B的第2层输出作为网络B的第8层输入，第8至11层卷积单元中卷积核的形状均为3×3，卷积核个数均为256，第11层的结果利用4个3×3的卷积核进行卷积处理后，从4个通道得到张量O的结果；

(3)网络C的结构

张量R和张量S作为输入，尺度分别为α×o×p×3和α×o×p×6，张量T作为输出，尺度为α×1，α为批次数量；

网络C作为判别器，由5层卷积单元组成，首先，将张量R和张量S按照最后一个通道进行拼接，获得尺度为α×o×p×9的张量，对于所述张量经过以下5层卷积单元处理后获得输出张量T；

网络C的第1至5层卷积单元分别包含1次卷积处理，1次激活处理，卷积核的尺度均为3×3，1至5层卷积单元中卷积核的个数分别为64、128、256、512、1，从第5层结果获得张量T；

步骤3：神经网络的训练

分别将数据集D、数据集E和数据集G中样本按照9:1划分为训练集和测试集，训练集中数据用于训练，测试集数据用于测试，在下列各步骤训练时，分别从对应的数据集中获取训练数据，统一缩放到分辨率p×o，输入到对应网络中，迭代优化，通过不断修改网络模型参数，使得每批次的损失达到最小；

在训练过程中，各损失的计算方法：

内部参数监督合成损失：在自然图像的网络模型参数训练中，将网络A的输出张量I作为深度，将网络B的输出结果L与训练数据的内部参数标签e_t(t＝1，2，3，4)分别作为位姿参数和摄像机内部参数，根据计算机视觉算法，利用图像b和图像d分别合成图像c视点处的两个图像，利用图像c分别与所述的两个图像，按照逐像素、逐颜色通道强度差之和计算得到；

无监督合成损失：在超声或者CT影像的网络模型参数训练中，将网络A的输出张量I作为深度，将网络B的输出张量L和张量O分别作为位姿参数和摄像机内部参数，根据计算机视觉算法，利用目标影像的两相邻影像分别合成目标影像视点处的影像，利用目标影像分别与所述目标影像视点处的影像，按照逐像素、逐颜色通道强度差之和计算得到；

内部参数误差损失：利用网络B的输出结果O与训练数据的内部参数标签e_t(t＝1，2，3，4)按照各分量差的绝对值之和计算得到；

空间结构误差损失：在超声或者CT影像的网络模型参数训练中，将网络A的输出张量I作为深度，将网络B的输出张量L和张量O分别作为位姿参数和摄像机内部参数，根据计算机视觉算法，以目标影像的视点为摄像机坐标系的原点，对目标影像进行重建，采用RANSAC算法对重建点进行空间结构的拟合，利用目标影像的每个重建点到空间几何结构的欧式距离计算得到；

(1)在数据集D上，分别对网络A及网络B的模块P训练80000次

每次从数据集D中取出训练数据，统一缩放到分辨率p×o，将图像c输入网络A，将图像c及图像τ输入网络B，对网络B的模块P进行训练，每批次的训练损失由内部参数监督合成损失计算得到；

(2)在数据集D上，对网络B的模块Q训练80000次

每次从数据集D中取出训练数据，统一缩放到分辨率p×o，将图像c输入网络A，将图像c及图像τ输入网络B，对网络B的模块Q进行训练，每批次的训练损失由内部参数监督合成损失和内部参数误差损失之和计算得到；

(3)在数据集E上，按照如下步骤对网络B的模块Q训练80000次

每次从数据集E中取出超声训练数据，统一缩放到分辨率p×o，将影像j输入网络A，将影像j及影像π输入网络B，对网络B的模块Q进行训练，每批次的训练损失由无监督合成损失和空间结构误差损失之和计算得到；

(4)在数据集E上，对网络C训练80000次得到模型参数ρ

每次从数据集E中取出超声训练数据，统一缩放到分辨率p×o，将影像j输入网络A，将影像j及影像π输入网络B，然后将网络A的输出作为深度，将网络B的输出作为位姿参数和摄像机内部参数，分别根据影像i和影像k合成影像j视点处的两张影像，将所述的两张影像按照颜色通道拼接，将影像j与所述拼接影像输入到网络C，网络C训练时，通过不断修改网络C的参数，迭代后得到最优的网络模型参数ρ；每批次的损失Ψ的计算方法为：

Ψ＝δ[log(T)]+δ[log(1-T)] (1)

其中，δ为真伪分布函数，log(T)是对T进行取自然对数，T是网络C的输出；

(5)在数据集G上，对网络C训练80000次得到模型参数ρ′

每次从数据集G中取出CT训练数据，统一缩放到分辨率p×o，将影像m输入网络A，将影像m及影像σ输入网络B，然后将网络A的输出作为深度，将网络B的输出作为位姿参数和摄像机内部参数，分别根据影像l和影像n合成影像m视点处的两张影像，将所述的两张影像按照颜色通道拼接，将影像m与所述拼接影像输入到网络C，网络C训练时，通过不断修改网络C的参数，迭代后得到最优的网络模型参数ρ′；每批次的损失为无监督合成损失、空间结构误差损失和摄象机平移运动损失Y之和，其中，Y由网络B的输出位姿参数，根据摄像机平移运动的约束计算得到；

步骤4：超声或者CT影像三维重建

利用自采样的一个超声或者CT序列影像，将每一帧影像统一缩放到分辨率p×o，使用模型参数ρ或者模型参数ρ′进行预测，对于超声序列影像，将影像j输入网络A，将影像j及影像π输入网络B，对于CT序列影像，将影像m输入网络A，将影像m及影像σ输入网络B，将网络A的输出作为深度，将网络B的输出作为位姿参数和摄像机内部参数，按照下列步骤选取关键帧，序列中第一帧作为当前关键帧，依次将序列影像中的每一帧作为目标帧，根据当前关键帧，利用摄像机位姿参数和内部参数，合成目标帧视点处的影像，利用所述合成影像与目标帧之间逐像素逐颜色通道强度差之和的大小计算误差λ，再根据目标帧的相邻帧，利用摄像机位姿参数和内部参数，合成目标帧视点处的影像，利用所述合成影像与目标帧之间逐像素、逐颜色通道强度差之和的大小计算误差γ，进一步利用公式(2)计算合成误差比Z，当Z大于阈值η时，1η2，将当前关键帧更新为此时的目标帧；

对任一目标帧，将其分辨率缩放到M×N，根据摄像机内部参数，依据计算机视觉的重建算法，计算每帧影像每个像素的摄象机坐标系中的三维坐标，进一步，将第一帧的视点作为世界坐标系的原点，再结合所有关键帧的位姿参数，利用三维空间几何变换，计算得到该序列每一帧影像每个像素的世界坐标系中的三维坐标。