Well drilling seismic data processing method and system and storage medium
By collecting drilling data at the wellhead and performing signal processing technology, the problem of explosives and excitation depth control in seismic exploration is solved, the seismic signal quality and acquisition efficiency are improved, and the cost and risk are reduced.
Patent Information
- Application Number
- CN202411829957.7
- Authority / Receiving Office
- CN · China
- Patent Type
- Applications(China)
- Current Assignee / Owner
- Filing Date
- 2024-12-12
- Publication Date
- 2025-04-18
- Estimated Expiration
- 2044-12-12
AI Technical Summary
In seismic exploration, it is difficult for the existing technology to accurately control the explosive quantity and excitation depth, resulting in low seismic signal quality, high collection cost, and different lithologies affect signal quality, making it difficult to accurately identify formation lithologies.
By laying a detector at the wellhead, collecting drilling data, and using Fourier transform, conjugated complex multiplication, phase adjustment and filtering technologies to process the signals, eliminate noise, extract effective signals, and accurately identify lithologies.
It improves the resolution and quality of seismic signals, reduces the amount of explosives, reduces the acquisition cost, and improves the safety of the acquisition process.
Smart Images

Figure CN119846700A_ABST
Abstract
Description
Technical Field
[0001] The present invention relates to the technical field of data processing, and in particular to a method, system and storage medium for processing seismic data while drilling. Background Art
[0002] The earthquake sources used in seismic exploration are mainly explosive earthquake sources. The amount of explosives is closely related to the location of the explosives, the lithology of the well excitation point and other factors. The excitation energy and excitation frequency are crucial to the signal-to-noise ratio and stratigraphic resolution of seismic waves. In the conventional field production and collection process, the amount of explosives is usually determined at one time based on the depth position of the standard well depth (between 8-20 meters) experiment to determine the data collection requirements of the entire work area. However, the excitation conditions in the field are ever-changing, and it is impossible to scientifically determine the amount of explosives (kilograms of explosives) and the stratum (especially the excitation of sandstone areas at appropriate depths). The explosives used for one excitation are 4-15 kilograms of TNT. In the process of using explosive earthquake sources, due to the non-uniformity of the rock formations at each well location, it is impossible to accurately give the amount of explosives used, and they can only be arranged based on daily experience. If too much explosive is used, the exploration cost will increase, and a large fracture zone will be formed, which will weaken the elasticity of the formation, move the seismic signal frequency band to the low frequency end, and lose the high-frequency signal; if too little explosive is used, the received signal will be weak, and effective high-quality seismic exploration will not be possible. Affected by the plasticity of mudstone and the elasticity of sandstone, when the same mass of explosive is used, there is a significant difference in the quality of the signal received when exciting sandstone and mudstone. This series of problems will lead to high data acquisition costs and low quality of the collected seismic data. Therefore, during the acquisition process, it is necessary to add the influence of stratum lithology on the amount of explosive to the control of the amount of explosive, so as to provide reliable data support for the next step of intelligent explosive research.
[0003] Lithology identification is the basis for studying geological reservoir characteristics, calculating reserves and conducting geological modeling. It can determine the shape and range of the spatial distribution of different lithologies and provide specific geological information. In the actual production process, lithology identification is mainly carried out through methods such as hand specimen identification, thin section identification, element identification, and mineral identification. These methods can stably obtain underground rock information and accurately determine lithology, but they require rich experience and knowledge, are inefficient, and require coring during drilling. In the process of drilling and coring, there are problems such as rock crushing and disordered arrangement order, and after the rock is brought out of the surface, its properties may change significantly. At the same time, it is costly and difficult to take out the core of each well. In addition to the above-mentioned lithology identification technologies, gravity and magnetic technology, seismic technology, and logging technology can also assist in lithology identification, but because the parameter ranges of gravity, magnetization, elasticity, etc. of rocks of different lithologies often overlap to a certain extent, it is difficult to clearly divide them in a single parameter. Therefore, high-resolution data is required to better distinguish lithology, but the cost of high-resolution data acquisition is high, and it is difficult to cover all drilling areas.
[0004] In the process of drilling production, there are tens of thousands of wells, and the vibration of drilling will generate seismic wave energy, which can be collected by high-precision instruments to identify the lithology of the formation. However, the seismic wave signals collected in this way have a large amount of data, a high sampling rate, and complex components (including effective seismic wave signals, vibration of drilling rigs, electrical interference, human interference, etc.), and a large amount of noise interference. The effective signal is submerged in the noise. Conventional seismic data processing methods have great difficulty in extracting effective signals from drilling signals.
[0005] In view of this, there is an urgent need to provide a method for processing while-drilling seismic data for lithology identification in oil and gas exploration. Summary of the invention
[0006] The present invention provides a method, system and storage medium for processing seismic data while drilling to solve one or more of the above problems.
[0007] This specification discloses a method for processing seismic data while drilling, comprising:
[0008] Arrange several geophones in a straight line at the wellhead to collect data while drilling;
[0009] Based on the drilling data, the seismic data signal is extracted and segmented according to its period to obtain m groups of segmented signals;
[0010] Perform Fourier transform on the first group and the mth group of split signals, convert from the time domain to the frequency domain, and perform conjugate complex multiplication on the frequency domains corresponding to the first group and the mth group of split signals to obtain the frequency domain correlation coefficient;
[0011] Perform inverse Fourier transform on the frequency domain correlation coefficient, convert the frequency domain correlation coefficient to the time domain, and obtain the time domain correlation coefficient;
[0012] Based on the time domain correlation coefficient, the largest group of correlation coefficients is extracted to obtain the phase difference of the segmented signal;
[0013] Based on the phase difference of the split signals, the phase of the mth group of split signals is adjusted to be consistent with the phase of the first group of split signals;
[0014] The mth group of split signals after phase adjustment is subjected to superposition processing to eliminate the influence of random noise on the signals and obtain the superpositioned signals;
[0015] A Blackman window is formed according to the filter, and the FIR filter coefficients are designed with the Blackman window so that the weighted square sum of the error between the actual filter response of the filter and the specified ideal frequency response is minimized, and the optimal sampling point weight is obtained;
[0016] Multiply the weight of the optimal sampling point by the Blackman window to obtain the optimal filter;
[0017] Using an optimal filter to filter the superimposed signal to obtain a filtered signal;
[0018] Take a differential of the filtered signal, calculate the change of the signal, and obtain a differential signal;
[0019] The secondary difference of the primary differential signal is obtained to obtain its variation, and finally the processed while drilling data is obtained.
[0020] In this specification, the expression of the drilling data is formula (1);
[0021] S = S0 + N (16);
[0022] In the formula, S is the acquired seismic data, S0 is the acquired effective seismic record, and N is the noise generated during the acquisition process;
[0023] Take out a 2s seismic data signal and divide it according to its period. The expression is formula (2);
[0024] x m (t) = S (T·(m-1) + t) (17);
[0025] In the formula, m represents the group number of the segmentation, x m (t) represents the mth group of signals after segmentation, T represents the signal period, and t represents the number of sampling points of the signal.
[0026] In this specification, the first group and the mth group of segmented signals are subjected to Fourier transform, and converted from the time domain to the frequency domain, and the expressions thereof are formula (3) and formula (4);
[0027]
[0028] Where X1(ω) represents the frequency domain of x1(t), X m (ω) represents x m (t) is the frequency domain, ω represents the amplitude value in the frequency domain, i is the imaginary unit, e is the natural constant, t represents the number of sampling points of the signal, It is the accumulation in the time domain.
[0029] In this specification, X1(ω) and X m The frequency domain correlation coefficient is obtained by multiplying the conjugate complex number of (ω), and its expression is formula (5);
[0030]
[0031] Where C(ω) represents the frequency domain correlation coefficient, For X m The complex conjugate of (ω);
[0032] Perform inverse Fourier transform on C(ω) to convert the frequency domain correlation coefficient into the time domain, which is expressed as formula (6);
[0033]
[0034] Where c(t) represents the time domain correlation coefficient.
[0035] In this specification, the largest group of correlation coefficients calculated in formula (6) is extracted to obtain the phase difference value of the split signal, which is expressed as formula (7);
[0036] (r,p)=max(c)(22);
[0037] In the formula, r represents the maximum value of c, p represents the value of t corresponding to the maximum value of c, max() is the function for finding the maximum value, and c(t) represents the time domain correlation coefficient;
[0038] Adjust x according to the calculation result of formula (7) m The phase of the signal is adjusted to be consistent with the phase of x1, and its expression is formula (8);
[0039]
[0040] In the formula, represents x after phase adjustment m (t).
[0041] In this specification, the signal obtained in formula (8) Perform superposition processing to eliminate the influence of random noise on the signal, and its expression is formula (9);
[0042]
[0043] In the formula, x(t) represents the signal after superposition, Represents x m (t) Data after phase adjustment, n is the total number of groups after data segmentation.
[0044] In this specification, a Blackman window is formed according to the filter, and its expression is formula (10);
[0045]
[0046] Where k represents the number of filter points, u represents the filter length, and w(k) represents the Blackman window function;
[0047] The FIR filter coefficients are designed using the generated Blackman window w(k) so that the weighted sum of squares of the errors between the actual filter response and the specified ideal frequency response is minimized. The objective function is expressed as formula (11);
[0048]
[0049] In the formula, E represents the error sum of squares, h(k) represents the weight of the kth sampling point, and H(e iωk ) represents the actual frequency response at frequency ω, and Dk represents the ideal amplitude response at frequency ω;
[0050] Multiply the optimal h(k) obtained from formula (11) by the Blackman window w(k) to obtain the optimal filter, which is expressed as formula (12);
[0051]
[0052] In the formula, represents the designed filter;
[0053] The superimposed signal is filtered according to the Blackman window in formula (12), and its expression is formula (13);
[0054]
[0055] In the formula, Represents the filtered superimposed signal at time t, that is, the filtered signal.
[0056] In this specification, a difference is taken for the signal in formula (13) to obtain the change in the signal, which is expressed as formula (14);
[0057]
[0058] In the formula, d(t) represents the first difference of the filtered signal at time t, that is, the first difference signal; represents the filtered superposition signal at time t; Represents the filtered superposition signal at time t+1.
[0059] The signal in formula (14) is taken to obtain the second difference, and its variation is obtained, and finally the processed while drilling data is obtained, and its expression is formula (15);
[0060] dd(t)=|d(t+1)||-|d(t)| (30);
[0061] Wherein, dd(t) represents the secondary difference of the filtered signal, i.e. the processed while drilling data; | | is the absolute value symbol; d(t) represents the primary difference of the filtered signal at time t, and d(t+1) represents the primary difference of the filtered signal at time t+1.
[0062] The while-drilling signal calculated in formula (15), i.e. the processed while-drilling data, is compared and analyzed with the core histogram.
[0063] This specification also discloses a drilling-while-drilling seismic data processing system, which is used to implement any one of the drilling-while-drilling seismic data processing methods described above, and the drilling-while-drilling seismic data processing system includes:
[0064] Arrangement acquisition module, used to arrange a number of detectors in a straight line at the wellhead to collect drilling data while drilling;
[0065] A segmentation module is used to extract seismic data signals based on the drilling data, segment them according to their periods, and obtain m groups of segmented signals;
[0066] Transformation module, used to:
[0067] Perform Fourier transform on the first group and the mth group of split signals, convert from the time domain to the frequency domain, and perform conjugate complex multiplication on the frequency domains corresponding to the first group and the mth group of split signals to obtain the frequency domain correlation coefficient;
[0068] Perform inverse Fourier transform on the frequency domain correlation coefficient, convert the frequency domain correlation coefficient to the time domain, and obtain the time domain correlation coefficient;
[0069] Based on the time domain correlation coefficient, the largest group of correlation coefficients is extracted to obtain the phase difference of the segmented signal;
[0070] Phase and superposition modules for:
[0071] Based on the phase difference of the split signals, the phase of the mth group of split signals is adjusted to be consistent with the phase of the first group of split signals;
[0072] The mth group of split signals after phase adjustment is subjected to superposition processing to eliminate the influence of random noise on the signals and obtain the superpositioned signals;
[0073] Filter module for:
[0074] A Blackman window is formed according to the filter, and the FIR filter coefficients are designed with the Blackman window so that the weighted square sum of the error between the actual filter response of the filter and the specified ideal frequency response is minimized, and the optimal sampling point weight is obtained;
[0075] Multiply the weight of the optimal sampling point by the Blackman window to obtain the optimal filter;
[0076] Using an optimal filter to filter the superimposed signal to obtain a filtered signal;
[0077] Differential modules for:
[0078] Take a differential of the filtered signal, calculate the change of the signal, and obtain a differential signal;
[0079] The secondary difference of the primary differential signal is obtained to obtain its variation, and finally the processed while drilling data is obtained.
[0080] This specification also discloses a computer-readable storage medium, wherein the storage medium stores computer instructions. When a computer reads the computer instructions, the computer executes any one of the above-mentioned methods for processing seismic data while drilling.
[0081] The embodiments of this specification can at least achieve the following beneficial effects:
[0082] Taking advantage of the high drilling density in seismic acquisition production, during the implementation of drilling, the energy generated by the drilling rig is used as the excitation source to optimize the ground observation method, and the drilling data is collected with high precision and high density to obtain seismic waves passing through different lithology strata. After processing the obtained seismic waves, the lithology interface is judged according to its morphological characteristics. By comparing the core column chart, this method can accurately describe the lithology changes of the strata during the drilling process, and use this as a guide for seismic data acquisition. The explosives used are reduced by 30%, and the energy efficiency ratio is significantly improved. The main frequency of the collected seismic signals is increased from 20-30Hz to 50-60Hz, the resolution of the seismic signals is significantly improved, and the safety of on-site management of explosives is also effectively improved. BRIEF DESCRIPTION OF THE DRAWINGS
[0083] In order to more clearly illustrate the technical solutions of the embodiments of the present invention, the accompanying drawings required for use in the description of the embodiments will be briefly introduced below. Obviously, the accompanying drawings described below are only some embodiments of the present invention. For ordinary technicians in this field, other accompanying drawings can be obtained based on these accompanying drawings without paying creative work.
[0084] Figure 1 It is a schematic diagram of the method for processing seismic data while drilling involved in the present invention.
[0085] Figure 2 It is a schematic diagram of the original data of the drilling data involved in the present invention.
[0086] Figure 3 This is a schematic diagram of phase-consistent data involved in the present invention.
[0087] Figure 4 It is a schematic diagram of superimposed data of phase-consistent data involved in the present invention.
[0088] Figure 5 It is a schematic diagram of the filtering signal involved in the present invention.
[0089] Figure 6 It is a schematic diagram of a primary differential signal involved in the present invention.
[0090] Figure 7 It is a schematic diagram of the processed while drilling data (ie, the signal after secondary difference) involved in the present invention.
[0091] Figure 8 It is a schematic diagram comparing the core column diagram of the drilling involved in the present invention and the processed while-drilling data. DETAILED DESCRIPTION
[0092] In the following, only some exemplary embodiments are briefly described. As those skilled in the art will appreciate, the described embodiments may be modified in various ways without departing from the spirit or scope of the embodiments of the present invention. Therefore, the drawings and descriptions are considered to be exemplary and non-restrictive in nature.
[0093] The disclosure below provides many different embodiments or examples to implement different structures of the embodiments of the present invention. In order to simplify the disclosure of the embodiments of the present invention, the parts and settings of specific examples are described below. Of course, they are only examples, and the purpose is not to limit the embodiments of the present invention. In addition, the embodiments of the present invention can repeat reference numbers and / or reference letters in different examples, and this repetition is for the purpose of simplification and clarity, which itself does not indicate the relationship between the various embodiments and / or settings discussed.
[0094] The embodiments of the present invention are described in detail below with reference to the accompanying drawings.
[0095] like Figure 1 As shown, this embodiment provides a method for processing seismic data while drilling, comprising:
[0096] Arrange several geophones in a straight line at the wellhead to collect data while drilling;
[0097] Based on the drilling data, the seismic data signal is extracted and segmented according to its period to obtain m groups of segmented signals;
[0098] Perform Fourier transform on the first group and the mth group of split signals, convert from the time domain to the frequency domain, and perform conjugate complex multiplication on the frequency domains corresponding to the first group and the mth group of split signals to obtain the frequency domain correlation coefficient;
[0099] Perform inverse Fourier transform on the frequency domain correlation coefficient, convert the frequency domain correlation coefficient to the time domain, and obtain the time domain correlation coefficient;
[0100] Based on the time domain correlation coefficient, the largest group of correlation coefficients is extracted to obtain the phase difference of the segmented signal;
[0101] Based on the phase difference of the split signals, the phase of the mth group of split signals is adjusted to be consistent with the phase of the first group of split signals;
[0102] The mth group of split signals after phase adjustment is subjected to superposition processing to eliminate the influence of random noise on the signals and obtain the superpositioned signals;
[0103] A Blackman window is formed according to the filter, and the FIR filter coefficients are designed with the Blackman window so that the weighted square sum of the error between the actual filter response of the filter and the specified ideal frequency response is minimized, and the optimal sampling point weight is obtained;
[0104] Multiply the weight of the optimal sampling point by the Blackman window to obtain the optimal filter;
[0105] Using an optimal filter to filter the superimposed signal to obtain a filtered signal;
[0106] Take a differential of the filtered signal, calculate the change of the signal, and obtain a differential signal;
[0107] The secondary difference of the primary differential signal is obtained to obtain its variation, and finally the processed while drilling data is obtained.
[0108] In some embodiments, the expression of the while drilling data is formula (1);
[0109] S = S0 + N (31);
[0110] In the formula, S is the acquired seismic data, S0 is the acquired effective seismic record, and N is the noise generated during the acquisition process;
[0111] Take out a 2s seismic data signal and divide it according to its period. The expression is formula (2);
[0112] x m (t) = S (T·(m-1) + t) (32);
[0113] In the formula, m represents the group number of the segmentation, x m (t) represents the mth group of signals after segmentation, T represents the signal period, and t represents the number of sampling points of the signal.
[0114] In some embodiments, the first group and the mth group of segmented signals are subjected to Fourier transform, and converted from the time domain to the frequency domain, and the expressions thereof are formula (3) and formula (4);
[0115]
[0116] Where X1(ω) represents the frequency domain of x1(t), X m (ω) represents x m (t) is the frequency domain, ω represents the amplitude value in the frequency domain, i is the imaginary unit, e is the natural constant, t represents the number of sampling points of the signal, It is the accumulation in the time domain.
[0117] In some embodiments, X1(ω) and X m The frequency domain correlation coefficient is obtained by multiplying the conjugate complex number of (ω), and its expression is formula (5);
[0118]
[0119] Where C(ω) represents the frequency domain correlation coefficient, For X m The complex conjugate of (ω);
[0120] Perform inverse Fourier transform on C(ω) to convert the frequency domain correlation coefficient into the time domain, which is expressed as formula (6);
[0121]
[0122] Where c(t) represents the time domain correlation coefficient.
[0123] In some embodiments, the largest group of correlation coefficients calculated in formula (6) is extracted to obtain the phase difference value of the split signal, which is expressed as formula (7);
[0124] (r,p)=max(c) (37);
[0125] In the formula, r represents the maximum value of c, p represents the value of t corresponding to the maximum value of c, max() is the function for finding the maximum value, and c(t) represents the time domain correlation coefficient;
[0126] Adjust x according to the calculation result of formula (7) m The phase of the signal is adjusted to be consistent with the phase of x1, and its expression is formula (8);
[0127]
[0128] In the formula, represents x after phase adjustment m (t).
[0129] In some embodiments, the signal obtained in formula (8) Perform superposition processing to eliminate the influence of random noise on the signal, and its expression is formula (9);
[0130]
[0131] In the formula, x(t) represents the signal after superposition, Represents x m (t) Data after phase adjustment, n is the total number of groups after data segmentation.
[0132] In some embodiments, a Blackman window is formed according to the filter, and its expression is formula (10);
[0133]
[0134] Where k represents the number of filter points, u represents the filter length, and w(k) represents the Blackman window function;
[0135] The FIR filter coefficients are designed using the generated Blackman window w(k) so that the weighted sum of squares of the errors between the actual filter response and the specified ideal frequency response is minimized. The objective function is expressed as formula (11);
[0136]
[0137] In the formula, E represents the error sum of squares, h(k) represents the weight of the kth sampling point, and H(e iωk ) represents the actual frequency response at frequency ω, D k represents the ideal amplitude response at frequency ω;
[0138] Multiply the optimal h(k) obtained from formula (11) by the Blackman window w(k) to obtain the optimal filter, which is expressed as formula (12);
[0139]
[0140] In the formula, represents the designed filter;
[0141] The superimposed signal is filtered according to the Blackman window in formula (12), and its expression is formula (13);
[0142]
[0143] In the formula, represents the filtered superimposed signal at time t, that is, the filtered signal.
[0144] In some embodiments, a difference is taken for the signal in formula (13) to obtain the change in the signal, which is expressed as formula (14);
[0145]
[0146] In the formula, d(t) represents the first difference of the filtered signal, that is, the first difference signal; represents the filtered superposition signal at time t; Represents the filtered superposition signal at time t+1.
[0147] The signal in formula (14) is taken to obtain the second difference, and its variation is obtained, and finally the processed while drilling data is obtained, and its expression is formula (15);
[0148] dd(t)=|d(t+1)|-|d(t)| (45);
[0149] Wherein, dd(t) represents the secondary difference of the filtered signal, i.e. the processed while drilling data; | | is the absolute value symbol; d(t) represents the primary difference of the filtered signal at time t; d(t+1) represents the primary difference of the filtered signal at time t+1;
[0150] The while-drilling signal calculated in formula (15), i.e. the processed while-drilling data, is compared and analyzed with the core histogram.
[0151] The embodiment of this specification also discloses a drilling-while-drilling seismic data processing system, which is used to implement any one of the drilling-while-drilling seismic data processing methods described above, and the drilling-while-drilling seismic data processing system includes:
[0152] Arrangement acquisition module, used to arrange a number of detectors in a straight line at the wellhead to collect drilling data while drilling;
[0153] A segmentation module is used to extract seismic data signals based on the drilling data, segment them according to their periods, and obtain m groups of segmented signals;
[0154] Transformation module, used to:
[0155] Perform Fourier transform on the first group and the mth group of split signals, convert from the time domain to the frequency domain, and perform conjugate complex multiplication on the frequency domains corresponding to the first group and the mth group of split signals to obtain the frequency domain correlation coefficient;
[0156] Perform inverse Fourier transform on the frequency domain correlation coefficient, convert the frequency domain correlation coefficient to the time domain, and obtain the time domain correlation coefficient;
[0157] Based on the time domain correlation coefficient, the largest group of correlation coefficients is extracted to obtain the phase difference of the segmented signal;
[0158] Phase and superposition modules for:
[0159] Based on the phase difference of the split signals, the phase of the mth group of split signals is adjusted to be consistent with the phase of the first group of split signals;
[0160] The mth group of split signals after phase adjustment is subjected to superposition processing to eliminate the influence of random noise on the signals and obtain the superpositioned signals;
[0161] Filter module for:
[0162] A Blackman window is formed according to the filter, and the FIR filter coefficients are designed with the Blackman window so that the weighted square sum of the error between the actual filter response of the filter and the specified ideal frequency response is minimized, and the optimal sampling point weight is obtained;
[0163] Multiply the weight of the optimal sampling point by the Blackman window to obtain the optimal filter;
[0164] Using an optimal filter to filter the superimposed signal to obtain a filtered signal;
[0165] Differential modules for:
[0166] Take a differential of the filtered signal, calculate the change of the signal, and obtain a differential signal;
[0167] The secondary difference of the primary differential signal is obtained to obtain its variation, and finally the processed while drilling data is obtained.
[0168] The embodiments of this specification also disclose a computer-readable storage medium, wherein the storage medium stores computer instructions. When a computer reads the computer instructions, the computer executes any one of the above-mentioned methods for processing seismic data while drilling.
[0169] In order to further illustrate the technical solution of the present invention, the technical concept of the present invention is as follows:
[0170] Seismic exploration is widely used in the exploration and development of geological resources, but the traditional method of using explosives cannot efficiently obtain high-quality seismic data. It is necessary to determine the depth and lithology of the explosives during the drilling process to optimize the amount of explosives. In order to find the appropriate location of explosives, improve the quality of seismic data acquisition, and reduce the cost of seismic acquisition, it is necessary to accurately determine the distribution of underground rock lithology.
[0171] The specific steps are as follows:
[0172] (1) Several geophones are arranged in a straight line at the wellhead to collect data while drilling. The expression is:
[0173] S=S0+N (46)
[0174] Where S is the acquired seismic data, S0 is the acquired effective seismic record, and N is the noise generated during the acquisition process.
[0175] Take out a 2s seismic data signal and divide it according to its period. The expression is:
[0176] x m (t) = S (T·(m-1) + t) (47)
[0177] In the formula, S represents the collected seismic data, m represents the group number, and x m (t) represents the mth group of signals after segmentation, T represents the signal period, and t represents the number of sampling points of the signal.
[0178] (2) For x1 and x m Perform Fourier transform from time domain to frequency domain, the expression is as follows:
[0179]
[0180] Where X1(ω) represents the frequency domain of x1(t), X m (ω) represents x m (t) is the frequency domain, t represents the number of sampling points of the signal, ω represents the amplitude value in the frequency domain, and i is the imaginary unit.
[0181] X1(ω) and X m Multiplying the conjugate complex number of (ω) gives the frequency domain correlation coefficient, which is expressed as follows:
[0182]
[0183] Where C(ω) represents the frequency domain correlation coefficient, X1(ω) represents the frequency domain of x1(t), and X m (ω) represents x m (t) in the frequency domain, For X m The conjugate complex number of (ω), ω represents the amplitude value in the frequency domain.
[0184] Perform inverse Fourier transform on C(ω) to convert the frequency domain correlation coefficient to the time domain. The expression is as follows:
[0185]
[0186] Where c(t) represents the time domain correlation coefficient, C(ω) represents the frequency domain correlation coefficient, t represents the number of sampling points of the signal, ω represents the amplitude value in the frequency domain, t represents the number of sampling points of the signal, and i is an imaginary unit.
[0187] (3) Extract the largest group of correlation coefficients calculated in formula (6) and calculate the phase difference of the split signal. The expression is:
[0188] (r,p)=max(c) (52)
[0189] In the formula, r represents the maximum value of y, p represents the value of t corresponding to the maximum value of y, and max() is the function for finding the maximum value.
[0190] Adjust x according to the calculation result of formula (7) m The phase of the signal is adjusted to match the phase of x1. The expression is:
[0191]
[0192] In the formula, represents x after phase adjustment m (t), x m (t) represents the mth group of signals after segmentation, p represents the value of t corresponding to the maximum value of y, and t represents the number of sampling points of the signal.
[0193] (4) The signal obtained in formula (8) Perform superposition processing to eliminate the influence of random noise on the signal. The expression is:
[0194]
[0195] In the formula, x(t) represents the signal after superposition, Represents x m (t) Data after phase adjustment, n is the total number of groups after data segmentation.
[0196] (5) A Blackman window is formed according to the filter, and its expression is:
[0197]
[0198] Where k represents the number of filter points, u represents the filter length, and w(k) represents the Blackman window function.
[0199] The FIR filter coefficients are designed with the generated Blackman window w(k) so that the weighted sum of squares of the error between the actual filter response and the specified ideal frequency response is minimized. The objective function is expressed as:
[0200]
[0201] In the formula, E represents the error sum of squares, e represents the filter length, h(k) represents the weight of the kth sampling point, ω represents the amplitude value in the frequency domain, and H(e iωk ) represents the actual frequency response at frequency ω, D k represents the ideal amplitude response at frequency ω, and k represents the number of filter points.
[0202] The optimal h(k) obtained from formula (11) is multiplied by the Blackman window w(k) to obtain the filter, which is expressed as:
[0203]
[0204] In the formula, represents the designed filter, h(k) represents the weight of the kth sampling point, w(k) represents the Blackman window function, and k represents the number of filter points.
[0205] The superimposed signal is filtered according to the Blackman window in formula (12), and its expression is:
[0206]
[0207] In the formula, represents the superimposed signal after filtering, represents the designed filter, t represents the number of sampling points of the signal, and k represents the number of filter points.
[0208] (6) Take a differential of the signal in formula (13) to obtain the change in the signal
[0209]
[0210] Where d(t) represents the first difference of the filtered signal, represents the superimposed signal after filtering, and t represents the number of sampling points of the signal.
[0211] (7) Take the second difference of the signal in formula (14) to find its change, and finally obtain the processed while drilling data.
[0212] dd(t)=|d(t+1)|-|d(t)| (60)
[0213] Where dd(t) represents the second difference of the filtered signal, d(t) represents the first difference of the filtered signal, || is the absolute value symbol, and t represents the number of sampling points of the signal.
[0214] (8) Compare and analyze the drilling signal calculated in formula (15) with the core column chart.
[0215] In a specific embodiment, during the drilling process of a well in an actual work area, instruments are deployed to collect data while drilling. The collected signals are as follows: Figure 2 As shown in the figure, there is a lot of irregular noise, but there is obvious periodicity in the data. The downhole data is cut into cycles to form multiple data files with the same cycle, and the phase-consistent data is obtained according to their mathematical relationship ( Figure 3 ). Figure 4 The superposition data of the phase-consistent data has shown the overall waveform characteristics, but there is still some residual random noise. The superimposed signal is filtered ( Figure 5 ). Take the difference of the filtered data to get the rate of change of the signal ( Figure 6 ). Further calculate the change of the signal to obtain the final processing result ( Figure 7 In order to verify the application effect of this technology, we compared the core column chart of the drilling with the processed near-surface while-drilling signal ( Figure 8 ), it can be seen that the processed LWD signal has obvious changes in the wave impedance interface, which corresponds well to the core column diagram, indicating that the processed LWD signal can better reflect the formation interface.
[0216] The above-described embodiments are used to illustrate the present invention, not to limit the present invention, so changes in the exemplified values or replacement of equivalent elements should still fall within the scope of the present invention.
[0217] From the above detailed description, it can be understood by those skilled in the art that the present invention can indeed achieve the aforementioned objectives and is in compliance with the provisions of the Patent Law.
[0218] Although the preferred embodiments of the present invention have been described, those skilled in the art may make additional changes and modifications to these embodiments once they are aware of the basic creative concepts. Therefore, the appended claims are intended to be interpreted as including the preferred embodiments and all changes and modifications that fall within the scope of the present invention. The above description is only a preferred embodiment of the present invention and is not intended to limit the present invention. It should be pointed out that any modifications, equivalent substitutions and improvements made within the spirit and principles of the present invention should be included in the scope of protection of the present invention.
[0219] It should be noted that the above description of the relevant process is only for example and explanation, and does not limit the scope of application of this specification. For those skilled in the art, various modifications and changes can be made to the process under the guidance of this specification. However, these modifications and changes are still within the scope of this specification.
[0220] The basic concepts have been described above. Obviously, for those of ordinary skill in the art who have read this application, the above invention disclosure is only for example and does not constitute a limitation of this application. Although not explicitly stated here, those of ordinary skill in the art may make various modifications, improvements and amendments to this application. Such modifications, improvements and amendments are suggested in this application, so such modifications, improvements and amendments still belong to the spirit and scope of the exemplary embodiments of this application.
[0221] At the same time, the present application uses specific words to describe the embodiments of the present application. For example, "one embodiment", "an embodiment", and / or "some embodiments" refer to a certain feature, structure or characteristic related to at least one embodiment of the present application. Therefore, it should be emphasized and noted that "one embodiment" or "an embodiment" or "an alternative embodiment" mentioned twice or more in different positions in this specification does not necessarily refer to the same embodiment. In addition, some features, structures or characteristics in one or more embodiments of the present application can be appropriately combined.
[0222] In addition, it will be appreciated by those skilled in the art that various aspects of the present application may be illustrated and described by a number of patentable categories or situations, including any new and useful combination of processes, machines, products or substances, or any new and useful improvements thereto. Therefore, various aspects of the present application may be implemented entirely by hardware, entirely by software (including firmware, resident software, microcode, etc.), or by a combination of hardware and software. The above hardware or software may all be referred to as "units", "modules" or "systems". In addition, various aspects of the present application may take the form of a computer program product embodied in one or more computer-readable media, wherein computer-readable program code is contained therein.
[0223] The computer program code required for the operation of each part of the application can be written in any one or more programming languages, including object-oriented programming languages such as Java, Scala, Smalltalk, Eiffel, JADE, Emerald, C++, C#, VB.NET, Python, etc., such as C programming language, VisualBasic, Fortran2103, Perl, COBOL2102, PHP, ABAP, such as Python, Ruby and Groovy dynamic programming language or other programming languages. The program code can be run completely on the user's computer, or run on the user's computer as an independent software package, or run partly on the user's computer and partly on the remote computer, or run completely on the remote computer or server. In the latter case, the remote computer can be connected to the user's computer through any network form, such as a local area network (LAN) or a wide area network (WAN), or connected to an external computer (e.g., via the Internet), or in a cloud computing environment, or used as a service such as software as a service (SaaS).
[0224] In addition, unless explicitly stated in the claims, the order of the processing elements and sequences described in this application, the use of alphanumeric characters, or the use of other names are not intended to limit the order of the processes and methods of this application. Although the above disclosure discusses some invention embodiments that are currently considered useful through various examples, it should be understood that such details are only for illustrative purposes, and the attached claims are not limited to the disclosed embodiments. On the contrary, the claims are intended to cover all modifications and equivalent combinations that are consistent with the essence and scope of the embodiments of this application. For example, although the implementation of the various components described above can be embodied in a hardware device, it can also be implemented as a pure software solution, for example, installation on an existing server or mobile device.
[0225] Similarly, it should be noted that in order to simplify the description disclosed in this application and thus help understand one or more embodiments of the invention, in the above description of the embodiments of the application, multiple features are sometimes combined into one embodiment, figure or description thereof. However, this method of the application should not be interpreted as reflecting the intention that the claimed object requires more features than those explicitly stated in each claim. On the contrary, the subject of the invention should have fewer features than the above single embodiment.
Claims
1. A method for processing seismic data while drilling, characterized in that: include: Arrange several geophones in a straight line at the wellhead to collect data while drilling; Based on the drilling data, the seismic data signal is extracted and segmented according to its period to obtain m groups of segmented signals; Perform Fourier transform on the first group and the mth group of split signals, convert from the time domain to the frequency domain, and perform conjugate complex multiplication on the frequency domains corresponding to the first group and the mth group of split signals to obtain the frequency domain correlation coefficient; Perform inverse Fourier transform on the frequency domain correlation coefficient, convert the frequency domain correlation coefficient to the time domain, and obtain the time domain correlation coefficient; Based on the time domain correlation coefficient, the largest group of correlation coefficients is extracted to obtain the phase difference of the segmented signal; Based on the phase difference of the split signals, the phase of the mth group of split signals is adjusted to be consistent with the phase of the first group of split signals; The mth group of split signals after phase adjustment is subjected to superposition processing to eliminate the influence of random noise on the signals and obtain the superpositioned signals; A Blackman window is formed according to the filter, and the FIR filter coefficients are designed with the Blackman window so that the weighted square sum of the error between the actual filter response of the filter and the specified ideal frequency response is minimized, and the optimal sampling point weight is obtained; Multiply the weight of the optimal sampling point by the Blackman window to obtain the optimal filter; Using an optimal filter to filter the superimposed signal to obtain a filtered signal; Take a differential of the filtered signal, calculate the change of the signal, and obtain a differential signal; The secondary difference of the primary differential signal is obtained to obtain its variation, and finally the processed while drilling data is obtained.
2. The method for processing seismic data while drilling according to claim 1, characterized in that: The expression of the drilling data is formula (1); S = S0 + N(1); In the formula, S is the acquired seismic data, S0 is the acquired effective seismic record, and N is the noise generated during the acquisition process; Take out a 2s seismic data signal and divide it according to its period. The expression is formula (2); x m (t)=S(T·(m-1)+t)(2); In the formula, m represents the group number of the segmentation, x m (t) represents the mth group of signals after segmentation, T represents the signal period, and t represents the number of sampling points of the signal.
3. The method for processing seismic data while drilling according to claim 2, characterized in that: Perform Fourier transform on the first group and the mth group of segmented signals, converting them from the time domain to the frequency domain, and the expressions are formula (3) and formula (4); Where X1(ω) represents the frequency domain of x1(t), X m (ω) represents x m (t) is the frequency domain, ω represents the amplitude value in the frequency domain, i is the imaginary unit, e is the natural constant, t represents the number of sampling points of the signal, It is the accumulation in the time domain.
4. The method for processing seismic data while drilling according to claim 3, characterized in that: X1(ω) and X m The frequency domain correlation coefficient is obtained by multiplying the conjugate complex number of (ω), and its expression is formula (5); Where C(ω) represents the frequency domain correlation coefficient, For X m The complex conjugate of (ω); Perform inverse Fourier transform on C(ω) to convert the frequency domain correlation coefficient into the time domain, which is expressed as formula (6); Where c(t) represents the time domain correlation coefficient.
5. The method for processing seismic data while drilling according to claim 4, characterized in that: The largest group of correlation coefficients calculated in formula (6) is extracted to obtain the phase difference of the split signal, which is expressed as formula (7); (r,p)=max(c)(7); In the formula, r represents the maximum value of c, p represents the value of t corresponding to the maximum value of c, max() is the function for finding the maximum value, and c(t) represents the time domain correlation coefficient; Adjust x according to the calculation result of formula (7) m The phase of the signal is adjusted to be consistent with the phase of x1, and its expression is formula (8); In the formula, represents x after phase adjustment m (t).
6. The method for processing seismic data while drilling according to claim 5, characterized in that: For the signal obtained in formula (8) Perform superposition processing to eliminate the influence of random noise on the signal, and its expression is formula (9); In the formula, x(t) represents the signal after superposition, Represents x m (t) Data after phase adjustment, n is the total number of groups after data segmentation.
7. The method for processing seismic data while drilling according to claim 6, characterized in that: A Blackman window is formed according to the filter, and its expression is formula (10); Where k represents the number of filter points, u represents the filter length, and w(k) represents the Blackman window function; The FIR filter coefficients are designed using the generated Blackman window w(k) so that the weighted sum of squares of the errors between the actual filter response and the specified ideal frequency response is minimized. The objective function is expressed as formula (11); In the formula, E represents the error sum of squares, h(k) represents the weight of the kth sampling point, and H(e iωk ) represents the actual frequency response at frequency ω, D k represents the ideal amplitude response at frequency ω; Multiply the optimal h(k) obtained from formula (11) by the Blackman window w(k) to obtain the optimal filter, which is expressed as formula (12); In the formula, represents the designed filter; The superimposed signal is filtered according to the Blackman window in formula (12), and its expression is formula (13); In the formula, Represents the filtered superimposed signal at time t, that is, the filtered signal.
8. The method for processing seismic data while drilling according to claim 7, characterized in that: Take a differential of the signal in formula (13) to obtain the change in the signal, which is expressed as formula (14); In the formula, d(t) represents the first difference of the filtered signal, that is, the first difference signal; represents the filtered superposition signal at time t, Represents the filtered superposition signal at time t+1. The signal in formula (14) is taken to obtain the second difference, and its variation is obtained, and finally the processed while drilling data is obtained, and its expression is formula (15); dd(t)=|d(t+1)|-|d(t)| (15); Wherein, dd(t) represents the second difference of the filtered signal at time t, i.e. the processed while drilling data; | | is the absolute value symbol; d(t) represents the first difference of the filtered signal at time t; d(t+1) represents the first difference of the filtered signal at time t+1; The while-drilling signal calculated in formula (15), i.e. the processed while-drilling data, is compared and analyzed with the core histogram.
9. A seismic data processing system while drilling, characterized in that: Used to implement the drilling-while-drilling seismic data processing method according to any one of claims 1 to 8, the drilling-while-drilling seismic data processing system comprises: Arrangement acquisition module, used to arrange a number of detectors in a straight line at the wellhead to collect drilling data while drilling; A segmentation module is used to extract seismic data signals based on the drilling data, segment them according to their periods, and obtain m groups of segmented signals; Transformation module, used to: Perform Fourier transform on the first group and the mth group of split signals, convert from the time domain to the frequency domain, and perform conjugate complex multiplication on the frequency domains corresponding to the first group and the mth group of split signals to obtain the frequency domain correlation coefficient; Perform inverse Fourier transform on the frequency domain correlation coefficient, convert the frequency domain correlation coefficient to the time domain, and obtain the time domain correlation coefficient; Based on the time domain correlation coefficient, the largest group of correlation coefficients is extracted to obtain the phase difference of the segmented signal; Phase and superposition modules for: Based on the phase difference of the split signals, the phase of the mth group of split signals is adjusted to be consistent with the phase of the first group of split signals; The mth group of split signals after phase adjustment is subjected to superposition processing to eliminate the influence of random noise on the signals and obtain the superpositioned signals; Filter module for: A Blackman window is formed according to the filter, and the FIR filter coefficients are designed with the Blackman window so that the weighted square sum of the error between the actual filter response of the filter and the specified ideal frequency response is minimized, and the optimal sampling point weight is obtained; Multiply the weight of the optimal sampling point by the Blackman window to obtain the optimal filter; Using an optimal filter to filter the superimposed signal to obtain a filtered signal; Differential modules for: Take a differential of the filtered signal, calculate the change of the signal, and obtain a differential signal; The secondary difference of the primary differential signal is obtained to obtain its variation, and finally the processed while drilling data is obtained.
10. A computer-readable storage medium, characterized in that: The storage medium stores computer instructions. When a computer reads the computer instructions, the computer executes the method for processing seismic data while drilling according to any one of claims 1 to 8.
Citation Information
Patent Citations
Method for extracting time difference in real time for logging-while-drilling acoustic underground signal acquisition and processing system
CN103742131A
Process and device for detecting seismic signals in order to obtain vertical seismic profiles during bore drilling operations
CN1101428A
Earthquake guiding method and system for rock drilling
CN111722276A
Tunnel vibroseis seismic data time-frequency detection and filtering method
CN117647835A
Method of obtaining seismic while drilling signal
US20230161060A1