Seismic section guided CPT data mixing field interpolation algorithm, device, medium and equipment

Abstract

The invention relates to a seismic section guided CPT data mixing field interpolation algorithm, device, medium and equipment, and the method comprises the steps: combining two-dimensional seismic data with hole site data to carry out the modeling and recognition of geological lithology, building a geological three-dimensional model to make a regional soil property parameter distribution frame, and carrying out the construction of a three-dimensional geologic model. And the purpose of predicting the soil property parameter distribution of the stratum at any point is achieved.

Description

Technical Field

[0001] The present invention relates to the field of engineering investigation and exploration technology, and in particular to a CPT data mixed field interpolation algorithm, device, medium and equipment guided by a seismic profile. Background Art

[0002] Marine geological survey technology uses survey tools to obtain data on the topography, stratum distribution, and rock and soil characteristics of the seabed. By integrating multi-source geological data, combining geological modeling with three-dimensional interpolation algorithms, key information such as seabed stratum structure and lithology distribution can be more accurately identified, providing a reliable basis for the site selection of offshore wind farms.

[0003] Prediction and interpolation of geological parameters is an important part of geological modeling. The seabed geological environment is complex, and various geophysical data are needed in the process of seabed geological interpolation. Geophysical data generally include seismic data and borehole data. Seismic data use strong and weak reflection axes to represent different underground lithologies, structures or stratigraphic interfaces. However, there are also great limitations in the identification of lithology by reflection axes, such as the difficulty in identifying the boundaries of lithology changes, and the difficulty in timely detection of lithology changes in the same layer by the naked eye.

[0004] CPT can determine the engineering mechanical properties of soil by measuring the stress in the stratum during the process of the probe penetrating the ground at a uniform speed, and can obtain high-precision, continuous and reliable data. In practical applications, although the borehole data and seismic data are usually combined for analysis, the degree of integration between the two is low, and generally only stays at the mutual reference in the interpretation of seismic data and geological stratification. Summary of the invention

[0005] In view of the above problems, the purpose of the present invention is to provide a seismic profile-guided CPT data mixed field interpolation algorithm, device, medium and equipment, which has the advantages of fusing seismic data and borehole data, ensuring the accuracy of stratum rock and soil parameter identification, and reflecting the real situation of seabed geology as much as possible, thereby providing reliable reference opinions for engineering.

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

[0007] The seismic profile-guided CPT data mixed domain interpolation algorithm of the present invention comprises the following steps:

[0008] Process the two-dimensional seismic images to obtain attribute features;

[0009] Based on the attribute characteristics, the spatial attributes and geological attributes of the points to be interpolated are calculated;

[0010] Based on the spatial attributes and geological attributes, CPT attribute parameters are extracted, and the CPT attribute parameter data are combined and merged with the features;

[0011] The combined and merged data are sliced ​​horizontally, and the slices are processed along the direction indicated by the two-dimensional seismic data to obtain the horizontal interpolation points and interpolation parameters to form multiple slice data sets; for each slice data set, it is processed according to the acquired two-dimensional seismic data;

[0012] In the transverse slice data set, the interpolation model is constructed using the 2D seismic data as the interpolation parameter and the cone tip resistance obtained from the CPT hole position as the interpolation result;

[0013] On the basis of the interpolation model, the nearest neighbor interpolation method and the inverse distance interpolation method are used for mixed field interpolation. The nearest neighbor interpolation and the inverse distance weighted interpolation are combined to fully consider the influence of spatial attributes and geological attributes to achieve smooth and accurate interpolation. A three-dimensional model is then established to form the final interpolation result.

[0014] The CPT data mixed domain interpolation algorithm preferably processes the two-dimensional seismic image to obtain attribute features, specifically comprising the following steps:

[0015] The two-dimensional seismic image is in the segy format file, which is a standard seismic data format used to store and exchange seismic reflection data. Seismic signal data is usually expressed as a time series, and the value of each sample point is the amplitude, that is:

[0016] Amplitude(t)=x(t) (1)

[0017] Where Amplitude(t) refers to the amplitude of the seismic signal at time t; x(t) is the function of the seismic signal over time, which represents the recorded seismic wave signal;

[0018] Compute phase and envelope characteristics:

[0019] The calculation of phase and envelope usually uses Hilbert transform, which is a tool in signal processing. Hilbert transform converts a real-valued signal into a complex-valued signal. The real part of the complex-valued signal is the original signal, and the imaginary part is its Hilbert transform. The formula of Hilbert transform is as follows:

[0020]

[0021] Where H(x(t)) is the signal after Hilbert transformation in time; x(τ) is the seismic signal value at time τ; t is the current time point; τ is the integral variable, representing all points on the time axis; PV is the principal value integral;

[0022] Through the Hilbert transform, we can get the complex valued signal, specifically:

[0023] z(t)=x(t)+jH(x(t)) (3)

[0024] Where z(t) is the complex value signal at time t; x(t) is the original seismic signal value at time t; j is the imaginary unit;

[0025] The instantaneous phase is the phase angle of a complex-valued signal and is expressed as:

[0026]

[0027] Where φ(t) is the instantaneous phase at time t; Im(z(t)) is the result of Hilbert transform, which is also the imaginary part of the complex-valued signal at time t; Re(z(t)) is the original seismic signal, which is also the real part of the complex-valued signal at time t;

[0028] The envelope is the modulus of the replicated signal and is calculated as:

[0029]

[0030] Where E(t) is the envelope amplitude at time t; |z(t)| is the modulus of the complex-valued signal at time t.

[0031] The CPT data mixed domain interpolation algorithm preferably requires the spatial attributes and geological attributes of the interpolation points to be calculated, and includes the following steps:

[0032] When calculating the depth of a point in a layer, the depth needs to be normalized. The formula for calculating the normalization is:

[0033]

[0034] In the formula, d is the original depth of the point; min is the minimum thickness of the formation, that is, the depth of the top layer of the formation; max is the maximum thickness of the formation, that is, the depth of the bottom layer of the formation; d NORM The value of the normalized depth of the point.

[0035] The CPT data mixed domain interpolation algorithm preferably extracts CPT attribute parameters and combines and merges CPT attribute parameter data with features, specifically including the following steps:

[0036] The data in the area where the CPT hole is located is selected as the training data, and the data obtained from the CPT is used as the test data. The parameters of the cone tip resistance are obtained to model and analyze the geological formation. c The calculation formula is:

[0037]

[0038] In the formula, F c is the cone tip penetration force; A c is the cross-sectional area of ​​the cone tip, calculated using a standard conical probe, using the formula:

[0039]

[0040] Where D c is the diameter of the cone tip;

[0041] Obtain the 2D seismic characteristics of the CPT hole location, including the amplitude, phase and envelope, merge and combine with the spatial attributes and geological attributes of the point, and use it together with the cone tip resistance data as a training data set;

[0042] The two-dimensional seismic characteristics of the points to be interpolated, including amplitude, phase and envelope, are merged and combined with the spatial attributes and geological attributes of the points, and the overall interpolation results are obtained after the interpolation model is trained and input.

[0043] The CPT data mixed field interpolation algorithm preferably comprises the following steps: in the transverse slice data set, the two-dimensional seismic data is used as the interpolation parameter, the cone tip resistance obtained at the CPT hole position is used as the interpolation result, and an interpolation model is constructed.

[0044] The cone tip resistance data set for the sample points is established as follows:

[0045]

[0046] In the formula, A collection of cone tip resistance data sets for sample points; f1 represents the cone tip resistance of the first sample point. K represents the cone tip resistance at the Kth sample point;

[0047] Its spatial attributes and two-dimensional seismic attributes are established as follows:

[0048]

[0049] Among them, x1 represents the amplitude, phase, envelope and spatial attributes of a CPT hole position;

[0050] Solve the diagonal equation, which is as follows:

[0051]

[0052] Where D(x) represents the metric tensor field; t(x) refers to the distance from the interpolation point to the nearest sample point x. k the time required; Represents points other than sample points, that is, points to be interpolated; Represents the gradient of the function t(x), that is, the spatial rate of change of t(x) at point x;

[0053] The calculation of the metric tensor field D(x) provides the anisotropy and spatial variation coefficient of the equation. By using the metric tensor field D(x) to measure the characteristic distance, it reflects the influence of geological structure on the characteristic distance, making two points in the same geological structure very close, while two points in different structures are much farther away;

[0054] Using spatial attributes to calculate t(x), the calculation formula of t(x) is as follows:

[0055]

[0056] In the formula, x mainly contains the spatial attributes of the interpolation points and sample points, that is, the spatial position; x k is the spatial attribute of the kth sample point; the formula for calculating t(x) for the sample point is as follows:

[0057]

[0058] The nearest neighbor interpolation method is used for the data. For each point to be interpolated, the nearest known data point is found and the value of the known data point is assigned to the point to be interpolated. For the known sample points, the cone tip resistance is used as the characteristic value;

[0059] The inverse distance interpolation method is used for the data. The inverse distance interpolation method makes the closest point have a greater impact on the interpolation result, and the farthest point has a smaller impact on the interpolation result. The calculation formula of the inverse distance interpolation method is:

[0060]

[0061] In the formula, q c (x) is the interpolation value of the point to be interpolated; N is the number of sample points that need to be provided in the vicinity.

[0062] The CPT data mixed domain interpolation algorithm, preferably, the formula of the mixed neighborhood interpolation is as follows:

[0063] q c (P) = w NN q c (Q nearest )+w IDW q c (IDW) (15)

[0064] In the formula, w NN is the weight of the nearest neighbor interpolation; q c (Q nearest ) is the result of nearest neighbor interpolation; w IDW is the weight of the inverse distance interpolation, qc (IDW) is the result of inverse distance interpolation; q c (P) is the final interpolation result.

[0065] The present invention also provides a seismic profile guided CPT data mixed field interpolation method, comprising:

[0066] A first processing unit is used to process the two-dimensional seismic image to obtain attribute features;

[0067] The second processing unit is used to calculate the spatial attributes and geological attributes of the points to be interpolated based on the attribute characteristics;

[0068] A third processing unit is used to extract CPT attribute parameters based on the spatial attributes and geological attributes, and to combine and merge the CPT attribute parameter data with the features;

[0069] The fourth processing unit is used to perform transverse slicing on the combined and merged data, and the slices are processed along the direction indicated by the two-dimensional seismic data to obtain transverse interpolation points and interpolation parameters to form multiple slice data sets; for each slice data set, processing is performed according to the acquired two-dimensional seismic data, including analyzing and extracting seismic features in the slice data;

[0070] The fifth processing unit is used to construct an interpolation model in the transverse slice data set using the two-dimensional seismic data as the interpolation parameter and the cone tip resistance obtained at the CPT hole position as the interpolation result;

[0071] The sixth processing unit is used to perform mixed field interpolation based on the interpolation model using the nearest neighbor interpolation method and the inverse distance interpolation method, combining the nearest neighbor interpolation and the inverse distance weighted interpolation, fully considering the influence of spatial attributes and geological attributes, and achieving smooth and accurate interpolation.

[0072] The present invention also provides a computer storage medium on which a computer program is stored. When the computer program is executed by a processor, the steps of the seismic profile-guided CPT data mixed field interpolation algorithm are implemented.

[0073] The present invention also provides a computer device, comprising a memory, a processor, and a computer program stored in the memory and executable on the processor, wherein the processor implements the steps of the seismic profile-guided CPT data mixed field interpolation algorithm when executing the computer program.

[0074] The present invention adopts the above technical solution, which has the following advantages:

[0075] The present invention combines two-dimensional seismic data and borehole location data to model and identify geological lithology, and creates a regional soil parameter distribution framework by establishing a geological three-dimensional model, so as to achieve the purpose of predicting the distribution of soil parameters at any point in the stratum;

[0076] The present invention has the advantage of fusing seismic data and borehole location data, ensuring the accuracy of stratum rock and soil parameter identification, reflecting the true situation of seabed geology as much as possible, and thus providing reliable reference opinions for engineering. BRIEF DESCRIPTION OF THE DRAWINGS

[0077] Various other advantages and benefits will become apparent to those of ordinary skill in the art by reading the detailed description of the preferred embodiments below. The accompanying drawings are only for the purpose of illustrating the preferred embodiments and are not to be considered as limiting the present invention. Throughout the accompanying drawings, the same reference numerals are used to represent the same components. In the accompanying drawings:

[0078] Figure 1 It is the original image of the 2D seafloor seismic data;

[0079] Figure 2 It is the result obtained by using the inverse distance interpolation method for the CPT cone tip resistance based on 2D seismic data;

[0080] Figure 3 It is the profile result obtained by using the mixed neighborhood interpolation method for the CPT cone tip resistance based on the two-dimensional seismic data;

[0081] Figure 4 The model top view and cross-sectional view are obtained by using the mixed neighborhood interpolation method of the method of the present invention to CPT cone tip resistance based on two-dimensional seismic data. DETAILED DESCRIPTION

[0082] The exemplary embodiments of the present invention will be described in more detail below with reference to the accompanying drawings. Although the exemplary embodiments of the present invention are shown in the accompanying drawings, it should be understood that the present invention can be implemented in various forms and should not be limited by the embodiments described herein. On the contrary, these embodiments are provided in order to enable a more thorough understanding of the present invention and to fully convey the scope of the present invention to those skilled in the art.

[0083] The present invention provides a seismic profile-guided CPT data hybrid field interpolation algorithm, which combines two-dimensional seismic data and borehole location data to model and identify geological lithology, and creates a regional soil parameter distribution framework by establishing a geological three-dimensional model, thereby achieving the purpose of predicting the distribution of soil parameters at any point in the stratum.

[0084] The present invention will be further described below in conjunction with the embodiments shown in the accompanying drawings.

[0085] The present invention provides a seismic profile-guided CPT data mixed domain interpolation algorithm, comprising the following steps:

[0086] Step 1: Process the two-dimensional seismic image to obtain attribute features;

[0087] Step 2: Based on the attribute characteristics, calculate the spatial attributes and geological attributes of the points to be interpolated. The spatial attributes include the horizontal and vertical coordinates of the points, the layer thickness, and the depth of the points in the layer. The stratigraphic data preparation stage of the region has been divided into 7 regions. The geological attributes include the region where the points are located and the stratigraphic attributes.

[0088] Step 3: Based on the spatial attributes and geological attributes, CPT attribute parameters are extracted, and the CPT attribute parameter data is combined and merged with the features;

[0089] Step 4: Slice the combined and merged data horizontally, and process the slices along the direction indicated by the two-dimensional seismic data to obtain horizontal interpolation points and interpolation parameters to form multiple slice data sets; for each slice data set, process it according to the acquired two-dimensional seismic data, including analyzing and extracting the seismic features in the slice data (including amplitude, phase, envelope, etc.);

[0090] Step 5: In the transverse slice data set, the interpolation model is constructed using the two-dimensional seismic data as the interpolation parameter and the cone tip resistance obtained at the CPT hole position as the interpolation result;

[0091] Step 6: Based on the interpolation model, the nearest neighbor interpolation method and the inverse distance interpolation method are used for mixed field interpolation. The nearest neighbor interpolation method and the inverse distance weighted interpolation method are combined to fully consider the influence of spatial attributes and geological attributes to achieve smooth and accurate interpolation (such as Figure 3 ); and establish a three-dimensional model to form the final interpolation result (as shown in Figure 4 shown).

[0092] In the above embodiment, preferably, the two-dimensional seismic image is processed in step 1 to obtain attribute features, which specifically includes the following steps:

[0093] Step 1.1, the two-dimensional seismic image is in segy format, which is a standard seismic data format used to store and exchange seismic reflection data (such as Figure 1 As shown in Figure 2), seismic signal data is usually expressed as a time series, and the value of each sample point is the amplitude, which is:

[0094] Amplitude(t)=x(t) (1)

[0095] Where Amplitude(t) refers to the amplitude of the seismic signal at time t; x(t) is the function of the seismic signal over time, which represents the recorded seismic wave signal;

[0096] Step 1.2, calculate the phase and envelope characteristics:

[0097] The calculation of phase and envelope usually uses Hilbert transform, which is a tool in signal processing. Hilbert transform converts a real-valued signal into a complex-valued signal. The real part of the complex-valued signal is the original signal, and the imaginary part is its Hilbert transform. The formula of Hilbert transform is as follows:

[0098]

[0099] Where H(x(t)) is the signal after Hilbert transformation in time; x(τ) is the seismic signal value at time τ; t is the current time point; τ is the integral variable, representing all points on the time axis; PV is the principal value integral;

[0100] Through the Hilbert transform, we can get the complex valued signal, specifically:

[0101] z(t)=x(t)+jH(x(t)) (3)

[0102] Where z(t) is the complex value signal at time t; x(t) is the original seismic signal value at time t; j is the imaginary unit;

[0103] The instantaneous phase is the phase angle of a complex-valued signal and is expressed as:

[0104]

[0105] Where φ(t) is the instantaneous phase at time t; Im(z(t)) is the result of Hilbert transform, which is also the imaginary part of the complex-valued signal at time t; Re(z(t)) is the original seismic signal, which is also the real part of the complex-valued signal at time t;

[0106] The envelope is the modulus of the replicated signal and is calculated as:

[0107]

[0108] Where E(t) is the envelope amplitude at time t; |z(t)| is the modulus of the complex-valued signal at time t.

[0109] In the above embodiment, preferably, the step 2 of calculating the spatial attributes and geological attributes of the points requiring interpolation includes the following steps:

[0110] Step 2.1: When calculating the depth of a point in a layer, the depth needs to be normalized. Normalizing the point depth is to normalize the depth of the seismic data to a standard range (usually 0 to 1). This makes it easier and more consistent to compare data between different strata thicknesses, and improves the accuracy of interpolation. Normalization first determines the thickness range of the strata, and then calculates the normalized value of the point. The formula for calculating normalization is:

[0111]

[0112] In the formula, d is the original depth of the point; min is the minimum thickness of the formation, that is, the depth of the top layer of the formation; max is the maximum thickness of the formation, that is, the depth of the bottom layer of the formation; d NORM The value of the normalized depth of the point.

[0113] In the above embodiment, preferably, the extraction of CPT attribute parameters in step 3 and combining and merging the CPT attribute parameter data with the features specifically include the following steps:

[0114] Step 3.1, select the data in the area where the CPT hole is located as training data, and use the data obtained from the CPT as test data to model and analyze the geological formation by obtaining the parameters of the cone tip resistance. c The calculation formula is:

[0115]

[0116] In the formula, F c A is the cone tip penetration force, which is the vertical force exerted on the cone tip during its penetration into the soil layer. This is directly measured by the sensor of the CPT device. c is the cross-sectional area of ​​the cone tip, calculated using a standard conical probe, using the formula:

[0117]

[0118] Where D c is the diameter of the cone tip;

[0119] Step 3.2, obtain the two-dimensional seismic characteristics of the CPT hole location, including the amplitude, phase and envelope, merge and combine the spatial attributes and geological attributes of the point, and use them together with the cone tip resistance data as a training data set;

[0120] Step 3.3: Merge and combine the two-dimensional seismic features of the points to be interpolated, including amplitude, phase and envelope, with the spatial attributes and geological attributes of the points, and wait for the interpolation model to be trained and input to obtain the overall interpolation result.

[0121] In the above embodiment, preferably, in step 5, in the transverse slice data set, the two-dimensional seismic data is used as the interpolation parameter, and the cone tip resistance obtained at the CPT hole position is used as the interpolation result to construct an interpolation model, which specifically includes the following steps:

[0122] Step 5.1: Create a collection of the cone tip resistance data set for the sample points as follows:

[0123]

[0124] In the formula, A collection of cone tip resistance data sets for sample points; f1 represents the cone tip resistance of the first sample point. K represents the cone tip resistance at the Kth sample point;

[0125] Its spatial attributes and two-dimensional seismic attributes are established as follows:

[0126]

[0127] Among them, x1 represents the amplitude, phase, envelope and spatial attributes of a CPT hole position;

[0128] Step 5.2, solve the diagonal equation, which is as follows:

[0129]

[0130] Where D(x) represents the metric tensor field; t(x) refers to the distance from the interpolation point to the nearest sample point x. k The time required; here, time is just an abbreviation for non-Euclidean distance, that is, characteristic distance; Represents points other than sample points, that is, points to be interpolated; Represents the gradient of the function t(x), that is, the spatial rate of change of t(x) at point x;

[0131] The calculation of the metric tensor field D(x) provides the anisotropy and spatial variation coefficient of the equation. By using the metric tensor field D(x) to measure the characteristic distance, it reflects the influence of geological structure on the characteristic distance, making two points in the same geological structure very close, while two points in different structures are much farther away;

[0132] Using spatial attributes to calculate t(x), the calculation formula of t(x) is as follows:

[0133]

[0134] In the formula, x mainly contains the spatial attributes of the interpolation points and sample points, that is, the spatial position; x k is the spatial attribute of the kth sample point; the formula for calculating t(x) for the sample point is as follows:

[0135]

[0136] Step 5.3: Use the nearest neighbor interpolation method to find the nearest known data point for each interpolation point, and assign the value of the known data point to the interpolation point. For known sample points, use the cone tip resistance as the characteristic value;

[0137] Step 5.4: Use the inverse distance interpolation method on the data. The inverse distance interpolation method makes the closest point have a greater impact on the interpolation result, while the farthest point has a smaller impact on the interpolation result (such as Figure 2 The calculation formula of the inverse distance interpolation method is:

[0138]

[0139] In the formula, q c (x) is the interpolation value of the point to be interpolated; N is the number of sample points that need to be provided in the vicinity.

[0140] In the above embodiment, preferably, the formula for the mixed neighborhood interpolation in step 6 is as follows:

[0141] q c (P) = w NN q c (Q nearest )+w IDW q c (IDW) (15)

[0142] In the formula, w NN is the weight of the nearest neighbor interpolation; q c (Q nearest ) is the result of nearest neighbor interpolation; w IDW is the weight of the inverse distance interpolation, q c (IDW) is the result of inverse distance interpolation; q c (P) is the final interpolation result.

[0143] The present invention also provides a seismic profile guided CPT data mixed field interpolation, which is characterized by comprising:

[0144] A first processing unit is used to process the two-dimensional seismic image to obtain attribute features;

[0145] The second processing unit is used to calculate the spatial attributes and geological attributes of the points to be interpolated based on the attribute characteristics;

[0146] A third processing unit is used to extract CPT attribute parameters based on the spatial attributes and geological attributes, and to combine and merge the CPT attribute parameter data with the features;

[0147] The fourth processing unit is used to perform transverse slicing on the combined and merged data, and the slices are processed along the direction indicated by the two-dimensional seismic data to obtain transverse interpolation points and interpolation parameters to form multiple slice data sets; for each slice data set, processing is performed according to the acquired two-dimensional seismic data, including analyzing and extracting seismic features in the slice data;

[0148] The fifth processing unit is used to construct an interpolation model in the transverse slice data set using the two-dimensional seismic data as the interpolation parameter and the cone tip resistance obtained at the CPT hole position as the interpolation result;

[0149] The sixth processing unit is used to perform mixed field interpolation based on the interpolation model using the nearest neighbor interpolation method and the inverse distance interpolation method, combining the nearest neighbor interpolation and the inverse distance weighted interpolation, fully considering the influence of spatial attributes and geological attributes, and achieving smooth and accurate interpolation.

[0150] The present invention also provides a computer storage medium on which a computer program is stored. When the computer program is executed by a processor, the steps of the seismic profile-guided CPT data mixed field interpolation algorithm are implemented.

[0151] The present invention also provides a computer device, comprising a memory, a processor, and a computer program stored in the memory and executable on the processor, wherein the processor implements the steps of the seismic profile-guided CPT data mixed field interpolation algorithm when executing the computer program.

[0152] Finally, it should be noted that the above embodiments are only used to illustrate the technical solutions of the present invention, rather than to limit it. Although the present invention has been described in detail with reference to the aforementioned embodiments, those skilled in the art should understand that they can still modify the technical solutions described in the aforementioned embodiments, or make equivalent replacements for some of the technical features therein. However, these modifications or replacements do not deviate the essence of the corresponding technical solutions from the spirit and scope of the technical solutions of the embodiments of the present invention.

Claims

1. A seismic profile-guided CPT data mixed domain interpolation algorithm, characterized in that: The steps include: Process the two-dimensional seismic images to obtain attribute features; Based on the attribute characteristics, the spatial attributes and geological attributes of the points to be interpolated are calculated; Based on the spatial attributes and geological attributes, CPT attribute parameters are extracted, and the CPT attribute parameter data are combined and merged with the features; The combined and merged data are sliced ​​horizontally, and the slices are processed along the direction indicated by the two-dimensional seismic data to obtain the horizontal interpolation points and interpolation parameters to form multiple slice data sets; for each slice data set, it is processed according to the acquired two-dimensional seismic data; In the transverse slice data set, the interpolation model is constructed using the 2D seismic data as the interpolation parameter and the cone tip resistance obtained from the CPT hole position as the interpolation result; On the basis of the interpolation model, the nearest neighbor interpolation method and the inverse distance interpolation method are used for mixed field interpolation. The nearest neighbor interpolation and the inverse distance weighted interpolation are combined to fully consider the influence of spatial attributes and geological attributes to achieve smooth and accurate interpolation. A three-dimensional model is then established to form the final interpolation result.

2. The CPT data mixed domain interpolation algorithm according to claim 1, characterized in that: The processing of the two-dimensional seismic image to obtain the attribute features specifically includes the following steps: The two-dimensional seismic image is in the segy format file, which is a standard seismic data format used to store and exchange seismic reflection data. Seismic signal data is usually expressed as a time series, and the value of each sample point is the amplitude, that is: Amplitude(t)=x(t) (1) Where Amplitude(t) refers to the amplitude of the seismic signal at time t; x(t) is the function of the seismic signal over time, which represents the recorded seismic wave signal; Compute phase and envelope characteristics: The calculation of phase and envelope usually uses Hilbert transform, which is a tool in signal processing. Hilbert transform converts a real-valued signal into a complex-valued signal. The real part of the complex-valued signal is the original signal, and the imaginary part is its Hilbert transform. The formula of Hilbert transform is as follows: Where H(x(t)) is the signal after Hilbert transformation in time; x(τ) is the seismic signal value at time τ; t is the current time point; τ is the integral variable, representing all points on the time axis; PV is the principal value integral; Through the Hilbert transform, we can get the complex valued signal, specifically: z(t)=x(t)+jH(x(t)) (3) Where z(t) is the complex value signal at time t; x(t) is the original seismic signal value at time t; j is the imaginary unit; The instantaneous phase is the phase angle of a complex-valued signal and is expressed as: Where φ(t) is the instantaneous phase at time t; Im(z(t)) is the result of Hilbert transform, which is also the imaginary part of the complex-valued signal at time t; Re(z(t)) is the original seismic signal, which is also the real part of the complex-valued signal at time t; The envelope is the modulus of the replicated signal and is calculated as: Where E(t) is the envelope amplitude at time t; |z(t)| is the modulus of the complex-valued signal at time t.

3. The CPT data mixed domain interpolation algorithm according to claim 1, characterized in that: The calculation requires the spatial attributes and geological attributes of the interpolation points, and includes the following steps: When calculating the depth of a point in a layer, the depth needs to be normalized. The formula for calculating the normalization is: In the formula, d is the original depth of the point; min is the minimum thickness of the formation, that is, the depth of the top layer of the formation; max is the maximum thickness of the formation, that is, the depth of the bottom layer of the formation; d NORM The value of the normalized depth of the point.

4. The CPT data mixed domain interpolation algorithm according to claim 1, characterized in that: The extracting of CPT attribute parameters and combining and merging CPT attribute parameter data with features specifically includes the following steps: The data in the area where the CPT hole is located is selected as the training data, and the data obtained from the CPT is used as the test data. The parameters of the cone tip resistance are obtained to model and analyze the geological formation. c The calculation formula is: In the formula, F c is the cone tip penetration force; A c is the cross-sectional area of ​​the cone tip, calculated using a standard conical probe, using the formula: Where D c is the diameter of the cone tip; Obtain the 2D seismic characteristics of the CPT hole location, including the amplitude, phase and envelope, merge and combine with the spatial attributes and geological attributes of the point, and use it together with the cone tip resistance data as a training data set; The two-dimensional seismic characteristics of the points to be interpolated, including amplitude, phase and envelope, are merged and combined with the spatial attributes and geological attributes of the points, and the overall interpolation results are obtained after the interpolation model is trained and input.

5. The CPT data mixed domain interpolation algorithm according to claim 1, characterized in that: In the transverse slice data set, the two-dimensional seismic data is used as the interpolation parameter, and the cone tip resistance obtained at the CPT hole position is used as the interpolation result to construct an interpolation model, which specifically includes the following steps: The cone tip resistance data set for the sample points is established as follows: In the formula, is a set of cone tip resistance data sets for sample points; f1 represents the cone tip resistance of the first sample point; f K represents the cone tip resistance at the Kth sample point; Its spatial attributes and two-dimensional seismic attributes are established as follows: Among them, x1 represents the amplitude, phase, envelope and spatial attributes of a CPT hole position; Solve the diagonal equation, which is as follows: Where D(x) represents the metric tensor field; t(x) refers to the distance from the interpolation point to the nearest sample point x. k the time required; Represents points other than sample points, that is, points to be interpolated; Represents the gradient of the function t(x), that is, the spatial rate of change of t(x) at point x; The calculation of the metric tensor field D(x) provides the anisotropy and spatial variation coefficient of the equation. By using the metric tensor field D(x) to measure the characteristic distance, it reflects the influence of geological structure on the characteristic distance, making two points in the same geological structure very close, while two points in different structures are much farther away; Using spatial attributes to calculate t(x), the calculation formula of t(x) is as follows: In the formula, x mainly contains the spatial attributes of the interpolation points and sample points, that is, the spatial position; x k is the spatial attribute of the kth sample point; the formula for calculating t(x) for the sample point is as follows: The nearest neighbor interpolation method is used for the data. For each point to be interpolated, the nearest known data point is found and the value of the known data point is assigned to the point to be interpolated. For the known sample points, the cone tip resistance is used as the characteristic value. The inverse distance interpolation method is used for the data. The inverse distance interpolation method makes the closest point have a greater impact on the interpolation result, and the farthest point has a smaller impact on the interpolation result. The calculation formula of the inverse distance interpolation method is: In the formula, q c (x) is the interpolation value of the point to be interpolated; N is the number of sample points that need to be provided in the vicinity.

6. The CPT data mixed domain interpolation algorithm according to claim 1, characterized in that: The formula for the mixed neighborhood interpolation is as follows: q c (P)=w NN q c (Q nearest )+w IDW q c (IDW) (15) In the formula, w NN is the weight of the nearest neighbor interpolation; q c (Q nearest ) is the result of nearest neighbor interpolation; w IDW is the weight of the inverse distance interpolation, q c (IDW) is the result of inverse distance interpolation; q c (P) is the final interpolation result.

7. A seismic profile guided CPT data mixed field interpolation, characterized in that: include: A first processing unit is used to process the two-dimensional seismic image to obtain attribute features; The second processing unit is used to calculate the spatial attributes and geological attributes of the points to be interpolated based on the attribute characteristics; A third processing unit is used to extract CPT attribute parameters based on the spatial attributes and geological attributes, and to combine and merge the CPT attribute parameter data with the features; The fourth processing unit is used to perform transverse slicing on the combined and merged data, and the slices are processed along the direction indicated by the two-dimensional seismic data to obtain transverse interpolation points and interpolation parameters to form multiple slice data sets; for each slice data set, processing is performed according to the acquired two-dimensional seismic data, including analyzing and extracting seismic features in the slice data; The fifth processing unit is used to construct an interpolation model in the transverse slice data set using the two-dimensional seismic data as the interpolation parameter and the cone tip resistance obtained at the CPT hole position as the interpolation result; The sixth processing unit is used to perform mixed field interpolation based on the interpolation model using the nearest neighbor interpolation method and the inverse distance interpolation method, combining the nearest neighbor interpolation and the inverse distance weighted interpolation, fully considering the influence of spatial attributes and geological attributes, and achieving smooth and accurate interpolation.

8. A computer storage medium having a computer program stored thereon, characterized in that: When the computer program is executed by a processor, the steps of the seismic profile-guided CPT data hybrid field interpolation algorithm described in any one of claims 1 to 6 are implemented.

9. A computer device comprising a memory, a processor and a computer program stored in the memory and executable on the processor, characterized in that: When the processor executes the computer program, the seismic profile-guided CPT data mixed field interpolation algorithm steps described in any one of claims 1-6 are implemented.

Images