CN1319492C - 一种变尺度的生物组织位移估计方法 - Google Patents

Abstract

本发明涉及一种变尺度的生物组织位移估计方法，属于超声弹性成像技术领域。本发明的方法为：第k条扫描线数据采用的尺度为T(k)；且如果M为奇数，T(k)满足T(1)＝T(M)＝T_max，T((M+1)/2)＝T_min，当1≤k≤(M+1)/2时，T(k)随着k的增大而减小，当(M+1)/2≤k≤M时，T(k)随着k的增大而增大；如果M为偶数，T(k)满足T(1)＝T(M)＝T_max，T(M/2)＝T(M/2+1)＝T_min，当1≤k≤M/2时，T(k)随着k的增大而减小，M/2+1≤k≤M时，T(k)随着k的增大而增大；其中，T_max和T_min分别为变尺度方法的最大尺度和最小尺度。本发明可在保证了足够的组织位移估计精度的前提下，提高了纵向分辨率；同时减小组织横向位移对组织纵向位移估计的影响。

Description

一种变尺度的生物组织位移估计方法

技术领域

本发明属于超声弹性成像技术领域，特别涉及生物组织位移估计方法。

背景技术

生物组织弹性模量的变化通常与其病理现象有关。例如，恶性的病理损害，例如乳房硬癌、前列腺癌、甲状腺癌及肝转移等，通常表现为硬的小结。乳房硬癌是乳腺癌的最常见形式，大约占乳腺癌总数的四分之三，由于其基质密度增大而表现为致密的硬块。而其他类型的乳腺癌如导管内癌和乳头状瘤则表现为柔软的组织，良性的乳腺纤维囊性病也很少表现为硬块。

生物组织的弹性模量信息对于疾病的诊断过程具有重要的参考价值。然而，包括X射线成像、超声成像、计算机断层成像(CT)和磁共振成像(MRI)等在内的传统医学成像模态都不能直接提供关于弹性模量这一组织的基本力学属性的信息。1991年，J.Ophir提出超声弹性成像(ultrasound elastography)的方法，对组织的弹性模量分布进行定量估计、成像。目前，超声弹性模量已经成为医学超声成像的一个研究热点，广泛应用于乳房、前列腺、动脉粥样斑块、心肌动力学以及高强度聚焦超声与射频消融引起的损害(lesion)的检测与评估。

超声弹性成像的基本原理为：将超声探头嵌于一块挤压平板中，沿着探头的纵向压缩组织，分别采集组织压缩前、后的射频信号；组织被压缩时，组织内将会产生一个沿压缩方向的应变，如果组织内部弹性模量分布不均匀，组织内的应变分布也会有所差异；弹性模量较大的区域，引起的应变比较小：反之，弹性模量较小的区域，相应的应变比较大。通过一些方法估计出组织内部不同位置的位移，从而计算出组织内部的应变分布情况，用来间接描述组织内部的弹性模量分布，从而描述组织的生理、病理状态。

对于二维超声弹性成像，一般采用线阵的B型超声探头，采集组织压缩前、后的探头每一条扫描线的射频信号，分别进行上面描述的位移估计，从而计算出每一条扫描线对应组织的一维应变分布。最后把所有扫描线对应的一维应变分布按扫描线顺序组成一个二维应变分布，以灰度图或者伪彩图的形式表示，用来间接描述组织内部的弹性模量分布。

一般的超声弹性成像方法包括以下步骤：

1.利用商用B型超声仪器(一般采用线阵探头)得到待测生物组织(一般为人体组织，也可以为动物组织，以下简称组织)压缩前的一幅数字化的二维射频信号(可以采用模拟射频信号输出端接信号放大器，再接高速数据采集卡，得到数字化的二维射频信号；也可以在数字化B型超声仪器上直接得到数字化的二维射频信号)；

2.手持该B型超声仪器的探头或者利用步进电机或者螺旋装置驱动该探头，沿着探头的纵向对该组织施加一个微小的挤压(组织压缩量一般控制在为1％的数量级)，得到组织压缩后的一幅数字化的二维射频信号；

3.从步骤1和2的得到的组织压缩前、后的二维射频信号中分别取出第一条扫描线的数据，设为s₁(n)和s₂(n)，n表示该两条扫描线上的数据序号，1≤n≤n_max，n的最大值n_max由该B型超声仪器的探查深度、发射的超声波在组织中的传播速度以及射频信号的采样频率决定；

4.从该扫描线数据s₁(n)中取一小段尺度为T的数据d₁，其数据个数为U，U＝round(T×U₁)，其中，T的单位为mm，U₁代表1mm的组织对应的数据个数，由发射的超声波在组织中的传播速度以及射频信号的采样频率决定，round(·)代表四舍五入的取整操作，该数据d₁的序号从n₁到n₁+U-1，n₁可在1≤n₁≤U的范围内选择；在τ₁到τ₂确定的搜索范围内求该小段数据与扫描线数据s₂(n)的互相关函数R(τ)，计算公式如下

$R\left(\mathrm{\τ}\right)=\frac{\underset{m=n_1}{\overset{n_1+U-1}{\mathrm{\Σ}}}{s}_1\left(m\right)s_2(m-\mathrm{\τ})}{\sqrt{\underset{m=n_1}{\overset{n_1+U-1}{\mathrm{\Σ}}}{s}_1^2\left(m\right)\·\underset{m=n_1}{\overset{n_1+U-1}{\mathrm{\Σ}}}{s}_2^2(m-\mathrm{\τ})}}({\mathrm{\τ}}_1\≤\mathrm{\τ}\≤{\mathrm{\τ}}_2)$

其中m为计算过程中表示数据序号的循环变量，τ₁为0，τ₂为对组织施加的压缩量，以采样数据的个数表示；(为了提高位移估计的精度，一般还需要对计算得到互相关函数进行插值，如抛物线插值)；

5.确定步骤4得到的互相关函数R(τ)的最大值对应的位置t₁，t₁就是数据d₁在组织压缩后的位移(即s₁(n)中的序号从n₁到n₁+U-1的小段数据d₁的在组织压缩后移动到s₂(n)中的序号从n₁-t₁到n₁+U-1-t₁的位置)；

6.依次从扫描线数据s₁(n)中取一小段尺度为T即数据个数为U的数据d₂、d₃、…、d_N，每段数据的序号依次错开V个采样数据(如V＝round(0.4×T×U₀)、V＝round(0.5×T×U₀)等)，直到再错开V个采样数据将超出s₁(n)的范围，按步骤4、5相同的方法依次得到各段数据对应的位移t₂、t₃、…、t_N，其中N为小段数据的总数；则位移序列t₁、t₂、…、t_N为第一条扫描线数据s₁(n)对应的组织的位移估计；

7.利用与步骤3-6相同的方法，依次得到第2、3、…、M条扫描线数据对应的组织的位移估计，其中M为表示探头的扫描线总数，由探头决定；

8.对第一条扫描线数据s₁(n)对应的组织的位移估计序列t₁、t₂、…、t_N求差分，得到组织第一条扫描线s₁(n)对应组织的应变分布，计算公式如下，

${\mathrm{\ϵ}}_1=\frac{t_2-t_1}{V},$ ${\mathrm{\ϵ}}_2=\frac{t_3-t_2}{V},\·\·\·,$ ${\mathrm{\ϵ}}_{N-1}=\frac{t_N-t_{N-1}}{V}$

其中，ε₁、ε₂、…、ε_N-1分别为d₁、d₂、…、d_N-1对应的组织应变；

9.利用与步骤8相同的方法，依次得到组织第2、3、…、M条扫描线数据对应的组织的应变分布；

10.将步骤9得到的M条扫描线数据对应的应变分布，按照扫描线的顺序组合成一个二维数据，并以灰度图或者伪彩图的形式表示出来，就得到组织的二维应变分布图。

在超声弹性成像中，关键的问题在于对组织的位移分布进行估计，也就是上面描述的方法的步骤3-7。互相关函数的值越大，说明压缩前、后的小段数据吻合得越好，互相关函数的最大值位置代表了压缩前的小段数据在压缩后对应的位置，从而可以求出该小段数据的位移，也就是该小段数据对应的组织的位移。

超声弹性成像中，对组织施加一个小的压缩量，利用互相关分析估计的组织位移是纵向位移，即沿压缩方向的位移。但是，对组织施加一个小的压缩量的时候，组织的运动是很复杂的，受到组织内部弹性模量分布、组织的几何形状、边界条件等因素的影响。组织不仅沿着压缩方向(即线阵探头的纵向)有一个压缩，沿着垂直于压缩方向的方向(包括线阵探头的横向以及垂直于探头扫描平面的方向)也会膨胀。研究表明，沿着垂直于压缩方向的位移会引起压缩前、后信号的互相关函数的幅度降低，也就是说，压缩前、后信号的相似性降低。而超声弹性成像正是利用压缩前、后信号的相似性来跟踪组织位移的，所以，沿着垂直于压缩方向的位移会使得组织位移估计的精度降低。并且，垂直于探头扫描平面方向的位移的影响比横向位移的影响小，所以，减少沿着垂直于压缩方向的位移的影响，主要是减小垂直于探头扫描平面的方向的位移影响。

为了减小垂直于探头扫描平面的方向的位移影响，有人提出用平板在探头横向限制组织位移的方法，避免组织在探头横向上的移动，虽然同时增大了垂直于探头扫描平面方向的位移，但是由于该方向的位移对组织位移估计的影响较小，因此该方法可以减少沿着垂直于压缩方向的位移的影响。但是，该方法主要应用于实验研究，而在实际应用中，很难对组织的横向移动加以限制。另外，从理论分析可以知道，减小组织压缩量、利用距离探头中心扫描线附近的超声射频数据进行弹性成像、增大扫描线宽度、提高探头中心频率或者增大带宽，都可以在一定程度上减少沿着垂直于压缩方向的位移的影响。然而，这些方法在实际上也不可行。

在超声弹性成像中，选择的压缩前信号中的小段数据用来跟踪对应的小段组织的位移，称为跟踪波段。其长度称为跟踪波段长度，或者尺度。在上述方法中，尺度T为常数，即对所有的扫描线数据的位移估计中，尺度T是恒定不变的。选取合适的尺度在超声弹性成像的位移估计中非常重要。当组织压缩量比较小的时候，尺度越大，包含的信息越多，位移估计的精度越高，对随机噪声的干扰也越不敏感。然而，增大尺度将降低超声弹性成像的纵向分辨率。因为超声弹性成像中，一般假定跟踪波段内的应变均匀。因此，随尺度的增大，对跟踪波段内的应变起到了平均作用。

发明内容

本发明为了减小垂直于探头扫描平面的方向的位移影响，提出一种变尺度的组织位移估计方法，可在保证了足够的组织位移估计精度的前提下，提高了纵向分辨率；同时减小组织横向位移对组织纵向位移估计的影响。

本发明提出的一种变尺度的组织位移估计方法，包括以下步骤：

1.从组织压缩前、后的二维射频信号中分别取出第一条扫描线的数据，设为s₁(n)和s₂(n)，n表示该两条扫描线上的数据序号，1≤n≤n_max，n的最大值n_max由该B型超声仪器的探查深度、发射的超声波在组织中的传播速度以及射频信号的采样频率决定；

2.从该扫描描线数据s₁(n)中取一小段尺度为T(1)的数据d₁，其数据个数为U，U＝round(T(1)×U₁)，其中，T的单位为mm，U₁代表1mm的组织对应的数据个数，由发射的超声波在组织中的传播速度以及射频信号的采样频率决定，round(·)代表四舍五入的取整操作，T(1)为尺度随扫描线序号k变化的函数T(k)在k＝1时的取值，该数据d₁的序号从n₁到n₁+U-1，n₁可在1≤n₁≤U的范围内选择；在τ₁到τ₂确定的搜索范围内计算该小段数据与扫描线数据s₂(n)的互相关函数R(τ)，计算公式如下

$R\left(\mathrm{\τ}\right)=\frac{\underset{m=n_1}{\overset{n_1+U-1}{\mathrm{\Σ}}}{s}_1\left(m\right)s_2(m-\mathrm{\τ})}{\sqrt{\underset{m=n_1}{\overset{n_1+U-1}{\mathrm{\Σ}}}{s}_1^2\left(m\right)\·\underset{m=n_1}{\overset{n_1+U-1}{\mathrm{\Σ}}}{s}_2^2(m-\mathrm{\τ})}}({\mathrm{\τ}}_1\≤\mathrm{\τ}\≤{\mathrm{\τ}}_2)$

其中m为计算过程中表示数据序号的循环变量，τ₁为0，τ₂为对组织施加的压缩量，以采样数据的个数表示；(为了提高位移估计的精度，一般还需要对计算得到互相关函数进行插值，如抛物线插值)；

3.确定该互相关函数R(τ)的最大值对应的位置t₁，t₁就是数据d₁在组织压缩后的位移(即s₁(n)中序号从n₁到n₁+U-1的小段数据d₁在组织压缩后移动到s₂(n)中的序号从n₁-t₁到n₁+U-1-t₁的位置)；

4.依次从该扫描线数据s₁(n)中取一小段尺度为T(1)即数据个数为U的数据d₂、d₃、…、d_N，每段数据的序号依次错开V个采样数据(如V＝round(0.4×T(1)×U₀)、V＝round(0.5×T(1)×U₀)等)，直到再错开V个采样数据将超出s₁(n)的范围，按步骤2、3相同的方法依次得到各段数据对应的位移t₂、t₃、…、t_N，其中N为小段数据的总数；则位移序列t₁、t₂、…、t_N为第一条扫描线数据s₁(n)对应的组织的位移估计；

5.利用与步骤1-4相同的方法，依次得到第2、3、…、M条扫描线数据对应的组织的位移估计，其中M为由探头决定的扫描线总数；

其特征在于，所述的第1、2、…、M条扫描线数据采用的尺度分别为T(1)，T(2)，…，T(M)，即第k条扫描线数据采用的尺度为T(k)；且如果M为奇数，T(k)满足T(1)＝T(M)＝Tmax， $T\left(\frac{M+1}{2}\right)=T\min ,$ 当1≤k≤(M+1)/2时，T(k)随着k的增大而减小，当(M+1)/2≤k≤M时，T(k)随着k的增大而增大；如果M为偶数，T(k)满足T(1)＝T(M)＝Tmax，T(M/2)＝T(M/2+1)＝Tmin，当1≤k≤M/2时，T(k)随着k的增大而减小，M/2+1≤k≤M时，T(k)随着k的增大而增大；其中，T_max和T_min分别为变尺度方法的最大尺度和最小尺度。

在上述方法中，探头中央部位的扫描线数据(即M为奇数时的第(M+1)/2条扫描线数据，M为偶数时的第M/2条和第M/2+1条扫描线数据)采用最小尺度T_min，探头两边缘部位的扫描线数据(即第1条扫描线和第M条扫描线)采用最大尺度T_max；通常情况下，可以选择取0.5-2mm作为最小尺度的范围(该范围的尺度在一般的超声弹性成像组织位移估计方法中经常使用)；而最大尺度还与探头的宽度和对组织压缩量关，例如探头宽度为60mm，对组织压缩量为1％的时候，可以选择最大尺度为8-10mm。如果探头宽度和对组织压缩量变化，最小尺度可以不变，而最大尺度做适当调整。

上述不同扫描线数据选择的尺度随扫描线序号k变化的方式T(k)，可为简单的线性变化方式(如果M为奇数，1≤k≤(M+1)/2时，T(k)随着k的增大而线性地减小，(M+1)/2≤k≤M时，T(k)随着k的增大而线性地增大；如果M为偶数，1≤k≤M/2时，T(k)随着k的增大而线性地减小，M/2+1≤k≤M时，T(k)随着k的增大而线性地增大)，该线性的尺度变化方式表示如下：

如果M为奇数，则

$T\left(k\right)=\left\{\begin{array}{cc}T\max -\frac{T\max -T\min }{(M+1)/2-1}\&times;(k-1)& (1\&le;k\&le;(M+1)/2)\\ T\min +\frac{T\max -T\min }{M-(M+1)/2}\&times;(k-(M+1)/2)& ((M+1)/2<k\&le;M)\end{array}\right.$

如果M为偶数，则

$T\left(k\right)=\left\{\begin{array}{cc}T\max -\frac{T\max -T\min }{M/2-1}\&times;(k-1)& (1\&le;k\&le;M/2)\\ T\min +\frac{T\max -T\min }{M-(M/2+1)}\&times;(k-(M/2+1))& (M/2+1<k\&le;M)\end{array}\right.$

不同扫描线数据选择的尺度随扫描线序号k变化的方式T(k)，也可为其他的变化方式，只要该变化方式满足本发明的特征，即如果M为奇数，T(k)满足T(1)＝T(M)＝Tmax， $T\left(\frac{M+1}{2}\right)=T\min ,$ 并且1≤k≤(M+1)/2时，T(k)随着k的增大而减小，(M+1)/2≤k≤M时，T(k)随着k的增大而增大；如果M为偶数，T(k)满足T(1)＝T(M)＝Tmax，T(M/2)＝T(M/2+1)＝Tmin，并且1≤k≤M/2时，T(k)随着k的增大而减小，M/2+1≤k≤M时，T(k)随着k的增大而增大。

本发明的原理：

本发明提出的一种变尺度的组织位移估计方法，在不同扫描线的数据(即组织不同的横向位置对应的数据)采用不同的跟踪波段长度即尺度的方法。在探头中央部位的扫描线数据采用最小尺度，估计组织的位移分布，从而得到应变分布，在保证了足够的组织位移估计精度的前提下，提高了纵向分辨率。而在探头两边缘部位的扫描线数据采用最大尺度，估计组织的位移分布以及应变分布，从而减小组织横向位移对组织纵向位移估计的影响。在从探头中央部位到探头各边缘部位过渡的扫描线数据，尺度逐渐增大，从最小尺度变化到最大尺度。可以采用简单的线性变化方式，也可以采用其他的变化方式。

虽然采用较小的尺度误差较大，但是在探头中央部位，横向位移的影响较小，所以即使采用较小的尺度也能够获得足够精确的位移估计，同时又提高了纵向分辨率。而在探头两边缘部位，则是牺牲一定的纵向分辨率来换取组织纵向位移的精度，因为在探头两边缘部位，组织横向位移对纵向位移估计的影响较大。

组织横向位移的影响主要取决于组织压缩量大小、扫描线宽度、扫描线位置以及探头的中心频率和带宽等参数。因此，本发明提出的变尺度的生物组织位移估计方法在实际应用中，最大尺度和最小尺度的选择也可能有所不同，这与探头的宽度(或探头两边缘部位与探头中央部位的距离)、组织压缩量大小、以及探头的中心频率和带宽等有关。

在特定的参数条件下，最大尺度和最小尺度的选择主要取决于在探头两边缘部位和探头中央部位分别要取得什么样的成像效果，需要在保证精度和保证分辨率中间综合考虑。因为探头中央部位的组织位移估计精度较高，因此主要考虑提高分辨率；而在探头两边缘部位主要考虑提高估计精度。

关于组织压缩量和其他参数条件下，如何选择合适的尺度，已经有丰富的研究结果。这些资料可以用来作为最小尺度和最大尺度选择的参考。另外，具体实施的时候可以先用计算机仿真来研究在特定的参数条件下，不同尺度进行组织位移估计的效果，包括估计精度和分辨率。通常情况下，超声弹性成像中可以选择取尺度0.5-2mm作为最小尺度的范围。而最大尺度还与探头的宽度和组织压缩量有关。例如探头宽度为60mm，组织压缩量为1％的时候，可以选择最大尺度范围为8-10mm。如果探头宽度和组织压缩量变化，最小尺度可以不变，而最大尺度做适当调整)，即如果探头宽度增大到两倍即120mm，或者组织压缩量增大到两倍即2％，最大尺度可以选择为原来的两倍即16-20mm。不同扫描线数据选择的尺度随扫描线序号变化的方式，可为简单的线性变化方式，也可为其他满足发明特征的变化方式。

本发明的特点：

在不同扫描线的数据(即组织不同的横向位置对应的数据)采用不同的跟踪波段长度即尺度的；

对横向位移引入误差较小的探头中央部位的扫描线数据，采用最小尺度来提高分辨率；

对横向位移引入误差较大的探头两边缘部位的扫描线数据，则采用最大尺度以减小组织位移估计误差。

从探头中央部位到探头各边缘部位过渡的扫描线数据，相应的尺度逐渐增大，从最小尺度变化到最大尺度。

探头中央部位到探头各边缘部位的扫描线数据选取的尺度的变化方式，可以采用简单的线性变化方式，也可以采用其他的变化方式。

附图说明

图1为本实施例的计算机仿真的组织模型；

图2为利用有限元分析计算得到的组织应变分布的理想结果；

图3为本实施例采用的尺度变化的简单的线性方式的示意图；

图4为本实施例采用的第二种尺度变化方式的示意图；

图5为本实施例采用的第三种尺度变化方式的示意图；

图6为尺度为1mm的一般方法估计得到的组织应变分布的计算机仿真结果；

图7为尺度为3mm的一般方法估计得到的组织应变分布的计算机仿真结果；

图8为尺度为9mm的一般方法估计得到的组织应变分布的计算机仿真结果；

图9为本实施例的变尺度的生物组织位移估计方法估计得到的组织应变分布的计算机仿真结果；

具体实施方式

本发明提出的变尺度的组织位移估计方法结合具体实施例及附图详细说明如下：

本实施例利用计算机程序和一般的超声散射模型仿真得到一块模拟的组织在压缩前和压缩后的二维射频信号。模拟的组织结构如图1，组织大小60×60mm²，组织内分布有5块弹性模量较大的组织层11、12、13、14和15，它们的弹性模量是周围组织的2倍，它们的厚度分别为10mm、5mm、2.5mm、1.25mm和0.625mm，它们之间的间隔16、17、18和19分别为5mm、2.5mm、1.25mm和0.625mm；组织压缩量为1％，即对组织施加的压缩量为0.6mm；探头中心频率为3.5MHz，-3dB带宽为2.0MHz，探头扫描线宽度和间隔分别为2mm和0.4mm，探头宽度和组织宽度一致，也是60mm，因此总共有151条扫描线，即M＝151，并且探头中央部位和探头边缘部位分别对应组织中央部位和组织边缘部位；射频信号的采样频率为20MHz，假设超声波在组织内的传播速度为1540m/s，因此1mm的组织长度对应1mm/(1540×10³mm/s×1/20×10⁶/2)Hz≈26个数据，因为组织深度为60mm，所以，每一条扫描线的数据为60×26＝1560个，对组织施加的压缩量以采样数据的个数表示为60×1％×26≈16个采样数据。

图2表示利用美国MSC公司的MARC软件进行有限元分析，计算得到的该实施例采用的组织模型的理想的应变分布。横向和纵向分别表示超声探头的横向位置(即扫描线的位置)和纵向位置(即组织深度)，灰度表示计算出来的理想应变的大小，灰度值越大(即颜色越亮或越白)，表示应变越大，灰度值越小(即颜色越暗或越黑)，表示应变越小，21为灰度值与应变大小的对照关系。图2中，较暗的区域(即22-26)与弹性模量较大的组织层(即图1中的11-15)对应，说明弹性模量较大的区域应变较小。

本实施例采用的最小尺度确定为1mm，采用的最大尺度确定为9mm。

本实施例的具体步骤如下：

1.从组织压缩前、后的二维射频信号(计算机仿真得到)中分别取出第一条扫描线的数据，设为s₁(n)和s₂(n)，n表示该两条扫描线上的数据序号，1≤n≤1560；

2.从该扫描描线数据s₁(n)中取一小段尺度为T(1)的数据d₁，T(1)＝9mm，其数据个数为U，U＝234，该数据的序号从13到246；在0到16的搜索范围内求该小段数据与扫描线数据s₂(n)的互相关函数R(τ)，计算公式如下

$R\left(\mathrm{\&tau;}\right)=\frac{\underset{m=n_1}{\overset{n_1+U-1}{\mathrm{\&Sigma;}}}{s}_1\left(m\right)s_2(m-\mathrm{\&tau;})}{\sqrt{\underset{m=n_1}{\overset{n_1+U-1}{\mathrm{\&Sigma;}}}{s}_1^2\left(m\right)\&CenterDot;\underset{m=n_1}{\overset{n_1+U-1}{\mathrm{\&Sigma;}}}{s}_2^2(m-\mathrm{\&tau;})}}(0\&le;\mathrm{\&tau;}\&le;16)$

3.确定该互相关函数R(τ)的最大值对应的位置t₁，t₁就是数据d₁在组织压缩后的位移(即s₁(n)中序号从13到246的小段数据d₁在组织压缩后移动到s₂(n)中的序号从13-t₁到246-t₁的位置)；

4.依次从该扫描线数据s₁(n)中取一小段尺度为T(1)即数据个数为U的数据d₂、d₃、…、d_N，每段数据的序号依次错开94个采样数据(V＝round(0.4×T(1)×U₀))，直到再错开94个采样数据将超出s₁(n)的范围，按步骤2、3相同的方法依次得到各段数据对应的位移t₂、t₃、…、t_N，其中N为小段数据的总数，N＝17；则位移序列t₁、t₂、…、t_N为第一条扫描线数据s₁(n)对应的组织的位移估计；

5.利用与步骤1-4相同的方法，依次得到第2、3、…、151条扫描线数据对应的组织的位移估计；采用的尺度从探头中央部位对应的第76条扫描线数据的1mm增大到探头两边缘部位对应的第1条和第151条扫描线数据的9mm。

本实施例采用的尺度变化选择简单的线性变化方式，如图3描述，探头中央部位的第76条扫描线(31)的数据对应的尺度为1mm，探头两边缘部位的第1条和第151条扫描线(32、33)的数据对应的尺度为9mm；从探头中央部位逐渐到探头各边缘部位过渡的扫描线(34、35)，相应的尺度从1mm线性地增大到9mm。该线性的尺度变化方式表示如下

$T\left(k\right)=\left\{\begin{array}{cc}-\frac{8}{75}\&times;k+\frac{683}{75}& (1\&le;k\&le;76)\\ \frac{8}{75}\&times;k-\frac{533}{75}& (76<k\&le;151)\end{array}\right.$

其中，1≤k≤76时的尺度变化对应图3中的34，76＜k≤151时的尺度变化对应图3中的35，34和35都是直线。

图4和图5是本实施例中采用的另外两种尺度变化方式。

图4描述的尺度变化方式表示如下

$T\left(k\right)=\left\{\begin{array}{c}0.00096284667713\&times;k^2-0.18080586080586\&times;k+9.17984301412872\\ (1\&le;k\&le;76)\\ 0.00096284667713\&times;k^2-0.11189952904239\&times;k+3.94296180010464\\ (76<k\&le;151)\end{array}\right.$

其中，1≤k≤76时的尺度变化对应图4中的41，76＜k≤151时的尺度变化对应图4中的42，41和42都是二次多项式曲线(抛物线)。

图5描述的尺度变化方式表示如下

$T\left(k\right)=\left\{\begin{array}{cc}-0.00001912960306\&times;k^3-0.00129124820660\&times;k^2+& \\ 0.09412721186035\&times;k+9.09285509325682& (1\&le;k\&le;76)\\ 0.00001912960306\&times;k^3-0.00743185078910\&times;k^2+& \\ 1.02749880439971\&times;k-42.56097560975238& (76<k\&le;151)\end{array}\right.$

其中，1≤k≤76时的尺度变化对应图5中的51，76＜k≤151时的尺度变化对应图5中的52，51和52都是三次多项式曲线。

尺度变化的方式并不限于图3-5的三种，还可为其他满足发明特征的变化方式。

如果组织压缩量增大，则探头两边缘部位的组织横向位移增大，需要提高最大尺度的大小；而探头中央部位的横向位移影响较小，所以最小尺度可以不变。例如，如果组织压缩量为2.0％，即压缩量增大一倍则最大尺度可以简单地增大一倍，为18mm；而最小尺度仍为1mm。

探头的宽度增大，则探头两边缘部位与探头中央部位的距离增大，也就是探头两边缘部位的扫描线位置增大，因此，探头两边缘部位的组织横向位移也变大，也需要提高最大尺度的大小；相同位置的扫描线对应的组织横向位移不变，对应的尺度不变，因此最小尺度也不变。例如，如果探头的宽度变为120mm(对应的组织宽度也变为120mm)，则最大尺度和最小尺度可分别设为18mm和1mm。

本实施例与一般方法的位移估计效果比较如下：

图6、图7和图8分别表示尺度分别为1mm、3mm和9mm的一般方法估计得到应变分布，图9为本实施例的变尺度的位移估计得到的应变分布。图6-9中，横向和纵向分别表示超声探头的横向位置(即扫描线位置)和纵向位置(即组织深度)，灰度表示估计出来的应变大小，灰度值越大(即颜色越亮或越白)，表示应变越大，灰度值越小(即颜色越暗或越黑)，表示应变越小，61、71、81和91分别为图6、图7、图8和图9的灰度值与应变大小的对照关系。

比较图6-8与图2可见，由于探头两边缘部位组织的横向位移较大，因此，一般的方法在探头两边缘部位(62、63、72、73、82、83)误差较大；比较图7-9还可以看出，增大尺度能够减小组织横向位移对组织位移估计和应变估计的影响，但是也降低了成像的分辨率；而从图9中可以看出，在探头的两边缘部位(92、93)，本实施例采用多尺度的方法可以减小组织横向位移的影响，从而减小两边缘部位(92、93)误差，而在探头的中央部位(94)，横向位移的影响小，采用多尺度的方法保证了足够的纵向分辨率。

Claims (2)

1、一种变尺度的组织位移估计方法，包括以下步骤：

1)从组织压缩前、后的二维射频信号中分别取出第一条扫描线的数据，设为s₁(n)和s₂(n)，n表示该两条扫描线上的数据序号，1≤n≤n_max，n_max表示n的最大值，由该B型超声仪器的探查深度、发射的超声波在组织中的传播速度以及射频信号的采样频率决定；

2)从该扫描描线数据s₁(n)中取一小段尺度为T(1)的数据d₁，其数据个数为U，U＝round(T(1)×U₁)，其中，T(1)为尺度随扫描线序号k变化的函数T(k)在k＝1时的取值，单位为mm，U₁代表1mm的组织对应的数据个数，round(·)代表四舍五入的取整操作，该数据d₁的序号n₁在1≤n₁≤U的范围内选择；在τ₁到τ₂确定的搜索范围内计算该小段数据与扫描线数据s₂(n)的互相关函数R(τ)如下

$R\left(\mathrm{\&tau;}\right)=\frac{\underset{m=n_1}{\overset{n_1+U-1}{\mathrm{\&Sigma;}}}{s}_1\left(m\right)s_2(m-\mathrm{\&tau;})}{\sqrt{\underset{m=n_1}{\overset{n_1+U-1}{\mathrm{\&Sigma;}}}{s}_1^2\left(m\right)\&CenterDot;\underset{m=n_1}{\overset{n_1+U-1}{\mathrm{\&Sigma;}}}{s}_2^2(m-\mathrm{\&tau;})}}({\mathrm{\&tau;}}_1\&le;\mathrm{\&tau;}\&le;{\mathrm{\&tau;}}_2)$

其中m为计算过程中表示数据序号的循环变量，τ₁为0，τ₂为对组织施加的压缩量，以采样数据的个数表示；

3)确定该互相关函数R(τ)的最大值对应的位置t₁，t₁就是数据d₁在组织压缩后的位移；

4)依次从该扫描线数据s₁(n)中取一小段尺度为T(1)即数据个数为U的数据d₂、d₃、…、d_N，每段数据的序号依次错开V个采样数据，直到再错开V个采样数据将超出s₁(n)的范围，按步骤2、3相同的方法依次得到各段数据对应的位移t₂、t₃、…、t_N，其中N为小段数据的总数；则位移序列t₁、t₂、…、t_N为第一条扫描线数据s₁(n)对应的组织的位移估计；

5)利用与步骤1-4相同的方法，依次得到第2、3、…、M条扫描线数据对应的组织的位移估计，其中M为表示探头的扫描线总数；

其特征在于，第k条扫描线数据采用的尺度为T(k)，且如果M为奇数，T(k)满足 $T\left(1\right)=T\left(M\right)=T\max ,T\left(\frac{M+1}{2}\right)=T\min ,$ 当1≤k≤(M+1)/2时，T(k)随着k的增大而减小，当(M+1)/2≤k≤M时，T(k)随着k的增大而增大；如果M为偶数，T(k)满足T(1)＝T(M)＝Tmax，T(M/2)＝T(M/2+1)＝Tmin，当1≤k≤M/2时，T(k)随着k的增大而减小，M/2+1≤k≤M时，T(k)随着k的增大而增大；其中，T_max和T_min分别为变尺度方法的最大尺度和最小尺度。

2、如权利要求1所述的变尺度的组织位移估计方法，其特征在于，

所述尺度随扫描线序号k变化的方式T(k)，为简单的线性变化方式：

如果M为奇数，则

$T\left(k\right)=\left\{\begin{array}{c}T\max -\frac{T\max -T\min }{(m+1)/2-1}\&times;(k-1),(1\&le;k\&le;(m+1)/2)\\ T\min +\frac{T\max -T\min }{M-(M+1)/2}\&times;(k-(m+1)/2),((M+1)/2<k\&le;M)\end{array}\right.$

如果M为偶数，则

$T\left(k\right)=\left\{\begin{array}{c}T\max -\frac{T\max -T\min }{M/2-1}\&times;(k-1),(1\&le;k\&le;M/2)\\ T\min +\frac{T\max -T\min }{M-(M/2+1)}\&times;(k-(M/2+1)),(M/2+1\&le;k\&le;M)\end{array}\right.$

Images