A method for removing strong shielding from hidden river channels based on adaptive hybrid L0-L1 norm

Abstract

This invention relates to the field of seismic data processing technology, specifically disclosing a method for removing strong reflection shielding in hidden river channels based on adaptive hybrid L0-L1 norm. The method includes: first, establishing an optimization objective function based on Bayes' theorem using hybrid L0-L1 norm; second, dynamically adjusting the L0-L1 norm weights using an adaptive weight function driven by the seismic signal, combined with reflection coefficient amplitude and residual information, enhancing sparsity constraints in strong reflection zones and reducing constraint strength in weak reflection zones; then, constructing a convex upper bound for the objective function using a minimization framework, and solving it iteratively in stages using an accelerated rapid iterative threshold shrinkage algorithm, while incorporating prior knowledge of seismic wave propagation laws and river channel deposition patterns to ensure the geological rationality of the solution. This invention solves the technical problems of traditional sparse processing methods, such as fixed parameters, lack of geological constraints, and inability to simultaneously address strong reflection suppression and weak signal protection, significantly improving the separation accuracy of strong reflections and the recovery rate of weak signals, effectively overcoming the strong reflection shielding effect.

Description

Technical Field

[0001] This invention belongs to the field of seismic data processing technology, specifically relating to a method for removing strong shielding from hidden river channels based on adaptive hybrid L0-L1 norm, which is suitable for extracting weak reflection signals and overcoming the strong reflection shielding effect. Background Technology

[0002] The core objective of seismic strong reflection separation is to accurately separate the strong reflection components from the original seismic record. The strong reflection components mainly correspond to the strong impedance difference response at stratigraphic interfaces, while the weak reflection components correspond to the response of target bodies such as reservoirs. The separation process must balance separation accuracy and the fidelity of the weak reflection signal. To achieve this goal, researchers have proposed various seismic strong reflection separation methods.

[0003] With the application of compressed sensing theory in the field of seismic signal processing, norm-constrained optimization methods have been gradually introduced into the separation of strong seismic reflections. Among them, the L0 norm can accurately characterize the sparsity of the signal and effectively locate the sparse distribution characteristics of strong reflection signals, thus achieving accurate extraction of strong reflection components. The L1 norm has a stable convergence process and a certain degree of noise resistance, and is widely used in sparse signal separation. However, the L1 norm is prone to causing amplitude distortion of weak reflection signals in the process of strong reflection separation, resulting in limited separation accuracy.

[0004] Therefore, developing a strong reflection separation method that can adaptively adapt to different seismic data characteristics, balance separation accuracy and signal fidelity, and does not require complex empirical parameter tuning has become an urgent technical problem to be solved in the field of seismic exploration data processing. Summary of the Invention

[0005] The purpose of this invention is to address the shortcomings of the prior art by providing a method for de-shielding strong reflections in concealed river channels based on adaptive hybrid L0-L1 norm. This invention can significantly reduce the shielding effect of strong reflections, enabling accurate identification and detailed characterization of weak signals against a strong reflection background.

[0006] To achieve the above-mentioned objectives, the technical solution adopted by this invention is as follows:

[0007] Step 1: To unify the modeling of sparse constraints and data fitting, this invention constructs an optimization objective function based on Bayes' theorem, incorporating adaptive hybrid L0-L1 norms. This design uses seismic wave propagation patterns and riverbed depositional models as prior embedded models, explicitly guaranteeing the geological rationality of the solution and avoiding deviations from actual geological laws by purely mathematical optimization. The optimization objective function is:

[0008] (1)

[0009] in, For observation data, For the target reflection coefficient signal, For system operators, For transform domain operators, An adaptive weighting function that depends on signal strength. This is the hyperparameter for regularization strength. It is the L0 norm. It is an L1 norm. It is the L2 norm. This indicates the search for the variable that minimizes the objective function. ;

[0010] Step 2: Design a data-driven adaptive weight update function. To achieve differentiated constraints for strong and weak reflection regions, this invention constructs a dynamic weight function based on the reflection coefficient amplitude and iterative residual. First, in the... In this iteration, initial adaptive weights are constructed based on the ratio of the Sigmoid function to the normalized residual:

[0011] (2)

[0012] in, For the first In the nth iteration One reflection coefficient component For the first The residual vector of the next iteration The sigmoid steepness coefficient. The threshold value for strong reflection amplitude. For residual adjustment parameters, Let be the infinite norm of the residual vector;

[0013] Based on the constructed initial adaptive weights, an exponentially weighted smoothing strategy is used to stabilize the initial adaptive weights:

[0014] (3)

[0015] in, For smoothing coefficients, For the first The weight values ​​for the next iteration. For the first The maximum value of the initial weight in the next iteration;

[0016] Based on the above adaptive weighting function, the L0 norm sparse constraint is automatically enhanced in the strong reflection region to accurately suppress strong energy interference, and the L1 norm stable constraint is automatically enhanced in the weak reflection region to protect weak signal details, thus solving the defect that fixed weights cannot adapt to different reflection intensities.

[0017] Step 3: To overcome the difficulty of solving the problem caused by the non-convexity and non-differentiability of the L0 norm, this invention first accurately locks the strong reflection position through the dual conditions of amplitude and residual. When the reflection coefficient components simultaneously satisfy:

[0018] and (4)

[0019] The location of a signal is determined to be a strong reflection location when it satisfies the above formula. The residual sensitivity threshold;

[0020] Then, the convex upper bound is constructed using the maximization (MM) framework, while the upper bound is corrected by introducing a geological prior constraint coefficient to avoid the weak reservoir signal being over-suppressed and to enhance the constraint strength of the strong reflection zone.

[0021] Step 4: To improve the efficiency and stability of the solution, this invention transforms the optimization problem of adaptive hybrid L0-L1 norm into a weighted L1 regularized sparse reconstruction problem, and uses the Accelerated Fast Iterative Shrinking Threshold (FISTA) algorithm to solve the convex upper bound function in stages to achieve efficient solution.

[0022] Step 5: When the iteration results meet the convergence accuracy or reach the maximum number of iterations, the separated strong reflection signal and the recovered weak reflection signal are output. The final result still maintains high robustness in noisy environments and can be directly used for fine characterization of oil and gas reservoirs.

[0023] Preferably, the strong reflection amplitude threshold mentioned in step 2 Determined adaptively through the statistical distribution of the reflection coefficient amplitude:

[0024] (5)

[0025] in, This represents the average amplitude of the reflection coefficient in the current iteration. The standard deviation of the reflection coefficient amplitude in the current iteration This is the threshold control coefficient.

[0026] Preferably, the specific process of constructing the convex upper bound function in step 3 is as follows:

[0027] First, we construct an element-level upper bound for the L0 norm, which is represented as:

[0028] (6)

[0029] in, For indicator functions, in the first... Next iteration point Construct an upper bound for the quadratic function using the indicator function:

[0030] (7)

[0031] in, The value should be a very small positive number to avoid the denominator approaching 0;

[0032] Secondly, physical constraint coefficients are introduced. The upper bound of the element level was revised to set the reservoir development depth range and the small fault prediction region. Set up for areas with high reflectivity and noise interference The revised upper bound is:

[0033] (8)

[0034] Finally, by combining the data fitting term and the L1 norm term, a convex upper bound function for the overall optimization objective function is constructed:

[0035] (9)

[0036] in, For the constructed element-level upper bound, Data-driven weights , These are the regularization coefficients for the L0 and L1 norms, respectively.

[0037] Preferably, the specific steps of iteratively solving using the FISTA algorithm in step 4 are as follows:

[0038] First, define the decomposition terms of the objective function:

[0039] (10)

[0040] Solve based on the decomposition terms of the objective function. gradient:

[0041] (11)

[0042] in, For system operators Transpose of;

[0043] Secondly, determine the gradient Lipschitz constant. ,in For matrix The largest eigenvalue. In the th... In this iteration, the gradient descent step is calculated. ,in Let be the intermediate variable in the t-th iteration; solve for the iterative update value using the proximal operator. ,in For proximal operators, For the soft threshold operator; repeat the above steps until the convex optimization subproblem of a single MM framework is solved, and update the target signal. .

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

[0045] This invention proposes a method for removing strong reflections from hidden channels based on adaptive hybrid L0-L1 norms. The method incorporates prior seismic knowledge through Bayesian probabilistic modeling to impose geological constraints on the physical model, explicitly ensuring the geological rationality of the learned information during optimization. Furthermore, it incorporates L0-L1 norms to impose sparsity constraints on the seismic data, and uses a dual mechanism of amplitude and residual discrimination to adaptively adjust the L0-L1 norm weights. In strong reflection zones, the sparsity constraint of the L0 norm is enhanced to efficiently suppress strong reflections, while in weak signal zones, the stability constraint of the L1 norm is enhanced to protect the details of weak channel reflections. Then, a convex upper bound function is constructed for the optimization objective function, and physical constraint coefficients are introduced to correct the element-level upper bound, preventing excessive suppression of weak reservoir signals and enhancing the constraint strength in strong reflection zones. The goal of separating strong seismic reflections is achieved by setting iterative convergence conditions. This invention can accurately separate strong reflection interference while effectively preserving effective weak reflection signals, significantly improving the resolution and signal-to-noise ratio of seismic data, and providing reliable technical support for fine description and reservoir characterization of oil and gas reservoirs. Attached Figure Description

[0046] To more clearly illustrate the embodiments of the present invention or the advantages of the prior art, the accompanying drawings used in the description of the embodiments or the prior art are briefly introduced below:

[0047] Figure 1 Flowchart for invention;

[0048] Figure 2 For synthesizing seismic record models;

[0049] Figure 3 The processing results of different methods are compared in a noise-free environment and a 5dB noise environment: (a) Processing result of compressed sensing method in a noise-free environment; (b) Processing result of compressed sensing method in a 5dB noise environment; (c) Processing result of compressed sensing combined with L1 norm regularization method in a noise-free environment; (d) Processing result of compressed sensing combined with L1 norm regularization method in a 5dB noise environment; (e) Processing result of the present method in a noise-free environment; (f) Processing result of the present method in a 5dB noise environment. Detailed Implementation

[0050] The invention will now be further described with reference to the accompanying drawings. See also: Figure 1-3 :

[0051] I. Constructing a synthetic seismic record

[0052] The following synthetic seismic record was constructed:

[0053] Ricker wavelet with a main frequency of 35Hz ( Figure 2 (b) is used as a seismic wavelet, and is compared with a preset reflection coefficient sequence ( Figure 2 (a) Perform convolution to generate a single-track synthetic seismic record ( Figure 2 (c) The reflection coefficient sequence includes three strong reflection interfaces and multiple sets of weak reflection segments, which are used to simulate the reservoir response under a strong shielding background.

[0054] Let the wavelet be The reflection coefficient is Then synthesize the seismic record Calculate using the following formula:

[0055]

[0056] in This indicates a convolution operation, with a sampling interval set to 1ms to ensure coverage of both strong and weak reflection target segments.

[0057] II. Initialization Parameters and Data Preparation

[0058] First, input the observed earthquake data. Constructing operators for forward modeling systems With transform domain operators Calculate the mean value based on the amplitude distribution of the initial reflection coefficient. with standard deviation The threshold for strong reflection amplitude is adaptively determined by the following formula. :

[0059]

[0060] The mean Take the statistical average and standard deviation of the reflection coefficient amplitude in the current iteration. Reflecting the degree of amplitude fluctuation, the threshold control coefficient is set based on the relative amplitude of the synthetic seismic record; in this example, it is set to 0.4, making... It can effectively distinguish strong reflections from background signals.

[0061] III. Adaptive Weight Calculation

[0062] In the In the nth iteration, for the th The initial adaptive weights for the reflection coefficient components are calculated using the following formula:

[0063]

[0064] The slackness coefficient of the Sigmoid function controls the steepness of the curve; the larger the slackness coefficient, the steeper the curve. In this example, the slackness coefficient of the Sigmoid function is set to 10, so that the weights are within the threshold range. Quickly switch to nearby locations; the residual adjustment parameter can be set to a positive number. Its value in the range (0, 1) indicates that the residual is reduced. When its value is 1, there is no residual adjustment by default. When its value is greater than 1, the residual is amplified. In this example, its value is set to 0.2 to enhance the guiding effect of the residual on the weight.

[0065] To avoid iterative oscillations, exponential weighted smoothing is used to stabilize the initial weights:

[0066]

[0067] The smoothing coefficient typically ranges from (0, 1). In this example, it is set to 0.8 to make the weight changes smoother and improve the stability of the iteration.

[0068] IV. Strong Reflection Dual Discrimination

[0069] A strong reflection dual discrimination mechanism is implemented, combining amplitude and residual dynamic weight adjustment. For each signal point, the following two conditions are simultaneously determined:

[0070]

[0071] When both conditions are met, the point is determined to be a strong reflection location. The residual sensitivity threshold is usually set in the range of (0, 1). In this example, it is set to 0.3, which can effectively eliminate false strong reflections caused by noise.

[0072] V. Structural and Geological Constraint Correction of the Convex Upper Boundary

[0073] The indicator function for the L0 norm at the current iteration point The second upper bound of the structure is:

[0074]

[0075] Here, 1e-6 is a very small positive number to avoid the denominator approaching 0.

[0076] Introducing geological prior constraint coefficients The upper limit is modified, and the reservoir development depth range and weak reflection target area are taken as... For areas with high reflection density and noise interference, take The revised upper bound is:

[0077]

[0078] By combining the corrected upper bound with the data fitting term and the L1 norm term, a convex upper bound function is obtained, which enables the original non-convex problem to be solved efficiently.

[0079] VI. Accelerating the FISTA Iterative Solution

[0080] Figure 1 In step 4, the sparse reconstruction process uses the accelerated FISTA algorithm to iteratively solve the convex upper bound function in stages. First, the objective function decomposition terms are defined:

[0081]

[0082]

[0083] calculate gradient:

[0084]

[0085] According to the matrix The largest eigenvalue determines the Lipschitz constant L. Gradient descent is performed in the t-th iteration:

[0086]

[0087] Iterative updates are achieved using the soft threshold operator:

[0088]

[0089] Repeat the above steps, update the reflection coefficient estimate after completing one minimization iteration, and enter the next loop.

[0090] VII. Convergence Judgment and Result Output

[0091] Figure 1 Step 5: The iteration convergence condition is set to the relative change of the results of two adjacent iterations being less than the relative change precision, or reaching the maximum number of iterations. The relative change precision and the maximum number of iterations can be adjusted according to actual needs. In this example, the relative change precision is set to 1e-6, and the maximum number of iterations is set to 200. The iteration stops when either condition is met. The final converged reflection coefficient is taken as the optimal solution. Based on this, the strong reflection component is separated and the weak reflection signal is recovered. The above parameter settings can simultaneously take into account both the noise-free case and the case with 5dB noise in the implementation case to achieve the optimal result.

[0092] Let's take a synthetic seismic record signal as an example. Figure 2 The convolution process of the synthetic record is shown, where the simulated seismic signal is synthesized by convolving a 35Hz Ricker wavelet with the reflection coefficient sequence. The reflection coefficient sequence indicates that times with higher reflection coefficients correspond to larger amplitudes in the synthetic record, and thus stronger seismic activity.

[0093] Figure 3 This figure compares the results of removing strong reflections from synthetic seismic records using different methods under noise-free and noise-containing conditions (with 5 dB of noise). The black curve represents the original signal, the blue curve represents the extracted strong reflection signal, and the red curve represents the remaining signal after separating the strong reflections. The horizontal axis represents time (ms), and the vertical axis represents amplitude (Amp).

[0094] from Figure 3 From the noise-free processing results of (a), (c), and (e), all three methods can effectively extract the strong reflection component from the original signal. The blue curve highly coincides with the strong reflection peak of the original signal, while the red curve shows a lower overall amplitude of the remaining signal, verifying the effectiveness of each method in separating strong reflection signals in noise-free scenarios. In comparison, the method of this invention ( Figure 3 (e) The waveform matching is better in the strong reflection region, and the fluctuation suppression effect of the remaining signal is more significant.

[0095] In a 5dB noise interference scenario ( Figure 3 (b)(d)(f)) The noise resistance of compressed sensing methods and L1 norm regularization methods has obvious limitations: Figure 3 (b) The blue extraction curve is significantly distorted due to noise interference, and the noise floor of the red residual signal is raised; Figure 3 (d) The peak broadening effect of the strong reflection signal is significant, and the signal details are severely lost. However, the method of this invention ( Figure 3 (f) demonstrates excellent noise robustness. The blue extraction curve can still accurately track the temporal and amplitude characteristics of the strong reflection peak of the original signal. The noise fluctuation of the red residual signal is effectively suppressed, retaining only a weak residual energy, which reflects the advantages of this method in balancing sparse constraints and residual fitting.

[0096] A comprehensive comparison shows that, in the task of separating strong reflection signals, the method of this invention has higher extraction accuracy and noise suppression capability than traditional compressed sensing and single L1 regularization methods in both noise-free and low signal-to-noise ratio environments, and is more suitable for strong reflection signal separation scenarios in complex noise backgrounds.

[0097] Although specific embodiments of the invention have been described in detail with reference to the accompanying drawings, this should not be construed as limiting the scope of protection of this patent. Various modifications and variations that can be made by a person skilled in the art without inventive effort within the scope described in the claims still fall within the scope of protection of this patent.

Claims

1. A method for removing strong shielding from hidden river channels based on adaptive hybrid L0-L1 norm, characterized in that, Includes the following steps: Step 1: Based on Bayes' theorem, establish an optimization objective function for de-strong reflection using adaptive hybrid L0-L1 norm regularization. The optimization objective function is: (1) in, For observation data, The target reflectance coefficient, This represents the intermediate value of the reflection coefficient during the iteration process. For system operators, For transform domain operators, An adaptive weighting function that depends on signal strength. This is the hyperparameter for regularization strength. It is the L0 norm. It is an L1 norm. It is the L2 norm. This indicates the search for the variable that minimizes the objective function. ; Step 2: Design an adaptive weight update function that dynamically adjusts the regularization weights of strong and weak reflection regions by combining amplitude and iteration residual information; Step 3: Establish a strong reflection dual discrimination mechanism, and at the same time, use the maximization framework to construct a convex upper bound function for the objective function to address the non-convexity and non-differentiability of the L0 norm. Meanwhile, introduce physical constraint coefficients to correct the element-level upper bound. Step 4: The optimization problem of adaptive hybrid L0-L1 norm is transformed into a sparse reconstruction problem with weighted L1 regularization. The convex upper bound function is solved in stages using an accelerated fast iterative shrinking threshold algorithm. Step 5: Set the iterative convergence condition, when the reflection coefficient reaches the intermediate value of the iteration. The iteration is terminated when the iteration result meets the convergence accuracy or reaches the maximum number of iterations, and the converged result is... Determined as the target reflectance coefficient and based on The output consists of the separated strong reflection signal and the recovered weak reflection signal.

2. The method for removing strong shielding from hidden river channels based on adaptive hybrid L0-L1 norm according to claim 1, characterized in that, The construction process of the adaptive weight update function described in step 2 is as follows: In the In this iteration, initial adaptive weights are constructed based on the ratio of the activation function to the normalized residual: (2)； in, For the first In the nth iteration One reflection coefficient component For the first The residual vector of the next iteration The steepness coefficient of the activation function. The threshold value for strong reflection amplitude. For residual adjustment parameters, Let the infinite norm of the residual vector be . It is an exponential function; An exponentially weighted smoothing strategy is used to stabilize the initial adaptive weights. (3)； in, For smoothing coefficients, For the first The weight values ​​for the next iteration. For the first The maximum value of the initial weight in the next iteration.

3. The method according to claim 2, characterized in that, The strong reflection amplitude threshold Determined adaptively through the statistical distribution of the reflection coefficient amplitude: (4) in, This represents the average amplitude of the reflection coefficient in the current iteration. The standard deviation of the reflection coefficient amplitude in the current iteration This is the threshold control coefficient.

4. The method according to claim 2, characterized in that, The criterion for establishing the strong reflection dual discrimination mechanism in step 3 is: when the reflection coefficient components simultaneously satisfy: (5)； Among them, the signal position is determined to be a strong reflection position when it satisfies formula (5). This is the residual sensitivity threshold.

5. The method for removing strong shielding from hidden river channels based on adaptive hybrid L0-L1 norm according to claim 1, characterized in that, The specific process of constructing the convex upper bound function in step 3 is as follows: Element-wise upper bounds are constructed for the L0 norm, and the L0 norm is expressed as: (6)； in, For indicator functions, in the first... Next iteration point Construct an upper bound for the quadratic function using the indicator function: (7)； in, The value should be a very small positive number to avoid the denominator approaching 0; Introducing physical constraint coefficients The upper bound of the element level was revised to set the reservoir development depth range and the small fault prediction region. Set up for areas with high reflectivity and noise interference The revised upper bound is: (8)； By combining the data fitting term and the L1 norm term, a convex upper bound function is constructed for the overall optimization objective function: (9)； in, For the constructed element-level upper bound, Data-driven weights , These are the regularization coefficients for the L0 and L1 norms, respectively.

6. The method for removing strong shielding from hidden river channels based on adaptive hybrid L0-L1 norm according to claim 1, characterized in that, The specific steps for iteratively solving the threshold using the accelerated fast iterative shrinkage algorithm in step 4 are as follows: First, define the decomposition terms of the objective function: (10)； Further, solve gradient: (11)； in, For system operators transpose; Then, determine the gradient Lipschitz constant. ,in For matrix The largest eigenvalue; In the In this iteration, the gradient descent step is calculated. ,in Let be an intermediate variable in the t-th iteration; Iterative update values ​​are solved using proximal operators. ,in For proximal operators, This is a soft threshold operator; Finally, repeat the above steps until the convex optimization subproblem of the single-step minimization framework is solved, and update the intermediate values ​​of the reflection coefficients during iteration. .

7. The method for removing strong shielding from hidden river channels based on adaptive hybrid L0-L1 norm according to claim 1, characterized in that, In step 5, the iteration convergence accuracy and the maximum number of iterations are set. The iteration can be stopped when the following conditions are met: (12)； in, For the first The intermediate value of the reflection coefficient in the next iteration. For the first The intermediate value of the reflection coefficient in the next iteration. This represents the maximum number of iterations.

Images