A method and system for compressing seismic multiple wave data

Abstract

The present invention discloses a method and system for compressing seismic multiple wave data, belonging to the field of data processing. The data compression method includes determining actual reflection waves based on raw seismic wave data; estimating each order of multiple waves based on the actual reflection waves; calculating forward multiple wave supergathers and reverse multiple wave supergathers based on the multiple waves; and calculating forward and reverse multiple wave plane waves after linearly delaying the forward and reverse multiple wave supergathers. The present invention can achieve secondary compression of seismic multiple wave data, significantly reduce the amount of migration data, and significantly improve the computational efficiency of multiple wave migration.

Description

Technical Field

[0001] The present invention relates to the field of geophysical technology, and more particularly to a method and system for compressing seismic multiple wave data. Background Art

[0002] Oil and natural gas will continue to dominate my country's energy mix for a long time to come. my country's offshore oil and gas reserves are extremely abundant, but the level of proven reserves is far below the global average. High-precision seismic migration imaging technology is key to the successful exploration and development of offshore oil and gas reservoirs.

[0003] Marine seismic data exhibit a significant abundance of multiple waves, which conventional reflection-based migration algorithms treat as data noise and suppress. However, multiple waves, as true reflections from structural interfaces, propagate back and forth underground, exhibiting smaller reflection angles, longer propagation paths, and a wider coverage area than reflection waves. Therefore, fully utilizing the unique advantages of multiple waves and performing imaging can effectively illuminate shadowed areas in reflection wave imaging and improve accuracy.

[0004] In recent years, multiple imaging has become a research hotspot in exploration geophysics. The core theory behind multiple imaging is the area-shot concept, where the original seismic data serves as a secondary area source for multiples. Unlike reflection imaging, which uses point-source propagation, multiple imaging uses complex seismic data. When incoherent multiples cross-correlate, this can lead to severe crosstalk noise, obscuring effective imaging events. Effectively suppressing this crosstalk noise is crucial for the practical application of multiple imaging.

[0005] Currently, the mainstream noise suppression strategies include the following two: 1) Based on the idea of ​​multiple wave ordering, each order of multiple waves is migrated separately and weighted superposition is performed. However, this method requires the pre-extraction of multiple waves of each order and imaging them separately, and requires processing data that is several times more than conventional reflection waves or multiple wave imaging, resulting in an exponential increase in the amount of migration calculations; 2) The least squares migration method is used to treat multiple wave imaging as a linear inverse problem, and the crosstalk noise is effectively suppressed by iteratively updating the approximate Hessian matrix. However, when performing least squares migration, each iteration requires migration and de-migration, which is computationally expensive. In addition, the complex formation mechanism of crosstalk noise will lead to serious inversion ambiguity, which slows down the convergence speed and poor suppression effect. In turn, more iterations are required to obtain multiple wave imaging results with better accuracy, further reducing computational efficiency.

[0006] Therefore, in order to solve the above problems, it is urgent to expand the theory and methods of multiple wave data compression. Summary of the Invention

[0007] In view of this, the present invention provides a seismic multiple wave data compression method and system, aiming to achieve a significant improvement in computational efficiency in compression, thereby enhancing the practicality of multiple wave imaging technology.

[0008] In order to achieve the above object, the present invention adopts the following technical solutions:

[0009] In a first aspect, the present application discloses a method for compressing seismic multiple wave data, the steps comprising:

[0010] determining actual reflected waves based on original seismic wave data;

[0011] Estimate multiple waves of each order based on actual reflected waves;

[0012] Calculate forward multiple wave super gathers and reverse multiple wave super gathers based on the multiple waves;

[0013] The forward and reverse multiple plane waves are calculated after linearly delaying the forward and reverse multiple super gathers.

[0014] Preferably, determining the actual reflected wave based on the original seismic wave data includes:

[0015] Predict multiple waves by multiplying the original seismic wave data by itself in the frequency domain;

[0016] Estimate true multiple waves based on hybrid norm matching algorithm;

[0017] Separate the true multiple waves from the original seismic wave data to obtain the actual reflected waves.

[0018] Preferably, based on the hybrid norm matching algorithm, the true multiple waves are estimated according to the following formula:

[0019]

[0020] Where, represents the adaptive filter to be obtained, Represents the adaptive filter The objective function, represents the estimated true multiples, and denote the L2 and L1 norms respectively, represents the adaptive weighting coefficient, represents the original seismic wave data, Represents the predicted multiples.

[0021] Preferably, estimating each order of multiple waves based on the actual reflected waves comprises:

[0022] Predict multiple waves by multiplying the actual reflected wave by itself in the frequency domain;

[0023] Based on the hybrid norm matching algorithm, the true multiple waves of each order are estimated.

[0024] Preferably, calculating the forward propagation multiple wave super gather and the reverse propagation multiple wave super gather based on the multiple waves includes:

[0025]

[0026] Where, Represents the orthodox multi-wave super-dao collection, Represents the reverse transmission of multiple wave super-dao collection, represents the location of the receiver point, t represents the time of the earthquake record, k represents the shot number, i represents the multiple wave index in the forward multiple wave super-trace set, j represents the multiple wave index in the reverse multiple wave super-trace set, and N represents the total order of multiple waves. represents the amplitude gain probability distribution function, It represents the amplitude gain value of the i-th order multiple wave in the k-th shot, and the random value range is -0.5~0.5. It represents the amplitude gain value of the j-order multiple wave in the k-th shot, sgn represents the random polarity reversal function, which takes the value of 1 or -1. Represents the discrete time series probability distribution function, which follows the standard normal distribution, where represents the time delay value of the i-th order forward multiple wave in the k-th shot, represents the time delay value of the j-order backpropagating multiple wave in the k-th shot, and They represent the forward i-order and reverse j-order multiple waves in the time domain respectively.

[0027] Preferably, linearly delaying the forward propagation multiple wave super gather and the reverse propagation multiple wave super gather comprises:

[0028] The time series linear distribution function is used to time delay the multi-shot multiple wave super gather according to the following formula;

[0029]

[0030] Where, represents the linearly delayed forward multiple super gather, represents the linearly delayed backpropagation multiple supergather, represents the location of the detection point, t represents the time of earthquake recording, k represents the shot number, Indicates the source location, linear time delay value varies linearly with the source position, where the slope , represents the surface incidence angle of the multiple plane wave, and v represents the surface velocity.

[0031] Preferably, the forward and reverse propagation multiple plane wave records are calculated according to the following formula:

[0032]

[0033] Where, represents the multiple plane wave record of the forward pass, represents the back-propagating multiple plane wave record, represents the location of the detection point, t represents the time of earthquake recording, is the slope of the linear delay, P is the number of slopes, S is the number of forward multiple wave super gathers / reverse multiple wave super gathers compressed once, k is the super gather index, represents the amplitude gain probability distribution function, with a random value range of -0.5~0.5, where The slope is The amplitude gain value of the forward and reverse propagation multiple wave super gather of the kth shot at time sgn represents the random polarity reversal function, which takes the value of 1 or -1. Represents the discrete time series probability distribution function, which follows the standard normal distribution. Represents the slope of The time delay value of the forward and reverse propagation multiple wave super gathers of the kth shot, represents the linearly delayed forward multiple super gather, represents the back-propagated multiple supergather after linear time delay.

[0034] In a second aspect, the present application further provides a seismic multiple wave data compression system, which applies any of the above-described seismic multiple wave data compression methods, including:

[0035] A reflection wave extraction module is used to determine the actual reflection wave based on the original seismic wave data;

[0036] A multiple wave determination module is used to estimate each order of multiple waves based on the actual reflected waves;

[0037] A primary compression module is used to calculate forward multiple wave super gathers and reverse multiple wave super gathers based on the multiple waves;

[0038] The secondary compression module is used to calculate the forward and reverse multiple wave plane waves after linearly delaying the forward multiple wave super gather and the reverse multiple wave super gather.

[0039] It can be seen from the above technical solutions that, compared with the prior art, the present invention provides a method and system for compressing seismic multiple wave data, which are mainly used for compressing seismic multiple wave data;

[0040] The compression method provided in this application can significantly reduce the amount of migration data and drastically cut the migration calculation cost, thereby significantly improving the calculation efficiency and the practicality of multiple wave imaging technology. BRIEF DESCRIPTION OF THE DRAWINGS

[0041] In order to more clearly illustrate the embodiments of the present invention or the technical solutions in the prior art, the following briefly introduces the drawings required for use in the embodiments or the description of the prior art. Obviously, the drawings described below are merely embodiments of the present invention. For ordinary technicians in this field, other drawings can be obtained based on the provided drawings without paying any creative work.

[0042] Figure 1 A flow chart of the seismic multiple wave data compression method provided by the present invention;

[0043] Figure 2 A flow chart of an embodiment of the seismic multiple wave data compression method provided by the present invention;

[0044] Figure 3 A schematic diagram of the velocity field of the complex model provided by the present invention;

[0045] Figure 4 Various types of data for the complex model provided by the present invention: (a) original data, (b) multiple wave data, (c) reflected wave data, (d) first-order multiple wave data, (e) second-order multiple wave data, and (f) third-order multiple wave data.

[0046] Figure 5 The complex model provided by the present invention, (a) is the forward propagation multiple wave super gather; (b) is the reverse propagation multiple wave super gather;

[0047] Figure 6 The multiple wave plane wave data of the complex model provided by the present invention when the incident angle is 0, wherein (a) is the forward propagation multiple wave plane wave data; (b) is the reverse propagation multiple wave plane wave data. DETAILED DESCRIPTION

[0048] The following will clearly and completely describe the technical solutions in the embodiments of the present invention in conjunction with the accompanying drawings. Obviously, the described embodiments are only part of the embodiments of the present invention, not all of the embodiments. Based on the embodiments of the present invention, all other embodiments obtained by ordinary technicians in this field without making creative efforts are within the scope of protection of the present invention.

[0049] The following description sets forth many specific details to facilitate a full understanding of the present invention. However, the present invention may also be implemented in other ways different from those described herein, and those skilled in the art may make similar generalizations without violating the scope of the present invention. Therefore, the present invention is not limited to the specific embodiments disclosed below.

[0050] The embodiment of the present invention first discloses a method for compressing seismic multiple wave data. Figure 1 , the steps include,:

[0051] determining actual reflected waves based on original seismic wave data;

[0052] Estimate multiple waves of each order based on actual reflected waves;

[0053] Calculate forward multiple wave super gathers and reverse multiple wave super gathers based on the multiple waves;

[0054] The forward and reverse multiple plane waves are calculated after linearly delaying the forward and reverse multiple super gathers.

[0055] In a specific embodiment, referring to Figure 2 , Figure 2 This is a flow chart of an embodiment of seismic multiple wave data compression. In the figure, S is the number of shot gathers after primary compression, P is the number of multiple wave plane wave records, and P is much smaller than S. In specific implementation, the execution process of each step is as follows:

[0056] Step 1: Determine the actual reflected wave based on the original seismic wave data. In this embodiment, the steps include:

[0057] 1.1 Predict multiple waves by multiplying the original seismic wave data by itself in the frequency domain; in the frequency domain, multiple waves can be expressed as the product of the original data by itself, that is,

[0058]

[0059] Where, Represents the original data; The bands represent predicted multiples;

[0060] 1.2 Based on the hybrid norm matching algorithm, the true multiple waves are estimated, including the true amplitude and actual phase of the multiple waves. After obtaining the predicted multiple waves, this application uses the hybrid L1-L2 norm matching algorithm to estimate the true multiple waves. The formula is expressed as:

[0061]

[0062] Where, represents the adaptive filter to be obtained, and denote the L2 and L1 norms respectively, Represents the adaptive weighting coefficient.

[0063] Furthermore, the least squares optimization algorithm is used to solve the above equation, so that The modulus of is the smallest, and the optimal adaptive matched filter is obtained. Based on this, the estimated multiple waves are expressed as:

[0064]

[0065] Where, Represents the estimated multiple waves. The matching process is first performed in the full spatiotemporal domain and then in the local spatiotemporal window.

[0066] 1.3 Separate the true multiple waves from the original seismic wave data to obtain the actual reflection wave. That is, the reflection wave data is represented by the difference between the original seismic wave data and the estimated multiple wave data:

[0067]

[0068] Where, Represents actual reflected wave data.

[0069] Step 2: Estimate the multiple waves of each order based on the actual reflected waves;

[0070] 2.1 This embodiment predicts 1st to Nth order multiples based on the frequency domain self-product of the reflected wave; where N is an integer greater than 1, and the 0th order multiple is the reflected wave;

[0071] The reflected wave is taken as the 0th-order multiple wave, the first self-product, and the multiple order is increased by 1 (the first self-product obtains the 1st-order multiple wave, the second self-product obtains the 2nd-order multiple wave, and so on). The nth self-product of the reflected wave corresponds to the nth-order water layer multiple wave, where n is an integer greater than 0 and less than N. The calculation formula is:

[0072]

[0073] in, represents the predicted i-order multiple wave; Represents the (i+1)th power of the reflected wave.

[0074] 2.2 Based on the hybrid norm matching algorithm, the true multiple waves of each order are estimated.

[0075] This embodiment uses the optimal adaptive matched filter obtained in step 1.2 to obtain the estimated multiple wave data of each order, namely:

[0076]

[0077] Where, is the estimated i-order multiple wave.

[0078] Step 3: Calculate forward multiple super gathers and reverse multiple super gathers based on the multiple waves to achieve primary compression; wherein the 0th to N-1th order multiple waves correspond to the forward multiple super gathers, and the 1st to Nth order multiple waves correspond to the reverse multiple super gathers;

[0079] In one embodiment, the forward propagation multiple super gathers corresponding to the 0th to N-1th order multiple waves and the reverse propagation multiple super gathers corresponding to the 1st to Nth order multiple waves are respectively calculated using a discrete time series probability distribution function, a random polarity reversal function, and an amplitude gain probability distribution function.

[0080] In the time domain, the calculation formula for the forward multiple wave super gather domain and the reverse multiple wave super gather domain is:

[0081]

[0082] Where, Represents the orthodox multi-wave super-dao collection, Represents the reverse transmission of multiple wave super-dao collection, represents the location of the detection point, t represents the time of the earthquake record, k represents the shot number, i represents the multiple wave index in the forward multiple wave super-trace set, j represents the multiple wave index in the reverse multiple wave super-trace set, represents the amplitude gain probability distribution function, It represents the amplitude gain value of the i-th order multiple wave in the k-th shot, and the random value range is -0.5~0.5. It represents the amplitude gain value of the j-order multiple wave in the k-th shot, sgn represents the random polarity reversal function, which takes the value of 1 or -1. Represents the discrete time series probability distribution function, which follows the standard normal distribution, where represents the time delay value of the i-th order forward multiple wave in the k-th shot, represents the time delay value of the j-order backpropagating multiple wave in the k-th shot, and They represent the forward i-order and reverse j-order multiple waves in the time domain respectively.

[0083] Through the above method, the amount of data required for multiple wave imaging can be compressed to 1 / N of the original amount, realizing the first compression of seismic multiple wave data.

[0084] Step 4, performing linear time delay on the forward propagation multiple wave super gather and the reverse propagation multiple wave super gather to calculate the forward propagation multiple wave plane waves;

[0085] The above method is limited to the single-shot level, meaning the number of compressed multiple supergathers equals the number of original shot gathers, and the migration computational efficiency is comparable to conventional algorithms. Due to the high computational cost of extracting multiples of various orders, further compression of the migrated data is necessary to improve the practicality of multiple imaging. However, data compression achieved by randomly adjusting time delay values, polarity, and amplitude is limited to the single-shot level of multiple imaging. When processing multi-shot multiple data, this can lead to a significant mismatch between the simulated and observed data during the inversion process.

[0086] Therefore, in a preferred embodiment, the present application performs a time delay operation on the multi-shot multiple wave super gathers by using a time series linear distribution function, which is expressed as follows:

[0087]

[0088] Where, represents the linearly delayed forward multiple super gather, represents the linearly delayed backpropagation multiple supergather, represents the location of the detection point, t represents the time of earthquake recording, k represents the shot number, Indicates the source location, linear time delay value varies linearly with the source position, where the slope , represents the surface incidence angle of the multiple plane wave, and v represents the surface velocity.

[0089] To further optimize the above technical solution, the multiple wave super gathers after linear time delay are superimposed and combined through the discrete time series probability distribution function, random polarity reversal function, and amplitude gain probability distribution function to generate multiple wave plane wave records;

[0090] That is, in this embodiment, P forward propagation multiple plane wave records and P reverse propagation multiple plane wave records are calculated according to the following formulas, where P is the number of multiple plane wave records;

[0091]

[0092] Where, represents the multiple plane wave record of the forward pass, represents the back-propagating multiple plane wave record, represents the location of the detection point, t represents the time of earthquake recording, is the slope of the linear delay, P is the number of slopes, that is, P multiple wave super gathers will be generated, which is usually dozens, much smaller than the number of original shot gathers S. S is the number of forward multiple wave super gathers / reverse multiple wave super gathers compressed once. It is usually hundreds or even thousands, k is the super gather index, represents the amplitude gain probability distribution function, with a random value range of -0.5~0.5, where The slope is The amplitude gain value of the forward and reverse propagation multiple wave super gather of the kth shot at time sgn represents the random polarity reversal function, which takes the value of 1 or -1. Represents the discrete time series probability distribution function, which follows the standard normal distribution. Represents the slope of The time delay value of the forward and reverse propagation multiple wave super gathers of the kth shot, represents the linearly delayed forward multiple super gather, represents the back-propagated multiple supergather after linear time delay.

[0093] Using the above formula, S shot multiple-wave supergathers (hundreds or even thousands) can be compressed into P (generally dozens) multiple-wave plane wave records, thereby further compressing the migration data to two orders of magnitude, greatly reducing the amount of migration data and significantly improving the efficiency of migration calculations.

[0094] This application adds an amplitude gain probability distribution function based on the randomly set time delay values ​​and polarities of multiple wave data of different orders, further enhancing the randomness of the superposition of multiple waves of each order.

[0095] In another embodiment, the present application provides a seismic multiple wave data compression system, which applies any of the above-described seismic multiple wave data compression methods, including:

[0096] A reflection wave extraction module is used to determine the actual reflection wave based on the original seismic wave data;

[0097] A multiple wave determination module is used to estimate each order of multiple waves based on the actual reflected waves;

[0098] A primary compression module is used to calculate forward multiple wave super gathers and reverse multiple wave super gathers based on the multiple waves;

[0099] The secondary compression module is used to calculate the forward and reverse multiple wave plane waves after linearly delaying the forward multiple wave super gather and the reverse multiple wave super gather.

[0100] The specific implementation process of each module is the same as the method steps in the above embodiment and will not be repeated here.

[0101] In order to further demonstrate the superiority of the present invention, the seismic multiple wave data compression method of the present invention is applied to Figure 3 The data of the complex Pluto 1.5 model is compressed. The original published data consists of 1387 shots, which are thinned out every 6 shots, resulting in a total of 232 shots for processing. The 76th shot is selected for various data types. Figure 4 As shown, Figure 4 The original data (a) is subjected to frequency domain self-convolution to predict the multiple wave data, and the application of the hybrid norm matching algorithm can obtain the estimated multiple wave ( Figure 4 middle (b); reflected wave data ( Figure 4 (c) is obtained by subtracting the multiple wave data from the original data; the present invention can predict the n-order multiple wave data through the n-order frequency domain self-convolution of the reflected wave data, and then use the mixed norm matching algorithm to obtain the estimated n-order multiple wave, such as Figure 4 The first-order multiple waves shown in (d) are Figure 4The second-order multiple waves shown in (e) are Figure 4 The third-order multiple waves shown in (f).

[0102] Conventional multiple imaging uses the original data as a surface source and multiples as boundary conditions for reverse continuation of the wavefield. When incoherent multiples cross-correlate, complex crosstalk noise is generated, severely degrading imaging quality. Extracting and migrating multiples of various orders effectively suppresses crosstalk noise, but the computational cost scales linearly with the order of the multiples considered for imaging—that is, the amount of migrated data. Considering third-order multiples for imaging triples the computational cost. To address this, compression of multiples of different orders is necessary.

[0103] Figure 5 This paper presents forward and reverse multiple supergathers calculated using the discrete time series probability distribution function, random polarity reversal function, and amplitude gain probability distribution function designed in this invention. This invention effectively processes first-order multiples, resulting in different arrival times, amplitudes, and polarities compared to the first-order multiples in the original data. This also applies to second- and third-order multiples, which will not be further detailed here.

[0104] Furthermore, this project applied the time series linear distribution function to the 232-shot multiple wave super gather data. On this basis, the discrete time series probability distribution function, random polarity reversal function, and amplitude gain probability distribution function were used to set the random delay value, random polarity, and random amplitude of the multiple wave super gather, and the random time delay value, random polarity, and random amplitude of the multiple wave super gather were obtained. Figure 6 The forward and reverse multiple plane wave records for an incident angle of 0 are shown in Figure 1 , where (a) is the forward plane wave and (b) is the reverse plane wave. There are 31 multiple plane wave records in total. Compared to the original 1387-shot data and the thinned 232-shot data, the amount of data required for migration is reduced by two and one orders of magnitude, respectively, significantly improving migration calculation efficiency. Figure 6 The data in the data are randomly distributed, which is mainly due to the fusion of time series linear distribution function, discrete time series probability distribution function, random polarity inversion function, amplitude gain probability distribution function, etc., which is consistent with the data compression theory of this patent.

[0105] The various embodiments in this specification are described in a progressive manner, with each embodiment focusing on the differences from other embodiments. Reference can be made to the common and similar parts between the various embodiments. For the devices disclosed in the embodiments, since they correspond to the methods disclosed in the embodiments, the description is relatively simple, and the relevant parts can be referred to the method description.

[0106] The above description of the disclosed embodiments is intended to enable one skilled in the art to implement or use the present invention. Various modifications to these embodiments will be readily apparent to one skilled in the art, and the general principles defined herein may be implemented in other embodiments without departing from the spirit or scope of the present invention. Therefore, the present invention is not limited to the embodiments shown herein but is intended to conform to the widest scope consistent with the principles and novel features disclosed herein.

Claims

1. A method for compressing seismic multiple wave data, characterized in that: determining actual reflected waves based on original seismic wave data; Estimate multiple waves of each order based on actual reflected waves; Calculate forward multiple wave super gathers and reverse multiple wave super gathers based on the multiple waves; Calculate the forward and reverse multiple plane waves after linearly delaying the forward and reverse multiple super gathers; Linear time delay of forward multiples super gathers and reverse multiples super gathers, including: The time series linear distribution function is used to time delay the multi-shot multiple wave super gather according to the following formula; Where, represents the linearly delayed forward multiple super gather, represents the linearly delayed backpropagation multiple wave supergather, x g represents the location of the detection point, t represents the time of earthquake recording, k represents the shot number, x k Indicates the source location, linear time delay value lx k It varies linearly with the focal position, where the slope l = sinβ / v, β represents the surface incidence angle of the multiple plane wave, and v represents the surface velocity; Calculate the forward and reverse propagation multiple plane wave records according to the following formula: Where, d pf represents the positive multiple plane wave record, d pb represents the back-propagating multiple plane wave record, x g represents the location of the detection point, t represents the time of earthquake recording, l P is the slope of the linear delay, P is the number of slopes, S is the number of forward multiple wave super gathers / reverse multiple wave super gathers compressed once, k is the super gather index, ψ represents the amplitude gain probability distribution function, where ψ(l P ; k) represents the slope of l P The amplitude gain value of the forward and reverse propagation multiple wave super gather of the kth shot at time sgn represents the random polarity reversal function, φ represents the discrete time series probability distribution function, φ(l P ; k) represents the slope of l P The time delay value of the forward and reverse propagation multiple wave super gathers of the kth shot, represents the linearly delayed forward multiple super gather, represents the back-propagated multiple super gather after linear time delay.

2. The seismic multiple wave data compression method according to claim 1, characterized in that: Determine the actual reflected wave based on the original seismic wave data, including: Predict multiple waves by multiplying the original seismic wave data by itself in the frequency domain; Estimate true multiple waves based on hybrid norm matching algorithm; Separate the true multiple waves from the original seismic wave data to obtain the actual reflected waves.

3. The seismic multiple wave data compression method according to claim 2, characterized in that: Based on the hybrid norm matching algorithm, the true multiple waves are estimated according to the following formula; J(κ)=γ||D0-κD m ′||2+(1-γ)||D0-κD m ′||1 In the formula, κ represents the adaptive filter to be obtained, J(κ) represents the objective function for obtaining the adaptive filter κ, and κD m ′ represents the estimated true multiple wave, ||||2 and ||||1 represent the L2 and L1 norms respectively, γ represents the adaptive weighting coefficient, D0 represents the original seismic wave data, D m ′ represents the predicted multiple waves.

4. The seismic multiple wave data compression method according to claim 1, characterized in that: Estimates multiples of each order based on actual reflections, including: Predict multiple waves by multiplying the actual reflected wave by itself in the frequency domain; Based on the hybrid norm matching algorithm, the true multiple waves of each order are estimated.

5. The seismic multiple wave data compression method according to claim 1, characterized in that: Calculate forward multiple supergathers and reverse multiple supergathers based on multiple waves, including: Where, d f Represents the positive transmission multiple wave super gather, d b represents the back-propagated multiple wave super gather, x g represents the location of the receiver point, t represents the time of the earthquake record, k represents the shot number, i is the multiple wave index in the forward multiple wave super-channel set, j is the multiple wave index in the reverse multiple wave super-channel set, N is the total order of multiple waves, ψ represents the amplitude gain probability distribution function, ψ i (k) represents the amplitude gain of the i-th order multiple wave in the k-th shot, ψ j-1 (k) represents the amplitude gain value of the back-propagated j-order multiple wave in the k-th shot, sgn represents the random polarity reversal function, ψ represents the discrete time series probability distribution function, which complies with the standard normal distribution, where φ i (k) represents the time delay of the i-th order forward multiple wave in the k-th shot, φ j-1 (k) represents the time delay of the jth order backpropagating multiple wave in the kth shot, d m,i with d m,j They represent the forward i-order and reverse j-order multiple waves in the time domain respectively.

6. A seismic multiple wave data compression system, characterized in that: The system using the seismic multiple wave data compression method according to any one of claims 1 to 5 comprises: A reflection wave extraction module is used to determine the actual reflection wave based on the original seismic wave data; A multiple wave determination module is used to estimate each order of multiple waves based on the actual reflected waves; A primary compression module is used to calculate forward multiple wave super gathers and reverse multiple wave super gathers based on the multiple waves; The secondary compression module is used to calculate the forward and reverse multiple wave plane waves after linearly delaying the forward multiple wave super gather and the reverse multiple wave super gather.

Images