CN112764040B - 一种基于射线理论相位修正的合成孔径波束形成方法 - Google Patents
一种基于射线理论相位修正的合成孔径波束形成方法 Download PDFInfo
- Publication number
- CN112764040B CN112764040B CN201911060860.3A CN201911060860A CN112764040B CN 112764040 B CN112764040 B CN 112764040B CN 201911060860 A CN201911060860 A CN 201911060860A CN 112764040 B CN112764040 B CN 112764040B
- Authority
- CN
- China
- Prior art keywords
- sound barrier
- ultrasonic
- imaging
- sound
- grid
- 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.)
- Active
Links
Images
Classifications
-
- G—PHYSICS
- G01—MEASURING; TESTING
- G01S—RADIO DIRECTION-FINDING; RADIO NAVIGATION; DETERMINING DISTANCE OR VELOCITY BY USE OF RADIO WAVES; LOCATING OR PRESENCE-DETECTING BY USE OF THE REFLECTION OR RERADIATION OF RADIO WAVES; ANALOGOUS ARRANGEMENTS USING OTHER WAVES
- G01S15/00—Systems using the reflection or reradiation of acoustic waves, e.g. sonar systems
- G01S15/88—Sonar systems specially adapted for specific applications
- G01S15/89—Sonar systems specially adapted for specific applications for mapping or imaging
- G01S15/8902—Side-looking sonar
- G01S15/8904—Side-looking sonar using synthetic aperture techniques
-
- G—PHYSICS
- G01—MEASURING; TESTING
- G01S—RADIO DIRECTION-FINDING; RADIO NAVIGATION; DETERMINING DISTANCE OR VELOCITY BY USE OF RADIO WAVES; LOCATING OR PRESENCE-DETECTING BY USE OF THE REFLECTION OR RERADIATION OF RADIO WAVES; ANALOGOUS ARRANGEMENTS USING OTHER WAVES
- G01S7/00—Details of systems according to groups G01S13/00, G01S15/00, G01S17/00
- G01S7/52—Details of systems according to groups G01S13/00, G01S15/00, G01S17/00 of systems according to group G01S15/00
- G01S7/52003—Techniques for enhancing spatial resolution of targets
Landscapes
- Engineering & Computer Science (AREA)
- Radar, Positioning & Navigation (AREA)
- Remote Sensing (AREA)
- Physics & Mathematics (AREA)
- Computer Networks & Wireless Communication (AREA)
- General Physics & Mathematics (AREA)
- Acoustics & Sound (AREA)
- Ultra Sonic Daignosis Equipment (AREA)
Abstract
本发明提供了一种基于射线理论相位修正的合成孔径波束形成方法,用于解决复杂介质超声成像问题,包括如下步骤:步骤S1,以预定的采样频率依次发射超声波信号并获取采集信号;步骤S2,获取材料的结构与其中的声速;步骤S3,成像区域空间离散化为二维方形网格集群;步骤S4,根据离散化网格,计算阵元到每个网格端点的初至时间,得到超声初至时间矩阵;步骤S5,根据超声初至时间矩阵,替代常规波数形成所需的时间矩阵,进行图像重建处理。本发明的方法有效地提高了图像重建的分辨率,且能够保持不同探测深度分辨率的一致性。
Description
技术领域
本发明属于超声聚焦成像领域,具体涉及一种基于射线理论相位修正的合成孔径波束形成方法。
背景技术
近年来,由于超声技术的诸多优点,广泛应用于无损检测、医学成像等领域。透过复杂介质的声屏障,进行超声聚焦成像以获取准确的结构信息较困难。例如声屏障与软组织的声速存在差异,导致超声波束合成的聚焦点偏移,进而导致颅骨内部软组织成像或脑功能成像出现误差。(可参考文献:Marsac L,Chauvet D,La Greca R,et al..Ex vivooptimisation of a heterogeneous speed of sound model of the human skull fornon-invasive transcranial focused ultrasound at 1MHz[J].Int J Hyperther,2017,33(1):635-645)
近年来,学者提出若干方法来解决相位畸变导致的超声图形质量下降的问题。例如在经颅成像领域,对超声接收信号做互相关处理,根据最小二乘法实现相位修正(可参考文献:Lindsey B D,Smith S W.Pitch-catch phase aberration correction ofmultipleisoplanatic patches for3-D transcranial ultrasound imaging.IEEE T UltrasonFerr,2013,60(3):463-480),或者根据颅骨的CT或MRI数据进行全波仿真以预估相位误差(可参考文献:Almquist S,Parker D L,Christensen D A.Rapid full-wave phaseaberration correction method for transcranial high-intensity focusedultrasound therapies.J Thera Ultrasound,2016,4(1):30-41)。
现有的方法具有以下缺点:1:互相关法无法分辨成像区域中距离较小的散射体;2:全波仿真方法消耗时间长,不利于高效的时延修正。
发明内容
本发明是为了解决上述问题而进行的,目的在于提供一种用于解决复杂介质超声成像问题的基于射线理论相位修正的合成孔径波束形成方法。
本发明提供了一种基于射线理论相位修正的合成孔径波束形成方法,具有这样的特征,包括如下步骤:
步骤S1,采用发射阵元总数为L的线性阵列,单个发射阵元以预定的采样频率向声屏障内部发射超声波信号,并通过全部L个接收阵元获取与发射阵元相对应的采集信号,L个发射阵元依次循环发射,然后进入步骤S2;
步骤S2,获取声屏障的结构与声速,然后进入步骤S3;
步骤S3,成像区域空间离散化为二维方形网格集群M×N,并将网格的端点作为成像的离散点,在成像区域内确定声屏障及其他介质所在区域Ω,根据声屏障以及其他复杂介质在空间结构及声速c(m,n)|m<M,n<N,计算对应的空间慢度分布S(m,n)|m<M,n<N=1/c(m,n)|m<M,n<N,然后进入步骤S4;
步骤S4,将发射阵元与接收阵元设定为点源,确定每个阵元在二维方形网格集群中的位置,分别以这些点源为起点,对程函方程进行离散化求解,计算其他网格端点处的超声初至时间t(m,n)|m<M,n<N,而后得到超声初至时间矩阵,然后进入步骤S5;
步骤S5,将步骤S4得到的超声初至时间矩阵替代常规波数形成所需的时间矩阵,进行图像重建处理。
在本发明提供的基于射线理论相位修正的合成孔径波束形成方法中,还可以具有这样的特征:其中,步骤S1中接收阵元的数量为L个,在接收采集信号时,每个接收阵元均为有效阵元,从而形成L次发射接收循环。
在本发明提供的基于射线理论相位修正的合成孔径波束形成方法中,还可以具有这样的特征:其中,步骤S2中声屏障的结构从X射线CT图像或者核磁共振MRI图像中获取,声屏障的声速通过超声测量获得。
在本发明提供的基于射线理论相位修正的合成孔径波束形成方法中,还可以具有这样的特征:其中,步骤S2中声屏障的结构和声速均通过地震层析反演技术估计获取。
在本发明提供的基于射线理论相位修正的合成孔径波束形成方法中,还可以具有这样的特征:其中,步骤S2中,在获取声速数据的前提下,通过如下步骤获取声屏障的结构:
步骤S2-1,利用成像的方法获取声屏障的上表面的形状及空间位置,假定上表面以下是声屏障的速度,进入步骤S2-2;
步骤S2-2,在S2-1的基础上,重复步骤S3至S5,可以获得声屏障的下表层的形状及空间位置。
在本发明提供的基于射线理论相位修正的合成孔径波束形成方法中,还可以具有这样的特征:其中,步骤S3中的网格是正方形网格,该网格的边长为h,S表示速度慢度空间分布。
在本发明提供的基于射线理论相位修正的合成孔径波束形成方法中,还可以具有这样的特征:其中,步骤S4中,离散化求解所述程函方程的方法为波前外推法、快速行进算法或快速扫掠算法,
当采用快速扫掠算法求解所述程函方时具体包括如下子步骤:
步骤S4-1,对整个所述成像区域的初至时间矩阵(M+2)×(N+2)预先设置大的初始时间,并将所述阵元所在空间点的初始时间设置为零,然后进行从左下角向右上角,即i=2:M+1,j=2:N+1的迭代循环;
步骤S4-2,为求解t(i,j),预先对比该点两侧时间最小值xmin=min[t(i-1,j),t(i+1,j)]和ymin=min[t(i,j-1),t(i,j+1)],假如t(i-1,j)=t(i+1,j)或t(i,j-1)=t(i,j+1),此处选择t(i-1,j)或者t(i,j-1);
步骤S4-3,在步骤S4-2的基础上,寻找t(i,j)对角点的局部最小值tmin,该最小值tmin为:
tmin=min[t(i-1,j-1),t(i+1,j+1),t(i-1,j+1),t(i+1,j-1)];
步骤S4-4,在步骤S4-2与步骤S4-3的基础上,当abs(xmin-ymin)>h/S(i,j)时,则t1(i,j)=min(xmin,ymin)+h/S(i,j),当abs(xmin-ymin)≤h/S(i,j)时,则
步骤S4-6,在步骤S4-4与步骤S4-5的基础上,获得此时该点的最终走时t(i,j)=min[t1(i,j),t2(i,j),t_last(i,j)],而后重复步骤S4-1,直至空间所有的初至走时被计算到,而后继续进行迭代直至t(i,j)-t_last(i,j)<Threshold,则结束迭代,
其中,i表示笛卡尔坐标系下水平方向离散点的空间坐标,j表示笛卡尔坐标系下垂直方向离散点的空间坐标,h表示正方形网格的边长,S表示速度慢度空间分布,t_last(i,j)表示上一次迭代计算的t(i,j)。
在本发明提供的基于射线理论相位修正的合成孔径波束形成方法中,还可以具有这样的特征:其中,迭代循环的方式还包括从右下角向左上角,即i=M+1:2,j=2:N+1,从左上角向右下角,即i=2:M+1,j=N+1:2,以及从右上角向左下角,即i=M+1:2,j=N+1:2。
在本发明提供的基于射线理论相位修正的合成孔径波束形成方法中,还可以具有这样的特征:其中,声屏障的厚度划分为至少两个网格宽度,即2h<H,H表示声屏障的最小厚度。
发明的作用与效果
根据本发明所涉及的基于射线理论相位修正的合成孔径波束形成方法,结合了合成孔径成像技术的发射接收模式及其相应的时域重建方法即延迟聚焦波束形成法,可以获取修正后的成像图形,并且由于采用了合成孔径的超声波发射-接收方式,利用其发射-接收具有动态聚焦的特点,有效地提高了图像重建的分辨率,且能够保持不同探测深度分辨率的一致性。
附图说明
图1是本发明的实施例中颅骨等声屏障及颅内软组织结构示意图;
图2是本发明的实施例中合成孔径实际焦点与虚拟焦点对比图;
图3是本发明的实施例中基于射线理论相位修正的合成孔径波束形成方法的流程框图;
图4是本发明的实施例中阵列探头第一单元(a)与第128单元(b)的初至时间等高图;
图5是本发明的实施例中应用于经颅超声成像修正仿真结果的示意图;
图6是本发明的实施例中实验结果的示意图。
具体实施方式
为了使本发明实现的技术手段与功效易于明白了解,以下结合实施例及附图对本发明作具体阐述。
本发明的一种基于射线理论相位修正的合成孔径波束形成方法,包括如下步骤:
步骤S1,采用发射阵元总数为L的线性阵列,单个发射阵元以预定的采样频率向声屏障内部发射超声波信号,并通过全部L个接收阵元获取与发射阵元相对应的采集信号,L个发射阵元依次循环发射,然后进入步骤S2。
步骤S1中接收阵元的数量为L个,在接收采集信号时,每个接收阵元均为有效阵元,从而形成L次发射接收循环。
步骤S2,获取声屏障的结构与声速,然后进入步骤S3。
为了保证离散化求解程函方程的准确性,声屏障的厚度划分为至少两个网格宽度,即2h<H,H表示声屏障的最小厚度。
本发明中,步骤S2中声屏障的结构从X射线CT图像或者核磁共振MRI图像中获取,声屏障的声速通过独立的声屏障超声测量获得、文献数据或者其他技术手段获取。
本发明中,步骤S2中声屏障的结构和声速还可以均通过地震层析反演技术估计获取,可以参考文献:Jerome R.Krebs1,John E.Anderson et al.Fast full-wavefieldseismic inversion using encoded sources.Geophysics,2009.74(6),177-188。
步骤S2中,在获取声速数据的前提下,通过如下步骤获取声屏障的结构:
步骤S2-1,利用成像的方法获取声屏障的上表面的形状及空间位置,假定上表面以下是声屏障的速度,进入步骤S2-2;
步骤S2-2,在S2-1的基础上,重复步骤S3至S5,可以获得声屏障的下表层的形状及空间位置。
步骤S3,成像区域空间离散化为二维方形网格集群M×N,并将网格的端点作为成像的离散点,在成像区域内确定声屏障及其他介质所在区域Ω,根据声屏障以及其他复杂介质在空间结构及声速c(m,n)|m<M,n<N,计算对应的空间慢度分布S(m,n)|m<M,n<N=1/c(m,n)|m<M,n<N,然后进入步骤S4。
步骤S3中的网格是正方形网格,该网格的边长为h,S表示速度慢度空间分布。
步骤S4,将发射阵元与接收阵元设定为点源,确定每个阵元在二维方形网格集群中的位置,分别以这些点源为起点,对程函方程进行离散化求解,计算其他网格端点处的超声初至时间t(m,n)|m<M,n<N,而后得到超声初至时间矩阵,然后进入步骤S5。
本发明中,离散化求解程函方程是射线理论求解的核心。
步骤S4中,离散化求解所述程函方程的方法为波前外推法、快速行进算法或快速扫掠算法。
当采用快速扫掠算法求解所述程函方时具体包括如下子步骤:
步骤S4-1,对整个所述成像区域的初至时间矩阵(M+2)×(N+2)预先设置大的初始时间,并将所述阵元所在空间点的初始时间设置为零,然后进行从左下角向右上角,即i=2:M+1,j=2:N+1的迭代循环;
步骤S4-2,为求解t(i,j),预先对比该点两侧时间最小值xmin=min[t(i-1,j),t(i+1,j)]和ymin=min[t(i,j-1),t(i,j+1)],假如t(i-1,j)=t(i+1,j)或t(i,j-1)=t(i,j+1),此处优先选择t(i-1,j)或者t(i,j-1);
步骤S4-3,在步骤S4-2的基础上,寻找t(i,j)对角点的局部最小值tmin,该最小值tmin为:
tmin=min[t(i-1,j-1),t(i+1,j+1),t(i-1,j+1),t(i+1,j-1)];
步骤S4-4,在步骤S4-2与步骤S4-3的基础上,当abs(xmin,ymin)>h/S(i,j)时,则t1(i,j)=min(xmin,ymin)+h/S(i,j),当abs(xmin,ymin)≤h/S(i,j)时,则
步骤S4-6,在步骤S4-4与步骤S4-5的基础上,获得此时该点的最终走时t(i,j)=min[t1(i,j),t2(i,j),t_last(i,j)],而后重复步骤S4-1,直至空间所有的初至走时被计算到,而后继续进行迭代直至t(i,j)-t_last(i,j)<Threshold,则结束迭代。
其中,i表示笛卡尔坐标系下水平方向离散点的空间坐标,j表示笛卡尔坐标系下垂直方向离散点的空间坐标,h表示网格的边长,S表示速度慢度空间分布,t_last(i,j)表示上一次迭代计算的t(i,j)。
步骤S4-1中的迭代循环的方式还包括从右下角向左上角,即i=M+1:2,j=2:N+1,从左上角向右下角,即i=2:M+1,j=N+1:2以及从右上角向左下角,即i=M+1:2,j=N+1:2。
此外,步骤S4-1与步骤S3中的参数有如下关系:i=m+1,j=n+1。
步骤S5,将步骤S4得到的超声初至时间矩阵替代常规波数形成所需的时间矩阵,进行图像重建处理。
实施例:
图1是本发明的实施例中颅骨等声屏障及颅内软组织结构示意图。
在本实施例中,首先通过时域伪谱法(Pseudospectral time-domain method,PSTD)进行仿真,根据仿真结果判断本实施例的基于射线理论相位修正的合成孔径波束形成方法是具有合理性的,然后再进行仿体的实验。
仿真实验的声速模型为水-骨-软组织多介质模型,PSTD仿真相关参数设置如下:如图1所示,颅骨部分的纵波速度为cl=2462m/s,横波速度为cs=1129m/s,密度为1182kg/m3。软组织的材料参数在空间中遵循均匀随机分布,其中声速为c=1505±75m/s,密度分布为ρ=1000±5kg/m3。软组织是直径分别为4.5mm,9.0mm,12.5mm和15.0mm的14个圆。颅骨及软组织外的其他部分全部设置为水,其声速为1540m/s,密度为1182kg/m3。
离体实验采用光敏树脂材料打印3D颅骨模型,软组织材料是由水,琼脂,碳粉按照96∶3∶1的比例混合而成。经测量可知颅骨板的平均厚度约为4mm,光敏树脂材料的横纵波速度为cl=2462m/s,cs=1129m/s,损耗约为α=4dB/MHz/cm。
图2是本发明的实施例中合成孔径实际焦点与虚拟焦点对比图。
如图2所示,本实施例中超声波的发信元为线性超声换能器阵列,由128个发射阵元和128个接收阵元组成,每个发射阵元和接收阵元之间的间距均为0.675mm,发射脉冲为高斯包络正弦波,中心频率f0=1MHz,采样频率fs=4MHz。
采用传统波束形成方法,时延为Δti和Δtj时,相应的焦点位于A,而实际超声波经过折射后聚焦在位置B,证明在使用传统波束形成方法时,B处的软组织信息会错误的显示在A处,可以预料的是:对于实际聚焦点B,所有阵列振源到该聚焦点的时延都需要进行修正,才能利用波束形成获得正确的图像。
合成孔径技术最早应用于雷达系统,20世纪70年代起应用于超声成像。合成孔径超声作为一种新颖的医学超声成像方法,采用动态聚焦,提升整体探测区域的分辨率及对比度,已广泛用于血管内成像(可参考文献:M.Yu et al.IntravascularUltrasoundImaging With Virtual Source Synthetic Aperture Focusing and Coherence FactorWeighting[J].IEEE Trans Med Imaging,2017,36(10):2171-2178.)等软组织检测。合成孔径成像方法因其发射接收模式,适用于本发明提及的相位修正方法,并且其发射-接收具有动态聚焦特点,能提高重建图像的分辨率。
图3是本发明的实施例中基于射线理论相位修正的合成孔径波束形成方法的流程框图。
如图3所示,本实施例的基于射线理论相位修正的合成孔径波束形成方法,用于对上述穿透颅骨等声屏障进行成像,包括如下步骤:
步骤S1,采用发射阵元总数为128的线性阵列,单个发射阵元以预定的采样频率向颅骨或其他声屏障内部发射超声波信号,并通过全部128个接收阵元获取与发射阵元相对应的采集信号,128个发射阵元依次循环发射,然后进入步骤S2。
本实施例中,在发射超声波时,每个发射阵元依次作为有效发射阵元,超声波能够发射并穿透颅骨等声屏障进入颅内,与颅内软组织碰撞反射形成回波信号。
在接收采集信号时,每个接收阵元均为有效接收阵元,从而接收从颅骨及颅内软组织的反射的回波信号作为采集信号,形成N次发射接收循环。
步骤S2,获取颅骨的结构与声速,然后进入步骤S3。
如图1所示,在本实施例中,选取颅骨颞窗部位的CT图像,并使用光敏树脂材料3D打印了颅骨仿体,此外我们预先测量光敏树脂材料的声学参数,如横纵波速度、声损耗和密度等,因此在仿真与实验中进行相位修正必须的颅骨结构与声学参数均可得到。
步骤S3,成像区域空间离散化为二维方形网格集群M×N。
本实施例中的发射阵元与接收阵元均采用阵列探头单元,并且由于本实施例中的阵列探头单元之间的间距为0.675mm,考虑颅骨的平均厚度,因此设定网格的间距为阵列探头单元之间的间距的1/3,即0.225mm,考虑成像深度大约为15cm。因此,本实施例中网格集群中N=382(127*3+1),M设定为400。
步骤S4,将发射阵元与接收阵元设定为点源,确定每个阵元在网格集群中的位置,分别以这些点源为起点,根据波前外推法,快速行进算法或快速扫掠算法等方法,离散化求解程函方程,计算其他网格端点处的超声初至时间t(m,n)|m<M,n<N。
本实施例中,以第k个阵列探头单元为例,设定该阵列探头单元的空间位置的初至时间为0,然后按照离散快速扫掠算法,计算整个空间M×N点的修正初至时间tk(m,n),如此依次计算1至128点的初至时间矩阵,并存储该矩阵。
步骤S5,将步骤S4得到的初至时间矩阵替代常规波数形成所需的时间阵列,进行图像重建处理。
图4是本发明实施例中阵列探头第一单元(a)与第128单元(b)的初至时间等高图。
图5是本发明实施例中仿真结果的示意图,图5(a)是仿真中未考虑颅骨等声屏障内横波的前提下,使用传统波束形成重建的图像,图5(b)未考虑颅骨等声屏障的横波的前提下,使用本发明中相位修正后重建的图像,图5(c)是仿真中考虑颅骨等声屏障的横纵波的前提下,使用传统波束形成重建的图像,图5(d)仿真中考虑颅骨等声屏障的横纵波的前提下,使用本发明中相位修正后重建的图像。
如图5(a)图5(b)所示,若只考虑颅骨的纵波信息进行时延修正,使用液体模型进行仿真时,颅骨及颅内软组织显示清晰,颅骨形状、大小、空间位置等与实际情况基本一致。当使用固体模型进行仿真时,贴近颅骨下表面的软组织1号和2号被伪影掩盖,且成像区域左侧软组织3号、6号和11号的左侧部分出现伪影。这是由于经颅超声的纵-横-纵模式转换,不适用于基于纵波速度的时延修正成像。
如图5(c)和图5(d)所示,入射波与骨表面夹角较小时,颅骨内横波能量较低,因此纵-横-纵模式转换只能掩盖住较小深度的软组织,如1号软组织何2号软组织;入射波与骨表面夹角较大时,颅骨内纵波能量较低,因此不能实现较深区域的软组织成像,如3号软组织、6号软组织以及11号软组织。
图6是本发明实施例中实验结果的示意图,图6(a)是实验信号经过传统波束形成重建的图像,图6(b)是实验信号经过传统波束形成经过维纳滤波后重建的图像,图6(c)是实验信号经过相位修正后重建的图像,图6(d)是实验信号经过相位修正且维纳滤波后重建的图像。
如图6所示,受限于低频超声相控阵的低成像分辨率以及噪声信号的影响,图像只能显示软组织的上下面,内部的成像细节较差。贴近颅骨下表面的软组织被伪影掩盖,且左侧部分的软组织伪影的形状相似;传统波束形成成像的软组织畸变以及时延修正波束形成成像的软组织信息与仿真相似。
图6(a)和图6(c)分别对可见软组织与真实软组织的空间信息进行了对比,并结合实验及仿真结果计算了偏移位移和图像的椭圆化程度。通过时延修正波束形成,颅内可见软组织成像的圆心位置平均误差从1.98mm下降至0.21mm,椭圆的平均离心率从0.63下降到0.19,最高最低亮度之差的平均值从37.36dB上升到42.41dB。
实施例的作用与效果
由实施例可知,本实施例的基于射线理论相位修正的合成孔径波束形成方法,结合了合成孔径成像方法,并且颅内软组织的成像质量获得提高,其主要表现在:软组织的空间位置得到纠正,即软组织成像的圆心位置平均误差从2.06mm下降至0.21mm,软组织成像的图像畸变减弱,即椭圆的平均离心率从0.62下降到0.19以及软组织与背景环境的可识别度增加,即最高最低亮度之差的平均值从37.36dB上升到42.40dB。
此外,由于成像阵列探头单元之间的间距与其相应波长的量级相似,并且成像区域距离阵列探头单元的距离约等于五倍波长,阵列探头单元可以当作点源,因此,通过射线追踪技术,可以获取从相控阵任意阵列探头单元到达空间所有空间成像点即步骤S4中的网格端点的时延修正量t。
因此,该方法结合了合成孔径成像技术的发射接收模式及其相应的时域重建方法即延迟聚焦波束形成法,可以获取修正后的成像图形,并且由于采用了合成孔径的超声波发射-接收方式,利用其发射-接收具有动态聚焦的特点,有效地提高了图像重建的分辨率,且能够保持不同探测深度分辨率的一致性。
上述实施方式为本发明的优选案例,并不用来限制本发明的保护范围。
Claims (9)
1.一种基于射线理论相位修正的合成孔径波束形成方法,其特征在于,包括如下步骤:
步骤S1,采用发射阵元总数为L的线性阵列,单个所述发射阵元以预定的采样频率向声屏障内部发射超声波信号,并通过全部L个接收阵元获取与所述发射阵元相对应的采集信号,L个所述发射阵元依次循环发射,然后进入步骤S2;
步骤S2,获取声屏障的结构与声速,然后进入步骤S3;
步骤S3,成像区域空间离散化为二维方形网格集群M×N,并将网格的端点作为成像的离散点,在成像区域内确定所述声屏障及其他介质所在区域Ω,根据所述声屏障以及其他复杂介质在空间结构及声速c(m,n)|m<M,n<N,计算对应的空间慢度分布S(m,n)|m<M,n<N=1/c(m,n)|m<M,n<N,然后进入步骤S4;
步骤S4,将每个所述发射阵元与所述接收阵元均设定为点源,确定每个阵元在所述二维方形网格集群中的位置,分别以这些点源为起点,对程函方程进行离散化求解,计算其他网格端点处的超声初至时间t(m,n)|m<M,n<N,而后得到超声初至时间矩阵,然后进入步骤S5;
步骤S5,将步骤S4得到的所述超声初至时间矩阵替代常规波数形成所需的时间矩阵,进行图像重建处理。
2.根据权利要求1所述的基于射线理论相位修正的合成孔径波束形成方法,其特征在于:
其中,所述步骤S1中所述接收阵元的数量为L个,
在接收所述采集信号时,每个所述接收阵元均为有效阵元,从而形成L次发射接收循环。
3.根据权利要求1所述的一种基于射线理论相位修正的合成孔径波束形成方法,其特征在于:
其中,所述步骤S2中所述声屏障的结构从X射线CT图像或者核磁共振MRI图像中获取,
所述声屏障的声速通过超声测量获得。
4.根据权利要求1所述的一种基于射线理论相位修正的合成孔径波束形成方法,其特征在于:
其中,所述步骤S2中所述声屏障的结构和声速均通过地震层析反演技术估计获取。
5.根据权利要求1所述的一种基于射线理论相位修正的合成孔径波束形成方法,其特征在于:
其中,所述步骤S2中,在获取声速数据的前提下,通过如下步骤获取所述声屏障的结构:
步骤S2-1,利用成像的方法获取所述声屏障的上表面的形状及空间位置,假定所述上表面以下是所述声屏障的速度,进入步骤S2-2;
步骤S2-2,在S2-1的基础上,重复步骤S3至S5,可以获得所述声屏障的下表层的形状及空间位置。
6.根据权利要求1所述的基于射线理论相位修正的合成孔径波束形成方法,其特征在于:
其中,所述步骤S3中的网格是正方形网格,该网格的边长为h,S表示速度慢度空间分布。
7.根据权利要求1所述的基于射线理论相位修正的合成孔径波束形成方法,其特征在于:
其中,所述步骤S4中,离散化求解所述程函方程的方法为波前外推法、快速行进算法或快速扫掠算法,
当采用快速扫掠算法求解所述程函方程时具体包括如下子步骤:
步骤S4-1,对于二维情形,对整个所述成像区域的初至时间矩阵(M+2)×(N+2)预先设置大的初始时间,并将所述阵元所在空间点的初始时间设置为零,然后进行从左下角向右上角,即i=2:M+1,j=2:N+1的迭代循环;
步骤S4-2,为求解t(i,j),预先对比该点两侧时间最小值xmin=min[t(i-1,j),t(i+1,j)]和ymin=min[t(i,j-1),t(i,j+1)],假如t(i-1,j)=t(i+1,j)或t(i,j-1)=t(i,j+1),此处选择t(i-1,j)或者t(i,j-1);
步骤S4-3,在步骤S4-2的基础上,寻找t(i,j)对角点的局部最小值tmin,该最小值tmin为:
tmin=min[t(i-1,j-1),t(i+1,j+1),t(i-1,j+1),t(i+1,j-1)];
步骤S4-4,在步骤S4-2与步骤S4-3的基础上,当abs(xmin-ymin)>h/S(i,j)时,则t 1(i,j)=min(xmin-ymin)+h/S(i,j),当abs(xmin-ymin)≤h/S(i,j)时,则
步骤S4-5,在步骤S4-2与步骤S4-3的基础上,若则
步骤S4-6,在步骤S4-4与步骤S4-5的基础上,获得此时该点的最终走时t(i,j)=min[t1(i,j),t 2(i,j),t_last(i,j)],而后重复步骤S4-1,直至空间所有的初至时间被计算到,而后继续进行迭代直至t(i,j)-t_last(i,j)<Threshold,则结束迭代,
其中,i表示笛卡尔坐标系下水平方向离散点的空间坐标,j表示笛卡尔坐标系下垂直方向离散点的空间坐标,h表示正方形网格的边长,S表示速度慢度空间分布,t_last(i,j)表示上一次迭代计算的t(i,j)。
8.根据权利要求7所述的基于射线理论相位修正的合成孔径波束形成方法,其特征在于:
其中,所述迭代循环的方式还包括:
从右下角向左上角,即i=M+1:2,j=2:N+1,
从左上角向右下角,即i=2:M+1,j=N+1:2,以及
从右上角向左下角,即i=M+1:2,j=N+1:2。
9.根据权利要求1所述的基于射线理论相位修正的合成孔径波束形成方法,其特征在于:
其中,所述声屏障的厚度划分为至少两个网格宽度,即2h≤H,H表示所述声屏障的最小厚度。
Priority Applications (1)
Application Number | Priority Date | Filing Date | Title |
---|---|---|---|
CN201911060860.3A CN112764040B (zh) | 2019-11-01 | 2019-11-01 | 一种基于射线理论相位修正的合成孔径波束形成方法 |
Applications Claiming Priority (1)
Application Number | Priority Date | Filing Date | Title |
---|---|---|---|
CN201911060860.3A CN112764040B (zh) | 2019-11-01 | 2019-11-01 | 一种基于射线理论相位修正的合成孔径波束形成方法 |
Publications (2)
Publication Number | Publication Date |
---|---|
CN112764040A CN112764040A (zh) | 2021-05-07 |
CN112764040B true CN112764040B (zh) | 2022-06-14 |
Family
ID=75692149
Family Applications (1)
Application Number | Title | Priority Date | Filing Date |
---|---|---|---|
CN201911060860.3A Active CN112764040B (zh) | 2019-11-01 | 2019-11-01 | 一种基于射线理论相位修正的合成孔径波束形成方法 |
Country Status (1)
Country | Link |
---|---|
CN (1) | CN112764040B (zh) |
Families Citing this family (1)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
CN113552573B (zh) * | 2021-06-29 | 2022-07-29 | 复旦大学 | 一种基于超声环阵合成孔径接收的快速成像算法 |
Citations (9)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
WO2009088845A1 (en) * | 2007-12-31 | 2009-07-16 | Brigham And Women's Hospital, Inc. | System and method for accelerated focused ultrasound imaging |
CN101658434A (zh) * | 2009-09-10 | 2010-03-03 | 复旦大学 | 用于定征松质骨微结构的超声频谱偏移量参数成像方法 |
WO2015042386A1 (en) * | 2013-09-20 | 2015-03-26 | Westerngeco Llc | Eikonal solver for quasi p-waves in anisotropic media |
WO2016073976A1 (en) * | 2014-11-07 | 2016-05-12 | Tessonics Corporation | An ultrasonic adaptive beamforming method and its application for transcranial imaging |
EP3187113A1 (en) * | 2015-03-27 | 2017-07-05 | Tsinghua University | Method and device for ultrasonic imaging by synthetic focusing |
CN109031292A (zh) * | 2018-07-13 | 2018-12-18 | 山东大学 | 一种sar与介电反演结合的成像方法 |
CN109475755A (zh) * | 2016-07-25 | 2019-03-15 | 医视特有限公司 | 使用反射的超声自动聚焦 |
CN109875606A (zh) * | 2019-01-22 | 2019-06-14 | 北京航空航天大学 | 基于先验反射成像的超声ct声速成像的方法 |
CN110051387A (zh) * | 2019-04-11 | 2019-07-26 | 华中科技大学 | 一种基于射线理论的超声ct图像重建方法及系统 |
Family Cites Families (3)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
US6880943B2 (en) * | 2003-06-04 | 2005-04-19 | Itt Manufacturing Enterprises, Inc. | High resolution image formation from a multiple aperture imaging system |
US11397167B2 (en) * | 2016-11-22 | 2022-07-26 | The Board Of Trustees Of The Leland Stanford Junior University | Local speed of sound estimation method for medical ultrasound |
JP6820791B2 (ja) * | 2017-04-26 | 2021-01-27 | 株式会社日立製作所 | 超音波撮像装置、および、超音波画像生成プログラム |
-
2019
- 2019-11-01 CN CN201911060860.3A patent/CN112764040B/zh active Active
Patent Citations (10)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
WO2009088845A1 (en) * | 2007-12-31 | 2009-07-16 | Brigham And Women's Hospital, Inc. | System and method for accelerated focused ultrasound imaging |
CN101658434A (zh) * | 2009-09-10 | 2010-03-03 | 复旦大学 | 用于定征松质骨微结构的超声频谱偏移量参数成像方法 |
WO2015042386A1 (en) * | 2013-09-20 | 2015-03-26 | Westerngeco Llc | Eikonal solver for quasi p-waves in anisotropic media |
WO2016073976A1 (en) * | 2014-11-07 | 2016-05-12 | Tessonics Corporation | An ultrasonic adaptive beamforming method and its application for transcranial imaging |
CN107205720A (zh) * | 2014-11-07 | 2017-09-26 | 艾因蒂克公司 | 一种超声适应性波束形成方法及其对经颅成像的应用 |
EP3187113A1 (en) * | 2015-03-27 | 2017-07-05 | Tsinghua University | Method and device for ultrasonic imaging by synthetic focusing |
CN109475755A (zh) * | 2016-07-25 | 2019-03-15 | 医视特有限公司 | 使用反射的超声自动聚焦 |
CN109031292A (zh) * | 2018-07-13 | 2018-12-18 | 山东大学 | 一种sar与介电反演结合的成像方法 |
CN109875606A (zh) * | 2019-01-22 | 2019-06-14 | 北京航空航天大学 | 基于先验反射成像的超声ct声速成像的方法 |
CN110051387A (zh) * | 2019-04-11 | 2019-07-26 | 华中科技大学 | 一种基于射线理论的超声ct图像重建方法及系统 |
Non-Patent Citations (9)
Title |
---|
Automated 3D ultrasound image segmentation to aid breast cancer image interpretation;PengGu 等;《Ultrasonics》;20160228;第51-58页 * |
Reconstructive tomography and applications to ultrasonics;R.K. Mueller 等;《Proceedings of the IEEE》;19790430;第567-587页 * |
Relationships of Ultrasonic Backscatter With Bone Densities and Microstructure in Bovine Cancellous Bone;Chengcheng Liu 等;《IEEE Transactions on Ultrasonics, Ferroelectrics, and Frequency Control》;20180924;第2311-2321页 * |
合成孔径雷达图像重建的断层投影技术;孙斌;《电子测量与仪器学报》;20040930;第18卷(第3期);第85-88页 * |
基于多层声速模型的合成孔径超声皮质骨成像;李云清 等;《物理学报》;20190923;第68卷(第18期);第1-10页 * |
基于非凸正则化项的合成孔径雷达图像分割新算法;尚晓清 等;《光子学报》;20120930;第41卷(第9期);第1124-1129页 * |
声屏障对城区敏感点降噪效果的研究;田欣 等;《环境工程》;20140831;第1005-1008、1012页 * |
经颅脑肿瘤HIFU治疗技术的研究现状;韩珍珍 等;《生命科学仪器》;20150430;第7-13、38页 * |
被动合成孔径的相位补偿自校正算法;李博 等;《声学技术》;20140831;第33卷(第4期);第392-395页 * |
Also Published As
Publication number | Publication date |
---|---|
CN112764040A (zh) | 2021-05-07 |
Similar Documents
Publication | Publication Date | Title |
---|---|---|
US10959703B2 (en) | Methods, systems and computer program products for single track location shear wave elasticity imaging | |
US20110144495A1 (en) | Perfusion Imaging of a Volume in Medical Diagnostic Ultrasound | |
US20150025388A1 (en) | Ultrasound waveform tomography with wave-energy-based preconditioning | |
US10324063B2 (en) | Methods and systems for measuring properties with ultrasound | |
Peralta et al. | Coherent multi-transducer ultrasound imaging | |
EP3577491B1 (en) | System and method for speed and attenuation reconstruction in ultrasound imaging | |
JP2015042344A (ja) | 変位計測装置、並びに、超音波診断装置 | |
US20120163691A1 (en) | Time-domain estimator for image reconstruction | |
JP7340868B2 (ja) | 超音波方法及び装置 | |
CN105411624B (zh) | 一种超声三维流体成像与测速方法 | |
Filoux et al. | Characterization of the spatial resolution of different high-frequency imaging systems using a novel anechoic-sphere phantom | |
Skjelvareid et al. | Synthetic aperture focusing of outwardly directed cylindrical ultrasound scans | |
Rohling | 3D freehand ultrasound: reconstruction and spatial compounding | |
Jensen et al. | Estimation of velocity vectors in synthetic aperture ultrasound imaging | |
CN112764040B (zh) | 一种基于射线理论相位修正的合成孔径波束形成方法 | |
Li et al. | 3D full-waveform inversion in ultrasound computed tomography employing a ring-array | |
Andresen et al. | Synthetic aperture focusing for a single-element transducer undergoing helical motion | |
Pinton | Ultrasound imaging of the human body with three dimensional full-wave nonlinear acoustics. Part 1: Simulations methods | |
Luo et al. | Fundamental performance assessment of 2-D myocardial elastography in a phased-array configuration | |
Albulayli et al. | Phase-shift depth migration for plane-wave ultrasound imaging | |
Theis et al. | Seismic imaging of medical ultrasound data: Towards in vivo applications | |
Plotnick et al. | Circular synthetic aperture sonar imaging of simple objects illuminated by an evanescent wavefield | |
CN114601496A (zh) | 一种基于线性阵列的三维超声剪切波弹性成像方法 | |
Li et al. | Ray theory-based compounded plane wave ultrasound imaging for aberration corrected transcranial imaging: Phantom experiments and simulations | |
Bezek et al. | Global speed-of-sound prediction using transmission geometry |
Legal Events
Date | Code | Title | Description |
---|---|---|---|
PB01 | Publication | ||
PB01 | Publication | ||
SE01 | Entry into force of request for substantive examination | ||
SE01 | Entry into force of request for substantive examination | ||
GR01 | Patent grant | ||
GR01 | Patent grant |