A method, system and storage medium for processing seismic data while drilling

Abstract

The present invention provides a method for processing seismic data while drilling, comprising: collecting data while drilling; segmenting the data while drilling to obtain m groups of segmented signals; performing Fourier transform on the segmented signals to obtain time domain correlation coefficients; obtaining phase differences of the segmented signals based on the time domain correlation coefficients; adjusting the phases of the segmented signals based on the phase differences of the segmented signals to make them consistent; superimposing the segmented signals after phase adjustment to obtain a superimposed signal; forming a Blackman window based on a filter, designing FIR filter coefficients to obtain an optimal filter; filtering the superimposed signal using the optimal filter to obtain a filtered signal; and obtaining first and second differences of the filtered signal to obtain processed data while drilling. The present invention can accurately describe lithologic changes in formations during drilling, and use this as a guide for seismic data acquisition.

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 primary source used in seismic exploration is an explosive source. The charge is closely related to factors such as the location of the explosives and the lithology of the well excitation point. The excitation energy and frequency are crucial to the signal-to-noise ratio and stratigraphic resolution of seismic waves. During conventional field production acquisition, the charge is typically determined based on the depth of the standard well depth (between 8 and 20 meters) to determine the data collection requirements for the entire work area. However, field excitation conditions vary widely, making it impossible to scientifically determine the charge (kilograms of explosives) and the stratum (especially for excitation of sandstone areas at appropriate depths). The explosive used for a single excitation is 4 to 15 kilograms of TNT. Due to the heterogeneity of the rock formations at each well location, the charge cannot be precisely determined, and placement can only be based on daily experience. Using too much explosive charge increases exploration costs and creates large fracture zones, weakening the elasticity of the formation and shifting the seismic signal band toward the lower end, resulting in loss of high-frequency signals. Using too little explosive charge results in a weak received signal, making effective, high-quality seismic exploration impossible. Due to the plasticity of mudstone and the elasticity of sandstone, even with the same explosive charge, the signal quality received when stimulating sandstone and mudstone differs significantly. This series of issues results in high data acquisition costs and low-quality seismic data. Therefore, the influence of formation lithology on explosive charge needs to be factored into explosive charge control during the acquisition process to provide reliable data support for further research on intelligent explosives.

[0003] Lithology identification is fundamental to studying geological reservoir characteristics, calculating reserves, and conducting geological modeling. It can determine the spatial distribution and range of different lithologies, providing detailed geological information. In practice, lithology identification is primarily performed through methods such as hand specimen identification, thin section analysis, element identification, and mineralogy. These methods can consistently obtain information about underground rock formations and accurately determine lithology, but they require extensive experience and knowledge, are inefficient, and require coring during drilling. The coring process can lead to rock fragmentation and disorganization, and once brought to the surface, its properties can undergo significant changes. Furthermore, extracting cores from every well is costly and challenging. In addition to the aforementioned lithology identification techniques, gravity and magnetic, seismic, and well logging can also assist in lithology identification. However, since the ranges of parameters such as gravity, magnetization, and elasticity for different lithologies often overlap, it can be difficult to clearly distinguish them using a single parameter. Therefore, high-resolution data is required to accurately distinguish lithologies. However, high-resolution data acquisition is expensive and difficult to cover all drilling areas.

[0004] During the drilling process, tens of thousands of wells are drilled. The vibrations from drilling generate seismic wave energy, which can be collected using high-precision instruments to identify formation lithology. However, the seismic wave signals collected in this way are large in volume, have high sampling rates, and are complex in composition (including valid seismic wave signals, drilling rig vibrations, electrical interference, and human interference). These signals are often submerged in the noise, making it difficult to extract valid signals from the drilling signals using conventional seismic data processing methods.

[0005] In view of this, there is an urgent need to provide a while-drilling seismic data processing method 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 while-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 and mth group of split signals, converting them from the time domain to the frequency domain, and perform conjugate complex multiplication on the frequency domains corresponding to the first and 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 split signal;

[0013] Based on the phase difference of the split signals, adjust the phase of the mth group of split signals to make it consistent with the phase of the first group of split signals;

[0014] The mth group of split signals after phase adjustment is superimposed to eliminate the influence of random noise on the signal and obtain the superimposed signal;

[0015] A Blackman window is formed according to the filter, and the FIR filter coefficients are designed using the Blackman window so that the weighted sum of squares of the errors between the actual filter response and the specified ideal frequency response of the filter is minimized, thereby obtaining the optimal sampling point weights;

[0016] Multiply the weight of the optimal sampling point by the Blackman window to obtain the optimal filter;

[0017] Filter the superimposed signal using an optimal filter 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] Where 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) 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 Multiplying the conjugate complex number of (ω) yields the frequency domain correlation coefficient, which is expressed as 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 to 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 of the split signal, which is expressed as formula (7);

[0036] (r,p)=max(c)(22);

[0037] Where 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] Where, 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] Where 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] Where E represents the sum of squared errors, 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] Where, 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] Where, Represents the filtered superimposed signal at time t, that is, the filtered signal.

[0056] In this specification, the signal in formula (13) is differentiated once to obtain the change in the signal, which is expressed as formula (14);

[0057]

[0058] Where d(t) represents the primary difference of the filtered signal at time t, that is, the primary 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, which is expressed as formula (15);

[0060] dd(t)=|d(t+1)||-|d(t)| (30);

[0061] Where dd(t) represents the second difference of the filtered signal, 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, and d(t+1) represents the first difference of the filtered signal at time t+1.

[0062] The LWD signal calculated in formula (15), i.e. the processed LWD data, is compared and analyzed with the core histogram.

[0063] This specification also discloses a drilling-while-drilling seismic data processing system for implementing any one of the above-mentioned drilling-while-drilling seismic data processing methods. 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 while-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 and mth group of split signals, converting them from the time domain to the frequency domain, and perform conjugate complex multiplication on the frequency domains corresponding to the first and 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 split signal;

[0070] Phase and Superposition Module for:

[0071] Based on the phase difference of the split signals, adjust the phase of the mth group of split signals to make it consistent with the phase of the first group of split signals;

[0072] The mth group of split signals after phase adjustment is superimposed to eliminate the influence of random noise on the signal and obtain the superimposed signal;

[0073] Filter module for:

[0074] A Blackman window is formed according to the filter, and the FIR filter coefficients are designed using the Blackman window so that the weighted sum of squares of the errors between the actual filter response and the specified ideal frequency response of the filter is minimized, thereby obtaining the optimal sampling point weights;

[0075] Multiply the weight of the optimal sampling point by the Blackman window to obtain the optimal filter;

[0076] Filter the superimposed signal using an optimal filter 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, which 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 achieve at least the following beneficial effects:

[0082] Taking advantage of the high drilling density in seismic acquisition production, during the drilling process, the energy generated by the drilling rig is used as the excitation source to optimize the ground observation method, and high-precision and high-density acquisition of drilling data is carried out to obtain seismic waves passing through different lithological strata. After processing the obtained seismic waves, the lithological interface is judged according to its morphological characteristics. By comparing the core histogram, this method can accurately describe the lithological changes of the strata during the drilling process and use this as a guide for seismic data acquisition. The amount of explosives used is reduced by 30%, and its 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 following briefly introduces the drawings required for use in the description of the embodiments. Obviously, the drawings described below are only some embodiments of the present invention. For ordinary technicians in this field, other drawings can be obtained based on these drawings without paying any creative work.

[0084] Figure 1 Schematic diagram of the method for processing seismic data while drilling involved in the present invention.

[0085] Figure 2 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-aligned data involved in the present invention.

[0087] Figure 4 Schematic diagram of superimposed data of phase-consistent data involved in the present invention.

[0088] Figure 5 Schematic diagram of the filtering signal involved in the present invention.

[0089] Figure 6 Schematic diagram of a primary differential signal involved in the present invention.

[0090] Figure 7 Schematic diagram of the processed while-drilling data (ie, the signal after secondary difference) involved in the present invention.

[0091] Figure 8 Schematic diagram comparing the core histogram of the well involved in the present invention and the processed while-drilling data. DETAILED DESCRIPTION

[0092] Hereinafter, only certain exemplary embodiments are briefly described. As will be appreciated by those skilled in the art, 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 description are to be regarded as illustrative in nature and not restrictive.

[0093] The disclosure below provides many different embodiments or examples for implementing different structures of the embodiments of the present invention. In order to simplify the disclosure of the embodiments of the present invention, the components and configurations of specific examples are described below. Of course, these are merely examples and are not intended to limit the embodiments of the present invention. In addition, the embodiments of the present invention may repeat reference numerals and / or reference letters in different examples. Such repetition is for the purpose of simplicity and clarity and does not in itself indicate the relationship between the various embodiments and / or configurations 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 while-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 and mth group of split signals, converting them from the time domain to the frequency domain, and perform conjugate complex multiplication on the frequency domains corresponding to the first and 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 split signal;

[0101] Based on the phase difference of the split signals, adjust the phase of the mth group of split signals to make it consistent with the phase of the first group of split signals;

[0102] The mth group of split signals after phase adjustment is superimposed to eliminate the influence of random noise on the signal and obtain the superimposed signal;

[0103] A Blackman window is formed according to the filter, and the FIR filter coefficients are designed using the Blackman window so that the weighted sum of squares of the errors between the actual filter response and the specified ideal frequency response of the filter is minimized, thereby obtaining the optimal sampling point weights;

[0104] Multiply the weight of the optimal sampling point by the Blackman window to obtain the optimal filter;

[0105] Filter the superimposed signal using an optimal filter 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 for the while drilling data is formula (1);

[0109] S=S0+N(31);

[0110] Where 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, which are expressed as formula (3) and formula (4);

[0115]

[0116] Where X1(ω) represents the frequency domain of x1(t), X m (ω) represents x m (t) 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 Multiplying the conjugate complex number of (ω) yields the frequency domain correlation coefficient, which is expressed as 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 to 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 of the split signal, which is expressed as formula (7);

[0124] (r,p)=max(c) (37);

[0125] Where 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] Where, 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] Where 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 based on the filter, which is expressed as 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] Where E represents the sum of squared errors, 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] Where, 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] Where, represents the filtered superposition 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] Where d(t) represents the first differential of the filtered signal, i.e., the first differential 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, which is expressed as formula (15);

[0148] dd(t)=|d(t+1)|-|d(t)| (45);

[0149] Where dd(t) represents the second difference of the filtered signal, 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;

[0150] The LWD signal calculated in formula (15), i.e. the processed LWD data, is compared and analyzed with the core histogram.

[0151] The embodiments of this specification further disclose a drilling-while-drilling seismic data processing system for implementing any one of the above-mentioned drilling-while-drilling seismic data processing methods. 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 while-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 and mth group of split signals, converting them from the time domain to the frequency domain, and perform conjugate complex multiplication on the frequency domains corresponding to the first and 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 split signal;

[0158] Phase and Superposition Module for:

[0159] Based on the phase difference of the split signals, adjust the phase of the mth group of split signals to make it consistent with the phase of the first group of split signals;

[0160] The mth group of split signals after phase adjustment is superimposed to eliminate the influence of random noise on the signal and obtain the superimposed signal;

[0161] Filter module for:

[0162] A Blackman window is formed according to the filter, and the FIR filter coefficients are designed using the Blackman window so that the weighted sum of squares of the errors between the actual filter response and the specified ideal frequency response of the filter is minimized, thereby obtaining the optimal sampling point weights;

[0163] Multiply the weight of the optimal sampling point by the Blackman window to obtain the optimal filter;

[0164] Filter the superimposed signal using an optimal filter 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 further 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. However, traditional methods using explosives cannot efficiently obtain high-quality seismic data. The depth and lithology of the explosives must be determined during the drilling process to optimize the explosive charge. To find the right explosives placement, improve seismic data acquisition quality, and reduce seismic acquisition costs, 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 a 2s seismic data signal and split 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, and its 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 complex conjugate of (ω), where ω 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. Its 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 of the frequency domain, t represents the number of sampling points of the signal, and i is the 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] Where 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] Where, 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] Where 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 based on 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 generated Blackman window w(k) is used to design the FIR filter coefficients 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:

[0200]

[0201] Where 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] Multiply the optimal h(k) obtained from formula (11) by the Blackman window w(k) to obtain the filter, which is expressed as:

[0203]

[0204] Where, Represents the designed filter, h(k) represents the weight of the k-th 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] Where, 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 variation, 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 histogram.

[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 drilling data is cut into cycles to form multiple data files with the same cycle. According to their mathematical relationship, the phase-consistent data is obtained ( 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 superposition 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 To verify the application effect of this technology, we compared the core histogram of the drilling well 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 histogram, indicating that the processed LWD signal can better reflect the formation interface.

[0216] The above embodiments are intended to illustrate the present invention, not to limit the present invention. Therefore, changes in illustrative values ​​or substitutions of equivalent components should still fall within the scope of the present invention.

[0217] From the above detailed description, it will be clear to 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 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 become aware of the basic inventive concepts. Therefore, the appended claims are intended to be interpreted as covering the preferred embodiments and all changes and modifications that fall within the scope of the invention. The foregoing description is merely a preferred embodiment of the present invention and is not intended to limit the present invention. It should be noted that any modifications, equivalent substitutions, and improvements made within the spirit and principles of the present invention are intended to be included within the scope of protection of the present invention.

[0219] It should be noted that the above description of the relevant processes is for illustration and purpose only 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 processes under the guidance of this specification. However, such modifications and changes are still within the scope of this specification.

[0220] The basic concepts have been described above. It will be apparent to those skilled in the art after reading this application that the above disclosures are merely illustrative and do not constitute limitations on this application. Although not explicitly stated herein, those skilled in the art may make various modifications, improvements, and amendments to this application. Such modifications, improvements, and amendments are suggested in this application and remain within the spirit and scope of the exemplary embodiments of this application.

[0221] At the same time, this application uses specific terms to describe the embodiments of this 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 this application. Therefore, it should be emphasized and noted that "one embodiment," "an embodiment," or "an alternative embodiment" mentioned twice or more in different places in this specification does not necessarily refer to the same embodiment. In addition, certain features, structures, or characteristics in one or more embodiments of this application may be appropriately combined.

[0222] In addition, it will be understood by those skilled in the art that various aspects of the present application can 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 can 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 can all be referred to as "units", "modules" or "systems". In addition, various aspects of the present application can 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, conventional procedural programming languages ​​such as C programming language, VisualBasic, Fortran2103, Perl, COBOL2102, PHP, ABAP, dynamic programming languages ​​such as Python, Ruby and Groovy 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 a remote computer, or run completely on a 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 be connected to an external computer (such as by the Internet), or in a cloud computing environment, or used as a service such as software as a service (SaaS).

[0224] In addition, unless expressly 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 embodiments of the invention 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 presentation of this disclosure and thereby facilitate understanding of one or more of the invention's embodiments, the foregoing descriptions of the embodiments of this disclosure sometimes combine multiple features into a single embodiment, figure, or description thereof. However, this approach should not be interpreted as reflecting an intention that the claimed subject matter requires more features than expressly recited in each claim. Rather, the subject matter of the invention may possess fewer features than the single embodiment described above.

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 while-drilling data, the seismic data signal is extracted and segmented according to its period to obtain m groups of segmented signals; Performing Fourier transform on the first and mth groups of split signals, converting from the time domain to the frequency domain, and performing conjugate complex multiplication on the frequency domains corresponding to the first and mth groups 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 split signal; Based on the phase difference of the split signals, adjust the phase of the mth group of split signals to make it consistent with the phase of the first group of split signals; The mth group of split signals after phase adjustment is superimposed to eliminate the influence of random noise on the signal and obtain the superimposed signal; A Blackman window is formed according to the filter, and the FIR filter coefficients are designed using the Blackman window so that the weighted sum of squares of the errors between the actual filter response and the specified ideal frequency response of the filter is minimized, thereby obtaining the optimal sampling point weights; Multiply the weight of the optimal sampling point by the Blackman window to obtain the optimal filter; Filter the superimposed signal using an optimal filter 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, wherein: The expression of the drilling data is formula (1); S=S0+N(1); Where 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, wherein: Perform Fourier transform on the first and mth group of segmented signals, converting them from the time domain to the frequency domain. The expressions are formula (3) and formula (4); Where X1(ω) represents the frequency domain of x1(t), X m (ω) represents x m (t) 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, wherein: X1(ω) and X m Multiplying the conjugate complex number of (ω) yields the frequency domain correlation coefficient, which is expressed as 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 to 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: Extract the largest group of correlation coefficients calculated in formula (6) and calculate the phase difference of the split signal, which is expressed as formula (7); (r,p)=max(c)(7); Where 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); Where, 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); Where 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, wherein: 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); Where E represents the sum of squared errors, 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); Where, represents the designed filter; The superimposed signal is filtered according to the Blackman window in formula (12), and its expression is formula (13); Where, 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, wherein: Take a differential of the signal in formula (13) to obtain the change in the signal, which is expressed as formula (14); Where d(t) represents the first differential of the filtered signal, i.e., the first differential 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, which is expressed as formula (15); dd(t)=|d(t+1)|-|d(t)| (15); Where 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 LWD signal calculated in formula (15), i.e. the processed LWD data, is compared and analyzed with the core histogram.

9. A seismic data processing system while drilling, characterized in that: For implementing the method for processing seismic data while drilling according to any one of claims 1 to 8, the seismic data processing system for processing seismic data while drilling 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 while-drilling data, segment them according to their periods, and obtain m groups of segmented signals; Transformation module, used to: Performing Fourier transform on the first and mth groups of split signals, converting from the time domain to the frequency domain, and performing conjugate complex multiplication on the frequency domains corresponding to the first and mth groups 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 split signal; Phase and Superposition Module for: Based on the phase difference of the split signals, adjust the phase of the mth group of split signals to make it consistent with the phase of the first group of split signals; The mth group of split signals after phase adjustment is superimposed to eliminate the influence of random noise on the signal and obtain the superimposed signal; Filter module for: A Blackman window is formed according to the filter, and the FIR filter coefficients are designed using the Blackman window so that the weighted sum of squares of the errors between the actual filter response and the specified ideal frequency response of the filter is minimized, thereby obtaining the optimal sampling point weights; Multiply the weight of the optimal sampling point by the Blackman window to obtain the optimal filter; Filter the superimposed signal using an optimal filter 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.

Images