Carbonate karst cave identification method based on spectral decomposition and attribute analysis

Abstract

The invention discloses a carbonate karst cave identification method based on spectral decomposition and attribute analysis, and the method comprises the steps: collecting well data and three-dimensional post-stack seismic data, determining a target layer interface through well-seismic calibration, and tracking a horizon; constructing a two-dimensional stratigraphic framework model containing a karst cave complexity parameter, defining karst cave seismic response characteristics, and performing frequency domain S transform frequency division processing by using a frequency weight coefficient; reconstructing single-frequency body features, combining F-X domain filtering and denoising, extracting various attributes such as intrinsic coherence, gradient mode and the like, introducing a boundary enhancement coefficient, improving a gray-level co-occurrence matrix and other innovative models to optimize attribute extraction; attribute slices are optimized through attribute comprehensive evaluation parameters, and karst cave distribution is highlighted through RGB fusion. According to the method, the defects of the traditional technology in data processing and karst cave feature recognition and attribute analysis are overcome, accurate recognition of the carbonate karst cave is achieved, and powerful technical support is provided for related exploration and development.

Description

Technical Field

[0001] The present invention relates to the field of carbonate cave identification, and in particular to a carbonate cave identification method based on spectrum decomposition and attribute analysis. Background Art

[0002] In the field of energy exploration and resource development, the precise identification of carbonate caves, often serving as important reservoirs for resources such as oil and gas, is crucial for exploration and development. With the continued growth of global energy demand, the efficient and accurate exploration of carbonate cave resources not only helps improve energy extraction efficiency but also holds strategic importance for ensuring energy supply security. However, traditional cave identification techniques have increasingly exposed their limitations when applied to carbonate caves in complex geological conditions, necessitating the urgent need for more advanced and effective identification methods.

[0003] Existing technologies have significant flaws in data processing. Traditional spectrum analysis methods lack targeted processing of data frequency differences, making it impossible to fully tap the potential value of cave information in data of different frequencies. When processing seismic data, due to the lack of appropriate weighting of different frequencies, data features within the effective frequency band are not prominent, making it difficult to clearly present cave response characteristics, thus affecting subsequent assessments of cave location and morphology. Furthermore, conventional data denoising methods, while removing noise, tend to over-smooth the data, resulting in a loss of cave feature details and reducing data usability.

[0004] Traditional methods also have significant shortcomings in cave feature identification and attribute analysis. On the one hand, the construction of cave models is overly simplified, failing to fully consider the impact of complex factors such as cave pore area and fracture length on cave morphology and distribution. This results in significant deviations between the simulated cave models and actual conditions, making it difficult to accurately determine the seismic response characteristics of the caves. On the other hand, the models and algorithms used in the attribute extraction process are relatively simple and traditional, failing to comprehensively and accurately extract various attribute information of the caves, such as boundary characteristics and secondary statistical attributes. This results in a less detailed and accurate depiction of the caves, severely restricting the efficiency and accuracy of carbonate cave exploration and development. Summary of the Invention

[0005] In order to overcome the shortcomings and deficiencies of the existing technology, the present invention provides a carbonate cave identification method based on spectrum decomposition and attribute analysis.

[0006] The technical solution adopted by the present invention is a carbonate cave identification method based on spectrum decomposition and attribute analysis, comprising the following steps: Step S1: collecting well data, wherein the well data includes well layer data and well logging curves, wherein the well logging curves include acoustic wave time difference curves DTC, volume density curves DEN, and gamma curves GR; and simultaneously collecting 3D post-stack seismic data; Step S2: Based on the collected 3D post-stack seismic data and well logging curves, the top interface Ttop and bottom interface Tbot of the target layer in the study area are determined by well-seismic calibration, and the horizon tracking operation is carried out in the entire area; Step S3: Based on the parameters of the target carbonate rock formation in the study area, namely velocity, density, and porosity, a two-dimensional stratigraphic framework model is constructed. A cave model with different characteristics is embedded in the model to generate a two-dimensional seismic forward profile. The forward profile is compared with the original profile to finally determine that the seismic response characteristics of the caves in the study area are characterized by beaded bright spot reflections. Step S4: Performing spectrum analysis on the 3D post-stack seismic data to obtain an effective frequency band of 10-60 Hz for the seismic data. Within this effective frequency band, frequency domain S-transformation is used to perform frequency division processing on the 3D post-stack seismic data, with a 2 Hz interval to obtain frequency division data volumes with discrete frequencies of 10 Hz, 12 Hz, 14 Hz, 16 Hz, 18 Hz, and 20 Hz, respectively. Step S5: Based on the response characteristics of the cave locations in each discrete frequency data volume, single frequency volumes with better responses to caves are selected from the obtained discrete frequency data volumes, namely, single frequency volumes of 16 Hz, 18 Hz, and 20 Hz, and feature reconstruction is performed on the selected single frequency volumes; Step S6: Perform FX domain filtering on the data volume after feature reconstruction to remove random noise and linear interference, thereby further improving the signal-to-noise ratio of the data volume; Step S7: extracting intrinsic coherence attributes from the data volume after denoising in S6. Based on the obtained intrinsic coherence attributes, the gradient modulus attributes are extracted using the intrinsic coherence attributes as input data, thereby characterizing the boundary characteristics of the cave. Step S8: performing GST dip and azimuth attribute extraction on the data volume after denoising in step S6, using the extracted dip and azimuth attributes as input data, and using the gray level co-occurrence matrix method to extract secondary statistical seismic attributes, namely, difference and entropy seismic attributes; Step S9: Decomposing the data volume after denoising in step S6 based on the structural gradient tensor method to obtain eigenvalues, i.e., Lamda-mean and Lamda-max seismic attributes; Step S10: extracting root mean square attributes and instantaneous frequency seismic attributes from the data volume after denoising in step S6; Step S11: extracting gradient modulus, difference, entropy, Lamda-max, Lamda-mean, root mean square and instantaneous frequency seismic attribute slices respectively, and selecting Lamda-max, difference and gradient modulus attribute slices that better describe the planar characteristics of the cave from the extracted attribute slices; Step S12: performing plane attribute RGB fusion on the gradient modulus, difference, and Lamda-max attribute slices that are selected to better depict the cave plane, so as to further highlight the distribution characteristics of the cave in the study area.

[0007] Furthermore, in step S3, when constructing the two-dimensional stratigraphic framework model, the cave complexity parameter is introduced. , through the formula Calculate, where is the total area of the cave pores, is the total area of the model, is the total length of the crack, is the total length of the model boundary, is the weight coefficient, and ; According to the calculated cave complexity parameters , adjust the morphology and distribution of the embedded cave model to simulate the actual situation of the caves in the study area.

[0008] Furthermore, in step S4, when performing the frequency domain S transform method frequency division processing, a frequency weight coefficient is introduced , for different frequencies The data body, the processed data value By formula Calculate, where The frequency before treatment is The data value of the data body, satisfy , , and the coefficient is used to perform weighted processing on data bodies of different frequencies.

[0009] Furthermore, in step S5, when reconstructing the features of the single-frequency volume, a feature fusion model is constructed to , the reconstructed data By formula Calculate, where For the The sub-data of a single frequency body at different spatial positions, is the number of sub-data, is the fusion weight coefficient, and According to the intensity and stability of the cave response characteristics in the single-frequency body, the model is used to reconstruct the characteristics of the single-frequency body and highlight the cave characteristics.

[0010] Furthermore, in step S7, when extracting the gradient modulus attribute, a boundary enhancement coefficient is introduced , for a point in the data body The gradient modulus , the enhanced gradient modulus By formula Calculate, where , is the standard deviation of the local data at this point, is the standard deviation of the entire data volume. This coefficient is used to enhance the gradient modulus at the cave boundary and characterize the cave boundary characteristics.

[0011] Furthermore, in step S8, when extracting secondary statistical earthquake attributes using the gray level co-occurrence matrix method, an improved gray level co-occurrence matrix model is constructed to , through the formula Calculate, where is the grayscale level, is the gray-level co-occurrence matrix element, is the weight coefficient, Pixel The distance between the two; the entropy of earthquake attributes , through the formula Calculate, where is the weight coefficient, ,The secondary statistical seismic attributes of caves are extracted by improving the model.

[0012] Furthermore, in step S9, when decomposing and obtaining eigenvalues based on the structural gradient tensor method, a tensor correction coefficient is introduced. , calculate the eigenvalue of the structure gradient tensor, and the corrected eigenvalue ( Represents the eigenvalue number) through the formula Calculate, where is the eigenvalue before correction, is the difference between the local area and the average density, is the average density of the study area, and the characteristic value is corrected by this coefficient to reflect the structural characteristics of the cave.

[0013] Furthermore, in step S10, when extracting the root mean square attribute, a weighted root mean square model is constructed to extract the data within a certain time window in the data body. , weighted root mean square attribute value By formula Calculate, where is the time window length, is the weight coefficient, for The frequency corresponding to the moment is used to highlight the amplitude characteristics of the cave at different frequencies through this model.

[0014] Furthermore, in step S11, when optimizing attribute slices, attribute comprehensive evaluation parameters are introduced. , for each attribute slice, comprehensive evaluation parameters By formula Calculate, where is the contrast between the cave area and the background area in the attribute slice, is the clarity of the cave boundary in the attribute slice, The integrity of the cave morphology in the attribute slice, is the weight coefficient, and , according to the calculated The value is optimized to slice the characteristic attributes of the cave plane.

[0015] Furthermore, in step S12, when performing plane attribute RGB fusion, a color intensity adjustment model is constructed to adjust the color intensity of a certain pixel in the fused image. , through the formula Calculate, where 、 、 are the color values of the red, green, and blue channels respectively, is the adjustment factor, and These are the attribute values corresponding to the three attribute slices. The model is used to highlight the distribution characteristics of caves in the study area.

[0016] Beneficial effects: The present invention proposes a carbonate cave identification method based on spectral decomposition and attribute analysis. In the data foundation construction stage, well data and three-dimensional post-stack seismic data are collected, and the target layer interface and tracking layer are determined through well-seismic calibration, laying a solid foundation for subsequent precise analysis. In terms of cave feature identification, the cave complexity parameter is introduced when constructing the two-dimensional stratigraphic framework model, and the cave hole area and fracture length are comprehensively considered to simulate the cave morphology and distribution more realistically, helping to clarify the cave seismic response characteristics; the frequency domain S transform frequency division processing uses the frequency weight coefficient to enhance the differences between different frequency data bodies, making the cave response characteristics easier to capture; the fusion model of single-frequency body feature reconstruction, the boundary enhancement coefficient of gradient modulus attribute extraction and other technical means improve the accuracy and clarity of cave identification from the perspectives of highlighting cave characteristics and clearly depicting cave boundaries. At the attribute analysis level, a variety of innovative models such as the improved grayscale co-occurrence matrix model and the weighted root mean square model are optimized for the extraction of different seismic attributes to accurately obtain cave-related attribute information. During the results integration stage, attribute slices are optimized through comprehensive attribute evaluation parameters, and then the color intensity adjustment model is adjusted through RGB fusion to highlight the distribution characteristics of caves in the study area, realizing full-process optimization from data processing to intuitive presentation of cave distribution, providing more accurate and efficient technical support for the exploration and development of carbonate caves, and helping to improve exploration efficiency and the accuracy and reliability of resource development. BRIEF DESCRIPTION OF THE DRAWINGS

[0017] Figure 1 The figure is a flow chart of the method steps of the present invention. DETAILED DESCRIPTION

[0018] It should be noted that, unless there is a conflict, the embodiments in this application and the features described in the embodiments can be combined with each other. The application is further described in detail below with reference to the accompanying drawings and specific embodiments.

[0019] like Figure 1 As shown in FIG, a carbonate cave identification method based on spectrum decomposition and attribute analysis includes the following steps: Step S1: collecting well data, wherein the well data includes well layer data and well logging curves, wherein the well logging curves include acoustic wave time difference curves DTC, volume density curves DEN, and gamma curves GR; and simultaneously collecting 3D post-stack seismic data; Specifically, in step S1, the first task is to collect basic data. Well data collection covers several key aspects. Well stratification data can clarify the distribution and boundaries of different strata, providing fundamental stratigraphic information for subsequent research. The DTC (Dialectic Time of Flight) curve, DEN (Dense Density Emission) curve, and GR (Gamma Ray Ray) curve in the well log reflect the physical properties of the subsurface rock from different dimensions. The DTC (Dialectic Time of Flight) curve can be used to determine rock porosity and lithology, the DEN (Dense Density Emission) curve provides information on rock density, and the GR (Gamma Ray Ray Ray) curve reflects the content of radioactive elements in the rock, assisting in distinguishing different rock types. Simultaneously, 3D post-stack seismic data is collected. This data, acquired through seismic exploration techniques, records the response of the subsurface geologic body to seismic waves and serves as a crucial foundation for subsequent cave identification and analysis. The collection of these data lays a solid foundation for the entire carbonate cave identification process. Only with accurate and comprehensive basic data can subsequent analysis be carried out smoothly and the reliability of the results be guaranteed.

[0020] When collecting well data, it is necessary to ensure the accuracy and completeness of the data, and to collect and record it in accordance with strict industry standards and specifications. For well logging curve data, the accuracy and stability of the measuring instruments must be guaranteed. During the measurement process, the instruments must be regularly calibrated and maintained to ensure that the acquired data truly reflects the physical properties of the underground rock. When collecting three-dimensional post-stack seismic data, it is necessary to determine appropriate acquisition parameters, such as the excitation energy of the seismic waves and the spacing between receiving channels, to ensure that the collected data can clearly reflect the underground geological structure. Reasonable acquisition parameter settings can enable seismic waves to better penetrate the underground medium and accurately receive reflected seismic signals, thereby providing high-quality data for subsequent analysis.

[0021] Step S2: Based on the collected 3D post-stack seismic data and well logging curves, the top interface Ttop and bottom interface Tbot of the target layer in the study area are determined by well-seismic calibration, and the horizon tracking operation is carried out in the entire area; Specifically, step S2 carries out well-seismic calibration based on the three-dimensional post-stack seismic data and well logging curves collected by S1. The core purpose of well-seismic calibration is to establish a corresponding relationship between well logging data and seismic data, so as to accurately determine the top interface Ttop and bottom interface Tbottom of the target layer in the study area. In the specific implementation process, first, according to the lithological and physical characteristics of the well logging curve, the corresponding seismic reflection phase axis is found on the seismic profile. Since the well logging curve reflects the detailed information at the wellbore, and the seismic data reflects the underground geological structure information of a larger range, by comparing the relationship between the two in time and depth, the depth information of the well logging curve is converted into time information, so that it corresponds to the time domain of the seismic data. In this way, the top and bottom interface positions of the target layer on the seismic profile can be accurately determined. This process is crucial for accurately defining the scope of the research target.

[0022] After the top and bottom interfaces of the target layer are determined, the layer tracking operation is immediately implemented throughout the entire area. Layer tracking uses the continuity and similarity of the reflection phase axis in the seismic data to extend and connect the top and bottom interfaces of the target layer determined at a single well throughout the entire study area. During the tracking process, it is necessary to comprehensively consider the attribute characteristics of the seismic data such as amplitude, frequency, phase, and the changing trend of the geological structure. For abnormal conditions such as interruption and distortion of the seismic reflection phase axis, it is necessary to analyze and judge them in combination with geological knowledge, eliminate interference caused by geological factors such as faults and formation pinch-outs, and ensure that the tracked layer can truly reflect the spatial distribution of the target layer. Accurate layer tracking can provide an accurate stratigraphic framework for subsequent geological modeling, attribute analysis, and other work.

[0023] Step S3: Based on the parameters of the target carbonate rock formation in the study area, namely velocity, density, and porosity, a two-dimensional stratigraphic framework model is constructed. A cave model with different characteristics is embedded in the model to generate a two-dimensional seismic forward profile. The forward profile is compared with the original profile to finally determine that the seismic response characteristics of the caves in the study area are characterized by beaded bright spot reflections. Specifically, step S3 constructs a two-dimensional stratigraphic framework model based on the parameters of the target carbonate formation in the study area. These stratigraphic parameters, including velocity, density, and porosity, are key indicators describing the physical properties of carbonate formations. Velocity reflects the speed at which seismic waves propagate through rock, density affects the reflection and transmission of seismic waves, and porosity is closely related to the rock's reservoir properties. When constructing the model, different stratigraphic parameters are appropriately set according to the actual geological distribution, forming a two-dimensional framework that reflects the stratigraphic structure of the study area. This framework model is then embedded with cave models with different characteristics, simulating the morphology, size, and distribution of actual caves. Seismic forward modeling is performed on the model to generate a two-dimensional seismic forward profile. Comparing the forward profile with the original seismic profile reveals that caves exhibit a characteristic beaded-like bright spot reflection on the seismic profile. The identification of this characteristic provides an important basis and marker for the subsequent identification of caves in actual seismic data.

[0024] When constructing a two-dimensional stratigraphic framework model and embedding a cave model, the geological background and actual exploration data of the study area must be fully considered. The selection of stratigraphic parameters should be based on regional geological data, exploration results in adjacent areas, and actual well logging data to ensure parameter accuracy and rationality. When setting the characteristics of the cave model, the geological conditions of cave development, such as rock type and tectonic movement, must be comprehensively considered. By continuously adjusting and optimizing model parameters, the generated seismic forward profiles are made as similar as possible to the original profiles, thereby improving the accuracy of the seismic response characteristics of the caves. Accurately understanding the seismic response characteristics of the caves enables the rapid and accurate identification of their location and extent in subsequent seismic data analysis.

[0025] Step S4: Performing spectrum analysis on the 3D post-stack seismic data to obtain an effective frequency band of 10-60 Hz for the seismic data. Within this effective frequency band, frequency domain S-transformation is used to perform frequency division processing on the 3D post-stack seismic data, with a 2 Hz interval to obtain frequency division data volumes with discrete frequencies of 10 Hz, 12 Hz, 14 Hz, 16 Hz, 18 Hz, and 20 Hz, respectively. Specifically, step S4 performs spectral analysis on the 3D post-stack seismic data to determine the effective frequency band range of the data. Spectral analysis can reveal the energy distribution of different frequency components in the seismic data. In practice, specialized spectral analysis software and algorithms are used to process the seismic data and generate spectral distribution curves. Analysis revealed that the effective frequency band of the seismic data in the study area is 10-60 Hz. Within this frequency range, the seismic data contains a wealth of information about the underground geological structure, particularly information related to karst caves. After determining the effective frequency band, the 3D post-stack seismic data is subjected to frequency-domain S-transformation processing. The frequency-domain S-transform is a method that allows for joint analysis of seismic data in the time and frequency domains. It decomposes the raw seismic data into sub-data volumes of varying frequencies. In this processing, discrete frequency-domain data volumes of 10 Hz, 12 Hz, 14 Hz, 16 Hz, 18 Hz, and 20 Hz are generated, using a 2 Hz interval. This frequency-domain processing allows for more detailed observation of the characteristic variations in the seismic data at different frequencies, providing a foundation for subsequent screening of frequency data volumes with a strong response to karst caves.

[0026] When performing spectral analysis and frequency division processing, it is necessary to rationally select analysis parameters and algorithms based on the geological characteristics of the study area and the actual seismic data. For spectral analysis, it is important to ensure that the time window length is appropriate, encompassing sufficient seismic information while avoiding interference between adjacent windows. When using the frequency domain S-transform method for frequency division processing, the appropriate frequency interval must be determined based on the effective frequency band and the required accuracy. Smaller frequency intervals provide more detailed frequency information but increase data processing volume and computation time; larger frequency intervals may miss some important frequency features. Therefore, choosing a 2Hz frequency interval represents a reasonable trade-off between ensuring research accuracy and improving processing efficiency. The discrete frequency data volume generated through frequency division processing can reveal underground geological structures from different frequency perspectives, facilitating more accurate cave identification.

[0027] Step S5: Based on the response characteristics of the cave locations in each discrete frequency data volume, single frequency volumes with better responses to caves are selected from the obtained discrete frequency data volumes, namely, single frequency volumes of 16 Hz, 18 Hz, and 20 Hz, and feature reconstruction is performed on the selected single frequency volumes; Specifically, step S5 selects single-frequency bodies that have a good response to caves from the obtained discrete frequency data volumes based on the response characteristics of the cave locations in each discrete frequency data volume. In actual operation, it is necessary to conduct a detailed analysis and comparison of each discrete frequency data volume to observe the manifestation of the cave in the data volumes of different frequencies, such as amplitude and morphological characteristics. Through analysis, it was found that the single-frequency bodies of 16Hz, 18Hz, and 20Hz have a more obvious response to caves. The caves show relatively prominent characteristics in these frequency data volumes, which can more clearly reflect the location and morphology of the caves. Therefore, the single-frequency bodies of these three frequencies are selected as the key data for subsequent analysis. After the single-frequency bodies are selected, feature reconstruction is performed on them. Feature reconstruction is intended to further enhance the characteristic manifestation of the cave in the single-frequency body, making it more prominent and obvious. By using specific signal processing algorithms and techniques to process the single-frequency body data, the characteristic parameters such as amplitude and phase of the data are adjusted, thereby highlighting the abnormal response of the cave and suppressing the influence of other interfering factors.

[0028] When optimizing single-frequency volumes and reconstructing features, it's necessary to combine geological knowledge and practical exploration experience. When assessing the response characteristics of caves in data volumes of different frequencies, it's important to consider the impact of factors such as cave size and burial depth on seismic wave propagation. Generally speaking, smaller caves may respond more clearly in higher-frequency seismic data, while larger caves can also perform well in lower-frequency data. During feature reconstruction, appropriate reconstruction algorithms and parameters should be selected based on the characteristics of the optimized single-frequency volume data. Different algorithms and parameter settings will enhance cave features to varying degrees. Through repeated trial and error, the optimal reconstruction scheme is determined to maximize the prominence of cave features in the single-frequency volume, providing clearer and more accurate data for subsequent cave identification and analysis.

[0029] Step S6: Perform FX domain filtering on the data volume after feature reconstruction to remove random noise and linear interference, thereby further improving the signal-to-noise ratio of the data volume; Specifically, step S6 performs FX domain filtering on the data volume after feature reconstruction. FX domain filtering is a filtering method based on the frequency-space domain that can effectively remove random noise and linear interference in the data volume. During the acquisition and processing of seismic data, various noises and interferences are inevitably introduced, which can affect the identification and analysis of cave features. Random noise has irregular characteristics, which can reduce the signal-to-noise ratio of seismic data and mask the effective signals of caves; linear interference usually manifests as some regular phase axes, which can be confused with the actual geological reflection signals and interfere with the judgment of the location and morphology of caves. FX domain filtering analyzes the characteristics of the data in the frequency-space domain, establishes a suitable filtering model, and identifies and removes noise and interference. In practical applications, according to the characteristics of the data volume and the type of noise interference, the filtering parameters, such as the filter window size and filter coefficient, are adjusted to ensure that while removing noise and interference, the effective signals of the cave are retained to the maximum extent, thereby further improving the signal-to-noise ratio of the data volume.

[0030] When performing FX-domain filtering, the data volume requires detailed analysis and evaluation to understand the distribution and characteristics of noise and interference. Spectral analysis and coherence analysis can be used to determine the frequency range and spatial distribution of noise and interference. Based on these analysis results, FX-domain filtering parameters can be appropriately set, employing different filtering strategies for different types of noise and interference. For example, high-frequency random noise can be suppressed using low-pass filtering; linear interference can be removed by designing targeted filtering operators leveraging its specific characteristics in the frequency-space domain. Through precise parameter settings and appropriate filtering strategies, noise and interference can be effectively removed from the data volume, improving data quality and providing a more reliable foundation for subsequent attribute extraction and cave identification.

[0031] Step S7: extracting intrinsic coherence attributes from the data volume after denoising in S6. Based on the obtained intrinsic coherence attributes, the gradient modulus attributes are extracted using the intrinsic coherence attributes as input data, thereby characterizing the boundary characteristics of the cave. Specifically, step S7 extracts the intrinsic coherence attribute of the data body after denoising in S6. The intrinsic coherence attribute is an attribute used to describe the similarity of seismic data, which can reflect the continuity and consistency of underground geological bodies. In carbonate rock formations, the presence of caves will destroy the continuity of the formation, resulting in changes in the similarity of seismic data at the boundaries of the caves. By extracting the intrinsic coherence attribute, this change can be highlighted, thereby identifying the boundary position of the cave. In the actual extraction process, a special algorithm is used to calculate the denoised data body. The algorithm is based on the waveform characteristics of the seismic data and calculates the intrinsic coherence value of each sampling point by comparing the similarity between adjacent seismic traces. These intrinsic coherence values constitute an intrinsic coherence attribute data body. In this data body, the intrinsic coherence value at the boundary of the cave is usually low, while the intrinsic coherence value in the continuous area of the formation is high, thereby forming a clear difference, which is convenient for identifying the boundary of the cave.

[0032] Based on the intrinsic coherence attributes, the gradient modulus attributes are used as input data to extract the intrinsic coherence attributes. The gradient modulus attributes are a further processing of the intrinsic coherence attributes, which can more clearly characterize the detailed features of the cave boundary. The gradient modulus attributes highlight the gradient changes at the boundary by calculating the spatial rate of change of the intrinsic coherence attributes. During the calculation process, the gradient operator is applied to the intrinsic coherence attribute data volume to obtain the gradient modulus value of each sampling point. These gradient modulus values can more accurately reflect the location and morphology of the cave boundary, capturing even subtle boundary changes. By extracting and analyzing the intrinsic coherence attributes and gradient modulus attributes, the cave boundary characteristics can be comprehensively and meticulously characterized, providing an important basis for subsequent cave morphological analysis and quantitative evaluation.

[0033] Step S8: Perform GST dip and azimuth attribute extraction on the data volume after denoising in step S6, use the extracted dip and azimuth attributes as input data, and then use the gray level co-occurrence matrix method to extract its secondary statistical seismic attributes, namely, difference and entropy seismic attributes; Specifically, step S8 performs GST dip and azimuth attribute extraction on the data body after denoising in S6. The GST (General Structure Tensor) method is a technology that can extract stratigraphic structural information from seismic data. This method can be used to calculate the dip and azimuth attributes of the stratum. In the identification of carbonate caves, the dip and azimuth attributes can reflect the spatial distribution morphology and structural characteristics of the stratum, which is of great significance for determining the relationship between the cave and the stratigraphic structure. In the actual extraction process, based on the principle of the GST method, tensor calculation is performed on the denoised data body to obtain the structural tensor information of each sampling point. Then, based on the eigenvalues and eigenvectors of the structural tensor, the dip and azimuth of the stratum are calculated. These dip and azimuth attribute data bodies can intuitively display the tilt direction and angle changes of the stratum in space.

[0034] The extracted dip and azimuth attributes are used as input data, and the gray-level co-occurrence matrix method is used to extract their secondary statistical seismic attributes, namely, the difference and entropy seismic attributes. The gray-level co-occurrence matrix method is a method for analyzing image texture features. Here, seismic attribute data is treated as image data. The difference seismic attribute reflects the intensity of grayscale changes in the image. In seismic attribute data, it can reflect the degree of stratigraphic changes. The presence of caves can cause stratigraphic changes, thereby changing the difference attribute value. The entropy seismic attribute reflects the uniformity of the grayscale distribution in the image. In seismic attribute data, it can reflect the complexity of the stratigraphic structure. The entropy attribute value will also vary when the stratigraphic structure around a cave is relatively complex. By extracting these two secondary statistical seismic attributes, the characteristics of the cave can be further characterized from different perspectives, providing more information for cave identification and analysis.

[0035] Step S9: Decomposing the data volume after denoising in step S6 based on the structural gradient tensor method to obtain its eigenvalues, namely Lamda-mean and Lamda-max seismic attributes; Specifically, step S9 decomposes the denoised data volume in S6 using the structural gradient tensor method to obtain its eigenvalues, namely the Lamda-mean and Lamda-max seismic attributes. The structural gradient tensor method is an analysis method based on the gradient information of seismic data, which can describe the spatial variation characteristics of seismic data. In carbonate cave identification, the presence of caves causes spatial gradient variations in seismic data. By performing structural gradient tensor decomposition on the data volume, eigenvalues reflecting these variations can be extracted. In practice, the gradient information of the data volume in three spatial directions (x, y, and z) is first calculated, and then a structural gradient tensor matrix is constructed. Eigenvalue decomposition is performed on this tensor matrix to obtain eigenvalues such as the Lamda-mean (average eigenvalue) and Lamda-max (maximum eigenvalue). These eigenvalues can reflect the changes in the stratigraphic structure surrounding the cave. For example, near the cave boundary, due to the sudden change in the stratigraphic structure, the Lamda-max value will be relatively large, and the Lamda-mean value will also change accordingly.

[0036] The Lamda-mean and Lamda-max seismic attributes can describe the structural characteristics of caves from different perspectives. The Lamda-mean attribute reflects the average degree of stratigraphic variation and can be used to assess the overall structural complexity of a cave region. The Lamda-max attribute more prominently reflects the maximum value of stratigraphic variation, clearly demonstrating cave boundaries and areas of intense local structural variation. Analysis of these two eigenvalue attributes provides a more comprehensive understanding of the structural morphology and distribution of caves, providing a crucial basis for their accurate identification and evaluation.

[0037] Step S10: extracting root mean square attributes and instantaneous frequency seismic attributes from the data volume after denoising in step S6; Specifically, step S10 extracts the root mean square attribute and instantaneous frequency seismic attribute of the data body after denoising in S6. The root mean square attribute is a commonly used seismic attribute, which reflects the average intensity of seismic wave energy by calculating the root mean square value of the amplitude of seismic data in a certain time window. In the identification of carbonate caves, the presence of caves will affect the propagation and reflection of seismic waves, resulting in changes in seismic wave energy. By extracting the root mean square attribute, this energy change can be captured, thereby identifying the location of the cave. In the actual extraction process, according to the research needs and data characteristics, the appropriate time window length is determined, and the root mean square value of the denoised data body is calculated in each time window. The root mean square value at different positions reflects the size of the seismic wave energy at that position. Due to its special geological structure, the cave area will cause the seismic wave energy to be abnormal, thereby showing characteristics different from the surrounding strata in the root mean square attribute data body.

[0038] Instantaneous frequency seismic attributes reflect the temporal frequency variations of seismic waves. In carbonate formations, the presence of caves alters the propagation path and velocity of seismic waves, thereby affecting their frequency characteristics. By extracting these instantaneous frequency attributes, it is possible to analyze the frequency variations of seismic waves and identify the location and extent of caves. To extract these attributes, a specialized algorithm is used to process the denoised data volume. This algorithm calculates the instantaneous frequency value at each sampling point based on the phase variations of the seismic waves. These instantaneous frequency values constitute the instantaneous frequency attribute data volume. Analysis of this data volume reveals that the instantaneous frequency of the cave area differs from that of the surrounding strata, providing important information for cave identification.

[0039] Step S11: extracting gradient modulus, difference, entropy, Lamda-max, Lamda-mean, root mean square and instantaneous frequency seismic attribute slices respectively, and selecting Lamda-max, difference and gradient modulus attribute slices that better describe the planar characteristics of the cave from the extracted attribute slices; Specifically, step S11 extracts gradient modulus, difference, entropy, Lamda-max, Lamda-mean, root mean square and instantaneous frequency seismic attribute slices respectively. Attribute slices are two-dimensional data extracted from a three-dimensional attribute data body along a specific plane, which can intuitively display the distribution characteristics of the attributes on the plane. In this operation, according to the geological conditions of the study area and the needs of cave identification, appropriate slice positions and directions are selected, and slices are extracted from the corresponding attribute data body respectively. These attribute slices reflect the characteristics of the cave from different angles. For example, the gradient modulus attribute slice can show the distribution of the cave boundary on the plane, the difference and entropy attribute slices can reflect the changes and complexity of the stratigraphic structure in the cave area, the Lamda-max and Lamda-mean attribute slices help to understand the structural morphology of the cave, and the root mean square and instantaneous frequency attribute slices can reflect the influence of the cave on the energy and frequency of seismic waves.

[0040] From the extracted attribute slices, select the Lamda-max, dissimilarity, and gradient norm attribute slices that best depict the planar features of the cave. During this selection process, each attribute slice is carefully analyzed and compared to evaluate its ability to depict the planar features of the cave. The quality of each attribute slice is comprehensively judged by observing the clarity of the cave, the accuracy of its boundaries, and the degree of consistency with the actual geological conditions.

[0041] Step S12: performing plane attribute RGB fusion on the gradient modulus, difference, and Lamda-max attribute slices that are selected to better depict the cave plane, so as to further highlight the distribution characteristics of the cave in the study area.

[0042] Specifically, step S12 is based on the Lamda-max, difference and gradient mode attribute slices selected in S11, and aims to achieve intuitive visualization of the distribution characteristics of the cave through the plane attribute RGB fusion technology. During implementation, the three attribute slices are first normalized, and the attribute values are uniformly mapped to the same numerical interval to eliminate the fusion interference caused by differences in attribute dimensions and numerical ranges. Subsequently, according to specific mapping rules, the Lamda-max attribute slice is mapped to the red channel. This attribute is highly sensitive to the structural differences between the cave and the surrounding strata and can highlight the unique structure of the cave; the difference attribute slice is assigned to the green channel, which can effectively reflect the degree of attribute changes inside and around the cave area; the gradient mode attribute slice is mapped to the blue channel, which accurately outlines the cave boundary contour with its high resolution of boundary changes. During the fusion process, the weight of each channel can be dynamically adjusted according to actual research needs. For example, when focusing on cave boundary identification, the weight of the blue channel is increased to enhance boundary contrast. The final RGB fusion image presents the three attribute information in the form of color overlay. By observing the color distribution and changes of the image, geological researchers can quickly and accurately grasp the distribution range, morphological direction and relationship of the caves in the study area with the surrounding strata, providing intuitive and critical decision-making basis for the exploration and development of carbonate caves.

[0043] Preferably, in step S3, when constructing the two-dimensional stratigraphic framework model, the cave complexity parameter is introduced. , through the formula Calculate, where is the total area of the cave pores, is the total area of the model, is the total length of the crack, is the total length of the model boundary, is the weight coefficient, and ; According to the calculated cave complexity parameters , the morphology and distribution of the embedded cave model were adjusted to more accurately simulate the actual conditions of the caves in the study area.

[0044] Specifically, a cave complexity parameter is introduced when constructing a two-dimensional stratigraphic framework model. This parameter comprehensively considers factors such as the total cave pore area, the total model area, the total fracture length, and the total length of the model boundary. By setting two weight coefficients, whose sum equals 1, the influence of cave pore area and fracture length in the complexity calculation can be flexibly adjusted according to the actual geological conditions. In implementation, the cave complexity parameter is first calculated based on detailed geological data of the study area, accurately obtaining data such as cave pore area and fracture length. The morphology and distribution of the caves embedded in the model are then dynamically adjusted based on the value of this parameter. A high parameter value indicates a high degree of pore and fracture development in the cave area. In this case, the cave model is adjusted to a more complex and irregular morphology, and its distribution density is increased within the model. Conversely, a low value simplifies the cave model morphology and reduces its distribution. This results in a model that better reflects the actual geological conditions, facilitates more accurate determination of the seismic response characteristics of caves, and provides a more reliable reference standard for subsequent cave identification based on seismic data.

[0045] Preferably, in step S4, when performing the frequency domain S transform method frequency division processing, a frequency weight coefficient is introduced , for different frequencies The data body, its processed data value By formula Calculate, where The frequency before treatment is The data value of the data body, satisfy , The coefficient is used to perform weighted processing on data bodies of different frequencies, thereby enhancing the differences of frequency data bodies within the effective frequency band and improving the cave identification effect.

[0046] Specifically, a frequency weighting coefficient is introduced during the frequency division process of the frequency-domain S-transform method. This coefficient is determined based on the varying contributions of different frequencies within the effective frequency band to cave identification, and its value gradually increases as the frequency moves from the lower limit to the upper limit of the effective frequency band. In practice, when processing each discrete frequency data volume, a corresponding weighting coefficient is first determined based on its corresponding frequency. This coefficient is then multiplied by the pre-processed data value of the corresponding frequency volume to obtain the processed data value. This weighting method gives low-frequency data volumes relatively low weights, while high-frequency data volumes relatively high weights, effectively enhancing the differences between data volumes of different frequencies. For example, in cave identification, high-frequency data volumes are more sensitive to cave details, making their features more prominent after weighting. Low-frequency data volumes, on the other hand, preserve the overall cave outline without excessively interfering with high-frequency features. This significantly improves the effectiveness of cave identification across the entire effective frequency band, making the response characteristics of caves at different frequencies easier to analyze and identify.

[0047] Preferably, in step S5, when reconstructing the features of the single-frequency volume, a feature fusion model is constructed to reconstruct the features of each single-frequency volume data. , the reconstructed data By formula Calculate, where For the The sub-data of a single frequency body at different spatial positions, is the number of sub-data, is the fusion weight coefficient, and According to the intensity and stability of the cave response characteristics in the single-frequency body, the model is used to reconstruct the characteristics of the single-frequency body and highlight the cave characteristics.

[0048] Specifically, when reconstructing features from a monofrequency volume, a feature fusion model is constructed that fully considers the sub-data at different spatial locations within the monofrequency volume and the corresponding fusion weight coefficients. The fusion weight coefficients are determined based on the strength and stability of the cave response characteristics within the monofrequency volume. In implementation, each monofrequency volume is first spatially partitioned to generate multiple sub-data regions. Next, the cave response characteristics within each sub-data region, such as amplitude and morphological integrity, are analyzed to assess their strength and stability. Sub-data regions with distinct and stable cave response characteristics are assigned higher fusion weight coefficients; conversely, sub-data regions with distinct and stable cave response characteristics are assigned lower weights. Each sub-data region is then multiplied by the corresponding fusion weight coefficient and summed to produce the reconstructed data. This model enhances the portions of the monofrequency volume with prominent cave characteristics and suppresses interfering information, making the reconstructed monofrequency volume more prominent in cave features. This facilitates a clearer and more accurate capture of cave representations in the single-frequency data during subsequent cave identification, providing a higher-quality data foundation for refined cave identification.

[0049] Preferably, in step S7, when extracting the gradient modulus attribute, a boundary enhancement coefficient is introduced , for a point in the data body The gradient modulus , the enhanced gradient modulus By formula Calculate, where , is the standard deviation of the local data at this point, is the standard deviation of the entire data volume. This coefficient is used to enhance the gradient modulus at the cave boundary and more clearly characterize the cave boundary characteristics.

[0050] Specifically, the process of extracting the gradient modulus attribute in step S7 is optimized and upgraded, introducing a boundary enhancement coefficient. This coefficient is obtained by calculating the ratio of the standard deviation of the local data at a certain point to the standard deviation of the entire data volume. In actual implementation, for each point in the data volume, the standard deviation of the data in the local area (for example, a neighborhood of a certain size centered on the point) is first calculated, and the standard deviation of the entire data volume is simultaneously obtained to obtain the boundary enhancement coefficient. This coefficient is then multiplied by the original gradient modulus value of the point to obtain the enhanced gradient modulus value. At the cave boundary, due to the drastic changes in the stratum structure, the local data standard deviation is relatively large, and the calculated boundary enhancement coefficient is also large, resulting in a significant enhancement of the gradient modulus value. In non-boundary areas, the local data standard deviation is smaller, and the gradient modulus value is enhanced less. In this way, the characteristics of the cave boundary in the data volume can be more clearly highlighted, making the depiction of the cave boundary more precise, helping geological researchers to more accurately determine the scope and morphology of the cave.

[0051] Preferably, in step S8, when extracting secondary statistical earthquake attributes using the gray level co-occurrence matrix method, an improved gray level co-occurrence matrix model is constructed to , through the formula Calculate, where is the grayscale level, is the gray-level co-occurrence matrix element, is the weight coefficient, Pixel The distance between the two; the entropy of earthquake attributes , through the formula Calculate, where is the weight coefficient, , by improving the model to more accurately extract the secondary statistical seismic attributes of caves.

[0052] Specifically, when calculating differential seismic attributes, a weight coefficient related to the distance between pixels is introduced. The closer the distance, the greater the weight, thereby emphasizing the relationship between adjacent pixels. During implementation, the grayscale level is first determined, the grayscale co-occurrence matrix is constructed, and the matrix elements are calculated. Then, the weight coefficient is calculated based on the position of each pixel pair, and each parameter is substituted into the formula to calculate the differential seismic attributes, so that the calculation results can more truly reflect the grayscale value changes in the cave area. When calculating the entropy seismic attributes, a weight coefficient related to the grayscale value difference of the pixel points is introduced. The smaller the difference, the greater the weight, highlighting the distribution of pixels with similar grayscale values. Similarly, after determining the grayscale level and the grayscale co-occurrence matrix elements, the corresponding weight coefficient is calculated and substituted into the formula to calculate the entropy seismic attributes. Through the improved model, the secondary statistical seismic attributes of the cave can be more accurately extracted, the cave characteristics can be more meticulously portrayed from the perspective of grayscale changes and distribution, and richer and more accurate attribute information can be provided for cave identification.

[0053] Preferably, in step S9, when decomposing and obtaining eigenvalues based on the structural gradient tensor method, a tensor correction coefficient is introduced. , calculate the eigenvalue of the structure gradient tensor, and the corrected eigenvalue ( Represents the eigenvalue number) through the formula Calculate, where is the eigenvalue before correction, is the difference between the local area and the average density, is the average density of the study area. The characteristic value is corrected by this coefficient to make it reflect the structural characteristics of the cave more accurately.

[0054] Specifically, this coefficient is associated with the difference between the local area and the average density, as well as the average density of the study area. In the actual implementation process, the data volume is first spatially divided, the average density of each local area is calculated, and the difference is obtained by comparing it with the average density of the study area. Then, the tensor correction coefficient is calculated based on the correlation relationship, and it is combined with the eigenvalue before correction to obtain the corrected eigenvalue. In the cave area, due to the difference between its density and the surrounding strata, the calculated difference is not zero. The eigenvalue is adjusted by the tensor correction coefficient so that the eigenvalue can more accurately reflect the changes in the structural characteristics of the cave. For example, near the boundary of the cave, the density changes significantly. The corrected eigenvalue will highlight the structural differences in the area, thereby improving the accuracy of identifying the structural characteristics of the cave and providing more reliable eigenvalue data support for the precise analysis of the cave.

[0055] Preferably, in step S10, when extracting the root mean square attribute, a weighted root mean square model is constructed to extract the data within a certain time window in the data body. , weighted root mean square attribute value By formula Calculate, where is the time window length, is the weight coefficient, for The frequency corresponding to the moment is used to highlight the amplitude characteristics of the cave at different frequencies through this model.

[0056] Specifically, the model introduces a weight coefficient associated with the frequency corresponding to the time instant. As the frequency changes within the effective frequency band, the weight coefficient also changes accordingly. During implementation, the appropriate time window length is first determined. For each time window in the data volume, the weight coefficient is calculated based on its corresponding frequency. The data value and weight coefficient are then substituted into the weighted RMS formula to calculate the RMS attribute value. This model ensures that when calculating RMS attributes, data of different frequencies contribute differently due to their weighting. In cave identification, high-frequency data is sensitive to cave details. When assigned a higher weight, it contributes more significantly to the RMS attribute calculation, making the cave's amplitude characteristics at high frequencies more prominent. Low-frequency data, while retaining the cave's overall energy information, participates in the calculation appropriately, thereby more comprehensively reflecting the cave's amplitude characteristics at different frequencies and effectively improving the effectiveness of RMS attributes for cave identification.

[0057] Preferably, in step S11, when optimizing attribute slices, attribute comprehensive evaluation parameters are introduced , for each attribute slice, its comprehensive evaluation parameter By formula Calculate, where is the contrast between the cave area and the background area in the attribute slice, is the clarity of the cave boundary in the attribute slice, The integrity of the cave morphology in the attribute slice, is the weight coefficient, and , according to the calculated The attribute slices that better describe the plane characteristics of the cave are selected.

[0058] Specifically, this parameter comprehensively considers key indicators within an attribute slice, including the contrast between the cave and background areas, the clarity of the cave boundary, and the integrity of the cave morphology. The influence of each indicator is adjusted using three weighting coefficients, the sum of which equals 1. During implementation, image segmentation and feature extraction techniques are first used to capture the cave and background areas, calculating the grayscale difference between them to determine contrast. Edge detection algorithms are then used to assess the clarity of the cave boundary, and morphological analysis methods are used to determine the integrity of the cave morphology. These three indicators are then calculated for each attribute slice. Appropriate weighting coefficients are then set based on research needs. Each indicator value is multiplied by the corresponding weighting coefficient and summed to obtain the attribute comprehensive evaluation parameter. Finally, attribute slices are sorted and selected based on the value of this parameter. Attribute slices with higher parameter values are more capable of depicting the planar features of the cave. This allows the selection of the most suitable attribute slices for cave planar feature analysis, providing high-quality data for subsequent accurate identification of the planar distribution of the cave.

[0059] Preferably, in step S12, when performing plane attribute RGB fusion, a color intensity adjustment model is constructed to adjust the color intensity of a pixel in the fused image. , through the formula Calculate, where 、 、 are the color values of the red, green, and blue channels respectively, is the adjustment factor, and These are the attribute values corresponding to the three attribute slices. The model is used to highlight the distribution characteristics of caves in the study area.

[0060] Specifically, the model optimizes the color intensity of the fused image using three adjustment coefficients, which sum to 3. These coefficients are calculated based on the maximum attribute values corresponding to the three attribute slices. During implementation, the maximum attribute values for the Lamda-max, dissimilarity, and gradient modulus attribute slices are first obtained, and the adjustment coefficients are calculated according to relevant rules. During RGB fusion, the color intensity of each pixel in the fused image is calculated based on its attribute values for the three attribute slices and the adjustment coefficients. This model dynamically adjusts the color intensity of the fused image based on the characteristics of the different attribute slices, making the colors of karst caves more prominent and distinct in the image. For example, when karst caves are prominent in the Lamda-max attribute slice, the corresponding red channel adjustment coefficient is larger. This increases the red component's contribution to the overall color intensity when calculating color intensity, thereby highlighting the karst caves in the image and more intuitively demonstrating the distribution of karst caves in the study area, providing geological researchers with clearer and more accurate visualization results.

[0061] In the description of the present invention, it should be noted that, unless otherwise expressly specified or limited, the terms "disposed," "installed," "connected," "connected," and "fixed" should be understood in a broad sense. For example, they may refer to fixed connections, detachable connections, or integral connections; they may refer to mechanical connections or electrical connections; they may refer to direct connections or indirect connections through an intermediate medium; and they may refer to internal communication between two components. Those skilled in the art will understand the specific meanings of the above terms in the present invention based on specific circumstances.

[0062] While embodiments of the present invention have been shown and described, it will be understood by those skilled in the art that various equivalent changes, modifications, substitutions and variations may be made to the embodiments without departing from the principles and spirit of the invention, and that the scope of the invention is defined by the appended claims and their equivalents.

Claims

1. A carbonate cave identification method based on spectrum decomposition and attribute analysis, characterized in that: The following steps are involved: Step S1: collecting well data, wherein the well data includes well layer data and well logging curves, wherein the well logging curves include acoustic wave time difference curves DTC, volume density curves DEN, and gamma curves GR; and simultaneously collecting 3D post-stack seismic data; Step S2: Based on the collected 3D post-stack seismic data and well logging curves, the top interface Ttop and bottom interface Tbot of the target layer in the study area are determined by well-seismic calibration, and the horizon tracking operation is carried out in the entire area; Step S3: Based on the parameters of the target carbonate rock formation in the study area, namely velocity, density, and porosity, a two-dimensional stratigraphic framework model is constructed. A cave model with different characteristics is embedded in the model to generate a two-dimensional seismic forward profile. The forward profile is compared with the original profile to finally determine that the seismic response characteristics of the caves in the study area are characterized by beaded bright spot reflections. Step S4: Performing spectrum analysis on the 3D post-stack seismic data to obtain an effective frequency band of 10-60 Hz for the seismic data. Within this effective frequency band, frequency domain S-transformation is used to perform frequency division processing on the 3D post-stack seismic data. With an interval of 2 Hz, frequency division data volumes with discrete frequencies of 10 Hz, 12 Hz, 14 Hz, 16 Hz, 18 Hz, and 20 Hz are obtained. Step S5: Based on the response characteristics of the cave locations in each discrete frequency data volume, single frequency volumes with better responses to caves are selected from the obtained discrete frequency data volumes, namely, single frequency volumes of 16 Hz, 18 Hz, and 20 Hz, and feature reconstruction is performed on the selected single frequency volumes; Step S6: Perform FX domain filtering on the data volume after feature reconstruction to remove random noise and linear interference, thereby improving the signal-to-noise ratio of the data volume; Step S7: extracting intrinsic coherence attributes from the data volume after denoising in S6. Based on the obtained intrinsic coherence attributes, the gradient modulus attributes are extracted using the intrinsic coherence attributes as input data, thereby characterizing the boundary characteristics of the cave. Step S8: performing GST dip and azimuth attribute extraction on the data volume after denoising in step S6, using the extracted dip and azimuth attributes as input data, and using the gray level co-occurrence matrix method to extract secondary statistical seismic attributes, namely, difference and entropy seismic attributes; Step S9: Decomposing the data volume after denoising in step S6 based on the structural gradient tensor method to obtain eigenvalues, i.e., Lamda-mean and Lamda-max seismic attributes; Step S10: extracting root mean square attributes and instantaneous frequency seismic attributes from the data volume after denoising in step S6; Step S11: extracting gradient modulus, difference, entropy, Lamda-max, Lamda-mean, root mean square and instantaneous frequency seismic attribute slices respectively, and selecting Lamda-max, difference and gradient modulus attribute slices that better describe the planar characteristics of the cave from the extracted attribute slices; Step S12: performing plane attribute RGB fusion on the gradient modulus, difference, and Lamda-max attribute slices that are selected to better depict the cave plane, so as to further highlight the distribution characteristics of the cave in the study area.

2. The method for identifying carbonate caves based on spectrum decomposition and attribute analysis according to claim 1, characterized in that: In step S3, when constructing the two-dimensional stratigraphic framework model, the cave complexity parameter is introduced. , through the formula Calculate, where is the total area of the cave pores, is the total area of the model, is the total length of the crack, is the total length of the model boundary, is the weight coefficient, and ; According to the calculated cave complexity parameters , adjust the morphology and distribution of the embedded cave model to simulate the actual situation of the caves in the study area.

3. The method for identifying carbonate caves based on spectrum decomposition and attribute analysis according to claim 1, characterized in that: In step S4, when performing the frequency domain S transform method frequency division processing, the frequency weight coefficient is introduced. , for different frequencies The data body, the processed data value By formula Calculate, where The frequency before treatment is The data value of the data body, satisfy , , and the coefficient is used to perform weighted processing on data bodies of different frequencies.

4. The method for identifying carbonate caves based on spectrum decomposition and attribute analysis according to claim 1, characterized in that: In step S5, when reconstructing the features of the single-frequency volume, a feature fusion model is constructed to reconstruct the features of each single-frequency volume data. , the reconstructed data By formula Calculate, where For the The sub-data of a single frequency body at different spatial positions, is the number of sub-data, is the fusion weight coefficient, and According to the intensity and stability of the cave response characteristics in the single-frequency body, the model is used to reconstruct the characteristics of the single-frequency body and highlight the cave characteristics.

5. The method for identifying carbonate caves based on spectrum decomposition and attribute analysis according to claim 1, characterized in that: In step S7, when extracting the gradient modulus attribute, a boundary enhancement coefficient is introduced. , for a point in the data body The gradient modulus , the enhanced gradient modulus By formula Calculate, where , is the standard deviation of the local data at this point, is the standard deviation of the entire data volume. This coefficient is used to enhance the gradient modulus at the cave boundary and characterize the cave boundary characteristics.

6. The method for identifying carbonate caves based on spectrum decomposition and attribute analysis according to claim 1, characterized in that: In step S8, when extracting secondary statistical earthquake attributes using the gray level co-occurrence matrix method, an improved gray level co-occurrence matrix model is constructed to analyze the difference in earthquake attributes. , through the formula Calculate, where is the grayscale level, is the gray-level co-occurrence matrix element, is the weight coefficient, Pixel The distance between the two; the entropy of earthquake attributes , through the formula Calculate, where is the weight coefficient, ,The secondary statistical seismic attributes of caves are extracted by improving the model.

7. The method for identifying carbonate caves based on spectrum decomposition and attribute analysis according to claim 1, characterized in that: In step S9, when decomposing and obtaining eigenvalues based on the structural gradient tensor method, a tensor correction coefficient is introduced. , calculate the eigenvalue of the structure gradient tensor, and the corrected eigenvalue ( Represents the eigenvalue number) through the formula Calculate, where is the eigenvalue before correction, is the difference between the local area and the average density, is the average density of the study area, and the characteristic value is corrected by this coefficient to reflect the structural characteristics of the cave.

8. The method for identifying carbonate caves based on spectrum decomposition and attribute analysis according to claim 1, characterized in that: In step S10, when extracting the root mean square attribute, a weighted root mean square model is constructed to extract the data within a certain time window in the data body. , weighted root mean square attribute value By formula Calculate, where is the time window length, is the weight coefficient, for The frequency corresponding to the moment is used to highlight the amplitude characteristics of the cave at different frequencies through this model.

9. The method for identifying carbonate caves based on spectrum decomposition and attribute analysis according to claim 1, characterized in that: In step S11, when optimizing attribute slices, attribute comprehensive evaluation parameters are introduced. , for each attribute slice, comprehensive evaluation parameters By formula Calculate, where is the contrast between the cave area and the background area in the attribute slice, is the clarity of the cave boundary in the attribute slice, The integrity of the cave morphology in the attribute slice, is the weight coefficient, and , according to the calculated The value is optimized to slice the characteristic attributes of the cave plane.

10. The method for identifying carbonate caves based on spectrum decomposition and attribute analysis according to claim 1, characterized in that: In step S12, when performing plane attribute RGB fusion, a color intensity adjustment model is constructed to adjust the color intensity of a pixel in the fused image. , through the formula Calculate, where 、 、 are the color values of the red, green, and blue channels respectively, is the adjustment factor, and These are the attribute values corresponding to the three attribute slices. The model is used to highlight the distribution characteristics of caves in the study area.

Images