Full waveform inversion method and device based on local fast matching decomposition

Abstract

This application proposes a full-waveform inversion method and apparatus based on local rapid matching decomposition. According to this application, accurate subscaled seismic data can be obtained through local rapid matching seismic data decomposition, with the advantages of high time-frequency resolution, no aliasing, low computational effort, and low storage requirements. Furthermore, the subscaled seismic data is used for step-by-step full-waveform inversion, enabling the inversion of velocity models from large to small scales, avoiding the problem of full-waveform inversion easily falling into local extremes, and obtaining a high-precision velocity model.

Description

Technical Field

[0001] The present invention relates to the field of seismic velocity modeling and seismic imaging in oil and gas exploration and development, and more particularly to a full waveform inversion method and device based on local fast matching decomposition. Background Art

[0002] Full waveform inversion is a high-precision seismic velocity modeling method, but it has strong nonlinearity. When the initial model is poor, it is easy to fall into the local extreme value problem during the inversion process, making it difficult to obtain good inversion results.

[0003] Existing approaches to this problem generally involve improving the accuracy of the initial model through other modeling methods. However, in practice, under complex geological conditions, ray-based modeling methods struggle to produce a high-precision initial model, and direct full-waveform inversion can still lead to local extrema. Summary of the Invention

[0004] In view of this, the present application proposes a full waveform inversion method based on local fast matching decomposition, and the present application also proposes corresponding devices, electronic devices, and computer-readable storage media.

[0005] According to one aspect of the present application, a full waveform inversion method based on local fast matching decomposition is proposed, the method comprising:

[0006] Step 1: Construct a seismic wavelet library, which includes multiple basic wavelets with different frequencies and time shifts. i represents the scale number, ω i refers to the frequency of the scale numbered i, and t represents time;

[0007] Step 2: Perform Hilbert transform on the seismic data d(t) to obtain the envelope of the seismic data d(t), and determine the time T corresponding to the peak of the envelope. h , calculate the earthquake data d(t) at T h The instantaneous frequency f h , and determine T h and f h The corresponding basic wavelet matching range is h=1, 2, ..., H, where H represents the number of peaks of the envelope;

[0008] Step 3: For each basic wavelet matching range, the seismic data d(t) is matched with each basic wavelet within the basic wavelet matching range to obtain a matching wavelet, and then the matching wavelet is screened according to different frequency scales to perform signal reconstruction to obtain seismic data with different frequency scales.

[0009] Step 4: Seismic data based on different frequency scales Substitute the seismic data matching decomposition function into the full waveform inversion target functional to obtain the step-by-step full waveform inversion target functional J(m (i) ) as the objective function:

[0010]

[0011] Among them, m (i) represents the velocity model of scale number i, Represents frequency as ω i Wavelet library, α i Represents the matching coefficient, G(m (i) ) represents the forward simulation data, Indicates finding the square of the L2 norm;

[0012] Step 5: For each scale numbered i, the objective function J(m (i) ) is iteratively solved to obtain the velocity model m of scale numbered i (i) ;

[0013] Step 6, select ω i , and according to ω i The speed model m is adjusted step by step from small to large (i) Perform full waveform inversion to obtain the final inversion result.

[0014] In some embodiments, in step 1, constructing a seismic wavelet library specifically includes:

[0015] The basic wavelets with different frequencies and time shifts are generated based on the following formula

[0016]

[0017] In some embodiments, in step 2, T is determined h and f h The corresponding basic wavelet matching range specifically includes:

[0018] T h The time window with the center being the preset time length and the width being the preset time length is the time window of the basic wavelet matching range;

[0019] f h The frequency band with the center as the preset frequency range and the width as the preset frequency range is the frequency band of the basic sub-wave matching range.

[0020] In some embodiments, in step 3, seismic data with different frequency scales are obtained. Specifically include:

[0021] Assume T h and f hThe corresponding basic wavelet matching range includes L basic wavelets For each basic wavelet Calculated based on the following formula:

[0022]

[0023] Among them, t j Refers to the instantaneous frequency of all H peaks f h The time corresponding to the peak value, J represents the instantaneous frequency of H peaks f h The number of peaks;

[0024] Select the one with the largest amplitude As

[0025] In some implementations, in step 5, the objective function is iteratively solved based on the following formula until the objective function is satisfied:

[0026] m (i) k =m (i) k-1 +Δm (i) k-1

[0027]

[0028] Where k represents the number of iterations, represents the result of the velocity model of scale number i after the kth iteration, a k-1 represents the iteration step, U (i) Indicates the forward wave field, G * Indicates TRON backpropagation.

[0029] In some embodiments, in step 6, ω is selected specifically by the following formula: i :

[0030]

[0031] Among them, z depth is the model depth, x offset is the offset distance.

[0032] According to another aspect of the present application, a full waveform inversion device based on local fast matching decomposition is also proposed, the device comprising:

[0033] The wavelet library construction unit is used to construct a seismic wavelet library, which includes multiple basic wavelets with different frequencies and time shifts. i represents the scale number, ω i refers to the frequency of the scale numbered i, and t represents time;

[0034] The matching range determination unit is used to perform Hilbert transform on the seismic data d(t) to obtain the envelope of the seismic data d(t) and determine the time T corresponding to the peak of the envelope. h , calculate the earthquake data d(t) at T h The instantaneous frequency f h , and determine T h and f h The corresponding basic wavelet matching range is h=1, 2, ..., H, where H represents the number of peaks of the envelope;

[0035] The scaled data acquisition unit is used to match the seismic data d(t) with each of the basic wavelets within the basic wavelet matching range to obtain matching wavelets, and then reconstruct the matching wavelets according to different frequency scales to obtain seismic data with different frequency scales.

[0036] Objective function determination unit for seismic data based on different frequency scales Substitute the seismic data matching decomposition function into the full waveform inversion target functional to obtain the step-by-step full waveform inversion target functional J(m (i) ) as the objective function:

[0037]

[0038] Among them, m (i) represents the velocity model of scale number i, Represents frequency as ω i Wavelet library, α i Represents the matching coefficient, G(m (i) ) represents the forward simulation data, Indicates finding the square of the L2 norm;

[0039] The scale iteration unit is used to calculate the objective function J(m (i) ) is iteratively solved to obtain the velocity model m of scale numbered i (i) ;

[0040] Full waveform inversion unit, select ω i , and according to ω i The speed model m is adjusted step by step from small to large (i) Perform full waveform inversion to obtain the final inversion result.

[0041] In some embodiments, the scaled data acquisition unit is specifically configured to:

[0042] Assume T h and fh The corresponding basic wavelet matching range includes L basic wavelets For each basic wavelet Calculated based on the following formula:

[0043]

[0044] Among them, t j Refers to the instantaneous frequency of all H peaks f h The time corresponding to the peak value, J represents the instantaneous frequency of H peaks f h The number of peaks;

[0045] Select the one with the largest amplitude As

[0046] According to another aspect of the present application, an electronic device is further provided, comprising:

[0047] a memory storing executable instructions;

[0048] A processor runs the executable instructions in the memory to implement the full waveform inversion device method based on local fast matching decomposition as described above.

[0049] According to another aspect of the present application, a computer-readable storage medium is provided, which stores a computer program. When the computer program is executed by a processor, the full waveform inversion method based on local fast matching decomposition as described above is implemented.

[0050] The technical solution proposed in this application obtains accurate sub-scale seismic data through local fast matching seismic data decomposition, which has the advantages of high time-frequency resolution and no false frequency generation, and requires small amount of calculation and low storage requirements; at the same time, the sub-scale seismic data is used for step-by-step full-wave inversion, which can invert the velocity model step by step from large scale to small scale, avoiding the problem that full waveform inversion is prone to falling into local extreme values, and obtaining a high-precision velocity model. BRIEF DESCRIPTION OF THE DRAWINGS

[0051] The above and other objects, features and advantages of the present application will become more apparent through a more detailed description of exemplary embodiments of the present application in conjunction with the accompanying drawings, wherein the same reference numerals generally represent the same components in the exemplary embodiments of the present application.

[0052] Figure 1 A flowchart of a full waveform inversion method based on local fast matching decomposition according to an embodiment of the present application is shown.

[0053] Figure 2A structural block diagram of a full waveform inversion device based on local fast matching decomposition according to an embodiment of the present application is shown.

[0054] Figure 3 A certain original velocity model is shown.

[0055] Figure 4 Shown is the initial velocity model used for full waveform inversion.

[0056] Figure 5 Show the use Figure 3 The raw velocity model shown is forward modeled to obtain a shot record from the seismic record.

[0057] Figure 6 Shown are subscaled seismic data obtained according to the present application.

[0058] Figure 7 The full waveform inversion results obtained according to the present application are shown. DETAILED DESCRIPTION

[0059] The preferred embodiments of the present application will be described in more detail below with reference to the accompanying drawings. Although the preferred embodiments of the present application are shown in the accompanying drawings, it should be understood that the present application can be implemented in various forms and should not be limited by the embodiments described herein. Instead, these embodiments are provided to make the present application more thorough and complete, and to fully convey the scope of the present application to those skilled in the art.

[0060] According to one aspect of the present application, a full waveform inversion method based on local fast matching decomposition is provided.

[0061] See Figure 1 . Figure 1 The flowchart of the full waveform inversion method based on local fast matching decomposition according to one embodiment of the present application is shown. As shown in the figure, the method includes steps 1 to 5.

[0062] Step 1: Construct a seismic wavelet library, which includes multiple basic wavelets with different frequencies and time shifts. i represents the scale number, ω i refers to the frequency of the scale numbered i, and t represents time.

[0063] The seismic wavelet library is composed of basic wavelets with different frequencies and time shifts. In some embodiments, the basic wavelets can be generated according to the following formula:

[0064]

[0065] Step 2: Perform Hilbert transform on the seismic data d(t) to obtain the envelope of the seismic data d(t), and determine the time T corresponding to the peak of the envelope.h , calculate the earthquake data d(t) at T h The instantaneous frequency f h , and determine T h and f h The corresponding basic wavelet matching range is h=1, 2, ..., H, where H represents the number of peaks of the envelope.

[0066] In some embodiments, determining T h and f h The corresponding basic wavelet matching range may specifically include:

[0067] T h The time window with the center being the preset time length and the width being the preset time length is the time window of the basic wavelet matching range;

[0068] f h The frequency band with the center as the preset frequency range and the width as the preset frequency range is the frequency band of the basic sub-wave matching range.

[0069] Step 3: For each basic wavelet matching range, the seismic data d(t) is matched with each basic wavelet within the basic wavelet matching range to obtain a matching wavelet, and then the matching wavelet is screened according to different frequency scales to perform signal reconstruction to obtain seismic data with different frequency scales.

[0070] Since seismic data is only matched with basic wavelets within the basic wavelet matching range, the amount of calculation is significantly reduced and the storage requirement is lowered.

[0071] In some embodiments, seismic data with different frequency scales can be obtained according to the following method: Specifically, they may include:

[0072] Assume T h and f h The corresponding basic wavelet matching range includes L basic wavelets For each basic wavelet Calculated based on the following formula:

[0073]

[0074] Among them, t j Refers to the instantaneous frequency of all H peaks f h The time corresponding to the peak value, J represents the instantaneous frequency of H peaks f h The number of peaks;

[0075] Select the one with the largest amplitude As

[0076] Step 4: Seismic data based on different frequency scales Substitute the seismic data matching decomposition function into the full waveform inversion target functional to obtain the step-by-step full waveform inversion target functional J(m (i) ) as the objective function:

[0077]

[0078] Among them, m (i) represents the velocity model of scale number i, Represents frequency as ω i Wavelet library, α i Represents the matching coefficient, G(m (i) ) represents the forward simulation data, It means to find the square of L2 norm.

[0079] Among them, the seismic data matching decomposition function and the full waveform inversion objective function are common knowledge in this field and will not be described in detail in this application.

[0080] Step 5: For each scale numbered i, the objective function J(m (i) ) is iteratively solved to obtain the velocity model m of scale numbered i (i) .

[0081] In some embodiments, the objective function can be iteratively solved based on the following formula until the objective function is satisfied, that is, the objective function J(m (i) )Minimum velocity model

[0082] m (i) k =m (i) k-1 +Δm (i) k-1

[0083]

[0084] Where k represents the number of iterations, represents the result of the velocity model of scale number i after the kth iteration, a k-1 represents the iteration step, U (i) Indicates the forward wave field, G * Indicates TRON backpropagation.

[0085] Step 6, select ω i , and according to ω i The speed model m is adjusted step by step from small to large (i) Perform full waveform inversion to obtain the final inversion result.

[0086] In some embodiments, ω can be selected by the following formula i :

[0087]

[0088] Among them, z depth is the model depth, x offset is the offset distance.

[0089] In the step-by-step inversion, the iteration result of the previous scale This is the initial velocity model for the next scale Where K represents the final iteration number of the velocity model of the scale numbered (i-1).

[0090] In the above embodiment, accurate sub-scale seismic data is obtained by local fast matching seismic data decomposition, which has the advantages of high time-frequency resolution and no false frequency generation, and requires small amount of calculation and low storage requirements; at the same time, the sub-scale seismic data is used to perform step-by-step full-waveform inversion, which can invert the velocity model step by step from large scale to small scale, avoiding the problem that full waveform inversion is prone to falling into local extreme values, and obtaining a high-precision velocity model.

[0091] According to another aspect of the present application, a full waveform inversion device based on local fast matching decomposition is also provided. Figure 2 . Figure 2 A structural block diagram of a full waveform inversion device based on local fast matching decomposition according to an embodiment of the present application is shown.

[0092] like Figure 2 As shown, the apparatus includes a wavelet library construction unit 202 , a matching range determination unit 204 , a sub-scale data acquisition unit 206 , an objective function determination unit 208 , a sub-scale iteration unit 210 , and a full waveform inversion unit 212 .

[0093] The wavelet library construction unit 202 is used to construct a seismic wavelet library, which includes a plurality of basic wavelets with different frequencies and time shifts. i represents the scale number, ω i refers to the frequency of the scale numbered i, and t represents time;

[0094] The matching range determination unit 204 is used to perform Hilbert transform on the seismic data d(t) to obtain the envelope of the seismic data d(t) and determine the time T corresponding to the peak of the envelope. h , calculate the earthquake data d(t) at T h The instantaneous frequency f h , and determine T h and f hThe corresponding basic wavelet matching range is h=1, 2, ..., H, where H represents the number of peaks of the envelope;

[0095] The scaled data acquisition unit 206 is used to match the seismic data d(t) with each of the basic wavelets within the basic wavelet matching range to obtain matching wavelets, and then reconstruct the matching wavelets according to different frequency scales to obtain seismic data with different frequency scales.

[0096] The objective function determination unit 208 is used to determine the seismic data based on different frequency scales. Substitute the seismic data matching decomposition function into the full waveform inversion target functional to obtain the step-by-step full waveform inversion target functional J(m (i) ) as the objective function:

[0097]

[0098] Among them, m (i) represents the velocity model of scale number i, Represents frequency as ω i Wavelet library, α i Represents the matching coefficient, G(m (i) ) represents the forward simulation data, Indicates finding the square of the L2 norm;

[0099] The scale iteration unit 210 is used to perform an iteration on the objective function J(m (i) ) is iteratively solved to obtain the velocity model m of scale numbered i (i) ;

[0100] Full waveform inversion unit 212, used to select ω i , and according to ω i The speed model m is adjusted step by step from small to large (i) Perform full waveform inversion to obtain the final inversion result.

[0101] In some implementations, the scaled data acquisition unit 206 may be specifically configured to:

[0102] Assume T h and f h The corresponding basic wavelet matching range includes L basic wavelets For each basic wavelet Calculated based on the following formula:

[0103]

[0104] Among them, t jRefers to the instantaneous frequency of all H peaks f h The time corresponding to the peak value, J represents the instantaneous frequency of H peaks f h The number of peaks;

[0105] Select the one with the largest amplitude As

[0106] For other detailed descriptions of this embodiment, please refer to the corresponding descriptions in the aforementioned embodiments, which will not be repeated here.

[0107] According to another aspect of the present application, an electronic device is provided. The electronic device includes:

[0108] a memory storing executable instructions;

[0109] A processor runs the executable instructions in the memory to implement the full waveform inversion method based on local fast matching decomposition as described above.

[0110] Specifically, the memory may include one or more computer program products, which may include various forms of computer-readable storage media, such as volatile memory and / or non-volatile memory. The volatile memory may include, for example, random access memory (RAM) and / or cache memory. The non-volatile memory may include, for example, read-only memory (ROM), a hard disk, a flash memory, etc.

[0111] The processor may be a central processing unit (CPU) or other forms of processing units having data processing capabilities and / or instruction execution capabilities, and may control other components in the electronic device to perform desired functions. In one embodiment of the present application, the processor is used to run the computer-readable instructions stored in the memory.

[0112] For detailed description of this embodiment, please refer to the corresponding description in the aforementioned embodiments, which will not be repeated here.

[0113] According to another aspect of the present application, a computer-readable storage medium is provided, which stores a computer program. When the computer program is executed by a processor, it implements the full waveform inversion method based on local fast matching decomposition as described above.

[0114] According to the computer-readable storage medium of the embodiment of the present application, non-transitory computer-readable instructions are stored thereon. When the non-transitory computer-readable instructions are executed by a processor, all or part of the steps of the aforementioned methods of the embodiments of the present application are executed.

[0115] The above-mentioned computer-readable storage media include, but are not limited to, optical storage media (e.g., CD-ROMs and DVDs), magneto-optical storage media (e.g., MOs), magnetic storage media (e.g., magnetic tapes or mobile hard disks), media with built-in rewritable non-volatile memory (e.g., memory cards), and media with built-in ROM (e.g., ROM cartridges).

[0116] Those skilled in the art should understand that in order to solve the technical problem of how to obtain a good user experience, this embodiment may also include well-known structures such as a communication bus and an interface, and these well-known structures should also be included in the scope of protection of this application.

[0117] For detailed description of this embodiment, please refer to the corresponding description in the aforementioned embodiments, which will not be repeated here.

[0118] The following examples can verify the accuracy of the technical solution disclosed in this application. Figure 3 A certain original velocity model is shown. Figure 4 The initial velocity model used for full-waveform inversion is shown. The initial velocity model has a grid cell size of 201 horizontally and 116 vertically, with a grid spacing of 20 meters horizontally and 20 meters vertically, respectively. Both the source and receivers are located at the surface, with a total of 40 shots and full coverage. The source is a Ricker wavelet with a dominant frequency of 20 Hz, and the time sampling interval is 0.8 ms.

[0119] The earthquake records obtained by forward modeling using the original velocity model are used as observed earthquake records. Figure 5 One of the shots is shown.

[0120] According to this application, the seismic data decomposition of the seismic records is performed locally and quickly, and the obtained sub-scale seismic data is as follows Figure 6 Then, according to this application, full waveform inversion is performed on each scale data step by step, and the velocity model obtained by inversion is as follows Figure 7 shown.

[0121] It can be seen that the velocity model obtained by full waveform inversion according to the present application is in good agreement with the original velocity model, which proves the accuracy of the present application.

[0122] The embodiments of the present application have been described above. The above description is illustrative and not exhaustive, and is not limited to the disclosed embodiments. Many modifications and variations will be apparent to those skilled in the art without departing from the scope and spirit of the described embodiments. The terminology used herein is selected to best explain the principles of the embodiments, practical applications, or technical improvements to existing technologies, or to enable other persons skilled in the art to understand the embodiments disclosed herein.

Claims

1. A full waveform inversion method based on local fast matching decomposition, characterized in that: The method comprises: Step 1: Construct a seismic wavelet library, which includes multiple basic wavelets with different frequencies and time shifts. i represents the scale number, ω i refers to the frequency of the scale numbered i, and t represents time; Step 2: Perform Hilbert transform on the seismic data d(t) to obtain the envelope of the seismic data d(t), and determine the time T corresponding to the peak of the envelope. h , calculate the earthquake data d(t) at T h The instantaneous frequency f h , and determine T h and f h The corresponding basic wavelet matching range, h = 1, 2, ..., H, H represents the number of peaks of the envelope; Step 3: For each basic wavelet matching range, the seismic data d(t) is matched with each basic wavelet within the basic wavelet matching range to obtain a matching wavelet, and then the matching wavelet is screened according to different frequency scales to perform signal reconstruction to obtain seismic data with different frequency scales. Step 4: Seismic data based on different frequency scales Substitute the seismic data matching decomposition function into the full waveform inversion target functional to obtain the step-by-step full waveform inversion target functional J(m (i) ) as the objective function: Among them, m (i) represents the velocity model of scale number i, Represents frequency as ω i Wavelet library, α i Represents the matching coefficient, G(m (i) ) represents the forward simulation data, Indicates the square of the L2 norm; Step 5: For each scale numbered i, the objective function J(m (i) ) is iteratively solved to obtain the velocity model m of scale numbered i (i) ; Step 6, select ω i , and according to ω i The speed model m is adjusted step by step from small to large (i) Perform full waveform inversion to obtain the final inversion result.

2. The method according to claim 1, characterized in that In step 1, building a seismic wavelet library specifically includes: The basic wavelets with different frequencies and time shifts are generated based on the following formula 3. The method according to claim 1, characterized in that In step 2, determine T h and f h The corresponding basic wavelet matching range specifically includes: T h The time window with the center being the preset time length and the width being the preset time length is the time window of the basic wavelet matching range; f h The frequency band with the center as the preset frequency range and the width as the preset frequency range is the frequency band of the basic sub-wave matching range.

4. The method according to claim 1, wherein In step 3, seismic data with different frequency scales are obtained Specifically include: Assume T h and f h The corresponding basic wavelet matching range includes L basic wavelets For each basic wavelet Calculated based on the following formula: Among them, t j Refers to the instantaneous frequency of all H peaks f h The time corresponding to the peak value, J represents the instantaneous frequency of H peaks f h The number of peaks; Select the one with the largest amplitude As 5. The method according to claim 1, characterized in that In step 5, the objective function is iteratively solved based on the following formula until the objective function is satisfied: Where k represents the number of iterations, represents the result of the velocity model of scale number i after the kth iteration, a k-1 represents the iteration step, U (i) Indicates the forward wave field, G * Indicates TRON backpropagation.

6. The method according to claim 1, characterized in that In step 6, ω is selected by the following formula: i : Among them, z depth is the model depth, x offset is the offset distance.

7. An electronic device, characterized in that: The electronic device comprises: a memory storing executable instructions; A processor, wherein the processor runs the executable instructions in the memory to implement the method according to any one of claims 1 to 6.

8. A computer-readable storage medium, characterized in that The computer-readable storage medium stores a computer program, and when the computer program is executed by a processor, the method according to any one of claims 1 to 6 is implemented.

9. A full waveform inversion device based on local fast matching decomposition, characterized in that: The device comprises: The wavelet library construction unit is used to construct a seismic wavelet library, which includes multiple basic wavelets with different frequencies and time shifts. i represents the scale number, ω i refers to the frequency of the scale numbered i, and t represents time; The matching range determination unit is used to perform Hilbert transform on the seismic data d(t) to obtain the envelope of the seismic data d(t) and determine the time T corresponding to the peak of the envelope. h , calculate the earthquake data d(t) at T h The instantaneous frequency f h , and determine T h and f h The corresponding basic wavelet matching range, h = 1, 2, ..., H, H represents the number of peaks of the envelope; The scaled data acquisition unit is used to match the seismic data d(t) with each of the basic wavelets within the basic wavelet matching range to obtain matching wavelets, and then reconstruct the matching wavelets according to different frequency scales to obtain seismic data with different frequency scales. Objective function determination unit for seismic data based on different frequency scales Substitute the seismic data matching decomposition function into the full waveform inversion target functional to obtain the step-by-step full waveform inversion target functional J(m (i) ) as the objective function: Among them, m (i) represents the velocity model of scale number i, Represents frequency as ω i Wavelet library, α i Represents the matching coefficient, G(m (i) ) represents the forward simulation data, Indicates the square of the L2 norm; The scale iteration unit is used to calculate the objective function J(m (i) ) is iteratively solved to obtain the velocity model m of scale numbered i (i) ; Full waveform inversion unit, select ω i , and according to ω i The speed model m is adjusted step by step from small to large (i) Perform full waveform inversion to obtain the final inversion result.

10. The device according to claim 9, characterized in that The sub-scale data acquisition unit is specifically used for: Assume T h and f h The corresponding basic wavelet matching range includes L basic wavelets For each basic wavelet Calculated based on the following formula: Among them, t j Refers to the instantaneous frequency of all H peaks f h The time corresponding to the peak value, J represents the instantaneous frequency of H peaks f h The number of peaks; Select the one with the largest amplitude As

Images