Sigmoid modulation phase correlation surface wave full waveform inversion method and system

Abstract

The invention provides a Sigmoid modulation phase correlation surface wave full waveform inversion method and system. The method comprises the following steps: acquiring a surface wave seismic record and constructing an initial S wave velocity model; performing forward modeling to generate a synthetic record; calculating instantaneous phase correlation of observation and synthetic records to construct an amplitude-independent target function; modulating the target function by using a Sigmoid function, and mapping the value of the target function to a predetermined interval; and calculating a gradient based on the modulated target function and iteratively updating the speed model until a termination condition is met, and outputting a final model. According to the method, the smooth saturation and bounded output characteristics of the Sigmoid function are utilized, the target function response is converted into a self-adaptive gradient scaling mechanism, a self-stabilization optimizer is introduced for the inversion process, and the convergence domain is widened by applying smooth constraint, so that the dependence on the initial model precision is reduced, and the first success rate and the result reliability of inversion under complex conditions are improved.

Description

Technical Field

[0001] This application relates to the field of seismic wave imaging technology, specifically to a Sigmoid-modulated phase-correlated full waveform inversion method and system for surface waves such as Rayleigh waves and Love waves. Background Technology

[0002] Surface wave full waveform inversion has become an important technique for obtaining high-resolution shear wave velocity structures at shallow surfaces. This technique fully utilizes the high sensitivity of surface waves to subsurface shear wave velocity distribution. By iteratively updating the subsurface medium model, it achieves optimal waveform characteristics matching between numerically simulated synthetic seismic records and field observation records, thereby inverting and deducing the physical parameters and structure of the subsurface medium. It has broad application prospects in the field of high-resolution geophysical imaging.

[0003] However, traditional least-squares norm-based full-waveform inversion methods still face two major challenges under complex shallow geological conditions: periodic jumps and amplitude errors. The period jump problem stems from the strong dispersion effect of surface waves, causing the phase of the seismic wave field to change rapidly with the propagation distance. When the initial velocity model is inaccurate, a time shift error of more than half a period may occur between the synthesized and observed waveforms, causing the inversion to get trapped in local minima and unable to converge to the true model; Amplitude error is caused by inherent uncertainties in field data acquisition, such as inconsistent coupling between the source and detector, differences in excitation energy, environmental noise interference, and inaccurate geometric diffusion and attenuation correction due to simplification from 3D to 2D. These errors reduce amplitude reliability and make traditional inversion results that rely on amplitude matching unstable or even distorted.

[0004] To address these bottlenecks, constructing an objective function insensitive to amplitude variations has become a key approach for robust full-waveform inversion. In recent years, various technical approaches have emerged in academia: for example, incorporating the minimization of the instantaneous phase objective function into the full-waveform inversion framework; proposing inversion methods based on the global correlation norm to reduce dependence on amplitude; utilizing instantaneous phase information to achieve time-domain full-waveform inversion and alleviate local minima; and systematically analyzing the advantages and limitations of phase information in adjoint state tomography. Of particular note is the development, inspired by the phase-weighted superposition method, of an amplitude-independent instantaneous phase correlation metric as the objective function, with its potential to address amplitude errors verified in the characterization of shallow surface structures at archaeological sites.

[0005] Although phase-based methods enhance amplitude robustness, the shape of their objective function can still produce spurious local minima under complex wavefield interference and strong noise, requiring high accuracy of the initial model and iterative stability. Therefore, while retaining the advantages of phase metrics, further optimizing the mathematical shape of the objective function to suppress spurious extrema and improve convergence stability remains a core scientific problem in current surface wave full waveform inversion.

[0006] It is worth noting that the sigmoid function is widely used in fields such as deep learning because its bounded output, smooth monotonicity, and saturation properties effectively suppress outliers and stabilize gradient propagation. Previous research has utilized the sigmoid function to define iterative weights in weighted envelope correlation inversion, achieving a smooth transition from background model to detailed updates. This inspires a new approach that combines the modulation capability of the sigmoid function with the robustness of instantaneous phase correlation, potentially leading to the construction of a novel objective function framework suitable for complex noisy environments.

[0007] In summary, existing full-waveform inversion techniques suffer from drawbacks such as sensitivity to amplitude errors, susceptibility to periodic jumps, and poor objective function shape. Developing a method that balances wavefield kinematic phase matching accuracy with robustness of the inversion process is of significant scientific and engineering value for promoting the practical application of surface wave technology in near-surface imaging. Summary of the Invention

[0008] This application provides a method and system for inverting the full waveform of a surface wave with Sigmoid modulation and phase correlation, which can solve the technical problems in the prior art where the full waveform inversion of surface waves is easily affected by period jumps and amplitude errors, and the instability of the inversion is caused by the poor shape of the objective function.

[0009] In a first aspect, this application provides a method for inverting the full waveform of a surface wave with Sigmoid modulation phase correlation, comprising the following steps: Acquire field surface wave seismic records and construct an initial subsurface S-wave velocity model; Forward modeling was performed based on the initial S-wave velocity model to generate synthetic seismic records; Calculate the instantaneous phase correlation between the observation record and the synthetic record, and construct a phase correlation objective function; The phase correlation target function is modulated using the Sigmoid function, and the target function value is mapped to a predetermined interval to obtain the modulated target function. Based on the modulated objective function, its gradient with respect to the underground model parameters is calculated, and the S-wave velocity model is iteratively updated using this gradient. When the preset inversion termination condition is met, the final inverted S-wave velocity structure model is output. A further optimization scheme involves constructing the instantaneous phase correlation based on the complex seismic traces of observational and synthetic data.

[0010] A further optimization scheme is as follows: the instantaneous phase correlation objective function Represented as: in, and To correspond to the observation data respectively and synthetic data instantaneous amplitude, This represents the summation over all sources and detectors. To record time.

[0011] A further optimization is that the Sigmoid function is the standard Sigmoid function, and its expression is: The modulated objective function Represented as: in, Let be the objective function for the instantaneous phase correlation. A further optimization scheme involves using the chain rule when calculating the gradient with respect to the model parameters based on the modulated objective function. The gradient expression is as follows: Among them, among them, The gradient of the modulated objective function; This represents the gradient of the original objective function.

[0012] A further optimization involves preprocessing the gradient before iteratively updating the S-wave velocity model, including: Apply a spatial taper function to the single-shot gradient; Apply a Gaussian filter to the superimposed gradients. A further optimization scheme involves adopting a multi-scale inversion strategy, including: The inversion process starts from the low-frequency band and gradually expands to higher frequency bands. The velocity model obtained from the previous frequency band inversion is used as the initial model for the next frequency band inversion.

[0013] A further optimization scheme is that the method is used to process field surface wave seismic data containing amplitude errors or random noise. Secondly, this application provides a sigmoid-modulated phase-correlated surface wave full waveform inversion system, comprising: The data acquisition module is used to acquire field surface wave seismic records and construct an initial subsurface S-wave velocity model; The forward modeling module is communicatively connected to the data acquisition module and is used to perform forward modeling based on the initial S-wave velocity model to generate synthetic seismic records. The objective function construction module is communicatively connected to the forward modeling module and is used to calculate the instantaneous phase correlation between the observation record and the synthetic record, and to construct the phase correlation objective function. The gradient calculation and inversion update module is communicatively connected to the objective function construction module. It is used to modulate the phase correlation objective function using the Sigmoid function, map the objective function value to a predetermined interval, and obtain the modulated objective function. The inversion control module is communicatively connected to the gradient calculation and inversion update module. It is used to calculate the gradient of the modulated objective function with respect to the underground model parameters, and to iteratively update the S-wave velocity model using the gradient. When the preset inversion termination condition is met, the final inverted S-wave velocity structure model is output. Thirdly, this application provides a computer-readable storage medium storing a Sigmoid-modulated phase-dependent surface wave full waveform inversion program, wherein when the Sigmoid-modulated phase-dependent surface wave full waveform inversion program is executed by a processor, the steps of the Sigmoid-modulated phase-dependent surface wave full waveform inversion method as described above are implemented.

[0014] The beneficial effects of the technical solutions provided in this application include at least the following: The instantaneous phase correlation objective function is nonlinearly modulated using the sigmoid function, leveraging its smooth saturation and bounded output characteristics to transform the objective function response into an adaptive gradient scaling mechanism. This mechanism automatically decays the gradient update step size when the phase difference is extremely large or small, effectively suppressing anomalous updates caused by periodic jumps and noise interference; while providing stronger gradient excitation when the phase difference is moderate. This intelligent modulation introduces a self-stabilizing optimizer into the inversion process, significantly widening the convergence region by applying smoothing constraints, thereby reducing the dependence on the accuracy of the initial model and improving the first-time success rate and reliability of inversion under complex wavefields and low signal-to-noise ratio conditions. Attached Figure Description

[0015] Figure 1 This is a flowchart of the Sigmoid modulation phase-dependent surface wave full waveform inversion method provided in this application embodiment; Figure 2 This is a schematic diagram illustrating the characteristics of the Sigmoid function provided in an embodiment of this application; Figure 3 A schematic diagram of the observed Ricker wavelet and its time-shifted replica for objective function characteristic analysis provided in this application embodiment; Figure 4A comparison graph of the relationship curves of time shift changes of different objective functions provided in the embodiments of this application; Figure 5 The morphological diagram of the cost function for sigmoid function parameter optimization provided in the embodiments of this application; Figure 6 This is a schematic diagram of the real S-wave velocity model for the synthetic experiment provided in the embodiments of this application; Figure 7 A schematic diagram of common shot point gathers based on a real model simulation provided for embodiments of this application; Figure 8 A comparison diagram of the inversion S-wave velocity models of the two objective functions provided in the embodiments of this application under clean data; Figure 9 Comparison of S-wave velocity inversion results for clean and noisy data provided by the three full waveform inversion methods in the embodiments of this application; Figure 10 This is a comparison chart of the fitting effects between the observation gather and the synthesized gather from the inversion model provided in the embodiments of this application; Figure 11 A schematic diagram of the initial S-wave velocity model and its smoothed version for a field research area provided in this application embodiment; Figure 12 This is a schematic diagram of the field-measured common shot point gather provided in the embodiments of this application; Figure 13 Comparison of the three objective functions provided in this application for the inversion S-wave velocity models under field data; Figure 14 A comparison of the fitting effects between field observation gathers and synthesized gathers from inversion models provided in the embodiments of this application. Detailed Implementation

[0016] To make the objectives, technical solutions, and advantages of the present invention clearer, the present invention will be further described in detail below with reference to the embodiments and accompanying drawings. The illustrative embodiments and descriptions of the present invention are only used to explain the present invention and are not intended to limit the present invention.

[0017] First, some of the technical terms used in this application will be explained to help those skilled in the art understand this application.

[0018] FWI: Full-waveform inversion; Vs: S-wave velocity; IPC: instantaneous-phase coherency; SH: SH-wave, SH wave (shear horizontal wave); PWS: phase-weighted stacks; GCN: global correlation norm; LSTM: Long Short-Term Memory; GRU: gated recurrent unit; WECI: Weighted envelope correlation inversion; C-PML: convolutional perfectly matched layer; RPE: relative percentage error; FATT: First-arrival traveltime tomography.

[0019] First, some of the technical terms used in this application will be explained to help those skilled in the art understand this application.

[0020] Firstly, such as Figure 1 As shown, this application provides a method for full waveform inversion of sigmoid-modulated phase-correlated surface waves, including the following steps: Step S1: Obtain surface wave seismic records and construct an initial subsurface S-wave velocity model; Step S2: Perform forward modeling based on the initial S-wave velocity model to generate a synthetic seismic record; Step S3: Calculate the instantaneous phase correlation between the observation record and the synthetic record, and construct the phase correlation objective function; Step S4: Modulate the phase correlation target function using the Sigmoid function, map the target function value to a predetermined interval, and obtain the modulated target function; Step S5: Based on the modulated objective function, calculate its gradient with respect to the underground model parameters, and use the gradient to iteratively update the S-wave velocity model; Step S6: When the preset inversion termination condition is met, output the final inverted S-wave velocity structure model.

[0021] The Sigmoid-modulated phase-correlated surface wave full-waveform inversion method provided in this application effectively overcomes the problems of traditional full-waveform inversion, such as sensitivity to amplitude errors and susceptibility to local extrema (cycle skipping), by introducing an instantaneous phase correlation objective function and combining it with Sigmoid function modulation. This method constructs an amplitude-independent phase correlation objective function, focusing the inversion on phase information, which is more critical to the velocity model. Furthermore, by utilizing the nonlinear mapping characteristics of the Sigmoid function, the objective function value is compressed to a stable range, enhancing the robustness of the inversion process in cases of inaccurate initial models or missing low-frequency data. In particular, it improves the imaging resolution and inversion accuracy for shallowly buried small-scale geological bodies or laterally inhomogeneous structures. In this application, S-wave is a concept in seismology, also known as shear wave or transverse wave.

[0022] In one embodiment, step S1: acquiring field surface wave seismic records and constructing an initial subsurface S-wave velocity model specifically includes the following steps: Step S11: Acquire surface wave seismic data recorded by geophones deployed in the work area and preprocess the data. Preprocessing includes: bad path removal, first arrival removal, bandpass filtering (preserving the effective frequency band of surface waves, such as 5-30Hz), amplitude recovery to compensate for geometric spread, etc., and finally obtain the preprocessed surface wave seismic record; Step S12: Based on the preprocessed field surface wave seismic data, construct a one-dimensional or two-dimensional initial S-wave velocity model (i.e., Vs model) for the work area; this initial S-wave velocity model can be obtained through one or more of the following methods: A. Obtain the shallow Vs structure by inverting the surface wave dispersion curve; B. Provide a layered model based on regional geological data or well logging data; C. Set a simple model with uniform or smooth gradient change as the starting point. This embodiment effectively eliminates noise interference and enhances the signal-to-noise ratio of the effective surface wave signal through preprocessing, providing high-quality data input for inversion. Meanwhile, the flexible and diverse initial model construction strategies (such as using dispersion curve inversion, combining prior geological information, or setting a smoothing model) significantly improve the convergence stability and adaptability of the inversion process, especially reducing the over-reliance on the initial model, and providing a robust starting point for subsequent Sigmoid modulation phase inversion.

[0023] In one embodiment, step S2: performing forward modeling based on the initial S-wave velocity model to generate a synthetic seismic record specifically includes the following steps: Step S21: Based on the initial S-wave velocity model, set the forward modeling parameters, including: spatial grid discretization, time step, source wavelet (such as the Ricker wavelet with a dominant frequency of 30Hz), and absorbing boundary conditions (such as convolutionally perfectly matched layer C-PML) to simulate the propagation of the wave in an infinite medium. Step S22: Using numerical methods such as the finite difference method or the spectral element method, based on the set forward modeling parameters, solve the elastic wave equation, simulate the propagation process of seismic waves in the given model, calculate the synthetic seismic waveform at each receiver point, and obtain the original wavefield data. Step S23: Extract and format the generated raw wavefield data, extract the wavefield values ​​of all detector locations and all time sampling points, and form a synthetic seismic record that is consistent with the field observation data in terms of trace number, point number, and sampling rate. This embodiment, by setting forward modeling parameters that include absorbing boundary conditions (such as C-PML) and using numerical methods such as finite difference to solve the elastic wave equation, can accurately simulate the propagation process of seismic waves in complex media and effectively suppress false boundary reflections. At the same time, by generating synthetic seismic records that are completely consistent with the format of field observation data, a high-fidelity waveform comparison basis is provided for subsequent inversion, thereby significantly improving the accuracy and numerical stability of forward modeling simulation in full waveform inversion.

[0024] In one embodiment, step S3: calculating the instantaneous phase correlation between the observed record and the synthesized record, and constructing the phase correlation objective function, specifically includes the following steps: Step S31: Process the observation data channels and synthetic data channels Perform complex seismic trace analysis separately and calculate the Hilbert transform for each trace. By seismic channels As the real part, its Hilbert transformation As a virtual part, construct complex seismic channels ;in The imaginary unit, Represents the Hilbert transform operator; Step S32: Calculate the instantaneous amplitude of each data point based on the complex seismic traces. To avoid traditional arctangent calculation To avoid potential phase entanglement problems, we instead use an exponential form to implicitly estimate the instantaneous phase: This method expresses phase information through a complex exponent, without directly calculating the phase angle, effectively avoiding phase transitions and improving robustness in noisy data; Step S33: Instantaneous phase based on exponential form and Calculate the instantaneous phase correlation measure between observed and synthetic data at each moment. The metric is defined as follows: Equation (1) Using trigonometric identities, this expression can be simplified to: Equation (2) in, The closer this metric is to -1, the more similar the instantaneous phases of the two signals are. Its core idea is to quantify the difference in instantaneous phases and set the inversion objective as minimizing this difference. Step S34: Integrate the instantaneous phase correlation measures of all sources s, all detectors r, and the entire recording time T to construct a total amplitude-independent phase correlation objective function. As shown in the following formula: Equation (3) This objective function transforms the local phase matching problem into a global optimization problem by integrating the spatiotemporal data over the entire observation system, providing a robust optimization objective that is insensitive to amplitude errors for full waveform inversion.

[0025] In this application, the observation data channel refers to the field seismic records actually collected by detectors (seismographs) deployed on the surface or in wells; the synthetic data channel refers to the seismic records obtained through computer numerical simulation.

[0026] This embodiment extracts instantaneous phase information from observed and synthesized data through complex seismic trace analysis and constructs an exponential-based instantaneous phase correlation metric, achieving accurate quantification of seismic waveform phase similarity; utilizing... By characterizing the instantaneous phase, the entanglement problem in traditional phase difference calculation is effectively avoided. At the same time, by designing a normalized phase correlation metric formula, the similarity evaluation results are compressed to the range of [-1, 0], so that the closer the metric value is to -1, the higher the phase matching degree. Thus, a robust objective function that is insensitive to amplitude changes and focuses only on phase information is constructed, providing a stable and physically meaningful error metric basis for subsequent Sigmoid modulation inversion. In one embodiment, step S4: modulates the phase correlation objective function using the Sigmoid function, mapping the objective function value to a predetermined interval to obtain the modulated objective function, thereby enhancing the inversion robustness. This specifically includes the following steps: Step S41: Use the standard Sigmoid function As a modulation function and to determine the optimal parameters: First, the Sigmoid function is chosen as the modulation tool, and its general form is: Equation (5) in, The upper boundary, Control the slope and direction. The central position. When hour, Smoothly increasing from 0 to ;when ,hour, Decreasing. In Place, Its derivative form is concise: Equation (6) The Sigmoid function possesses bounded output and smooth monotonicity, effectively suppressing outliers and stabilizing gradient propagation. To optimize the modulation effect, the optimal parameters are determined by minimizing the squared difference between the least-squares objective function and the IPC objective function after the Sigmoid transformation on a discrete grid. Equation (7) in, and These are the least squares and IPC objective function values, respectively.

[0027] The parameters are determined by least-squares optimization on the discrete grid (Equation (7)). Figure 5 As shown, when the optimization result is in To obtain the global minimum, the standard Sigmoid function (L=1) is used: This function maps the input to the (0, 1) interval, and its output is bounded and smooth, making it suitable for target function modulation, such as... Figure 2 As shown; Figure 2 (a) is a graph showing the effect of different slopes k on the shape of the Sigmoid function when the center position x0=0; Figure 2 (b) is a graph showing the influence of different center positions x0 on the shape of the Sigmoid function when the slope k=1. Step S42: Obtain the phase correlation objective function As input to the Sigmoid function, nonlinear modulation is performed, and the target function after modulation is defined as: Equation (9) This transformation will change the original The range of values ​​is mapped from approximately [-1, 1] to the interval (0, 1), effectively compressing the influence range of extreme values ​​while smoothing the overall shape of the objective function. For example, when ≈ At point 1 (high phase coherence), ≈0.268; when When ≈1 (low phase coherence), ≈0.732. This mapping reduces the impact of outliers on the inversion. Modulated objective function It possesses the characteristics of being bounded and well-formed. By comparing the relationship between different objective functions and time shifts, such as... Figure 4 As shown, the Sigmoid-modulated IPC objective function (blue curve) effectively compresses the influence of extreme tail values ​​compared to the original IPC objective function (red curve), and forms a smoother curvature near the minimum. This shape helps to achieve a more stable optimization process, especially in scenarios where the time shift may be large in the early stages of inversion, suppressing spurious local minima and improving convergence reliability. In summary, Sigmoid modulation, through an adaptive scaling mechanism, balances the inversion process, providing stable updates when phase differences are extreme and strong momentum when they are moderate, thereby enhancing overall performance.

[0028] In one embodiment, step S5: Based on the modulated objective function, calculate its gradient with respect to the subsurface model parameters, and use this gradient to iteratively update the S-wave velocity model, specifically including the following steps: Step S51: Based on the modulated objective function The gradient of the underground S-wave velocity model parameter m was calculated using the adjoint state method. According to the chain rule, this gradient can be expressed as: Equation (10) in, The derivative of the Sigmoid function As an adaptive scaling factor, when the phase is highly coherent ( ≈0.268) or highly incoherent ( When the phase difference is approximately 0.732, the scaling factor is approximately 0.196, which plays a stabilizing role; when the phase difference is moderate ( When the scaling factor reaches 0.5, it reaches a peak of 0.25, providing stronger update motivation. The gradient of the original objective function is expressed as follows: Equation (4) in, This represents the synthetic wavefield perturbation caused by model perturbation. This gradient guides the iterative update of the subsurface model to minimize the phase-based objective function. Step S52: Preprocess the calculated gradients to improve stability, including applying a spatial taper function to the single-shot gradient to suppress the acquisition boundary effect, and performing Gaussian smoothing filtering on the result of the superposition of gradients from all shot points to suppress high-frequency noise and preserve the geologically significant structure. This preprocessing step effectively eliminates false gradients caused by the boundary of the acquisition system, smooths the gradient field, avoids model oscillations caused by high-frequency noise, and significantly improves the stability and rationality of model updates. Step S53: Using optimization algorithms such as gradient descent or quasi-Newton method, the S-wave velocity model is iteratively updated using the preprocessed gradient direction. The update formula is: in, Let be the step size for the k-th iteration.

[0029] This embodiment accurately calculates the gradient of the modulated objective function with respect to the S-wave velocity model parameters using the adjoint state method, and introduces the derivative of the Sigmoid function as an adaptive scaling factor to achieve intelligent adjustment of the gradient amplitude. Furthermore, the gradient is preprocessed using a spatial taper function and Gaussian smoothing filtering, effectively suppressing boundary effects and high-frequency noise, significantly improving the gradient's signal-to-noise ratio and geological rationality. Finally, iterative updates are performed using optimization algorithms such as gradient descent, ensuring the numerical stability and convergence efficiency of the model update process, thus providing a solid optimization foundation for obtaining high-precision and high-reliability underground S-wave velocity structures.

[0030] In one embodiment, step S6: When the preset inversion termination condition is met, the final inverted S-wave velocity structure model is output, specifically including the following steps: Step S61: After each iteration update, the system performs a comprehensive evaluation of the current inversion state to determine whether the preset termination conditions are met. The termination conditions are set with multiple judgment criteria, including monitoring the objective function value. The mechanism checks whether the rate of decline is continuously below a set threshold (e.g., the relative decline rate is less than 1% for two consecutive iterations) to capture the convergence trend of the inversion; it also checks whether the preset maximum number of iterations has been reached to ensure computational efficiency; it calculates whether the current gradient norm is below a given tolerance value to verify the necessity of model updates from an optimization perspective; and it evaluates the degree of improvement of the inverted model relative to the initial model, judging whether it tends to stabilize by the trend of the relative percentage error (RPE). This multi-condition joint judgment mechanism ensures the scientific nature of the stopping decision. Step S62: When any of the above termination conditions is triggered, the inversion iteration process stops immediately. At this time, the optimal S-wave velocity model obtained in the current iteration cycle is output as the final inversion result. This model is a high-resolution underground S-wave velocity structure model obtained through full waveform inversion, which accurately characterizes the spatial distribution of shear wave velocity in the underground medium of the work area, providing key physical property parameters for subsequent geological interpretation, engineering exploration, or resource exploration.

[0031] This embodiment provides a scientific and reliable stopping criterion for the full waveform inversion process by setting multiple inversion termination conditions (including the objective function descent rate, the maximum number of iterations, and the gradient norm threshold). This ensures that the algorithm terminates in a timely manner when the model converges to satisfactory accuracy, avoiding model distortion or waste of computational resources caused by excessive iteration. Its core advantage lies in its ability to sensitively capture the convergence trend of the inversion by monitoring the relative change in the objective function descent rate; the gradient norm threshold verifies the necessity of model updates from an optimization perspective; and the setting of the maximum number of iterations ensures computational efficiency. This multi-condition joint judgment mechanism not only guarantees the stability and reliability of the inversion results but also improves the practicality and automation of the algorithm, providing a crucial guarantee for the final output of a high-resolution, high-precision underground S-wave velocity structure model.

[0032] In a preferred embodiment, the method further includes employing a multi-scale inversion strategy: The inversion begins with data from a lower frequency band, such as 0-10Hz. After successful convergence, the inversion results are used as the initial model for the next higher frequency band inversion, gradually expanding to 0-20Hz, 0-30Hz, and so on, up to 0-60Hz. By progressively incorporating higher frequency information, the method can recover more detailed subsurface structures, effectively avoiding period jump problems and improving inversion stability and resolution.

[0033] It is worth noting that sigmoid modulation may slightly reduce sensitivity to high-contrast features in noise-free data, but significantly improves robustness in noisy data.

[0034] In this example, due to the use of a high-quality initial model, the advantages of sigmoid modulation were not very significant in improving the final model, but it still showed a positive effect on data fitting (lower RPE) and stability of the inversion process. It can be expected that the advantages of sigmoid modulation will be more pronounced when the initial model accuracy is lower or the data quality is worse.

[0035] Figure 7 (a) is a diagram illustrating the common shot point gather of the horizontal components under clean data. Figure 7 (b) is a diagram illustrating the gathers of the same common shot point after adding random noise; The following evaluation assesses the full waveform inversion performance based on the Sigmoid modulation IPC objective function through synthetic experiments and field data applications. The real S-wave velocity vs. the model from the synthetic test is shown below. Figure 6 As shown, the data consists of a uniform background with a velocity (Vs) of 300 m / s and a density (ρ) of 2000 kg / m³, and four anomalies of equal size. The anomaly velocities deviate from the background velocity by ±10%, or ±30 m / s. The acquisition geometry includes 40 horizontal component geophones in the 5-45 m range at 1 m intervals, and an SH wave source in the 10-40 m range at a 2 m step size. The top boundary is a free surface, while the lateral and bottom boundaries employ a fully matched convolutional layer (C-PML), referencing the work of Komatitsch and Martin (2007). The computational domain is a 200×40 grid with a spatial step size of 0.25 m. The finite difference method is used to solve the two-dimensional elastic wave equations, referencing the work of Bohlen (2002), to simulate the wave propagation process of a 30 Hz band-limited sharp pulse source with a duration of 0.2 s and a time step of 0.05 ms, generating 16 common shot gathers. To assess noise resistance, random noise was added to the clean synthetic data with a signal-to-noise ratio (S / N) of 3, and a bandpass filter of 5-130 Hz was applied to obtain the second dataset. Figure 7 The representative horizontal component common shot gathers at a focal depth of 10m are presented, including clean and noisy data.

[0036] Figure 8 (a) shows the inversion results of the standard IPC method without Sigmoid modulation; Figure 8 (b) shows the inversion result of the Sigmoid modulation IPC method; When assessing the sensitivity to the periodic jump problem, only Vs is updated during inversion, with the density fixed at ρ = 2000 kg / m³. The initial model uses a uniform background of Vs = 300 m / s, employing a single-scale inversion strategy with an inversion bandwidth limited to 0-40 Hz. Optimization terminates when the relative decrease in the objective function is less than 1% twice consecutively. Figure 8 As shown, the unmodulated IPC inversion only recovers shallow anomalous bodies but introduces sidelobe artifacts, while the Sigmoid-modulated IPC successfully recovers all four anomalous bodies with clearer boundaries and Vs values ​​closer to the true model, confirming the effect of modulation in suppressing spurious local minima. Figure 8The inversion results of two objective functions, unmodulated standard IPC and Sigmoid-modulated IPC, were compared. In the unmodulated case, the inversion successfully recovered two shallow anomalies, but introduced sidelobe artifacts that masked deeper structures and distorted the Vs value, indicating its susceptibility to local minima and periodic jumps. In contrast, Sigmoid-modulated IPC successfully recovered all four anomalies with clearer boundaries, fewer artifacts, and a Vs value closer to the true model. This demonstrates that the Sigmoid function effectively reshapes the objective function, suppresses spurious local minima, and yields more stable, high-resolution inversion results.

[0037] Figure 9 (a) shows the results of traditional least squares full waveform inversion under clean data. Figure 9 (b) shows the results of IPC full waveform inversion under clean data. Figure 9 (c) shows the results of Sigmoid modulation IPC full waveform inversion under clean data. Figure 9 (d) shows the results of traditional least-squares full waveform inversion on noisy data; Figure 9 (e) shows the results of IPC full waveform inversion under noisy data; Figure 9 (f) shows the results of full waveform inversion of Sigmoid-modulated IPC under noisy data; When analyzing robustness to random noise, the relative percentage error (RPE) is used to quantify the difference between the inversion model and the true model, as well as between the observed gather and the synthesized gather. RPE is defined as RPE = 100 × ||AB|| / ||B||, where a lower RPE value indicates better consistency. Under noisy data, the RPE of traditional least-squares FWI increases significantly. Both IPC and sigmoid-modulated IPC can suppress noise, but the modulation method achieves the lowest RPE. Figure 9 Middle (e) Figure 9 As shown in (f), this is due to the sigmoid function's suppression of bounded residuals. Figure 9 The image shows the Vs reconstruction results of three full waveform inversion methods: traditional least squares full waveform inversion, IPC-based full waveform inversion, and Sigmoid-modulated IPC-based full waveform inversion. The first row shows the inversion results for clean data: all three methods can recover the four-anomaly structure and approximate velocity contrast, but there are significant differences in RPE—traditional least squares full waveform inversion. Figure 9In (a), the RPE is the lowest, followed by IPC full waveform inversion, with the highest RPE for Sigmoid-modulated IPC. This trend is in line with expectations: Sigmoid modulation enhances robustness by reducing the impact of amplitude outliers, but in noise-free data, this modulation also slightly reduces sensitivity to high-contrast features, resulting in a slight decrease in resolution and a slightly higher RPE compared to unmodulated IPC. Figure 9 The second row shows the inversion results for noisy data: Traditional least-squares full-waveform inversion completely fails, overfitting the noise and producing spurious oscillations and block artifacts, leading to a significant increase in RPE and obscuring the true underground geometry. In contrast, both IPC and Sigmoid-modulated IPC suppress noise-driven updates, generating models highly similar to the real structure. Notably, Sigmoid-modulated IPC achieves the lowest RPE in noisy data inversion. This improvement is attributed to the bounded influence of the Sigmoid function—the large residuals caused by noise are further weighted and suppressed, preventing them from dominating gradient calculations, thus ensuring the inversion proceeds stably around a coherent, phase-consistent signal.

[0038] Figure 10 (a) is a comparison of the synthetic gather and the observation gather in clean data using traditional least squares full waveform inversion. Figure 10 (b) is a comparison of the synthetic gather and the observation gather in IPC full waveform inversion under clean data; Figure 10 (c) is a comparison of the synthesized gather and the observed gather in clean data for full waveform inversion of Sigmoid-modulated IPC. Figure 10 (d) is a comparison of the synthetic gather and the observation gather in the noisy data using traditional least squares full waveform inversion. Figure 10 (e) is a comparison of the synthetic gather and the observation gather in noisy data using IPC full waveform inversion; Figure 10 (f) is a comparison of the synthesized gather and the observed gather in noisy data for full waveform inversion of Sigmoid-modulated IPC; Figure 10 The synthetic gathers obtained from the forward modeling of the inversion model were compared with the observed data: the first row shows the simulated gathers and observed gathers of clean data superimposed, and the second row shows the comparison results of noisy data. Under noisy conditions, the Sigmoid modulation IPC has a lower RPE and better waveform alignment, compared with... Figure 9 (e) and Figure 9The Vs reconstruction results shown in (f) are consistent; however, in the clean data, the RPE of the gather generated by the Sigmoid modulated IPC inversion is slightly higher than that of the standard IPC. This may be because the noise-suppressing modulation weakens the large residuals that are truly valuable for information, such as the residuals related to strong impedance contrast, which leads to a slight reduction in the effective information content of the objective function in scenarios where robustness guarantees are not required.

[0039] Figure 11 (a) is an illustration of the initial S-wave velocity model obtained from First Arrival Travel Time Tomography (FATT) at the Fossa Carolina site; Figure 11 (b) is Figure 11 The smoothed version of (a) is used as an illustration of the initial model for full waveform inversion; Figure 12 (a) is a diagram illustrating the original preprocessed common shot gather; Figure 12 Figure (b) illustrates the same gather after adding random noise and performing bandpass filtering; In terms of field data application, the proposed method was applied to a publicly available field dataset of a medieval silted canal in Fossa, Carolina, southern Germany. The survey employed a horizontal hammer-driven seismic source and a fixed geophone array of 48 horizontal component geophones spaced 0.75 m apart, with the seismic source moving in 0.75 m increments. Köhn et al. (2019) initially constructed a Vs model using first-arrival travel-time tomography (FATT), successfully identifying low-velocity canal fill and surrounding high-velocity embankments—features subsequently validated by archaeological excavations. This model was used as a benchmark for subsequent full-waveform inversion studies, referencing the work of Tohti et al. (2002). Figure 12 Figure (a) shows a representative shot gather after preprocessing, where surface wave energy is significant. To test noise immunity, random noise with a signal-to-noise ratio (S / N) of 10 was added to the gather, and a 10-130 Hz bandpass filter was applied. Figure 12 (b) Simulates typical field data challenges. The inversion uses a smoothed version of the FATT model as the initial Vs model. Figure 11 In (b), only Vs is updated, while the density remains constant. A multi-scale strategy is adopted, referencing the study by Bunks et al. (195): starting from 0-10Hz, the bandwidth is gradually expanded to 20, 30, 40, 50, and 60Hz, with the inversion results of the previous stage serving as the initial model for the next stage. To enhance stability and mitigate the impact of acquisition artifacts, the gradient is preprocessed: a semi-circular tape is applied to the contribution of each shot point before stacking, and then a two-dimensional Gaussian filter is used with an effective radius approximately half that of the dominant wavelength. Referring to the work of Ravaut et al. (2004), the stacked gradient is smoothed to suppress subwavelength noise and preserve geologically significant contrast.

[0040] Figure 13 (a) shows the S-wave velocity inversion model results of the traditional least squares full waveform inversion method under the condition of preprocessed field data (clean data); Figure 13 Figure (b) shows the S-wave velocity inversion model results of the full waveform inversion method based on instantaneous phase correlation (IPC) under the condition of preprocessed field data (clean data); Figure 13 (c) shows the S-wave velocity inversion model results of the full waveform inversion method based on Sigmoid modulation IPC under the condition of preprocessed field data (clean data); Figure 13 The middle (d) figure shows the S-wave velocity inversion model results of the traditional least squares full waveform inversion method under the condition of noisy field data (signal-to-noise ratio S / N=10, 10-130Hz bandpass filter); Figure 13 The middle (e) figure shows the S-wave velocity inversion model results of the full waveform inversion method based on instantaneous phase correlation (IPC) under the condition of noisy field data (signal-to-noise ratio S / N=10, 10-130Hz bandpass filter); Figure 13 (f) shows the S-wave velocity inversion model results of the full waveform inversion method based on Sigmoid modulation IPC under noisy field data (signal-to-noise ratio S / N=10, 10-130Hz bandpass filter); Figure 13 Three inversion Vs models with different objective functions are presented, including clean field data and noisy field data: In the case of clean data, traditional least-squares full-waveform inversion fails due to numerical instability, such as... Figure 13 As shown in (a), a significant artifact occurs—this stems from the inconsistency in amplitude recorded in the field, which may be caused by changes in detector coupling, imperfect amplitude compensation for geometric diffusion and attenuation, or source coupling errors; while the two IPC methods Figure 13 (b) Figure 13 In (c), geologically plausible models were generated, successfully reconstructing the shallow low-velocity zone of approximately 50-150 m / s, corresponding to the bottom of the canal, and the asymmetric high-velocity body at a depth of approximately 1 m and a lateral depth of 15-27 m. These characteristics are consistent with previous archaeological findings, confirming that phase-based objective functions can effectively stabilize full-waveform inversion even when amplitude is unreliable. (In the case of noisy data...) Figure 13 (Lower row) Least square full waveform inversion further deteriorates, generating severely distorted or invalid models; IPC-type inversion remains robust, recovering the main anomalies and their approximate depths, but due to the reduced signal-to-noise ratio, there are slight fuzziness and positioning errors.

[0041] Figure 14(a) shows a comparison between synthetic gathers and observation gathers in traditional least-squares full waveform inversion using preprocessed field data (clean data); Figure 14 (b) is a comparison of synthetic gathers and observation gathers obtained from IPC full waveform inversion using preprocessed field data (clean data); Figure 14 (c) is a comparison of the synthetic gather and the observation gather in the preprocessed field data (clean data) using the full waveform inversion of Sigmoid-modulated IPC. Figure 14 (d) is a comparison of the synthesized gather and the observed gather in the traditional least squares full waveform inversion using noisy field data (S / N=10, 10-130Hz bandpass filter); Figure 14 (e) is a comparison of the synthesized gather and the observed gather based on IPC full waveform inversion using noisy field data (S / N=10, 10-130Hz bandpass filter); Figure 14 (f) is a comparison of the synthesized gather and the observed gather based on the full waveform inversion of Sigmoid-modulated IPC in noisy field data (S / N=10, 10-130Hz bandpass filter); like Figure 14 As shown, the data fitting effect is evaluated by overlaying the synthetic seismogram simulated by the inversion model with the corresponding observation gathers. The top row represents clean field data, and the bottom row represents noisy data. In both cases, the RPE of the Sigmoid-modulated IPC is lower than that of the standard IPC, confirming that Sigmoid modulation enhances robustness to outliers. These results indicate that the proposed method can be directly applied to actual field data, and under given acquisition geometry and initial model conditions, its performance is slightly better than that of the unmodulated IPC. Although the high-quality initial model limits the observability advantage of modulation in this specific case, combined with the noise resistance verified in the synthetic test, the Sigmoid-modulated IPC remains a practical and reliable full-waveform inversion scheme for field-scale surface waves.

[0042] This application proposes an enhanced full-waveform inversion FWI objective function—integrating Sigmoid modulation into the instantaneous phase-correlation (IPC) framework. Theoretical analysis and numerical experiments confirm that this improvement effectively reshapes the objective function, suppresses spurious local minima, and retains the inherent robustness of IPC to amplitude errors. Synthetic inversion results show that, compared with standard IPC full-waveform inversion, the Sigmoid-modulated IPC objective function significantly improves convergence characteristics and depth resolution, especially in the recovery of deep subsurface anomalies. Furthermore, this method exhibits superior robustness to random noise, achieving lower relative percentage error (RPE) for both the inversion model under noisy conditions and the synthetic data fitting. Applying this method to a benchmark field dataset from the Fossa Carolina archaeological site in southern Germany, both modulated and unmodulated IPC inversions generated S-wave velocity-vs structures that highly match independently verified archaeological features. Although the high-quality initial model limited the observable advantages of Sigmoid modulation in this specific case, the modulation method consistently achieved comparable or slightly improved waveform fitting results, demonstrating low sensitivity to the initial model and high stability under noisy conditions. The results show that Sigmoid modulation IPC is a practical and reliable alternative to surface wave full waveform inversion, especially suitable for complex shallow surface scenarios with high risk of period jumps, inconsistent amplitudes, or low signal-to-noise ratios.

[0043] Secondly, embodiments of this application provide a Sigmoid-modulated phase-correlated surface wave full waveform inversion system, including: The data acquisition module is used to acquire field surface wave seismic records and construct an initial subsurface S-wave velocity model; The forward modeling module is communicatively connected to the data acquisition module and is used to perform forward modeling based on the initial S-wave velocity model to generate synthetic seismic records. The objective function construction module is communicatively connected to the forward modeling module and is used to calculate the instantaneous phase correlation between the observation record and the synthetic record, and to construct the phase correlation objective function. The gradient calculation and inversion update module is communicatively connected to the objective function construction module. It is used to modulate the phase correlation objective function using the Sigmoid function, map the objective function value to a predetermined interval, and obtain the modulated objective function. The inversion control module is communicatively connected to the gradient calculation and inversion update module. It is used to calculate the gradient of the modulated objective function with respect to the underground model parameters, and to iteratively update the S-wave velocity model using the gradient. When the preset inversion termination condition is met, the final inverted S-wave velocity structure model is output.

[0044] The functions of each module in the above-mentioned Sigmoid-modulated phase-dependent surface wave full waveform inversion system correspond to the steps in the above-mentioned Sigmoid-modulated phase-dependent surface wave full waveform inversion method embodiment, and their functions and implementation processes will not be described in detail here.

[0045] Thirdly, embodiments of this application provide a Sigmoid-modulated phase-dependent surface wave full waveform inversion device, which can be a personal computer (PC), laptop computer, server, or other device with data processing capabilities.

[0046] In this embodiment of the application, the sigmoid-modulated phase-dependent surface wave full waveform inversion device may include a processor, a memory, a communication interface, and a communication bus.

[0047] The communication bus can be of any type and is used to interconnect the processor, memory, and communication interface.

[0048] The communication interface includes input / output (I / O) interfaces, physical interfaces, and logical interfaces for interconnecting devices within the Sigmoid-modulated phase-correlated surface wave full waveform inversion device, as well as interfaces for interconnecting the Sigmoid-modulated phase-correlated surface wave full waveform inversion device with other devices (such as other computing devices or user equipment). Physical interfaces can be Ethernet interfaces, fiber optic interfaces, ATM interfaces, etc.; user equipment can be displays, keyboards, etc.

[0049] Memory can be various types of storage media, such as random access memory (RAM), read-only memory (ROM), non-volatile RAM (NVRAM), flash memory, optical storage, hard disk, programmable ROM (PROM), erasable PROM (EPROM), electrically erasable PROM (EEPROM), etc.

[0050] The processor can be a general-purpose processor, which can call the Sigmoid modulation phase-dependent surface wave full waveform inversion program stored in memory and execute the Sigmoid modulation phase-dependent surface wave full waveform inversion method provided in the embodiments of this application. For example, the general-purpose processor can be a central processing unit (CPU). The method executed when the Sigmoid modulation phase-dependent surface wave full waveform inversion program is called can refer to the various embodiments of the Sigmoid modulation phase-dependent surface wave full waveform inversion method of this application, and will not be repeated here.

[0051] Fourthly, embodiments of this application also provide a readable storage medium.

[0052] The present application stores a Sigmoid-modulated phase-dependent surface wave full waveform inversion program on a readable storage medium, wherein when the Sigmoid-modulated phase-dependent surface wave full waveform inversion program is executed by a processor, the steps of the Sigmoid-modulated phase-dependent surface wave full waveform inversion method described above are implemented.

[0053] The method implemented when the Sigmoid-modulated phase-dependent surface wave full waveform inversion procedure is executed can be referred to in various embodiments of the Sigmoid-modulated phase-dependent surface wave full waveform inversion method of this application, and will not be repeated here.

[0054] It should be noted that the sequence numbers of the embodiments in this application are for descriptive purposes only and do not represent the superiority or inferiority of the embodiments.

[0055] Through the above description of the embodiments, those skilled in the art can clearly understand that the methods of the above embodiments can be implemented by means of software plus necessary general-purpose hardware platforms. Of course, they can also be implemented by hardware, but in many cases the former is a better implementation method. Based on this understanding, the technical solution of this application, in essence, or the part that contributes to the prior art, can be embodied in the form of a software product. This computer software product is stored in a storage medium (such as ROM / RAM, magnetic disk, optical disk) as described above, and includes several instructions to cause a terminal device to execute the Sigmoid modulation phase-related surface wave full waveform inversion method described in the various embodiments of this application.

[0056] The above are merely preferred embodiments of this application and do not limit the patent scope of this application. Any equivalent structural or procedural transformations made using the content of this application's specification and drawings, or direct or indirect applications in other related technical fields, are similarly included within the patent protection scope of this application.

Claims

1. A method for full waveform inversion of sigmoid-modulated phase-correlated surface waves, characterized in that, Includes the following steps: Acquire field surface wave seismic records and construct an initial subsurface S-wave velocity model; Forward modeling was performed based on the initial S-wave velocity model to generate synthetic seismic records; Calculate the instantaneous phase correlation between the observation record and the synthetic record, and construct a phase correlation objective function; The phase correlation target function is modulated using the Sigmoid function, and the target function value is mapped to a predetermined interval to obtain the modulated target function. Based on the modulated objective function, its gradient with respect to the underground model parameters is calculated, and the S-wave velocity model is iteratively updated using this gradient. When the preset inversion termination condition is met, the final inverted S-wave velocity structure model is output.

2. The sigmoid-modulated phase-correlated surface wave full waveform inversion method according to claim 1, characterized in that, The instantaneous phase correlation is constructed based on the complex seismic traces of observational and synthetic data.

3. The sigmoid-modulated phase-correlated surface wave full waveform inversion method according to claim 2, characterized in that, The instantaneous phase correlation objective function Represented as: ； in, and To correspond to the observation data respectively and synthetic data The instantaneous amplitude, This represents the summation over all sources and detectors. To record time.

4. The sigmoid-modulated phase-correlated surface wave full waveform inversion method according to claim 1, characterized in that, The Sigmoid function mentioned is the standard Sigmoid function, and its expression is: ； The modulated objective function Represented as: ； in, Let be the objective function for the instantaneous phase correlation.

5. The sigmoid-modulated phase-correlated surface wave full waveform inversion method according to claim 1, characterized in that, When calculating the gradient with respect to the model parameters based on the modulated objective function, the chain rule is used, and the gradient expression is: ； in, The gradient of the modulated objective function; This represents the gradient of the original objective function.

6. The sigmoid-modulated phase-correlated surface wave full waveform inversion method according to claim 1, characterized in that, Before iteratively updating the S-wave velocity model, the gradient is preprocessed, including: Apply a spatial taper function to the single-shot gradient; Apply a Gaussian filter to the superimposed gradients.

7. The sigmoid-modulated phase-correlated surface wave full waveform inversion method according to claim 1, characterized in that, A multi-scale inversion strategy is employed, including: The inversion process starts from the low-frequency band and gradually expands to higher frequency bands. The velocity model obtained from the previous frequency band inversion is used as the initial model for the next frequency band inversion.

8. The sigmoid-modulated phase-correlated surface wave full waveform inversion method according to claim 1, characterized in that, The method is used to process field surface wave seismic data that includes amplitude errors or random noise.

9. A sigmoid-modulated phase-correlated surface wave full waveform inversion system, characterized in that, include: The data acquisition module is used to acquire field surface wave seismic records and construct an initial subsurface S-wave velocity model; The forward modeling module is communicatively connected to the data acquisition module and is used to perform forward modeling based on the initial S-wave velocity model to generate synthetic seismic records. The objective function construction module is communicatively connected to the forward modeling module and is used to calculate the instantaneous phase correlation between the observation record and the synthetic record, and to construct the phase correlation objective function. The gradient calculation and inversion update module is communicatively connected to the objective function construction module. It is used to modulate the phase correlation objective function using the Sigmoid function, map the objective function value to a predetermined interval, and obtain the modulated objective function. The inversion control module is communicatively connected to the gradient calculation and inversion update module. It is used to calculate the gradient of the modulated objective function with respect to the underground model parameters, and to iteratively update the S-wave velocity model using the gradient. When the preset inversion termination condition is met, the final inverted S-wave velocity structure model is output.

10. A computer-readable storage medium, characterized in that, The computer-readable storage medium stores a Sigmoid-modulated phase-dependent surface wave full waveform inversion program, wherein when the Sigmoid-modulated phase-dependent surface wave full waveform inversion program is executed by a processor, it implements the steps of the Sigmoid-modulated phase-dependent surface wave full waveform inversion method as described in any one of claims 1 to 8.

Images