Sound velocity non-uniform spiral wave spectrum passive cavitation imaging method, device, equipment and medium

Abstract

The application provides a sound speed non-uniform spiral wave spectrum passive cavitation imaging method, device, equipment and medium, wherein the method comprises the following steps: a convex array ultrasonic transducer passively receives an ultrasonic time domain signal p(r, φ, t) in a sound speed non-uniform medium, the elements of the convex array ultrasonic transducer are distributed on an arc of a circle with a radius R, wherein r is the extrapolated radius of the current layer, φ is the circumferential angle in the coordinate system, t is time, and R is the distribution plane of the elements of the convex array ultrasonic transducer; based on the ultrasonic time domain signal p(r, φ, t), the sound speed distribution of the propagation medium is determined, and the aligned sound speed map is generated by aligning the probe of the convex array ultrasonic transducer; the ultrasonic time domain signal p(r, φ, t) is signal converted in time t and the circumferential direction φ to obtain a time frequency domain wave spectrum and a spatial frequency domain wave spectrum P n (r, ω), wherein ω represents frequency, and n is the spatial frequency corresponding to the Fourier transform with φ as the independent variable; based on the aligned sound speed map, the passive cavitation image I(r, φ) in the coordinate system is obtained by selecting the sound field intensity superposition in a specific frequency band Ω.

Description

Technical Field

[0001] This invention belongs to the field of ultrasonic imaging technology, and in particular relates to a passive cavitation imaging method, device, equipment and medium for sound velocity non-uniform spiral spectrum. Background Technology

[0002] Focused ultrasound therapy (UHEP), based on the acoustic cavitation effect, can improve the permeability of biological barriers, deliver therapeutic genes or drugs to lesions, and produce the desired biological effects to achieve therapeutic goals. It holds significant application potential in oncology and central nervous system diseases. However, acoustic cavitation activity is dynamically influenced by factors such as sound pressure and cavitation threshold in the treatment area, making it difficult to predict and control. Therefore, UHEP carries the risk of off-target effects, requiring efficient methods to monitor microbubble cavitation in real time and provide precise localization to ensure the safety and reliability of the UHEP procedure. Summary of the Invention

[0003] In view of the shortcomings of the prior art described above, the purpose of this invention is to provide a passive cavitation imaging method, apparatus, device and medium for sound velocity non-uniform spiral spectrum, in order to solve the above problems.

[0004] The passive cavitation imaging method for non-uniform sound velocity spiral spectrum provided by the present invention includes:

[0005] A convex array ultrasonic transducer passively receives ultrasonic time-domain signals p(r,φ,t) in a non-uniform sound velocity medium. The array elements of the convex array ultrasonic transducer are distributed on an arc of radius R, where r is the extrapolation radius of the current layer, φ is the circumferential angle in the coordinate system, t is time, and R is the distribution surface of the array elements of the convex array ultrasonic transducer.

[0006] Based on the ultrasonic time-domain signal p(r,φ,t), the sound velocity distribution in the propagation medium is determined and aligned with the probe of the convex array ultrasonic transducer to generate an aligned sound velocity map.

[0007] The ultrasonic time-domain signal p(r,φ,t) is transformed in time t and circumferential direction φ to obtain the time-frequency domain spectrum. and spatial frequency domain spectrum P n (r,ω), where ω represents the frequency and n is the spatial frequency corresponding to the Fourier transform with φ as the independent variable;

[0008] Based on the aligned sound velocity map, and by superimposing the sound field intensity within a specific frequency band Ω, a passive cavitation image I(r,φ) in the coordinate system is obtained.

[0009] In one embodiment of the present invention, the ultrasonic time-domain signal is converted in both time and circumferential directions, including:

[0010] Perform a Fourier transform on the ultrasonic time-domain signal p(r,φ,t) with respect to time t to convert the ultrasonic time-domain signal p(r,φ,t) into a time-frequency domain spectrum.

[0011] Time-frequency domain spectrum Perform a Fourier transform along the circumferential φ to convert the time-frequency spectrum. Converted to spatial frequency domain spectrum P n (r,ω).

[0012] In one embodiment of the present invention, a passive cavitation image I(r,φ) in a coordinate system is obtained by superimposing the sound field intensity within a specific frequency band Ω based on the aligned sound velocity map, including:

[0013] Based on the aligned sound velocity map, the slowness u in the medium is determined;

[0014] If r = r m The corresponding spatial frequency domain spectrum P n (r m ,ω), probe distribution surface r m+1 =r m +Δr, [r] m ,r m+1 The medium slowness u within the interval is decomposed into the interval average slowness u0 and the deviation Δu between the slowness in the interval and the average slowness.

[0015] Based on the average slowness u0 of the interval and the deviation Δu between the slowness in the interval and the average slowness, the sound field intensity in a specific frequency band Ω is selected and superimposed to obtain the passive cavitation image I(r,φ) in the coordinate system.

[0016] In one embodiment of the present invention, a passive cavitation image I(r,φ) in a coordinate system is obtained by superimposing the sound field intensity within a specific frequency band Ω based on the aligned sound velocity map, including:

[0017] Calculation of spiral spectrum using first-order direct integration method For the helical spectrum P n (r m+1 Perform an inverse Fourier transform on ω to obtain the frequency domain sound field.

[0018] or

[0019] Based on the first-order direct integration method, r m The distributed surface wave field propagates to r m+1 Distribution surface, determine r m+1 Spiral spectrum P of the distribution surface n (r m+1 ,ω), for the spiral spectrum P n (r m+1Perform an inverse Fourier transform on ω to obtain the frequency domain sound field.

[0020] or

[0021] Determining the frequency domain sound field using the step-by-step spiral spectroscopy method

[0022] In one embodiment of the present invention, the spiral spectrum is calculated by the first-order direct integration method. include:

[0023]

[0024] in,

[0025] k0 = ωu0;

[0026]

[0027] This is a Hankel equation of the second kind, order n.

[0028] This is a Hankel equation of the first kind, order n.

[0029] In one embodiment of the present invention, based on the first-order direct integration method, r m The distributed surface wave field propagates to r m+1 Distribution surface, determine r m+1 Spiral spectrum P of the distribution surface n (r m+1 ,ω), including:

[0030]

[0031] in,

[0032]

[0033] In one embodiment of the present invention, the frequency domain sound field is determined by the stepwise spiral spectral method. include:

[0034]

[0035] The passive cavitation imaging device for non-uniform sound velocity spiral spectrum provided by the present invention includes:

[0036] A convex array ultrasonic transducer is used to passively receive ultrasonic time-domain signals p(r,φ,t) in a non-uniform sound velocity medium. The array elements of the convex array ultrasonic transducer are distributed on an arc of radius R, where r is the extrapolation radius of the current layer, φ is the circumferential angle in the coordinate system, t is time, and R is the distribution surface of the array elements of the convex array ultrasonic transducer.

[0037] The aligned sound velocity map generation module is used to determine the sound velocity distribution of the propagation medium based on the ultrasonic time-domain signal p(r,φ,t), and align it with the probe of the convex array ultrasonic transducer to generate an aligned sound velocity map.

[0038] The processing module is used to perform signal conversion on the ultrasonic time-domain signal p(r,φ,t) in time t and circumferential φ to obtain the time-frequency domain spectrum. and spatial frequency domain spectrum P n (r,ω), where ω represents the frequency;

[0039] An imaging module is used to obtain a passive cavitation image I(r,φ) in a coordinate system by superimposing the sound field intensity within a specific frequency band Ω based on the aligned sound velocity map.

[0040] The electronic device provided by the present invention includes:

[0041] One or more processors;

[0042] A storage device for storing one or more programs, which, when executed by one or more processors, enable the electronic device to implement the passive cavitation imaging method for non-uniform helical spectra of sound speeds.

[0043] The present invention provides a computer-readable storage medium storing a computer program thereon, which, when executed by a computer processor, causes the computer to perform the aforementioned passive cavitation imaging method for non-uniform sound velocity spiral spectrum.

[0044] The beneficial effects of this invention are as follows: Compared with existing spiral spectroscopy methods, the passive cavitation imaging method of the sound velocity non-uniform spiral spectrum (i.e., the sound velocity heterogeneous spiral spectrum method) can accurately locate the sound source in a sound velocity non-uniform medium and more effectively suppress distortion errors and artifacts. Compared with the time delay accumulation integration method, it has the advantages of low algorithm complexity and low computation time.

[0045] It should be understood that the above general description and the following detailed description are exemplary and explanatory only, and do not limit this application. Attached Figure Description

[0046] The accompanying drawings, which are incorporated in and form part of this specification, illustrate embodiments consistent with this application and, together with the description, serve to explain the principles of this application. It is obvious that the drawings described below are merely some embodiments of this application, and those skilled in the art can obtain other drawings based on these drawings without any inventive effort. In the drawings:

[0047] Figure 1 This is a reference coordinate system shown in an exemplary embodiment of this application.

[0048] Figure 2 This is a flowchart illustrating a passive cavitation imaging method for sound velocity non-uniform spiral spectrum, as shown in an exemplary embodiment of this application.

[0049] Figure 3 This is a schematic diagram of the coordinate system for imaging a convex array ultrasonic transducer in a non-uniform medium, illustrating an exemplary embodiment of this application.

[0050] Figure 4 This is an exemplary embodiment of the present application illustrating an experimental diagram showing a sound source distributed in the shape of blood vessels.

[0051] Figure 5 This is an exemplary embodiment of the present application illustrating an experimental diagram of a sound source distributed in a grid pattern.

[0052] Figure 6 This is an exemplary embodiment of the present application illustrating the use of a translation stage to move the focusing probe to acquire data at different locations.

[0053] Figure 7 This is an exemplary embodiment of the present application illustrating a comparison between HHWS image results corrected for distortion using sound speed information and HWS results using average sound speed.

[0054] Figure 8 This is a PAM image of HWS and HHWS from a multi-cavitation source experiment, as illustrated in an exemplary embodiment of this application.

[0055] Figure 9 This is an exemplary embodiment of the present application, showing a comparison between HHWS image results using the exact Hankel function and HWS results using the average speed of sound without corrected distortion.

[0056] Figure 10 This is an exemplary embodiment of the present application showing a comparison of HHWS image results using the exact Hankel function and the approximate function. Detailed Implementation

[0057] The embodiments of the present invention will be described below with reference to the accompanying drawings and preferred embodiments. Those skilled in the art can easily understand other advantages and effects of the present invention from the content disclosed in this specification. The present invention can also be implemented or applied through other different specific embodiments, and various details in this specification can also be modified or changed based on different viewpoints and applications without departing from the spirit of the present invention. It should be understood that the preferred embodiments are only for illustrating the present invention and not for limiting the scope of protection of the present invention.

[0058] It should be noted that the illustrations provided in the following embodiments are only schematic representations of the basic concept of the present invention. Therefore, the drawings only show the components related to the present invention and are not drawn according to the actual number, shape and size of the components in the actual implementation. In the actual implementation, the form, quantity and proportion of each component can be arbitrarily changed, and the layout of the components may also be more complex.

[0059] It is worth noting that passive acoustic mapping (PAM) reconstructs images by passively receiving acoustic cavitation signals using an ultrasound transducer, visually monitoring cavitation activity generated during ultrasound therapy. Abdominal organs are an important application scenario for ultrasound cavitation therapy, and abdominal ultrasound imaging mostly employs convex arrays to obtain a wider imaging field. However, the types of tissues within the abdominal cavity are complex; the abdominal wall contains heterogeneous acoustic components such as skin, fat, muscle, and other organs, and the distortion of sound waves propagating within the abdominal cavity affects the image reconstruction quality. Currently, the main method suitable for convex arrays in non-uniform acoustic media is the passive cavitation imaging algorithm based on "time-delay accumulation integration." In non-uniform acoustic media, the time-delay accumulation integration method calculates the time delay by integrating along the sound propagation path, adjusts the waveform based on the time delay, and then coherently accumulates the adjusted signal and calculates the signal energy as image intensity information. This method has high algorithm complexity and is computationally time-consuming. Frequency-wavenumber domain methods, such as the angular spectral method, have lower computational complexity. In Cartesian coordinates, the heterogeneous angular spectral method has been applied to linear transcranial passive cavitation imaging to overcome the distortion problem caused by the heterogeneity of sound velocity. However, the heterogeneous angular spectral method cannot be directly applied to convex arrays. In recent years, the Helical Wave Spectrum Method (HWS) has been proposed as an extension of the angular spectral method in cylindrical coordinates. It is suitable for passive cavitation imaging based on convex arrays in homogeneous media, but not for media with non-homogeneous sound velocity.

[0060] In the following description, numerous details are explored to provide a more thorough explanation of embodiments of the invention. However, it will be apparent to those skilled in the art that embodiments of the invention may be practiced without these specific details. In other embodiments, well-known structures and devices are shown in block diagram form rather than in detail to avoid obscuring embodiments of the invention.

[0061] For example, the array elements of a convex array ultrasonic transducer can be distributed on a cylindrical or circumferential surface. For a convex array probe with one-dimensional array elements, y is always 0, the distribution surface is circumferential, and the passive cavitation imaging method with non-uniform sound velocity spiral spectrum can be regarded as operating in polar coordinates. Below is the establishment (see attached document) Figure 1 (As shown). For a convex array probe (i.e., a probe of a convex array ultrasonic transducer) with array elements arranged in two dimensions on a cylindrical surface, the y-value is variable and the distribution surface is cylindrical, which can be regarded as being established in a cylindrical coordinate system.

[0062] This application uses a one-dimensional array element arrangement as an example for illustration. The algorithm flow of the proposed convex array acoustic velocity heterogeneous passive cavitation imaging method is as follows: Figure 2 As shown, given the sound field distribution on the circumferential surface r, the sound field on the circumferential surface r+Δr is calculated, where ω is the angular frequency corresponding to the Fourier transform with respect to time t.

[0063] Step S100, data acquisition.

[0064] Specifically, see the attached document. Figure 3 As shown, to monitor acoustic cavitation activity, a convex array ultrasonic transducer is used to passively receive the ultrasonic time-domain signal p(r,φ,t) in a non-uniform sound velocity medium. The array elements of the convex array ultrasonic transducer are distributed on an arc of radius R. Here, r is the extrapolation radius of the current layer when acquiring images at different imaging depths, φ is the circumferential angle in polar coordinates, and t is time.

[0065] Step S200: Generate an aligned sound velocity map.

[0066] Specifically, an ultrasound image is generated based on the ultrasound time-domain signal p(r,φ,t). The sound velocity distribution in the propagation medium is obtained from the ultrasound image or sound velocity reconstruction, and aligned with the position of the ultrasound probe (referring to the probe of a convex array ultrasound transducer, or simply an ultrasound probe) to obtain an aligned sound velocity map. The sound velocity at each sampling point is obtained based on the aligned sound velocity map to determine the slowness of the sound wave in the medium.

[0067] Step S300: Sound field distribution calculation.

[0068] Step S310: Perform a Fourier transform on the ultrasonic time-domain signal p(r,φ,t) received by the ultrasonic probe with respect to time t. To convert the ultrasonic time-domain signal p(r,φ,t) into a time-frequency domain spectrum. Recorded as:

[0069]

[0070] Step S320, perform time-frequency domain spectrum analysis. Perform a Fourier transform along the circumferential φ to convert the time-frequency spectrum. Converted to spatial frequency domain spectrum P n (r,ω), the process is represented as:

[0071]

[0072] Where n is the spatial frequency corresponding to the Fourier transform with φ as the independent variable, i.e., the circumferential wave number.

[0073] For example, suppose the circumferential surface r = r m The corresponding spatial frequency domain spectrum is P n (r m Since the convex array can only receive sound waves emitted inward from the cavitation source located on the outer side, we only consider the inward incoming waves. From the aligned sound velocity diagram, we know that the slowness in the medium is u(l) = 1 / c(l), where c(l) is the sound velocity in the medium at spatial coordinates l = (r, φ). Let r... m+1 =r m +Δr, for [r] m ,r m+1 The interval is defined as follows: the slowness of the medium within the interval is decomposed into two parts, u0 and Δu, denoted as u(r,φ)=u0+Δu(r,φ), where u0 is the average slowness of the interval (constant), which is independent of the position coordinates φ and r, and Δu is the deviation of the slowness in the interval from the average slowness.

[0074] For example, the acoustic velocity heterogeneous spiral spectral method (i.e., the acoustic velocity non-uniform spiral spectral passive cavitation imaging method) applicable to convex arrays can be divided into three types according to different sound field propagation calculation processes: first-order and second-order direct integration methods, and step-by-step spiral spectral methods. Among them, the step-by-step spiral spectral method has the simplest algorithm and the highest robustness, but it has a slight error compared to the direct integration method. Therefore, it is recommended to use the step-by-step spiral spectral method for initial verification of imaging effects, and then replace it with the direct integration method when it is necessary to improve the imaging effect. The second-order direct integration method should be used when the axial interval Δr < λ / 2, and the first-order direct integration method should be used when Δr < λ / 2, where λ is the wavelength of the sound wave.

[0075] Method 1: The first-order direct integration method uses the first-order forward Euler method for calculation, as follows:

[0076]

[0077] in,

[0078] k0=ωu0, This is a Hankel equation of the second kind, nth order. This is a Hankel equation of the first kind, order n.

[0079] Method 2: If the result of the first-order direct integration method will be used in subsequent calculations, for ease of distinction, the result of the first-order direct integration method will be denoted as... The calculation can be performed using the second-order direct integration method.

[0080] The second-order direct integration method, based on the first-order direct integration method, uses a second-order prediction-correction method to convert r... m Peripheral wave field propagates to r m+1 For the perimeter, the trapezoidal numerical integration method is used to calculate r. m+1 Circumferential spiral spectrum P n (r m+1 The specific process is as follows:

[0081]

[0082] in,

[0083] This is the inverse Fourier transform.

[0084] For the spiral spectrum P obtained by the first-order direct integration method or the second-order direct integration method n (r m+1 Perform an inverse Fourier transform on the frequency domain sound field (ω) to obtain the sound field in the frequency domain. Represented as:

[0085]

[0086] Method 3: Determining the frequency domain sound field using the step-by-step spiral spectral method. Its calculation is expressed as:

[0087]

[0088] For the three methods mentioned above, the Hankel equation in the first-order direct integration method, the second-order direct integration method, or the stepwise spiral spectroscopy method can be simplified by using asymptotic analytical solutions (i.e., using asymptotic formulas), specifically including:

[0089]

[0090] in,

[0091] for . conjugate.

[0092] The calculation involving Hankel's equations simplifies to:

[0093]

[0094] in,

[0095]

[0096] For the non-evanescent wave portion, the asymptotic formula can accurately approximate the wave, reducing computational costs.

[0097] Starting from the circumferential surface of the ultrasonic probe (r = R), the spectral distribution of the entire imaging space is obtained by iterative calculation using the above method. The sound field intensity within a specific frequency band Ω is then superimposed to obtain the passive cavitation image I(r,φ) in polar coordinates, as shown below:

[0098]

[0099] Step S400, Coordinate System Transformation

[0100] Specifically, an interpolation algorithm is used to transform the polar coordinate sampling points into an image in a Cartesian coordinate system for easier observation and evaluation. The origin is the intersection of the arc where the probe of the convex array ultrasonic transducer is located and the axis of symmetry, with the axial direction as the x-axis and the transverse direction as the z-axis.

[0101] This embodiment derives a sound velocity heterogeneous spiral wave spectrum propagation method in the frequency-wavenumber domain, starting from the weakly heterogeneous wave equation of sound velocity in polar coordinates, and applies it to passive cavitation imaging of convex arrays. The sound velocity heterogeneous spiral wave spectrum method (HHWS) can be divided into first-order direct integration, second-order direct integration, and step-by-step spiral wave spectrum methods, depending on the different sound field propagation calculation processes. Compared with the spiral wave spectrum method, the sound velocity heterogeneous spiral wave spectrum passive cavitation imaging method can accurately locate the sound source in a non-uniform sound velocity medium and more effectively suppress distortion errors and artifacts. Compared with the time-delay accumulation integration method, it has the advantages of lower algorithm complexity and lower computation time. The method proposed in this embodiment can correct the phase distortion caused by sound velocity heterogeneity, providing key technical support for real-time, accurate monitoring and guidance of abdominal ultrasound cavitation therapy. The sound velocity heterogeneous spiral wave propagation method derived in this embodiment is a general solution to the sound velocity weak heterogeneity wave equation. As a fast solution method for the circumferential wave field, it can also be used to reconstruct B-mode ultrasound images or as a forward model for solving the circumferential sound velocity field.

[0102] Experiment 1: Simulation Experiment

[0103] Step S100: Data Acquisition

[0104] A 128-channel convex array ultrasonic transducer (aperture radius R = 49.6 mm, scanning angle 74.6°) was simulated using k-Wave technology to receive point source acoustic signals in a non-uniform sound velocity medium. The signal sampling frequency was 13.89 MHz. The acoustic characteristics of the medium were simulated to resemble human soft tissue (1540 m / s, 1050 kg / m³). 3 0.54dB / cm / MHz), fat (1430m / s, 928kg / m 3 (0.60 dB / cm / MHz) and human abdominal wall muscles (1580 m / s, 1041 kg / m 3 (0.57 dB / cm / MHz). In soft tissue, 200 point sources were randomly placed in a vascular-shaped region approximately 50 mm deep from the convex array (Fig. 4(a)). In addition, 9×13 point sources were placed in a rectangular grid with lateral and axial positions ranging from -30 to 30 mm and 40 to 80 mm, respectively (Fig. 5(a)), and these point sources sequentially emitted a 100-cycle sine wave signal at 2 MHz.

[0105] Step S200: Obtain the aligned sound velocity map

[0106] The sound velocity diagram shows the sound velocity distribution set in the simulation.

[0107] Step S300: Sound field distribution calculation

[0108] Taking the second-order direct integration method as an example, according to step S310, the received microbubble cavitation signal is subjected to a Fourier transform in the time t direction. Then, according to step S320, a Fourier transform is performed in the circumferential direction to obtain the spatial frequency domain spectrum P. n (r,ω). At depths from r = R to R + 100 mm, the spatial frequency domain spectrum P is derived by iteratively applying first-order and second-order direct integration methods at intervals of Δr = 0.2 mm. n (r,ω). Inverse Fourier transform to time-frequency domain spectrum. Then, the sound field intensity of the single frequency component corresponding to the selected sound source is superimposed to obtain the passive cavitation image in polar coordinates, with a sampling grid of 384×501(φ×r).

[0109] Step S400: Coordinate system transformation

[0110] After converting to a Cartesian coordinate system, the location of the simulated sound source is used as the positioning reference (l r Compare the HHWS image results using sound velocity information to correct distortion with those using average sound velocity. For the HWS method, the localization errors of 2MHz, 3MHz, and 4MHz signal sources (localization error at the sound intensity peak position l: |ll) are also compared. rThe actual values ​​were 4.8±2.1mm, 4.2±2.0mm, and 3.4±1.3mm, respectively. After using the HHWS method, the errors were reduced to 0.5±0.3mm, 0.7±0.3mm, and 0.8±0.3mm, respectively. After calibration, the positioning signal source and the shape of the blood vessel can be matched well. See the attached diagram for details. Figure 4 (b)-(c).

[0111] Appendix Figure 4 (a) is a schematic diagram of a simulation experiment showing the sound source distributed in the shape of blood vessels. (b) compares the localized (orange dots) and actual (gray dots) sound source locations. The color of the localized points is normalized according to the maximum intensity at that frequency (2MHz, 3MHz, and 4MHz). (c) is the PAM image of the location indicated by the arrow (the red cross marks the actual sound source location).

[0112] For experiments with sound sources distributed in a grid, compared to the HWS method, the average localization error of the HHWS method decreased from 5.6±3.0 mm to 0.7±0.5 mm. The average energy diffusion area (ESA, the area in the image with intensity greater than half of the maximum value) decreased from 20.9±18.0 mm. 2 Reduced to 8.8±4.0mm 2 See attached document Figure 5 (b)-(c), attached Figure 5 (a) Schematic diagram of the simulation experiment with the sound source distributed in a grid. (b) Localization error and (c) ESA data of the PAM image reconstructed using the HWS and HHWS methods.

[0113] Experiment 2: Single-source voidification experiment

[0114] Step S100: Data Acquisition

[0115] A 128-channel convex array ultrasonic transducer with a radius of curvature of 49.57 mm and an imaging field of view of 74.6° was used to receive acoustic cavitation signals generated by excited ultrasonic microbubbles for reconstructing passive cavitation images. The transducer's center frequency was 3.57 MHz, and the sampling frequency was 13.89 MHz. An elliptical cylindrical phantom with a sound velocity similar to that of muscle (1580 m / s) was placed below the convex array ultrasonic transducer in water (1493 m / s). A 3D-printed bracket was used to fix the relative position of a focused ultrasonic probe with a center frequency of 1 MHz and the rectangular wallless phantom, which was connected to a translation stage and placed below the elliptical phantom. The focused ultrasonic probe excited microbubbles in the channels at a repetition frequency of 30 Hz to generate cavitation activity. The focused sound pressure amplitude was set to 800 kPa, and the focused ultrasonic probe was moved using a translation stage to acquire data at 9 different positions centered at [0 mm, 60 mm]. Figure 6 ).

[0116] Step S200: Obtain the aligned sound velocity map

[0117] The elliptical phantom is fabricated using a 3D-printed mold, with a fixed semi-axis length. The sound velocity map is reconstructed by identifying the position of the elliptical phantom using B-mode image recognition. Figure 4 As shown.

[0118] Step S300: Sound field distribution calculation

[0119] Taking the second-order direct integration method as an example, according to step S310, the received microbubble cavitation signal is subjected to a Fourier transform in the time direction. According to step S320, a Fourier transform is performed in the circumferential direction to obtain the spatial frequency domain spectrum P. n (r,ω). At depths from r = R to R + 100 mm, the spectrum P is derived by iteratively applying the first-order direct integration method and the second-order direct integration method at intervals of Δr = 0.2 mm. n (r,ω). Inverse Fourier transform to time-frequency domain spectrum. Then, the sound field intensities of each frequency within the selected frequency band Ω are superimposed to obtain a passive cavitation image in polar coordinates, with a sampling grid of 384×501 (φ×r). A panoramic image is reconstructed by selecting all frequency components (1~6MHz) covering the working bandwidth of the ultrasound probe.

[0120] Step 4: Coordinate System Transformation

[0121] After transforming to a Cartesian coordinate system, the HWS method reconstruction result under a homogeneous medium (water) with the elliptical hypersonic phantom removed is used as the positioning reference (l r Data from sound source locations at (0.8mm, 48.8mm) and (-8.8mm, 58.4mm) are selected. The HHWS image results corrected for distortion using sound velocity information are compared and presented with the HWS results using average sound velocity. Figure 7 As shown. (Attached) Figure 7 Example images of a single cavitation source experiment: (a) PAM image of a sound source located at (0.8 mm, 48.8 mm). (b) PAM image of a sound source located at (-8.8 mm, 58.4 mm). The crosshairs indicate the reference position of the sound source.

[0122] The localization errors and ESA data for the nine sound source locations are listed in Table (1). In the uncorrected (HWS) case, the average localization error (localization error at the peak sound intensity location l: |ll) is... r The energy diffusion area (ESA, the region in the image with an intensity greater than half of the maximum value) ranges from 6.4 to 12.8 mm and from 32.3 to 95.5 mm, respectively. 2 Between. After correction (HHWS), these values ​​range from 0.8-1.9 mm and 10.3-22.8 mm, respectively. 2The image quality and positioning accuracy were improved after correction.

[0123] Table 1. Experimental positioning error and ESA data for single cavitation source

[0124]

[0125]

[0126] Experiment 3: Multi-voidification Source Experiment

[0127] Step S100: Data Acquisition

[0128] A convex array ultrasonic transducer (128 channels, radius of curvature 49.57 mm, imaging field of view 74.6°) was used to passively receive microbubble cavitation signals (i.e., ultrasonic time-domain signals) and reconstruct passive cavitation images. The transducer (short for convex array ultrasonic transducer) had a center frequency of 3.57 MHz and a sampling frequency of 13.89 MHz. Two elliptical phantoms with a sound velocity of 1580 m / s were placed below the transducer in water (23°, 1493 m / s). A focused ultrasonic probe with a center frequency of 1 MHz was fixed to a dual-channel rectangular wallless phantom under the elliptical phantoms, and microbubbles were excited in the channels at a repetition frequency of 30 Hz and a sound pressure of 400 kPa to generate cavitation activity.

[0129] Step S200: Obtain the aligned sound velocity map

[0130] Specifically, the sound velocity map is reconstructed by identifying the position of the elliptical phantom through B-mode image recognition.

[0131] Step S300: Sound field distribution calculation

[0132] Taking the second-order direct integration method as an example, according to step S310, the received microbubble cavitation signal is subjected to a Fourier transform in the time direction. According to step S320, a Fourier transform is performed in the circumferential direction to obtain the spatial frequency domain spectrum P. n (r,ω). At depths from 0 to 100 mm, the formulas from the first-order and second-order direct integration methods are iterated layer by layer at intervals of Δr = 0.2 mm to derive the spectral P. n (r,ω). Inverse Fourier transform to time-frequency domain spectrum. Then, the sound field intensities of each frequency within the selected frequency band Ω are superimposed to obtain a passive cavitation image in polar coordinates, with a sampling grid of 384×501 (φ×r). A panoramic image is reconstructed by selecting all frequency components (1~6MHz) covering the ultrasonic probe's frequency band.

[0133] Step 4: Coordinate System Transformation

[0134] After transforming to a Cartesian coordinate system, the HWS method reconstruction result under a homogeneous medium (water) with the elliptical hypersonic phantom removed is used as the positioning reference (l r The results of HHWS image distortion correction using sound velocity information are compared with those using average sound velocity, as shown in the figure. Figure 8 As shown. Compared to the uncorrected case, the positioning error (positioning at the peak sound intensity location l, positioning error: |ll) is... r The area of ​​incision decreased from 4.3 ± 2.3 mm (8.6 ± 9.1 mm) at the left (right) site to 1.6 ± 3.4 mm (1.9 ± 2.6 mm), and the ESA decreased from 47.0 ± 30.1 mm at the left (right) site. 2 Reduced to 14.7±8.5mm 2 (from 79.8±48.8mm) 2 Reduced to 20.0±16.0mm 2 ).

[0135] Experiment 4: Dual-modal Experiment

[0136] Step S100: Data Acquisition

[0137] A convex array ultrasonic transducer (128 channels, radius of curvature 49.57 mm, imaging field of view 74.6°) was used to passively receive microbubble cavitation signals and reconstruct passive cavitation images. The transducer's center frequency was 3.57 MHz, and the sampling frequency was 13.89 MHz. A fan-shaped, wallless phantom was embedded in a hypersonic elliptical phantom and placed below the transducer. A focused ultrasonic probe excited microbubbles in the channels with a repetition frequency of 30 Hz and a focused sound pressure amplitude of 800 kPa, generating cavitation activity.

[0138] Step S200: Obtain the aligned sound velocity map

[0139] The sound velocity map is reconstructed by identifying the position of the elliptical phantom using B-mode images.

[0140] Step S300: Sound field distribution calculation

[0141] Taking the second-order direct integration method as an example, according to step S310, the received microbubble cavitation signal is subjected to a Fourier transform in the time direction. According to step S320, a Fourier transform is performed in the circumferential direction to obtain the spatial frequency domain spectrum P. n (r,ω). Within the imaging range from r = R to r = R + 100 mm, iterate layer by layer at intervals of Δr = 0.2 mm. Derive the spatial frequency domain spectrum P using the formulas in the second-order direct integration method and the formulas in the second-order direct integration method after replacing the asymptotic formula. n (r,ω). Inverse Fourier transform to time-frequency domain spectrum. Then, the sound field intensity of each frequency (1-6MHz) within the selected frequency band Ω is superimposed to obtain the passive cavitation image in polar coordinates, with a sampling grid of 384×501(φ×r).

[0142] Step 4: Coordinate System Transformation

[0143] After converting to a Cartesian coordinate system, the HHWS image results using the exact Hankel function (Matlab's "besselh" function) are compared with the HWS results using the average sound velocity without distortion correction. (Example image follows) Figure 9 As shown. Compared to the uncorrected case, the passively cavitation image is correctly focused on the channel. A comparison of HHWS image results using the exact Hankel function (formula in the second-order direct integral method) and the approximate function (formula in the second-order direct integral method after replacing the asymptotic formula) is shown below. Figure 10 As shown. The difference between the PAM image results of the approximate function and the accurate Hankel function results is within 10%. -3 Order of magnitude. The HHWS method enables real-time imaging. For a frequency window of 1 (739), the computation time using the accurate Hankel function is 102.1 ± 3.1 ms (74.4 ± 1.2 s) on the CPU and 18.9 ± 0.5 ms (66.9 ± 0.7 ms) on the GPU. After using the approximate Hankel function, the computation time on the CPU and GPU for the 739 frequency window is reduced to 27.8 ± 0.4 s and 46.7 ± 0.8 ms, respectively.

[0144] Reference Appendix Figure 9 As shown, attached Figure 9 The PAM image and B-mode image are overlaid to reconstruct the sound source at (a) and (b) positions shifted 10 mm to the right using the accurate Hankel function. The orange arrows in the image indicate the flow direction of the microbubble (MB) suspension.

[0145] Reference Appendix Figure 10 As shown, attached Figure 10 PAM images reconstructed using (a) the exact Hankel function and (b) the approximation function. (c) The difference between PAM images reconstructed using the exact Hankel function and the approximation function. All images were evaluated by normalizing to the maximum value in the image formed using the exact Hankel function.

[0146] In view of the above-mentioned passive cavitation imaging method for non-uniform sound velocity spiral spectrum, this application also provides a passive cavitation imaging device for non-uniform sound velocity spiral spectrum, comprising:

[0147] A convex array ultrasonic transducer is used to passively receive ultrasonic time-domain signals p(r,φ,t) in a non-uniform sound velocity medium. The array elements of the convex array ultrasonic transducer are distributed on an arc of radius R, where r is the extrapolation radius of the current layer, φ is the circumferential angle in polar coordinates, and t is time.

[0148] The aligned sound velocity map generation module is used to determine the sound velocity distribution of the propagation medium based on the ultrasonic time domain signal p(r,φ,t), and align it with the ultrasonic probe of the convex array ultrasonic transducer to generate an aligned sound velocity map.

[0149] The processing module is used to perform signal conversion on the ultrasonic time-domain signal p(r,φ,t) in time t and circumferential φ to obtain the time-frequency domain spectrum. and spatial frequency domain spectrum P n (r,ω), where ω represents the frequency;

[0150] The imaging module is used to derive the spatial frequency domain spectrum P of the imaging space based on the aligned sound velocity map. n (r,ω), and the passive cavitation image I(r,φ) in polar coordinates is obtained by superimposing the sound field intensity within a specific frequency band Ω.

[0151] This application also provides an apparatus, comprising:

[0152] One or more processors and memory,

[0153] The memory stores a computer program, which, when executed by one or more processors, causes the device to perform the passive cavitation imaging method for non-uniform sound velocity spiral spectrum.

[0154] Another aspect of this application provides a computer-readable storage medium storing a computer program thereon, which, when executed by a computer's processor, causes the computer to perform the passive cavitation imaging method for non-uniform sound velocity spiral spectra as described above. This computer-readable storage medium may be included in the electronic device described in the above embodiments, or it may exist independently and not assembled into the electronic device.

[0155] Another aspect of this application provides a computer program product or computer program including computer instructions stored in a computer-readable storage medium. A processor of a computer device reads the computer instructions from the computer-readable storage medium and executes the computer instructions, causing the computer device to perform the passive cavitation imaging method for non-uniform helical spectra of sound speeds provided in the various embodiments above.

[0156] The above embodiments are merely illustrative of the principles and effects of the present invention and are not intended to limit the invention. Any person skilled in the art can modify or alter the above embodiments without departing from the spirit and scope of the present invention. Therefore, all equivalent modifications or alterations made by those skilled in the art without departing from the spirit and technical concept disclosed in the present invention should still be covered by the claims of the present invention.

Claims

1. A method of passive cavitation imaging with a non-uniform spectrum of acoustic velocity spirals, characterized in that, comprises: Passive receiving of ultrasonic time domain signals by a convex array ultrasonic transducer in a non-uniform medium , where r is the extrapolated radius of the current layer, φ is the circumferential angle in the coordinate system, t is time, is the distribution plane of the elements of the convex array ultrasonic transducer determining a propagation medium sound speed distribution based on the time domain signal , and aligning the probe with the convex array ultrasonic transducer to generate an aligned sound speed map; The ultrasonic time domain signal The signal is converted in time t and circumferentially φ to obtain a time frequency domain spectrum And a spatial frequency domain spectrum Where ω represents the frequency, The Fourier transform corresponding to the spatial frequency is taken as the independent variable The Fourier transform corresponding to the spatial frequency is taken as the independent variable Based on the aligned sound velocity map, and selecting a specific frequency band The sound field intensity in the coordinate system is superimposed to obtain a passive cavitation image ; Wherein, based on the alignment sound velocity map, and selecting a specific frequency band The sound field intensity superposition in the coordinate system obtains the passive cavitation image , comprising: determining a slowness in the medium based on the aligned velocity maps ; like The corresponding spatial frequency domain spectrum is Probe distribution surface ,Will Medium slowness within the range Decomposed into interval average slowness And the deviation between the slowness in the interval and the average slowness ; based on the interval average slowness and the deviation of the interval slowness from the average slowness , the sound field intensity in a certain frequency band is superimposed to obtain a passive cavitation image in the coordinate system ; Determining a frequency domain sound field , comprising: The helicon wave spectrum is calculated by a first order direct integration method The helicon wave spectrum is calculated by a first order direct integration method The inverse Fourier transform is applied to the frequency domain sound field ; or On the basis of the first-order direct integration method, the distributed surface wave field is propagated to the distributed surface, the spiral wave spectrum of the distributed surface is determined ; the inverse Fourier transform is performed on the spiral wave spectrum to obtain the frequency domain sound field ; or Determining a frequency domain acoustic field by step-spiral spectroscopy .

2. The acoustic velocity non-uniform spiral wave spectrum passive cavitation imaging method of claim 1, wherein, signal transforming the ultrasound time domain signal in time and circumferential direction, comprising: ultrasound time-domain signal a Fourier transform is performed with respect to time t to convert the ultrasound time-domain signal into a time-frequency domain spectrum ; time-frequency spectrum A Fourier transform is performed in the circumferential direction φ to convert the time-frequency spectrum into a spatial-frequency spectrum .

3. The acoustic velocity non-uniform spiral wave spectrum passive cavitation imaging method of claim 1, wherein, The helical wave spectrum is calculated by a first order direct integration method comprising: wherein, ， ； ； for the second category Hankel equation of order For the first class order Hankel equation.

4. The acoustic velocity non-uniform spiral wave spectrum passive cavitation imaging method of claim 1, wherein, On the basis of the first-order direct integration method, the distributed surface wave field propagates to the distributed surface, determines the spiral wave spectrum of the distributed surface , comprising: wherein, 。 5. The acoustic velocity non-uniform spiral wave spectrum passive cavitation imaging method of claim 1, wherein, Determining a frequency domain acoustic field by step-spiral spectroscopy comprising: 。 6. A passive cavitation imaging device of acoustic velocity non-uniform spiral wave spectrum, characterized in that, comprises: A convex array ultrasonic transducer for passive receiving ultrasonic time domain signals in a non-uniform medium with respect to sound velocity , the elements of the convex array ultrasonic transducer are distributed on an arc of a circle with radius R, where r is the extrapolated radius of the current layer, φ is the circumferential angle in the coordinate system, t is time, is the distribution plane of the elements of the convex array ultrasonic transducer; an aligned velocity map generation module configured to determine a propagation medium velocity distribution based on the ultrasound time domain signals and aligned with a probe of the convex array ultrasound transducer to generate an aligned velocity map; a processing module for transforming the ultrasonic time domain signal a signal transformation in time t and in the circumferential direction φ to obtain a time frequency domain spectrum and a spatial frequency domain spectrum where ω denotes the frequency; An imaging module is used to obtain passive cavitation images in a coordinate system based on the aligned sound velocity map and selected frequency band and the sound field intensity superposition in the inner part of the selected frequency band ; Wherein, wherein, based on the alignment sound velocity map, and select a specific frequency band The sound field intensity superposition in the selected frequency band is obtained to obtain the passive cavitation image in the coordinate system , comprising: determining a slowness in the medium based on the aligned velocity maps ; like The corresponding spatial frequency domain spectrum is Probe distribution surface ,Will Medium slowness within the range Decomposed into interval average slowness And the deviation between the slowness in the interval and the average slowness ; based on the interval average slowness and the deviation of the interval slowness from the average slowness , the sound field intensity in a certain frequency band is superimposed to obtain a passive cavitation image in the coordinate system ; Determining a frequency domain sound field , comprising: The helicon wave spectrum is calculated by a first order direct integration method The helicon wave spectrum is calculated by a first order direct integration method The inverse Fourier transform is performed to obtain the frequency domain acoustic field ; or On the basis of the first-order direct integration method, the distributed surface wave field is propagated to the distributed surface, the spiral wave spectrum of the distributed surface is determined ; the inverse Fourier transform is performed on the spiral wave spectrum , and the frequency domain sound field is obtained ; or Determining a frequency domain acoustic field by step-spiral spectroscopy .

7. An apparatus, comprising: comprises: one or more processors and a memory, the memory having stored thereon a computer program which, when executed by the one or more processors, causes the device to perform the method of any one of claims 1-5.

8. A computer-readable storage medium, characterized in that, having stored thereon a computer program which, when executed by one or more processors, causes a device to perform the method of any one of claims 1-5.

Images