CN105748100A - 准静态超声弹性成像位移计算方法和装置 - Google Patents

准静态超声弹性成像位移计算方法和装置 Download PDF

Info

Publication number
CN105748100A
CN105748100A CN201410800079.6A CN201410800079A CN105748100A CN 105748100 A CN105748100 A CN 105748100A CN 201410800079 A CN201410800079 A CN 201410800079A CN 105748100 A CN105748100 A CN 105748100A
Authority
CN
China
Prior art keywords
cross
correlation
displacement
time shift
tau
Prior art date
Legal status (The legal status is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the status listed.)
Granted
Application number
CN201410800079.6A
Other languages
English (en)
Other versions
CN105748100B (zh
Inventor
刘德清
冯乃章
董颖
Current Assignee (The listed assignees may be inaccurate. Google has not performed a legal analysis and makes no representation or warranty as to the accuracy of the list.)
Sonoscape Medical Corp
Original Assignee
Sonoscape Medical Corp
Priority date (The priority date is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the date listed.)
Filing date
Publication date
Application filed by Sonoscape Medical Corp filed Critical Sonoscape Medical Corp
Priority to CN201410800079.6A priority Critical patent/CN105748100B/zh
Publication of CN105748100A publication Critical patent/CN105748100A/zh
Application granted granted Critical
Publication of CN105748100B publication Critical patent/CN105748100B/zh
Active legal-status Critical Current
Anticipated expiration legal-status Critical

Links

Abstract

本发明公开了一种准静态超声弹性成像位移计算方法,以解决现有技术中存在的容易因迭代导致连锁错误,二维互相关搜索的计算量太大等技术问题。本发明实施例还提供相应的装置。本发明一些可行的实施方式中,方法可包括:对于任一个位移计算中心点,利用先验时移建立基于邻域窗的迭代互相关函数,进行迭代互相关计算;获取所述迭代互相关计算的互相关可靠性指数值,判断是否小于第一可信阈值;如果大于或等于,根据所述先验时移确定横向时移,根据所述先验时移和迭代互相关结果确定纵向时移;如果小于,使用邻域内的二维综合互相关搜索,得到所述位移计算中心点的纵向时移和横向时移;遍历全部超声射频信号,生成所述被测生物组织的位移场分布。

Description

准静态超声弹性成像位移计算方法和装置
技术领域
本发明涉及超声成像技术领域,具体涉及一种准静态超声弹性成像位移计算方法和装置。
背景技术
生物组织的软硬程度变化通常都和组织的病理变化直接相关。比如乳腺的恶性癌变,四分之三以上的都表现为硬块,而良性的增生性病变绝大部分都表现为和正常组织相似的软硬程度。超声(Ultrasound)弹性成像(“Elastography”or“StrainImaging”)是近十多年来出现的医学超声成像的新模式,主要是利用超声对组织力学特性参数进行检测并进行成像,为临床疾病诊断提供更加丰富的功能信息。超声弹性成像根据激励组织方式的不同,可以分为准静态超声弹性成像,瞬时超声弹性成像和基于声辐射力的剪切波弹性成像等。
准静态超声弹性成像的基本原理是,操作者手持超声探头,使超声探头表面和被测量组织表面轻轻的接触,然后按照一定的频率垂直按压组织,超声仪器会分别采集组织被压缩前和被压缩后的射频信号,计算出轴向方向上组织各点的位移分布,然后沿轴向方向对位移求梯度,得到组织的应变分布。弹性模量较大的区域,引起的应变比较小;反之弹性模量比较小的区域,相应的应变比较大。因此就得到了被检测组织的软硬程度的分布。
申请号为201210363616.6的中国专利申请公开了一种实时准静态超声弹性成像方法,该方法通过建立两帧超声射频信号的复数互相关函数,求取对应的位移场分布图,进而根据位移场分布图求取超声射频信号对应的生物组织应变分布图。
实践发现,上述方法中,每一次互相关计算的结果都作为下一次互相关计算的先验时移,采用依次迭代的方式进行位移的计算,所以,当探测计算过程中某一次信号片段互相关计算出现错误时,这种错误会在后续的迭代过程中引发连锁反应,导致后续的互相关计算全部错误,造成较大片的错误区域。
另外,由于生物组织中存在血管、包膜和纤维等复杂成分,导致其受到探头垂直向下的按压力时,组织发生的位移非常复杂,而远远不仅仅是纵向上发生位移,而是形成包括横向、纵向多个方向的位移场。所以只在纵向方向进行迭代跟踪的方法难以在人体组织的弹性测量中取得比较好的效果。
有一些算法采用二维互相关搜索的方法能够同时检测纵向方向上和横向方向上发生的位移,能够在一定程度上增加弹性成像位移跟踪计算的准确性。但是由于这种方法对于每一个点都需要做两维方向的搜索,然后找寻该点的最佳结果,所以导致其计算量相对于常规迭代跟踪算法的呈数十倍增长,会造成弹性成像实时性实现比较困难。
发明内容
本发明实施例提供一种准静态超声弹性成像位移计算方法和装置,以解决现有技术中存在的容易因迭代导致连锁错误,二维互相关搜索的计算量太大等技术问题。
本发明第一方面提供一种准静态超声弹性成像位移计算方法,包括:
采集得到连续的多帧超声射频信号;
对于两帧超声射频信号上的任一个位移计算中心点,利用先验时移建立基于邻域窗的迭代互相关函数,进行迭代互相关计算;
获取所述迭代互相关计算的互相关可靠性指数值,判断所述互相关可靠性指数值是否小于第一可信阈值;
如果所述互相关可靠性指数值大于或等于第一可信阈值,根据所述先验时移确定所述位移计算中心点的横向时移,根据所述先验时移和所述迭代互相关计算的迭代互相关结果确定所述位移计算中心点的纵向时移;
如果所述互相关可靠性指数值小于第一可信阈值,对所述位移计算中心点使用邻域内的二维综合互相关搜索,根据邻域内的二维综合互相关系数最大时的二维综合互相关结果,得到所述位移计算中心点的纵向时移和横向时移;
根据遍历全部超声射频信号得到的所有位移计算中心点的纵向时移和横向时移,生成所述被测生物组织的位移场分布。
结合第一方面,在第一种可能的实现方式中,所述对于两帧超声射频信号上的任一个位移计算中心点,利用先验时移建立基于邻域窗的迭代互相关函数包括:
记先验时移为[tp,tq],建立如下的迭代互相关函数:
R ( t w 0 , t h 0 ) = Σ t h = t h 0 - h / 2 t h 0 + h / 2 Σ t w = t w 0 - w / 2 t w 0 + w / 2 x ( t w , t h ) y ( t w - t p , t h - t q ) * = Ae jω τ w , τ h ,
其中,x()和y()分别表示所述两帧超声射频信号,tw和th是所述超声射频信号在邻域窗内的横向坐标和纵向坐标,[tw0,th0]是第一帧超声射频信号x()上的所述位移计算中心点的坐标,[tw0-tp,th0-tq]是第二帧超声射频信号y()上对应的位移计算中心点的坐标,[w,h]为邻域窗的宽度和高度,R(tw0,th0)表示[tw0,th0]的邻域窗内的数据与[tw0-tp,th0-tq]的邻域窗内的数据的互相关函数,“*”表示共轭运算,A表示迭代互相关结果的幅度,为迭代互相关结果的相位。
结合第一方面的第一种可能的实现方式,在第二种可能的实现方式中,所述获取所述迭代互相关计算的互相关可靠性指数值包括:
采用公式计算所述互相关可靠性指数值Qc,其中,
R sum = Σ t h = t h 0 - h / 2 t h 0 + h / 2 Σ t w = t w 0 - w / 2 t w 0 + w / 2 | | x ( t w , t h ) | + | y ( t w - t p , t h - t q ) | | ,
R sub = Σ t h = t h 0 - h / 2 t h 0 + h / 2 Σ t w = t w 0 - w / 2 t w 0 + w / 2 | | x ( t w , t h ) | - | y ( t w - t p , t h - t q ) | | .
结合第一方面的第一种可能的实现方式,在第三种可能的实现方式中,所述获取所述迭代互相关计算的互相关可靠性指数值包括:
采用以下公式计算所述互相关可靠性指数值Qc:
Qc = | R ( t w 0 , t h 0 ) | R ( t w , t h ) R ( t w - t p , t h - t q ) ,
其中,
R ( t w , t h ) = Σ t h = t h 0 - h / 2 t h 0 + h / 2 Σ t w = t w 0 - w / 2 t w 0 + w / 2 x ( t w , t h ) 2 ,
R ( t w - t p , t h - t q ) = Σ t h = t h 0 - h / 2 t h 0 + h / 2 Σ t w = t w 0 - w / 2 t w 0 + w / 2 y ( t w - t p , t h - t q ) 2 .
结合第一方面的第一种至第三种可能的实现方式中的任一种,在第四种可能的实现方式中,所述根据所述先验时移确定所述位移计算中心点的横向时移,根据所述先验时移和所述迭代互相关计算的迭代互相关结果确定所述位移计算中心点的纵向时移包括:
根据所述迭代互相关结果的相位计算则,所述位移计算中心点的时移为:Uh=Ps+th,Uw=tw
其中,Uh为纵向时移,Uw为横向时移,λ是频率为发射脉冲中心频率的超声波在所述被测生物组织中传播的波长。
结合第一方面或第一方面的第一种至第三种可能的实现方式中的任一种,在第五种可能的实现方式中,所述对所述位移计算中心点使用邻域内的二维综合互相关搜索,根据邻域内的二维综合互相关系数最大时的二维综合互相关结果,得到所述位移计算中心点的纵向时移和横向时移包括:
记二维综合互相关搜索的范围是[W,H],二维综合互相关搜索的公式为:
R ( t w 0 , t h 0 , τ w , τ h ) = Σ t h = t h 0 - h / 2 t h 0 + h / 2 Σ t w = t w 0 - w / 2 t w 0 + w / 2 x ( t w 0 , t h 0 ) y ( t w 0 - τ w , t h 0 - τ h ) * = Ae j ω τ w , τ h ,
其中,x()和y()分别表示所述的两帧超声射频信号,[tw0,th0]是第一帧超声射频信号x()上的所述位移计算中心点的坐标,[τwh]是第二帧超声射频信号y()的搜索范围内偏离所述位移计算中心点的横向偏移量和纵向偏移量,*表示共轭运算,A表示互相关结果的幅度,为二维综合互相关结果的相位,τw在[-W/2,W/2]范围内取值,τh在[-H/2,H/2]范围内取值;
二维综合互相关系数为:
Cr ( t w 0 , t h 0 , τ w , τ h ) = | R ( t w 0 , t h 0 , τ w , τ h ) | R ( t w 0 , t h 0 , τ w , τ h ) R ( t w 0 - τ w , t h 0 - τ h ) ,
其中, R ( t w 0 , t h 0 ) = Σ t h = t h 0 - h / 2 t h 0 + h / 2 Σ t w = t w 0 - w / 2 t w 0 + w / 2 x ( t w , t h ) 2 ,
R ( t w 0 - τ w , t h 0 - τ h ) = Σ t h = t h 0 - h / 2 t h 0 + h / 2 Σ t w = t w 0 - w / 2 t w 0 + w / 2 y ( t w - τ w , t h - τ h ) 2 ,
其中,“||”表示复数取模或实数取绝对值运算;
记二维综合互相关系数最大时,对应的偏移量是[τwMhM],对应的相角为则,
所述位移计算中心点的时移如下:
Ps = ω τ wM , τ hM 2 π λ , U h = Ps + τ hM , U w = τ wM ,
其中,Uh为二维综合互相关计算结果的纵向时移,Uw为二维综合互相关计算结果的横向时移,λ是频率为发射脉冲中心频率的超声波在所述被测生物组织中传播的波长。
结合第一方面,在第六种可能的实现方式中,每一帧所述多帧超声射频信号包括多条扫描线,所述探测步骤中,纵向方向上每一条扫描线的数据的采样率是所述超声探头发射脉冲的中心频率的至少2倍。
结合第一方面,在第七种可能的实现方式中,以同一根扫描线上、纵向方向上的前一个位移计算中心点的时移,作为所述先验时移。
结合第一方面,在第八种可能的实现方式中所述方法还包括:
判断所述位移场分布中每一个位移计算中心点的互相关可靠性指数是否小于第二可信阈值,如果不小于,判断所述位移计算中心点的位移为可信位移;
对所有可信位移区域进行平滑微分梯度计算,得到所述被测生物组织的二维应变分布结果。
本发明第二方面提供一种准静态超声弹性成像位移计算装置,包括:
采集模块,用于采集得到连续的多帧超声射频信号;
迭代互相关模块,用于对于两帧超声射频信号上的任一个位移计算中心点,利用先验时移建立基于邻域窗的迭代互相关函数,进行迭代互相关计算;
可靠性判断模块,用于获取所述迭代互相关计算的互相关可靠性指数值,判断所述互相关可靠性指数值是否小于第一可信阈值;
第一时移计算模块,用于如果所述互相关可靠性指数值大于或等于第一可信阈值,根据所述先验时移确定所述位移计算中心点的横向时移,根据所述先验时移和所述迭代互相关计算的迭代互相关结果确定所述位移计算中心点的纵向时移;
第二时移计算模块,用于如果所述互相关可靠性指数值小于第一可信阈值,对所述位移计算中心点使用邻域内的二维综合互相关搜索,根据邻域内的二维综合互相关系数最大时的二维综合互相关结果,得到所述位移计算中心点的纵向时移和横向时移;
位移场生成模块,用于根据遍历全部超声射频信号得到的所有位移计算中心点的纵向时移和横向时移,生成所述被测生物组织的位移场分布。
结合第二方面,在第一种可能的实现方式中,所述迭代互相关模块,具体用于:
记先验时移为[tp,tq],建立如下的迭代互相关函数:
R ( t w 0 , t h 0 ) = Σ t h = t h 0 - h / 2 t h 0 + h / 2 Σ t w = t w 0 - w / 2 t w 0 + w / 2 x ( t w , t h ) y ( t w - t p , t h - t q ) * = Ae jω τ w , τ h ,
其中,x()和y()分别表示所述两帧超声射频信号,tw和th是所述超声射频信号在邻域窗内的横向坐标和纵向坐标,[tw0,th0]是第一帧超声射频信号x()上的所述位移计算中心点的坐标,[tw0-tp,th0-tq]是第二帧超声射频信号y()上对应的位移计算中心点的坐标,[w,h]为邻域窗的宽度和高度,R(tw0,th0)表示[tw0,th0]的邻域窗内的数据与[tw0-tp,th0-tq]的邻域窗内的数据的互相关函数,“*”表示共轭运算,A表示迭代互相关结果的幅度,为迭代互相关结果的相位。
结合第二方面的第一种可能的实现方式,在第二种可能的实现方式中,所述可靠性判断模块具体用于:
采用公式计算所述互相关可靠性指数值Qc,其中,
R sum = Σ t h = t h 0 - h / 2 t h 0 + h / 2 Σ t w = t w 0 - w / 2 t w 0 + w / 2 | | x ( t w , t h ) | + | y ( t w - t p , t h - t q ) | | ,
R sub = Σ t h = t h 0 - h / 2 t h 0 + h / 2 Σ t w = t w 0 - w / 2 t w 0 + w / 2 | | x ( t w , t h ) | - | y ( t w - t p , t h - t q ) | | .
结合第二方面的第一种可能的实现方式,在第三种可能的实现方式中,所述可靠性判断模块具体用于:
采用下述公式计算所述互相关可靠性指数值Qc:
Qc = | R ( t w 0 , t h 0 ) | R ( t w , t h ) R ( t w - t p , t h - t q ) ,
其中,
R ( t w , t h ) = Σ t h = t h 0 - h / 2 t h 0 + h / 2 Σ t w = t w 0 - w / 2 t w 0 + w / 2 x ( t w , t h ) 2 ,
R ( t w - t p , t h - t q ) = Σ t h = t h 0 - h / 2 t h 0 + h / 2 Σ t w = t w 0 - w / 2 t w 0 + w / 2 y ( t w - t p , t h - t q ) 2 .
结合第二方面的第一种至第三种可能的实现方式中的任一种,在第四种可能的实现方式中,所述第一时移计算模块具体用于:
根据所述迭代互相关结果的相位计算则,所述位移计算中心点的时移为:Uh=Ps+th,Uw=tw
其中,Uh为纵向时移,Uw为横向时移,λ是频率为发射脉冲中心频率的超声波在所述被测生物组织中传播的波长。
结合第二方面或第二方面的第一种至第三种可能的实现方式中的任一种,在第五种可能的实现方式中,所述第二时移计算模块具体用于:
记二维综合互相关搜索的范围是[W,H],二维综合互相关搜索的公式为:
R ( t w 0 , t h 0 , τ w , τ h ) = Σ t h = t h 0 - h / 2 t h 0 + h / 2 Σ t w = t w 0 - w / 2 t w 0 + w / 2 x ( t w 0 , t h 0 ) y ( t w 0 - τ w , t h 0 - τ h ) * = Ae j ω τ w , τ h ,
其中,x()和y()分别表示所述的两帧超声射频信号,[tw0,th0]是第一帧超声射频信号x()上的所述位移计算中心点的坐标,[τwh]是第二帧超声射频信号y()的搜索范围内偏离所述位移计算中心点的横向偏移量和纵向偏移量,*表示共轭运算,A表示互相关结果的幅度,为二维综合互相关结果的相位,τw在[-W/2,W/2]范围内取值,τh在[-H/2,H/2]范围内取值;
二维综合互相关系数为:
Cr ( t w 0 , t h 0 , τ w , τ h ) = | R ( t w 0 , t h 0 , τ w , τ h ) | R ( t w 0 , t h 0 , τ w , τ h ) R ( t w 0 - τ w , t h 0 - τ h ) ,
其中, R ( t w 0 , t h 0 ) = Σ t h = t h 0 - h / 2 t h 0 + h / 2 Σ t w = t w 0 - w / 2 t w 0 + w / 2 x ( t w , t h ) 2 ,
R ( t w 0 - τ w , t h 0 - τ h ) = Σ t h = t h 0 - h / 2 t h 0 + h / 2 Σ t w = t w 0 - w / 2 t w 0 + w / 2 y ( t w - τ w , t h - τ h ) 2 ,
其中,“||”表示复数取模或实数取绝对值运算;
记二维综合互相关系数最大时,对应的偏移量是[τwMhM],对应的相角为则,
所述位移计算中心点的时移如下:
Ps = ω τ wM , τ hM 2 π λ , U h = Ps + τ hM , U w = τ wM ,
其中,Uh为二维综合互相关计算结果的纵向时移,Uw为二维综合互相关计算结果的横向时移,λ是频率为发射脉冲中心频率的超声波在所述被测生物组织中传播的波长。
结合第二方面,在第六种可能的实现方式中,所述可靠性判断模块,还用于判断所述位移场分布中每一个位移计算中心点的互相关可靠性指数是否小于第二可信阈值,如果不小于,判断所述位移计算中心点的位移为可信位移;
所述装置还包括:组织应变生成模块,用于对所有可信位移区域进行平滑微分梯度计算,得到所述被测生物组织的二维应变分布结果。
由上可见,本发明实施例通过采用上述技术方案,取得了以下技术效果:
1、本发明技术方案中,可通过互相关可靠性指数值来判断迭代互相关计算的结果是否可信,根据迭代互相关结果是否可信,自适应的决定:是通过迭代互相关结果和先验位移来计算位移计算中心点的时移,还是进一步通过二维综合互相关搜索计算位移计算中心点的时移。
2、在判断迭代互相关结果可信,利用迭代互相关结果计算纵向时移时,由于迭代互相关函数是基于邻域窗建立的,因而能够快速检测纵向方向上发生的时移,同时,能够降低横向时移对纵向时移计算的影响,从而可以提高纵向时移的准确度。
3、在判断迭代互相关结果不可信时,能够自适应的使用二维综合互相关搜索,能够准确检测纵向方向上和横向方向上发生的位移。
4、在判断迭代互相关结果为不可信时,使用邻域内的二维综合互相关搜索,而不是根据迭代互相关结果,来得到位移计算中心点的时移,可以将不可靠的计算结果剔除,可避免后续连锁错误;
5、在判断迭代互相关结果为可信时,可以直接根据先验时移确定位移计算中心点的横向时移,可根据先验时移和迭代互相关结果计算得到位移计算中心点的纵向时移,不必进行二维综合互相关搜索,因而大大减少了计算量,有利于实现实时成像。
附图说明
为了更清楚地说明本发明实施例技术方案,下面将对实施例和现有技术描述中所需要使用的附图作简单地介绍,显而易见地,下面描述中的附图仅仅是本发明的一些实施例,对于本领域普通技术人员来讲,在不付出创造性劳动的前提下,还可以根据这些附图获得其它的附图。
图1是本发明实施例提供的一种准静态超声弹性成像位移计算方法的流程示意图;
图2a是多帧超声射频信号组织成的矩阵结构的示意图;
图2b是矩阵结构上邻域窗的示意图;
图2c是矩阵结构上邻域窗与二维综合互相关搜索范围的示意图;
图3是本发明实施例提供的另一种准静态超声弹性成像位移计算方法的流程示意图;
图4是本发明实施例提供的一种准静态超声弹性成像位移计算装置的结构示意图。
具体实施方式
本发明实施例提供一种准静态超声弹性成像位移计算方法,以解决现有技术中存在的容易因迭代导致连锁错误,二维互相关搜索的计算量太大等技术问题。本发明实施例还提供相应的装置。
为了使本技术领域的人员更好地理解本发明方案,下面将结合本发明实施例中的附图,对本发明实施例中的技术方案进行清楚、完整地描述,显然,所描述的实施例仅仅是本发明一部分的实施例,而不是全部的实施例。基于本发明中的实施例,本领域普通技术人员在没有做出创造性劳动前提下所获得的所有其他实施例,都应当属于本发明保护的范围。
下面通过具体实施例,分别进行详细的说明。
实施例一、
请参考图1,本发明实施例提供一种准静态超声弹性成像位移计算方法,可包括:
110、采集得到连续的多帧超声射频信号。
在使用超声测量仪的超声探头对被测生物组织例如人体皮肤进行按压探测时,超声探头会周期性发送超声波,并接收被测生物组织的超声回波,从而采集得到多帧超声射频信号。采集得到的每一帧超声射频信号可包括多条扫描线,所述扫描线可以是垂直于被测生物组织的表面的探测线,每条扫描线可包括多个采样点。优选的,纵向方向上每一条扫描线的数据的采样率可以是超声探头发射脉冲的中心频率的至少2倍。
本步骤中,采集得到的多帧超声射频信号可被组织成一个矩阵结构,如图2a所示,图中每个小圆圈表示一个采样点,每一列的多个采样点可表示一条扫描线。探测过程中,使用者采用超声探头按压被测生物组织,采集任意相邻的两帧超声射频信号的时刻,被测生物组织受力不同,产生不同的形变,采集得到的超声射频信号可反映该形变。
120、对于两帧超声射频信号上的任一个位移计算中心点,利用先验时移建立基于邻域窗的迭代互相关函数,进行迭代互相关计算。
本发明实施例中,依次将每一帧超声射频信号上的每一个采样点,或者按照预设规则选中的每一个采样点,作为位移计算中心点,针对每个位移计算中心点,计算该点在两次采集超声射频信号时的位移。由于超声信号在被测生物组织中的传播速度不变,因此,位移与超声信号的传送时间延迟成线性关系,可以用超声信号的传送时间延迟来表示位移,可将超声信号的传送时间延迟称为时移。
本发明实施例中,对于两帧超声射频信号,可以将前一个位移计算中心点的时移,作为当前的位移计算中心点的先验时移。优选的,所说的前一个位移计算中心点,是指同一根扫描线上、纵向方向上的前一个位移计算中心点。即,可以以同一根扫描线上、纵向方向上的前一个位移计算中心点的时移,作为所述先验时移。
所说的先验时移包括纵向(扫描线方向)时移和横向(垂直于扫描线方向)时移,例如可以用[tp,tq]表示。利用该先验时移,可以建立基于邻域窗的迭代互相关函数。所说的邻域窗是一个如2a所示的矩阵结构上的二维窗口,纵向(扫描线方向)和横向(垂直于扫描线方向)方向上均包括多个采样点。
所说的迭代互相关函数,可以是第一邻域窗内所有采样点的数据与第二邻域窗内所有采样点的数据的复数互相关函数,第一邻域窗是指以第一帧超声射频信号上当前的位移计算中心点为中心的邻域窗,第二邻域窗是指以第二帧超声射频信号上以第一邻域窗的中心坐标减去先验时移后的坐标为中心的邻域窗,第一帧超声射频信号和第二帧超声射频信号是组织形变不同的两帧超声射频信号。
利用迭代互相关函数进行互相关计算,可得到一个复数的互相关结果。本发明实施例中,该互相关结果可用于作为计算当前的位移计算中心点的时移的依据之一。
130、获取所述迭代互相关计算的互相关可靠性指数值,判断所述互相关可靠性指数值是否小于第一可信阈值。
本发明实施例中,可计算获取上一步骤迭代互相关计算的互相关可靠性指数值,利用该互相关可靠性指数值来判断先验时移(即前一位移计算中心点的时移)是否可信,根据可信与否,来决定后续采用不同的方式计算获取当前的位移计算中心点的时移。判断方法可以是,判断所述互相关可靠性指数值是否小于第一可信阈值,如果不小于(即大于或等于),则判断为先验时移可信,如果小于,则判断为先验时移不可信。
140、如果所述互相关可靠性指数值大于或等于第一可信阈值,根据所述先验时移确定所述位移计算中心点的横向时移,根据所述先验时移和所述迭代互相关计算的迭代互相关结果确定所述位移计算中心点的纵向时移。
本发明实施例中,如果互相关可靠性指数值大于或等于第一可信阈值,即先验时移可信,则可以直接利用先验时移和所述迭代互相关计算的互相关结果来确定当前的位移计算中心点的时移,例如,将先验时移和互相关结果的和,作为当前的位移计算中心点的时移。由于直接利用先验时移确定当前的时移,本步骤计算复杂度被大大简化,计算量大大减少,可提高实时性。
150、如果所述互相关可靠性指数值小于第一可信阈值,对所述位移计算中心点使用邻域内的二维综合互相关搜索,根据邻域内的二维综合互相关系数最大时的二维综合互相关结果,得到所述位移计算中心点的纵向时移和横向时移。
本发明实施例中,如果互相关可靠性指数值小于第一可信阈值,即先验时移不可信,则,为了避免当前的位移计算错误,以及为了避免后续的连锁错误,此时,不再根据先验实时移来计算当前的位移计算中心点的时移,而是,对所述位移计算中心点使用邻域内的二维综合互相关搜索,根据邻域内的二维综合互相关系数最大时的二维综合互相关结果,得到所述位移计算中心点的时移。这种通过二维综合互相关搜索,计算得到时移的方式,可以准确的得到当前的时移,不被先验时移错误影响,并且,可同时得到纵向和横向上的时移。
160、根据遍历全部超声射频信号得到的所有位移计算中心点的纵向时移和横向时移,生成所述被测生物组织的位移场分布。
本发明实施例中,通过遍历采集得到的全部超声射频信号,对每一个采样点,或者,按照预设规则选中的每一个采样点,按照上述步骤120-150进行计算处理,可得到各个采样点的时移。然后,根据得到的所有时移,将时移乘以超声波的速度转换为位移,可生成被测生物组织的位移场分布。
然后,对位移场分布做纵向和横向上的微分运算,即可得到被测生物组织的二维应变分布结果。本发明实施例中,为了提高应变分布结果的准确性,可选的,上述方法还可以包括以下步骤:
判断所述位移场分布中每一个位移计算中心点的互相关可靠性指数是否小于第二可信阈值,如果不小于,判断所述位移计算中心点的位移为可信位移;对所有可信位移区域进行平滑微分梯度计算,得到所述被测生物组织的二维应变分布结果。
本步骤中,通过将每一个位移计算中心点的互相关可靠性指数与第二可信阈值进行比较,将不可信的位移排除,后续,仅对可信位移区域进行平滑微分梯度计算,可得到准确度更高、更可信的被测生物组织的二维应变分布结果。
由上可见,本发明实施例提供了一种准静态超声弹性成像位移计算方法,该方法通过采用上述技术方案,取得了以下技术效果:
1、本发明技术方案中,可通过互相关可靠性指数值来判断迭代互相关计算的结果是否可信,根据迭代互相关结果是否可信,自适应的决定:是通过迭代互相关结果和先验位移来计算位移计算中心点的时移,还是进一步通过二维综合互相关搜索计算位移计算中心点的时移。
2、在判断迭代互相关结果可信,利用迭代互相关结果计算纵向位移时,由于迭代互相关函数是基于邻域窗建立的,因而能够快速检测纵向方向上发生的位移,同时,能够降低横向位移对纵向位移计算的影响;从而可以提高纵向位移的准确度。
3、在判断迭代互相关结果不可信时,能够自适应的使用二维综合互相关搜索,能够准确检测纵向方向上和横向方向上发生的位移。
4、在判断迭代互相关结果为不可信时,使用邻域内的二维综合互相关搜索,而不是根据迭代互相关结果,来得到位移计算中心点的时移,可以将不可靠的计算结果剔除,可避免后续连锁错误。
5、在判断迭代互相关结果为可信时,可以直接根据先验时移和迭代互相关结果计算得到位移计算中心点的时移,不必进行二维综合互相关搜索,因而大大减少了计算量,有利于实现实时成像。
为便于更好的理解本发明实施例提供的技术方案,下面通过一个具体场景下的实施方式为例进行介绍。
请参考图3,本发明实施例的准静态超声弹性成像位移计算方法,可包括:
301、获取压缩前后两帧超声射频信号。
本步骤中,获取连续的多帧超声射频信号中的、压缩前后(或者说,组织形变不同)的两帧超声射频信号。其中,两帧超声射频信号可以分别用x()和y()表示。x和y都是时间或者说时移t的函数。时移t可包括纵向时移th和横向时移tw。x()和y()可以连续的两帧超声射频信号,也可以是不连续的两帧超声射频信号。
302、针对当前的位移计算中心点,建立邻域窗。
例如,可假定当前的位移计算中心点的坐标是[tw0,th0],假设邻域窗的尺寸为[w,h],其中w和h分别为邻域窗的宽度和高度,邻域窗以[tw0,th0]为中心点。如图2b所示,是一种邻域窗的示意图,其中,实心圆圈表示当前的位移计算中心点,小方框表示邻域窗。本文中,为便于识别,可称当前的位移计算中心点的邻域窗为第一邻域窗。
303、利用先验时移建立基于邻域窗的迭代互相关函数,进行迭代互相关计算。
假设前一个位移计算中心点的时移,即先验时移是[tp,tq],则第二帧超声射频信号上对应的位移计算中心点的坐标是[tw0-tp,th0-tq],同时,可建立以[tw0-tp,th0-tq]为中心,尺寸为[w,h]的邻域窗。本文中,为便于识别,可称第二帧超声射频信号上对应的位移计算中心点的邻域窗为第二邻域窗。利用先验时移建立基于邻域窗的迭代互相关函数,实际上是[tw0,th0]的邻域窗即第一邻域窗内的数据与[tw0-tp,th0-tq]的邻域窗即第二邻域窗内的数据的互相关函数,可用R(tw0,th0)表示。
本发明一些实施例中,迭代互相关函数具体可以是如下公式(1):
R ( t w 0 , t h 0 ) = Σ t h = t h 0 - h / 2 t h 0 + h / 2 Σ t w = t w 0 - w / 2 t w 0 + w / 2 x ( t w , t h ) y ( t w - t p , t h - t q ) * = Ae jω τ w , τ h . . . ( 1 )
其中,x()和y()分别表示两帧超声射频信号,tw和th分别是所述超声射频信号在邻域窗内的横向坐标和纵向坐标,[tw0,th0]是第一帧超声射频信号x()上的所述位移计算中心点的坐标,[tw0-tp,th0-tq]是第二帧超声射频信号y()上对应的位移计算中心点的坐标,[w,h]为邻域窗的宽度和高度,R(tw0,th0)表示[tw0,th0]的邻域窗内的数据与[tw0-tp,th0-tq]的邻域窗内的数据的互相关函数,“*”表示共轭运算,A表示迭代互相关结果的幅度,为迭代互相关结果的相位。等效于两个邻域窗内的平均相位差。
304、计算迭代互相关计算的互相关可靠性指数值,判断互相关可靠性指数值是否小于第一可信阈值。
互相关可靠性指数值是反映迭代互相关计算结果的可靠性的一个数值。
本发明一些实施例中,可采用以下公式计算互相关可靠性指数值:
R sum = Σ t h = t h 0 - h / 2 t h 0 + h / 2 Σ t w = t w 0 - w / 2 t w 0 + w / 2 | | x ( t w , t h ) | + | y ( t w - t p , t h - t q ) | | . . . ( 2 )
R sub = Σ t h = t h 0 - h / 2 t h 0 + h / 2 Σ t w = t w 0 - w / 2 t w 0 + w / 2 | | x ( t w , t h ) | - | y ( t w - t p , t h - t q ) | | . . . ( 3 )
Qc = R sub R sum . . . ( 4 )
本发明另一些实施例中,也可以采用以下公式计算互相关可靠性指数:
R ( t w , t h ) = Σ t h = t h 0 - h / 2 t h 0 + h / 2 Σ t w = t w 0 - w / 2 t w 0 + w / 2 x ( t w , t h ) 2 . . . ( 5 )
R ( t w - t p , t h - t q ) = Σ t h = t h 0 - h / 2 t h 0 + h / 2 Σ t w = t w 0 - w / 2 t w 0 + w / 2 y ( t w - t p , t h - t q ) 2 . . . ( 6 )
Qc = | R ( t w 0 , t h 0 ) | R ( t w , t h ) R ( t w - t p , t h - t q ) . . . ( 7 )
本发明其它实施例中,还可以采用其它公式计算互相关可靠性指数,本文对此不作限制。
其中,Qc是所述互相关可靠性指数值,“||”表示复数取模或者实数取绝对值运算。当互相关可靠性指数值大于或等于第一可信阈值时,表示可信;小于第一可信阈值时,表示不可信。
所说的第一可信阈值,可以根据超声平台(例如超声测量仪)的计算处理能力,以及对于非可靠性弹性结果的容忍程度,优化的选择。例如,可以根据经验值确定。
305、进行二维综合互相关搜索,计算位移计算中心点的时移。
如果所述互相关可靠性指数值小于第一可信阈值,表示先验时移不可信。则,可以对所述位移计算中心点使用邻域内的二维综合互相关搜索,根据邻域内的二维综合互相关系数最大时的二维综合互相关结果,和所述二维综合互相关搜索的偏移量,得到所述位移计算中心点的时移。
本实施例中,可设立二维综合互相关搜索的范围[W,H],W为互相关搜索范围宽度,H为互相关搜索范围高度,W和H要分别大于迭代互相关邻域窗的宽度w和高度h。如图2c所示,图中较小的方框表示迭代互相关的邻域窗,较大的虚线方框表示二维综合互相关搜索的范围。
本发明一些实施例中,二维综合互相关搜索可以采用如下公式(8):
R ( t w 0 , t h 0 , τ w , τ h ) = Σ t h = t h 0 - h / 2 t h 0 + h / 2 Σ t w = t w 0 - w / 2 t w 0 + w / 2 x ( t w 0 , t h 0 ) y ( t w 0 - τ w , t h 0 - τ h ) * = Ae j ω τ w , τ h . . . ( 8 )
其中,x()和y()分别表示两帧超声射频信号,[tw0,th0]是第一帧超声射频信号x()上的所述位移计算中心点的坐标,[τwh]是第二帧超声射频信号y()的搜索范围内偏离所述位移计算中心点[tw0,th0]的横向偏移量和纵向偏移量,*表示共轭运算,A表示互相关结果的幅度,为互相关结果的相位。等效于两个邻域窗内的数据的平均相位差。
本发明一些实施例中,二维综合互相关系数可用如下公式表示:
R ( t w 0 , t h 0 ) = Σ t h = t h 0 - h / 2 t h 0 + h / 2 Σ t w = t w 0 - w / 2 t w 0 + w / 2 x ( t w , t h ) 2 . . . ( 9 )
R ( t w 0 - τ w , t h 0 - τ h ) = Σ t h = t h 0 - h / 2 t h 0 + h / 2 Σ t w = t w 0 - w / 2 t w 0 + w / 2 y ( t w - τ w , t h - τ h ) 2 . . . ( 10 )
Cr ( t w 0 , t h 0 , τ w , τ h ) = | R ( t w 0 , t h 0 , τ w , τ h ) | R ( t w 0 , t h 0 , τ w , τ h ) R ( t w 0 - τ w , t h 0 - τ h ) . . . ( 11 )
其中,Cr(tw0,th0wh)即是位移计算中心点[tw0,th0]在第二帧超声射频信号上偏移量为[τwh]的偏移互相关计算的二维综合互相关系数,“||”表示复数取模或实数取绝对值运算,τw在[-W/2,W/2]范围内取值,τh在[-H/2,H/2]范围内取值。
同时,本步骤中,可按照公式(2)-(4)或者(5)-(7)计算每一个位移计算中心点的互相关可靠性指数值。
按照上述计算方法遍历当前的位移计算中心点的整个二维综合互相关搜索的范围,可得到多个例如W×H个复数二维综合互相关结果、二维综合互相关系数和互相关可靠性指数。本步骤中,可在这W×H个二维综合互相关系数中找寻最大值,和最大值对应的偏移量[τwMhM],以及对应的二维综合互相关互相关函数R(tw0,th0wMhM)的相角
本步骤中,可根据得到的最大二维综合互相关系数值对应的相角和偏移量,计算两帧超声射频信号在当前的位移计算中心点的邻域窗内的平均时移,公式如下:
Ps = ω τ wM , τ hM 2 π λ , . . . ( 12 )
Uh=Ps+τhM,…………………(13)
Uw=τwM,……………………(14)
上述的Uh和Uw分别是所述位移计算中心点的纵向和横向时移。
其中,Uh为二维综合互相关计算结果的纵向时移,Uw为二维综合互相关计算结果的横向时移,即,Uh和Uw分别是位移计算中心点的纵向和横向时移;λ是频率为发射脉冲中心频率的超声波在所述被测生物组织中传播的波长。
306、根据先验时移和迭代互相关计算的互相关结果,确定位移计算中心点的时移。
如果所述互相关可靠性指数值大于或等于第一可信阈值,表示先验时移可信。则,可以根据先验时移和迭代互相关计算的互相关结果,确定位移计算中心点的时移,即,将先验时移和迭代互相关计算的互相关结果的和,作为当前的位移计算中心点的时移,公式如下:
Ps = ω τ w τ h 2 π λ . . . ( 15 )
Uh=Ps+th…………………(16)
Uw=tw……………………………(17)
其中,是迭代互相关计算的迭代互相关结果的相位,λ是频率为发射脉冲中心频率的超声波在所述被测生物组织中传播的波长;Uh为计算得到的纵向时移,Uw为计算得到的横向时移。
综上,通过步骤305和306,可得到当前的位移计算中心点的纵向时移和横向时移。
307、遍历采集得到的全部超声射频信号。
采用上述方法,遍历采集得到的全部超声射频信号,每一次计算之后,判断是否已经遍历完全部超声射频信号,如果还没有遍历完,则,将本次计算得到的时移作为下一次计算的先验时移。遍历完成后,得到所有计算的采样点的时移,进而可生成被测生物组织的位移场分布。
308、利用互相关可靠性指数对位移场分布进行分类。
本步骤中,可设一个第二可信阈值,判断所述位移场分布中每一个位移计算中心点的互相关可靠性指数是否小于第二可信阈值,如果不小于,判断所述位移计算中心点的位移为可信位移。第二可信阈值可以和第一可信阈值相同,也可以不相同,具体可以根据对非可靠弹性结果的容忍程度优化的选择。例如,可以根据经验值确定。
309、对所有可信位移区域进行微分求应变计算。
例如,可以对所有可信位移区域进行平滑微分梯度计算,后者其它算法计算,得到被测生物组织的二维应变分布结果。所有微分求应变计算在本质上可以理解于邻域窗内在受力方向上做平均差分运算。
由上可见,本发明实施例提供了一种准静态超声弹性成像位移计算方法。
该方法通过采用可信阈值判断的方法,自适应决策是进行迭代互相关计算位移,还是二维综合互相关搜索计算位移,能够做到按需要决策,所以能够做到既不增加很多运算量,又能够实现二维的位移跟踪计算。
该方法使用互相关可靠性指数值来判断位移跟踪结果的可靠性,将不可信的计算结果剔除,可防止异常的计算结果影响正常的弹性色彩映射显示。互相关可靠性指数的实现计算方法可以有很多种,可以是前文所述的计算方法,也可以是互相关系数,还可以是适用于一切反映两邻域数据窗相关性的计算以及运算结果。
综上,本发明实施例提供的准静态超声弹性成像位移计算方法,通过采用上述技术方案,取得了以下技术效果:
1、本发明技术方案中,可通过互相关可靠性指数值来判断迭代互相关计算的结果是否可信,根据迭代互相关结果是否可信,自适应的决定:是通过迭代互相关结果和先验位移来计算位移计算中心点的时移,还是进一步通过二维综合互相关搜索计算位移计算中心点的时移。
2、在判断迭代互相关结果可信,利用迭代互相关结果计算纵向时移时,由于迭代互相关函数是基于邻域窗建立的,因而能够快速检测纵向方向上发生的时移,同时,能够降低横向时移对纵向时移计算的影响,从而可以提高纵向时移的准确度。
3、在判断迭代互相关结果不可信时,能够自适应的使用二维综合互相关搜索,能够准确检测纵向方向上和横向方向上发生的位移。
4、在判断迭代互相关结果为不可信时,使用邻域内的二维综合互相关搜索,而不是根据迭代互相关结果,来得到位移计算中心点的时移,可以将不可靠的计算结果剔除,可避免后续连锁错误;
5、在判断迭代互相关结果为可信时,可以直接根据先验时移确定位移计算中心点的横向时移,可根据先验时移和迭代互相关结果计算得到位移计算中心点的纵向时移,不必进行二维综合互相关搜索,因而大大减少了计算量,有利于实现实时成像。
请参考图4,本发明实施例提供准静态超声弹性成像位移计算装置,该装置可包括:
采集模块410,用于采集得到连续的多帧超声射频信号;
迭代互相关模块420,用于对于两帧超声射频信号上的任一个位移计算中心点,利用先验时移建立基于邻域窗的迭代互相关函数,进行迭代互相关计算;
可靠性判断模块430,用于获取所述迭代互相关计算的互相关可靠性指数值,判断所述互相关可靠性指数值是否小于第一可信阈值;
第一时移计算模块440,用于如果所述互相关可靠性指数值大于或等于第一可信阈值,根据所述先验时移确定所述位移计算中心点的横向时移,根据所述先验时移和所述迭代互相关计算的迭代互相关结果确定所述位移计算中心点的纵向时移;
第二时移计算模块450,用于如果所述互相关可靠性指数值小于第一可信阈值,对所述位移计算中心点使用邻域内的二维综合互相关搜索,根据邻域内的二维综合互相关系数最大时的二维综合互相关结果,得到所述位移计算中心点的纵向时移和横向时移;
位移场生成模块460,用于根据遍历全部超声射频信号得到的所有位移计算中心点的纵向时移和横向时移,生成所述被测生物组织的位移场分布。
本发明一些实施例中,所述迭代互相关模块420,具体可用于:
记先验时移为[tp,tq],建立如下的迭代互相关函数:
R ( t w 0 , t h 0 ) = Σ t h = t h 0 - h / 2 t h 0 + h / 2 Σ t w = t w 0 - w / 2 t w 0 + w / 2 x ( t w , t h ) y ( t w - t p , t h - t q ) * = Ae jω τ w , τ h
其中,x()和y()分别表示所述两帧超声射频信号,tw和th分别是所述超声射频信号在邻域窗内的横向坐标和纵向坐标,[tw0,th0]是第一帧超声射频信号x()上的所述位移计算中心点的坐标,[tw0-tp,th0-tq]是第二帧超声射频信号y()上对应的位移计算中心点的坐标,[w,h]为邻域窗的宽度和高度,R(tw0,th0)表示[tw0,th0]的邻域窗内的数据与[tw0-tp,th0-tq]的邻域窗内的数据的互相关函数,“*”表示共轭运算,A表示迭代互相关结果的幅度,为迭代互相关结果的相位。
本发明一些实施例中,所述可靠性判断模块430具体可用于:
采用公式计算所述互相关可靠性指数值Qc,其中,
R sum = Σ t h = t h 0 - h / 2 t h 0 + h / 2 Σ t w = t w 0 - w / 2 t w 0 + w / 2 | | x ( t w , t h ) | + | y ( t w - t p , t h - t q ) | | ,
R sub = Σ t h = t h 0 - h / 2 t h 0 + h / 2 Σ t w = t w 0 - w / 2 t w 0 + w / 2 | | x ( t w , t h ) | - | y ( t w - t p , t h - t q ) | | .
本发明一些实施例中,所述可靠性判断模块430具体可用于:
采用以下公式计算所述互相关可靠性指数值Qc:
Qc = | R ( t w 0 , t h 0 ) | R ( t w , t h ) R ( t w - t p , t h - t q )
其中,
R ( t w , t h ) = Σ t h = t h 0 - h / 2 t h 0 + h / 2 Σ t w = t w 0 - w / 2 t w 0 + w / 2 x ( t w , t h ) 2 ,
R ( t w - t p , t h - t q ) = Σ t h = t h 0 - h / 2 t h 0 + h / 2 Σ t w = t w 0 - w / 2 t w 0 + w / 2 y ( t w - t p , t h - t q ) 2
本发明一些实施例中,所述第一时移计算模块240具体可用于:
根据所述迭代互相关结果的相位计算则,所述位移计算中心点的时移为:Uh=Ps+th,Uw=tw
其中,Uh为纵向时移,Uw为横向时移,λ是频率为发射脉冲中心频率的超声波在所述被测生物组织中传播的波长。
本发明一些实施例中,所述第二时移计算模块450具体可用于:
记二维综合互相关搜索的范围是[W,H],二维综合互相关搜索的公式为:
R ( t w 0 , t h 0 , τ w , τ h ) = Σ t h = t h 0 - h / 2 t h 0 + h / 2 Σ t w = t w 0 - w / 2 t w 0 + w / 2 x ( t w 0 , t h 0 ) y ( t w 0 - τ w , t h 0 - τ h ) * = Ae j ω τ w , τ h
其中,x()和y()分别表示所述的两帧超声射频信号,[tw0,th0]是第一帧超声射频信号x()上的所述位移计算中心点的坐标,[τwh]是第二帧超声射频信号y()的搜索范围内偏离所述位移计算中心点的横向偏移量和纵向偏移量,*表示共轭运算,A表示互相关结果的幅度,为二维综合互相关结果的相位,τw在[-W/2,W/2]范围内取值,τh在[-H/2,H/2]范围内取值;
二维综合互相关系数为:
Cr ( t w 0 , t h 0 , τ w , τ h ) = | R ( t w 0 , t h 0 , τ w , τ h ) | R ( t w 0 , t h 0 , τ w , τ h ) R ( t w 0 - τ w , t h 0 - τ h ) ,
其中, R ( t w 0 , t h 0 ) = Σ t h = t h 0 - h / 2 t h 0 + h / 2 Σ t w = t w 0 - w / 2 t w 0 + w / 2 x ( t w , t h ) 2 ,
R ( t w 0 - τ w , t h 0 - τ h ) = Σ t h = t h 0 - h / 2 t h 0 + h / 2 Σ t w = t w 0 - w / 2 t w 0 + w / 2 y ( t w - τ w , t h - τ h ) 2 ,
其中,“||”表示复数取模或实数取绝对值运算;
记二维综合互相关系数最大时,对应的偏移量是[τwMhM],对应的相角为则,
所述位移计算中心点的时移如下:
Ps = ω τ wM , τ hM 2 π λ , U h = Ps + τ hM , U w = τ wM ,
其中,Uh为二维综合互相关计算结果的纵向时移,Uw为二维综合互相关计算结果的横向时移,λ是频率为发射脉冲中心频率的超声波在所述被测生物组织中传播的波长。
本发明一些实施例中,所述装置还包括:组织应变生成模块;
所述可靠性判断模块430,还用于判断所述位移场分布中每一个位移计算中心点的互相关可靠性指数是否小于第二可信阈值,如果不小于,判断所述位移计算中心点的位移为可信位移;
所述组织应变生成模块,用于对所有可信位移区域进行平滑微分梯度计算,得到所述被测生物组织的二维应变分布结果。
综上,本发明实施例提供的准静态超声弹性成像位移计算装置,关于该装置的更详细的说明,请参考上述方法实施例,该装置通过采用上述技术方案,取得了以下技术效果:
1、本发明技术方案中,可通过互相关可靠性指数值来判断迭代互相关计算的结果是否可信,根据迭代互相关结果是否可信,自适应的决定:是通过迭代互相关结果和先验位移来计算位移计算中心点的时移,还是进一步通过二维综合互相关搜索计算位移计算中心点的时移。
2、在判断迭代互相关结果可信,利用迭代互相关结果计算纵向时移时,由于迭代互相关函数是基于邻域窗建立的,因而能够快速检测纵向方向上发生的时移,同时,能够降低横向时移对纵向时移计算的影响,从而可以提高纵向时移的准确度。
3、在判断迭代互相关结果不可信时,能够自适应的使用二维综合互相关搜索,能够准确检测纵向方向上和横向方向上发生的位移。
4、在判断迭代互相关结果为不可信时,使用邻域内的二维综合互相关搜索,而不是根据迭代互相关结果,来得到位移计算中心点的时移,可以将不可靠的计算结果剔除,可避免后续连锁错误;
5、在判断迭代互相关结果为可信时,可以直接根据先验时移确定位移计算中心点的横向时移,可根据先验时移和迭代互相关结果计算得到位移计算中心点的纵向时移,不必进行二维综合互相关搜索,因而大大减少了计算量,有利于实现实时成像。
在上述实施例中,对各个实施例的描述都各有侧重,某个实施例中没有详细描述的部分,可以参见其它实施例的相关描述。
需要说明的是,对于前述的各方法实施例,为了简单描述,故将其都表述为一系列的动作组合,但是本领域技术人员应该知悉,本发明并不受所描述动作顺序的限制,因为依据本发明,某些步骤可以采用其它顺序或者同时进行。其次,本领域技术人员也应该知悉,说明书中所描述的实施例均属于优选实施例,所涉及的动作和模块并不一定是本发明所必须的。
以上对本发明实施例所提供的准静态超声弹性成像位移计算方法和装置进行了详细介绍,但以上实施例的说明只是用于帮助理解本发明的方法及其核心思想,不应理解为对本发明的限制。本技术领域的技术人员,依据本发明的思想,在本发明揭露的技术范围内,可轻易想到的变化或替换,都应涵盖在本发明的保护范围之内。

Claims (16)

1.一种准静态超声弹性成像位移计算方法,其特征在于,包括:
采集得到连续的多帧超声射频信号;
对于两帧超声射频信号上的任一个位移计算中心点,利用先验时移建立基于邻域窗的迭代互相关函数,进行迭代互相关计算;
获取所述迭代互相关计算的互相关可靠性指数值,判断所述互相关可靠性指数值是否小于第一可信阈值;
如果所述互相关可靠性指数值大于或等于第一可信阈值,根据所述先验时移确定所述位移计算中心点的横向时移,根据所述先验时移和所述迭代互相关计算的迭代互相关结果确定所述位移计算中心点的纵向时移;
如果所述互相关可靠性指数值小于第一可信阈值,对所述位移计算中心点使用邻域内的二维综合互相关搜索,根据邻域内的二维综合互相关系数最大时的二维综合互相关结果,得到所述位移计算中心点的纵向时移和横向时移;
根据遍历全部超声射频信号得到的所有位移计算中心点的纵向时移和横向时移,生成所述被测生物组织的位移场分布。
2.根据权利要求1所述的方法,其特征在于,所述对于两帧超声射频信号上的任一个位移计算中心点,利用先验时移建立基于邻域窗的迭代互相关函数包括:
记先验时移为[tp,tq],建立如下的迭代互相关函数:
R ( t w 0 , t h 0 ) = Σ t h = t h 0 - h / 2 t h 0 + h / 2 Σ t w = t w 0 - w / 2 t w 0 + w / 2 x ( t w , t h ) y ( t w - t p , t h - t q ) * = Ae jω τ w , τ h ,
其中,x()和y()分别表示所述两帧超声射频信号,tw和th是所述超声射频信号在邻域窗内的横向坐标和纵向坐标,[tw0,th0]是第一帧超声射频信号x()上的所述位移计算中心点的坐标,[tw0-tp,th0-tq]是第二帧超声射频信号y()上对应的位移计算中心点的坐标,[w,h]为邻域窗的宽度和高度,R(tw0,th0)表示[tw0,th0]的邻域窗内的数据与[tw0-tp,th0-tq]的邻域窗内的数据的互相关函数,“*”表示共轭运算,A表示迭代互相关结果的幅度,为迭代互相关结果的相位。
3.根据权利要求2所述的方法,其特征在于,所述获取所述迭代互相关计算的互相关可靠性指数值包括:
采用公式计算所述互相关可靠性指数值Qc,其中,
R sum = Σ t h = t h 0 - h / 2 t h 0 + h / 2 Σ t w = t w 0 - w / 2 t w 0 + w / 2 | | x ( t w , t h ) | + | y ( t w - t p , t h - t q ) | | ,
R sub = Σ t h = t h 0 - h / 2 t h 0 + h / 2 Σ t w = t w 0 - w / 2 t w 0 + w / 2 | | x ( t w , t h ) | - | y ( t w - t p , t h - t q ) | | .
4.根据权利要求2所述的方法,其特征在于,所述获取所述迭代互相关计算的互相关可靠性指数值包括:
采用以下公式计算所述互相关可靠性指数值Qc:
Qc = | R ( t w 0 , t h 0 ) | R ( t w , t h ) R ( t w - t p , t h - t q ) ,
其中,
R ( t w , t h ) = Σ t h = t h 0 - h / 2 t h 0 + h / 2 Σ t w = t w 0 - w / 2 t w 0 + w / 2 x ( t w , t h ) 2 ,
R ( t w - t p , t h - t q ) = Σ t h = t h 0 - h / 2 t h 0 + h / 2 Σ t w = t w 0 - w / 2 t w 0 + w / 2 y ( t w - t p , t h - t q ) 2 .
5.根据权利要求2至4中任一所述的方法,其特征在于,所述根据所述先验时移确定所述位移计算中心点的横向时移,根据所述先验时移和所述迭代互相关计算的迭代互相关结果确定所述位移计算中心点的纵向时移包括:
根据所述迭代互相关结果的相位计算则,所述位移计算中心点的时移为:Uh=Ps+th,Uw=tw
其中,Uh为纵向时移,Uw为横向时移,λ是频率为发射脉冲中心频率的超声波在所述被测生物组织中传播的波长。
6.根据权利要求1至4中任一所述的方法,其特征在于,所述对所述位移计算中心点使用邻域内的二维综合互相关搜索,根据邻域内的二维综合互相关系数最大时的二维综合互相关结果,得到所述位移计算中心点的纵向时移和横向时移包括:
记二维综合互相关搜索的范围是[W,H],二维综合互相关搜索的公式为:
R ( t w 0 , t h 0 , τ w , τ h ) = Σ t h = t h 0 - h / 2 t h 0 + h / 2 Σ t w = t w 0 - w / 2 t w 0 + w / 2 x ( t w 0 , t h 0 ) y ( t w 0 - τ w , t h 0 - τ h ) * = Ae jω τ w , τ h ,
其中,x()和y()分别表示所述的两帧超声射频信号,[tw0,th0]是第一帧超声射频信号x()上的所述位移计算中心点的坐标,[τwh]是第二帧超声射频信号y()的搜索范围内偏离所述位移计算中心点的横向偏移量和纵向偏移量,*表示共轭运算,A表示互相关结果的幅度,为二维综合互相关结果的相位,τw在[-W/2,W/2]范围内取值,τh在[-H/2,H/2]范围内取值;
二维综合互相关系数为:
Cr ( t w 0 , t h 0 , τ w , τ h ) = | R ( t w 0 , t h 0 , τ w , τ h ) | R ( t w 0 , t h 0 , τ w , τ h ) R ( t w 0 - τ w , t h 0 - τ h ) ,
其中, R ( t w 0 , t h 0 ) = Σ t h = t h 0 - h / 2 t h 0 + h / 2 Σ t w = t w 0 - w / 2 t w 0 + w / 2 x ( t w , t h ) 2 ,
R ( t w 0 - τ w , t h 0 - τ h ) = Σ t h = t h 0 - h / 2 t h 0 + h / 2 Σ t w = t w 0 - w / 2 t w 0 + w / 2 y ( t w - τ w , t h - τ h ) 2 ,
其中,“||”表示复数取模或实数取绝对值运算;
记二维综合互相关系数最大时,对应的偏移量是[τwMhM],对应的相角为则,
所述位移计算中心点的时移如下:
Uh=Ps+τhM,Uw=τwM
其中,Uh为二维综合互相关计算结果的纵向时移,Uw为二维综合互相关计算结果的横向时移,λ是频率为发射脉冲中心频率的超声波在所述被测生物组织中传播的波长。
7.根据权利要求1所述的方法,其特征在于,每一帧所述多帧超声射频信号包括多条扫描线,所述探测步骤中,纵向方向上每一条扫描线的数据的采样率是所述超声探头发射脉冲的中心频率的至少2倍。
8.根据权利要求1所述的方法,其特征在于,
以同一根扫描线上、纵向方向上的前一个位移计算中心点的时移,作为所述先验时移。
9.根据权利要求1所述的方法,其特征在于,所述方法还包括:
判断所述位移场分布中每一个位移计算中心点的互相关可靠性指数是否小于第二可信阈值,如果不小于,判断所述位移计算中心点的位移为可信位移;
对所有可信位移区域进行平滑微分梯度计算,得到所述被测生物组织的二维应变分布结果。
10.一种准静态超声弹性成像位移计算装置,其特征在于,包括:
采集模块,用于采集得到连续的多帧超声射频信号;
迭代互相关模块,用于对于两帧超声射频信号上的任一个位移计算中心点,利用先验时移建立基于邻域窗的迭代互相关函数,进行迭代互相关计算;
可靠性判断模块,用于获取所述迭代互相关计算的互相关可靠性指数值,判断所述互相关可靠性指数值是否小于第一可信阈值;
第一时移计算模块,用于如果所述互相关可靠性指数值大于或等于第一可信阈值,根据所述先验时移确定所述位移计算中心点的横向时移,根据所述先验时移和所述迭代互相关计算的迭代互相关结果确定所述位移计算中心点的纵向时移;
第二时移计算模块,用于如果所述互相关可靠性指数值小于第一可信阈值,对所述位移计算中心点使用邻域内的二维综合互相关搜索,根据邻域内的二维综合互相关系数最大时的二维综合互相关结果,得到所述位移计算中心点的纵向时移和横向时移;
位移场生成模块,用于根据遍历全部超声射频信号得到的所有位移计算中心点的纵向时移和横向时移,生成所述被测生物组织的位移场分布。
11.根据权利要求10所述的装置,其特征在于,所述迭代互相关模块,具体用于:
记先验时移为[tp,tq],建立如下的迭代互相关函数:
R ( t w 0 , t h 0 ) = Σ t h = t h 0 - h / 2 t h 0 + h / 2 Σ t w = t w 0 - w / 2 t w 0 + w / 2 x ( t w , t h ) y ( t w - t p , t h - t q ) * = Ae jω τ w , τ h ,
其中,x()和y()分别表示所述两帧超声射频信号,tw和th是所述超声射频信号在邻域窗内的横向坐标和纵向坐标,[tw0,th0]是第一帧超声射频信号x()上的所述位移计算中心点的坐标,[tw0-tp,th0-tq]是第二帧超声射频信号y()上对应的位移计算中心点的坐标,[w,h]为邻域窗的宽度和高度,R(tw0,th0)表示[tw0,th0]的邻域窗内的数据与[tw0-tp,th0-tq]的邻域窗内的数据的互相关函数,“*”表示共轭运算,A表示迭代互相关结果的幅度,为迭代互相关结果的相位。
12.根据权利要求11所述的方法,其特征在于,所述可靠性判断模块具体用于:
采用公式计算所述互相关可靠性指数值Qc,其中,
R sum = Σ t h = t h 0 - h / 2 t h 0 + h / 2 Σ t w = t w 0 - w / 2 t w 0 + w / 2 | | x ( t w , t h ) | + | y ( t w - t p , t h - t q ) | | ,
R sub = Σ t h = t h 0 - h / 2 t h 0 + h / 2 Σ t w = t w 0 - w / 2 t w 0 + w / 2 | | x ( t w , t h ) | - | y ( t w - t p , t h - t q ) | | .
13.根据权利要求11所述的方法,其特征在于,所述可靠性判断模块具体用于:
采用下述公式计算所述互相关可靠性指数值Qc:
Qc = | R ( t w 0 , t h 0 ) | R ( t w , t h ) R ( t w - t p , t h - t q ) ,
其中,
R ( t w , t h ) = Σ t h = t h 0 - h / 2 t h 0 + h / 2 Σ t w = t w 0 - w / 2 t w 0 + w / 2 x ( t w , t h ) 2 ,
R ( t w - t p , t h - t q ) = Σ t h = t h 0 - h / 2 t h 0 + h / 2 Σ t w = t w 0 - w / 2 t w 0 + w / 2 y ( t w - t p , t h - t q ) 2 .
14.根据权利要求11至13中任一所述的装置,其特征在于,所述第一时移计算模块具体用于:
根据所述迭代互相关结果的相位计算则,所述位移计算中心点的时移为:Uh=Ps+th,Uw=tw
其中,Uh为纵向时移,Uw为横向时移,λ是频率为发射脉冲中心频率的超声波在所述被测生物组织中传播的波长。
15.根据权利要求10至13中任一所述的装置,其特征在于,所述第二时移计算模块具体用于:
记二维综合互相关搜索的范围是[W,H],二维综合互相关搜索的公式为:
R ( t w 0 , t h 0 , τ w , τ h ) = Σ t h = t h 0 - h / 2 t h 0 + h / 2 Σ t w = t w 0 - w / 2 t w 0 + w / 2 x ( t w 0 , t h 0 ) y ( t w 0 - τ w , t h 0 - τ h ) * = Ae jω τ w , τ h ,
其中,x()和y()分别表示所述的两帧超声射频信号,[tw0,th0]是第一帧超声射频信号x()上的所述位移计算中心点的坐标,[τwh]是第二帧超声射频信号y()的搜索范围内偏离所述位移计算中心点的横向偏移量和纵向偏移量,*表示共轭运算,A表示互相关结果的幅度,为二维综合互相关结果的相位,τw在[-W/2,W/2]范围内取值,τh在[-H/2,H/2]范围内取值;
二维综合互相关系数为:
Cr ( t w 0 , t h 0 , τ w , τ h ) = | R ( t w 0 , t h 0 , τ w , τ h ) | R ( t w 0 , t h 0 , τ w , τ h ) R ( t w 0 - τ w , t h 0 - τ h ) ,
其中, R ( t w 0 , t h 0 ) = Σ t h = t h 0 - h / 2 t h 0 + h / 2 Σ t w = t w 0 - w / 2 t w 0 + w / 2 x ( t w , t h ) 2 ,
R ( t w 0 - τ w , t h 0 - τ h ) = Σ t h = t h 0 - h / 2 t h 0 + h / 2 Σ t w = t w 0 - w / 2 t w 0 + w / 2 y ( t w - τ w , t h - τ h ) 2 ,
其中,“||”表示复数取模或实数取绝对值运算;
记二维综合互相关系数最大时,对应的偏移量是[τwMhM],对应的相角为则,
所述位移计算中心点的时移如下:
Uh=Ps+τhM,Uw=τwM
其中,Uh为二维综合互相关计算结果的纵向时移,Uw为二维综合互相关计算结果的横向时移,λ是频率为发射脉冲中心频率的超声波在所述被测生物组织中传播的波长。
16.根据权利要求10所述的装置,其特征在于,
所述可靠性判断模块,还用于判断所述位移场分布中每一个位移计算中心点的互相关可靠性指数是否小于第二可信阈值,如果不小于,判断所述位移计算中心点的位移为可信位移;
所述装置还包括:组织应变生成模块,用于对所有可信位移区域进行平滑微分梯度计算,得到所述被测生物组织的二维应变分布结果。
CN201410800079.6A 2014-12-19 2014-12-19 准静态超声弹性成像位移计算方法和装置 Active CN105748100B (zh)

Priority Applications (1)

Application Number Priority Date Filing Date Title
CN201410800079.6A CN105748100B (zh) 2014-12-19 2014-12-19 准静态超声弹性成像位移计算方法和装置

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
CN201410800079.6A CN105748100B (zh) 2014-12-19 2014-12-19 准静态超声弹性成像位移计算方法和装置

Publications (2)

Publication Number Publication Date
CN105748100A true CN105748100A (zh) 2016-07-13
CN105748100B CN105748100B (zh) 2019-04-16

Family

ID=56340998

Family Applications (1)

Application Number Title Priority Date Filing Date
CN201410800079.6A Active CN105748100B (zh) 2014-12-19 2014-12-19 准静态超声弹性成像位移计算方法和装置

Country Status (1)

Country Link
CN (1) CN105748100B (zh)

Cited By (5)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN107970043A (zh) * 2017-12-28 2018-05-01 深圳开立生物医疗科技股份有限公司 一种剪切波的检测方法及装置
CN109512463A (zh) * 2018-10-16 2019-03-26 深圳大学 超声弹性成像位移估计方法、系统、终端及可读存储介质
CN110327077A (zh) * 2019-07-09 2019-10-15 深圳开立生物医疗科技股份有限公司 一种血流显示方法、装置及超声设备和存储介质
JP2020003247A (ja) * 2018-06-26 2020-01-09 公益財団法人鉄道総合技術研究所 波形データの高精度位置補正方法及びシステム
WO2021103493A1 (zh) * 2019-11-27 2021-06-03 深圳开立生物医疗科技股份有限公司 一种基于剪切波的成像方法、系统及装置

Citations (5)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
US20050074153A1 (en) * 2003-09-30 2005-04-07 Gianni Pedrizzetti Method of tracking position and velocity of objects' borders in two or three dimensional digital images, particularly in echographic images
CN102860842A (zh) * 2012-09-26 2013-01-09 浙江大学 一种实时准静态超声弹性成像方法
CN103735287A (zh) * 2013-12-05 2014-04-23 中国科学院苏州生物医学工程技术研究所 一种血管内超声弹性成像二维多级混合位移估计方法
CN103845081A (zh) * 2012-11-28 2014-06-11 深圳迈瑞生物医疗电子股份有限公司 超声弹性成像系统和方法、实时动态帧间处理方法
CN104188689A (zh) * 2013-10-25 2014-12-10 中国科学院深圳先进技术研究院 基于超声回波射频信号的组织位移估算方法和系统

Patent Citations (5)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
US20050074153A1 (en) * 2003-09-30 2005-04-07 Gianni Pedrizzetti Method of tracking position and velocity of objects' borders in two or three dimensional digital images, particularly in echographic images
CN102860842A (zh) * 2012-09-26 2013-01-09 浙江大学 一种实时准静态超声弹性成像方法
CN103845081A (zh) * 2012-11-28 2014-06-11 深圳迈瑞生物医疗电子股份有限公司 超声弹性成像系统和方法、实时动态帧间处理方法
CN104188689A (zh) * 2013-10-25 2014-12-10 中国科学院深圳先进技术研究院 基于超声回波射频信号的组织位移估算方法和系统
CN103735287A (zh) * 2013-12-05 2014-04-23 中国科学院苏州生物医学工程技术研究所 一种血管内超声弹性成像二维多级混合位移估计方法

Non-Patent Citations (2)

* Cited by examiner, † Cited by third party
Title
SHAOGUO CUI ET AL: "A Robust Phase Zero Estimator for Ultrasonic Elastography Using Quality-", 《 ICALIP2010 》 *
崔少国 等: "基于先验估计和相位相关的实时超声弹性成像", 《中国生物医学工程学报》 *

Cited By (8)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN107970043A (zh) * 2017-12-28 2018-05-01 深圳开立生物医疗科技股份有限公司 一种剪切波的检测方法及装置
JP2020003247A (ja) * 2018-06-26 2020-01-09 公益財団法人鉄道総合技術研究所 波形データの高精度位置補正方法及びシステム
JP2023022130A (ja) * 2018-06-26 2023-02-14 公益財団法人鉄道総合技術研究所 波形データの高精度位置補正方法及びシステム
JP7446698B2 (ja) 2018-06-26 2024-03-11 公益財団法人鉄道総合技術研究所 波形データの高精度位置補正方法及びシステム
CN109512463A (zh) * 2018-10-16 2019-03-26 深圳大学 超声弹性成像位移估计方法、系统、终端及可读存储介质
CN109512463B (zh) * 2018-10-16 2021-12-03 深圳大学 超声弹性成像位移估计方法、系统、终端及可读存储介质
CN110327077A (zh) * 2019-07-09 2019-10-15 深圳开立生物医疗科技股份有限公司 一种血流显示方法、装置及超声设备和存储介质
WO2021103493A1 (zh) * 2019-11-27 2021-06-03 深圳开立生物医疗科技股份有限公司 一种基于剪切波的成像方法、系统及装置

Also Published As

Publication number Publication date
CN105748100B (zh) 2019-04-16

Similar Documents

Publication Publication Date Title
Feigin et al. A deep learning framework for single-sided sound speed inversion in medical ultrasound
CN105748100A (zh) 准静态超声弹性成像位移计算方法和装置
CN102469980B (zh) 空间上精细的剪切波分散超声振动测定采样
CN103327904B (zh) 超声波摄像装置和超声波摄像方法
JP6063553B2 (ja) 超音波イメージング方法及び超音波イメージング装置
US20110313291A1 (en) Medical image processing device, medical image processing method, medical image diagnostic apparatus, operation method of medical image diagnostic apparatus, and medical image display method
Vogt et al. Development and evaluation of a high-frequency ultrasound-based system for in vivo strain imaging of the skin
CN106680825B (zh) 一种声学阵列成像系统与方法
CN102695457A (zh) 超声波诊断装置及测量内中膜的厚度的方法
CN110432926B (zh) 弹性测量检测方法及系统
Tong et al. Deep learning inversion with supervision: A rapid and cascaded imaging technique
EP3459082B1 (en) Method of, and apparatus for, non-invasive medical imaging using waveform inversion
CN102764141B (zh) 弹性成像方法和系统及其中的生物组织位移估计方法和系统
JP5692079B2 (ja) 変位推定方法、変位推定装置
CN109561883A (zh) 超声波诊断装置
CN107616814A (zh) 一种生物组织剪切波波速测量方法及医用超声波设备
CN103735287B (zh) 一种血管内超声弹性成像二维多级混合位移估计方法
CN112022215B (zh) 超声弹性成像角膜检测方法、装置、系统和存储介质
CN104739442A (zh) 压力弹性成像位移检测方法、装置和超声成像设备
CN102860842A (zh) 一种实时准静态超声弹性成像方法
CN104622507A (zh) 弹性模量测量方法和系统
Sun et al. Trajectory-based deformation correction in ultrasound images
CN112674791B (zh) 肌肉超声弹性成像的优化方法及系统
JP4761999B2 (ja) 超音波診断装置およびその画像処理方法、その画像処理プログラム
CN109589138A (zh) 一种剪切波波速计算方法及弹性成像设备

Legal Events

Date Code Title Description
C06 Publication
PB01 Publication
C10 Entry into substantive examination
SE01 Entry into force of request for substantive examination
GR01 Patent grant
GR01 Patent grant