A uniform medium shear wave velocity detection method and system based on frequency domain analysis

Abstract

The application provides a uniform medium shear wave velocity detection method and system based on frequency domain analysis, and belongs to the technical field of ultrasonic elastography, which comprises the following steps: (1) in a phantom, a shear wave is excited at a set focus depth to extract time domain information of shear wave propagation; (2) the time domain information in step (1) is converted into frequency domain information; (3) a straight line is detected in the frequency domain by using an improved radon algorithm; and (4) shear wave velocity is obtained according to parameters of the detected straight line. The method can greatly reduce noise in the peak displacement algorithm in the time domain, can more accurately reflect real data, and can obtain relatively accurate viscoelastic information of tissues and high-quality shear wave velocity images.

Description

Technical Field

[0001] This invention relates to the field of ultrasonic elastography technology, specifically to a method and system for detecting shear wave velocity in a uniform medium based on frequency domain analysis. Background Technology

[0002] Viscoelastic information of biological tissues can provide valuable reference information for the diagnosis of tissue function and the degree of lesions. For example, domestic patent 201880097002.3 discloses a shear wave elastography method, apparatus, and computer storage medium. This method involves emitting an ultrasonic beam towards a target to detect motion interference; generating a shear wave in the target and emitting a probe beam to detect the shear wave; receiving a first echo signal corresponding to the ultrasonic beam and a second echo signal corresponding to the probe beam; determining the displacement parameter corresponding to the motion interference based on the first echo signal, and determining the intensity parameter and signal-to-noise ratio corresponding to the shear wave based on the second echo signal, thereby determining the image quality of the shear wave elastography. Therefore, it can be seen that determining the shear wave characteristics is crucial for the imaging results in ultrasonic shear wave elastography.

[0003] Therefore, shear wave-based elastography often uses pulses to generate shear waves within the body, directly calculating the tissue's elasticity by measuring their propagation velocity. The propagation velocity of shear waves in a medium is related to the longitudinal modulus of the biological tissue; thus, the tissue's elastic modulus can be estimated by measuring its velocity. The core of ultrasonic shear wave elastography lies in calculating the wave velocity of the generated shear waves and reconstructing images using this data. A commonly used method for calculating shear wave velocity is the peak-time algorithm. The basic idea of ​​this algorithm is to calculate the shear wave velocity by detecting the time between the peak position of the impact acoustic radiation force generated outside the excitation region and the peak position after the shear wave is generated. This calculation method is robust even with low displacement signal-to-noise ratios.

[0004] However, the peak position and time data obtained using this algorithm are affected by noise in the time domain, which may interfere with the arrival time of the peak, or the data may be unstable when the peak is small, affecting the accuracy of the data and thus causing some errors in the calculated shear wave velocity. Summary of the Invention

[0005] To address the problems existing in the prior art, this invention provides a method and system for detecting shear wave velocity in a uniform medium based on frequency domain analysis. This method can improve upon the shortcomings of peak time algorithms, is less affected by environmental interference, has good linearity, and can improve the accuracy of shear wave velocity calculation, thereby facilitating the acquisition of high-quality shear wave velocity image information.

[0006] To achieve the above objectives, the present invention provides the following technical solution:

[0007] A method for detecting shear wave velocity in a uniform medium based on frequency domain analysis, comprising the following steps:

[0008] (1) In the phantom, a shear wave is generated at the set focusing depth, and the time domain information of the shear wave propagation is extracted.

[0009] (2) Convert the time domain information from step (1) into frequency domain information;

[0010] (3) Detecting straight lines in frequency domain information using the improved Radon algorithm;

[0011] (4) Based on the parameters of the detection line, the shear wave velocity is obtained.

[0012] Further specifying, step (1) specifically includes:

[0013] (1.1) In the phantom, acoustic radiation force is generated at the set focusing depth to generate shear wave propagation data U(z,x,t) in the entire phantom, where z represents the axial dimension variable, x represents the transverse dimension variable, and t represents the slow time dimension variable.

[0014] (1.2) Take the data u(x,t) at the depth z0 where the focal point is located, and obtain the spatiotemporal shear wave propagation image represented by the particle displacement signal, and obtain the time domain information of the spatiotemporal shear wave propagation image.

[0015] Further specifying, the phantom refers to a homogeneous medium that is isotropic, linearly elastic, and viscous.

[0016] Further specifying, step (2) specifically includes:

[0017] (2.1) The time-domain information from step (1) is subjected to a two-dimensional discrete Fourier transform using the FFT function to obtain the spectral information of the shear wave propagating in a homogeneous medium.

[0018] (2.2) The spectrum information is shifted so that the zero frequency is in the center of the frequency domain image;

[0019] (2.3) The amplitude of the spectrum is obtained using the abs function, and the frequency domain information and frequency domain information diagram of the shear wave propagating in the homogeneous medium are obtained, as shown below:

[0020] S(k x f) = FFT[FFT] x,t (u(x,t))]

[0021] Where, k x The space wave number is represented by 1 / m; f represents the frequency by Hz.

[0022] Further specifying, step (2.2) specifically involves using the fftshift function to shift the spectral information so that the zero frequency is located in the center of the frequency domain image.

[0023] Further specifying, step (3) specifically includes:

[0024] In the frequency domain information graph, the data corresponding to the x and y coordinates of the line y = ax + b in the XYD reference frame will be transformed into the intercept and slope parameter space. Below, using the cyclic displacement method, in the intercept slope parameter space... The lower detection line transforms y = ax + b into the intercept and slope parameter space. The corresponding point (a, b) is used to complete the improved Radon algorithm transformation of the frequency domain information.

[0025] Further specifying, the cyclic displacement method is as follows:

[0026] The intercept b of the line y = ax + b has M discretized values, representing b1 to b2. M The slope k has N discretized values, represented as k1 to k2. N ;

[0027] Establish parameter space Parameter space The quantization scales of the two axes are the same as those of the original space XYD, both being matrices of size M×N.

[0028] Parameter space The value of the data in the i-th row and j-th column is k i Let b represent the i-th value of k, where i = 1, 2, ..., N. j Let b be the j-th value, where j = 1, 2, ..., M; parameter space. The value of the j-th row and i-th column It should be: the sum of the values ​​corresponding to all two-dimensional coordinates (X,Y) that satisfy the equation Y=kX+b, thereby determining the parameter space. The i-th column can be obtained by cyclically shifting and accumulating each column of the parameter space XYD;

[0029] When k i X L When the value is less than 0, the data in the Lth column of the parameter space XYD is cyclically shifted upwards by -k. i X L Bit;

[0030] When k i X LWhen the value is greater than 0, the data in the Lth column of the parameter space XYD is shifted downwards by k in a circular manner. i X L Bit;

[0031] When k i X L When the value is 0, no cyclic displacement operation is performed;

[0032] Complete the cyclic displacement.

[0033] Further specifying, step (4) specifically includes:

[0034] After transforming the point (a,b) obtained by the improved Radon algorithm, we can convert it into the parameter space of energy accumulation. The point with the highest energy accumulation is then determined by the slope and intercept of the corresponding line. This can be calculated using the following formula:

[0035]

[0036] V g ω is the group velocity, in m / s; ω is the angular frequency, in rad / s; k is the space wave number, in rad / m.

[0037] The slope of the straight line is the shear wave velocity.

[0038] A system for detecting shear wave velocity in a homogeneous medium based on frequency domain analysis, which implements the aforementioned method for detecting shear wave velocity in a homogeneous medium based on frequency domain analysis; specifically including:

[0039] The shear wave excitation module is used to generate shear waves in the phantom at a set focusing depth.

[0040] The shear wave extraction module is used to extract the time-domain information of shear wave propagation;

[0041] The shear wave conversion module is used to convert time-domain information into frequency-domain information.

[0042] The shear wave velocity detection module is used to detect straight lines in the frequency domain using an improved Radon algorithm, and to detect the shear wave velocity based on the parameters of the detected straight lines.

[0043] Further specifying, the shear wave velocity detection module includes:

[0044] The line detection module is used to transform the data corresponding to the horizontal and vertical coordinates of the line y = ax + b in the XYD reference frame into the intercept and slope parameter space in the frequency domain information graph. Below, using the cyclic displacement method, in the intercept slope parameter space... The lower detection line transforms y = ax + b into the intercept and slope parameter space. The next corresponding point is (a, b);

[0045] The wave velocity calculation module is used to transform the point (a,b) into the parameter space of energy accumulation, and to obtain the slope and intercept of the corresponding straight line at the point with the highest energy accumulation. The shear wave velocity is then calculated based on the slope of the straight line.

[0046] Compared with the prior art, the present invention has the following beneficial effects:

[0047] This invention provides a method and system for calculating shear wave velocity in a homogeneous medium based on frequency domain analysis. Shear wave velocity calculation is a crucial parameter reflecting tissue viscoelasticity and a core component of ultrasonic shear wave elastography. Therefore, this invention processes data in the frequency domain. First, it uses a two-dimensional discrete Fourier transform to convert time-domain information into frequency-domain information, then performs translation processing. An improved Radon algorithm is used to detect straight lines in the frequency domain, and the shear wave velocity is determined based on the detected lines. This method significantly reduces noise in peak displacement algorithms in the time domain, more accurately reflecting real data. Furthermore, the improved Radon algorithm can fit straight lines with reduced computational load, significantly improving computational efficiency. Calculating the shear wave velocity in parameter space improves calculation accuracy and allows for more intuitive observation of the required data. It also improves the accuracy and noise resistance of reconstructed shear wave velocity images, overcomes the shortcomings of peak time algorithms, is less affected by environmental interference, and has good linearity, thus facilitating the acquisition of high-quality shear wave velocity image information. Ultimately, it obtains more accurate tissue viscoelasticity information and high-quality shear wave velocity images. Attached Figure Description

[0048] Figure 1 This is a schematic diagram of a shear wave velocity detection system for a uniform medium based on frequency domain analysis.

[0049] Figure 2 This is a flowchart illustrating the steps of the shear wave velocity calculation method of the present invention.

[0050] Figure 3 A time-domain information diagram of shear wave propagation;

[0051] Figure 4 A frequency domain information diagram for shear wave propagation;

[0052] Figure 5 Detecting the parameter space graph of a line using an improved Radon algorithm;

[0053] Figure 6 A schematic diagram of the traditional Radon algorithm;

[0054] Figure 7 This is a schematic diagram of the transformation of a straight line in parameter space;

[0055] Figure 8 Example diagram of the improved Radon algorithm;

[0056] Figure 9 Flowchart for prosthetic body manufacturing;

[0057] Figure 10 This is a diagram showing the trajectory of a shear wave propagating in a homogeneous medium.

[0058] Figure 11 This is a graph showing the information of shear waves in the frequency domain.

[0059] Figure 12 The calculation results of the improved Radon algorithm are shown in the figure. Detailed Implementation

[0060] To enable those skilled in the art to better understand the present invention, the technical solutions of the present invention will be clearly and completely described below with reference to the accompanying drawings of the embodiments. Obviously, the described embodiments are merely some embodiments of the present invention, and not all embodiments. Based on the embodiments of the present invention, all other embodiments obtained by those skilled in the art without creative effort should fall within the scope of protection of the present invention.

[0061] It should be noted that the terms "first," "second," etc., in the specification, claims, and accompanying drawings of this invention are used to distinguish similar objects and are not necessarily used to describe a specific order or sequence. It should be understood that such data can be used interchangeably where appropriate so that the embodiments of the invention described herein can be implemented in orders other than those illustrated or described herein. Furthermore, the terms "comprising" and "having," and any variations thereof, are intended to cover a non-exclusive inclusion; for example, a process, method, system, product, or apparatus that comprises a series of steps or units is not necessarily limited to those steps or units explicitly listed, but may include other steps or units not explicitly listed or inherent to such processes, methods, products, or apparatus.

[0062] The present invention will now be described in further detail with reference to the accompanying drawings;

[0063] See Figure 1 This invention provides a shear wave velocity detection system for a homogeneous medium based on frequency domain analysis, comprising:

[0064] The shear wave excitation module is used to generate shear waves in the phantom at a set focusing depth.

[0065] The shear wave extraction module is used to extract the time-domain information of shear wave propagation;

[0066] The shear wave conversion module is used to convert time-domain information into frequency-domain information and frequency-domain information graphs.

[0067] The shear wave velocity detection module is used to detect straight lines in frequency domain information and frequency domain information graphs using an improved Radon algorithm, and to detect the shear wave velocity based on the parameters of the detected straight lines.

[0068] To further explain, the shear wave velocity detection module includes:

[0069] The line detection module is used to transform the data corresponding to the horizontal and vertical coordinates of the line y = ax + b in the XYD reference frame into the intercept and slope parameter space in the frequency domain information graph. Below, using the cyclic displacement method, in the intercept slope parameter space... The lower detection line transforms y = ax + b into the intercept and slope parameter space. The next corresponding point is (a, b);

[0070] The wave velocity calculation module is used to transform the point (a,b) into the parameter space of energy accumulation, and to obtain the slope and intercept of the corresponding straight line at the point with the highest energy accumulation. The shear wave velocity is then calculated based on the slope of the straight line.

[0071] like Figure 2 As shown, the frequency domain analysis-based method for detecting shear wave velocity in a uniform medium, implemented by the above-mentioned frequency domain analysis-based system, includes the following steps:

[0072] (1) In the phantom, a shear wave is generated at a set focusing depth, and the time domain information of the shear wave propagation is extracted; the phantom refers to a homogeneous medium that is isotropic, linearly elastic and viscous, such as a biomimetic tissue medium.

[0073] (1.1) In the phantom, acoustic radiation force is generated at the set focusing depth to generate shear wave propagation data U(z,x,t) in the entire phantom, where z represents the axial dimension variable, x represents the transverse dimension variable, and t represents the slow time dimension variable.

[0074] (1.2) Take the data u(x,t) at the depth z0 where the focal point is located to obtain the spatiotemporal shear wave propagation image represented by the particle displacement signal, as shown below. Figure 3 As shown.

[0075] (2) Convert the time domain information from step (1) into frequency domain information;

[0076] (2.1) The time-domain information from step (1) is subjected to a two-dimensional discrete Fourier transform using the FFT function to obtain the spectral information of the shear wave propagating in a homogeneous medium.

[0077] (2.2) The spectrum information is shifted so that the zero frequency is in the center of the frequency domain image;

[0078] (2.3) The amplitude of the spectrum is obtained using the abs function, and the frequency domain information and frequency domain information diagram of the shear wave propagating in the homogeneous medium are obtained, as shown below:

[0079] S(k x f) = FFT[FFT] x,t (u(x,t))]

[0080] Where, k x The space wave number is represented by 1 / m; f represents the frequency by Hz.

[0081] This embodiment was completed in the MATLAB R2019b environment. The time-domain data was converted into frequency-domain data, and the FFT function was used to perform a fast two-dimensional Fourier transform on the time-domain information of the shear wave propagation. However, the zero-frequency signal after the two-dimensional fast Fourier transform may not be distributed in the center of the spectrum. Therefore, the fftshift function is needed to shift the spectrum information, that is, to shift the spectrum to the center position. Then, after using the abs function to obtain the amplitude of the spectrum, the frequency domain information of the shear wave propagating in the homogeneous medium can be obtained. Figure 4 This is a frequency domain information diagram of shear wave propagation.

[0082] Therefore, it can be seen that at this time, the signal energy is distributed along a straight line.

[0083] In general, the formula for calculating the group velocity of a simple harmonic wave is as follows.

[0084]

[0085] Combination Figure 4 As can be seen from the frequency domain representation of the shear wave transmission in the figure, the slope of the linear distribution of signal energy in the figure is the group velocity of the shear wave.

[0086] Therefore, we propose to use an improved Radon algorithm to detect the straight line and obtain the slope of the linear distribution of signal energy, thereby reversing the propagation speed of the shear wave.

[0087] (3) Detecting straight lines in the frequency domain using the improved Radon algorithm;

[0088] To facilitate a deeper understanding of the improved Radon algorithm, let's first briefly introduce the traditional Radon algorithm: The traditional Radon algorithm is a line detection algorithm, and its specific implementation process can be found in [link to relevant documentation]. Figure 6This is a schematic diagram of the traditional Radon algorithm. Here, d and θ are parameters of line L, and the line at different positions is determined by changing their values; f(x,y) is the region of lines to be detected; the point on L is determined by performing a line integral on line L within f(x,y); by repeatedly calculating the integral values ​​on L under different parameters, the line L with the largest integral value is selected as the most reasonable detected line; point R in the diagram is the point where the integral value falls, which can be used to observe the magnitude of the integral value.

[0089] While the traditional Radon algorithm can detect straight lines, it requires transforming the line L into a parameter-controlled element before performing multiple integration operations, resulting in a large computational load and a complex process. Based on the above shortcomings, the following algorithm will be improved to obtain the improved Radon algorithm of this invention.

[0090] The improved Radon algorithm of this invention is also a line detection algorithm. It can transform the region to be detected into a parameter space, and then use a cyclic shift method to detect lines. Compared with the Radon algorithm, the improved Radon algorithm has less computation and lower complexity. Its principle is as follows:

[0091] Establish a parameter space XYD, which is an M×N matrix, where M and N are the row and column numbers of the input frequency domain data, respectively; the value of the data in the L-th row and K-th column of parameter space XYD is D. XL , YK X L This represents the Lth value of the x-axis, where L = 1, 2, ..., N, Y K This represents the Kth value of the ordinate, where K = 1, 2, ..., M.

[0092] The problem of estimating parameters k and b is actually part of the line detection problem; refer to Figure 8 This is a first principle diagram for line detection in an embodiment of the present invention; Figure 8 In (a), the horizontal axis represents parameter X, and the vertical axis represents parameter Y; Figure 8 In (b), the horizontal axis represents the parameter k, and the vertical axis represents the parameter b; Figure 8 In (a), the values ​​of parameters k and b can be obtained simply by finding the intercept and slope of the line Y = kX + b. Figure 8 In (a), the straight line represented by a series of discrete points is a straight line in the plane XY with a slope of k and an intercept of b. Each point of this straight line in the XY space corresponds to a straight line in the kb space. The straight lines corresponding to each point on the straight line Y = kX + b in the plane kb all intersect at the same point. That is to say, the coordinates (k, b) of the intersection point of these straight lines in the kb space are the slope and intercept of the straight line in the XY space.

[0093] The cyclic shifting method in this embodiment is specifically as follows:

[0094] The intercept b is pre-defined to have M discretized values, representing b1 to b2 respectively. M The slope k is pre-defined to have N discretized values, represented as k1 to k2. N ;

[0095] Establish parameter space Parameter space The quantization scales of the two axes are the same as those of the original space XYD, i.e., they are also matrices of size M×N; parameter space The value of the data in the i-th row and j-th column is k i Let b represent the i-th value of k, where i = 1, 2, ..., N. j Let b be the j-th value, where j = 1, 2, ..., M; parameter space. The value of the j-th row and i-th column It should be: the sum of the values ​​corresponding to all two-dimensional coordinates (X,Y) that satisfy the equation Y=kX+b. It can be seen that the parameter space... The i-th column can be obtained by cyclically shifting and accumulating each column of the parameter space XYD.

[0096] The following is a detailed explanation of cyclic shifting:

[0097] When k i X L When the value is less than 0, the data in the Lth column of the parameter space XYD is cyclically shifted upwards by -k. i X L Bit;

[0098] When k i X L When the value is greater than 0, the data in the Lth column of the parameter space XYD is shifted downwards by k in a circular manner. i X L Bit;

[0099] When k i X L When the value is 0, no cyclic displacement operation is performed.

[0100] The following example illustrates the cyclic displacement process:

[0101] like Figure 10This is an example diagram of the improved Radon algorithm; the red dot represents the line Y = -X + 8; the green dot represents the data cyclic shift when X = -1, the blue dot represents the data cyclic shift when X = 2, and the orange dot represents the data cyclic shift when X = 1; the red circle marks the position where the maximum value of the accumulated points after the cyclic shift occurs; when X = -1, the corresponding data in the parameter space XYD is [0,0,0,0,0,0,0,0,0,1,0]. T Where T represents transpose; at this time, in the parameter space XYD, when k = -3, kX = 3 > 0, and the data should be shifted down by 3 units. The corresponding data is: [0,0,0,0,0,0,1,0,0,0,0] T In the parameter space XYD, when k = -1, kX = 1 > 0. In this case, the data should be shifted downwards by 1 unit. The corresponding data is: [0,0,0,0,0,0,0,0,1,0,0] T Similarly, when X = 1, in the parameter space XYD, when k = -1, kX = -1 < 0. In this case, the data should be shifted upwards by 1 unit. The corresponding data is: [0,0,0,0,0,0,0,0,1,0,0] T When X = 2, in the parameter space XYD, when k = -1, kX = -2 < 0. In this case, the data should be shifted upwards by 2 units. The corresponding data is: [0,0,0,0,0,0,0,0,1,0,0] T When k = -1, the parameter space The corresponding data, when summed up, is: [0,0,0,0,0,0,0,0,3,0,0] T The accumulated data for other k values ​​are not as large as those for b=8; after taking multiple discrete points and then performing cyclic displacement accumulation, the values ​​in the parameter space are... This will result in a point with a very large value (-1, 8), which corresponds to the line Y = -X + 8 in the parameter space XYD. In other words, in the frequency domain information graph, the data corresponding to the x and y coordinates of the line y = ax + b in the XYD reference frame will be transformed into the intercept and slope parameter space. Below, using the cyclic displacement method, in the intercept slope parameter space... The lower detection line transforms y = ax + b into the intercept and slope parameter space. The corresponding point (a, b) is used to complete the improved Radon algorithm transformation of the frequency domain information.

[0102] Typically, the propagation of shear waves in the time domain can be observed in a time-displacement coordinate system. If the frequency domain information obtained after Fourier transform is represented in a two-dimensional plane as y = ax + b, then 'a' corresponds to the shear wave velocity. However, calculations in the time domain are affected by noise, impacting the accuracy of the data and causing slight errors in the calculated shear wave velocity. Therefore, to address the noise issue, the frequency domain can mitigate the noise's influence and improve the accuracy of shear wave velocity calculations. Thus, in the frequency domain, time and displacement are transformed into frequency and wave number. In this two-dimensional plane, since... (where V) g ω is the group velocity (m / s); k is the angular frequency (rad / s); k is the space wave number (rad / m). The slope of the detected straight line can represent the shear wave velocity.

[0103] Figure 5 The parameter space graph of the improved Radon algorithm is used to detect the line. At this point, it can be seen that there is a clear line in the graph.

[0104] (4) Based on the parameters of the detection line, the shear wave velocity is obtained.

[0105] The velocity of a shear wave refers to the ratio of the parameter space XYD to the parameter space... The transformation of corresponding points in the middle, if the line is y = ax + b, then the corresponding... In space, it would be (a, b);

[0106] Will Figure 5 The point (a,b) after the improved Radon transform is transformed to Figure 7 In the parameter space of energy accumulation, the point with the highest energy accumulation can be obtained, which corresponds to the slope and intercept of the straight line. This can be achieved using the following formula:

[0107]

[0108] V g ω is the group velocity, in m / s; ω is the angular frequency, in rad / s; k is the space wave number, in rad / m.

[0109] This allows us to determine the slope of the straight line, which is the shear wave velocity.

[0110] The calculated shear wave velocity was compared with the data from QIBA (Quantitative Imaging Biomarkers Alliance), and the results are shown in Table 1.

[0111] Table 1 compares the two results for determining the shear wave velocity in a homogeneous medium.

[0112] The method described in this invention will be verified by specific experiments below:

[0113] In this experiment, the instruments and materials used for phantom preparation included gelatin, nano-titanium dioxide, a digitally displayed constant-temperature magnetic stirring oil bath, and a vacuum drying oven; the phantom fabrication process was as follows: Figure 9 Flowchart of prosthetic body manufacturing process:

[0114] Dissolve the gelatin granules in deionized water (7% gelatin) according to the specified ratio. Dissolve in an oil bath at 60°C and stir for half an hour to one hour. Note that the beaker should be covered with plastic wrap during dissolution to prevent the absorption of moisture from the air. After complete dissolution, the solution should be free of bubbles and transparent.

[0115] Dissolve a certain amount of titanium dioxide powder in anhydrous ethanol (1.5% titanium dioxide) and stir for half an hour to ensure thorough mixing.

[0116] The titanium dioxide solution was injected into the gelatin solution using a pipette, and then stirred again in a 60°C oil bath for half an hour to one hour to obtain a uniform gel. The mixture was then vacuumed at 70°C using a vacuum drying oven and left at room temperature to form an elastomer.

[0117] The experiment used a Verasonics ultrasonic testing system for shear wave excitation, acquisition, and imaging. An L11-4V ultrasonic linear array probe with a center frequency of 6.25MHz was used, employing dynamic focusing excitation. A continuous sine wave was emitted as the excitation signal, and ultrafast plane wave imaging was used to capture the shear wave. The system's multifocal excitation and acquisition parameters are shown in Table 2.

[0118] Table 2 Multi-focus excitation acquisition parameters

[0119]

[0120]

[0121] Shear waves were generated on the phantom using multifocal excitation, and then the propagation trajectory of the shear waves at 1ms intervals was extracted as shown in the diagram. Figure 10 As shown;

[0122] The time-domain information at all surfaces along the direction of shear wave propagation will be converted into frequency-domain information, such as... Figure 11 As shown;

[0123] Figure 12 The shear wave velocity was calculated using an improved Radon algorithm.

[0124] The shear wave velocities calculated in the prepared gelatin phantoms are shown in Table 3.

[0125] Table 3 shows the calculated shear wave velocity.

[0126]

[0127] The TTP algorithm in the table above is a classic method for calculating shear wave velocity. The shear wave velocity calculated by the improved Radon algorithm is very close to the result of the TTP algorithm, which shows that the proposed shear wave velocity algorithm has high accuracy and practicality.

[0128] Experimental results show that the improved Radon algorithm can calculate the velocity of shear waves relatively accurately, providing a new approach to calculating shear wave velocity. Moreover, it has high overall detection efficiency, is less affected by environmental interference, and has good noise resistance.

[0129] Finally, it should be noted that the above embodiments are only used to illustrate the technical solutions of the present invention and not to limit it. Although the present invention has been described in detail with reference to the above examples, those skilled in the art should understand that modifications or equivalent substitutions can still be made to the specific implementation schemes of the present invention. Any modifications or equivalent substitutions that do not depart from the spirit and scope of the present invention should be covered within the scope of protection of the claims of the present invention.

Claims

1. A method for detecting shear wave velocity in a homogeneous medium based on frequency domain analysis, characterized in that, Includes the following steps: (1) In the phantom, a shear wave is generated at the set focusing depth, and the time domain information of the shear wave propagation is extracted; (2) Convert the time-domain information from step (1) into frequency-domain information; (3) The detection line is obtained from the frequency domain information using the improved Radon algorithm; Step (3) specifically involves: In the frequency domain information graph, the data corresponding to the horizontal and vertical coordinates of the line y=ax+b in the XYD reference system are transformed into the intercept slope parameter space kb-φ. Using the cyclic displacement method, the line is detected in the intercept slope parameter space kb-φ, so that y=ax+b is transformed into a point (a,b) corresponding to the intercept slope parameter space kb-φ, thus completing the improved Radon algorithm transformation of the frequency domain information. The cyclic displacement method is specifically as follows: The intercept b of the line y=ax+b has M discretized values, representing b1 to b2. M The slope k has N discretized values, represented as k1 to k2. N ; Establish a parameter space kb-φ, where the quantization scale of the two axes in parameter space kb-φ is the same as that of the original space XYD, both being matrices of size M×N; In the parameter space kb-φ, the data in the i-th row and j-th column takes the value φ. ki,bj k i Let b represent the i-th value of k, i = 1, 2, ..., N. j This represents the j-th value of b, where j = 1, 2, ..., M; and the value φ in the i-th row and j-th column of the parameter space kb-φ. ki,bj It should be: the sum of the values ​​corresponding to all two-dimensional coordinates (X,Y) that satisfy the equation Y=kX+b, thereby determining that the i-th column of the parameter space kb-φ can be obtained by cyclically displacing and accumulating each column of the parameter space XYD; When k i X L When the value is less than 0, the data in the Lth column of the parameter space XYD is cyclically shifted upwards by -k. i X L Bit; When k i X L When the value is greater than 0, the data in the Lth column of the parameter space XYD is shifted downwards by k in a circular manner. i X L Bit; When k i X L When the value is 0, no cyclic displacement operation is performed; Complete the loop; (4) Based on the parameters of the detection line, the shear wave velocity is obtained; Step (4) specifically involves: After transforming the point (a,b) obtained by the improved Radon algorithm, we can convert it into the parameter space of energy accumulation. The point with the highest energy accumulation is then determined by the slope and intercept of the corresponding line. This can be calculated using the following formula: Group velocity, unit: m / s; Angular frequency, unit: rad / s; Space wavenumber, unit: rad / m; The slope of the straight line is the shear wave velocity.

2. The method for detecting shear wave velocity in a uniform medium based on frequency domain analysis according to claim 1, characterized in that, The specific steps (1) are as follows: (1.1) In the phantom, acoustic radiation force is generated at the set focusing depth to generate shear wave propagation data U(z,x,t) in the entire phantom, where z represents the axial dimension variable, x represents the transverse dimension variable, and t represents the slow time dimension variable. (1.2) Take the data U(x,t) at the depth z0 where the focal point is located, and obtain the spatiotemporal shear wave propagation image represented by the particle displacement signal, and obtain the time domain information of the spatiotemporal shear wave propagation image.

3. The method for detecting shear wave velocity in a uniform medium based on frequency domain analysis according to claim 1, characterized in that, The phantom refers to a homogeneous medium that is isotropic, linearly elastic, and viscous.

4. The method for detecting shear wave velocity in a uniform medium based on frequency domain analysis according to claim 1, characterized in that, Step (2) specifically involves: (2.1) The time-domain information from step (1) is subjected to a two-dimensional discrete Fourier transform using the FFT function to obtain the spectral information of the shear wave propagating in a homogeneous medium; (2.2) The spectral information is shifted so that the zero frequency is in the center of the frequency domain image, as shown below: in, Space wavenumber, unit: 1 / m; Frequency is expressed in Hz. (2.3) The amplitude of the spectrum is obtained by using the abs function, and the frequency domain information and frequency domain information diagram of the shear wave propagating in the homogeneous medium are obtained.

5. The method for detecting shear wave velocity in a homogeneous medium based on frequency domain analysis according to claim 4, characterized in that, The specific step (2.2) is as follows: use the FFTSHIFT function to shift the spectral information so that the zero frequency is in the center of the frequency domain image.

6. A system for detecting shear wave velocity in a uniform medium based on frequency domain analysis, characterized in that, The method for detecting shear wave velocity in a uniform medium based on frequency domain analysis, as described in any one of claims 1 to 5, specifically includes: The shear wave excitation module is used to generate shear waves in the phantom at a set focusing depth. The shear wave extraction module is used to extract the time-domain information of shear wave propagation; The shear wave conversion module is used to convert time-domain information into frequency-domain information. The shear wave velocity detection module is used to obtain the detection line in the frequency domain information using the improved Radon algorithm, and to detect the shear wave velocity based on the parameters of the detection line.

7. The shear wave velocity detection system for a uniform medium based on frequency domain analysis according to claim 6, characterized in that, The shear wave velocity detection module includes: The line detection module is used to transform the data corresponding to the horizontal and vertical coordinates of the line y=ax+b in the XYD reference system in the frequency domain information graph to the intercept slope parameter space kb-φ. Using the cyclic displacement method, the line is detected in the intercept slope parameter space kb-φ, so that y=ax+b is transformed into a point (a,b) corresponding to the intercept slope parameter space kb-φ. The wave velocity calculation module is used to transform the point (a,b) into the parameter space of energy accumulation, and to obtain the slope and intercept of the corresponding straight line at the point with the highest energy accumulation. The shear wave velocity is then calculated based on the slope of the straight line.

Images