专利名称：超声测温二值图像模糊跟踪方法
技术领域：
本发明属于医疗影像技术领域。
背景技术：
目前，在医疗领域的许多方面，人们迫切需要实现人体内部温度的无损测量，比如超声聚焦、超低温冷冻肿瘤治疗、全身灌注热疗等。并且，即使是常规的人体深部温度无损监控也将对现有的医疗诊断技术做出重大的贡献。
由于超声波成像技术已相当成功地应用于临床实际，因此超声无损测温被寄予较大期望。超声无损测温技术中，最主要的模型是基于生物离散随机介质的时移模型。超声时移模型是基于检测来自组织区的后向散射超声的回波时移。回波时移是由于局部温度对声速的影响和加热区热膨胀的综合结果。
在超声时移模型下，有公式(1)δθ(z)=c0(z)2(1α(z)-β(z))∂∂z(δt(z))----(1)]]>其中，δθ(z)是深度z处的温度变化；c0(z)是深度z处的初始声速；δt(z)是深度z处的累积时移；k(z)=1α(z)-β(z)]]>是个媒质相关参数，可由实验确定。由此可见，只要能够测出累积时移δt(z)随深度z的变化，就可以依据公式(1)计算出任意深度z的温度变化δθ(z)。
到目前为止，研究人员主要使用互相关方法来得到相邻的两帧超声脉冲回波之间的时移，从而估计出温度变化。但是，这种方法显著的缺点是精度低(3到5个波长)，误差会累积，而且需要的计算量很大。
发明内容
本发明的目的是解决超声无损测温中互相关方法精度低，误差会累积，且计算量大的问题。
本发明提出了一种二值图像模糊跟踪(BIFT，Binary Image Fuzzy Tracking)方法，首先利用超声测温仪中的超声探头发射超声脉冲，经人体生物组织反馈后接收返回的超声回波，其特征在于该方法还包括如下步骤1)超声测温仪的采样频率为超声探头发射的超声脉冲中心频率的N倍，其中N≥8，据此设置跟踪窗的尺寸dx×cy，dx方向为帧内方向，dx取N/2～N，dy方向为帧间方向，dy至少为4；2)利用超声测温仪中的微处理器对采样得到的前三帧超声回波数据进行二值化处理，波峰波谷特征点置1，其它点置0；将所述前三帧的波峰波谷特征点直接设定为已评估特征点，其整数累加时移值设为0，并在波峰波谷特征点附近取5个点作4次多项式曲线拟合，所述4次多项式被微分以确定其最大值从而得到偏移量；3)利用超声测温仪中的微处理器对采样得到的下一帧超声回波数据进行二值化处理，波峰波谷特征点置1，其它点置0；每个波峰波谷特征点设定相应的VMARK值，以当前帧回波的最高峰值VM为标准，按黄金分割法分成5级，当峰值Vk≥VM·(1-α)i(i＝1，....，5)时，VMARK＝i，其中α＝0.618；剩下的波峰波谷特征点则在当前帧最高峰值VM的0.8％之下，令VMARK＝0；4)对当前帧中VMARK≠0的特征点，用跟踪窗进行遍历跟踪，即以当前帧中VMARK≠0的特征点为当前跟踪窗的基点，按VMARK值由低到高处理，在dx×dy方向按权从小到大的顺序搜索已评估特征点，直至搜索到第一个已评估特征点，获取所述基点的整数累加时移值，若搜索不到已评估特征点，则将所述基点的VMARK置0，完成评估过程；其中，在当前跟踪窗内，前一帧中与基点的边相邻的采样点的权为1，离所述权为1的采样点越远，权值越大；每评估完一个基点，就在该基点附近取5个点作4次多项式曲线拟合，所述4次多项式被微分以确定其最大值从而得到偏移量，该偏移量减去搜索到的已评估特征点的偏移量，然后加到该基点的整数累加时移值上，获得总的累加时移值；5)对VMARK＝0的特征点，以步骤4)中已评估特征点为基准，通过线性插值得到这些特征点的累加时移值；6)使用中值滤波器对当前帧中的已评估特征点进行帧内过滤；7)对当前帧中的非特征点，通过帧内线性插值来得到其累加时移值；8)对当前帧所有点的累加时移值进行曲线拟合，差分后得到该帧所有点的局部时移值；9)返回步骤3)进入下一帧超声回波数据的处理，直至超声测温结束。
在本发明所述的步骤4)中，在当前跟踪窗内，前一帧中与基点的边相邻的采样点的权为1，与所述权为1的采样点全相邻的采样点权为2，以此类推，离所述权为1的采样点越远，权值越大。
在本发明所述的步骤4)中，在当前跟踪窗内，前一帧中与基点的边相邻的采样点的权为1，与所述权为1的采样点边相邻的采样点权为2，以此类推，离所述权为1的采样点越远，权值越大。
在本发明所述的步骤6)中所述中值滤波器为3到5点中值滤波器。
在本发明所述的步骤8)中所述曲线拟合采用6～10阶多项式。
利用背景技术中所述的公式(1)，将上述方法测出的局部时移值代入，即可得到所测人体生物组织的温度。
本发明所述方法中提到的特征点，也可以是除波峰波谷之外的其他特征点。
本发明所述方法是一种实时方法，所述二值数据一帧帧连起来从总体上看就是二值图像，因此称为二值图像模糊跟踪(BIFT，Binary Image Fuzzy Tracking)方法，该方法精度可以达到1/4波长，误差不会累积，并且计算量也大幅度下降，优于互相关方法。本发明提出的BIFT方法是一种独立于超声测温的方法，对非超声测温的其他类似图像都能有效跟踪，即该方法可以应用在其他需要图像跟踪的领域。
图1为本发明所述方法的流程框图；图2为跟踪窗示意图；图3为波峰波谷特征点在跟踪窗内的三种情况分布图；图4为跟踪窗的权图；图5为4次多项式曲线拟合的示意图；图6为某一帧回波的时移累加值、其拟合曲线，以及差分后得到的局部时移图；图7为实现1/4波长分辨率的示意图。
具体实施例方式
下面结合具体实例和附图来进一步说明本发明。
如图1所示，超声测温二值图像模糊跟踪方法，首先利用超声测温仪中的超声探头发射超声脉冲，经人体生物组织反馈后接收返回的超声回波，其特征在于该方法还包括如下步骤1)超声测温仪的采样频率为超声探头发射的超声脉冲中心频率的N倍，其中N≥8，这里N取8，据此设置跟踪窗的尺寸dx×dy，dx方向为帧内方向，dx取N/2～N，这样各个跟踪窗在dx方向就不会重叠，既可以防止一些误判，也最大限度地减低了噪声的影响，这里dx取7。dy方向为帧间方向，dy取得越大，则计算量越大，但抗噪性越好，这里dy取4；如图2所示，在dx×dy的二维格子中，横向表示帧内相邻的dx个采样点，纵向表示相邻的dy帧。斜线格表示帧内的波峰波谷特征点，点格表示帧内其它点。
2)由于温度变化过程中的回波时移很直观的反映在其波峰波谷上，所以直接使用回波的波峰波谷作为特征点进行跟踪。利用超声测温仪中的微处理器对采样得到的前三帧超声回波数据进行二值化处理，波峰波谷特征点置1，其它点置0；将所述前三帧的波峰波谷特征点直接设定为已评估特征点，其整数累加时移值设为0，并在波峰波谷特征点附近取5个点作4次多项式曲线拟合，所述4次多项式被微分以确定其最大值从而得到偏移量，所述偏移量即为总的累加时移值。
3)利用超声测温仪中的微处理器对采样得到的下一帧超声回波数据进行二值化处理，波峰波谷特征点置1，其它点置0；每个波峰波谷特征点设定相应的VMARK值，VMARK表示该点的可跟踪度，以当前帧回波的最高峰值VM为标准，按黄金分割法分成5级，当峰值Vk≥VM·(1-α)i(i＝1，...，5)时，VMARK＝i，其中α＝0.618；剩下的波峰波谷特征点则在当前帧最高峰值VM的0.8％之下，不列入跟踪范围，令VMARK＝0。
4)对当前帧中VMARK≠0的特征点，用跟踪窗进行遍历跟踪，即以当前帧中VMARK≠0的特征点为当前跟踪窗的基点，按VMARK值由低到高处理，在dx×dy方向按权从小到大的顺序搜索已评估特征点，直至搜索到第一个已评估特征点，获取所述基点的整数累加时移值，若搜索不到已评估特征点，则将所述基点的VMARK置0，完成评估过程；每评估完一个基点得到的时移是整数，因为真正的波峰波谷特征点并不一定在采样点上，所以为了得到更精确的估计，需要在当前跟踪窗基点的附近取5个点作4次多项式曲线拟合，如图5所示，所述4次多项式被微分以确定其最大值从而得到偏移量，该偏移量减去搜索到的已评估特征点的偏移量，然后加到该基点的整数累加时移值上，获得总的累加时移值。
步骤4)得到整数累加时移值的具体分析如下波峰波谷特征点在dy方向的分布有三种基本情况。一种是一直都有，如图3(a)；第二种是从有到无，如图3(b)；最后一种是从无到有，如图3(c)；其他情况可由这三种情况组合得到。对于第一种情况，除了可以直接得到时移评估外，还可以作为断裂波(即第三种情况)插值的依据。第二种情况可以不考虑，因为当前没有需要评估的特征点。对于第三种情况，则有依据不同的情况有不同的解决方法。当当前跟踪窗的当前处理帧之前有至少一个已评估的波峰波谷特征点时，可以用跟踪窗解决；否则，将所述当前跟踪窗基点的VMARK置0。
步骤4)中所述跟踪窗内的权的计算方法如图4所示。在当前跟踪窗内，前一帧中与基点的边相邻的采样点的权为1，与所述权为1的采样点全相邻的采样点权为2，以此类推，离所述权为1的采样点越远，权值越大，见图4(a)。在当前跟踪窗内，前一帧中与基点的边相邻的采样点的权为1，与所述权为1的采样点边相邻的采样点权为2，以此类推，离所述权为1的采样点越远，权值越大，见图4(b)。
权的这种设计既能处理帧向连续点，也能处理帧向断裂点。
5)以步骤4)中已评估特征点为基准，通过线性插值得到VMARK＝0的特征点的累加时移值。
6)使用3到5点中值滤波器对当前帧中的已评估特征点进行帧内过滤。这个过滤器移去了可能出现在位移数据内的任何尖峰。
7)对当前帧中的非特征点，通过帧内线性插值来得到其累加时移值。
8)对当前帧所有点的累加时移值使用p阶(典型的，p＝6～10)多项式来对数据进行曲线拟合。这一步是必须的，因为位移评估中的小误差在差分中可能被转化为大误差。因为得到的时移是累计值，所以差分之后才能得到该帧的局部时移值，图6是单点加热时，某一帧回波的时移累加值、其拟合曲线，以及差分后得到的局部时移。
9)返回步骤3)进入下一帧超声回波数据的处理，直至超声测温结束。
如图7所示，提高特征点的密度可以提高本发明所述方法的分辨率。考虑波峰、波谷与零点作为图像匹配特征量，则每一个周期得到四个特征点，可以轻松的得到1/4波长的分辨率。
实验证明，本发明所述方法能有效地估计出回波累计时移和差分时移。与互相关方法相比较，本发明所述方法在空间、时间和温度的分辨率上都有改善。
首先，互相关方法的跟踪窗位置是固定的，而本发明所述方法中跟踪窗随特征点一起移动。因此，本发明所述方法具有很好的机动性，容易在局部取得更大的空间分辨率。理论上，使用本文提到的图像跟踪方法可以实现1/4波长的空间分辨率。而互相关方法的分辨率为3到5个波长。
其次，互相关方法涉及大量乘法，计算量很大。而本发明所述方法只有简单的逻辑判断及加减运算，所以计算量相对小很多，时间分辨率有很大提高。
最后，由于本发明所述方法不会导致误差积累，所以温度分辨率保持恒定，不会随时间的增加而降低。而互相关方法误差会累加，温度分辨率呈发散态势。
权利要求
1.超声测温二值图像模糊跟踪方法，首先利用超声测温仪中的超声探头发射超声脉冲，经人体生物组织反馈后接收返回的超声回波，其特征在于该方法还包括如下步骤1)超声测温仪的采样频率为超声探头发射的超声脉冲中心频率的N倍，其中N≥8，据此设置跟踪窗的尺寸dx×dy，dx方向为帧内方向，dx取N/2～N，dy方向为帧间方向，dy至少为4；2)利用超声测温仪中的微处理器对采样得到的前三帧超声回波数据进行二值化处理，波峰波谷特征点置1，其它点置0；将所述前三帧的波峰波谷特征点直接设定为已评估特征点，其整数累加时移值设为0，并在波峰波谷特征点附近取5个点作4次多项式曲线拟合，所述4次多项式被微分以确定其最大值从而得到偏移量；3)利用超声测温仪中的微处理器对采样得到的下一帧超声回波数据进行二值化处理，波峰波谷特征点置1，其它点置0；每个波峰波谷特征点设定相应的VMARK值，以当前帧回波的最高峰值VM为标准，按黄金分割法分成5级，当峰值Vk≥VM·(1-α)i(i＝1，...，5)时，VMARK＝i，其中α＝0.618；剩下的波峰波谷特征点则在当前帧最高峰值VM的0.8％之下，令VMARK＝0；4)对当前帧中VMARK≠0的特征点，用跟踪窗进行遍历跟踪，即以当前帧中VMARK≠0的特征点为当前跟踪窗的基点，按VMARK值由低到高处理，在dx×dy方向按权从小到大的顺序搜索已评估特征点，直至搜索到第一个已评估特征点，获取所述基点的整数累加时移值，若搜索不到已评估特征点，则将所述基点的VMARK置0，完成评估过程；其中，在当前跟踪窗内，前一帧中与基点的边相邻的采样点的权为1，离所述权为1的采样点越远，权值越大；每评估完一个基点，就在该基点附近取5个点作4次多项式曲线拟合，所述4次多项式被微分以确定其最大值从而得到偏移量，该偏移量减去搜索到的已评估特征点的偏移量，然后加到该基点的整数累加时移值上，获得总的累加时移值；5)对VMARK＝0的特征点，以步骤4)中已评估特征点为基准，通过线性插值得到这些特征点的累加时移值；6)使用中值滤波器对当前帧中的已评估特征点进行帧内过滤；7)对当前帧中的非特征点，通过帧内线性插值来得到其累加时移值；8)对当前帧所有点的累加时移值进行曲线拟合，差分后得到该帧所有点的局部时移值；9)返回步骤3)进入下一帧超声回波数据的处理，直至超声测温结束。
2.根据权利要求1所述的超声测温二值图像模糊跟踪方法，其特征在于步骤4)中，在当前跟踪窗内，前一帧中与基点的边相邻的采样点的权为1，与所述权为1的采样点全相邻的采样点权为2，以此类推，离所述权为1的采样点越远，权值越大。
3.根据权利要求1所述的超声测温二值图像模糊跟踪方法，其特征在于步骤4)中，在当前跟踪窗内，前一帧中与基点的边相邻的采样点的权为1，与所述权为1的采样点边相邻的采样点权为2，以此类推，离所述权为1的采样点越远，权值越大。
4.根据权利要求1所述的超声测温二值图像模糊跟踪方法，其特征在于步骤6)中所述中值滤波器为3到5点中值滤波器。
5.根据权利要求1所述的超声测温二值图像模糊跟踪方法，其特征在于步骤8)中所述曲线拟合采用6～10阶多项式。
全文摘要
超声测温二值图像模糊跟踪方法，属于医疗影像技术领域。为了解决超声无损测温中互相关方法精度低，误差会累积，且计算量大的问题，本发明公开了一种超声测温二值图像模糊跟踪方法，该方法利用超声测温仪中的超声探头发射超声脉冲，经人体生物组织反馈后接收返回的超声回波，对每帧超声回波的数据进行二值化处理，得到含有波峰波谷特征点的二值图像；然后，使用跟踪窗对二值图像进行遍历跟踪，得到超声回波累加时移值；最后，对累加时移值进行曲线拟合，差分后得到局部时移。本发明所述方法是一种实时方法，精度可以达到1/4波长，误差不会累积，并且计算量也大幅度下降，优于互相关方法。
文档编号G01S15/00GK1555760SQ200410000009
公开日2004年12月22日 申请日期2004年1月2日 优先权日2004年1月2日
发明者周兆英, 徐海东, 张毓笠 申请人:清华大学