A method for detecting and correcting time drift of node seismic data

Abstract

The present invention relates to the technical field of seismic data processing, and in particular to a method for detecting and correcting the time drift of node seismic data. The method comprises the following steps: sorting the seismic data collected by the node into a common shot point track set; calculating the theoretical first arrival time; opening a time window near the theoretical first arrival time; calculating the cross-correlation time difference of adjacent tracks within the time window; applying the cross-correlation time difference and calculating the cross-correlation coefficient of adjacent tracks; setting a threshold to filter out seismic tracks with time drift; and applying the cross-correlation time difference to correct the time drift seismic tracks. The method of the present invention can quickly and accurately detect node seismic data with time drift, and provides an effective technical means for realizing efficient and accurate node acquisition quality monitoring; at the same time, the method of the present invention can also automatically detect and correct the time drift of node seismic data, and has a wide range of application prospects in node acquisition projects.

Description

Technical Field

[0001] The present invention relates to the technical field of seismic data processing, and in particular to a method for detecting and correcting time drift of node seismic data. Background Art

[0002] When using cable-free node instruments for seismic exploration, the precise timing function of the node instruments relies on GPS timing. When the node instruments are flooded or blocked by thick clouds, the node instruments will not receive GPS signals and can only use internal clock timing. If the GPS signal is lost for a long time, the internal clock of the node instrument cannot be calibrated, resulting in a large time accumulation error, which eventually causes time drift in the recorded seismic data. Seismic data with large time drift will seriously affect data quality and subsequent processing. The conventional processing method for seismic data with large time drift is to delete it as a waste channel to prevent it from affecting normal data.

[0003] Chinese patent application CN111983684A discloses a time calibration method for four-dimensional seismic data, which includes: collecting the starting seismic data at the time of starting exploration and a monitoring seismic data group at a set time, wherein the monitoring seismic data group includes time-shifted seismic data corresponding to a number of time-shifted time points with different set time intervals from the set time; cross-correlating each time-shifted seismic data in the monitoring seismic data group with the starting seismic data, generating an average cross-correlation function with the time-shifted time point as the independent variable; determining the time drift of the monitoring seismic data group according to the maximum position of the average cross-correlation function; and time-calibrating the monitoring seismic data group according to the time drift to generate time-calibrated four-dimensional seismic data. This patent is for time calibration of seismic data collected multiple times (2 times or more) in different periods, and overall time shift correction of seismic data collected once, and its time drift is obtained by cross-correlating seismic data of different periods of the seismic channel of the same detection point. The time calibration of node seismic data is for seismic data collected at a single time in the same period. Time drift and correction only occur on some data channels, not on the entire data. The time drift of node seismic data can only be calculated through seismic channel data at different detection points. Therefore, an effective method for detecting and correcting the time drift of node seismic data is currently needed. Summary of the invention

[0004] The present invention provides a method for detecting and correcting the time drift of node seismic data. The method calculates the cross-correlation of multiple adjacent traces in a common shot point gather, detects seismic traces with time drift through the cross-correlation time difference and the cross-correlation coefficient, and then uses the cross-correlation time difference to perform time shift correction.

[0005] To achieve the above object, the present invention adopts the following technical solutions:

[0006] The present invention provides a method for detecting and correcting the time drift of node seismic data, which comprises the following steps: sorting the seismic data collected by the node into common shot point trace sets; calculating the theoretical first arrival time; opening a time window near the theoretical first arrival time; calculating the cross-correlation time difference of adjacent traces within the time window; applying the cross-correlation time difference and calculating the cross-correlation coefficient of adjacent traces; setting a threshold to screen out seismic traces with time drift; and applying the cross-correlation time difference to correct the seismic traces with time drift.

[0007] Furthermore, the raw data collected by the node are cut according to the excitation time of each shot to eliminate invalid data, then merged into normal channel data, and then sorted into common shot point data sets according to the shot point and detection point pile numbers, and the common shot point data sets are sorted by channel number.

[0008] Furthermore, when the near-surface structure is simple, the method for calculating the theoretical first arrival time is:

[0009] The theoretical first arrival time is calculated based on the offset using the following formula:

[0010] t=t0+d / v

[0011] Where t is the theoretical first arrival time, t0 is the zero offset time from the wellhead, d is the offset, and v is the refracted wave velocity;

[0012] The offset is calculated by the coordinates of the shot check point. The calculation formula of the offset is:

[0013]

[0014] Where d is the offset, (x s ,y s ) is the shot point coordinate, (x d ,y d ) are the coordinates of the detection point; or,

[0015] The theoretical first arrival time is calculated by fitting and interpolating the offset distance generated by manually picking the first arrival with the first arrival time table.

[0016] When the near-surface structure is complex, the near-surface velocity model is used to forward calculate the theoretical first arrival time.

[0017] Furthermore, a time window is opened near the theoretical first arrival time to automatically pick up the first arrival. The time window takes the theoretical first arrival time as the midpoint and extends upward and downward. The extension length of the time window should depend on the drift time of the seismic trace, and the time window length should be greater than the drift time of the seismic trace.

[0018] Further, calculate the cross-correlation function between the channel to be detected and multiple adjacent channels within a time window. The time corresponding to the maximum value of the cross-correlation function is the cross-correlation time difference. Calculate the mean value of the cross-correlation time differences as the final cross-correlation time difference;

[0019] The calculation formula for the cross-correlation function is:

[0020]

[0021] In the above formula, C is the cross-correlation function, t is the time shift amount, i is the time index of the samples within the time window, n is the number of samples within the time window, S is the seismic channel to be detected, and T is the adjacent seismic channel.

[0022] Further, first perform time shift on the seismic channel to be detected, then calculate the cross-correlation coefficient between the channel to be detected and multiple adjacent channels within the time window, and calculate the mean value of the cross-correlation coefficients as the final cross-correlation coefficient;

[0023] The time shift time is the calculated cross-correlation time difference.

[0024] Even further, the calculation formula for the cross-correlation coefficient is:

[0025]

[0026] In the above formula, r is the cross-correlation coefficient, i is the time index of the samples within the time window, n is the number of samples within the time window, S is the seismic channel to be detected, T is the adjacent seismic channel, is the mean value of the samples of the seismic channel detected within the time window, is the mean value of the samples of the adjacent seismic channel within the time window.

[0027] Further, the adjacent channels mentioned above are: the adjacent channels of the Mth channel are the (M - N)th channel to the (M + N)th channel, where N is a positive integer satisfying 1 < N < M;

[0028] The value of N is determined according to the signal-to-noise ratio of the data. When the signal-to-noise ratio is relatively high, N takes a smaller value.

[0029] Further, set a threshold through the product of the cross-correlation time difference and the cross-correlation coefficient.

[0030] Compared with the prior art, the present invention has the following advantages:

[0031] The method of the present invention can quickly and accurately detect the nodal seismic data with time drift, providing an effective technical means for realizing efficient and accurate quality monitoring of nodal acquisition; at the same time, the method of the present invention can also automatically detect and correct the time drift of the nodal seismic data, and has a wide application prospect in nodal acquisition projects. BRIEF DESCRIPTION OF THE DRAWINGS

[0032] The accompanying drawings in the specification, which constitute a part of the present invention, are used to provide a further understanding of the present invention. The exemplary embodiments of the present invention and their descriptions are used to explain the present invention and do not constitute improper limitations on the present invention.

[0033] Figure 1 It is a flow chart of a method for detecting and correcting time drift of node seismic data according to a specific embodiment of the present invention;

[0034] Figure 2 A common shot point gather collected by a node according to a specific embodiment of the present invention, wherein the seismic trace numbered 5773 has a time drift, as indicated by the arrow in the figure;

[0035] Figure 3 for Figure 2 The gathers in the figure and their theoretical first arrivals have time windows extending 500ms above and below: the white lines in the figure are the top and bottom of the time window;

[0036] Figure 4 A cross-plot of the channel number, the cross-correlation time difference and the product of the cross-correlation coefficient according to a specific embodiment of the present invention;

[0037] Figure 5 This is a gather after time drift correction according to a specific embodiment of the present invention. DETAILED DESCRIPTION

[0038] It should be noted that the following detailed descriptions are exemplary and are intended to provide further explanation of the present invention. Unless otherwise specified, all technical and scientific terms used herein have the same meanings as those commonly understood by those skilled in the art to which the present invention belongs.

[0039] It should be noted that the terms used herein are only for describing specific embodiments and are not intended to limit exemplary embodiments according to the present invention. As used herein, unless the context clearly indicates otherwise, the singular form is also intended to include the plural form. In addition, it should be understood that when the terms "comprising" and / or "including" are used in this specification, it indicates the presence of features, steps, operations and / or combinations thereof.

[0040] In order to enable those skilled in the art to more clearly understand the technical solution of the present invention, the technical solution of the present invention will be described in detail below in conjunction with specific embodiments.

[0041] Example 1

[0042] like Figure 1 As shown, the node seismic data time drift detection and correction method comprises the following steps:

[0043] 101) Sort the seismic data collected by the nodes into common shot point gathers:

[0044] The raw data collected by the node is cut according to the firing time of each shot to remove invalid data, and then merged into normal data. The data header contains information such as the shot checkpoint pile number, coordinates and elevation, and then sorted into a common shot point gather according to the shot point and detection point pile number, and the common shot point gather is sorted by the track number. Figure 2 shown.

[0045] 102) Calculate the theoretical first arrival time:

[0046] When the near-surface structure is simple, the method for calculating the theoretical first arrival time is:

[0047] The theoretical first arrival time is calculated based on the offset using the following formula:

[0048] t=t0+d / v

[0049] Where t is the theoretical first arrival time, t0 is the zero offset time from the wellhead, d is the offset, and v is the refracted wave velocity;

[0050] The offset is calculated by the coordinates of the shot check point. The calculation formula of the offset is:

[0051]

[0052] Where d is the offset, (x s ,y s ) is the shot point coordinate, (x d ,y d ) are the coordinates of the detection point; or,

[0053] The theoretical first arrival time is calculated by fitting and interpolating the offset distance generated by manually picking the first arrival with the first arrival time table.

[0054] When the near-surface structure is complex, the near-surface velocity model is used to forward calculate the theoretical first arrival time.

[0055] 103) Open a time window near the theoretical first arrival time:

[0056] A time window is opened near the theoretical first arrival time to automatically pick up the first arrival. The time window takes the theoretical first arrival time as the midpoint and extends upward and downward. The extension length of the time window should depend on the drift time of the seismic trace, and the time window length should be greater than the drift time of the seismic trace.

[0057] In this embodiment, the time window is taken as the midpoint of the theoretical first arrival time, extending upward by 500ms and downward by 500ms to open the time window. Figure 3 shown.

[0058] 104) Calculate the cross-correlation time difference between adjacent channels within the time window:

[0059] Calculate the cross-correlation function between the track to be detected and the adjacent first two tracks and the next two tracks in the time window. The time corresponding to the maximum value is the cross-correlation time difference. The mean of the cross-correlation time difference is taken as the final cross-correlation time difference.

[0060] The calculation formula of the cross-correlation function is:

[0061]

[0062] In the above formula, C is the cross-correlation function, t is the time shift, i is the time index of the sample point in the time window, n is the number of sample points in the time window, S is the seismic trace to be detected, and T is the adjacent seismic trace.

[0063] The time corresponding to the maximum value of the cross-correlation function is the cross-correlation time difference t c , the average of the cross-correlation time differences between the detected channel and the first two adjacent channels and the last two adjacent channels is

[0064] 105) Apply the cross-correlation time difference and calculate the cross-correlation coefficient of adjacent channels:

[0065] First, the seismic trace to be detected is time-shifted, and the time shift time is Then, the cross-correlation coefficients between the track to be detected and the first two tracks and the last two tracks adjacent to it are calculated within the time window, and the average of the cross-correlation coefficients is taken as the final cross-correlation coefficient.

[0066] The formula for calculating the mutual correlation coefficient is:

[0067]

[0068] In the above formula, r is the correlation coefficient, i is the time index of the sample point in the time window, n is the number of sample points in the time window, S is the seismic trace to be detected, T is the adjacent seismic trace, is the mean value of the seismic trace sample points detected in the time window, is the mean of adjacent seismic trace sample points in the time window.

[0069] The average of the correlation coefficients between the detected track and the first two tracks and the last two tracks is

[0070] 106) Set the threshold to filter out seismic traces with time drift:

[0071] Repeat steps 104 and 105 until all seismic traces are detected and the cross-correlation time difference of each trace is calculated. Mutual Correlation Coefficient The product of the channel number, the cross-correlation time difference and the product of the cross-correlation coefficient are plotted, such as Figure 4 As shown in Figure 2. According to the intersection diagram, set an appropriate threshold to filter out seismic traces with time drift. Figure 4It is shown that the seismic trace with the trace number 5773 in this embodiment has obvious anomalies, and the threshold range can be set as:

[0072] 107) Applying cross-correlation time difference to correct time drift seismic traces: according to the mean cross-correlation time difference calculated in step 104 The time shift correction is performed on the seismic channel with channel number 5773. Therefore, the data is shifted downward by 118ms; Figure 5 shown.

[0073] To sum up, the node seismic data time drift detection and correction method of the present invention can quickly and accurately detect node seismic data with time drift, providing an effective technical means for achieving efficient and accurate node acquisition quality monitoring; it can automatically detect and correct the time drift of node seismic data, and has broad application prospects in node acquisition projects.

[0074] The above embodiments are preferred implementation modes of the present invention, but the implementation modes of the present invention are not limited to the above embodiments. Any other changes, modifications, substitutions, combinations, and simplifications that do not deviate from the spirit and principles of the present invention should be equivalent replacement methods and are included in the protection scope of the present invention.

Claims

1. A method for detecting and correcting time drift of node seismic data, characterized in that: Including the following steps: Sort the seismic data collected by nodes into common shot gathers; Calculate the theoretical first arrival time; open a time window near the theoretical first arrival time; calculate the cross-correlation time difference between adjacent traces within the time window; apply the cross-correlation time difference and calculate the cross-correlation coefficient between adjacent traces; set a threshold to screen out seismic traces with time drift; apply the cross-correlation time difference to correct the time-drifted seismic traces.

2. The method according to claim 1, characterized in that: The original data collected by nodes is subjected to data cutting to remove invalid data according to the excitation time of each shot, then merged into normal trace data, and then sorted into common shot gathers according to the shot point and geophone stake numbers, and sorted within the common shot gather according to the trace number.

3. The method according to claim 1, characterized in that: When the near-surface structure is simple, the method for calculating the theoretical first arrival time is: Calculate the theoretical first arrival time according to the shot-receiver offset using the following formula: t = t0 + d / v Where t is the theoretical first arrival time, t0 is the zero-offset wellhead time, d is the shot-receiver offset, and v is the refraction wave velocity; Calculate the shot-receiver offset through the shot point and receiver coordinates, and the calculation formula for the shot-receiver offset is: Where d is the offset, (x s ,y s ) is the shot point coordinate, (x d ,y d ) are the coordinates of the detection point; or, Use the shot-receiver offset and first arrival time table generated by manual picking of the first arrival to perform fitting interpolation to calculate the theoretical first arrival time; When the near-surface structure is complex, use the near-surface velocity model to forward model and calculate the theoretical first arrival time.

4. The method according to claim 1, characterized in that: Automatically pick the first arrival by opening a time window near the theoretical first arrival time. The time window takes the theoretical first arrival time as the midpoint and extends upward and downward to open the time window. The extension length of the time window should depend on the drift time of the seismic trace, and the time window length should be greater than the drift time of the seismic trace.

5. The method according to claim 1, characterized in that: Calculate the cross-correlation function between the trace to be detected and multiple adjacent traces within the time window. The maximum value corresponds to the cross-correlation time difference, and find the mean value of the cross-correlation time differences as the final cross-correlation time difference; The calculation formula for the cross-correlation function is: In the above formula, C is the cross-correlation function, t is the time shift, i is the time index of the samples within the time window, n is the number of samples within the time window, S is the seismic trace to be detected, and T is the adjacent seismic trace.

6. The method according to claim 1, characterized in that: First shift the seismic trace to be detected in time, then calculate the cross-correlation coefficient between the trace to be detected and multiple adjacent traces within the time window, and find the mean value of the cross-correlation coefficients as the final cross-correlation coefficient; The time shift time is the calculated cross-correlation time difference.

7. The method according to claim 6, characterized in that: The calculation formula for the cross-correlation coefficient is: In the above formula, r is the correlation coefficient, i is the time index of the sample point in the time window, n is the number of sample points in the time window, S is the seismic trace to be detected, T is the adjacent seismic trace, is the mean value of the seismic trace sample points detected in the time window, is the mean of adjacent seismic trace sample points in the time window.

8. The method according to any one of claims 1, 5 and 6, characterized in that: The adjacent traces mentioned above are: the adjacent traces of the Mth trace are from the (M - N)th trace to the (M + N)th trace, where N is a positive integer satisfying 1 < N < M; The selection of N depends on the signal-to-noise ratio of the data. When the signal-to-noise ratio is relatively high, N takes a smaller value.

9. The method according to claim 1, characterized in that: Set a threshold through the product of the cross-correlation time difference and the cross-correlation coefficient.

Images