CN105023282A

CN105023282A - 一种基于压缩感知的稀疏投影超声ct图像重建方法

Info

Publication number: CN105023282A
Application number: CN201410180951.1A
Authority: CN
Inventors: 尉迟明; 丁明跃; 花少炎
Original assignee: Huazhong University of Science and Technology
Current assignee: Huazhong University of Science and Technology
Priority date: 2014-04-30
Filing date: 2014-04-30
Publication date: 2015-11-04

Abstract

本发明公开了一种基于压缩感知的稀疏投影超声CT图像重建方法，该方法包含以下步骤：首先使用超声CT以远低于Nyquist定理所需次数对目标进行扫描，获得有限次方向上投影数据。对得到的投影数据进行一维傅里叶变换，依据傅里叶衍射投影定理，得到目标断层截面空间频率信息。通过选择变换基，使得目标在变换域具有稀疏性，依照压缩感知原理，构建基于l₁范数优化的逆问题。最后通过共轭迭代算法，重建目标图像。本方法能大幅度降低投射式超声CT投影次数，降低系统的控制复杂度，降低设备的数据存储量，而且显著降低系统扫描时间，临床上能减少病人呼吸运动带来的伪影；通过迭代求解，获得更高的图像重建质量，在医学成像领域具有极大的应用价值。

Description

一种基于压缩感知的稀疏投影超声CT图像重建方法

技术领域

本发明属于医学成像领域，，更具体地，涉及一种基于压缩感知的稀疏投影超声CT图像重建方法。

背景技术

随着信息处理技术及加工工艺的进步，超声层析成像(ComputedTomography，简称CT)获得越来越多的关注。与传统的B超显示人体组织结构影像不同，它能够提供人体组织声学参量如折射率，衰减系数，密度的分布，是一种功能性成像，具有分辨率高，对致密组织敏感等特点，具有很高的临床应用价值。例如，在乳腺检查中，与X射线检查相比，它无辐射危害，而且可以区分良性和恶性肿瘤。早期超声CT借鉴X-CT技术，假设超声波在人体内沿直线传播，忽略超声波在介质中的折射衍射效应，然后利用X-CT算法重建图像。其所重建图像误差较大，不能获得目标内部细节信息。目前广泛采用的是超声衍射层析成像技术，其利用波动方程理论，将介质看成是声学参量连续变换的非均与体。在弱散射条件下，利用Bom或Rytov近似，可以建立傅里叶衍射投影定理，将目标的空间频率和投影散射场的傅里叶变换联系起来，进而可以通过傅里叶变换快速的重建目标。

依据Nyquist采样理论，对超声CT而言，为了抑制伪影，就必须增加投影次数和超声探头的接收通道数目。对实际成像系统而言，这就会增加扫描时间和设备成本，同时对系统的控制精度也提出了很高的，多次投影导致采样信息中包含大量的冗余信息。因此，用有限次投影高精度重建图像具重要的理论和实际意义。此外，传统的重建算法如，时域插值法，频率插值法，误差大，重建质量差。

发明内容

针对现有技术的以上缺陷或改进需求，本发明提供了一种基于压缩感知的稀疏投影超声CT图像重建方法，其目的在于，利用压缩感知理论，突破Nuquist采样理论限制，以有限次稀疏投影对目标扫描，以此解决现有超声CT系统扫描时间长、系统复杂度高的技术问题，同时利用迭代重建提高了图像质量。

为实现上述目的，按照本发明的一个方面，提供了一种基于压缩感知理论的超声CT图像重建方法，包括以下步骤：

(1)基于傅里叶衍射投影定理，获得超声CT目标的空间频率；具体包括以下子步骤：

(1-1)将发射阵列和接收阵列以相对位置放置在超声CT目标两侧，其中发射阵列用于发射平面波，接收阵列用于获得超声CT目标对平面波的散射声场；

(1-2)对接收到的散射声场，利用一维傅里叶变换，建立散射声场与超声CT目标的空间频率之间的映射关系：FFT_1D{u_s，θ(ζ)}＝FFT_2D{f(x，y)}；其中f(x，y)为超声CT目标截面在二维笛卡尔坐标系(x，y)}下的声学参量函数，u_s，θ(ζ)代表入射角为θ的平面波的散射声场，FFT_1D{u_s，θ(ζ)}表示散射声场的傅里叶变换，(ζ，η)为与入射角θ一致的笛卡尔坐标系，FFT_2D{f(x，y)}表示超声CT目标的目标空间频率，其轨迹为超声CT目标空间频率域内的一条圆弧线，其圆心为其中为平面波入射方向的单位向量，k₀为入射波波数。

(1-3)随机的改变发射阵列的平面波发射角度达10至30次并重复步骤(1-2)，以获得不同的目标空间频率FFT_2D{f(x，y)}和观测矩阵Φ，这些目标空间频率的轨迹位于半径为2k₀的圆周内，其中观测矩阵Φ的维数为m×n，其中m为投影角度改变的次数乘以探头的阵元数，n为声学参量函数f(x，y)离散的维度；

(2)确定声学参量函数f(x，y)的系数变换基Ψ和变换域系数S；

(3)利用步骤(1)测得的目标空间频率和步骤(2)确立的声学参量函数f(x，y)的系数变换基Ψ，并基于压缩感知方法建立超声CT目标重建的逆问题

\min_{s} {{| | s | |}_{1} : Θs = p}

或

\min_{s} {α {| | s | |}_{1} + {| | Θs - p | |}_{2}^{2}};

其中α为规则化参数，感知矩阵Θ＝ΦΨ，p＝FFT_1D{u_s，θ(ζ)}；

(4)利用包括非一致快速傅里叶变换的迭代算法对步骤(3)建立的超声CT目标重建的逆问题进行求解，以得到声学参量函数f(x，y)。

优选地，步骤(2)具体为，利用目标图像库函数对声学参量函数f(x，y)进行训练以得到系数变换基Ψ，然后利用f＝Ψs，计算声学参量函数f(x，y)在系数变换基Ψ下的变换域系数S，并使变换域系数S中的大部分为零或其幅度按指数衰减。

优选地，目标图像库函数包括超声CT目标的CT、MR、或B模式图像。

优选地，步骤(2)具体为，结合人体组织结构特点，利用TV方法刻画超声CT目标的梯度稀疏性，以获取声学参量函数f(x，y)的系数变换基Φ和变换域系数S。

优选地，迭代算法可以采用包括非一致快速傅里叶变换的共轭梯度法或内点法。

总体而言，通过本发明所构思的以上技术方案与现有技术相比，能够取得下列有益效果：

1、本发明利用目标的稀疏性，以压缩感知理论构建逆问题，突破Nyquist限制，降低投影次数，减少通道数目，降低系统设备成本，减少扫描时间。

2、本发明利用傅里叶衍射投影定理，充分考虑超声CT扫描时声波在介质内传播过程中的衍射效应，从而能获得较高质量的图像。

3、本发明将非一致快速傅里叶变换用于迭代重建，加快了算法重建速度。

附图说明

图1是本发明基于压缩感知的稀疏投影超声CT图像重建方法的流程图。

图2示出本发明中所用的超声CT原理图。

图3示出本发明中对投影所得散射场进行傅里叶变换所得目标空间频率分布。

图4(a)是传统插值算法的示意图，图4(b)是本发明方法的示意图。

图5(a)为16次稀疏投影下传统插值方法的重建结果，图5(b)为16次稀疏投影下本发明方法的重建结果。

图6(a)为96次稀疏投影下传统插值方法的重建结果，图6(b)为96次稀疏投影下本发明方法的重建结果。

具体实施方式

为了使本发明的目的、技术方案及优点更加清楚明白，以下结合附图及实施例，对本发明进行进一步详细说明。应当理解，此处所描述的具体实施例仅仅用以解释本发明，并不用于限定本发明。此外，下面所描述的本发明各个实施方式中所涉及到的技术特征只要彼此之间未构成冲突就可以相互组合。

在实际中，由于超声发射阵列也可以用于接收超声回波，因此，超声CT有反射式和透射式两种。无论是透射式还是反射式，根据傅里叶衍射投影定理，都可以建立散射场和目标空间频率信息的联系，其不同在于：反射模式中接收声场包含目标的高频信息，透射模式接收声场包含目标的低频信息。本发明中以透射模式说明以压缩感知理论为基础的超声CT图像重建，但本方法可以用于反射式超声CT。

如图1所示，本发明基于压缩感知的稀疏投影超声CT图像重建方法包括以下步骤：

(1)基于傅里叶衍射投影定理，获得超声CT目标的空间频率；本步骤具体包括以下子步骤：

(1-1)利用两个线性阵列超声探头，以相对位置放置在超声CT目标两侧，选定一个做发射阵列，另一个做接收阵列，其中发射阵列用于发射平面波，接收阵列用于获得超声CT目标对平面波的散射声场；

(1-2)对接收到的散射声场，利用一维傅里叶变换，建立散射声场与超声CT目标的空间频率之间的映射关系：FFT_1D{u_s，θ(ζ)}＝FFT_2D{f(x，y)}；其中f(x，y)为超声CT目标截面在二维笛卡尔坐标系(x，y)}下的声学参量函数，u_s，θ(ζ)代表入射角为θ(如图2所示)的平面波的散射声场，FFT_1D{u_s，θ(ζ)}表示散射声场的傅里叶变换，(ζ，η)为与入射角θ一致的笛卡尔坐标系，FFT_2D{f(x，y)}表示超声CT目标的目标空间频率，其轨迹为超声CT目标空间频率域内的一条圆弧线，如图3所示，其圆心为其中为平面波入射方向的单位向量，k₀为入射波波数，实部圆弧代表透射声场所对应的空间频率域轨迹，虚线圆弧代表反射声场所对应的空间频率域轨迹。

(1-3)随机的改变发射阵列的平面波发射角度达10至30次并重复步骤(1-2)，以获得不同的目标空间频率FFT_2D{f(x，y)}和观测矩阵Φ，这些目标空间频率的轨迹位于半径为2k₀的圆周内，其中观测矩阵Φ的维数为m×n，其中m为投影角度改变的次数乘以探头的阵元数，n为声学参量函数f(x，y)离散的维度；由于目标空间频率的轨迹为弧形，传统方法利用逆傅里叶变换，需要对空间网格进行插值，如图4(b)所示。为了保证插值精度，需要大量投影和增加阵列阵元个数，而本发明利用压缩感知，只需超声CT目标以随机的角度进行投影，获得稀疏的目标空间频率值，如图4(a)所示。

(2)确定声学参量函数f(x，y)的系数变换基Ψ和变换域系数S；具体而言，利用已有的目标图像库函数(如同一超声CT目标的CT、MR、或B模式图像)对声学参量函数f(x，y)进行训练以得到系数变换基Ψ，然后利用f＝Ψs，计算声学参量函数f(x，y)在系数变换基Ψ下的变换域系数S，并使变换域系数S中的大部分为零或其幅度按指数衰减。备选地，也可结合人体组织结构特点，利用总变分(Total variant，简称TV)方法刻画超声CT目标的梯度稀疏性，以获取声学参量函数f(x，y)的系数变换基Φ和变换域系数S。

(3)利用步骤(1)测得的目标空间频率和步骤(2)确立的声学参量函数f(x，y)的系数变换基Ψ，并基于压缩感知(Compressed sensing，简称CS)方法建立超声CT目标重建的逆问题或其中α为规则化参数，感知矩阵Θ＝ΦΨ，p＝FFT_1D{u_s，θ(ζ)}；

(4)利用包括非一致快速傅里叶变换的迭代算法对步骤(3)建立的超声CT目标重建的逆问题进行求解，以得到声学参量函数f(x，y)；具体而言，迭代算法可以采用包括非一致快速傅里叶变换的共轭梯度法、内点法等多种求解算法。

图5显示传统插值法和本发明方法在透射模式下重建结果的对比，其中步骤(1-3)中的次数设置为16次。图5(a)为插值法的重建结果，图5(b)为本发明的重建结果。从图中可以看出，在有限次投影条件下，本发明方法优于传统的插值算法。由于投影次数远低于Nyquist定理的要求，图5(a)存在严重的伪影，除了比较大的结构轮廓可以辨识外，细节几乎全被伪影淹盖。图5(b)除了部分外形轮廓出现失真，整体机构清晰，尤其是底部三个小椭圆，虽然不能完全分离，但整体轮廓可以分辨，表明本发明的方法具有更好的分辨率。

不仅在有限次投影下，压缩感知能较好的重建目标，在投影次较多的情况下，压缩感知也能提高图像质量。图6为插值法和本发明方法在96次投影下重建结果对比。图6(a)为插值法的重建结果，图6(b)为本发明的重建结果。从图中可以发现，在96次投影下，通过CS重建的图像质量高于插值法。对于插值法，随着投影次数的增加，虽然伪影尤其是背景部分逐渐减少，但伪影始终存在，图像内部细节无法清晰的识别出来。但通过本发明的方法，模型内部细节特征已经能够彻底的重建出来。

本领域的技术人员容易理解，以上所述仅为本发明的较佳实施例而已，并不用以限制本发明，凡在本发明的精神和原则之内所作的任何修改、等同替换和改进等，均应包含在本发明的保护范围之内。

Claims

1.一种基于压缩感知的稀疏投影超声CT图像重建方法，其特征在于，包括以下步骤：

(2)确定声学参量函数f(x，y)的系数变换基Ψ和变换域系数S；

\min_{s} {{| | s | |}_{1} : Θs = p}

或

\min_{s} {α {| | s | |}_{1} + {| | Θs - p | |}_{2}^{2}};

2.根据权利要求1所述的稀疏投影超声CT图像重建方法，其特征在于，步骤(2)具体为，利用目标图像库函数对声学参量函数f(x，y)进行训练以得到系数变换基Ψ，然后利用f＝Ψs，计算声学参量函数f(x，y)在系数变换基Ψ下的变换域系数S，并使变换域系数S中的大部分为零或其幅度按指数衰减。

3.根据权利要求2所述的稀疏投影超声CT图像重建方法，其特征在于，目标图像库函数包括超声CT目标的CT、MR、或B模式图像。

4.根据权利要求1所述的稀疏投影超声CT图像重建方法，其特征在于，步骤(2)具体为，结合人体组织结构特点，利用TV方法刻画超声CT目标的梯度稀疏性，以获取声学参量函数f(x，y)的系数变换基Φ和变换域系数S。

5.根据权利要求1所述的稀疏投影超声CT图像重建方法，其特征在于，迭代算法可以采用包括非一致快速傅里叶变换的共轭梯度法或内点法。