CN101926657B - 一种超声图像特征追踪方法及其系统 - Google Patents

Abstract

本发明公开了一种超声图像特征追踪方法及其系统，其包括：设定初始的感兴趣特征区域；计算区域灰度相似性，获得区域灰度相似性参量；读取所述感兴趣特征区域内的多普勒速度信息，建立约束条件；结合所述区域灰度相似性参量和约束条件，构造相似性度量，计算所述相似性度量的极值，获得追踪的轨迹点位置。本发明结合了彩色组织多普勒成像的特点，有效地提高了特征区域追踪的准确性和鲁棒性。

Description

一种超声图像特征追踪方法及其系统

技术领域

本发明涉及组织多普勒成像中图像特征区域的追踪技术领域，具体涉及一种图像特征追踪技术设计方法和系统。

背景技术

M模式(Motion Mode)超声主要用于心脏及大血管检查。M模式超声沿声束方向将心肌运动回波信号的强度值沿时间轴展开，构成距离时间曲线图。由于它所描记的是声束所经心脏各层组织结构的运动轨迹，因此又称为M模式超声心动图(M-Ultrasound Cardiogram&Echocardiogram)。根据瓣膜的形态、厚度、反射强弱、活动速度等改变，利用M模式超声心动图可确诊二尖瓣狭窄、瓣膜赘生物、腱索断裂、心肌肥厚等病变。它与心电图配合使用，是测定多项心功能指标的重要工具。

随着组织多普勒成像(TDI)技术的方展，将应用多普勒效应得到的心肌或血管壁等组织运动的多普勒频移信号进行彩色编码，代替心肌或血管壁的回波强度信号，得到心肌等组织运动随时间变化的曲线图，称为组织多普勒M模式成像。

由于在组织多普勒M模式成像中，当每一组新的时间采样数据产生后，就可以计算得到一组新的多普勒频移信息，因此M模式是所有多普勒组织成像模式中时间分辨率最高的一种方式，在心肌定量分析方面得到了广泛应用。另一方面，组织多普勒M模式成像可以显示不同时相上室壁运动的梯度和不同层次肌纤维的方向和速度，也可以确定心肌相位速度分布，同时，位移，加速度，能量，应变，应变率等多种运动参数都可以用彩色编码的形式显示出来，极大的扩展了应用的范围，使得组织多普勒M模式成像在心脏功能检测中起到了越来越重要的作用。

为了对心肌进行准确的定量分析，需要对M模式图像中某些特征点(心肌特定位置)进行追踪或者边界描迹，从而能够准确的描绘出对应特征点的速度、位移等运动参数随时间变化的曲线，为临床提供更好的诊断依据。这些特征点包括心肌室壁，二尖瓣等临床诊断比较关注的心肌位置。

由于普通M模式超声心动图一样，组织多普勒M模式成像仅仅是一条取样线上的成像，因此必须将取样线放置在有临床意义的部位。这就要求，在进行M模式成像时，B图像的质量应当比较高，使得心脏不同结构各部分轮廓比较清晰。

通常如图1所示，一种典型的彩色多普勒M模式超声系统的结构框图，其包括：超声探头101向人体相应检查部位(如心脏)发射超声波102，成像扫描可以通过ECG控制器103触发控制；接收到的回波信号经过波束合成105、正交解调106等信号处理环节，分别送入非多普勒信号处理模块107和多普勒信号处理模块108，通过非多普勒信号处理模块107计算，获得人体组织(如心脏)的解剖结构灰度图像(例如B图像和M图像)，通过多普勒信号处理模块108计算，获得人体组织(如心脏)的彩色多普勒运动图像，并保存到电影回放数据存储单元109。其中，图像计算单元111读取图像数据，包括灰度图像和多普勒速度信息，进行特征区域的追踪计算，然后根据得到的特征区域位置信息，读取对应的多普勒速度信息，如速度等，进行心肌各项功能参数计算分析。得到的位置信息和功能参数送到显示器进行显示112。此超声系统允许用户手工控制，即进入操作控制环节110，选择设定初始特征区域，调整追踪计算结果，选择最终计算显示的图像参数等。

图像计算单元111在进行图像特征区域的追踪时，常用的是基于图像灰度信息的相似性度量进行追踪，获得的心肌短轴的M图像。用户在设定感兴趣特征区域后，根据灰度特性来追踪感兴趣特征区域随时间变化的运动情况。在M模式图像特征点追踪中，单纯利用相邻区域灰度相似性进行目标追踪，容易受到灰度噪声干扰，使计算结果的鲁棒性也容易受到影响，并且，由于心肌运动情况比较复杂，容易出现错误，使得计算结果出现误差，易受到多普勒速度信息中噪声信息的干扰，从而影响特征区域追踪的准确性和鲁棒性。

可见，现有技术中存在一定的缺陷，需要进一步地改进。

发明内容

本发明的目的在于提供一种超声图像特征追踪方法及其系统，其结合了彩色组织多普勒成像的特点，有效地提高了特征区域追踪的准确性和鲁棒性。

为了实现上述目的，本发明采用如下技术方案：

本发明提供了一种超声图像特征追踪方法，其包括：

设定初始的感兴趣特征区域；

计算区域灰度相似性，获得区域灰度相似性参量；

读取所述感兴趣特征区域内的多普勒速度信息，建立约束条件；

结合所述区域灰度相似性参量和约束条件，构造相似性度量，计算所述相似性度量的极值，获得追踪的轨迹点位置。

基于上述方法，本发明还提供了一种超声系统，包括：图像特征计算单元和参数计算模块，所述图像特征计算单元包括：灰度相似性计算模块，用于对输入的图像数据，计算区域灰度相似性，获得区域灰度相似性参量；约束条件建立单元，用于读取所述感兴趣特征区域内的多普勒速度信息，建立约束条件；相似性度量构造单元，用于结合所述区域灰度相似性参量和约束条件，构造相似性度量；追踪计算单元，用于计算所述相似性度量的极值，获得追踪的轨迹点位置，依此结果沿时间轴逐线在深度方向上建立感兴趣特征区域的运动轨迹，并将运动轨迹的最终结果输入至所述参数计算模块。

另外，本发明还提供了一种彩色多普勒M模式感兴趣特征区域追踪装置，包括：用于存储或读写超声图像数据的存储及读写模块、追踪检测模块、以及参数计算模块；所述追踪检测模块包括：灰度相似性计算模块，用于对输入的图像数据，计算区域灰度相似性，获得区域灰度相似性参量；约束条件建立单元，用于读取所述感兴趣特征区域内的多普勒速度信息，建立约束条件；相似性度量构造单元，用于结合所述区域灰度相似性参量和约束条件，构造相似性度量；追踪计算单元，用于计算所述相似性度量的极值，获得追踪的轨迹点位置，依此结果沿时间轴逐线在深度方向上建立感兴趣特征区域的运动轨迹，并将运动轨迹的最终结果输入至所述参数计算模块。

可见，本发明提出一种新的彩色多普勒M模式图像追踪方法，与现有技术不同，本发明综合利用M模式图像的灰度相似性信息和多普勒速度信息进行相似性度量实现图像追踪，在考虑了灰度特性的基础上，同时利用多普勒速度信息构成另一个追踪条件，从而约束灰度相似性构成的相似性度量，有效地提高了追踪的准确性，具有更强的鲁棒性。

附图说明

图1为彩色多普勒M模式超声系统的结构框图；

图2为本发明彩色多普勒M模式图像特征追踪技术算法流程图；

图3为利用本发明追踪技术后的心肌长轴切面M模式图像及追踪结果示意图；

图4为本发明超声系统中图像计算单元的结构示意图；

图5为图像计算单元中追踪计算单元的结构示意图。

具体实施方式

本发明在现有彩色多普勒M模式超声系统的基础上，给出了一种新的超声图像特征追踪方法，其包括以下步骤：

设定初始的感兴趣特征区域；

计算区域灰度相似性，获得区域灰度相似性参量；

读取所述感兴趣特征区域内的多普勒速度信息，建立约束条件；

结合所述区域灰度相似性参量和约束条件，构造相似性度量，计算所述相似性度量的极值，获得追踪的轨迹点位置。

根据上述过程的描述可知，本发明同时利用区域灰度相似性和多普勒速度信息构造相似性度量，在区域灰度相似性的基础上考虑多普勒速度信息对相似性度量可能产生的影响，从而提高了追踪结果的精确性。

在上述方法的基础上，为了更加提高追踪的准确性和鲁棒性，本发明在上述相似性度量中增加的约束条件采用位移差，该位移差为利用多普勒速度信息产生的位移变量与待追踪区域的位移变量之差。本发明在区域灰度相似信息的基础上，增加多普勒速度信息产生的位移变量特征与待追踪区域的位移变量之间的相似性作为另一追踪条件，考了角度纠正等因素，更加提高了追踪结果的准确性，增强了鲁棒性，提高了计算过程的稳定性；相对于现有技术中，本发明将灰度信息与多普勒速度信息有机的结合在一起，分别保留了多普勒速度信息和灰度信息进行图像追踪的优点，而且还对彩色多普勒M型的时间分辨率要求不用很高，并且不容易受到多普勒速度信息中噪声信息的干扰，避免计算结果出现偏差。

在上述方法的基础上，区域灰度相似性参量可以包括：沿时间轴逐线在深度方向上进行相邻区域灰度相似性计算，获得的第一灰度相似性；那么，相似性度量可以采用以下形式：

第一灰度相似性与一单调递减函数之积，该单调递减函数以所述位移差的绝对值为变量。

为了提高追踪结果的准确性和鲁棒性，也还可以在此第一灰度相似性基础上增加第二灰度相似性，即沿时间轴逐线在深度方向上将待追踪区域和所述初始感兴趣特征区域进行区域灰度相似性计算。那么，相似性度量可以采用以下形式：

所述第一灰度相似性与第二灰度相似性的和，与一单调递减函数之积，该单调递减函数以所述位移差的绝对值为变量。

本发明采用第一灰度相似性和第二灰度相似性相结合的方式时，可同时利用追踪过程中图像相邻特征点特征信息的相似性与待追踪区域和初始图像特征点特征信息相似性两个相似信息来进行追踪，相对于现有技术而言，本发明考虑了初始区域的灰度特性，使计算结果的鲁棒性不容易受误差的影响，有效提高了追踪结果的准确性和鲁棒性。

以下将结合附图详细描述本发明的一最佳实现方式，此最佳实现方式，即考虑了初始区域的灰度特性、还考虑了多普勒速度信息的影响。

本发明采用图1所示的系统对心肌进行超声图像的采样，并且获得相应的图像数据，以下是图像数据的处理过程。如图2所示：

步骤200，输入目标超声图像，包含多普勒图像和非多普勒图像。在此步骤中，选择合适的心脏切面进行扫描，调整探头角度及图像参数获得清晰的彩色多普勒M模式图像；冻结图像，选取图像质量较好的图像进行计算。此目标超声图像可以是彩色多普勒M型图像，包括彩色血流多普勒M型图像和彩色组织多普勒M型图像。

步骤201，对图像进行相应的减噪滤波处理。此属于降噪滤波环节，用于消除噪声对于感兴趣区域内的干扰，剔除无效数据的影响。在此过程中，还可以通过进行低灰度抑制或者速度空间平滑，实现进一步地减噪处理。

步骤202，用户根据需求，设定感兴趣特征区域。此感兴趣区域可以为已包含心肌的矩形区域，或者用户根据需要设定的其他形状，比如线段、圆形或者平行四边形。并且感兴趣特征区域的选择可以手动，或者自动并且允许手动调整感兴趣特征区域的大小。

步骤203，沿时间轴逐线在深度方向上进行相邻区域灰度相似性计算，获得第一灰度相似性；

步骤204，沿时间轴逐线在深度方向上将待追踪区域和最初设定的感兴趣区域进行区域灰度相似性计算，获得第二灰度相似性；

步骤205，读取多普勒速度信息，添加速度约束条件，利用计算得到的位移进行相似性计算。这里的约束条件是关于位移差的，该位移差为利用多普勒速度信息产生的位移变量与待追踪区域的位移变量之差。上述步骤203、204以及205可以并行计算，不限制其先后顺序。

步骤206，结合所述第一灰度相似性、第二灰度相似性和速度约束条件，构造相似性度量，计算相似性度量的极值，获得追踪的轨迹点位置；

步骤207，依据获得的轨迹点位置，沿时间轴逐线在深度方向上建立感兴趣特征区域的运动轨迹。此过程中主要包括以下步骤：

根据计算获得轨迹点位置选取两个相邻的轨迹点，沿时间轴逐线搜索此相邻轨迹点之间的正向运动轨迹和反向运动轨迹；对所述正向运动轨迹和反向运动轨迹进行加权计算，获得所述感兴趣特征区域的运动轨迹。通过搜索正向运动轨迹和反向运动轨迹，利用它们的相关性获得最终的运动轨迹，可以提高计算准确度。

另外，为了减少计算误差，本发明在此步骤207的过程中添加计算检验步骤，比如：对相邻两次计算结果进行阈值判断，具体过程如下：

首先，计算相邻两次所述感兴趣特征区域运动轨迹的结果差异；

然后，判断所述结果差异是否超过预设阈值，若是，则修正本次运动轨迹的计算结果，利用前两次计算结果进行修正，用以调整计算；若否，则继续执行下一个步骤。

步骤208，调整计算过程，判断是否需要手工调整；若是，则对于计算不够理想的部分，进行手工调整，然后执行下一个步骤209；否则，直接执行步骤209。

步骤209，根据所述运动轨迹提取感兴趣特征区域的多普勒速度信息进行运动参数计算，包括提取速度，方差等多普勒速度信息进行位移，速度梯度、应变，应变率等运动参数计算；

步骤210，根据追踪结果分析生成参数曲线；

步骤211，显示参数曲线。

以上各部分技术细节如下所述：

用户根据需要，在TDI的M模式图像上设定一个或多个大小为N×T像素(N为深度方向，T为时间轴方向)的感兴趣特征区域R₀作为初始感兴趣特征区域位置，并且沿时间轴间断地点取目标位置(通常为心动周期内的特征时刻)，最后一个点取的位置为终止位置(至少点取两个点)。在初始位置与终止位置之间进行图像追踪和运动参数计算，并且每次计算在选择的两个相邻的轨迹点之间进行。

上述步骤203和步骤204中均涉及到区域灰度相似性计算，本发明在这里提供三种实现方式，但本发明决不限于此三种方式，其他相关性的计算方法在本发明中一样适用。

1、上述两个区域灰度分布间的相似性度量可以采用下述公式：

${\mathrm{Cor}}_{i,i+1}=\frac{\underset{k=0}{\overset{N}{\mathrm{\Σ}}}(x_k-\stackrel{\‾}{x})(y_k-\stackrel{\‾}{y})}{{\left(\left(\underset{k=0}{\overset{N}{\mathrm{\Σ}}}{(x_k-\stackrel{\‾}{x})}^2\right)\left(\underset{k=0}{\overset{N}{\mathrm{\Σ}}}{(y_k-\stackrel{\‾}{y})}^2\right)\right)}^{1/2}}---\left(1\right)$

其中，Cor_i，i+1表示沿时间轴逐线搜索时第i、i+1时刻对应得计算区域R_i、R_i+1区域内灰度的相似性；x_k、y_k表示第i、i+1时刻对应的计算区域R_i、R_i+1内对应点的灰度，x_k∈R_i，y_k∈R_i+1，k＝0…N；x，y分别表示区域R_i、R_i+1内灰度的平均值。上述公式(1)给出的灰度分布相似性度量去除了灰度差异对于相似性的影响，并且有一定的抗噪声能力。

2、上述两个区域灰度分布间的相似性度量可以采用下述公式：

${\mathrm{Cor}}_{i,i+1}=\underset{a,b}{\mathrm{\Σ}}{P}_{i,i+1}(a,b)\lg \left(\frac{P_{i,i+1}(a,b)}{P_i\left(a\right)P_{i+1}\left(b\right)}\right)---\left(2\right)$

其中，Cor_i，i+1表示沿时间轴逐线搜索时第i、i+1时刻对应得计算区域R_i，R_i+1区域内灰度的相似性；P_i(a)、P_i+1(b)分别表示计算区域R_i、R_i+1基于灰度直方图的图像像素灰度概率分布；P_i，i+1(a，b)表示基于计算区域R_i、R_i+1的二维灰度直方图的联合概率分布。

3、上述两个区域灰度分布间的相似性度量可以采用下述公式：

${\mathrm{Cor}}_{i,i+1}=\underset{k=0}{\overset{N}{\mathrm{\Σ}}}|x_k-y_k|---\left(3\right)$

其中，Cor_i，i+1表示沿时间轴逐线搜索时第i、i+1时刻对应得计算区域R_i、R_i+1区域内灰度的相似性；x_k、y_k表示第i、i+1时刻对应的计算区域R_i、R_i+1内对应点的灰度。

采用上述三个公式计算的第一灰度相似性即Cor_i，i+1，第二灰度相似性即Cor_0，i+1。上述步骤206中，结合所述第一灰度相似性Cor_i，i+1、第二灰度相似性Cor_0，i+1和约束条件，构造相似性度量，此相似性度量可以为：上述第一灰度相似性与第二灰度相似性的和，与以步骤205中位移差的绝对值为变量的单调递减函数之积。这里的单调递减函数可以为：倒数、指数函数或者幂函数等等。例如：

1、以所述位移差的绝对值与一系数的乘积加上1之和的倒数，作为相似性度量的约束条件。首先读取初始感兴趣区域R₀内的平均速度v₀，以及扫描线间时间隔常数T。在已知第i时刻位置p_i，搜索第i+1时刻位置p_i+1时，可基于如下相似性度量：

E＝(Cor_i，i+1+Cor_0，i+1)/(1+a·|p_i+v_i·T-p_i+1|)    (4)

其中，Cor_i，i+1表示沿时间轴逐线搜索时第i、i+1时刻对应得计算区域R_i，R_i+1区域内灰度的相似性；a≥0为权重系数；在公式(4)中，分子项中的Cor_i，i+1表示目标区域与已知参考区域灰度的相似性，能够约束相邻时刻位置的连续性，Cor_0，i+1表示目标区域与初始感兴趣区域R₀灰度分布的相似性。通过结合上述两个相似性度量项，可以同时保证搜索位置前后的连续性和当前目标与初始感兴趣区域的相似性，从而保证目标跟踪的正确性。公式(4)中分母项表示利用多普勒速度信息对于追踪结果的约束，利用速度信息对于位移判断的估计来约束追踪结果，系数a越大，速度信息对追踪结果的约束越强，即在相似性度量中引入了基于组织多普勒速度信息对目标运动的弹性约束，其中v_i·T表示利用多普勒速度信息产生的位移量。

2、以上述位移差的绝对值与一系数乘积的指数函数，作为相似性度量的约束条件。

$E=({\mathrm{Cor}}_{i,i+1}+{\mathrm{Cor}}_{0,i+1})e^{-a\·|p_i+v_i\·T-p_{i+1}|}---\left(5\right)$

3、以上述位移差的绝对值与一常数之和的幂函数，作为相似性度量的约束条件，其中α＝1时为公式(4)。

E＝(Cor_i，i+1+Cor_0，i+1)*(1+|p_i+v_i·T-p_i+1|)^-a    (6)

上述公式(4)至公式(6)可统一表示为：

E＝(Cor_i，i+1+Cor_0，i+1)f(|p_i+v_i·T-p_i+1|)，其中f(x)＞0，并且单调递减。

上述步骤207的运动轨迹计算过程中，在实际搜索时，为了提高计算的准确度，在手动选取的两个相邻的轨迹点之间沿时间轴分别搜索计算正向运动轨迹p_i，i＝0…M和反向运动轨迹q_i，i＝0…M，M表示两轨迹点之间的线数，即搜索时间范围。对于两次计算的结果进行加权平均，公式如(7)所示：

P_i＝[1-w(i)]p_i+w(i)q_i    (7)

将加权计算得到的结果P_i作为最终的运动轨迹。其中w(i)＝i/M。该加权方法，在靠近起始点时p_i的权重大，靠近结束点时q_i的权重大，即考虑跟踪所得位置靠近越初始位置可信度越高。

另外，为了消除计算结果因为噪声引入的偶发计算错误，对于第i次的结果P_i与第i-1次的结果P_i-1进行比较，如果两次结果的差值的绝对值超出设定阈值Th，则认为本次计算错误，需要利用前两次计算的结果对当前结果进行修正，修正公式如(4)所示加权计算：

P_i＝P_i-1+β*(P_i-1-P_i-2)    (8)

其中，β为权重系数。通过这样的阈值判断保证，可以消除偶发噪声造成的位置追踪计算错误，提高追踪计算的准确性和鲁棒性。

在上述步骤209中，对于自动追踪计算得到的区域，用户可以根据结果判断是否需要进行手动调整，对于个别不理想的结果，可以手动调整追踪结果，得到理想的感兴趣区域运动轨迹曲线。如图3所示，利用本发明可以在组织多普勒M模式图像中比较准确的获取心肌室壁和室间隔运动的轨迹曲线，其中，虚线部分表示追踪结果。

从图2所示的流程可以看出，本发明对心肌中各部分组织灰度分布的差异，利用TDI-M模式的高时间分辨率，基于多普勒速度信息和灰度信息一起构成相似性度量，来进行感兴趣特征区域的追踪，能有效克服噪声的影响。获取感兴趣区域的运动轨迹曲线以后，沿运动轨迹曲线提取对应位置上的多普勒速度信息，进行感兴趣区域运动参数的计算，获取速度、速度梯度、位移等精确的运动参数来进行心肌运动的定量分析测量，为临床心脏疾病的诊断提供便捷有力的工具。本发明也可以同时选择多个感兴趣特征区域，同时利用本发明进行运动轨迹追踪，从而得到进行多个感兴趣区域的运动参数曲线，比较不同位置心肌在不同时相的定量运动信息。

上述方法可以通过软件编程，或者进行硬件改进在原有超声系统的基础上实现相应功能，基于此点，本发明主要通过技术手段在图1所示的图像计算单元111实现内部性能和构成或功能上的改变。如图4所示，图像计算单元111主要包括：图像特征计算单元300，和参数计算模块400，其中图像特征计算单元300包括：

灰度相似性计算模块301，用于对输入的图像数据，计算区域灰度相似性，获得区域灰度相似性参量；

约束条件建立单元302，用于读取所述感兴趣特征区域内的多普勒速度信息，建立约束条件；

相似性度量构造单元304，用于结合所述区域灰度相似性参量和约束条件，构造相似性度量；

追踪计算单元305，用于计算所述相似性度量的极值，获得追踪的轨迹点位置，依此结果沿时间轴逐线在深度方向上建立感兴趣特征区域的运动轨迹，并将运动轨迹的最终结果输入至所述参数计算模块400。

如图4所示，基于上述结构，约束条件建立单元301包括：

速度信息计算单元311，用于计算多普勒速度信息产生的位移变量与待追踪区域的位移变量之差，构造关于位移差的约束条件，并将结果添加到所述相似性度量构造单元304中。

如图4所示，基于上述结构，灰度相似性计算模块302包括：

第一灰度相似性计算单元312，用于对输入的图像数据，沿时间轴逐线在深度方向上进行相邻区域灰度相似性计算，获得第一灰度相似性，并将结果添加到所述相似性度量构造单元304中。

为了考虑初始特征区域灰度特性的影响，灰度相似性计算模块302中还包括：

第二灰度相似性计算单元312，用于对输入的图像数据，沿时间轴逐线在深度方向上将待追踪区域和所述初始感兴趣特征区域进行区域灰度相似性计算，获得第二灰度相似性，并将结果添加到所述相似性度量构造单元304中。

如图5所示，上述追踪计算单元305包括：

追踪单元365，用于计算所述相似性度量的极值，获得追踪的轨迹点位置；

正向运动轨迹计算单元315，用于根据计算获得轨迹点位置选取两个相邻的轨迹点，沿时间轴逐线搜索此相邻轨迹点之间的正向运动轨迹；

反向运动轨迹计算单元325，根据计算获得轨迹点位置选取两个相邻的轨迹点，沿时间轴逐线搜索此相邻轨迹点之间的反向运动轨迹；

加权计算单元335，用于对所述正向运动轨迹和反向运动轨迹进行加权计算，获得所述感兴趣特征区域的运动轨迹；

阈值判断单元345，用于计算相邻两次所述感兴趣特征区域运动轨迹的结果差异，并判断所述结果差异是否超过预设阈值，在所述结果差异超过预设阈值时启动修正单元，在所述结果差异不超过预设阈值时直接控制将所述加权计算单元335输入至所述参数计算模块400；

修正单元355，用于利用先前的计算结果修正本次计算结果，并将该结果输入至所述参数计算模块400。

上述改进是基于现有的超声系统，本发明还针对产业链生产中环节，单独给出了一种彩色多普勒M模式感兴趣特征区域追踪装置，其包括：用于存储或读写超声图像数据的存储及读写模块、追踪检测模块、参数计算模块以及显示模块；并且，所述追踪检测模块与上述超声系统中图像特征计算单元300具有相同的结构组成，也包括灰度相似性计算模块、约束条件建立单元、相似性度量构造单元以及追踪计算单元，对应的相关功能模块内部的结构组成参见图4与图5的相关说明。这里的彩色多普勒M模式感兴趣特征区域追踪装置可以通过USB等数据端口与现有的超声系统相连，进行实时数据的采集和处理及显示，也可以从现有的超声系统获取相关数据后，进行单独的数据处理。本发明彩色多普勒M模式感兴趣特征区域追踪装置可以安装在超声系统内部，也可以外置，使用方便，结构简单，并且改装成本低，也不需要对现有的超声系统硬件进行改装。

综上所述，本发明主要应用于彩色多普勒M模式，根据彩色多普勒M模式超声图像的特点，提出了一种新的感兴趣特征区域运动追踪技术，对于用户选择的感兴趣特征区域，采用基于相邻区域灰度分布相似性和与初始位置区域灰度分布相似性组合构造相似性度量，同时利用多普勒速度信息对该相似性度量进行有效约束，来准确的检测感兴趣区域的运动轨迹。在检测过程中，利用正反向计算加权平均和阈值判断修正的措施，进一步保证了计算结果的鲁棒性。另外，本发明对手工控制点之间数据进行正反向同时追踪，对追踪结果进行加权处理。对搜索区域，相邻线之间的误差进行了限制和修正，进一步提高了追踪算法的稳定性和准确性。

上述各具体步骤的举例说明较为具体，并不能因此而认为是对本发明的专利保护范围的限制，本发明还可以拓展到彩色多普勒成像模式，如TDI，CFM等模式，将搜索范围从一维扩展到两维甚至三维空间上来，上述涉及的计算公式不需要改变，公式中的小标i，k等表示的是时间关系，并不是位置坐标。本发明的主要应用是一维彩色M型成像，在选择图像特征时，可以选择一线或者多线上一段区域来进行特征追踪，扩展到如TDI等两维成像时，选择一定区域内，如5×5的成像区域内的数据作为特征点来追踪即可；而对于三维TDI等心肌成像时，选择区域可以扩展到5×5×5的一个立体区域来作为特征区域进行追踪，在原理上只是从一维空间扩展到三维空间。本发明的专利保护范围应以所附权利要求为准。

Claims (14)

1.一种超声图像特征追踪方法，其特征在于，包括： 

设定初始的感兴趣特征区域； 

计算区域灰度相似性，获得区域灰度相似性参量； 

读取所述感兴趣特征区域内的多普勒速度信息，建立约束条件； 

结合所述区域灰度相似性参量和约束条件，构造相似性度量，计算所述相似性度量的极值，获得追踪的轨迹点位置； 

所述约束条件为：计算多普勒速度信息产生的位移变量与待追踪区域的位移变量之差，所获得的位移差； 

所述区域灰度相似性参量包括：沿时间轴逐线在深度方向上进行相邻区域灰度相似性计算，获得的第一灰度相似性； 

所述相似性度量为：所述第一灰度相似性与一单调递减函数之积，该单调递减函数以所述位移差的绝对值为变量。 

2.根据权利要求1所述的方法，其特征在于，所述方法中，所述区域灰度相似性参量还包括：沿时间轴逐线在深度方向上将待追踪区域和所述初始感兴趣特征区域进行区域灰度相似性计算，获得的第二灰度相似性。 

3.根据权利要求2所述的方法，其特征在于，所述相似性度量为：所述第一灰度相似性与第二灰度相似性的和，与一单调递减函数之积，该单调递减函数以所述位移差的绝对值为变量。 

4.根据权利要求1所述的方法，其特征在于，所述区域灰度相似性的计算采用下述公式： 

其中，Cor_i,i+1表示沿时间轴逐线搜索时第i、i+1时刻对应得计算区域R_i、R_i+1区域内灰度的相似性；x_k、y_k表示第i、i+1时刻对应的计算区域R_i、R_i+1 内对应点的灰度，x_k∈R_i,y_k∈R_i+1,k＝0...N；

分别表示区域R_i、R_i+1内灰度的平均值。 

5.根据权利要求1所述的方法，其特征在于，所述区域灰度相似性的计算采用下述公式： 

其中，Cor_i,i+1表示沿时间轴逐线搜索时第i、i+1时刻对应得计算区域R_i,R_i+1区域内灰度的相似性；P_i(a)、P_i+1(b)分别表示计算区域R_i、R_i+1基于灰度直方图的图像像素灰度概率分布；P_i,i+1(a,b)表示基于计算区域R_i、R_i+1的二维灰度直方图的联合概率分布。 

6.根据权利要求1所述的方法，其特征在于，所述区域灰度相似性的计算采用下述公式： 

其中，Cor_i,i+1表示沿时间轴逐线搜索时第i、i+1时刻对应得计算区域R_i、R_i+1区域内灰度的相似性；x_k、y_k表示第i、i+1时刻对应的计算区域R_i、R_i+1内对应点的灰度。 

7.根据权利要求1所述的方法，其特征在于，还包括： 

依据获得的轨迹点位置，沿时间轴逐线在深度方向上建立感兴趣特征区域的运动轨迹； 

根据所述运动轨迹提取感兴趣特征区域的多普勒速度信息进行运动参数计算。 

8.根据权利要求7所述的方法，其特征在于，所述建立运动轨迹的过程包括： 

根据计算获得轨迹点位置选取两个相邻的轨迹点，沿时间轴逐线搜索此相邻轨迹点之间的正向运动轨迹和反向运动轨迹； 

对所述正向运动轨迹和反向运动轨迹进行加权计算，获得所述感兴趣 特征区域的运动轨迹。 

9.根据权利要求8所述的方法，其特征在于，所述建立运动轨迹的过程中，还包括阈值判断步骤： 

计算相邻两次所述感兴趣特征区域运动轨迹的结果差异； 

判断所述结果差异是否超过预设阈值，若是，则修正本次所述感兴趣特征区域运动轨迹的结果；若否，则继续执行所述建立运动轨迹的过程的下一个步骤。 

10.一种彩色多普勒M模式感兴趣特征区域追踪装置，包括：用于存储或读写超声图像数据的存储及读写模块、追踪检测模块、以及参数计算模块；其特征在于，所述追踪检测模块包括： 

灰度相似性计算模块，用于对输入的图像数据，计算区域灰度相似性，获得区域灰度相似性参量； 

约束条件建立单元，用于读取所述感兴趣特征区域内的多普勒速度信息，建立约束条件； 

相似性度量构造单元，用于结合所述区域灰度相似性参量和约束条件，构造相似性度量； 

追踪计算单元，用于计算所述相似性度量的极值，获得追踪的轨迹点位置，依此结果沿时间轴逐线在深度方向上建立感兴趣特征区域的运动轨迹，并将运动轨迹的最终结果输入至所述参数计算模块； 

所述约束条件建立单元包括： 

速度信息计算单元，用于计算多普勒速度信息产生的位移变量与待追踪区域的位移变量之差，获得位移差，构造关于位移差的约束条件，并将结果添加到所述相似性度量构造单元中； 

所述灰度相似性计算模块包括： 

第一灰度相似性计算单元，用于对输入的图像数据，沿时间轴逐线在深度方向上进行相邻区域灰度相似性计算，获得第一灰度相似性，并将结果添加到所述相似性度量构造单元中； 

所述相似性度量构造单元构造的相似性度量为：所述第一灰度相似性与一单调递减函数之积，该单调递减函数以所述位移差的绝对值为变量。 

11.根据权利要求10所述的装置，其特征在于，所述灰度相似性计算模块还包括： 

第二灰度相似性计算单元，用于对输入的图像数据，沿时间轴逐线在深度方向上将待追踪区域和所述感兴趣特征区域进行区域灰度相似性计算，获得第二灰度相似性，并将结果添加到所述相似性度量构造单元中； 

所述相似性度量构造单元构造的相似性度量为：所述第一灰度相似性与第二灰度相似性的和，与一单调递减函数之积，该单调递减函数以所述位移差的绝对值为变量。 

12.根据权利要求11所述的装置，其特征在于，所述追踪计算单元包括： 

追踪单元，用于计算所述相似性度量的极值，获得追踪的轨迹点位置； 

正向运动轨迹计算单元，用于根据计算获得轨迹点位置选取两个相邻的轨迹点，沿时间轴逐线搜索此相邻轨迹点之间的正向运动轨迹； 

反向运动轨迹计算单元，根据计算获得轨迹点位置选取两个相邻的轨迹点，沿时间轴逐线搜索此相邻轨迹点之间的反向运动轨迹； 

加权计算单元，用于对所述正向运动轨迹和反向运动轨迹进行加权计算，获得所述感兴趣特征区域的运动轨迹； 

阈值判断单元，用于计算相邻两次所述感兴趣特征区域运动轨迹的结果差异，并判断所述结果差异是否超过预设阈值，在所述结果差异超过预设阈值时启动修正单元，在所述结果差异不超过预设阈值时直接控制将所述加权计算单元输入至所述参数计算模块； 

修正单元，用于利用先前的计算结果修正本次计算结果，并将该结果 输入至所述参数计算模块。 

13.一种超声系统，包括：图像特征计算单元和参数计算模块，其特征在于，所述图像特征计算单元包括： 

灰度相似性计算模块，用于对输入的图像数据，计算区域灰度相似性，获得区域灰度相似性参量；所述灰度相似性计算模块包括： 

第一灰度相似性计算单元，用于对输入的图像数据，沿时间轴逐线在深度方向上进行相邻区域灰度相似性计算，获得第一灰度相似性，并将结果添加到所述相似性度量构造单元中;

约束条件建立单元，用于读取感兴趣特征区域内的多普勒速度信息，建立约束条件，所述感兴趣特征区域为设定的特征区域； 

相似性度量构造单元，用于结合所述区域灰度相似性参量和约束条件，构造相似性度量； 

追踪计算单元，用于计算所述相似性度量的极值，获得追踪的轨迹点位置，依此结果沿时间轴逐线在深度方向上建立感兴趣特征区域的运动轨迹，并将运动轨迹的最终结果输入至所述参数计算模块； 

所述约束条件建立单元包括：速度信息计算单元，用于计算多普勒速度信息产生的位移变量与待追踪区域的位移变量之差，获得位移差，构造关于位移差的约束条件，并将结果添加到所述相似性度量构造单元中； 

所述灰度相似性计算模块包括：第一灰度相似性计算单元，用于对输入的图像数据，沿时间轴逐线在深度方向上进行相邻区域灰度相似性计算，获得第一灰度相似性，并将结果添加到所述相似性度量构造单元中； 

所述相似性度量构造单元构造的相似性度量为：所述第一灰度相似性与一单调递减函数之积，该单调递减函数以所述位移差的绝对值为变量。 

14.根据权利要求13所述的超声系统，其特征在于，所述灰度相似性计算模块还包括： 

第二灰度相似性计算单元，用于对输入的图像数据，沿时间轴逐线在深度方向上将待追踪区域和所述感兴趣特征区域进行区域灰度相似性计算，获得第二灰度相似性，并将结果添加到所述相似性度量构造单元中； 

所述相似性度量构造单元构造的相似性度量为：所述第一灰度相似性与第二灰度相似性的和，与一单调递减函数之积，该单调递减函数以所述位移差的绝对值为变量。 

Images