Underground engineering geological disaster source earthquake scattered wave high-resolution imaging method and system

Abstract

The present application relates to the technical field of seismic wave detection, and provides a high-resolution imaging method and system for seismic scattered waves of geological disaster sources of underground engineering, which comprises: obtaining multi-component seismic signals of a target area, establishing a velocity model of the target area through wave velocity inversion; decomposing the multi-component seismic signals into a primary scattered wave field and a multiple scattered wave field by iteratively solving a Malchenko focusing function; performing polarization analysis on the multi-component seismic signals to extract polarization characteristics; based on the velocity model, performing forward and reverse wave field continuation on the multi-component seismic signals, calculating the theoretical arrival time and propagation direction of different waveforms to obtain a propagation path; based on the polarization characteristics and the propagation path, optimizing the primary scattered wave field and the multiple scattered wave field; and performing joint imaging of the primary scattered wave field and the multiple scattered wave field by solving a least squares Malchenko imaging target function. The reliability and interpretation accuracy of scattered wave imaging are improved.

Description

Technical Field

[0001] This invention belongs to the field of seismic wave detection technology, and in particular relates to a high-resolution imaging method and system for seismic scattered waves of underground engineering geological disaster sources. Background Technology

[0002] The statements in this section are merely background information related to the present invention and do not necessarily constitute prior art.

[0003] Quaternary overburden generally exhibits low velocity and strong heterogeneity, making seismic wave fields susceptible to interference during propagation. High-frequency energy is strongly scattered and absorbed, limiting the usable frequency band. Simultaneously, surface waves and environmental noise dominate in shallow layers, further reducing the signal-to-noise ratio of reflected signals. Furthermore, in complex media with strong lateral variations or filled with scatterers, seismic waves undergo multiple mode conversions, generating numerous multiple scattered waves that overwhelm effective reflected signals. This significantly weakens the ability of conventional migration imaging methods to recover high wavenumber information and limits the imaging resolution of lower strata.

[0004] Scattered waves possess advantages such as rich high wavenumber components, wide-angle path propagation, and point-focusing imaging, enabling them to break through the limitations of conventional reflected wave imaging and achieve higher-resolution imaging of disaster sources. The Marchenko method is an effective approach for eliminating interlayer multiple scattering and accurately constructing subsurface wavefield information to achieve high-quality imaging. This method adaptively constructs a Green's function without multiple scattering using data and converts multiple scattering into effective illumination, overcoming the problems of traditional imaging relying on single scattering and the difficulty in accurately modeling multiple scattering paths in traditional scattering imaging methods. Existing Marchenko methods are mostly based on acoustic approximations, treating converted waves as noise. In complex media, this leads to problems such as inaccurate imaging position and poor wavefield focusing. Therefore, to more accurately reconstruct the wavefield and achieve imaging in complex media, it is necessary to establish an elastic wave Marchenko method that considers PS (P-S) mode conversion and elastic multiple scattering.

[0005] In summary, high-precision imaging of underground engineering geological hazard sources using scattered wave fields can leverage the advantages of scattered waves—rich high wavenumber components, wide-angle path propagation, and point-focusing imaging—to overcome the limitations of conventional reflected wave imaging and achieve high-precision imaging of complex underground media. However, due to the weak energy, complex wavefield, and unpredictable propagation paths of multiple scattering waves, there are currently no relevant methods or technologies. Therefore, using scattered waves to achieve high-resolution imaging of underground engineering geological hazard sources presents several challenges:

[0006] (1) Difficulty in separating single scattering and multiple scattering: Under complex medium conditions, scattered waves undergo multiple reflections, diffractions and mode conversions during propagation. Different orders of scattering have highly overlapping time, frequency and spatial characteristics, making it difficult to effectively distinguish between single scattering and multiple scattering. Existing separation methods are highly dependent on medium parameters and prior assumptions, and the separation effect is unstable in strongly non-uniform or high-contrast environments, which in turn affects the reliability of scattering imaging results.

[0007] (2) Existing scattering imaging methods have strong artifacts: Scattering imaging is often affected by limited aperture, insufficient sampling and multiple scattering interference. The imaging operator deviates from the real scattering process, which easily introduces false anomalies and strip artifacts into the imaging results. Especially in complex structures or strong noise backgrounds, the artifact energy may be comparable to the real scattering response, reducing the resolution and interpretation credibility of the imaging results.

[0008] (3) Scattered wave imaging methods are computationally intensive and slow: Scattered wave imaging usually involves multiple wavefield forward and inverse operations, frequency or time domain integration, and large-scale matrix operations, which require high computing resources. With the increasing demand for three-dimensional high-resolution imaging, the computational scale grows exponentially, resulting in low imaging efficiency and long computation cycle, which makes it difficult to meet the requirements of fast processing and real-time interpretation in engineering applications. Summary of the Invention

[0009] To address the technical problems mentioned above, this invention provides a high-resolution imaging method and system for seismic scattered waves from underground engineering geological hazard sources. Utilizing an accelerated three-dimensional multi-component elastic wave Marchenko method combined with multi-component polarization characteristic analysis, the observed wavefield is iteratively decomposed. Through focusing function constraints, propagation path matching, and polarization consistency discrimination, high-precision separation of primary and secondary scattering of P-waves, S-waves, and mode-converted waves is achieved. This effectively reduces the impact of multi-order scattering overlap on separation accuracy under complex medium conditions, improving the reliability and interpretation accuracy of scattered wave imaging.

[0010] To achieve the above objectives, the present invention adopts the following technical solution:

[0011] The first aspect of the present invention provides a high-resolution imaging method for seismic scattered waves from underground engineering geological hazard sources, comprising:

[0012] Acquire multi-component seismic signals of the target area and establish a velocity model of the target area through wave velocity inversion;

[0013] By iteratively solving the Marchenko focusing function, multi-component seismic signals are decomposed into primary and multiple scattered wave fields. Polarization analysis is performed on the multi-component seismic signals to extract polarization features. Based on the velocity model, forward and reverse wavefield extensions are performed on the multi-component seismic signals to calculate the theoretical arrival time and propagation direction of different waveforms and obtain the propagation path. Based on the polarization features and propagation path, the primary and multiple scattered wave fields are optimized.

[0014] By solving the least-squares Marchenko imaging objective function, joint imaging of the primary and multiple scattered wave fields is performed.

[0015] Furthermore, it also includes: accelerating the iterative solution of the Marchenko focusing function, and the acceleration adopts a combination of frequency segmentation and preconditional Krylov subspace iteration.

[0016] Furthermore, the preconditional Krylov subspace iteration is expressed as:

[0017] ;

[0018] in, For the frequency domain focusing function, Angular frequency, Indicates the number of frequency bands. For frequency segmentation and smoothing of time windows, m represents the velocity model. Here is the Krylov precondition matrix. For the Marchenko operator matrix, It serves as an intermediate auxiliary variable.

[0019] Furthermore, the least-squares Marchenko imaging objective function is: ;in, For the first scattering weight, For multiple scattering weights, Based on the first scattered wave field The first scattering imaging operator, For multiple scattering wave fields Multiple scattering imaging operator, The reflectivity model is updated based on the first scattered wave field. This is a reflectivity model updated from multiple scattered wave fields, where d represents the target matching number. For the regularization term of the reflectivity model of the first-order scattered wave field, For the regularization term of the reflectivity model of multiple scattered wave fields, and This is a regularization parameter that is adaptively adjusted based on the imaging residuals or convergence characteristics.

[0020] Furthermore, the multi-component seismic signal includes direct waves, reflected waves, scattered waves, and transmitted waves.

[0021] Furthermore, the polarization characteristics include the polarization direction, amplitude ratio, and phase relationship of the longitudinal wave, transverse wave, and longitudinal-transverse mode-converted wave.

[0022] A second aspect of the present invention provides a high-resolution imaging system for seismic scattered waves from underground engineering geological hazard sources, comprising:

[0023] The inversion module is configured to: acquire multi-component seismic signals of the target area and establish a velocity model of the target area through wave velocity inversion;

[0024] The decomposition module is configured to: decompose multi-component seismic signals into primary and multiple scattered wave fields by iteratively solving the Marchenko focusing function; perform polarization analysis on the multi-component seismic signals to extract polarization features; perform forward and reverse wave field extension on the multi-component seismic signals based on the velocity model, calculate the theoretical arrival time and propagation direction of different waveforms, and obtain the propagation path; and optimize the primary and multiple scattered wave fields based on the polarization features and propagation path.

[0025] The joint imaging module is configured to perform joint imaging of the primary and multiple scattered wave fields by solving the least-squares Marchenko imaging objective function.

[0026] Furthermore, it also includes an acceleration module configured to accelerate the iterative solution of the Marchenko focusing function by combining frequency segmentation with preconditioned Krylov subspace iteration.

[0027] Furthermore, the preconditional Krylov subspace iteration is expressed as:

[0028] ;

[0029] in, For the frequency domain focusing function, Angular frequency, Indicates the number of frequency bands. For frequency segmentation and smoothing of time windows, m represents the velocity model. Here is the Krylov precondition matrix. For the Marchenko operator matrix, It serves as an intermediate auxiliary variable.

[0030] Furthermore, the least-squares Marchenko imaging objective function is: ;in, For the first scattering weight, For multiple scattering weights, Based on the first scattered wave field The first scattering imaging operator, For multiple scattering wave fields Multiple scattering imaging operator, The reflectivity model is updated based on the first scattered wave field. This is a reflectivity model updated from multiple scattered wave fields, where d represents the target matching number. For the regularization term of the reflectivity model of the first-order scattered wave field, For the regularization term of the reflectivity model of multiple scattered wave fields, and This is a regularization parameter that is adaptively adjusted based on the imaging residuals or convergence characteristics.

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

[0032] This invention utilizes the accelerated three-dimensional multi-component elastic wave Marchenko method combined with multi-component polarization feature analysis to iteratively decompose the observed wavefield. Through focusing function constraints, propagation path matching, and polarization consistency discrimination, it achieves high-precision separation of primary and secondary scattering of longitudinal waves, transverse waves, and mode-converted waves. This effectively reduces the impact of multi-order scattering overlap on separation accuracy under complex medium conditions, and improves the reliability and interpretation accuracy of scattered wave imaging.

[0033] This invention introduces regularization constraints into the least-squares Marchenko imaging framework to suppress false anomalies and banded artifacts caused by finite aperture, undersampling, and multiple scattering. It adaptively adjusts the regularization parameters to enhance imaging stability and ensure that imaging results still have high resolution and high reliability in complex structures or high-noise backgrounds.

[0034] This invention employs an acceleration strategy to effectively reduce the computational load of three-dimensional high-resolution scattering imaging, supports rapid iterative solution of the focusing function, realizes joint imaging of single and multiple scattering, significantly shortens the computation cycle, and meets the engineering application requirements of rapid processing and real-time interpretation in underground engineering. Attached Figure Description

[0035] The accompanying drawings, which form part of this invention, are used to provide a further understanding of the invention. The illustrative embodiments of the invention and their descriptions are used to explain the invention and do not constitute an improper limitation of the invention.

[0036] Figure 1 This is a flowchart of a high-resolution imaging method for seismic scattered waves of underground engineering geological disaster sources according to Embodiment 1 of the present invention;

[0037] Figure 2 This is a flowchart of the joint imaging of primary and multiple scattered wave fields according to Embodiment 1 of the present invention;

[0038] Figure 3 This is a schematic diagram of the source and station setup according to Embodiment 1 of the present invention;

[0039] Figure 4 This is a schematic diagram of a high-resolution numerical simulation model for detecting seismic scattered waves from underground engineering geological hazard sources according to Embodiment 1 of the present invention;

[0040] Figure 5 This is a schematic diagram of the conventional reverse time migration result of Embodiment 1 of the present invention;

[0041] Figure 6 This is a schematic diagram of the migration results of a high-resolution imaging method for seismic scattered waves of underground engineering geological disaster sources according to Embodiment 1 of the present invention;

[0042] Among them, 1. the earthquake source, 2. the receiving station, and 3. the underground medium anomaly. Detailed Implementation

[0043] To make the objectives, technical solutions, and advantages of the embodiments of the present invention clearer, the technical solutions of the embodiments of the present invention will be clearly and completely described below with reference to the accompanying drawings.

[0044] It should be noted that the following detailed description is illustrative and intended to provide further explanation of the invention. Unless otherwise specified, all technical and scientific terms used herein have the same meaning as commonly understood by one of ordinary skill in the art to which this invention pertains.

[0045] Example 1

[0046] This embodiment provides a high-resolution imaging method for seismic scattered waves from underground engineering geological disaster sources.

[0047] This embodiment provides a high-resolution imaging method for seismic scattered waves of underground engineering geological disaster sources, which can obtain images of complex media and the distribution of adverse geological bodies within a specific depth range underground.

[0048] This embodiment provides a high-resolution imaging method for seismic scattered waves from underground engineering geological hazard sources, such as... Figure 1 As shown, it includes the following steps:

[0049] Step 1: Arrange the ground-based two-dimensional regular observation array and perform time synchronization and spatial positioning.

[0050] The specific process of arranging the ground-based two-dimensional regular observation array includes: on the surface above the detection area, the receiving stations 2 are arranged as regularly as possible according to the terrain to form a ground-based two-dimensional regular observation array, such as... Figure 3 As shown.

[0051] Among them, the ground two-dimensional regular observation array selects the channel spacing between the receiving stations 2 according to the resolution of the detection target, and selects the length and width of the array according to the detection depth.

[0052] Among them, receiving station 2 uses a three-component detector to simultaneously acquire longitudinal wave, transverse wave and mode-converted wave signals.

[0053] Step 2: The seismic source generates a signal at a specific excitation location, and the ground two-dimensional regular observation array synchronously acquires and stores multi-component seismic signals (three-dimensional multi-component seismic data), i.e., multi-component observation wavefield.

[0054] like Figure 3 As shown, source 1 can be selected as a pulse source or a controllable source according to the environment and target requirements. The controllable source must be able to achieve high-energy excitation of P-waves and S-waves, and different sources must meet the detection depth requirements.

[0055] The multi-component seismic signal includes information on direct waves, reflected waves, scattered waves, and transmitted waves.

[0056] Step 3: Based on the acquired multi-component seismic signals, perform preprocessing to obtain high signal-to-noise ratio seismic data.

[0057] The preprocessing of the acquired multi-component seismic signals includes the following specific steps:

[0058] Denoising, DC component removal, time correction, and gather quality control are performed on multi-component seismic signals.

[0059] Based on the noise type, filtering is performed to obtain filtered seismic data; the filtering includes frequency domain filtering, FK filtering, τp filtering and / or predictive deconvolution filtering;

[0060] Preliminary P-wave and S-wave separation was performed on the filtered seismic data.

[0061] Step 4: Based on high signal-to-noise ratio seismic data, establish a velocity model for the target area using the wave velocity inversion method.

[0062] Among them, wave velocity inversion methods include travel time inversion, waveform inversion, or joint inversion methods.

[0063] The velocity model includes both the P-wave velocity model (P-wave velocity profile) and the S-wave velocity model (S-wave velocity profile).

[0064] Step 5: Separation of the scattered wave field in the target area: Under the constraint of the velocity model, the multi-component observed wave field is decomposed by the Marchenko method of three-dimensional multi-component elastic waves to separate the primary scattered wave field and the multiple scattered wave field containing the P-wave and S-wave mode conversion.

[0065] The Marchenko method for three-dimensional multi-component elastic waves specifically includes:

[0066] Step 501: Based on the velocity model, according to the target area, uniformly set the focal points, take each focal point as the source position and the actual position of the receiving station 2 as the receiving point position, perform forward modeling to obtain the transmitted wave field at each focal point, where the direct transmitted wave is the initial value of the iterative solution of the Marchenko focusing function.

[0067] Step 502: Using the Marchenko focusing function, iteratively decompose the multi-component observed wavefield to obtain the primary and secondary scattered wavefields.

[0068] The focusing function includes focusing conditions for longitudinal waves, transverse waves, and mode-converting waves, and is used to characterize the propagation differences between single scattering and multiple scattering.

[0069] Specifically, the formula for calculating the focusing function is:

[0070] ;

[0071] ;

[0072] ;

[0073] in, These are the coordinates of the multi-component seismic signal receiving point. These are the coordinates of the earthquake's source excitation point. Angular frequency, For the focusing function, These are the downlink and uplink focus functions, respectively. It is a multi-component reflection response matrix. For time window functions, It is a unit vector. Let R be the conjugate vector of R.

[0074] In this embodiment, the separation of the primary and secondary scattered wave fields requires the calculation of the Green's function:

[0075] ;

[0076] ;

[0077] ;

[0078] Where G is the reconstructed Green's function, For the upward Green's function, For the descending Green's function, It is the conjugate function of the downlink focusing function.

[0079] The primary scattered wave field and the multiple scattered wave field are separated by the following formula:

[0080] ;

[0081] ;

[0082] in, For multiple scattered wave fields, This is a primary scattered wave field.

[0083] Step 503: Utilize multi-component polarization characteristic results to improve the separation accuracy of mode-converted scattered waves. The specific process includes:

[0084] Polarization analysis was performed on the multi-component observed wavefield to extract the polarization direction, amplitude ratio, and phase relationship of the longitudinal wave (P-wave), transverse wave (S-wave), and the P-wave / S-wave mode-converted wave, thus obtaining the polarization characteristics.

[0085] Based on the existing three-dimensional velocity model (i.e., the velocity model obtained in step 4), the multi-component observed wavefield is extended in the forward and reverse directions, and the theoretical arrival time and propagation direction information corresponding to different waveforms are calculated to obtain the propagation path.

[0086] Based on polarization characteristics and propagation paths, constraints and optimizations are performed on the primary and secondary scattering wave fields.

[0087] The polarization analysis employs a sliding time window approach to improve the identification capability of non-stationary scattered waves within the time window. The formula for calculating the polarization covariance matrix is: The sliding time window divides a period of time into multiple sub-windows, where t0 is the start time of each sub-window. This refers to the length of the sliding time window. For example, if a 10-second time interval is divided into 5 sub-windows, and the length of each sliding time window is 2 seconds, then t0 would be 0, 2, 4, 6, and so on. , The covariance matrix represents a multi-component seismic signal. , , and These represent the seismic signal components along the x, y, and z axes of the three-dimensional coordinate system; eigenvalue decomposition of the covariance matrix: Where E is the eigenvector matrix, The eigenvalue vector is the largest eigenvalue, which corresponds to the main polarization direction. The eigenvalue ratio is used to distinguish between P-wave / S-wave / P-wave / S-wave mode-converting waves.

[0088] In this embodiment, the Marchenko method for three-dimensional multi-component elastic waves is accelerated, that is, the process of iteratively solving the Marchenko focusing function is accelerated.

[0089] Among them, the accelerated computation strategies for the Marchenko method for three-dimensional multi-component elastic waves include one or more of the following: frequency segmentation, parallel computation, low-rank approximation, or graphics processing unit (GPU) acceleration.

[0090] Specifically, to accelerate the Marchenko method for three-dimensional multi-component elastic waves, a combination of frequency segmentation and preconditional Krylov subspace iteration is employed, with GPU acceleration used when necessary. The preconditional Krylov subspace iteration represents the Marchenko focusing function solution process as a linear least squares problem, and iteratively solves it using preconditional CGLS (conjugate gradient least squares) or PCG (preprocessed conjugate gradient method). In each iteration, historical residual information is comprehensively utilized to construct the optimal search direction, thereby significantly reducing the number of iterations, improving convergence stability, and lowering the computational complexity of three-dimensional multi-component scattering imaging.

[0091] The preconditional Krylov subspace iteration is expressed as:

[0092] ;

[0093] in, For the frequency domain focusing function, Angular frequency, Indicates the number of frequency bands. For frequency segmentation and smoothing time windows, , , , It is a multi-component reflection response matrix. for The conjugate vector, It is the identity matrix. Let m be the focusing function and the time window function, representing the velocity model. This is the Krylov precondition matrix. The Marchenko operator matrix, As an intermediate auxiliary variable, This is the upward focusing function.

[0094] Step 6: Construct a least-squares Marchenko imaging objective function, perform least-squares reverse time migration to achieve joint imaging of the primary and multiple scattered wave fields of the target region, thereby identifying features such as... Figure 3 The subsurface anomaly 3 shown is the source of underground engineering geological hazards, such as... Figure 4 , Figure 5 and Figure 6 As shown.

[0095] Introducing regularization constraints during joint imaging to improve imaging stability and suppress artifacts involves the following steps:

[0096] A regularization term is introduced into the least squares Marchenko imaging objective function to impose smoothness, sparsity, or energy constraints on the imaging results.

[0097] By adjusting the regularization weight parameter, a balance is struck between the data fitting term and the regularization term, suppressing imaging instability caused by noise, multiple scattering, and finite aperture.

[0098] During iterative imaging, the regularization parameter is adaptively adjusted based on the imaging residual or convergence characteristics to reduce imaging artifacts and improve the reliability of imaging results.

[0099] The least squares Marchenko imaging objective function is: ;in, For the first scattering weight, For multiple scattering weights, Based on the first scattered wave field The first scattering imaging operator, For multiple scattering wave fields Multiple scattering imaging operator, This is a reflectivity model updated from the first scattered wave field, including P-wave and S-wave reflectivity models. This is a reflectivity model updated from multiple scattered wave fields, including P-wave and S-wave reflectivity models. d represents the target matching data, here consisting of the separated primary and multiple scattered wave fields. For the regularization term of the reflectivity model of the first-order scattered wave field, This is the regularization term for the reflectivity model of the multiple scattered wave field, which can take the form of Tikhonov regularization, first-order smoothing (gradient) regularization, or sparse regularization. and This is the regularization parameter.

[0100] like Figure 2 As shown, the specific steps for joint imaging of primary and multiple scattered wave fields include:

[0101] (1) Input the primary and multiple scattered wave fields into the reflection coefficient model for imaging to obtain the initial reflection coefficient distribution, i.e. and ;

[0102] (2) By combining the velocity model and performing forward modeling, synthetic seismic data are obtained. ;

[0103] (3) The Marchenko scattering wave field separation algorithm based on frequency segmentation and precondition Krylov subspace linear iteration acceleration algorithm is combined to separate the synthetic seismic data and obtain the synthetic primary scattering wave field and multiple scattering wave field.

[0104] (4) Calculate the residuals between the synthesized and the actual observed primary and multiple scattering wavefield data to measure the fit between the current model and the observed data: .

[0105] (5) Substitute the data residuals into the gradient calculation formula to obtain The gradients are respectively:

[0106] ;

[0107] ;

[0108] in, for The conjugate operator, for The conjugate operator, The gradient of the single scattering model. The gradient of the multiple scattering model. , These are the regularization gradients for single-pass and multiple-pass scattering models, used for smoothing, sparsity, or energy constraints.

[0109] (6) Update the single scattering and multiple scattering models based on gradient information:

[0110] ;

[0111] ;

[0112] in, and This is the step size coefficient, which can be determined through line search or by preset step size; For the kth generation , For the kth generation , For the (k-1)th generation , For the (k-1)th generation , For the (k-1)th generation , For the (k-1)th generation .

[0113] (7) Determine whether the convergence condition is met. If yes, output the least square Marchenko imaging result; otherwise, return to step (2).

[0114] This embodiment provides a high-resolution imaging method for seismic scattered waves of underground engineering geological hazard sources. Based on the acquisition of high-energy three-component seismic data by a regular array on the surface, high signal-to-noise ratio seismic data is obtained by preprocessing the information of the acquired direct waves, reflected waves, scattered waves and transmitted waves. Then, a smooth velocity model of the study area is established by using travel-time inversion, waveform inversion or joint inversion methods. Under the constraint of the velocity model, the observed wavefield is decomposed based on the accelerated three-dimensional multi-component elastic wave Marchenko method, and the least squares Marchenko imaging operator is constructed to realize the joint imaging of single scattering and multiple scattering of the target area.

[0115] This embodiment provides a high-resolution imaging method for seismic scattered waves from underground engineering geological hazard sources, achieving high-precision separation of single and multiple scattering. It utilizes the accelerated three-dimensional multi-component elastic wave Marchenko method combined with multi-component polarization characteristic analysis to iteratively decompose the observed wavefield. Through focusing function constraints, propagation path matching, and polarization consistency discrimination, it achieves high-precision separation of single and multiple scattering of P-waves, S-waves, and mode-converted waves. This effectively reduces the impact of multi-order scattering overlap on separation accuracy under complex medium conditions, improving the reliability and interpretation accuracy of scattered wave imaging.

[0116] This embodiment provides a high-resolution imaging method for seismic scattered waves from underground engineering geological hazard sources, which significantly reduces scattering imaging artifacts. Regularization constraints (smoothness, sparsity, or energy constraints) are introduced into the least-squares Marchenko imaging framework to suppress false anomalies and banded artifacts caused by finite aperture, undersampling, and multiple scattering. Adaptive adjustment of the regularization parameters enhances imaging stability, ensuring high resolution and high reliability of imaging results even in complex structures or high-noise backgrounds.

[0117] This embodiment provides a high-resolution imaging method for seismic scattered waves from underground engineering geological hazard sources, which improves computational efficiency and imaging speed. It employs multiple acceleration strategies, including frequency segmentation, parallel computing, low-rank approximation, and GPU acceleration, effectively reducing the computational load of three-dimensional high-resolution scattering imaging. It supports rapid iterative solution of the focusing function, enabling joint imaging of single and multiple scattering waves, significantly shortening the computation cycle, and meeting the engineering application requirements for rapid processing and real-time interpretation of underground engineering data.

[0118] This embodiment provides a high-resolution imaging method for seismic scattered waves of underground engineering geological hazard sources, which enhances the ability to identify complex media and adverse geological bodies. Through high-resolution imaging of three-dimensional multi-component wavefields, the precise spatial distribution of complex media and adverse geological bodies within a specific depth range underground can be obtained; the combined use of P-waves, S-waves, and mode-converted waves improves the detection capability of small scatterers or weak anomalies under complex geological conditions.

[0119] This embodiment provides a high-resolution imaging method for seismic scattered waves from underground engineering geological hazard sources, balancing imaging accuracy and engineering applicability. It integrates observation array deployment, data preprocessing, velocity modeling, wavefield decomposition, and regularized imaging, achieving integrated processing from raw seismic data to high-resolution scattered imaging. It is highly adaptable, allowing flexible selection of source type, array size, and processing parameters according to different engineering scenarios, thus balancing imaging accuracy and engineering practicality.

[0120] Example 2

[0121] This embodiment provides a high-resolution imaging system for seismic scattered waves from underground engineering geological hazard sources, comprising:

[0122] The inversion module is configured to: acquire multi-component seismic signals of the target area and establish a velocity model of the target area through wave velocity inversion;

[0123] The decomposition module is configured to: decompose multi-component seismic signals into primary and multiple scattered wave fields by iteratively solving the Marchenko focusing function; perform polarization analysis on the multi-component seismic signals to extract polarization features; perform forward and reverse wave field extension on the multi-component seismic signals based on the velocity model, calculate the theoretical arrival time and propagation direction of different waveforms, and obtain the propagation path; and optimize the primary and multiple scattered wave fields based on the polarization features and propagation path.

[0124] The joint imaging module is configured to perform joint imaging of the primary and multiple scattered wave fields by solving the least-squares Marchenko imaging objective function.

[0125] Furthermore, it also includes an acceleration module configured to accelerate the iterative solution of the Marchenko focusing function by combining frequency segmentation with preconditioned Krylov subspace iteration.

[0126] Furthermore, the preconditional Krylov subspace iteration is expressed as:

[0127] ;

[0128] in, For the frequency domain focusing function, Angular frequency, Indicates the number of frequency bands. For frequency segmentation and smoothing of time windows, m represents the velocity model. Here is the Krylov precondition matrix. For the Marchenko operator matrix, It serves as an intermediate auxiliary variable.

[0129] Furthermore, the least-squares Marchenko imaging objective function is: ;in, For the first scattering weight, For multiple scattering weights, Based on the first scattered wave field The first scattering imaging operator, For multiple scattering wave fields Multiple scattering imaging operator, The reflectivity model is updated based on the first scattered wave field. This is a reflectivity model updated from multiple scattered wave fields, where d represents the target matching number. For the regularization term of the reflectivity model of the first-order scattered wave field, For the regularization term of the reflectivity model of multiple scattered wave fields, and This is a regularization parameter that is adaptively adjusted based on the imaging residuals or convergence characteristics.

[0130] It should be noted that each module in this embodiment corresponds one-to-one with each step in Embodiment 1, and their specific implementation processes are the same, so they will not be repeated here.

[0131] The above description is merely a preferred embodiment of the present invention and is not intended to limit the invention. Various modifications and variations can be made to the present invention by those skilled in the art. Any modifications, equivalent substitutions, improvements, etc., made within the spirit and principles of the present invention should be included within the scope of protection of the present invention.

Claims

1. A method for high-resolution imaging of seismic scattered waves from a source of geohazards in underground engineering, characterized in that, The method comprises the following steps: obtaining a multi-component seismic signal of a target area, and establishing a velocity model of the target area through wave velocity inversion; decomposing the multi-component seismic signal into a primary scattered wave field and a multiple scattered wave field by iteratively solving a Malchenko focusing function; performing polarization analysis on the multi-component seismic signal to extract polarization characteristics; based on the velocity model, performing forward and reverse wave field continuation on the multi-component seismic signal, calculating theoretical travel times and propagation directions of different waveforms, and obtaining propagation paths; and based on the polarization characteristics and the propagation paths, optimizing the primary scattered wave field and the multiple scattered wave field; performing joint imaging of the primary scattered wave field and the multiple scattered wave field by solving a least squares Malchenko imaging target function. The least-squares Marchenko imaging objective function is: ; wherein, is a single-scattering weight, is a multiple-scattering weight, is a single-scattering imaging operator based on the single-scattering wavefield , is a multiple-scattering imaging operator based on the multiple-scattering wavefield , is a reflectivity model updated by the single-scattering wavefield, is a reflectivity model updated by the multiple-scattering wavefield, d denotes the number of target matches, is a regularization term for the single-scattering wavefield reflectivity model, is a regularization term for the multiple-scattering wavefield reflectivity model, and are regularization parameters, which are adaptively adjusted according to imaging residuals or convergence characteristics. The polarization characteristics include polarization directions, amplitude ratios and phase relationships of P waves, S waves and P-S mode conversion waves.

2. The method according to claim 1, wherein, The method further comprises the following steps: accelerating the process of iteratively solving the Malchenko focusing function, and the acceleration adopts a combination of frequency segmentation and preconditioned Krylov subspace iteration.

3. The method according to claim 2, wherein the method is characterized by, The preconditioned Krylov subspace iteration is represented as: ； wherein is a frequency domain focusing function, is an angular frequency, denotes the number of divided frequency bands, is a frequency section smoothing time window, m denotes a velocity model, is a Krilov preconditioning matrix, is a Martenchko operator matrix, is an intermediate auxiliary variable.

4. The method of claim 1, wherein the method is characterized by, The multi-component seismic signal includes direct waves, reflected waves, scattered waves and transmitted waves.

5. An underground engineering geological disaster source seismic scattered wave high-resolution imaging system, characterized in that, The method comprises the following steps: an inversion module configured to obtain a multi-component seismic signal of a target area, and establish a velocity model of the target area through wave velocity inversion; a decomposition module configured to decompose the multi-component seismic signal into a primary scattered wave field and a multiple scattered wave field by iteratively solving a Malchenko focusing function; performing polarization analysis on the multi-component seismic signal to extract polarization characteristics; based on the velocity model, performing forward and reverse wave field continuation on the multi-component seismic signal, calculating theoretical travel times and propagation directions of different waveforms, and obtaining propagation paths; and based on the polarization characteristics and the propagation paths, optimizing the primary scattered wave field and the multiple scattered wave field; a joint imaging module configured to perform joint imaging of the primary scattered wave field and the multiple scattered wave field by solving a least squares Malchenko imaging target function. The least square Marchenko imaging objective function is: ; wherein, is a single scattering weight, is a multiple scattering weight, is a single scattering imaging operator based on the single scattering wavefield , is a multiple scattering imaging operator based on the multiple scattering wavefield , is a reflectivity model updated by the single scattering wavefield, is a reflectivity model updated by the multiple scattering wavefield, d denotes the number of target matches, is a regularization term of the single scattering wavefield reflectivity model, is a regularization term of the multiple scattering wavefield reflectivity model, and are regularization parameters, which are adaptively adjusted according to imaging residuals or convergence characteristics. The polarization characteristics include polarization directions, amplitude ratios and phase relationships of P waves, S waves and P-S mode conversion waves.

6. The system for high resolution imaging of seismic scattered waves from geological disaster sources in underground structures according to claim 5, characterized in that, The method further comprises an acceleration module configured to accelerate the process of iteratively solving the Malchenko focusing function, and the acceleration adopts a combination of frequency segmentation and preconditioned Krylov subspace iteration.

7. The system of claim 6, wherein the system is configured to perform the steps of: determining a plurality of seismic scattering wave data sets, each of the plurality of seismic scattering wave data sets corresponding to a different one of the plurality of different frequencies; and determining a plurality of scattering wave images, each of the plurality of scattering wave images corresponding to a different one of the plurality of different frequencies. The preconditioned Krylov subspace iteration is represented as: ； wherein is a frequency domain focusing function, is an angular frequency, denotes the number of divided frequency bands, is a frequency section smoothing time window, m denotes a velocity model, is a Krilov preconditioning matrix, is a Martenchko operator matrix, is an intermediate auxiliary variable.

Images