CN111213066B - 基于模型的图像重建方法 - Google Patents

Abstract

本发明涉及一种用于例如超声成像的脉搏波图像重建方法。所提出的方法基于图像测量模型及其伴随算子。所提出的方法介绍了测量模型及其伴随算子的无矩阵形式。所提出的方法具有重建图像质量高、重建速度快的优点。

Description

基于模型的图像重建方法

技术领域

本发明涉及一组适用于例如超声成像的基于模型的图像重建方法。本发明还涉及一种用于执行该方法的成像设备。

背景技术

脉搏波成像应用(例如，超声(US)成像应用)可以用于各种技术领域。在医学应用中，用于观察内部身体结构，例如，肌腱、肌肉、关节、血管和内部器官。然而，脉搏波成像技术也广泛应用于医疗领域之外。超声波的一个商业用途是无损检测(NDT)，即在不破坏零件或系统的适用性的情况下，检查、测试或估计材料、部件或组件的不连续性或特性差异的过程。换言之，当检查或测试完成后，零件仍可使用。无损检测的一个示例是油管。在这种情况下，超声波用于沿管道向下传播，并且以这种方式，可以检测是否存在作为脉冲反射的裂缝或缺陷。通过了解管道中的波的速度以及从裂缝反射到发送波之间的时间间隔，可以估计缺陷的位置。用于管道裂缝检测的相同理论可用于检测飞机结构/机翼面板等的裂缝，并检测例如机翼表面上的结冰。脉冲波成像当然还有许多其他可能的NDT应用。脉搏波成像因其安全性、便携性和实时性而广受欢迎。

传统的US成像使用多个聚焦的时间分离的传输光束沿着扫描线重建图像。因此，帧速率受到聚焦发射波束数量的限制，每秒不能超过几十个图像。虽然这种帧速率对于许多应用来说已经足够了，但是更复杂的现象的理解(例如，用于弹性成像的剪切波的传播)需要更高的帧速率。

克服这一限制的一种方法是尽可能减少发射波束的数量。合成孔径(SA)方法试图通过仅使用几个换能器元件依次声穿透整个介质来克服这个问题。另一种选择是超快US(UFUS)成像，利用平面波(PW)或发散波(DW)的概念来立即对整个视场进行声穿透，使US系统每秒钟能达到数千帧。

这种方法的一个重要限制是图像质量下降。实际上，与能量集中在感兴趣的特定区域的聚焦发射波束相比，PW或DW的能量分散在声穿透介质上，导致较低的信噪比和较差的空间分辨率。解决这个问题的一种方法是对多个图像进行平均，这些图像是在PW成像的情况下以不同的声穿透角获得的，是在UFUS的DW成像的情况下以不同的点源位置获得的，或者是在SA成像的情况下以不同的换能器元件组获得的。在UFUS中，这一过程称为相干合成。虽然这种技术的实现相当简单，但需要多次声穿透，从而降低了帧速率。

对多次声穿透求平均的另一种方法是使用比传统技术更有效的重建方法，传统技术主要围绕延迟和求和(DAS)波束形成。一组流行的方法依赖于使用迭代算法来解决由US图像重建引起的不适定图像形成问题。这些方法建立在问题的前向模型上。已知解决方案的主要问题在于其计算复杂性，通常转化为相应矩阵表示的存储要求。一些提出的模型要求二维(2D)矩阵系数存储数百千兆字节。鉴于传统方法，这个问题严重限制了迭代方法的吸引力。

发明内容

本发明的一个目的是克服至少一些与脉搏波成像中的图像重建相关的上述问题。

根据本发明的第一方面，提供了一种根据权利要求1所述的图像重建方法。

所提出的方法具有重建图像质量高、重建速度快的优点。换言之，与现有技术解决方案相比，图像对比度和分辨率可以显著提高。

根据本发明的第二方面，提供了一种成像设备，其被配置为执行根据权利要求15所述的根据本发明的第一方面的方法。

在所附的从属权利要求中叙述本发明的其他方面。

附图说明

参考附图，从下面对非限制性示例性实施例的描述中，本发明的其他特征和优点将变得显而易见，其中：

图1是示出根据本发明的示例的成像系统的一些元件的框图，这些元件有助于理解本发明的教导；

图2是示出根据本发明的一个示例的所提出的成像方法中使用的成像配置的示意图；

图3是示出根据本发明的一个示例的所提出的成像方法中使用的成像配置和连续域反射率超曲面的示意图；

图4是示出根据本发明的一个示例的所提出的成像方法中使用的成像配置和离散域反射率超曲面的示意图；

图5是概括根据本发明示例的图像重建方法的流程图；

图6是根据本发明的一个示例更详细地示出图像重建方法的流程图；

图7是根据本发明的另一示例更详细地示出图像重建方法的流程图；

图8是示出根据本发明的示例实现图像测量模型的过程的流程图；以及

图9是示出根据本发明的示例实现测量模型的伴随算子的过程的流程图。

具体实施方式

现在将参考附图详细描述本发明的实施例。在US成像的背景下描述该实施例，但是本发明的教导不限于这种环境。本发明的教导同样适用于可以使用脉搏波成像的其他领域。出现在不同附图中的相同或相应的功能和结构元件具有相同的附图标记。

本实施例中的本发明提出了一种基于解决以下形式的不适定线性逆问题的2D和3D超声成像环境中的图像形成方法和设备或系统：

m＝Hγ+n

其中，m是测量值，γ是受检查的对象，或者更具体地，是反射率函数，H是测量模型，n是噪声。所提出的图像形成方法依赖于两个主要支柱：快速且无矩阵的测量模型H和图像重建方法，该方法允许在给定测量值m的情况下检索受检查的对象γ的估计。该对象可以是较大结构的内部部分或元件。所提出的方法可以由标准的US系统执行，该系统包括US探头、图像形成模块以及后处理和显示模块，如图1所示。

图1是示出成像系统或设备1的一些元件的简化框图，这些元件有助于理解本发明的教导。系统1包括超声波探头3，在该示例中，超声波探头3包括压电元件或传感器的线性阵列，也称为换能器。这种探头适合于2D成像的目的，并且可以包括圆柱形聚焦方案，以在高度上将声能量集中在图像平面附近。然而，这种探头仅仅是一个示例。所提出的方法可以容易地适应不同的探头几何形状(例如，凸形探头、相控阵等)或技术(包括聚偏氟乙烯(PVDF)共聚物和电容性微机械超声换能器(CMUT)等)。类似地，所提出的方法可以扩展到使用2D矩阵探头的情况，所述探头被设计成传输声束并从声穿透的体积收集回波，以便提供体积成像，由此成像过程产生受检查的体积的3D表示。

US探头3连接到图像形成模块或设备5，其被配置为执行图像形成。在图像形成过程中使用在测量模型单元7中估计的测量模型。稍后详细解释的所提出的测量模型基于脉冲回波空间脉冲响应模型。本发明介绍了测量模型及其伴随的无矩阵公式，在连续算子的情况下，算子的伴随被定义为方阵转置的连续扩展，这表示脉搏波成像中的一种新的成像方法。测量模型的描述可以通过两个主要步骤来实现：

·导出了脉冲波往返飞行时间的参数方程。这将表明，模型可以重铸为超曲面上的一个积分，如后文所述。

·所提出的积分可随后离散，如后所述。

可以例如为PW、DW和SA成像导出该模型。该模型与射频、同相和正交(IQ)数据兼容。最后，导出了所提出的测量模型的伴随的无矩阵形式。将证明也可以被重铸为超曲面上的积分，在2D成像的情况下被定义为1D流形，在3D成像的情况下被定义为2D流形，进行参数化，并且可以使用为测量模型提议的相同的离散方案。

图像形成设备5还包括用于运行图像重建过程或方法的图像重建单元9。图像重建方法旨在检索受检查的图像的估计

估计值可以写成：

其中，f说明了具有超参数Ξ的图像重建过程。在本说明书中呈现了两种不同的图像重建方法：

·f(m,H,Ξ)具有显示表达式的分析方法(例如，使用H伴随算子的反向传播方法)。

·f(m,H,Ξ)是一种迭代算法的迭代方法，解决了以下问题：

其中，

表示目标函数，该目标函数包括测量Hγ和m之间误差的函数，也称为数据差异项，

表示非负函数(正则化)，表示正则化参数。在这种情况下，Ξ是与用于解决最小化问题的优化算法相关的超参数，例如，迭代次数、算法的停止标准等。

在图像重建之后，重建的数据发送到后处理和显示单元或模块11。后处理步骤涵盖了一系列应用，例如，B模式成像、彩色多普勒成像、矢量多普勒成像和弹性成像。

对于使用射频数据的B模式成像，包络检测应用到重建数据上。例如，包络检测可以通过希尔伯特变换实现，之后是幅度检测和可选的低通滤波。也可以通过对信号进行平方和低通滤波来实现。如果对IQ数据使用B模式成像，则提取信号的幅度。包络检测步骤之后是归一化和动态压缩步骤。对于多普勒成像和弹性成像，重建的射频或IQ数据直接使用，无需后处理。

在详细解释所提出的方法之前，首先解释本说明书中使用的一些符号。在本说明书中考虑的脉冲回波配置中，由传感器13检测以其反射率(函数)γ(r)(其中，r＝[x,y,z]^T∈Ω)为特征的来自介质

的回波。图2示出了一种用于PW成像的特定简化测量配置。在本说明书中使用了以下符号：

数学或类似的符号：

·矢量用粗体表示。

·

表示向量

的

-范数，其中，

表示实数空间。

·希尔伯特空间W的平方可积函数的空间表示为L₂(W)。两个函数f,g∈L₂(W)之间的内积定义为<f|g>_W＝∫_x∈Wf(x)g(x)dx。

超声相关符号：

·m(p,t)是位于p位置的传感器13在某一时刻t接收到的电信号。换言之，m(p,t)定义了在测量网格点(p,t)计算的一组连续测量组的一个样本(相当于一个元素原始数据的样本)。

·γ(r)定义在反射率网格点r计算的一组反射率值(相当于图像样本)中的一个值。

·Ω表示感兴趣的介质或对象，即反射率样本的一组连续位置。

·o(r,p)是位于p的传感器和位于r的点之间的元件方向性。

·a(r,p)定义位于p的点源和位于r的点之间的衰减值。

·v_pe(t)是接收到的脉冲形状(即脉冲回波波形)。

·t_Tx(r,q)表示传输传播延迟，即由位于q的传输装置或元件13传输的脉冲波到达位于r的点所用的时间。

·t_Rx(r,p)表示接收传播延迟，即由位于r的点传输的脉冲回波到达位于p的传感器所用的时间

测量网格相关符号：

·一组传输元件：

其中，N_elt指传输元件的数量，也称为发射机。

·一组接收传感器：

其中，N_elr表示接收传感器的数量，也称为接收器。

·传输元件内部的一组点源位置：

·接收传感器内部的一组点源位置：

·一组时间样本：

反射率网格相关符号：

·反射率值的一组位置：Ω_d＝{rⁿ∈Ω,n＝{1...N_γ}}，其中，N_γ指反射率值的数量。

所提出的方法定义了都位于3D空间中的第一网格和第二网格。第一网格称为反射率网格，其包括在两条交叉网格线的交点处的反射率网格点，在该示例中，这两条网格线相对于彼此形成大致90度的角度。反射率网格的每个反射率网格点都属于介质Ω。第二网格称为测量网格，其包括在两条交叉网格线的交点处的测量网格点，在该示例中，这两条网格线相对于彼此形成大致90度的角度。测量网格点在一个维度上与接收传感器的位置一致，在另一维度上与每个传感器的采样时刻一致。

接下来解释根据本发明在给定测量网格点的测量模型的估计、评估、计算或实现。在这个示例中，在每个测量网格点估计测量模型。在本公开中使用的测量模型是脉冲回波空间脉冲响应模型，例如，G.E.Tupholme的“Generation of acoustic pulses by baffledplane pistons”，Mathematika，第16卷，第209页，1969年以及P.R.Stepanishen的“TheTime-Dependent Force and Radiation Impedance on a Piston in a Rigid InfinitePlanar Baffle”，J.Acoust.Soc.Am.，第49卷，第76页，1971年，首先介绍了该模型。根据这个模型，可以写出下面的方程，该方程说明了在位置p和时间t记录的原始数据样本。

接下来解释一个传输点源元件和一个接收点源传感器的配置的测量模型的连续参数公式。

测量模型定义为：

上述方程可表示如下，并将测量模型定义为沿超曲面Γ(p,q,t)的积分：

其中，

Γ(p,q,t)＝{r∈Ω|g(r,p,q,t)＝0}

g(r,p,q,t)＝t-t_Tx(r,q)-t_Rx(r,p)

在上面的方程中，超曲面(或者更具体地，反射率超曲面)Γ(p,q,t)由方程g(r,p,q,t)＝0定义，其中，感兴趣的变量是r。当r位于2D平面(抛物线或椭圆)时，Γ(p,q,t)定义二次曲线，当位于3D体积中时定义二次曲线。相同的形状也适用于测量超曲面。图3示意性地示出了连续域中的测量配置，还显示了反射率超曲面。

反射率超曲面Γ(p,q,t)现在可以被参数化，以减少集成维数。由该方程g(r,p,q,t)＝0定义的反射率超曲面Γ(p,q,t)可以等效地由一组反射率参数方程定义。在数学中，参数方程将一组量定义为一个或多个自变量(即参数)的函数。参数方程通常用于表示构成几何对象(例如，曲面或曲线)的点的坐标，在这种情况下，这些方程统称为对象的参数表示或参数化。这组反射率参数方程现在可以如下定义：

等式右边的方程定义了r的反射率参数方程，因为r现在依赖于一个参数α∈A，其中，如果

则

并且当

则

需要注意的是，每个测量值或样本是否有一组反射率参数方程。此外，可以通过绘制一组反射率参数方程来获得超曲面。配备了这组反射率参数方程的上述公式，可以如下导出测量模型H的连续参数公式：

其中，

接下来将解释用于一个任意形状的传输有限长度元件和一个接收点源传感器的配置的连续测量模型。在这种情况下，传输元件具有有限的长度和任意的形状。其特征在于一组连续的点源位置X^j。

测量模型定义为

可以看出，与结合一个传输点源元件和一个接收点源传感器的配置所描述的模型的主要区别在于，还有一个积分，其对应于形成任意形状的有限长度传感器的所有点源的贡献的总和。

相应的参数公式如下：

接下来解释用于一个任意形状的传输有限长度元件和一个接收有限长度传感器的配置的连续测量模型。在这种情况下，接收传感器具有有限的长度和任意的形状。由这组点源位置Πⁱ定义。所得信号的坐标由m(ξⁱ,t)表示，其中，ξⁱ∈Πⁱ。例如，可以是Πⁱ中心的点坐标。

可以写成：

通过沿由Πⁱ定义的传感器表面集成获得所得测量值。当将定义m(p,t)的方程集成为上述积分时，以下关系成立：

相应的参数公式如下：

接下来解释用于任意形状的多个传输有限长度元件和任意形状的一个接收有限长度传感器的配置的连续测量模型。在这种情况下，一组点源位置X^j与传输元件j相关联，导致以下方程：

这个方程非常类似于结合任意形状的一个传输有限长度元件和任意形状的一个接收有限长度传感器的配置所描述的方程。唯一的区别在于介绍的总和，考虑了多个传输元件。

相应的参数公式如下：

接下来解释用于任意形状的多个传输有限长度元件和任意形状的多个接收有限长度传感器的配置的连续测量模型。

当使用一组多个接收传感器时，一组点源位置Πⁱ和坐标ξⁱ与接收传感器i相关联。因此，接收测量字段可以定义为

其仅针对离散坐标ξⁱ定义。在本说明书的剩余部分中，以下内容成立：

接下来离散测量模型的连续参数公式。该步骤包括离散结合任意形状的多个传输有限长度元件和任意形状的一个接收有限长度传感器的配置描述的测量模型的参数公式。然而，上面描述的任何模型都可以离散。离散是指我们选择参数α的一些离散值α^p，这些离散值定义了以下一组参数：A_d＝{α^p∈A,p＝{1...N_α}}。可以看出，参数α的离散值的选择导致超曲面Γ(p,q,t)的离散，现在仅通过N_α点

来描述。图4示意性示出了离散域中的测量配置，还显示了反射率超曲面。测量模型的参数公式的离散可以表示为：

其中，

是在测量网格的每个时间样本上估计的脉冲形状，

其中，

其中，w(q^k)、u(α^p)和z(pⁿ)是集成权重(与连续积分的离散相关)。可以注意到，对于传输元件的每个点q^k，上述方程包括在反射率超曲面的点上估计的一组反射率值

称为超曲面反射率样本或值。

超曲面反射率值接下来进行插值，以基本上与反射率网格点一致。超曲面反射率值通常是未知的，因为它们不在反射率网格点上。因此，超曲面反射率值需要通过插值内核从反射率网格上定义的反射率样本进行近似，如P.Thévenaz、T.Blu和M.Unser在“Interpolation revisited”，美国电气和电子工程师学会医学成像学报，第19卷，第739-758页，2000中所提出的。因此，让我们介绍插值内核

以便：

上述插值方程将超曲面反射率值γ(r(α^p,pⁿ,q^k,t^l))(通常不在反射率网格点上)与反射率样本γ^s相关联。

接下来对离散和插值的超曲面反射率样本或值进行求和。插值方程插入离散测量模型，并执行以下求和：

在该示例中，求和的输出最终与脉冲形状进行卷积，如在前面的方程中所描述的，以便在所考虑的测量网格点处获得测量值的估计，也称为测量值估计。对每个测量值执行上述过程。因此，测量模型的估计或处理的结果是各个测量网格点处的测量值的估计。

接下来在每个反射率网格点进行测量模型的伴随算子的估计、评估、计算或实现。首先，获得测量模型的伴随算子的连续参数公式。可以使用函数分析工具从测量模型中导出伴随算子的连续公式，即从上面给出的任何连续测量模型的示例中导出。从任意形状的多个有限长度传输元件和任意形状的多个有限长度接收传感器的情况导出下面的方程。实际上，如果我们定义算子

则

的伴随算子

定义为：

让我们考虑以下内积：

其中，u(t)＝v_pe(-t)是脉冲回波波形的匹配滤波器。我们现在有：

其中

配备测量模型的伴随算子，给定测量值的估计，我们将反射率网格点r处的反射率值的估计定义为：

这种情况比测量模型简单。实际上，参数公式是直接的，测量参数方程可以表示为简单的飞行时间方程：

t＝t_Tx(r,q)+t_Rx(r,p)

对应于上述方程的测量超曲面是：

接下来离散测量模型的伴随算子的参数公式。与测量模型的离散应用相同的过程。在这种情况下，离散过程更直接，因为直接应用于p上。一组离散的参数是

因此，我们具有以下方程：

其中，w(q^k)和z(pⁿ)是集成权重，

和

是在每个时间采样点估计的脉冲形状的匹配滤波器。可以注意到，上述方程包括位于测量超曲面上的一组卷积测量样本

称为超曲面测量样本或值。

接下来执行超曲面测量样本的插值，以基本上与测量网格点一致。超曲面测量样本的值通常是未知的，因为它们不位于测量网格点上。因此，这些值需要通过插值内核从测量网格上定义的测量样本中进行近似。我们介绍插值内核

使得：

上述方程将超曲面测量样本(m(ξⁱ)*u)(t_Tx(r,q^k)+t_Rx(r,pⁿ))与卷积的测量样本(m(ξⁱ)*u)(t^l)相关联。

之后，对离散和插值的超曲面测量样本求和。插值方程插入到测量模型的离散伴随算子中，并执行以下求和：

求和的输出是在所考虑的反射率网格点的反射率值的估计。

接下来解释图像重建过程。首先定义逆问题。测量值m可用于使用如上所述的测量模型的伴随算子来生成反射率值

的估计。类似地，使用前面描述的测量模型，可以使用反射率值

的估计来生成不等于m的测量值

的估计。因此，我们现在知道了测量值和反射率值之间的关系。

图像重建方法的目的是在给定一组测量值m的情况下，检索在反射率网格点估计的反射率值的估计

测量值和未知反射率值之间的关系由以下线性反问题定义：

其对应于在测量网格的每个点估计的测量模型。

此外，因为算子

相对于γ(可以从在给定网格点估计测量模型时描述的公式中推导出)是线性的，所以存在矩阵

(与线性算子

相关联)，使得

m＝Hγ

其定义了逆问题。

下面解释与逆问题相关的优化问题。为了解决上述反问题，将其表示为以下优化问题：

其中，

表示优化问题中所涉及的目标函数，

是用于数据差异项的下半连续函数，测量这些测量值Hγ的估计值和测量值m之间的距离，

是用于描述可选的先验项的下半连续函数，用于说明反射率值的额外信息，例如，特定的统计行为，λ＞0是正则化参数。

下面列出了用于数据差异术语的函数的一些示例：

欧几里德距离：

半径∈为的

球上的指示器函数：

范数：

先验项所使用的函数的一些示例如下：

幂p的

范数：

当p＜1，是γ的稀疏性的度量。

当p∈[1,2]，是数据与广义高斯分布(GGD)的拟合程度的度量。

给定模型Ψ中的幂p的

范数：

Ψ可以是一般的变换，例如，小波变换或傅立叶变换。

Ψ可以是学术字典。

接下来解释如何用一个图像先验来解决优化问题。根据

和

的属性，可以使用不同类型的算法。

接下来解释

和

是可微的情况。在这种情况下，通过寻找目标函数

相对于γ的导数的根来解决逆问题，其可以写成

可以看出，这个问题涉及到导数

的计算中H的计算。根据函数

导数的计算也可能涉及测量模型H的伴随算子的计算。这种问题的著名示例是Tikhonov正规化，其中，

and

在这种情况下，我们有如下的解：

其中，

是单位矩阵，N_γ×N_γ方阵在主对角线上有1，在其他地方有0。

接下来解释仅

可微且

是凸的情况。用于解决这类问题的一组流行方法称为投影梯度法。投影梯度法利用了最优化问题的解满足一个固定点方程这一事实。一种流行的算法是向前-向后拆分，涉及以下固定点方程系统，在这个特定的示例中仅限于一个方程：

其中，

表示与

相关的邻近算子，由P.Combettes和J.-C.Pesquet的“Proximal Splitting Methods in Signal Processing”，Fixed-Point Algorithms for Inverse Problems in Science and Engineering，第185-212页，2011年中介绍了邻近算子。上述方程被表示为固定点方程，因为

出现在方程的两边。因此，通过执行以下迭代来解决该问题：

直到达到收敛标准。

另一种流行的算法是快速迭代收缩阈值算法(FISTA)，由A.Beck和M.Teboulle的“A Fast Iterative Shrinkage-Thresholding Algorithm for Linear InverseProblems”，SIAM Journal on Imaging Science，第2卷，第183-202页，2009年中介绍，这是上述方法的变体。

向前-向后拆分的一个示例是

和

时。在这种情况下，执行以下迭代：

其中，

是一个超参数。

实际上，与

范数相关联的邻近算子被定义为

其中，软阈值算子被元素式地定义为：sign(x_i)max(|x_i|-τ,0)，其中，当x_i＞0时，sign(x_i)等于1，并且当x_i＝0以及-1时，sign(x_i)等于0。

下面解释

和

不可微的情况。在这种情况下，可以使用几种方法，例如，由S.Boyd、N.Parikh、E.Chu、B.Peleato和J.Eckstein的“Distributed Optimization andStatistical Learning via the Alternating Direction Method of Multipliers”，Foundations and Trends in Machine Learning，第3卷，第1-122页，2011年介绍的交替方向乘数法(ADMM)，或者由P.Combettes、L.Condat、J.-C.Pesquet和BC Vu的“A forward-backward view of some primal-dual optimization methods in image recovery”Proceedings of the 2014 IEEE International Conference on Image Processing，第4141-4145页，2014年介绍的原始对偶前向-后向(PDFB)算法。

接下来解释如何解决具有M个图像先验的优化问题。在这种情况下，可以认为先验项可以写成：

在这种情况下，可以使用几种算法，例如：由N.Pustelnik、C.Chaux和J.-C.Pesquet的“Parallel Proximal Algorithm for Image Restoration Using HybridRegularization”，IEEE Transactions on Image Processing，第20卷，第2450-2462页，2011年介绍的原始对偶前向-后向并行近似算法(PPXA)或由S.Setzer、G.Steidl和T.Teuber的“Deblurring Poissonian Images by Split Bregman Techniques”，Journalof Visual Communication and Image Representation，第21卷，第193-199页，2010年介绍的同时方向乘数法(SDMM)。这些算法的共同特性是需要计算

(或导数

)和M图像先验项的邻近算子。

的导数或邻近算子的计算可能涉及测量模型的估计和取决于函数

的测量的伴随算子。

上述图像重建方法可以通过图5至图9的流程图来概括。图5的流程图给出了该方法的概述。在步骤21中，一组脉冲波由一组传输元件传输。在这个示例中，传输元件具有有限的尺寸和任意的形状。这组传输脉冲可以只包括一个脉冲或者一个序列中传输的多于一个脉冲。在步骤23中，一组传感器接收一组反射回波波形，在这个示例中，这些传感器具有有限的尺寸和任意的形状。在步骤25中，定义了将这组测量值与这组反射率值相关联的逆问题。在步骤27中，逆问题的解被表示为涉及目标函数的优化问题，该目标函数包括数据差异项和一个或多个图像先验项。在步骤29中，通过使用投影获得反射率值的估计来解决优化问题。

图6的流程图更详细地概述了在函数

可微的情况下提出的方法。前四个步骤(即步骤31、33、35和37)与图6流程图中的前四个步骤相同。在步骤39中，作为初始化，通过在测量值m上估计测量模型的伴随算子来计算第一反射率值估计

在步骤41中，根据图8的流程图计算或估计变量

在步骤43中，计算导数

在某些实现示例中，在步骤45中，伴随算子用于计算数据差异项的导数。例如，当

时，

的导数表示为：

这涉及测量模型的伴随算子的计算。在步骤47中，计算图像先验项的邻近算子。在步骤49中，更新

在步骤50中，确定是否满足收敛标准。满足的话，这个过程结束了。如果不满足标准，则过程在步骤41继续。

图7的流程图更详细地概述了在数据差异项是欧几里德距离并且逆问题包含一个图像先验项的情况下所提出的方法。前五个步骤(即步骤51、53、55、57、59)基本上与图6的流程图中的前五个步骤相同，不同之处在于，现在逆问题包括一个图像先验项。在步骤61中，根据图8的流程图为每个测量网格点计算变量

在步骤63中，计算残差rⁿ＝vⁿ-m。在步骤65中，根据图10的流程图估计残差

上的测量模型的伴随算子。在步骤67中，计算该值

在步骤69中，计算图像先验项

上的投影。在步骤70，确定是否满足收敛标准。满足的话，这个过程结束了。如果不满足标准，则过程在步骤61继续。

图8的流程图解释了估计或处理测量模型的示例过程。在步骤71中，为每个测量网格点(定义任意形状的有限尺寸接收传感器)、对应于测量网格点的任意形状的有限尺寸接收器中的每个点源位置、任意形状的每个有限尺寸传输器以及任意形状的有限尺寸传输器中的每个点源位置，生成一组反射率参数方程。换言之，为每个传输器和每个接收器中的每对点源之间的每个传输器-接收器波路径生成一组反射率参数方程。反射率参数方程定义了一组反射率超曲面。在步骤73中，反射率超曲面离散，以获得一组超曲面反射率样本。在步骤75中，在反射率网格上内插超曲面反射率样本。在步骤77中，对插值和离散的超曲面反射率样本进行求和。在步骤79中，在步骤77中获得的值添加到测量网格点处的测量值估计。在步骤80中，将获得的累积测量值估计与脉冲形状进行卷积。

图9的流程图解释了估计或处理测量模型的伴随算子的示例过程。在步骤81中，测量值估计与脉冲形状的匹配滤波器进行卷积。在步骤83中，为反射率网格的每个点、为任意形状的每个有限尺寸传输器、为任意形状的每个有限尺寸接收器、为任意形状的有限尺寸传输器中的每个点源位置以及为任意形状的有限尺寸接收器中的每个点源位置，生成一组测量参数方程。换言之，为每个传输器和每个接收器中的每对点源之间的每个传输器-接收器路径生成一组测量参数方程。测量参数方程定义了一组测量超曲面。在步骤85中，测量超曲面离散，以获得一组超曲面测量样本。在步骤87中，在测量网格上内插超曲面测量样本。在步骤89中，对插值和离散的超曲面测量样本进行求和。在步骤90中，在步骤89中获得的值添加到反射率网格点处的反射率值估计。

传输的脉搏波可以由至少一个由电信号激励的机电换能装置产生。可以从至少一个互易机电换能装置产生的电信号中获得测量值。机电换能装置在空间上设置成：具有多元件结构的一个或多个线阵探头，所述多元件结构具有基本上一个波长间距，并且沿着直线、凸线或凹线排列；具有多元件结构的一个或多个相控阵探头，所述多元件结构具有基本上半波长间距，并且沿着直线、凸线或凹线排列；具有多元件结构的一个或多个矩阵阵列探头，所述多元件结构具有基本上一个波长间距，并且在平坦表面上排列；或者具有多元件结构的一个或多个矩阵阵列探头，所述多元件结构具有基本上半波长间距，并且在平坦表面上排列。

虽然已经在附图和前面的描述中详细示出和描述本发明，但是这种示出和描述被认为是说明性的或示例性的，而不是限制性的，本发明不限于所公开的实施例。基于对附图、公开内容和所附权利要求的研究，本领域技术人员在实施要求保护的发明时可以理解并实现其他实施例和变型。

例如，本发明的教导可以应用于光声成像，其中，传输声波被电磁辐射脉冲代替。在这种特殊情况下，通过将所有声传输传播延迟设置为零，可以简化本说明书中介绍的整个公式，因为电磁辐射以接近光速的速度传播，即通常比声速快五个数量级。因此，当考虑接收时的声学传播延迟时，可以忽略传输传播延迟。在本发明中考虑的另一修改是使用不同于压电材料的换能器技术，例如，电容式加工超声换能器(CMUT)。另一修改是将本发明应用于时分复用方案，由此传输和/或接收换能器元件的子集依次连续寻址，并且存储所得到的信号，用于随后的重组。所有上述应用也可以转化为本领域已知的基于微泡的对比度特定的非线性成像技术，例如，反向脉冲、振幅调制、振幅调制反向脉冲、谐波成像等。

本发明还提出了一种计算机程序产品，其包括存储在非暂时性介质上的指令，用于当在电子数据处理装置(例如，成像设备1)的计算装置上加载和运行时实现如上所述的方法的步骤。

在权利要求中，单词“包括”不排除其他元件或步骤，不定冠词“a”或“an”不排除多个。在相互不同的从属权利要求中引用不同特征这一事实并不表示不能有利地使用这些特征的组合。

Claims (14)

1.一种脉冲波图像重建方法，包括：

·通过一组传输元件将一组脉冲波传输(21)到要至少部分重建的对象；

·通过一组传感器从对象接收(23)一组反射回波波形，所述反射回波波形的特征在于，在包括反射率网格点的反射率网格上定义的一组反射率值，基于所述反射回波波形确定一组测量值，将测量值定义在包括测量网格点的测量网格上；

·定义(25)逆问题，用于通过测量模型从这组测量值获得一组反射率值估计；

·将所述逆问题表示(27)为涉及包括数据差异项的目标函数的优化问题；

·将优化问题表示(27)为涉及数据差异项贡献的方程组；

·通过在方程组中执行迭代，直到所述方程组迭代获得的解达到给定的收敛阈值，来获得(29)这组反射率值估计，至少一些迭代至少涉及一组测量网格点中每个测量网格点的测量值估计和一组反射率网格点中每个反射率网格点的反射率值估计；所述数据差异项被表示为正泛函，估计所述一组测量网格点中的每个测量网格点的测量值和测量值估计之间的距离；

所述方法的特征在于，获得方程组中的贡献包括：

·从测量模型获得这组测量网格点中的每个测量网格点的测量值估计，包括：

·生成(71)每个传输元件-传感器路径的一组反射率参数方程，其定义了包括超曲面反射率值的一组反射率超曲面；

·在这组反射率超曲面中的每个反射率超曲面上集成(77)超曲面反射率值，以获得各个测量网格点的测量值估计；

·从测量模型的伴随算子中获得这组反射率网格点中的每个反射率网格点的反射率值估计，包括：

·生成(83)每个传输元件-传感器路径的一组测量参数方程，其定义了包括超曲面测量值的一组测量超曲面；

·在这组测量超曲面中的每个测量超曲面上集成(89)超曲面测量值，以获得各个反射率网格点的反射率值估计。

2.根据权利要求1所述的方法，其中，所述一组传输元件中的每个传输元件被认为是一组传输点源，而所述一组传感器中的每个传感器(13)被认为是一组接收点源。

3.根据权利要求1或2所述的方法，其中，在集成超曲面反射率值之前，所述方法还包括：

·离散(73)这组反射率超曲面中的每个反射率超曲面，以获得超曲面反射率值；

·插值(75)离散的超曲面反射率值，以使插值后的所述离散的超曲面反射率值基本上与反射率网格点一致。

4.根据权利要求1所述的方法，其中，在集成超曲面测量值之前，所述方法还包括：

·离散(85)这组测量超曲面中的每个测量超曲面，以获得超曲面测量值；

·插值(87)离散的超曲面测量值，以使插值后的所述离散的超曲面测量值基本上与测量网格点一致。

5.根据权利要求1所述的方法，其中，所述正泛函是正可微泛函。

6.根据权利要求1所述的方法，其中，所述目标函数被表示为所述数据差异项和至少一个图像先验项的组合，所述图像先验项包括关于反射率值的先验信息。

7.根据权利要求6所述的方法，其中，所述至少一个图像先验项被表示为给定模型中的这组反射率值估计的表示的幂p的l_p范数，为向量

定义该l_p范数，作为

其中，

表示实数空间。

8.根据权利要求1所述的方法，其中，所述脉冲波是脉冲声波或电磁辐射波。

9.根据权利要求1所述的方法，其中，所述脉冲波包括一组转向平面波或一组发散波。

10.根据权利要求1所述的方法，其中，由电信号激励的至少一个机电换能装置产生所述脉冲波，并且从由至少一个互易机电换能装置产生的电信号中获得所述测量值。

11.根据权利要求10所述的方法，其中，所述机电换能装置在空间上设置成

·具有多元件结构的线阵探头，所述多元件结构具有基本上一个波长间距，并且沿着直线、凸线或凹线排列；

·具有多元件结构的相控阵探头，所述多元件结构具有基本上半波长间距，并且沿着直线、凸线或凹线排列；

·具有多元件结构的矩阵阵列探头，所述多元件结构具有基本上一个波长间距，并且在平坦表面上排列；或者

·具有多元件结构的矩阵阵列探头，所述多元件结构具有基本上半波长间距，并且在平坦表面上排列。

12.根据权利要求1所述的方法，其中，当所述一组反射率值位于2D空间中时，所述反射率超曲面和所述测量超曲面具有圆锥形状，而当所述一组反射率值位于3D空间中时，所述反射率和所述测量超曲面具有二次曲面形状。

13.根据权利要求1所述的方法，其中，通过数值集成来实现所述超曲面测量或反射率值的集成，并且其中，所述方法还包括：将获得的测量值估计与相应测量网格点处的先前测量值估计相加(79)，以获得累积测量值估计；并且将累积的测量值估计与传输的脉冲波的脉冲形状进行卷积(80)。

14.一种用于重建对象图像的成像设备(1)，所述设备包括用于以下操作的装置：

·通过一组传输元件将一组脉冲波传输到要至少部分重建的对象；

·通过一组传感器从对象接收一组反射回波波形，所述反射回波波形的特征在于在包括反射率网格点的反射率网格上定义的一组反射率值，基于所述反射回波波形确定一组测量值，将测量值定义在包括测量网格点的测量网格上；

·定义逆问题，用于通过测量模型从这组测量值获得一组反射率值估计；

·将所述逆问题表示为涉及包括数据差异项的目标函数的优化问题；

·将优化问题表示为涉及数据差异项贡献的方程组；

·通过在方程组中执行迭代，直到所述方程组迭代获得的解达到给定的收敛阈值，来获得这组反射率值估计，至少一些迭代至少涉及一组测量网格点中每个测量网格点的测量值估计和一组反射率网格点中每个反射率网格点的反射率值估计；所述数据差异项被表示为正泛函，估计所述一组测量网格点中的每个测量网格点的测量值和测量值估计之间的距离；

所述设备的特征在于，为了获得方程组中的贡献，所述设备包括用于以下操作的装置：

·从测量模型获得这组测量网格点中的每个测量网格点的测量值估计，所述测量值估计被配置为通过以下方式获得：

·生成(71)每个传输元件-传感器路径的一组反射率参数方程，其定义了包括超曲面反射率值的一组反射率超曲面；

·在这组反射率超曲面中的每个反射率超曲面上集成(77)超曲面反射率值，以获得各个测量网格点的测量值估计；

·从测量模型的伴随算子中获得这组反射率网格点中的每个反射率网格点的反射率值估计，所述反射率值估计被配置为通过以下方式获得：

·生成(83)每个传输元件-传感器路径的一组测量参数方程，其定义了包括超曲面测量值的一组测量超曲面；

·在这组测量超曲面中的每个测量超曲面上集成(89)超曲面测量值，以获得各个反射率网格点的反射率值估计。

Images