CN114176639A - 用于介质的超声表征的方法和系统 - Google Patents

用于介质的超声表征的方法和系统 Download PDF

Info

Publication number
CN114176639A
CN114176639A CN202111075417.0A CN202111075417A CN114176639A CN 114176639 A CN114176639 A CN 114176639A CN 202111075417 A CN202111075417 A CN 202111075417A CN 114176639 A CN114176639 A CN 114176639A
Authority
CN
China
Prior art keywords
transducer
image
virtual transducer
medium
wavefront
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.)
Pending
Application number
CN202111075417.0A
Other languages
English (en)
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.)
Sound Imaging Co ltd
Centre National de la Recherche Scientifique CNRS
Ecole Superieure de Physique et Chimie Industrielles de Ville Paris
Original Assignee
Sound Imaging Co ltd
Centre National de la Recherche Scientifique CNRS
Ecole Superieure de Physique et Chimie Industrielles de Ville Paris
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 Sound Imaging Co ltd, Centre National de la Recherche Scientifique CNRS, Ecole Superieure de Physique et Chimie Industrielles de Ville Paris filed Critical Sound Imaging Co ltd
Publication of CN114176639A publication Critical patent/CN114176639A/zh
Pending legal-status Critical Current

Links

Images

Classifications

    • AHUMAN NECESSITIES
    • A61MEDICAL OR VETERINARY SCIENCE; HYGIENE
    • A61BDIAGNOSIS; SURGERY; IDENTIFICATION
    • A61B8/00Diagnosis using ultrasonic, sonic or infrasonic waves
    • A61B8/08Detecting organic movements or changes, e.g. tumours, cysts, swellings
    • A61B8/0833Detecting organic movements or changes, e.g. tumours, cysts, swellings involving detecting or locating foreign bodies or organic structures
    • A61B8/085Detecting organic movements or changes, e.g. tumours, cysts, swellings involving detecting or locating foreign bodies or organic structures for locating body or organic structures, e.g. tumours, calculi, blood vessels, nodules
    • GPHYSICS
    • G10MUSICAL INSTRUMENTS; ACOUSTICS
    • G10KSOUND-PRODUCING DEVICES; METHODS OR DEVICES FOR PROTECTING AGAINST, OR FOR DAMPING, NOISE OR OTHER ACOUSTIC WAVES IN GENERAL; ACOUSTICS NOT OTHERWISE PROVIDED FOR
    • G10K11/00Methods or devices for transmitting, conducting or directing sound in general; Methods or devices for protecting against, or for damping, noise or other acoustic waves in general
    • G10K11/18Methods or devices for transmitting, conducting or directing sound
    • G10K11/26Sound-focusing or directing, e.g. scanning
    • G10K11/28Sound-focusing or directing, e.g. scanning using reflection, e.g. parabolic reflectors
    • GPHYSICS
    • G01MEASURING; TESTING
    • G01SRADIO 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/00Systems using the reflection or reradiation of acoustic waves, e.g. sonar systems
    • G01S15/88Sonar systems specially adapted for specific applications
    • G01S15/89Sonar systems specially adapted for specific applications for mapping or imaging
    • G01S15/8906Short-range imaging systems; Acoustic microscope systems using pulse-echo techniques
    • G01S15/8909Short-range imaging systems; Acoustic microscope systems using pulse-echo techniques using a static transducer configuration
    • G01S15/8915Short-range imaging systems; Acoustic microscope systems using pulse-echo techniques using a static transducer configuration using a transducer array
    • AHUMAN NECESSITIES
    • A61MEDICAL OR VETERINARY SCIENCE; HYGIENE
    • A61BDIAGNOSIS; SURGERY; IDENTIFICATION
    • A61B8/00Diagnosis using ultrasonic, sonic or infrasonic waves
    • A61B8/13Tomography
    • AHUMAN NECESSITIES
    • A61MEDICAL OR VETERINARY SCIENCE; HYGIENE
    • A61BDIAGNOSIS; SURGERY; IDENTIFICATION
    • A61B8/00Diagnosis using ultrasonic, sonic or infrasonic waves
    • A61B8/44Constructional features of the ultrasonic, sonic or infrasonic diagnostic device
    • A61B8/4483Constructional features of the ultrasonic, sonic or infrasonic diagnostic device characterised by features of the ultrasound transducer
    • AHUMAN NECESSITIES
    • A61MEDICAL OR VETERINARY SCIENCE; HYGIENE
    • A61BDIAGNOSIS; SURGERY; IDENTIFICATION
    • A61B8/00Diagnosis using ultrasonic, sonic or infrasonic waves
    • A61B8/44Constructional features of the ultrasonic, sonic or infrasonic diagnostic device
    • A61B8/4483Constructional features of the ultrasonic, sonic or infrasonic diagnostic device characterised by features of the ultrasound transducer
    • A61B8/4488Constructional features of the ultrasonic, sonic or infrasonic diagnostic device characterised by features of the ultrasound transducer the transducer being a phased array
    • AHUMAN NECESSITIES
    • A61MEDICAL OR VETERINARY SCIENCE; HYGIENE
    • A61BDIAGNOSIS; SURGERY; IDENTIFICATION
    • A61B8/00Diagnosis using ultrasonic, sonic or infrasonic waves
    • A61B8/52Devices using data or image processing specially adapted for diagnosis using ultrasonic, sonic or infrasonic waves
    • AHUMAN NECESSITIES
    • A61MEDICAL OR VETERINARY SCIENCE; HYGIENE
    • A61BDIAGNOSIS; SURGERY; IDENTIFICATION
    • A61B8/00Diagnosis using ultrasonic, sonic or infrasonic waves
    • A61B8/52Devices using data or image processing specially adapted for diagnosis using ultrasonic, sonic or infrasonic waves
    • A61B8/5215Devices using data or image processing specially adapted for diagnosis using ultrasonic, sonic or infrasonic waves involving processing of medical diagnostic data
    • GPHYSICS
    • G01MEASURING; TESTING
    • G01NINVESTIGATING OR ANALYSING MATERIALS BY DETERMINING THEIR CHEMICAL OR PHYSICAL PROPERTIES
    • G01N29/00Investigating or analysing materials by the use of ultrasonic, sonic or infrasonic waves; Visualisation of the interior of objects by transmitting ultrasonic or sonic waves through the object
    • G01N29/02Analysing fluids
    • G01N29/024Analysing fluids by measuring propagation velocity or propagation time of acoustic waves
    • GPHYSICS
    • G01MEASURING; TESTING
    • G01NINVESTIGATING OR ANALYSING MATERIALS BY DETERMINING THEIR CHEMICAL OR PHYSICAL PROPERTIES
    • G01N29/00Investigating or analysing materials by the use of ultrasonic, sonic or infrasonic waves; Visualisation of the interior of objects by transmitting ultrasonic or sonic waves through the object
    • G01N29/02Analysing fluids
    • G01N29/032Analysing fluids by measuring attenuation of acoustic waves
    • GPHYSICS
    • G01MEASURING; TESTING
    • G01NINVESTIGATING OR ANALYSING MATERIALS BY DETERMINING THEIR CHEMICAL OR PHYSICAL PROPERTIES
    • G01N29/00Investigating or analysing materials by the use of ultrasonic, sonic or infrasonic waves; Visualisation of the interior of objects by transmitting ultrasonic or sonic waves through the object
    • G01N29/04Analysing solids
    • G01N29/06Visualisation of the interior, e.g. acoustic microscopy
    • G01N29/0654Imaging
    • GPHYSICS
    • G01MEASURING; TESTING
    • G01NINVESTIGATING OR ANALYSING MATERIALS BY DETERMINING THEIR CHEMICAL OR PHYSICAL PROPERTIES
    • G01N29/00Investigating or analysing materials by the use of ultrasonic, sonic or infrasonic waves; Visualisation of the interior of objects by transmitting ultrasonic or sonic waves through the object
    • G01N29/04Analysing solids
    • G01N29/07Analysing solids by measuring propagation velocity or propagation time of acoustic waves
    • GPHYSICS
    • G01MEASURING; TESTING
    • G01NINVESTIGATING OR ANALYSING MATERIALS BY DETERMINING THEIR CHEMICAL OR PHYSICAL PROPERTIES
    • G01N29/00Investigating or analysing materials by the use of ultrasonic, sonic or infrasonic waves; Visualisation of the interior of objects by transmitting ultrasonic or sonic waves through the object
    • G01N29/22Details, e.g. general constructional or apparatus details
    • G01N29/26Arrangements for orientation or scanning by relative movement of the head and the sensor
    • G01N29/262Arrangements for orientation or scanning by relative movement of the head and the sensor by electronic orientation or focusing, e.g. with phased arrays
    • GPHYSICS
    • G01MEASURING; TESTING
    • G01NINVESTIGATING OR ANALYSING MATERIALS BY DETERMINING THEIR CHEMICAL OR PHYSICAL PROPERTIES
    • G01N29/00Investigating or analysing materials by the use of ultrasonic, sonic or infrasonic waves; Visualisation of the interior of objects by transmitting ultrasonic or sonic waves through the object
    • G01N29/44Processing the detected response signal, e.g. electronic circuits specially adapted therefor
    • GPHYSICS
    • G01MEASURING; TESTING
    • G01NINVESTIGATING OR ANALYSING MATERIALS BY DETERMINING THEIR CHEMICAL OR PHYSICAL PROPERTIES
    • G01N29/00Investigating or analysing materials by the use of ultrasonic, sonic or infrasonic waves; Visualisation of the interior of objects by transmitting ultrasonic or sonic waves through the object
    • G01N29/44Processing the detected response signal, e.g. electronic circuits specially adapted therefor
    • G01N29/4463Signal correction, e.g. distance amplitude correction [DAC], distance gain size [DGS], noise filtering
    • GPHYSICS
    • G01MEASURING; TESTING
    • G01NINVESTIGATING OR ANALYSING MATERIALS BY DETERMINING THEIR CHEMICAL OR PHYSICAL PROPERTIES
    • G01N29/00Investigating or analysing materials by the use of ultrasonic, sonic or infrasonic waves; Visualisation of the interior of objects by transmitting ultrasonic or sonic waves through the object
    • G01N29/44Processing the detected response signal, e.g. electronic circuits specially adapted therefor
    • G01N29/4472Mathematical theories or simulation
    • GPHYSICS
    • G01MEASURING; TESTING
    • G01NINVESTIGATING OR ANALYSING MATERIALS BY DETERMINING THEIR CHEMICAL OR PHYSICAL PROPERTIES
    • G01N29/00Investigating or analysing materials by the use of ultrasonic, sonic or infrasonic waves; Visualisation of the interior of objects by transmitting ultrasonic or sonic waves through the object
    • G01N29/44Processing the detected response signal, e.g. electronic circuits specially adapted therefor
    • G01N29/46Processing the detected response signal, e.g. electronic circuits specially adapted therefor by spectral analysis, e.g. Fourier analysis or wavelet analysis
    • GPHYSICS
    • G01MEASURING; TESTING
    • G01SRADIO 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/00Systems using the reflection or reradiation of acoustic waves, e.g. sonar systems
    • G01S15/006Theoretical aspects
    • GPHYSICS
    • G01MEASURING; TESTING
    • G01SRADIO 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/00Details of systems according to groups G01S13/00, G01S15/00, G01S17/00
    • G01S7/52Details of systems according to groups G01S13/00, G01S15/00, G01S17/00 of systems according to group G01S15/00
    • G01S7/52017Details of systems according to groups G01S13/00, G01S15/00, G01S17/00 of systems according to group G01S15/00 particularly adapted to short-range imaging
    • G01S7/52019Details of transmitters
    • G01S7/5202Details of transmitters for pulse systems
    • G01S7/52022Details of transmitters for pulse systems using a sequence of pulses, at least one pulse manipulating the transmissivity or reflexivity of the medium
    • GPHYSICS
    • G01MEASURING; TESTING
    • G01SRADIO 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/00Details of systems according to groups G01S13/00, G01S15/00, G01S17/00
    • G01S7/52Details of systems according to groups G01S13/00, G01S15/00, G01S17/00 of systems according to group G01S15/00
    • G01S7/52017Details of systems according to groups G01S13/00, G01S15/00, G01S17/00 of systems according to group G01S15/00 particularly adapted to short-range imaging
    • G01S7/52023Details of receivers
    • G01S7/52025Details of receivers for pulse systems
    • G01S7/52026Extracting wanted echo signals
    • GPHYSICS
    • G01MEASURING; TESTING
    • G01SRADIO 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/00Details of systems according to groups G01S13/00, G01S15/00, G01S17/00
    • G01S7/52Details of systems according to groups G01S13/00, G01S15/00, G01S17/00 of systems according to group G01S15/00
    • G01S7/52017Details of systems according to groups G01S13/00, G01S15/00, G01S17/00 of systems according to group G01S15/00 particularly adapted to short-range imaging
    • G01S7/52023Details of receivers
    • G01S7/52036Details of receivers using analysis of echo signal for target characterisation
    • GPHYSICS
    • G01MEASURING; TESTING
    • G01SRADIO 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/00Details of systems according to groups G01S13/00, G01S15/00, G01S17/00
    • G01S7/52Details of systems according to groups G01S13/00, G01S15/00, G01S17/00 of systems according to group G01S15/00
    • G01S7/52017Details of systems according to groups G01S13/00, G01S15/00, G01S17/00 of systems according to group G01S15/00 particularly adapted to short-range imaging
    • G01S7/52023Details of receivers
    • G01S7/52036Details of receivers using analysis of echo signal for target characterisation
    • G01S7/52042Details of receivers using analysis of echo signal for target characterisation determining elastic properties of the propagation medium or of the reflective target
    • GPHYSICS
    • G01MEASURING; TESTING
    • G01SRADIO 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/00Details of systems according to groups G01S13/00, G01S15/00, G01S17/00
    • G01S7/52Details of systems according to groups G01S13/00, G01S15/00, G01S17/00 of systems according to group G01S15/00
    • G01S7/52017Details of systems according to groups G01S13/00, G01S15/00, G01S17/00 of systems according to group G01S15/00 particularly adapted to short-range imaging
    • G01S7/52046Techniques for image enhancement involving transmitter or receiver
    • GPHYSICS
    • G01MEASURING; TESTING
    • G01SRADIO 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/00Details of systems according to groups G01S13/00, G01S15/00, G01S17/00
    • G01S7/52Details of systems according to groups G01S13/00, G01S15/00, G01S17/00 of systems according to group G01S15/00
    • G01S7/52017Details of systems according to groups G01S13/00, G01S15/00, G01S17/00 of systems according to group G01S15/00 particularly adapted to short-range imaging
    • G01S7/52046Techniques for image enhancement involving transmitter or receiver
    • G01S7/52049Techniques for image enhancement involving transmitter or receiver using correction of medium-induced phase aberration
    • GPHYSICS
    • G10MUSICAL INSTRUMENTS; ACOUSTICS
    • G10KSOUND-PRODUCING DEVICES; METHODS OR DEVICES FOR PROTECTING AGAINST, OR FOR DAMPING, NOISE OR OTHER ACOUSTIC WAVES IN GENERAL; ACOUSTICS NOT OTHERWISE PROVIDED FOR
    • G10K11/00Methods or devices for transmitting, conducting or directing sound in general; Methods or devices for protecting against, or for damping, noise or other acoustic waves in general
    • G10K11/18Methods or devices for transmitting, conducting or directing sound
    • G10K11/26Sound-focusing or directing, e.g. scanning
    • G10K11/34Sound-focusing or directing, e.g. scanning using electrical steering of transducer arrays, e.g. beam steering
    • G10K11/341Circuits therefor
    • G10K11/346Circuits therefor using phase variation
    • GPHYSICS
    • G01MEASURING; TESTING
    • G01NINVESTIGATING OR ANALYSING MATERIALS BY DETERMINING THEIR CHEMICAL OR PHYSICAL PROPERTIES
    • G01N2291/00Indexing codes associated with group G01N29/00
    • G01N2291/01Indexing codes associated with the measuring variable
    • G01N2291/011Velocity or travel time
    • GPHYSICS
    • G01MEASURING; TESTING
    • G01NINVESTIGATING OR ANALYSING MATERIALS BY DETERMINING THEIR CHEMICAL OR PHYSICAL PROPERTIES
    • G01N2291/00Indexing codes associated with group G01N29/00
    • G01N2291/02Indexing codes associated with the analysed material
    • G01N2291/024Mixtures
    • G01N2291/02475Tissue characterisation
    • GPHYSICS
    • G01MEASURING; TESTING
    • G01NINVESTIGATING OR ANALYSING MATERIALS BY DETERMINING THEIR CHEMICAL OR PHYSICAL PROPERTIES
    • G01N2291/00Indexing codes associated with group G01N29/00
    • G01N2291/10Number of transducers
    • G01N2291/106Number of transducers one or more transducer arrays

Abstract

用于介质的超声表征的方法,包括:生成一系列入射超声波的步骤;生成在作为输入的发射基(i)与作为输出的接收基(u)之间定义的实验反射矩阵Rui(t)的步骤;确定在基于作为所述实验反射矩阵的输入的聚焦计算的输入虚拟换能器(TVin)与基于来自所述实验反射矩阵的输出的聚焦计算的输出虚拟换能器(TVout)之间的所述介质的聚焦反射矩阵RFoc(rin,rout,δt)的步骤,在相对于所述输入虚拟换能器(TVin)的响应的时刻移位额外延迟δt的时刻获得所述输出虚拟换能器(TVout)的响应。

Description

用于介质的超声表征的方法和系统
技术领域
本说明书涉及用于介质的超声表征的方法和系统,且具体来说适用于医学成像或无损检测,且更一般来说适用于可以使用超声成像的所有领域。
背景技术
在声学成像的领域中,我们试图通过使用超声波主动探测来表征完全或部分未知的环境。这是用于医学成像的超声机器的原理。
声学成像系统的分辨率可以被定义为用于辨别对象的小细节的能力。原则上,声学成像系统受衍射限制,并且理论分辨率由λ/2(其中λ是介质中声音的波长)或由检测器的有限角孔径给出。然而,实际上,当传播介质为异质时,分辨率通常会因声速的变化而降低。
实际上,在声学成像中的大部分时间,介质被视为均质的且具有恒定声速c0。然而,均质环境的假设并不总是适用。例如,在肝脏超声的情况下,探头放置在患者的肋骨之间。声波在到达目标器官之前穿过脂肪和肌肉层。软组织各自具有不同的机械特性。因此,声速远非均质的,并且它可以在例如脂肪组织的1450m/s与肝脏的1600m/s之间变化。声速的变化会导致波的不同相移,这取决于它们传播的区域。这会导致声波波前的像差,从而导致所得超声图像失真并因此降低其分辨率和对比度。这些像差可能导致例如在医学检查期间无法重建可靠的图像,从而影响结果。
如图1A到1C中所说明,常规超声方法使用压电换能器11阵列10,压电换能器可以独立地发射和/或接收超声脉冲。每个换能器的位置由向量u识别。当此阵列面向人们希望研究的介质放置时,可以以各种方式对介质进行声穿透和成像。
生成要研究的介质的超声图像的第一方式是从阵列的换能器中的一个发射超声脉冲,所述阵列的位置由向量uin(图1A,左图)识别。这产生1D(或2D)换能器阵列的发散圆柱形(或球形)入射波。此波由介质20的散射体21反射并且反向散射场由每个换能器11记录为时间的函数(图1A、右图)。通过每个换能器连续地用作源来重复此操作,测量每个换能器之间的脉冲响应的集合R(uout,uin,t),其中向量uout表示检测器的位置。这些响应形成以换能器为基础表示的反射矩阵Ruu(t)。此测量的优点在于,此矩阵包含有关所分析介质的所有信息,然后可以对其应用一组矩阵运算,以便例如对介质进行成像。另一方面,这种采集假设介质在测量的持续时间内保持固定,这在体内使用的情况下可能非常困难。另外,单个压电元件发射的能量很低,这会导致较差的信噪比。
已知用于生成要分析的介质的图像的其它方法,其中使用波束成形技术执行聚焦发射。如图1B的左图中所示,这些方法包括基于均质速度模型向换能器11施加一组适当的延迟,以便校正波的传播时间,以便所有脉冲一起到达位置rin处的目标焦点。所采用的假设声速将表示为c0。由于衍射的物理限制,所发射的超声波集中在由超声探头的孔口定界的区域中。为了构建超声图像,还在接收时执行聚焦步骤。接着处理由阵列10中的元件11捕获的回波集合,以模拟在接收时的透镜效果,如在图1B的右图中所示。由换能器接收的信号进行时移以使其回到相位中。这些延迟与发射时施加的延迟相同。在发射阶段,所有信号在位置点rin处发生干扰。在接收时,来自此相同点rout=rin的信号通过对弹道时间t=(‖uout-rin‖+‖uin-rin‖)/C0的信号求和发生电子干扰。此总和提供在接收时的聚焦的最终结果。图1B中所说明的方法(称为在发射和接收时的共焦双聚焦方法)可以直接对介质的反射率成像,其具有受衍射限制的横向分辨率、仅受初始脉冲的持续时间限制的极佳轴向分辨率,以及极佳的对比度。然而,此方法是耗时的,因为它需要在每一行图像上在介质的每个点或至少在给定深度对发射进行物理聚焦。
最近开发的另一成像技术由通过用一系列平面波使声穿透介质来生成介质图像组成。图1C说明此所谓的平面波超声的原理,例如在G.Montaldo等人的文章“用于极高帧速率超声波检查和瞬态弹性成像的相干平面波复合(Coherent plane-wave compoundingfor very high frame rate ultrasonography and transient elastography)”(《IEEE超声学、铁电体技术与频率控制汇刊》第56卷,第489到506页,2009年)中描述。在发射时对每个信号施加延迟(图1C、左图)以形成相对于换能器10阵列以角度θin倾斜的波前。在接收(图1C、右图)时,由介质反向散射的场R(uoutin,t)通过用于一系列入射平面波的所有位置传感器uout测量,其中改变入射角θin。这些响应的集合形成在作为输入的空间傅里叶基(或平面波基(plane wave basis))与作为输出的换能器基之间定义的反射矩阵R(t)。一旦记录此矩阵,就在相干地求和之前对信号进行时移,以在发射和接收时对每个位置点rin的数据进行数字聚焦。因此,与标准超声(聚焦发射)相比,形成超声图像所需的数据捕获数目有利地减少,并且这对于超声图像的相同对比度和分辨率水平来说是正确的。
图2说明环境像差对常规超声成像方法(图1A到1C)的影响。当介质c(r)中的声速并不对应于具有恒定声速c0的均质介质的假设时,出现这些像差。最初基于此假设确定并在发射和接收时施加到阵列的每个换能器的延迟然后对于评估介质的图像来说不是最佳的。在图2中,像差层22引发入射波前的失真。在发射或激发时,在步骤25,所使用的延迟规则不允许声能集中在由衍射极限界定的区域中,所述区域通常称为焦斑。在接收时,在步骤26中,所使用的延迟规则不允许正确地选择源自介质焦点的超声信号,并且混合源自同样异常的焦斑的信号。这会导致图像构建过程中的双重像差,从而大大降低图像的分辨率。然后可以重新计算新的延迟规则以补偿像差层的影响,例如通过向波束成形中通常使用的延迟添加额外的延迟规则。
然而,这些像差校正并不完全校正这些像差或分辨率降低。需要更好地估计介质中的聚焦质量。
文档“脉冲回波测量中的冯西泰-塞尼克定理(The van Cittert-Zerniketheorem in pulse echo measurements)”(Raoul Mallart和Mathias Fink,美国声学学会杂志90(5),1991年11月)研究在简单的散射条件下由随机介质反射的场的统计特性。具体来说,已经证明,对于聚焦入射波,来自远场的反射场的空间协方差与发射孔径函数的傅里叶变换成比例。换句话说,此定理解释对远场中的反射场的统计特性的研究可以确定入射波在介质中的聚焦质量。
然而,这种方法仅提供超声图像的分辨率的整体平均估计,因为它需要在大量无序实施方案上,即在入射波的大量焦点上对反射场的相关性进行统计平均化。所述方法不允许对图像的每个点的聚焦质量进行精确的局部评估。此外,这种方法仅在简单的散射条件下有效。
因此,有必要提出一种克服每个上述缺点的方法。
发明内容
根据第一方面,本说明书涉及一种用于介质的超声表征以确定积分声速值的方法,所述方法包括:
-通过换能器(11)阵列(10)在所述介质的区域中生成一系列入射超声波(USin)的步骤,所述系列的入射超声波是发射基(i);以及
-生成在作为输入的发射基(i)与作为输出的接收基(u)之间定义的实验反射矩阵Rui(t)的步骤;
-确定聚焦反射矩阵RFoc(rin,rout,δt)的步骤,所述聚焦反射矩阵包括在空间位置rin的输入虚拟换能器(TVin)与空间位置rout的输出虚拟换能器(TVout)之间的介质的响应,在相对于输入虚拟换能器(TVin)的响应的时刻移位额外延迟δτ的时刻获得输出虚拟换能器(TVout)的响应,
-针对输入虚拟换能器(TVin)以及针对额外延迟间隔确定波前图像的步骤,所述波前图像根据介质中的声速c0确定,并且所述波前图像基于以下项确定:
-聚焦反射矩阵RFoc(rin,rout,δt),以及
-以下类型的弹道传播关系
δt(Δrout)=-sign(Δzout)·|Δrout|/c0,其可以从聚焦反射矩阵提取值以便构建波前图像,并且其中:
δt是额外延迟,
|Δrout|是输入虚拟换能器(TVin)与输出虚拟换能器(TVout)之间的向量的模数,其中Δrout=rout-rin
Δzout是沿着空间位置向量Δrout的深度轴Z的分量。
-确定波前图像中的焦斑的中心的深度位置Δz(0)(rin)的步骤,以及
-基于以下公式计算积分声速c(1)(rin)的步骤:
Figure BDA0003261953800000041
其中zin是沿着输入虚拟换能器(TVin)的空间位置向量rin的深度轴Z的分量。
通过这些布置,所述方法有利地使得可以在任何点和任何方向上并且以相对于超声波在介质中的弹道传播时间的任何时间滞后来局部地探测介质。
此方法随后可以确定在弹道参考系统内聚焦在波上的图像,此图像称为“波前图像”并跟随波向焦点的传播,从而可以估计聚焦质量。随后可以推导出此焦点的声速值(积分)。
积分声速的局部确定构成超声成像中的新对比度,与反射率的图像(“经典”超声图像)互补。在任何焦点处的声速的此确定随后还可以通过计算提高超声图像的质量,因此无需迭代新的发射和/或采集,从而例如在对人进行医学检查的情况下呈现工作人员、患者和设备的时间和可用性方面的经济优势和显着节省。还可以几乎实时地表征传播介质,尤其在体内测量期间这是一个重要优势。
在根据本公开的方法的各种实施例中,可以任选地依赖以下布置中的一个或多个。
根据一个变体,通过在波前图像中搜索具有最大值的点来确定焦斑的中心。
根据一个变体,仅对深度轴Z执行波前图像的确定。
根据一个变体:
-在确定波前图像的步骤与确定焦斑的深度位置Δz(0)(rin)的步骤之间,执行改进波前图像的步骤,其中执行对应于给定相干区域的一组波前图像的线性组合,在不同空间位置rin的选定输入虚拟换能器TVin与空间位置rout的输出虚拟换能器TVout之间获得所述集合的每个波前图像,使得rout=Δrout+rin,其中Δrout预定义且针对所述集合的所有波前图像相同,并且选定输入虚拟换能器彼此靠近,以便获得与参考输入虚拟换能器(TVin,ref)相关联的改进的波前图像,此参考输入虚拟换能器TVin,ref表示所使用的波前图像集合的输入虚拟换能器并且与所选择的相干区域ZC相关联,以及
-在确定深度位置Δz(0)(rin)的步骤中,代替波前图像使用改进的波前图像,焦斑的中心的深度位置与参考输入虚拟换能器TVin,ref的空间位置相关,并且焦斑的中心的此深度位置可以估计在参考输入虚拟换能器TVin,ref的空间位置处的积分声速c(1)(rin,ref)。
根据一个变体,通过计算波前图像的集合的奇异值分解(SVD)以获得与奇异值分解的最大绝对值的奇异值相关联的奇异向量(W1)来确定线性组合,此奇异向量(W1)随后是对应于所述参考输入虚拟换能器(TVin,ref)且针对相同额外延迟δt的改进的波前图像。
根据一个变体,通过计算积分声速并且对于波前图像的集合的线性组合,使用基本上覆盖介质中的整个所关注区域的与选定输入虚拟换能器(TVin)相对应的一组波前图像来确定介质的最佳声速。
根据一个变体:
-代替先前使用的声速或代替在第一步骤中使用的声速c0,使用在前一迭代中确定的积分声速c(n)迭代确定聚焦反射矩阵RFoc(rin,rout,δt)、确定波前图像、确定波前图像中的焦斑的中心的深度位置Δz(n)的步骤,以及
-在计算积分声速的步骤期间,使用以下递归公式,
Figure BDA0003261953800000051
以及
其中在对应于输入虚拟换能器的介质点处的积分声速值是在所述方法的步骤n中计算出的积分声速c(n)(rin),此步骤n通过预定迭代次数或通过用于两个连续积分声速值之间的差或两者的组合的停止阈值确定。
根据一个变体:
-输入虚拟换能器和输出虚拟换能器的角色反转,以确定相对于输出虚拟换能器的积分声速c(1)(rout),以及
-将参考输入虚拟换能器的积分声速c(1)(rin)与参考输出虚拟换能器的积分声速c(1)(rout)组合以获得改进的积分声速。
根据一个变体,所述方法进一步包括通过确定介质中的多个点的积分声速来确定积分声速图像的步骤,每个点对应于空间位置rin的输入虚拟换能器(TVin)或空间位置rin,ref的参考输入虚拟换能器(TVin,ref)。
根据一个变体,确定积分声速图像,其中从积分声速图像的值计算此声速图像的每个点处的值。
根据一个变体,在确定聚焦反射矩阵的步骤中:
输入虚拟换能器(TVin)的响应的计算对应于基于实验反射矩阵Rui(t)在输入处的聚焦过程,所述实验反射矩阵使用发射基与输入虚拟换能器(TVin)之间的波的向外飞行时间以在空间位置rin处产生输入焦斑,
输出虚拟换能器(TVout)的响应的计算对应于基于实验反射矩阵Rui(t)在输出处的聚焦过程,所述实验反射矩阵使用输出虚拟换能器(TVout)与接收基u的换能器之间的波的返回飞行时间以在空间位置rout处产生输出焦斑,
额外延迟δt是在聚焦过程期间添加到向外和返回飞行时间的时间滞后。
根据一个变体,通过以下公式计算聚焦反射矩阵:
Figure BDA0003261953800000061
其中
Nin是发射基(i)的元素的数目,
Nout是在输出处的接收基(u)的元素的数目,
Rui(t)是实验反射矩阵,其中Rui(uout,iin,τ(rin,rout,uout,iin,δt))是在发射基中的索引iin的发射之后且在时间τ处由空间位置uout的换能器记录的实验反射矩阵Rui(t)的元素,
τ是时间,所述时间是在发射基(i)的换能器与空间位置rin的输入虚拟换能器(TVin)之间的超声波的向外飞行时间τin,以及在空间位置rout的输出换能器(TVout)与接收基u的换能器之间的超声波的返回飞行时间τout,以及额外延迟δt的总和,如通过以下公式解释:
τ(rin,rout,uout,iin,δt)=τin(rin,iin)+τout(rout,uout)+δt
根据第二方面,本说明书涉及一种用于介质的超声表征以确定在介质中的任一点处的积分声速值的系统,并且所述系统被配置用于实施如上文所描述的用于超声表征的方法。根据第二方面的用于超声表征的系统包括:
-换能器阵列,其适合于在介质的区域中生成一系列入射超声波并且记录通过所述区域反向散射的超声波作为时间的函数;以及
-计算单元,其与换能器阵列相关联并且适合于实施根据第一方面的方法。
附图说明
通过阅读参考附图出于说明性目的以非限制性方式呈现的以下详细描述,上文呈现的技术的其它特征和优点将变得显而易见,其中:
图1A到1C(已经描述)说明用于超声成像和量化的已知发射/接收机构;
图2(已经描述)说明根据现有技术的超声成像中像差的影响;
图3说明根据本说明书的用于实施用于超声表征的方法的用于超声表征的系统的实例;
图4说明根据本说明书的在用于超声表征的方法中使用的定义;
图5示出超声图像,其中,选择与回声元件相关联的位置且针对若干额外延迟获得围绕此位置的传播图像;
图6示出超声图像,其中,选择与具有相当的反射率的一组子分辨率散射体相关联的位置并且针对若干额外延迟获得围绕此位置的传播图像;
图7示出超声图像,其中,选择与具有相当的反射率的子分辨率散射体相关联的一组位置并且针对若干额外延迟获得从与每个选定位置相关联的传播图像的组合产生的相干波传播图像;
图8A示出与对应于具有相当的反射率的子分辨率散射体的位置相关联的传播图像以及相干波传播图像的中心点的强度的时间变化曲线;
图8B示出图8A的曲线的频谱;
图9A示出与图5的位置相同的位置相关联的波前的图像振幅;
图9B示出与图9A中使用的波前图像相同的波前图像的实部;
图10示出与为图5选择且针对三个假设声速获得的位置相同的位置相关联的若干波前图像的振幅,以及这些波前图像纵轴Δz上的强度曲线;
图11示出超声图像以及积分声速的对应图像;
图12示出根据本公开的在不具有像差校正(图像A)、具有常见横向像差校正(图像B)以及具有使用波前图像中的测量的轴向像差校正的情况下获得的超声图像;
图13示出包括具有在各个方向上倾斜的肌纤维的若干区域的超声图像(A),以及对应于介质的这些各种区域的波前图像(B、C、D、E);
图14示出用于计算肌纤维在介质中的优选倾斜方向的角度的曲线;
图15示出叠加在具有肌纤维的介质的超声图像上的所确定优选方向的矩阵;
图16示出超声图像(A)、此超声图像(B)围绕特定点的区域的放大、本地时间信号(C)的实部和此特定点的所估计振幅(D),以及介质的此本地时间信号的频谱分析;
图17示出超声图像(A)以及作为此图的对应超声图像的深度Z的函数的平均频谱的估计值;
图18示出超声图像(A)以及针对此图的对应超声图像确定的频谱相关图像。
在将参考附图描述的各种实施例中,除非另有说明,否则相似或相同的元件具有相同的附图标记。
具体实施方式
在以下详细描述中,仅详细描述某些实施例以确保描述的清楚,但这些实例并不旨在限制从此描述中出现的原理的一般范围。
本说明书中描述的各种实施例和方面可以以多种方式组合或简化。具体来说,除非另有说明,否则各可以重复、反转和/或并行执行种方法的步骤。
本说明书涉及用于介质的超声表征的方法和系统,并且具体来说适用于活体或非活体组织的医学成像。介质是例如异质介质,人们试图表征所述异质介质以例如标识和/或表征异质性。任选地,这些方法和系统可以应用于无损探测,例如,金属部件等。因此,这些表征技术在随后保存的介质中是非侵入性的。
图3说明根据本说明书的用于实施用于介质(例如,异质介质20)的超声表征的方法的超声表征系统40的实例。系统40包括换能器11的至少一个阵列10,例如线性或二维或矩阵阵列;换能器是例如压电超声换能器,所述换能器可以采用与介质20直接或间接接触的刚性杆的传统形式。换能器阵列是例如探测装置41(通常称为探头)的一部分;换能器阵列连接到计算单元42,所述计算单元自身可以连接到显示装置43或与显示装置43相关联;计算单元从换能器11中的每一个发射电信号和/或将电信号记录到换能器11中的每一个。超声换能器然后将这些电信号转换成超声波,且反之亦然。探测装置41、计算单元42和显示装置43之间的“连接”或“链接”应理解成表示任何类型的有线连接、电连接或光学连接,或使用例如WiFiTM、BluetoothTM或其它协议的任何协议的任何类型的无线连接。这些连接或链接是单向的或双向的。
计算单元42被配置成实施计算或处理步骤,具体来说实施根据本说明书的方法的步骤。按照惯例,通过取第一轴X和垂直于第一轴的第二轴Z来定义介质20的空间参考系统。为简单起见,第一轴X对应于换能器11针对线性阵列对准的横向方向,并且第二轴Z对应于介质20相对于换能器11的此阵列10的深度。此定义可以适合于上下文,且因此例如在二维阵列10的情况下延伸到三轴空间参考系统。
在图3中,如在本说明书的其余部分中,涉及用于发射和接收的换能器阵列,应了解,在更一般情况下,可以同时使用换能器的若干阵列。类似地,阵列可以由相同类型或不同种类的一个(1)到N个换能器组成。换能器可以是发射器和接收器两者,或者仅用于一些的发射器而仅用于其它的接收器。
换能器阵列例如充当发射器和接收器两者,或由换能器的若干子阵列组成,一些专用于超声波的发射,其它专用于接收。术语“换能器阵列”应理解成表示至少一个换能器、一系列对准或非对准的换能器,或换能器的矩阵。
在本说明书中,当涉及计算或处理步骤时,特别是用于实施方法的步骤时,应理解每个计算或处理步骤可以由软件、硬件、固件、微码或这些技术或相关技术的任何合适组合来实施。当使用软件时,每个计算或处理步骤可以通过例如可以被解释或执行的计算机程序指令或代码来实施。这些指令可以存储在或传输到计算机(或计算单元)可读的存储介质和/或由计算机(或计算单元)执行以实施这些计算或处理步骤。
介质中的点通过聚焦反射矩阵的分析
本说明书描述用于介质的超声表征的方法和系统。在实际情况下,假设介质是异质的。这些方法和系统基于图4中所示的定义:
我们在介质中定义:
-在介质的空间参考系统中的空间位置rin的第一点P1,
-在介质的空间参考系统中的空间位置rout的第二点P2。
这些空间位置rin和rout以粗体指示,以表示这些元件是位置向量,即在介质的空间参考系统中获取的向量(X,Z)。点位置的其它表示和定义对于超声专业的任何专家都是可能的并且是可访问的。
在超声频率下将这两个点P1和P2选择为彼此相距较短距离,这表示彼此几毫米,以及例如二十(20)毫米或更小。
如图4中表示,通过系统40的计算单元42实施的用于超声表征的方法包括:
-通过换能器11阵列10在所述介质的区域中生成一系列入射超声波USin的步骤,所述系列的入射超声波是发射基i;以及
-生成在作为输入的发射基i与作为输出的接收基u之间定义的实验反射矩阵Rui(t)的步骤;
-确定聚焦反射矩阵RFoc(rin,rout,δt)的步骤,所述聚焦反射矩阵包括在空间位置rin的输入虚拟换能器TVin与空间位置rout的输出虚拟换能器TVout之间的介质的响应,在相对于输入虚拟换能器TVin的响应的时刻移位额外延迟δt的时刻获得输出虚拟换能器TVout的响应。
聚焦反射矩阵RFoc(rin,rout,δt)的响应对应于在介质中的任一点处计算的声压场。
作为输入的发射基i是例如各自由阵列10的换能器11中的一个生成的波的基,或相对于轴线X具有角倾斜θ的平面波的基,如上文在说明书中针对图1A到1C所描述。
接收基u例如是换能器11的基。任选地,可以在接收时使用另一接收基。
因此,生成超声波的步骤预期在发射基i与接收基u之间。因此,针对任何类型的聚焦或非聚焦的超声波(例如平面波)定义生成超声的此步骤。
在生成矩阵的步骤中,在作为输入的发射基i与作为输出的接收基u之间定义实验反射矩阵Rui(t)。此矩阵包含在时间t处由空间坐标uout的每个换能器11针对每次发射iin测量的介质的时间响应集合。应理解,以索引“in”命名的元素是指发射(即输入),并且以索引“out”命名的元素是指接收(即输出)。也可以将此实验矩阵记录和/或存储在例如计算单元的存储器中,或在可移动或不可移动的任何其它介质上,从而实现永久或临时存储。
更精确地说,在确定聚焦反射矩阵RFoc(rin,rout,δt)的步骤中,我们应用:
-基于实验反射矩阵Rui(t)在输入处的聚焦过程,所述实验反射矩阵使用发射基(i)与输入虚拟换能器TVin之间的波的向外飞行时间并且围绕空间位置rin的第一点P1产生所谓的输入焦斑,所述输入焦斑对应于输入虚拟换能器TVin
-基于实验反射矩阵Rui(t)在输出处的聚焦过程,所述实验反射矩阵使用输出虚拟换能器(TVout)与接收基(u)的换能器之间的波的返回飞行时间并且围绕空间位置rout的第二点P2产生所谓的输出焦斑,所述输出焦斑对应于输出虚拟换能器TVout
-额外延迟δt,所述额外延迟是在聚焦过程期间添加到向外和返回飞行时间的时间滞后。
在输入和输出处的这些聚焦过程实际上形成输入-输出处的聚焦过程,这在本说明书的其余部分中称为聚焦过程。
换句话说,在用于超声表征的此方法中,输入虚拟换能器TVin对应于位于介质中的空间位置rin处的超声“虚拟源”,并且输出虚拟换能器TVout对应于位于空间位置rout处的超声“虚拟传感器”。此虚拟源和此虚拟传感器通过其空间位置的差异Δr=rout-rin在空间上分离。它们还可以通过额外延迟δt在时间上分离,所述额外延迟是独立于空间距离|Δr|的任意且可调整延迟。因此,所述方法能够在空间上和/或在时间上探测围绕点P1和/或点P2的介质,这可以获得关于波传播的在这两个维度(空间和时间)上的新信息。
例如,通过输入和输出处的所述聚焦过程计算在输入虚拟换能器TVin与输出虚拟换能器TVout之间的介质的聚焦反射矩阵RFoc(rin,rout,δt)是可以通过以下简化公式表示的改进的波束成形方法:
Figure BDA0003261953800000111
其中
Nin是发射基i的元素的数目,
Nout是在输出处的接收基u的元素的数目,
Rui(t)是实验反射矩阵,其中Rui(uout,iin,τ(rin,rout,uout,iin,δt))是在时间τ处在发射iin之后由换能器uout记录的实验反射矩阵Rui(t)的元素。
时间τ是在发射基i的换能器与空间位置rin(第一点P1)的输入虚拟换能器TVin之间的超声波的向外飞行时间τin、在空间位置rout(第二点P2)的输出虚拟换能器TVout与接收基u的换能器之间的超声波的返回飞行时间τout,以及额外延迟δt的总和,如通过以下公式解释:
τ(rin,rout,uout,iin,δt)=τin(rin,iin)+τout(rout,uout)+δt
(方程式2)
从声速模型计算飞行时间τin和τout。最简单假设由假设具有恒定声速c0的均质介质组成。在这种情况下,基于探头的换能器与虚拟换能器之间的距离直接获得飞行时间。
发射基Nin的元素数目例如大于或等于一(1),并且有利地大于或等于二(2)。接收基Nout的元素数目例如大于或等于二(2)。
因此,此改进的波束成形公式是记录在实验反射矩阵Rui中的时间响应的两倍总和:根据发射基i的第一总和表示发射时的聚焦,以及根据接收基u的第二总和与接收时的聚焦相关,此计算针对(空间位置rin、rout的)两个点P1和P2的空间坐标执行。因此,此改进的波束成形公式的结果是这两个空间坐标(rin、rout)的时间信号,但也是输入与输出之间的额外延迟δt的函数,此额外延迟进行任意调整。
这种波束成形公式还可以通过通常称为接收和/或传输切趾法的输入和输出加权项来补充。因此,本领域技术人员可以用这些权重补充波束形成公式的其余部分。
所记录的实验反射矩阵Rui(t)可以是“实数”矩阵,即由时域中的实数系数组成,由换能器中的每一个记录的电信号是实数。或者,在用于同相和正交波束成形(“IQ波束成形”)的解调的情况下,此矩阵可以是“复数”矩阵,即由复数值组成。
我们因此获得包含时间信号的聚焦反射矩阵RFoc(rin,rout,δt)。此聚焦反射矩阵在线性探头的情况下具有五(5)个维度;用于空间位置rin和rout的两个空间,以及额外延迟δt,这与现有技术的聚焦反射矩阵非常不同并且信息丰富得多。
在此分析中,由于额外延迟δt,不在相同时刻定义输入虚拟换能器TVin和输出虚拟换能器TVout,这可以虚拟地突出输入虚拟换能器TVin的第一点P1与输出虚拟换能器TVout的第二点P2之间的超声波的传播。此额外延迟δt可以是正的或负的,这可以在超声波在介质中的路径的参考时刻之前和之后分别探测超声波在第二点P2处的聚焦。
此参考时刻称为弹道时间tb。此弹道时间是发射基i的换能器与输入虚拟换能器TVin之间、随后输出虚拟换能器TVout与接收基u的换能器之间的超声波的往返时间。
此弹道时间tb通过以下公式定义:
tb=(‖uout-rOut‖+‖uin-rin‖)/c0
(方程式3)
其中:
c0是介质的假设声速(超声波的传播速度)。
由于这些布置,所述方法可以相对于第一点P1在第二点P2处非常局部地探测介质,其中来自这两个点的信号之间具有额外延迟δt。这种局部信息完全包含在从介质的聚焦反射矩阵RFoc(rin,rout,δt)计算的时间响应的值中,并且可以在事后(且无任何新发射和/或采集)用于标准介质的每个点。
因此,可以通过考虑共焦信号的绝对值而从波束成形之后的此时间响应推导出介质的反射率的估计值,所述共焦信号的特征在于,在输入和输出rin=rout以及在零额外延迟δt=0时(即,在不具有此额外延迟的弹道时间)空间位置相等。介质的反射率的此估计值是介质的超声图像的像素值。因此,为了构建超声图像,一个人可以扫描或选择对应于超声图像中的一组像素位置的一组空间位置r=rin=rout
然后可以通过取r=rin=rout以及δt=0基于聚焦反射矩阵RFoc(rin,rout,δt)构建超声图像I(0)(r),这表示:
I(0)(r)=RFoc(rin,routin,δt=0)
(方程式4)
围绕介质点的传播的图像
然后可以通过基于聚焦反射矩阵RFoc(rin,rout,δt)构造一个或多个传播图像来补充由系统40的计算单元42实施的用于超声表征的方法,针对额外延迟δt的一个或多个值、针对一个输入虚拟换能器TVin(第一点P1)以及针对多个输出虚拟换能器TVout(第二点P2)确定此或这些传播图像,输出虚拟换能器TVout位于围绕空间位置rin的输入虚拟换能器TVin的空间位置rout处。
在单个传播图像的情况下,针对单个预定额外延迟δt从聚焦反射矩阵RFoc(rin,rout,δt)确定此传播图像。
此传播图像表示超声波在虚拟换能器之间且例如接近输入虚拟换能器在等于额外延迟的时刻(视为相对于弹道时间的时刻)传播的方式。
系统40随后任选地能够将一个或多个传播图像显示在显示装置43上。
计算单元42还可以针对若干额外的时间上连续延迟计算一系列传播图像,例如以在输入虚拟换能器TVin(第一点P1)周围构建超声波的传播膜。此传播膜可以任选地显示在显示装置43上或任何其它介质上。
在我们的实例中,在额外延迟间隔内应用为了构建此传播膜而应用的时间上连续的额外延迟。
例如,额外延迟间隔可以采用时间跨度的形式,所述时间跨度适合于从空间位置rin的输入虚拟换能器TVin延伸到空间位置rout的所有输出虚拟换能器TVout。此额外延迟间隔随后例如表示为[-δtmin,+δtmax],其中δtmin=zout max-zin/c0和δtmax=zout min-zin/c0,其中zin和zout分别是在空间位置rin的输入虚拟换能器TVin和空间位置rout的输出虚拟换能器TVout的第二轴Z的正方向上的深度。
例如,额外延迟间隔可以围绕零值(δt=0)对称且具有振幅δtmax,此额外延迟间隔表示为[-δtmax,+δtmax]。例如,对于用于传播图像的输出换能器TVout,这可以由δtmax=max(|Δr|)/c0定义。
图5中表示为A的图像示出分析的样本介质或体模的超声图像,其包括若干类型的预定异质性。在此介质中,我们考虑通过计算扫描的矩形分析区域ZAout(由输出虚拟换能器TVout的第二点P2组成),以便围绕空间位置rin的输入虚拟换能器TVin的第一点P1构建一个或多个传播图像,所述输入虚拟换能器此处位于分析区域ZAout的介质中。分析区域可以独立于输入虚拟换能器的位置而位于任何位置处。然而,具有围绕输入虚拟换能器的分析区域备受关注。
在此参考图像A上,输入虚拟换能器TVin(第一点P1)位于介质的反射元件(回声目标)上或附近。
图5中表示为B到F的图像是针对五(5)个额外延迟值δt的图5中的图像A的分析区域ZAout的传播图像。在我们的说明性实例中,这些额外延迟是-3.86μs、-1.93μs、0μs、1.93μs、3.86μs。每个传播图像由以下项组成:
-与分析区域ZAout中的一组点的聚焦反射矩阵的值的振幅相对应的索引1(例如B1)的第一图像,以及
-与分析区域ZAout中的同一组点的聚焦反射矩阵的值的实部相对应的索引2(例如B2)的第二图像。
在这些图像中,振幅的水平或实部的水平由灰度表示,其中图例出现在图5的图像B1和B2中。这些传播图像的点或像素对于空间位置具有Δr=rout-rin,这表示空间位置rout的输出虚拟换能器TVout关于输入虚拟换能器TVin的位置rin的相对位置。在说明此实例的图式中,在这些图像上的坐标在横坐标表示为Δx且在纵坐标上表示为Δz。
这些传播图像说明上面给出的关于使用额外延迟δt计算的聚焦反射矩阵的解释:它们可以使相干波的传播可视化。具体来说,对于趋向于零的负额外延迟,此相干波朝向输入虚拟换能器TVin的第一点P1会聚,并且理想地集中并聚焦在由零额外延迟(δt=0)的衍射极限定义的焦斑中。此相干波然后针对正的且增加的额外延迟发散。
此相干波由来自位于空间位置rin的输入虚拟换能器处的虚拟源的回波的数字时间反转过程产生,这些回波由探头的换能器测量。通过对于空间位置rin的输入虚拟换能器周围的一组空间位置rout在接收时以及在各种额外时间δt执行波束形成,这说明在共焦位置外接收时的聚焦(即rin=rout)。
由于针对位于介质的反射元件(回声目标)上或附近的输入虚拟换能器TVin的第一点P1获得这些传播图像,因此相干波在这些传播图像中可容易地识别并且与相邻信号相比呈现良好信噪比。
图6中表示为A的图像示出与图5相同的超声图像,但是其中考虑另一矩形分析区域ZA'out(在此实例中具有相同维度),所述另一矩形分析区域通过计算扫描,以围绕空间位置rin的输入虚拟换能器TV'in的另一第一点P1'构建传播图像。
输入虚拟换能器TV'in的此另一第一点P1'此处与分辨率单元相关联,所述分辨率单元包含任意地布置且具有相当的反射率的一组子分辨率散射体。在波长尺度下,这种介质称为“超声散斑”并且特征在于由每个子分辨率散射体之间的破坏性和建设性相互作用产生的随机反射率,从而造成B模式超声图像的颗粒效应。
图6中表示为B到F的图像是对于与图5中的图像B到F相同的五个额外延迟值δt,在图6中的图像A的此另一分析区域ZA'out的传播图像。
此另一分析区域ZA'out的一组第二点的聚焦反射矩阵的值的振幅和实部在此处以相同方式表示。
散射体的这些传播图像还示出会聚、集中在输入虚拟换能器TV'in的第一点P1'处,然后发散的相干超声波(由虚线环绕)。然而,由于通过位于焦平面上游或下游的散射体产生的回波具有与所分析的虚拟源的反射率相当的反射率,因此识别这种相干波更加困难。
此外,并且与对传播图像的先前定义相反,还可以在多个输入虚拟换能器TVin(第一点P1)与输出虚拟换能器TVout(第二点P2)之间构建一个或多个传播图像。因此,传播图像由聚焦反射矩阵RFoc(rin,rout,δt)构成,针对额外延迟δt的一个或多个值、针对一个输出虚拟换能器TVout(第二点P2)以及针对多个输入虚拟换能器TVin(第一点P1)确定此或这些传播图像,输入虚拟换能器TVin位于围绕空间位置rout的输出虚拟换能器TVout的空间位置rin处。
传播图像相对于输入和输出换能器实际上相反。由于波传播的互易性,产生的图像非常相似并且根据这些传播图像执行以及下文解释的各种计算和确定可以以类似的方式执行。为简单起见,本发明的详细描述将仅解释输入换能器与多个输出虚拟换能器之间的第一方向。然而,应当理解,在本文件中出现的每个定义中,可以互换具有“out”索引和“in”索引的元素以及术语“输入”和“输出”。
另外,还可以使用两种类型的传播图像(在第一和第二方向上)并且将它们组合或将这两个传播图像平均化以获得与介质中的波传播相比更具代表性和更具对比度的平均传播图像。还可以组合来自这两种类型的图像或从这两种类型的图像确定的结果,以获得通常更精确的结果。
如上文所定义的聚焦反射矩阵RFoc(rin,rout,δt)使用输入虚拟换能器TVin和输出虚拟换能器TVout的空间位置rin、rout。这些空间位置是空间参考系统内的绝对位置。然而,还可以使用单个绝对空间位置以及相对于此绝对空间位置的空间位置。例如,我们可以采用输入虚拟换能器的绝对空间位置rin以及输出虚拟换能器的相对空间位置Δrout,其中Δrout=rout-rin。相反地,我们可以采用输出虚拟换能器的绝对空间位置rout和输入虚拟换能器的相对空间位置Δrin,其中Δrin=rin-rout。本说明书中的每个计算和/或确定可以使用前述定义或任何其它类似和/或等效定义中的任一个来执行。
相干波提取
可以通过应用组合步骤来补充由系统40的计算单元42实施的用于超声表征的方法,其中执行一组传播膜的线性组合,所述集合的每个传播膜在不同空间位置rin的选定输入虚拟换能器TVin与空间位置rout的输出虚拟换能器TVout之间捕获,使得rout=Δrout+rin,其中Δrout预定义且针对所述集合的所有传播膜相同,并且选定输入虚拟换能器彼此靠近。
换句话说,选择选定输入虚拟换能器TVin的一组相邻空间位置,这组空间位置形成相关性的感兴趣区域,更简单地称为空间相关区域ZC,并且可以使这些输入虚拟换能器的传播膜相关。此空间相关区域例如是参考点周围的矩形区域。它也可以是整个图像,或者形状对称或不对称的任何区域。相邻空间位置例如是彼此靠近的空间位置。
通过一组若干传播膜的此组合,获得例如在相干性和对比度方面得到改进的改进相干波传播膜。针对相同额外延迟δt并且针对相同相对位置Δrout获得此新传播膜(称为相干波传播膜)的图像。
此新相干波传播膜随后可以与空间位置rin,ref的参考输入虚拟换能器TVin,ref相关联,所述参考输入虚拟换能器表示来自传播膜的集合的选定输入虚拟换能器(空间相关区域的输入虚拟换能器)。
根据第一实例,参考输入虚拟换能器TVin,ref是在与选定输入虚拟换能器的空间位置的平均值相对应的空间位置处的输入虚拟换能器。因此,在此以上变体中,可以通过以下公式表达参考输入虚拟换能器的空间位置:
Figure BDA0003261953800000161
Figure BDA0003261953800000162
是所选择的输入虚拟换能器,
Figure BDA0003261953800000163
是构成空间相关区域的选定输入虚拟换能器的数目。
根据另一实例,参考输入虚拟换能器TVin,ref是在与选定输入虚拟换能器的空间位置的加权平均值相对应的空间位置处的输入虚拟换能器,所述权重例如基于选定输入虚拟换能器的每个点的反射率值。因此,在此变体中,可以通过以下公式表达参考输入虚拟换能器的空间位置:
Figure BDA0003261953800000164
例如,通过表示为SVD的奇异值分解确定或执行此线性组合,在所述奇异值分解期间计算传播膜的集合的奇异值分解以获得与最大绝对值的奇异值相关联的奇异向量V1,此奇异向量V1随后是与所述参考输入虚拟换能器TVin,ref相关联且针对相同额外延迟δt的相干波传播膜。
所述集合的多个传播膜此处通过奇异值分解处理以便组合多个薄膜,这意味着在靠近输入虚拟换能器的区域中进行多次声学无序测量或实验,这可以提高传播膜对比度且因此有利地改进其使用。
为了执行此奇异值分解计算(特别是因为当前传统的奇异值分解工具使用二维矩阵),可以构建级联聚焦反射矩阵RFoc',其中此级联聚焦反射矩阵RFoc'的行是空间位置rin的输入虚拟换能器TVin的选定索引,并且此级联聚焦反射矩阵RFoc'的列是每个选定输入虚拟换能器TVin的级联传播膜{Δrout,δt}(图像集),这些传播膜对于相同时间序列的额外延迟δt获得。因此,此级联聚焦反射矩阵是再聚焦于输入rin处的聚焦点上的聚焦反射矩阵RFoc。
例如,写入此级联聚焦反射矩阵RFoc':
RFoc′=[RFoc(rin,{Δrout,δt})]=[RFoc(rin,{rin+Δrout,δt})]
奇异值分解SVD的此步骤随后提供奇异向量V1,所述奇异向量将选定输入虚拟换能器TVin的每个源之间的相关性最大化。奇异向量V1与来自奇异值分解的最大绝对值的奇异值相关联。奇异向量V1随后是与参考输入虚拟换能器TVin,ref且针对相同额外延迟δt的相干波传播膜。
因此,奇异值分解SVD的使用可以组合多个波传播膜,同时避免由散斑类型条件引入的随机反射率。由于相干波是每个传播膜的共同元素,因此它在组合过程中出现,而位于每个输入虚拟换能器TVin外部的散射体的贡献被相消干涉擦除。这相当于对传播膜应用滤波以提取相干波。
图7中表示为A的图像示出与图5和6中相同的超声图像。在此图中,考虑与从选定输入虚拟换能器TVin的集合中选择的输入虚拟换能器相关联的分析区域ZAout的实例,这些选定虚拟换能器在此图像A中由矩形点的网格表示。此网格的点表示用于执行传播膜的相干组合的选定输入虚拟换能器TVin的集合,称为相干区域ZC(即,相邻的选定输入虚拟换能器)。
图7中表示为B到F的图像是对于表示第一奇异向量V1的若干额外延迟值δt,图7的图像A的分析区域ZAout的相干波传播图像。此实例中使用与前面图中所示相同的表示,即第一图像为振幅且第二图像为实部。
图7中的图像示出也可以对于位于散斑中的一组第一点P1(输入虚拟换能器TVin)提取超声波的相干部分。实际上,在这些图像中,观察到从底部向顶部移动,同时集中在输入虚拟换能器TVin的位置处的单个相干波,而不通过本实验的奇异值分解过程SVD处理的传播图像将类似于图6中表示为B到F的图像中呈现的那些传播图像。
奇异值分解可以以非常可靠的方式从传播图像/膜提取相干波。例如,在图8A中,第一曲线A1对应于共焦点处的信号的振幅,所述振幅作为针对属于区域ZC的输入虚拟换能器TVin获得的一个传播膜的额外延迟δt(此处在-3μs与+3μs之间)的函数。共焦点是由Δx=Δz=|Δr|=0(rin=rout)定义的传播图像的点,其在我们的实例中由位于每个传播图像的中心处的x表示并且对应于输入虚拟换能器的位置。在所表示的情况下,此曲线A1非常混乱,因为共焦位置|r|与“散斑”类型区域重合。然后相干波由来自位于此共焦区域上游或下游的散射体的回波完全或部分隐藏。在此图中,曲线A2对应于从前述传播膜的奇异值分解产生且针对相同共焦点的相干波传播膜(第一奇异向量V1)的信号的振幅。此曲线A2示出以零额外延迟δt为中心的单个峰值,这表明即使对于涉及反射不良元件的这种特定情况也能很好地聚焦波。
图8B示出图8A的信号的频谱,曲线S1对应于曲线A1的信号的频谱并且曲线S2对应于曲线A2的信号的频谱。然而,观察到相干波的时间分辨率的损失(在图7中可见),这导致所研究的信号的谱宽减小。必要时,可以通过应用频谱均衡步骤来校正这种现象。
相干波传播图像类似于与回声散射体相关联的传播图像,但谱宽减小。
这些曲线A2、S2说明用于提取或过滤具有单个峰值(单个主波)的相干波传播膜的组合/奇异值分解步骤的效率。
弹道参考系统中的相干波
聚焦反射矩阵RFoc(rin,rout,δt)的计算假设超声波在介质中的速度模型(例如,恒定的声速c0)。实际上,波的向外飞行时间τin和返回飞行时间τout通常用几何公式计算,用于计算换能器11与介质中每个点之间的距离并假设声速恒定。
因此,上文所计算的传播图像、传播膜和相干波传播膜包含恒定声速c0的此假设。在这些图像和膜中,相干波是基于假设的声速模型的数字时间反转过程产生的。因此,此波以假设声速c0传播。在时间δt=0处,它位于输入虚拟换能器TVin的深度处(在这些图中的中心x),这意味着Δz=0因此,相干波的飞行时间遵循以下弹道传播关系:
δt(Δrout)=-sign(Δzout)·|Δrout|/c0
(方程式7)
其中:
c0是介质中的声速,
|Δrout|是输入虚拟换能器TVin与输出虚拟换能器TVout之间的向量的模数,Δrout=rout-rin
δt是额外延迟,
Δzout是沿着空间位置向量Δrout的第二轴Z的分量。
换句话说,在这些传播图像中,以声速c0传播的理论波形成以图像原点(即,空间位置rin的输入虚拟换能器TVin)为中心的圆弧。因此,弹道传播关系通过声速c0将相对位置Δrout与额外延迟δt联系起来。负号强调这是一个数字时间反转的过程。
然后可以从传播膜或从相干波传播膜中提取聚焦在弹道参考系内的波上的图像,此图像称为波前图像并以声速c0跟随此理论波:对于每个传播图像或相干波传播图像,在额外延迟δt处,我们提取位于此圆弧上的值(声压值)(即满足以上弹道传播关系)。因此构造称为波前图像的新图像,其表示弹道参考系统中的传播膜或相干波传播膜的演进。因此,此波前图像是弹道参考系统内的波前图像。
根据第一变体,通过计算传播膜或相干波传播膜并且如上文所解释通过从此膜提取适当数据来间接地确定波前图像,以便确定在额外延迟间隔期间的波前图像。
因此,可以通过应用针对输入虚拟换能器TVin或针对参考输入虚拟换能器TVin,ref以及针对额外延迟间隔确定波前图像的步骤来补充由系统40的计算单元42实施的用于超声表征的方法,所述波前图像从以下项确定:
-传播膜或相干波传播膜的图像,以及
-以下类型的弹道传播关系:
δt(Δrout)=-sign(Δzout)·|Δrout|/c0,其可以从膜的图像中的每一个提取值以便构建波前图像。
根据第二变体,通过施加以上弹道传播关系从实验反射矩阵Rui(t)直接确定波前图像。
因此,可以通过应用针对输入虚拟换能器TVin以及针对额外延迟间隔确定波前图像的步骤来补充由系统40的计算单元42实施的用于超声表征的方法,所述波前图像从以下项确定:
-聚焦反射矩阵RFoc(rin,rout,δt),以及
-以下类型的弹道传播关系
δt(Δrout)=-sign(Δzout)·|Δrout|/c0,其可以从聚焦反射矩阵提取值以便构建波前图像。
在所有这些变体中,波前图像可以基于通过探头的换能器测量的回波而估计由输入虚拟换能器TVin或参考输入虚拟换能器TVin,ref产生的压力场(在发射-接收期间的响应)。
应注意,包含在波前图像中的信号是聚焦反射矩阵的子矩阵。因此,对于计算,我们可以将自己限制于满足以上弹道传播关系的信号。在这种情况下,波前图像是聚焦反射矩阵RFoc(rin,rout,δt)。
这些波前图像的点或像素具有空间位置Δrout=rout-rin,这意味着相对于输入虚拟换能器TVin的位置rin的位置。因此,坐标在这些图像的横坐标上表示为Δx且在纵坐标上表示为Δz。还可以针对三维成像方法确定这些波前图像。其它坐标随后用于表示各种平面中的波前图像。
图9A示出此波前图像的振幅,并且图9B示出此波前图像的实部。在图9B中,注意在通过输入虚拟换能器TVin(在此图中坐标Δx=Δz=0的空间位置rin)的聚焦点的通道处π弧度的相移。此相移被称为戈伊相移。波前图像清楚地说明此现象。
与传播图像一样,可以反转输入虚拟换能器TVin和输出虚拟换能器TVout所扮演的角色。在这种情况下,获得通过聚焦产生的压力场的估计值作为输出。
积分声速的确定
根据本公开的且由系统40的计算单元42实施的用于介质的超声表征的方法和系统还能够确定在介质中的点处的积分声速。积分声速是探测装置41的换能器与介质的点之间的声速的平均值的估计。更精确地说,此积分声速将超声波的向外然后返回路径穿过的区域的所有局部声速积分。
在这种情况下,所述方法包括:
-针对输入虚拟换能器TVin以及针对额外延迟间隔确定波前图像的步骤,所述波前图像如上所述根据介质中的声速c0确定,
-确定输入虚拟换能器TVin的波前图像中焦斑中心的深度位置Δz(0)的步骤,以及
-基于以下公式计算积分声速c(1)的步骤:
Figure BDA0003261953800000201
其中zin是沿着输入虚拟换能器TVin的空间位置向量rin的第二轴Z的分量。
“波前图像中的焦斑的中心”应理解成表示例如波前图像中的焦斑的最大值的位置;表示在整个波前图像中具有最大值的像素的位置。应注意,在波前图像中可以仅观察到一个焦斑,因此其位置是唯一的。因此,焦斑的中心的位置也是唯一的,并且对于与输入虚拟换能器TVin的空间位置rin相对应的介质中的点,表示用于校正声速c0的深度位置Δz(0)(rin)。
例如,通过在波前图像中搜索具有最大值的点的空间位置来确定焦斑的中心,并且焦斑的中心的深度位置Δz(0)随后是在此最大值点的深度轴Z的方向上(对应于轴线Δz)的分量。
应注意,针对介质中获得的每个输入虚拟换能器TVin或相反地针对介质中获取的每个输出虚拟换能器TVout确定深度位置Δz(0)。更一般来说,此深度位置取决于空间位置r的每个考虑点并且可以表示为Δz(0)(r),其中r=rin或r=rout
实际上,在传播膜或相干波传播膜的图像中,仅当用于通过向外飞行时间和返回飞行时间的计算来计算聚焦反射矩阵RFoc(rin,rout,δt)以及通过弹道传播关系计算波前图像的声速c0是声速值时,在零额外延迟δt(δt=0)的时刻聚焦超声波,所述声速值对应于校正探测装置41的换能器11之间的实际介质以及对应于空间位置rin的输入虚拟换能器TVin的介质点的积分声速。
例如,图10说明此过程。在此图10中,表示为A、B和C的图像示出分别以1440m/s、1540m/s和1640m/s的预定义声速c0获得的波前图像。在这些波前图像中,我们观察到沿着纵轴Δz(表示深度(方向Z))移动的焦斑。表示为D的图形随后在此纵轴Δz(表示其中Δx=0的轴线)上示出这些波前图像的三个强度曲线CIA、CIB和CIC
例如,如图10中所说明获得焦斑的深度位置Δz(0)(rin),即通过确定波前图像的值中的最大值在纵轴Δz上的深度位置,使得Δx=0。通过在波前图像中搜索在纵轴Δz上具有此波前图像中的最大值的位置来获得波前图像中的焦斑的深度位置Δz(0)(rin),此纵轴Δz对应于波前图像中的零Δx横坐标。
例如,对于图形D的强度曲线CIA,深度位置Δz(0)(rin)基本上等于4.5mm,这将在选定输入虚拟换能器TVin的位置rin处产生大于初始假设的声速c(0)的积分声速c(1)(rin)的估计,使得沿着图像A的焦斑的Δz轴的垂直位置将向上且因此朝向输入虚拟换能器的原点(Δx=Δz=0)移动,这对应于通过计算在输入虚拟换能器TVin的介质中的此点的积分声速进行的调整。
因此,实际上,可以通过计算Δz轴上的波前图像的值来满足以确定声速或积分声速,其中Δx=0。
因此,用于介质的超声表征的方法包括用于确定积分声速的以下步骤:
-通过换能器11阵列10在所述介质的区域中生成一系列入射超声波USin的步骤,所述系列的入射超声波是发射基i;以及
-生成在作为输入的发射基i与作为输出的接收基u之间定义的实验反射矩阵Rui(t)的步骤;
-确定聚焦反射矩阵RFoc(rin,rout,δt)的步骤,所述聚焦反射矩阵包括在空间位置rin的输入虚拟换能器TVin与空间位置rout的输出虚拟换能器TVout之间的介质的响应,在相对于输入虚拟换能器TVin的响应的时刻移位额外延迟δt的时刻获得输出虚拟换能器TVout的响应,
-针对输入虚拟换能器TVin以及针对额外延迟间隔确定波前图像的步骤,所述波前图像根据介质中的声速c0确定,并且所述波前图像基于以下项确定:
-聚焦反射矩阵RFoc(rin,rout,δt),以及
-以下类型的弹道传播关系
δt(Δrout)=-sign(Δzout)·|Δrout|/c0,其可以从聚焦反射矩阵提取值以便构建波前图像,并且其中:
δt是额外延迟,
|Δrout|是输入虚拟换能器TVin与输出虚拟换能器TVout之间的向量的模数,其中Δrout=rout-rin
Δzout是沿着空间位置向量Δrout的深度轴Z的分量,
-确定波前图像中的焦斑的中心的深度位置Δz(0)的步骤,以及
-从以下公式计算积分声速c(1)的步骤:
Figure BDA0003261953800000221
其中zin是沿着输入虚拟换能器TVin的空间位置向量rin的深度轴Z的分量。
任选地,通过基于用先前的积分声速c(n)获得的波前图像的确定、焦斑的中心的深度位置Δz(n)的确定,以及新的积分声速c(n)的计算以相同迭代公式计算新的积分声速c(n +1),此方法可以如上文所定义迭代一次或多次:
Figure BDA0003261953800000222
实际上,此迭代过程极快地收敛到最佳积分声速,所述最佳积分声速对应于探测装置的换能器11和介质中的所选择点(输入虚拟换能器)的最佳积分声速。
此外,或者,用于确定积分声速的此方法可以通过在确定波前图像的步骤与确定焦斑的深度位置Δz(0)(rin)的步骤之间执行改进波前图像的步骤来改进,其中执行对应于给定相干区域ZC的一组波前图像的线性组合,在不同空间位置rin的选定输入虚拟换能器TVin与空间位置rout的输出虚拟换能器TVout之间获得所述集合的每个波前图像,使得rout=Δrout+rin,其中Δrout预定义且针对所述集合的所有波前图像相同,并且选定输入虚拟换能器彼此靠近。因此获得与参考输入虚拟换能器TVin,ref相关联的改进的波前图像或相干波前图像,此参考输入虚拟换能器TVin,ref表示所使用的波前图像集合的输入虚拟换能器,针对相同相对位置Δrout与所选择的相干区域ZC相关联。
例如,如上文已针对传播膜的情况所解释,参考输入虚拟换能器TVin,ref是与选定输入虚拟换能器的空间位置的平均值或选定输入虚拟换能器的空间位置的加权平均值相对应的空间位置的输入虚拟换能器。
总而言之,在本公开的方法中,添加以下步骤:
-在确定波前图像的步骤与确定焦斑的深度位置Δz(0)(rin)的步骤之间,执行改进波前图像的步骤,其中执行对应于相干区域的一组波前图像的线性组合,在不同空间位置rin的选定输入虚拟换能器(TVin)与空间位置rout的输出虚拟换能器(TVout)之间获得每个波前图像,使得rout=Δrout+rin,其中Δrout预定义且针对所述集合的所有波前图像相同,并且选定输入虚拟换能器彼此靠近,以便获得与参考输入虚拟换能器(TVin,ref)相关联的改进的波前图像,此参考输入虚拟换能器TVin,ref表示所使用的波前图像集合的输入虚拟换能器并且与相干区域ZC相关联,以及
-在确定深度位置Δz(0)(rin)的步骤中,代替波前图像使用改进的波前图像,焦斑的中心的深度位置与参考输入虚拟换能器TVin,ref的空间位置相关,并且焦斑的中心的此深度位置可以估计在参考输入虚拟换能器TVin,ref的空间位置处的积分声速c(1)(rin,ref)。
改进的波前图像(相干波前图像)随后用于(代替波前图像)确定焦斑的中心的轴向位置。此距离或深度位置Δz(0)(rin,ref)随后表示声速的不正确模型,并且可以用于估计与参考输入虚拟换能器TVin,ref的空间位置rin,ref相关联的积分声速c(1)(rin,ref)。
根据一个实施例,通过计算波前图像的集合的奇异值分解SVD以获得与奇异值分解的最大绝对值的奇异值相关联的奇异向量W1来确定线性组合,此奇异向量W1随后是对应于所述参考输入虚拟换能器TVin,ref且针对相同额外延迟δt的改进的波前图像。
所述集合的多个波前图像此处可以通过奇异值分解处理以便组合在靠近输入虚拟换能器的区域中的若干声学无序测量或实验,这可以避免与无序相关的波动且改进波前图像的对比度以及其使用。
另外,可以通过计算如上所述的积分声速以及通过针对波前图像的集合的线性组合使用与基本上覆盖介质中的整个感兴趣区域的选定输入虚拟换能器(TVin)相对应的一组波前图像来确定介质的最佳声速(对于整个介质而言是真实的)。具体来说,这些选定输入虚拟换能器可以以预定间隔规则地分布在介质的整个感兴趣区域上。例如,这些选定输入虚拟换能器可以表示例如用于构建覆盖待研究区域的介质的超声图像的输入虚拟换能器的数目的20%或更多。
应注意,由于在来自空间位置rin或rin,ref的回波的反向传播步骤中所经历的像差,深度距离或位置Δz(0)(rin)或Δz(0)(rin,ref)可以解释为输出处的聚焦误差。还可以通过探测在向外路径期间波前所经历的像差来确定积分声速测量。此测量通过反转上述方程式中的“输入”和“输出”指示,同时反转输入和输出虚拟换能器的角色来描述,以获得积分声速c(1) out的另一个估计。
另外,可以通过组合从在向外和/或返回行程上生成的像差获得积分声速的测量或估计,即积分声速c(1) in和c(1) out来改进积分声速的估计。
所述方法随后通过以下步骤补充:
-输入虚拟换能器和输出虚拟换能器的角色反转,以确定相对于输出虚拟换能器的积分声速c(1)(rout),以及
-将参考输入虚拟换能器的积分声速c(1)(rin)与参考输出虚拟换能器的积分声速c(1)(rout)组合以获得改进的积分声速。
积分声速图像
由系统40的计算单元42实施的用于超声表征的方法可以通过构建一个或多个 分声速图像来补充,通过上文针对对应于空间位置rin的输入虚拟换能器TVin(第一点P1)的介质中的多个点描述的积分声速的至少一个计算来确定此或这些积分声速图像。
图11说明此图像的两个实例。
在第一实例中,对应于图11的图像A1和A2,介质是具有基本上为cref=1542m/s给出的参考声速的体模类型。图像A1是标准超声图像,而图像A2是通过以上方法获得的积分声速图像。此积分声音图像A2允许估计介质中声速1544m/s的平均值,标准偏差为+/-3m/s,这与此介质的声速的参考值完全一致。
在第二实例中,对应于图11的图像B1和B2,介质是分层介质,第一层几乎20毫米厚,具有纤维结构和大约1570m/s的声速,位于相同的体模类型介质上,具有由构造预定的基本上1542m/s的声速。图像B1是此介质的标准超声图像,而图像B2是通过上文公开的方法获得的积分声速图像。此积分声音图像B2反映第一层中约1580m/s的较高声速,以及低于但不等于预期声速且对应于在第一实例中研究的此介质的声速的较低声速。这种效果归因于以下事实:通过此方法计算出的声速是积分声速,所述积分声速对应于换能器11与介质的点之间的波的整个向外和返回路径上的平均或积分声速。
轴向像差的校正
根据本公开的且由系统40的计算单元42实施的用于介质的超声表征的方法和系统还能够确定轴向校正
用于介质的超声表征以通过轴向校正确定超声聚焦的时间和局部标准的方法包括已针对获得聚焦反射矩阵说明的以下步骤:
-通过换能器11阵列10在所述介质的区域中生成一系列入射超声波USin的步骤,所述系列的入射超声波是发射基i;以及
-生成在作为输入的发射基i与作为输出的接收基u之间定义的实验反射矩阵Rui(t)的步骤;
-确定聚焦反射矩阵RFoc(rin,rout,δt)的步骤,所述聚焦反射矩阵包括在空间位置rin的输入虚拟换能器TVin与空间位置rout的输出虚拟换能器TVout之间的介质的响应,在相对于输入虚拟换能器TVin的响应的时刻移位额外延迟δt的时刻获得输出虚拟换能器TVout的响应,以及
-针对输入虚拟换能器TVin以及针对额外延迟间隔确定波前图像的步骤,所述波前图像如上所述根据介质中的声速c0确定,
-确定波前图像中的焦斑的中心的深度位置Δz(0)(rin)的步骤。
“波前图像中的焦斑的中心”应理解成表示例如波前图像中的焦斑的最大值的位置;表示在整个波前图像中具有最大值的像素的位置。可以根据上文已经描述的技术中的一个找到/确定焦斑的中心和深度位置。
此方法进一步包括通过在深度方向Z上对空间平移的聚焦反射矩阵RFoc(rin,rout,δt)的响应进行平移来确定校正聚焦反射矩阵RFoc(1)(rin,rout,δt)的步骤,所述空间平移是先前确定的深度位置Δz(0)(rin)的函数。
根据第一变体,通过具有等于2.Δz(0)(rin)的校正值Δzcorr(rin)的空间位置rout的输出虚拟换能器TVout的轴向分量的空间平移(沿着深度轴Z)执行空间平移,以获得校正的聚焦反射矩阵RFoc(1)(rin,rout,δt),使得:
RFoc(1)(rin,rout,δt)=RFoc(rin,{xout,zout+ΔzCorr(rout)},δt)
(方程式11)
校正的超声图像I(1)(rin)随后可以由校正的聚焦反射矩阵RFoc(1)(rin,rout,δt)(其特征在于,r=rin=rout和δt=0)构成以获得:
I(1)(r)=RFoc(1)(rin,rout=rin,δt=0)
(方程式12)
相反地,空间平移也可以对应于分量沿着具有等于2.Δz(0)(rin)的校正值Δzcorr(rin)的空间位置rin的输入虚拟换能器TVin的深度轴Z的空间平移,以获得随后的校正的聚焦反射矩阵RFoc(1)(rin,rout,δt)。
RFoc(1)(rin,rout,δt)=RFoc(xin,zin+ΔzCorr(rin),rOut,δt)
(方程式13)
应注意,深度位置Δz(0)针对介质中获取的每个输入虚拟换能器TVin确定并且表示返回行程期间所经历的像差。通过反转“输入”和“输出”符号,可以针对在介质中获取的每个输出虚拟换能器TVout确定表示在向外行程期间所经历的像差的位置Δz(0)(rout)。换句话说,更一般来说,此深度位置取决于空间位置r的每个考虑点并且可以表示为Δz(0)=Δz(0)(r),其中r=rin或r=rout
根据第二变体,通过以下操作执行空间平移:
-计算由以下公式确定的等于Δz(1)(r)的校正值Δzcorr(r):
Δz(1)(r)=z(1)(r)-zin
其中
Figure BDA0003261953800000261
其中此方程式适用于
r=rin和r=rout
z=zin和z=zout是沿着输入虚拟换能器TVin空间位置rin或输出虚拟换能器TVout的空间位置rout的深度轴Z的分量,
-通过分量沿着所述校正值Δzcorr(rin)的空间位置rin的输入虚拟换能器TVin以及所述校正值Δzcorr(rout)的空间位置rout的输出虚拟换能器TVout的深度轴Z进行空间平移来计算校正的聚焦反射矩阵RFoc(1)(rin,rout,δt),使得:
RFoc(1)(rin,rout,δt)=RFoc(xin,zin+ΔzCorr(rin)},{xouT,zouT+ΔzCorr(rout)},δt)
此计算还可以基于以下公式表达为积分声速c(1)(r)的函数:
Figure BDA0003261953800000262
其中zin是沿着输入虚拟换能器TVin的空间位置向量rin的第二轴Z的分量,以及
-平移计算Δz(1)(r)通过:
Δz(1)(r)=z(1)(r)-zin
其中
Figure BDA0003261953800000263
根据对上述两个变体的一些修改,可以通过计算空间傅里叶变换、通过相位斜坡相移(其中斜率取决于校正值)、然后空间逆傅里叶变换来实现平移。此实施方案的优点是可以为新的空间坐标组合平移和插值。
作为实例,因此实施的方法执行:
-根据方程式确定空间频率矩阵RFreqz(rin,xout,kzout,δt)的步骤,所述空间频率矩阵是在聚焦反射矩阵RFoc(rin,rout,δt)的深度方向上的空间傅里叶变换:
RFreqz(rin,xout,kz,out,δt)=TFzout[RFoc(rin,rout,δt)]
(方程式17)
其中
TFzout是深度方向Δzout上的空间傅里叶变换,
kzout是包括在间隔[ω-/c0+/c0]内的对应波数,其中脉冲ω-和ω+是界定超声波的带宽的脉冲,以及
xout是空间位置rout的每个输出虚拟换能器TVout的在X轴方向上的横向分量。
-根据以上变体通过深度校正值Δzcorr(即,等于2.Δz(0)(rin))的相位斜坡确定空间频率矩阵RFreqz(rin,xout,kzout,δt)的产品的对应于相同深度方向上的逆空间傅里叶变换的校正的聚焦反射矩阵RFoc(1)(rin,rout,δt),并且所述校正的聚焦反射矩阵针对输入虚拟换能器TVin的每个空间位置确定并且其中应用以下公式:
Figure BDA0003261953800000271
其中
e-ix是复指数函数,
Δzcorr是通过波前图像中的焦斑的中心的深度位置确定的校正值。
可以例如通过以下空间离散傅里叶变换公式描述在方向Δzout上的空间傅里叶变换:
Figure BDA0003261953800000272
存在傅里叶变换和空间傅里叶变换的其它公式。
随后可以通过以下反演公式说明在方向Δzout上的逆空间傅里叶变换:
Figure BDA0003261953800000273
Figure BDA0003261953800000281
根据第三变体,通过使用替代假设声速c0的新声速c1(r)计算或确定校正的聚焦反射矩阵RFoc(1)(rin,rout,δt)来轴向地平移响应。
此第三变体的方法因此进一步包括用于获得轴向校正的聚焦反射矩阵的以下步骤:
-从以下公式计算积分声速c(1)(r)的步骤:
Figure BDA0003261953800000282
其中zin是沿着输入虚拟换能器TVin的空间位置向量rin的第二轴Z的分量,以及
-计算校正的聚焦反射矩阵RFoc(1)(rin,rout,δt),所述校正的聚焦反射矩阵包括在空间位置rin的输入虚拟换能器TVin与空间位置rout的输出虚拟换能器TVout之间的介质的响应,所述响应各自通过取决于输入虚拟换能器的校正声速获得。
对于这些变体中的每一个,校正的聚焦反射矩阵RFoc(1)(rin,rout,δt)是聚焦反射矩阵,即已校正轴向像差的聚焦反射矩阵的轴向校正。由于此校正的聚焦反射矩阵,有利地可以构建具有减小的轴向像差的超声图像。因此,此校正的超声图像中的轴向方向上的距离更精确并且例如可以获得更好质量的图像。
通过空间平移获得校正的聚焦反射矩阵RFoc(1)(rin,rout,δt),所述空间平移是在一个或两个虚拟换能器(TVin和/或TVout)的轴向分量的Z方向上的空间位置的平移,或通过改变声速c的平移。这些替代方案可以改进波束形成步骤,这类似于将实验反射矩阵Rui(t)的实验信号(也通常称为RF信号)给出的时间信息通过关系t=z/c转换为空间信息的过程。因此,介质点的空间位置在深度方向Z上进行轴向校正,这可以获得具有更精确垂直定位的图像。
例如,图12说明此过程。在此图12中,表示为A的图像对应于通过c0=1540m/s的声速c0获得的超声图像,所述声速是体模中的声速,但不是体模上方水层中具有声速cwater=1480m/s=1480m/s的声速。因此,归因于所研究介质的异质性,通常通过降低的分辨率和对比度获得超声图像A。已知像差校正技术可以获得表示为B的图像,所述图像横向地改进且比常规技术提供更佳质量的图像。然而,在此图像中,不校正反射元件的深度位置(参见图像A与B之间的水平箭头)。
表示为C的图像对应于通过在上文呈现的方法中提出的轴向校正获得的超声图像。在此图像C中,反射元件略微向上移位(朝向外表面),这是处与体模相比水中的声速降低的影响。因此,归因于此轴向校正,图像的点的(深度)轴向位置更靠近观察到的介质的真实性质,并且在此图像中测量到的距离更靠近确切值。
还可以通过使用一组波前图像的组合来确定改进的波前图像以及在涉及确定积分声速的部分时如上文所解释确定输入和输出两者处的声速,以及例如通过奇异值分解技术来改进以上三个变体中的任一个的技术。
所述集合的多个波前图像此处通过奇异值分解处理以便组合在靠近输入虚拟换能器的区域中的若干声学无序测量或实验,这可以非常有利地改进波前图像的对比度以及其使用。
通过轴向像差校正的校正超声图像
由系统40的计算单元42实施的用于超声表征以确定轴向校正的方法随后可以通过构建一个或多个校正超声图像补充,通过以下操作确定校正超声图像:基于校正的聚焦反射矩阵RFoc(1)(rin,rout,δt)计算介质的多个点的超声强度值,每个点对应于空间位置rin的输入虚拟换能器TVin;以及施加与输入虚拟换能器TVin一致的输出虚拟换能器TVout,即rin=rout
介质中的散射体的各向异性的优选方向的确定
根据本公开的且由系统40的计算单元42实施的用于介质的超声表征的方法和系统还能够局部地确定介质中的散射体的各向异性的优选方向。
散射体各向异性表征任何散射体,当它在特定入射方向被声穿透时,所述散射体能够在优选方向上生成回波。因此,此各向异性涉及尺寸大于波长的任何散射体。具体来说,这在纤维、器官壁、手术器械(如活检针等)的医学成像中很重要。
在这种情况下,所述方法包括与上面已经解释的步骤相似或相同的步骤,直至:
-针对输入虚拟换能器TVin以及针对额外延迟间隔确定波前图像的步骤,所述波前图像如上所述根据介质中的声速c0确定。
所述方法还包括:
-通过对所述波前图像进行图像处理确定波前图像中的焦斑的优选方向的步骤。
例如,图13说明此过程。在此图13中,表示为A的图像对应于在组织的各向异性的方向上具有空间变化的超声图像。在此超声中成像的介质对应于患者的肌肉(在本此实例中为小腿),其中观察到几个区域的纤维在非常不同的方向上倾斜。这种对肌肉图像的应用仅是可以应用所述方法的各向异性介质的一个实例。然而,此常规超声图像提供失真的一般信息。现有技术的超声不允许可靠地观察作为介质局部特征的此种各向异性,因为超声波在这种介质中的传播既不是以恒定速度传播,也不是从探头的换能器在直线方向传播。
此图13中表示为B、C、D和E的图像对应于为超声图像的由箭头连接的小区域构建的波前的图像。此处,这些波前图像通过在这些区域中的每一个中的多个虚拟换能器的奇异值分解处理以捕获或探测此区域中的多个声学无序实验,并且因此改进所产生的波前图像的对比度以及其分析。
所有这些波前图像B、C、D和E示出在垂直方向(深度轴Δz的方向)上伸长,但具有不同倾斜的焦斑。波前图像中的焦斑的此倾斜是与所考虑区域中的肌纤维的倾斜的实际值高度相关的局部倾斜信息。焦斑的倾斜轴实际上基本上垂直于纤维的方向,特别是在图像的中心处,即入射波具有基本上在深度方向Z上的方向的位置。
因此,所述方法通过对此波前图像进行图像处理确定波前图像中的焦斑的优选方向。
根据第一变体,所述方法可以例如通过在低于此波前图像中的最大值的水平处,例如在此最大值的50%或70%处的阈值提取焦斑的轮廓。从这个轮廓,我们可以推导出优选方向或主方向(焦斑的最大尺寸方向)和次方向(最小尺寸方向)。然而,其它图像处理技术可以提取焦斑的优选方向。
根据第二变体,所述方法可以例如:
-将波前图像U(rin,Δxout,Δzout)从笛卡尔坐标中的参考系转换到类型U(rin,Δsout,Δφout)的极坐标中的参考系
-对所述波前图像在极坐标中的参考系中的值在多个径向距离偏差值Δsout上相加,以获得多个角度值Δφout的角灵敏度函数f(rin,Δφout),
-确定对应于角灵敏度函数的最大值的最佳角度值Δφmax out(rin),所述最佳角度值Δφmax out(rin)对应于与输入虚拟换能器TVin相关联的焦斑的优选方向。
因此我们具有:
Figure BDA0003261953800000301
图14示出使用对应于图13B的极坐标中的参考系的所述波前图像的角灵敏度函数f(rin,Δφout)的实例曲线,所述角灵敏度函数在此说明性实例中被归一化为具有等于一(1)的最大值。此曲线具有朝向Δφout=-11°的最大值,这是此示例的介质中相关点的局部优选方向的所估计角度。
任选地,校正应用于与从换能器观察到的介质中相关点的视角相对应的角度值Δφout以获得角各向异性值γout(rin),其表示位于输入虚拟换能器的空间位置处的散射体的各向异性方向。
此处,基于输出信号的相关性进行此估计。相反地,可以估计另一个角度各向异性值γin(rout),其表示位于输出虚拟换能器的空间位置处的散射体的各向异性方向。有利地,可以将两个角各向异性值γout(rin)和γin(rout)组合,以获得介质的各向异性方向的更佳局部表征。
根据一个实例,所述方法可以通过以下步骤补充:
-输入虚拟换能器的角色和输出虚拟换能器的角色反转,以确定相对于输出虚拟换能器的优选方向,以及
-将参考输入虚拟换能器的优选方向和参考输出虚拟换能器的优选方向组合以获得改进的优选方向。
根据另一实例,所述方法可以通过以下步骤补充:
-输入虚拟换能器和输出虚拟换能器的角色反转,以确定相对于输出虚拟换能器γout(rout)的角各向异性值,以及
-将参考输入虚拟换能器γout(rin)的角各向异性值以及参考输出虚拟换能器γout(rout)的角各向异性值组合以获得改进的角各向异性值。
可以通过以下公式给出角各向异性值的实例计算(在第一情况下,参考输入虚拟换能器γout(rin)的角各向异性值):
Figure BDA0003261953800000311
角各向异性值的此计算例如来自在以下文档中说明的计算:Rodriguez-Molares等人在《IEEE超声学、铁电体技术与频率控制汇刊》(第64卷第9期,2017年9月)中公开的“镜面波束成形”。
我们添加例如所述类型的输入虚拟换能器的空间位置rin点的视角的定义:
Figure BDA0003261953800000312
其中
Figure BDA0003261953800000313
其中:
Figure BDA0003261953800000321
是阵列中的换能器的空间位置的最大值和最小值。
本领域技术人员可以想到计算角各向异性值的其它公式,以便获得优选方向的角度的更现实的值。
图15示出超声图像,在所述超声图像上叠加与分布在超声图像上的一组点的优选方向的估计相对应的线。人们会注意到估计的优选方向与超声图像中可见的底层结构之间的高度一致性。所提出的方法有利地可以在整个超声图像上充分地估计优选方向。
此优选方向的测量,即焦斑最大尺寸的倾斜角是改进此区域中的超声图像质量的一个重要参数:了解这一点可以例如通过选择具有特定倾角的平面波或聚焦在特定位置的波而适应入射超声波USin的特性。这也使得能够在波束形成步骤期间适应在接收时选择的切趾。
此局部优选方向的测量可以分析较大区域的各向异性程度,由此确定组织中病灶的潜在存在且发现其位置。
因此,用于介质的超声表征以局部地确定各向异性的优选方向的方法包括以下步骤:
-通过换能器11阵列10在所述介质的区域中生成一系列入射超声波USin的步骤,所述系列的入射超声波是发射基i;以及
-生成在作为输入的发射基i与作为输出的接收基u之间定义的实验反射矩阵Rui(t)的步骤;
-确定聚焦反射矩阵RFoc(rin,rout,δt)的步骤,所述聚焦反射矩阵包括在空间位置rin的输入虚拟换能器TVin与空间位置rout的输出虚拟换能器TVout之间的介质的响应,在相对于输入虚拟换能器TVin的响应的时刻移位额外延迟δt的时刻获得输出虚拟换能器TVout的响应,
-针对输入虚拟换能器TVin以及针对额外延迟间隔确定波前图像的步骤,所述波前图像根据介质中的声速c0确定,并且所述波前图像基于以下项确定:
-聚焦反射矩阵RFoc(rin,rout,δt),以及
-以下类型的弹道传播关系
δt(Δrout)=-sign(Δzout)·|Δrout|/c0,其允许从聚焦反射矩阵提取值以便构建波前图像,并且其中:
δt是额外延迟,
|Δrout|是输入虚拟换能器TVin与输出虚拟换能器TVout之间的向量的模数,其中Δrout=rout-rin
Δzout是沿着空间位置向量Δrout的深度轴Z的分量,
-通过对所述波前图像进行图像处理确定波前图像中的焦斑的优选方向的步骤。
可以通过在涉及确定积分声速的部分中使用如上文所解释的一组波前图像的组合确定改进的波前图像,以及例如通过奇异值分解技术来改进所提出的技术。在这种情况下,从改进的波前图像获得的优选方向可以表征对应于所选择相干区域的介质的各向异性并且归因于参考虚拟换能器的空间位置rin,ref
此处,所述集合的多个波前图像通过奇异值分解处理以便组合在靠近输入虚拟换能器的区域中的若干声学无序测量或实验,这可以改进波前图像的对比度以及因此改进其使用。
在所述方法中,可以因此添加以下步骤:
-在确定波前图像的步骤与确定焦斑的优选方向的步骤之间,执行改进波前图像的步骤,其中执行对应于相干区域的一组波前图像的线性组合,在不同空间位置rin的选定输入虚拟换能器(TVin)与空间位置rout的输出虚拟换能器(TVout)之间获得每个波前图像,使得rout=Δrout+rin,其中Δrout预定义且针对所述集合的所有波前图像相同,并且选定输入虚拟换能器彼此靠近,以便获得与参考输入虚拟换能器(TVin,ref)相关联的改进的波前图 ,此参考输入虚拟换能器TVin,ref表示所使用的波前图像集合的输入虚拟换能器并且与所选择的相干区域ZC相关联,以及
-在确定焦斑的优选方向的步骤中,代替波前图像使用改进的波前图像;焦斑的优选方向与参考输入虚拟换能器TVin,ref的空间位置相关。
另外,通过反转输入虚拟换能器TVin和输出虚拟换能器TVout的角色,换句话说,通过反转“输入”和“输出”符号,可以确定与空间位置rout的输出虚拟换能器TVout相关联的焦斑的优选方向Δφmax in(rout)。将与位置r相关联的两个优选方向,即Δφmax in(r)和Δφmax out(r)组合可以改进散射体各向异性的测量。
由于优选方向和这些图像的此计算,我们可以表征介质的散射体的各向异性,或者表征例如介质中的各向异性结构,例如引入到组织的针,或分隔不同组织的壁。散射体的各向异性应理解成表示大于超声波的波长的任何元素。
用于共焦点的时间信号的分析
根据本公开的且由系统40的计算单元42实施的用于介质的超声表征的方法和系统还能够执行超声聚焦的局部频谱分析。
在此分析中,共焦响应备受关注,这意味着空间位置rin的输入虚拟换能器TVin叠加在空间位置rout的输出虚拟换能器TVout上;即其中rin=rout=r。
额外延迟δt随后用于探测由这些虚拟换能器选择的散射体的时间响应。
在这种情况下,所述方法包括已经解释以获得聚焦反射矩阵,但应用于同一空间位置,即共焦位置的以下步骤:
-通过换能器11阵列10在所述介质的区域中生成一系列入射超声波USin的步骤,所述系列的入射超声波是发射基i;以及
-生成在作为输入的发射基i与作为输出的接收基u之间定义的实验反射矩阵Rui(t)的步骤;
-确定聚焦反射矩阵RFoc(r,δt)的步骤,所述聚焦反射矩阵包括在空间位置rin的输入虚拟换能器TVin与空间位置rout的输出虚拟换能器TVout之间的介质的响应,输入和输出虚拟换能器叠加在相同空间位置r处,其中rin=rout=r,并且在相对于输入虚拟换能器TVin的响应的时刻移位额外延迟δt的时刻获得输出虚拟换能器TVout的响应。
所述方法随后包括可以执行局部频谱分析的以下步骤:
-确定频率矩阵RFreqt(r,ω)的步骤,所述频率矩阵是聚焦反射矩阵RFoc(r,δt)的时间傅里叶变换
RFreqt(r,ω)=TFt[RFoc(r,δt)]
(方程式25)
其中
TFt是时间傅里叶变换,以及
ω是脉冲,其中ω=2πf,f是对应于所述脉冲的频率。
可以例如通过以下离散时间傅里叶变换公式解释时间傅里叶变换:
Figure BDA0003261953800000341
其它傅里叶变换和时间傅里叶变换公式例如以离散或连续形式存在且归一化或没有归一化,并且也可以使用。
RFreqt(r,ω)随后含有由介质反向散射的回波的频谱的局部估计。更精确地说,这些回波来自包括在以位置r为中心的单色焦斑中的散射体。在不存在像差的情况下,这些尺寸因此由在介质反向散射的回波的中心频率处定义的衍射极限提供。
出于改进空间分辨率的目的,此方法因此可以基于反向散射回波的频率分析通过任何医学成像技术补充。具体来说,此方法允许在执行任何频谱分析之前在接收时针对每个频率进行空间波束成形。应注意,共焦配置有利地可以限制脉冲衍射现象。
例如,此方法可以通过滤波步骤补充,在所述滤波步骤期间执行频率矩阵RFreqt(r,ω)的元素的频率滤波。具体来说,可以执行低通、带通或高通频率滤波,以便根据目标应用在聚焦反射矩阵的响应中提取所需组件。例如,频率滤波可以任选地适合于提取入射超声波USin的基频的谐波分量。
例如,图16说明此过程。在此图16中,表示为A的图像说明包括气泡的介质的超声图像。这些气泡被破坏超声波图像的介质中的谐振结构,因为它们在入射波通过之后继续振荡。因此,它们生成以大于弹道飞行时间的飞行时间到达接收换能器的回波,这会在气泡下游的超声图像中产生伪影。图14中表示为B的图像示出图像A的放大,其中在空间位置r=[x,z]=[11,17]mm处观察到气泡的明亮回波,并且其下游伪影位于此位置下方(即,竖直深度)。表示为C1和C2的图像分别针对零的额外延迟δt对应于传播图像的振幅和实部。
所述方法的聚焦反射矩阵RFoc(r,δt)可以研究此气泡振荡的时间信号。图14中表示为D-E-F的图像分别对应于在对应于此气泡位置的空间位置r的点处响应RFoc(r,δt)的实部、振幅和频谱的图。在图像D和E上,在以δt=0为中心的主回波后约1.5μs处观察到第二回波。图像F示出不包括此第二回波的第一频谱图以及具有此第二回波的第二频谱图。此第二频谱图包括对应于入射波的频率的约6MHz的主频率,以及对应于气泡的共振频率(振荡)的约3MHz的另一频率。
因此,所述方法执行频谱分析,这可以例如识别观察到的介质中的气泡或任何其它谐振结构的谐振频率。
因此,可以例如通过预定带通滤波器对聚焦反射矩阵的响应进行滤波,并且然后使用这些滤波后的响应计算改进的超声图像。然后可以在超声图像中减弱或消除谐振的影响。
相反地,可以通过仅将这些谐振保持在聚焦反射矩阵的响应中来构建共振频率图像。应注意,气泡的谐振频率与其大小有关,并且可以用于估计介质中的局部压力。
在第二实例中,RFreqt(r,ω)可以用于研究介质的衰减。实际上,此现象取决于频率。由于高频比低频衰减得更多,因此可以例如通过比较来自观察到的介质中的两个不同深度的回波频谱来推导出衰减系数。因此,用于估计来自给定区域的回波的局部频谱的上述技术对于确定衰减是理想的。为实现这一点,所述方法可以例如通过确定深度S(z,ω)处的平均频谱的步骤补充,所述平均频谱由在介质中的预定深度z处的频率矩阵的频谱的平均值确定。
例如,在某一深度处的此平均频谱通过以下公式计算,所述平均频谱是归一化平均值且在预定间隔内包括的相同深度z和横向坐标x的一组空间位置上平均。
Figure BDA0003261953800000361
例如,图17说明通过构建超声图像的所有深度的频谱集合的图像进行的在深度方向上的平均频谱深度的此计算。在此图17中,表示为A的图像说明健康个体的小腿的体内超声图像并且表示为B的图像以灰度示出在深度方向上的平均频谱。深度频谱的此图像示出大深度的最强高频衰减。
使用此图像,我们可以通过使用全部频率内容通过用于在理论和/或实验模型与此图像之间调整的技术估计随深度而变的衰减的演变。
在第三实例中,所述方法还可以通过以下步骤补充:通过计算频率矩阵RFreqt(r,ω)的每个频谱的自相关的半高全宽,即通过以下公式确定空间位置r的点的频谱相关宽度δω(r):
Figure BDA0003261953800000362
其中
FWHM是用于计算半高全宽的函数
()*是共轭复数函数,
ω-和ω+是边界脉冲,Δω=ω+-是边界脉冲之间的间隔,即所涉及的超声波带宽。
归因于矩阵RFreqt(r,ω)的空间分辨率,频谱相关宽度δω(r)是可以用于标准以空间位置r为中心的单色焦斑中包含的散射体性质的局部值。如果焦斑包含单个非谐振散射体,则频谱相关宽度δω(r)具有超声波信号带宽的数量级。如果焦斑包含具有相同强度(超声散斑条件)的一组任意分布的散射体,则频谱相关宽度δω(r)的值变得比带宽Δω小得多。
所述方法还可以包括确定至少一个频谱相关图像的步骤,所述频谱相关图像通过确定介质的多个点的频谱宽度δω(r)获得,每个点对应于空间位置r的介质的点。
例如,图18说明此过程。在此图18中,表示为A的图像是包含若干不同元素的体模介质的超声图像:点目标和回声圆柱体。表示为B的对应图像是先前超声图像的频谱相关图像,其通过计算此介质的一组点的频谱相关宽度δω(r)获得。在此图像B中,圆柱体和点目标的边缘具有比由大量任意分布的子分辨率散射体组成的介质的其余部分大的频率相关宽度δω(r)。
通过频谱相关宽度和这些图像的此计算,可以表征介质中的目标的性质。例如,可以区分亮散斑点与单个散射体。例如,这可以有助于识别用于对比度成像的气泡,或表示存在肿瘤(尤其在乳腺癌中)的微钙化。

Claims (13)

1.用于介质的超声表征的方法,以确定所述介质中的积分声速,所述方法包括:
-通过换能器(11)阵列(10)在所述介质的区域中生成一系列入射超声波(USin)的步骤,所述系列的入射超声波是发射基(i);以及
-生成在作为输入的所述发射基(i)与作为输出的接收基(u)之间定义的实验反射矩阵Rui(t)的步骤;
-确定聚焦反射矩阵RFoc(rin,rout,δt)的步骤,所述聚焦反射矩阵包括在空间位置rin的输入虚拟换能器(TVin)与空间位置rout的输出虚拟换能器(TVout)之间的所述介质的响应,在相对于所述输入虚拟换能器(TVin)的响应的时刻移位额外延迟δt的时刻获得所述输出虚拟换能器(TVout)的响应,
-针对所述输入虚拟换能器(TVin)以及针对额外延迟间隔确定波前图像的步骤,所述波前图像根据所述介质中的声速c0确定,并且所述波前图像基于以下项确定:
-所述聚焦反射矩阵RFoc(rin,rout,δt),以及
-以下类型的弹道传播关系
δt(Δrout)=-sign(Δzout)·|Δrout|/c0,其可以从所述聚焦反射矩阵提取值以便构建所述波前图像,并且其中:
δt是所述额外延迟,
|Δrout|是所述输入虚拟换能器(TVin)与所述输出虚拟换能器(TVout)之间的向量的模数,其中Δrout=rout-rin
Δzout是沿着空间位置向量Δrout的深度轴Z的分量,
-确定所述波前图像中的焦斑的中心的深度位置Δz(0)(rin)的步骤,以及
-基于以下公式计算积分声速c(1)(rin)的步骤:
Figure FDA0003261953790000011
其中zin是沿着所述输入虚拟换能器(TVin)的空间位置向量rin的所述深度轴Z的分量。
2.根据权利要求1所述的方法,其特征在于,通过在所述波前图像中搜索具有最大值的点的空间位置来确定所述焦斑的中心。
3.根据权利要求1和2中任一项所述的方法,其特征在于,仅对所述深度轴Z执行所述波前图像的确定。
4.根据权利要求1至3中任一项所述的方法,其特征在于:
-在确定波前图像的步骤与确定焦斑的深度位置Δz(0)(rin)的步骤之间,执行改进所述波前图像的步骤,其中执行对应于给定相干区域ZC的一组波前图像的线性组合,在不同空间位置rin的选定输入虚拟换能器TVin与空间位置rout的输出虚拟换能器TVout之间获得所述集合的每个波前图像,使得rout=Δrout+rin,其中Δrout预定义且针对所述集合的所有波前图像相同,并且所述选定输入虚拟换能器彼此靠近,以便获得与参考输入虚拟换能器(TVin,ref)相关联的改进的波前图像,此参考输入虚拟换能器TVin,ref表示所使用的所述波前图像集合的所述输入虚拟换能器并且与所述所选择的相干区域ZC相关联,以及
-在确定深度位置Δz(0)(rin)的步骤中,代替所述波前图像使用所述改进的波前图像,所述焦斑的中心的深度位置与所述参考输入虚拟换能器TVin,ref的空间位置相关,并且所述焦斑的中心的此深度位置可以估计在所述参考输入虚拟换能器TVin,ref的空间位置处的积分声速c(1)(rin,ref)。
5.根据权利要求4所述的方法,其特征在于,通过计算所述波前图像集合的奇异值分解(SVD)以获得与所述奇异值分解的最大绝对值的奇异值相关联的奇异向量(W1)来确定所述线性组合,此奇异向量(W1)随后是对应于所述参考输入虚拟换能器(TVin,ref)且针对所述相同额外延迟δt的所述改进的波前图像。
6.根据权利要求4至5中任一项所述的方法,其特征在于,通过计算积分声速以及通过对于所述波前图像集合的所述线性组合,使用基本上覆盖所述介质中的整个所关注区域的与所述选定输入虚拟换能器(TVin)相对应的一组波前图像来确定所述介质的最佳声速。
7.根据权利要求1至6中任一项所述的方法,其特征在于:
-代替先前使用的声速或代替在第一步骤中使用的所述声速c0,使用在前一迭代中确定的所述积分声速c(n)迭代确定聚焦反射矩阵RFoc(rin,rout,δt)、确定波前图像,确定所述波前图像中的所述焦斑的中心的深度位置Δz(n)的步骤,以及
-在计算积分声速的步骤期间,使用以下递归公式,
Figure FDA0003261953790000021
以及
其中在对应于所述输入虚拟换能器的所述介质的点处的积分声速值是在所述方法的步骤n中计算出的所述积分声速c(n)(rin),此步骤n通过预定迭代次数或通过用于两个连续积分声速值之间的差或两者的组合的停止阈值确定。
8.根据权利要求1至7中任一项所述的方法,其特征在于:
-所述输入虚拟换能器和所述输出虚拟换能器的角色反转,以确定相对于输出虚拟换能器的积分声速c(1)(rout),以及
-将参考所述输入虚拟换能器的所述积分声速c(1)(rin)与参考所述输出虚拟换能器的所述积分声速c(1)(rout)组合以获得改进的积分声速。
9.根据权利要求1至8中任一项所述的方法,其进一步包括通过确定所述介质中的多个点的积分声速来确定积分声速图像的步骤,每个点对应于空间位置rin的输入虚拟换能器(TVin)或空间位置rin,ref的参考输入虚拟换能器(TVin,ref)。
10.根据权利要求9所述的方法,其特征在于,确定积分声速图像,其中从所述积分声速图像的值计算此声速图像的每个点处的值。
11.根据权利要求1至10中任一项所述的方法,其特征在于,在确定所述聚焦反射矩阵的步骤中:
所述输入虚拟换能器(TVin)的响应的计算对应于基于所述实验反射矩阵Rui(t)在输入处的聚焦过程,所述实验反射矩阵使用所述发射基与所述输入虚拟换能器(TVin)之间的所述波的向外飞行时间以在空间位置rin处产生输入焦斑,
所述输出虚拟换能器(TVout)的响应的计算对应于基于所述实验反射矩阵Rui(t)在输出处的聚焦过程,所述实验反射矩阵使用所述输出虚拟换能器(TVout)与所述接收基u的所述换能器之间的所述波的返回飞行时间以在空间位置rout处产生输出焦斑,
所述额外延迟δt是在所述聚焦过程期间添加到所述向外和返回飞行时间的时间滞后。
12.根据权利要求1至11中任一项所述的方法,其特征在于,通过以下公式计算所述聚焦反射矩阵:
Figure FDA0003261953790000031
其中
Nin是所述发射基(i)的元素的数目,
Nout是在输出处的所述接收基(u)的元素的数目,
Rui(t)是实验反射矩阵,其中Rui(uout,iin,τ(rin,rout,uout,iin,δt))是在所述发射基中的索引iin的发射之后且在时间τ处由空间位置uout的所述换能器记录的所述实验反射矩阵Rui(t)的所述元素,
τ是时间,所述时间是在所述发射基(i)的所述换能器与空间位置rin的所述输入虚拟换能器(TVin)之间的所述超声波的所述向外飞行时间τin,以及在空间位置rout的所述输出换能器(TVout)与所述接收基u的所述换能器之间的所述超声波的所述返回飞行时间τout,以及所述额外延迟δt的总和,如通过以下公式解释:
τ(rin,rout,uout,iin,δt)=τin(rin,iin)+τout(rout,uout)+δt 。
13.用于介质(20)的超声表征以确定所述介质中的局部积分声速值的系统(40),所述系统包括:
-换能器阵列(10),其适合于在所述介质的区域中生成一系列入射超声波并且记录通过所述区域反向散射的所述超声波作为时间的函数;以及
-计算单元(42),其连接到所述换能器阵列并且适合于实施根据前述权利要求中任一项所述的方法。
CN202111075417.0A 2020-09-15 2021-09-14 用于介质的超声表征的方法和系统 Pending CN114176639A (zh)

Applications Claiming Priority (2)

Application Number Priority Date Filing Date Title
FR2009312 2020-09-15
FR2009312A FR3114157B1 (fr) 2020-09-15 2020-09-15 Procédé et système de caractérisation ultrasonore d’un milieu

Publications (1)

Publication Number Publication Date
CN114176639A true CN114176639A (zh) 2022-03-15

Family

ID=74592039

Family Applications (1)

Application Number Title Priority Date Filing Date
CN202111075417.0A Pending CN114176639A (zh) 2020-09-15 2021-09-14 用于介质的超声表征的方法和系统

Country Status (8)

Country Link
US (1) US11776526B2 (zh)
EP (1) EP3967239A1 (zh)
JP (1) JP2022048982A (zh)
KR (1) KR20220036339A (zh)
CN (1) CN114176639A (zh)
AU (1) AU2021209266A1 (zh)
CA (1) CA3127900A1 (zh)
FR (1) FR3114157B1 (zh)

Families Citing this family (4)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
FR3114155B1 (fr) * 2020-09-15 2022-07-29 Supersonic Imagine Procédé et système de caractérisation ultrasonore d’un milieu
FR3114159A1 (fr) * 2020-09-15 2022-03-18 Supersonic Imagine Procédé et système de caractérisation ultrasonore d’un milieu
FR3114157B1 (fr) * 2020-09-15 2022-07-29 Supersonic Imagine Procédé et système de caractérisation ultrasonore d’un milieu
FR3114158B1 (fr) * 2020-09-15 2022-07-29 Supersonic Imagine Procédé et système de caractérisation ultrasonore d’un milieu

Family Cites Families (12)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
EP0256481B1 (de) 1986-08-20 1994-01-12 Siemens Aktiengesellschaft Verfahren und Einrichtung zur adaptiven Fokussierung bei einem medizinischen Ultraschall-Bildgabegerät
US5551433A (en) 1994-08-05 1996-09-03 Acuson Corporation Method and apparatus for a geometric aberration transform in an adaptive focusing ultrasound beamformer system
FR2815717B1 (fr) 2000-10-20 2003-01-10 Centre Nat Rech Scient Procede et dispositif non invasif de focalisation d'ondes acoustiques
FR2883982B1 (fr) * 2005-04-05 2009-05-29 Centre Nat Rech Scient Procede et dispositif d'imagerie utilisant des ondes de cisaillement
WO2008124923A1 (en) * 2007-04-13 2008-10-23 Centre Hospitalier De L'universite De Montreal Method and system of ultrasound scatterer characterization
FR2946753B1 (fr) * 2009-06-11 2011-07-22 Centre Nat Rech Scient Procede et dispositif ultrasonores pour caracteriser un milieu
US8915852B2 (en) * 2012-07-06 2014-12-23 Centre Hospitalier De L'universite De Montreal System and method for ultrasound scatterer characterization
CA2971676A1 (en) * 2014-12-24 2016-06-30 Super Sonic Imagine Shear wave elastography method and apparatus for imaging an anisotropic medium
FR3084166B1 (fr) * 2018-07-19 2020-10-16 Centre Nat Rech Scient Procedes et systemes de caracterisation ultrasonore non invasive d'un milieu heterogene
FR3114157B1 (fr) * 2020-09-15 2022-07-29 Supersonic Imagine Procédé et système de caractérisation ultrasonore d’un milieu
FR3114158B1 (fr) * 2020-09-15 2022-07-29 Supersonic Imagine Procédé et système de caractérisation ultrasonore d’un milieu
EP4006579A1 (en) * 2020-11-25 2022-06-01 Supersonic Imagine Method and system for compensating depth-dependent attenuation in ultrasonic signal data

Also Published As

Publication number Publication date
AU2021209266A1 (en) 2022-03-31
US20220084496A1 (en) 2022-03-17
FR3114157A1 (fr) 2022-03-18
EP3967239A1 (fr) 2022-03-16
KR20220036339A (ko) 2022-03-22
JP2022048982A (ja) 2022-03-28
FR3114157B1 (fr) 2022-07-29
US11776526B2 (en) 2023-10-03
CA3127900A1 (en) 2022-03-15

Similar Documents

Publication Publication Date Title
Vignon et al. Capon beamforming in medical ultrasound imaging with focused beams
Deffieux et al. The variance of quantitative estimates in shear wave imaging: theory and experiments
US11346819B2 (en) Methods and systems for non-invasively characterizing a heterogeneous medium using ultrasound
CN114176639A (zh) 用于介质的超声表征的方法和系统
US11761928B2 (en) Method and system for ultrasonic characterization of a medium
US11768181B2 (en) Method and system for ultrasonic characterization of a medium
EP3111206A1 (en) Methods and systems for measuring properties with ultrasound
Lambert et al. Ultrasound Matrix Imaging—Part I: The focused reflection matrix, the F-factor and the role of multiple scattering
US11275006B2 (en) Ultrasound estimation of nonlinear bulk elasticity of materials
US20220082693A1 (en) Method and system for ultrasonic characterization of a medium
US20240036004A1 (en) Method and system for ultrasonic characterization of a medium
US20240032889A1 (en) Method and system for non-invasively characterising a heterogeneous medium using ultrasound
Churchin History, fame, legacy: Malcolm Glazer: the life story and remembrance of Malcolm Glazer
Verma Beamforming strategies for plane wave vascular elastography

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