Three-dimensional ultrasonic tomographic image reconstruction method based on time-domain full waveform inversion

Abstract

The application relates to the technical field of medical image reconstruction, and discloses a three-dimensional ultrasonic tomographic image reconstruction method based on time-domain full waveform inversion, which comprises the following steps: establishing a time-domain wave equation of ultrasonic wave propagation in biological tissues, generating a simulated sound pressure distribution according to the time-domain wave equation; using a mean square error to measure the error between a real sound pressure distribution and the simulated sound pressure distribution, using an adjoint field method to calculate the functional gradient of the mean square error, and using a momentum gradient descent method to minimize the error; in the process of minimizing the error, reconstructing a sound velocity distribution according to the functional gradient of the error with respect to the sound velocity, then using the reconstructed sound velocity distribution to reconstruct a sound absorption coefficient distribution, and repeatedly and alternately reconstructing until a set condition is reached, so that the reconstruction of the sound velocity distribution and the sound absorption coefficient distribution is completed. The application simultaneously reconstructs the sound velocity and the sound absorption coefficient in the time domain, avoids the problem that a traditional full waveform inversion algorithm needs to reconstruct the sound velocity in the time domain and reconstruct the sound absorption coefficient in the frequency domain, and reduces the algorithm complexity.

Description

Technical Field

[0001] This invention relates to the field of medical image reconstruction technology, specifically to a method for reconstructing three-dimensional ultrasound tomographic images based on time-domain full waveform inversion. Background Technology

[0002] Computational ultrasound computed tomography (USCT) is a novel ultrasound imaging modality that combines ultrasound imaging with computed tomography. Traditional medical ultrasound relies on pulse-echo signals to locate ultrasound scatterers and image tissue structures; it can also measure blood flow velocity using the Doppler effect. As a new ultrasound imaging modality, USCT can quantitatively image the acoustic properties of tissues, such as sound velocity, sound absorption coefficient, and density, providing information unavailable with traditional medical ultrasound.

[0003] Imaging algorithms based on time-of-flight (TOF) are widely used in computational ultrasound tomography. The simplest acoustic model in these algorithms assumes that sound waves propagate in a straight ray, neglecting scattering, diffraction, and refraction during propagation. This leads to inaccurate TOF calculations, reconstructing only a rough distribution of sound velocity or absorption coefficient. Bent-ray models can account for refraction during sound propagation to some extent. Their calculation method involves a linear correlation between the measured TOF and the integral of the slowness along the ray (the reciprocal of the sound velocity). By constructing a system matrix, standard linear algebra methods are used for inversion. Bent-ray models can obtain more accurate arrival times, thus reconstructing higher resolution and more accurate sound velocity or absorption coefficient distributions. Full waveform inversion (FWI) algorithms, based on the acoustic wave equation, consider the complete physical characteristics of sound wave propagation and obtain the equation parameter distributions (sound velocity distribution, absorption coefficient distribution, and density distribution) through optimization methods, enabling the reconstruction of the highest resolution and most accurate images.

[0004] In recent years, various three-dimensional ultrasound computational tomography systems have emerged. For example, QT Imaging developed a system using a plane wave transmitter and a two-dimensional plane detector; Delphinus Medical Technologies introduced a scanning ring transducer array system; and the Karlsruhe Institute of Technology designed a system based on a hemispherical transducer array. Despite the progress in three-dimensional transducer array design, many full-waveform inversion algorithms remain limited to two-dimensional reconstruction, assuming that sound wave propagation is confined to the transducer plane. This assumption introduces errors because actual sound waves propagate in three dimensions, and reflected waves outside the plane may re-enter the reconstruction plane, introducing artifacts. Such simplification often incorrectly characterizes out-of-plane energy loss or gain, particularly affecting the accuracy of the reconstructed acoustic absorption coefficient distribution. Another challenge of full-waveform inversion is the need for different computational frameworks to effectively reconstruct different acoustic properties. Typically, sound velocity distribution reconstruction is performed in the time domain, while acoustic absorption coefficient distribution is usually reconstructed in the frequency domain. Existing techniques have proposed methods to simultaneously reconstruct sound velocity and acoustic absorption coefficient distributions in the time domain, but these methods are still limited to two-dimensional reconstruction. Summary of the Invention

[0005] To address the aforementioned technical problems, this invention provides a method for reconstructing three-dimensional ultrasonic tomographic images based on time-domain full waveform inversion, which can simultaneously reconstruct the three-dimensional distribution of sound velocity and sound absorption coefficient.

[0006] To solve the above-mentioned technical problems, the present invention adopts the following technical solution:

[0007] A method for reconstructing three-dimensional ultrasound tomographic images based on time-domain full waveform inversion includes the following steps:

[0008] Step 1: Establish the time-domain wave equation for the propagation of ultrasound in biological tissues, and generate a simulated sound pressure distribution ρ(r) based on the time-domain wave equation. m ,t);

[0009] Step two, use mean square error to measure the true sound pressure distribution p measured by the transducer. obs (r m ,t) and the simulated sound pressure distribution p(r) generated using the time-domain wave equation m The error ε between ,t) is calculated, and the functional gradients of ε with respect to the sound speed c and the sound absorption coefficient α0 are obtained based on the adjoint field method. The distributions of the sound speed c and the sound absorption coefficient α0 are updated by the momentum gradient descent method to minimize the error ε.

[0010] Step 3, Alternating Reconstruction: The sound velocity distribution is reconstructed based on the functional gradient of ε with respect to the sound velocity. When the error ε changes less than a set threshold between two iterations, the sound velocity distribution is considered to have converged. Then, the reconstructed sound velocity distribution is used to reconstruct the sound absorption coefficient distribution. Next, the reconstructed sound absorption coefficient distribution is used to iteratively update the sound velocity distribution. The alternating reconstruction process is repeated until the set stopping condition is reached, thus completing the reconstruction of the sound velocity distribution and the sound absorption coefficient distribution.

[0011] Furthermore, step one specifically includes:

[0012] The propagation of ultrasound waves in biological tissues is described by a set of first-order equations:

[0013]

[0014] Where v is the particle velocity of the sound wave transmission medium, r represents the spatial position, t represents time, p and ρ are the sound pressure distribution and sound density distribution, respectively, ρ0 represents the ambient density, m is the mass source term, and d represents the particle displacement of the sound wave transmission medium.

[0015] For lossy media, the sound absorption term follows a frequency power law:

[0016] α=α0ω y (4)

[0017] Where ω is the time frequency, y is the power law exponent, α represents the sound absorption term, and α0 represents the sound absorption coefficient;

[0018] Sound absorption and sound dispersion are defined using a fractional-order Laplace model:

[0019]

[0020] τ(r)=-2α0(r)c(r) y-1 (6)

[0021] η(r)=2α0(r)c(r) y tan(πy / 2); (7)

[0022] Where τ and η are the parts of equation (5) related to sound absorption and sound dispersion, respectively, and L represents the fractional-order Laplace operator. Represents the Nabla operator;

[0023] Integrating equations (1) to (7), we get:

[0024]

[0025] in, Known as the lossy d'Alembert operator:

[0026]

[0027] Where p(r,t) is an equation about spatial location r and time t, used to describe sound pressure; S is the sound source term; equations (1) to (3) are used to generate the simulated sound pressure distribution p(r). m ,t).

[0028] Furthermore, in step two, the mean square error is used to measure the true sound pressure distribution p measured by the transducer. obs (r m ,t) and the simulated sound pressure distribution p(r) generated using the time-domain wave equation m The error ε between t and t specifically includes:

[0029]

[0030] Where M represents the number of transmit-receive pairs, M = N(N-1), N represents the number of transducers; T represents the total duration of the transducer's signal detection.

[0031] Further, in step two, the functional gradients of ε with respect to the sound velocity c and the sound absorption coefficient α0, based on the adjoint field method, specifically include:

[0032] Partial derivative of ε with respect to the speed of sound for:

[0033]

[0034] Solve equation (15) using the adjoint field method:

[0035] The adjoint field wave equation is defined as:

[0036]

[0037] p * (r, Tt) is called p(r) m The accompanying field of t), r m The position of the transducer is indicated by T, and the total duration of the transducer's detection signal is indicated by T. Equation (15) can be written as:

[0038]

[0039]

[0040] c(r) represents the speed of sound at spatial location r, where Γ1, Γ2, and Γ3 are all intermediate variables;

[0041] In equation (17), Γ2 and Γ3 are much smaller than Γ1. Ignoring Γ2 and Γ3, equation (17) simplifies to:

[0042]

[0043] Solve equation (19) to calculate the functional gradient of ε with respect to the speed of sound c;

[0044] The functional gradient of ε with respect to the acoustic absorption coefficient α0 is calculated based on the adjoint field method:

[0045]

[0046] The spatial gradient is solved in the Fourier domain:

[0047]

[0048] in, and These are the Fourier transform and the inverse Fourier transform, respectively, where k is the wave number.

[0049] Furthermore, in step two, updating the distributions of the sound velocity c and the sound absorption coefficient α0 using the momentum gradient descent method to minimize the error ε specifically includes:

[0050]

[0051] θ n+1 =θ n -λ n v n (12)

[0052] Among them, v n It is the velocity vector of the nth iteration; γ is the momentum parameter; It is the gradient parameter; θ n λ is the parameter updated in the nth iteration. n This represents the step size of the nth iteration; Indicates J relative to θ n The partial derivatives; J represents the loss function;

[0053] The update formulas for sound velocity and sound absorption coefficient are as follows:

[0054]

[0055] c n This represents the speed of sound after the nth iteration. v represents the step size of the sound speed in the (n+1)th iteration. c，n-1 Let α represent the velocity vector in the momentum gradient descent method for the sound speed at the (n-1)th iteration. 0，n This represents the sound absorption coefficient after the nth iteration. v represents the step size of the acoustic absorption coefficient in the (n+1)th iteration. α，n-1 Let represent the velocity vector of the momentum gradient descent method when the acoustic absorption coefficient is obtained in the (n-1)th iteration.

[0056] Compared with the prior art, the beneficial technical effects of the present invention are:

[0057] Compared to traditional time-of-flight-based algorithms, the three-dimensional ultrasound tomography image reconstruction method based on time-domain full waveform inversion proposed in this invention considers all physical processes in ultrasound imaging, and can obtain reconstructed images with higher accuracy.

[0058] This invention performs sound field simulation and image reconstruction in three-dimensional space, avoiding the image artifacts and accuracy degradation problems inherent in traditional two-dimensional full-waveform inversion. Furthermore, this invention simultaneously reconstructs sound velocity and sound absorption coefficient in the time domain, avoiding the need for traditional full-waveform inversion to reconstruct sound velocity in the time domain and sound absorption coefficient in the frequency domain, thus reducing algorithm complexity. Attached Figure Description

[0059] Figure 1 A flowchart of the present invention is shown schematically.

[0060] Figure 2 A schematic diagram of a simulation setup according to an embodiment of the present invention is shown; wherein, (a) is a numerical phantom of a three-dimensional cardiovascular shape used in the simulation (maximum intensity projection); (b) is a 45° view projection of a spherical transducer used in the simulation; (c) is the relative position of the transducer and the phantom during the simulation process; and (d) is the maximum intensity projection of the phantom.

[0061] Figure 3 The diagram schematically illustrates the simulation results of a single acoustic property distribution reconstruction according to the present invention. In this diagram, (a) through (h) all use images processed by maximum intensity projection; (a) is the true value of the sound velocity distribution; (b) is the input value of the sound velocity distribution; (c) is the reconstructed sound velocity distribution; (d) is the absolute error between the reconstructed result and the true value; (e) is the true value of the sound absorption coefficient distribution; (f) is the input value of the sound absorption coefficient distribution; (g) is the reconstructed sound absorption coefficient distribution; (h) is the absolute error between the reconstructed sound absorption coefficient result and the true value; (i) through (l) are the cross-sections corresponding to (a) through (d); and (m) through (p) are the cross-sections corresponding to (e) through (h).

[0062] Figure 4 The diagram illustrates the variation of the root mean square error between the reconstructed result and the true value when reconstructing a single acoustic property distribution. In this diagram, (a) shows the variation corresponding to the sound velocity, and (b) shows the variation corresponding to the sound absorption coefficient.

[0063] Figure 5The diagram schematically illustrates the reconstruction results of simultaneously reconstructing the sound velocity and sound absorption coefficient distributions according to the present invention. In the diagram, (a) represents the reconstructed sound velocity distribution (maximum intensity projection); (b) represents a cross-section of the reconstructed sound velocity distribution; (c) represents the absolute error between the reconstructed sound velocity distribution and the true value; (d) represents a cross-section of the absolute error between the reconstructed sound velocity distribution and the true value; (e) represents the reconstructed sound absorption coefficient distribution (maximum intensity projection); (f) represents a cross-section of the reconstructed sound absorption coefficient distribution; (g) represents the absolute error between the reconstructed sound absorption coefficient distribution and the true value; and (h) represents a cross-section of the absolute error between the reconstructed sound absorption coefficient distribution and the true value.

[0064] Figure 6 The diagram illustrates the variation of the root mean square error between the reconstructed distribution map and the true distribution map when the present invention simultaneously reconstructs the distributions of sound velocity and sound absorption coefficient. In this diagram, (a) shows the variation corresponding to sound velocity, and (b) shows the variation corresponding to sound absorption coefficient. Detailed Implementation

[0065] A preferred embodiment of the present invention will now be described in detail with reference to the accompanying drawings.

[0066] like Figure 1 As shown, a three-dimensional ultrasound tomographic image reconstruction method based on time-domain full waveform inversion in this invention includes establishing a time-domain wave equation describing the ultrasound propagation process, optimization based on gradient descent method, functional gradient calculation based on adjoint field method, and alternating reconstruction.

[0067] 1. Establish the time-domain wave equation describing the propagation process of ultrasound:

[0068] The propagation of ultrasound waves in biological tissues can be described by a set of first-order equations:

[0069]

[0070] Where v is the particle velocity of the sound wave transmission medium, r represents the spatial position, p and ρ are the sound pressure distribution and sound density distribution, respectively, ρ0 represents the ambient density, m is the mass source term, and c represents the speed of sound.

[0071] For lossy media, the sound absorption term follows a frequency power law:

[0072] α=α0ω y (4)

[0073] Where α0 is the sound absorption coefficient, with units of Np (rad / s). -y m -1 ω is the time frequency, measured in rad / s; y is the power law exponent.

[0074] Sound absorption and sound dispersion are defined using a fractional-order Laplace model:

[0075]

[0076] τ(r)=-2α0(r)c(r) y-1 (6)

[0077] η(r)=2α0(r)c(r) y tan(πy / 2); (7)

[0078] Where τ and η represent the sound absorption and sound dispersion terms, respectively. Integrating the above equations, we get:

[0079]

[0080] in, This is called the lossy d'Alembert operator, defined as follows:

[0081]

[0082] Where p(r,t) is an equation describing sound pressure with respect to spatial location r and time t; S is a sound source term. This invention uses equations (1) to (3) to generate the simulated sound pressure distribution p(r). m ,t).

[0083] 2. Optimization based on gradient descent:

[0084] This invention uses mean square error to measure the true sound pressure distribution p measured by the transducer. obs (r m The simulated sound pressure distribution p(r) obtained using the time-domain wave equation is given by t) and t). m The error ε between ,t) is:

[0085]

[0086] Where M represents the number of transmit-receive pairs, M = N(N-1), and N represents the number of transducers. This invention uses the gradient descent method to minimize the error ε in equation (10). The functional gradients of ε with respect to the sound velocity c and the sound absorption coefficient α0 are calculated respectively. and By applying the gradient descent method, the distributions of the sound velocity c and the sound absorption coefficient α0 that result in a smaller ε can be updated.

[0087] To mitigate potential oscillations caused by the traditional gradient descent method, this invention employs the momentum gradient descent method to make the update process smoother. The parameter update formula is as follows:

[0088]

[0089] θ n+1 =θ n -λ n v n (12)

[0090] Among them, v n η is the velocity vector of the nth iteration; γ is the momentum parameter (between 0 and 1); η is the gradient parameter, usually set to (1-γ); θ n v represents the parameter updated in the nth iteration. n This represents the step size for the nth iteration. The update formulas for the sound velocity and sound absorption coefficient are as follows:

[0091]

[0092] 3. Functional gradient calculation based on the adjoint field method:

[0093] The partial derivative of ε with respect to the speed of sound can be written as:

[0094]

[0095] This invention uses the adjoint field method to solve this partial differential. The adjoint field wave equation is defined as:

[0096]

[0097] p * (r, Tt) is called p(r) m The adjoint field of t). The adjoint field method then allows equation (15) to be written as...

[0098]

[0099] In equation (17), Γ2 and Γ3 are much smaller than Γ1. Therefore, Γ2 and Γ3 can be ignored, and equation (17) can be simplified to

[0100]

[0101] Existing technology has proven that equation (19) can yield a relatively accurate solution.

[0102] The gradient of the acoustic absorption coefficient based on the adjoint field method is calculated as follows:

[0103]

[0104] The spatial gradient is solved in the Fourier domain:

[0105]

[0106] in and These are the Fourier transform and the inverse Fourier transform, respectively, where k is the wave number.

[0107] 4. Alternating Reconstruction:

[0108] As can be seen from equations (19) and (20), the gradient calculation of the sound absorption coefficient depends on the accurate sound velocity distribution, while the gradient calculation of the sound velocity depends less on the sound absorption coefficient. Therefore, this invention employs an alternating reconstruction method to reconstruct the distributions of the two acoustic properties. First, the sound velocity is reconstructed according to equation (19). When the value of ε changes by less than 1% between two iterations, it is considered to have converged. After the sound velocity reconstruction converges, the reconstructed sound velocity distribution is used to reconstruct the sound absorption coefficient. After the reconstruction of the sound absorption coefficient converges, the reconstructed sound absorption coefficient is used to iteratively update the sound velocity distribution again to make it more accurate. By repeating the above process, the accurate sound velocity and sound absorption coefficient distributions can finally be obtained.

[0109] The technical solution of the present invention is further illustrated by the following two embodiments. Both embodiments are based on MATLAB programming and run on the same computer. Their main configuration is: Intel 13700K CPU, 32GB of memory and an NVIDIA RTX 4070Ti GPU with 12GB of video memory.

[0110] Example 1

[0111] The numerical phantom of the heart and blood vessel shape is used for three-dimensional full waveform inversion simulation reconstruction. In this embodiment, the actual sound absorption coefficient distribution is used when reconstructing the sound velocity, and the actual sound velocity distribution is used when reconstructing the sound absorption coefficient. The alternating reconstruction method is not used.

[0112] This embodiment uses a numerical phantom with a heart-shaped blood vessel. The phantom is located at the exact center of the computational grid, such as... Figure 2 As shown in (a) in the figure. Figure 2 (b) in the figure is a 45° projection of the spherical transducer used in the simulation; Figure 2 (c) in the figure represents the relative position of the transducer and the phantom during the simulation process; Figure 2 In the diagram, (d) represents the maximum intensity projection of the phantom, where the white dashed line indicates the location of the cross-section shown later.

[0113] The acoustic simulation was implemented using the k-Wave Toolbox in MATLAB. The detector array is a spherical array with its center at the center of the computational grid. The computational grid size is 54mm × 54mm × 54mm, discretized into 108 × 108 × 108 grids. The time-domain sampling length is 740, and the sampling interval is 100ns. In this embodiment, Figure 3 (a) in the figure represents the true value of the sound speed distribution; Figure 3 In the diagram, (b) represents the input value for the sound speed distribution; Figure 3 (c) in the figure represents the reconstructed sound velocity distribution; Figure 3 In the figure, (d) represents the absolute error between the reconstruction result and the true value; Figure 3 In this context, (e) represents the true value of the sound absorption coefficient distribution; Figure 3 In this context, (f) represents the input value for the sound absorption coefficient distribution; Figure 3 In the diagram, (g) represents the reconstructed distribution of the sound absorption coefficient. Figure 3 In this context, (h) represents the absolute error between the reconstructed sound absorption coefficient and the true value. Figure 3 (i) to (l) are Figure 3 The cross sections corresponding to (a) to (d) in the text; Figure 3 (m) to (p) in the middle are Figure 3 The cross-sections corresponding to (e) to (h) in the diagram. From Figure 3 As can be seen, this invention can achieve high-precision reconstruction of acoustic properties when reconstructing only one acoustic property. Figure 4 Figure (a) shows the variation of the root mean square error between the reconstructed sound velocity and the true value in this embodiment. Figure 4 (b) in this example shows the variation of the root mean square error between the reconstructed acoustic absorption coefficient and the true value.

[0114] In this invention, the unit of sound speed is m / s, and the unit of sound absorption coefficient is dB / (cm·MHz). y ).

[0115] Example 2

[0116] This embodiment performs three-dimensional full waveform inversion simulation and reconstruction of the numerical phantom of the heart and blood vessel shape, and simultaneously reconstructs the distribution of sound velocity and sound absorption coefficient.

[0117] The numerical model of the heart and blood vessels is the same as in Example 1, such as... Figure 2 Unlike Example 1, which only considers sound velocity or sound absorption coefficient in reconstruction, Example 2 reconstructs the distribution of sound velocity and sound absorption coefficient simultaneously based on an alternating reconstruction method. Figure 5 The reconstruction results of this embodiment are shown. From Figure 5 As can be seen, while the present invention cannot achieve the accuracy level of Example 1 when simultaneously reconstructing two acoustic properties, it can still achieve good accuracy. Figure 6 Figure (a) shows the change in root mean square error between the reconstructed sound velocity distribution map and the actual map in this embodiment. Figure 6 (b) in this example shows the change in root mean square error between the reconstructed distribution map of the acoustic absorption coefficient and the actual map.

[0118] It will be apparent to those skilled in the art that the present invention is not limited to the details of the exemplary embodiments described above, and that the invention can be implemented in other specific forms without departing from its spirit or essential characteristics. Therefore, the embodiments should be considered in all respects as exemplary and non-limiting, and the scope of the invention is defined by the appended claims rather than the foregoing description. Thus, all variations falling within the meaning and scope of equivalents of the claims are intended to be included within the present invention, and no reference numerals in the claims should be construed as limiting the scope of the claims.

[0119] Furthermore, it should be understood that although this specification describes embodiments, not every embodiment contains only one independent technical solution. This narrative style is merely for clarity. Those skilled in the art should consider the specification as a whole, and the technical solutions in each embodiment can also be appropriately combined to form other embodiments that can be understood by those skilled in the art.

Claims

1. A method for reconstructing three-dimensional ultrasound tomographic images based on time-domain full waveform inversion, comprising the following steps: Step 1: Establish the time-domain wave equation for the propagation of ultrasound in biological tissues, and generate a simulated sound pressure distribution p(r) based on the time-domain wave equation. m ,t); Step two, use mean square error to measure the true sound pressure distribution p measured by the transducer. obs (r m ,t) and the simulated sound pressure distribution p(r) generated using the time-domain wave equation m The error ε between ,t) is calculated, and the functional gradients of ε with respect to the sound speed c and the sound absorption coefficient α0 are obtained based on the adjoint field method. The distributions of the sound speed c and the sound absorption coefficient α0 are updated by the momentum gradient descent method to minimize the error ε. Step 3, Alternating Reconstruction: The sound velocity distribution is reconstructed based on the functional gradient of ε with respect to the sound velocity. When the error ε changes less than a set threshold between two iterations, the sound velocity distribution is considered to have converged. Then, the reconstructed sound velocity distribution is used to reconstruct the sound absorption coefficient distribution. Next, the reconstructed sound absorption coefficient distribution is used to iteratively update the sound velocity distribution. The alternating reconstruction process is repeated until the set stopping condition is reached, thus completing the reconstruction of the sound velocity distribution and the sound absorption coefficient distribution.

2. The three-dimensional ultrasound tomography image reconstruction method based on time-domain full waveform inversion according to claim 1, characterized in that, Step one specifically includes: The propagation of ultrasound waves in biological tissues is described by a set of first-order equations: Where v is the particle velocity of the sound wave transmission medium, r represents the spatial position, t represents time, p and ρ are the sound pressure distribution and sound density distribution, respectively, ρ0 represents the ambient density, m is the mass source term, and d represents the particle displacement of the sound wave transmission medium. For lossy media, the sound absorption term follows a frequency power law: α=α0ω y (4) Where ω is the time frequency, y is the power law exponent, α represents the sound absorption term, and α0 represents the sound absorption coefficient; Sound absorption and sound dispersion are defined using a fractional-order Laplace model: τ(r)=-2α0(r)c(r) y-1 (6) η(r)=2α0(r)c(r) y tan(πy / 2); (7) Where τ and η are the parts of equation (5) related to sound absorption and sound dispersion, respectively, and L represents the fractional-order Laplace operator. Represents the Nabla operator; Integrating equations (1) to (7), we get: in, Known as the lossy d'Alembert operator: Where p(r,t) is an equation about spatial location r and time t, used to describe sound pressure; S is the sound source term; equations (1) to (3) are used to generate the simulated sound pressure distribution p(r). m ,t).

3. The three-dimensional ultrasound tomography image reconstruction method based on time-domain full waveform inversion according to claim 2, characterized in that, In step two, the mean square error is used to measure the true sound pressure distribution p measured by the transducer. obs (r m ,t) and the simulated sound pressure distribution p(r) generated using the time-domain wave equation m The error ε between ,t) specifically includes: Where M represents the number of transmit-receive pairs, M = N(N-1), N represents the number of transducers; T represents the total duration of the transducer's signal detection.

4. The three-dimensional ultrasound tomography image reconstruction method based on time-domain full waveform inversion according to claim 3, characterized in that, In step two, the functional gradients of ε with respect to the sound velocity c and the sound absorption coefficient α0, based on the adjoint field method, are calculated, specifically including: Partial derivative of ε with respect to the speed of sound for: Solve equation (15) using the adjoint field method: The adjoint field wave equation is defined as: p * (r,Tt) is called p(r) m The accompanying field of ,t), r m The position of the transducer is indicated by T, and the total duration of the transducer's detection signal is indicated by T. Equation (15) can be written as: c(r) represents the speed of sound at spatial location r, where Γ1, Γ2, and Γ3 are all intermediate variables; In equation (17), Γ2 and Γ3 are much smaller than Γ1. Ignoring Γ2 and Γ3, equation (17) simplifies to: Solve equation (19) to calculate the functional gradient of ε with respect to the speed of sound c; The functional gradient of ε with respect to the acoustic absorption coefficient α0 is calculated based on the adjoint field method: The spatial gradient is solved in the Fourier domain: in, and These are the Fourier transform and the inverse Fourier transform, respectively, where k is the wave number.

5. The three-dimensional ultrasound tomography image reconstruction method based on time-domain full waveform inversion according to claim 3, characterized in that, In step two, the method of updating the distribution of sound velocity c and sound absorption coefficient α0 using the momentum gradient descent method to minimize the error ε specifically includes: i n+1 =θ n -l n v n (12) Among them, v n It is the velocity vector of the nth iteration; γ is the momentum parameter; It is the gradient parameter; θ n λ is the parameter updated in the nth iteration. n This represents the step size of the nth iteration; Indicates J relative to θ n The partial derivatives; J represents the loss function; The update formulas for sound velocity and sound absorption coefficient are as follows: c n This represents the speed of sound after the nth iteration. v represents the step size of the sound speed in the (n+1)th iteration. c,n-1 Let α represent the velocity vector in the momentum gradient descent method for the sound speed at the (n-1)th iteration. 0,n This represents the sound absorption coefficient after the nth iteration. v represents the step size of the acoustic absorption coefficient in the (n+1)th iteration. α,n-1 Let represent the velocity vector of the momentum gradient descent method when the acoustic absorption coefficient is obtained in the (n-1)th iteration.

Images