Method, device and electronic equipment for automatic identification and occurrence tracking of electrical imaging logging fractures

Abstract

The present invention relates to the field of imaging logging fracture identification technology, and is a method, device and electronic equipment for automatic identification and occurrence tracking of electrical imaging logging fractures, including preprocessing of FMI dynamic logging image data to be analyzed; performing Hough transform on the preprocessed FMI dynamic logging image data to obtain Hough transform values, and obtaining Hough correction coefficients from an error function constructed by a deviation distribution relative to the Hough mean; and obtaining two types of quasi-Hough characteristic transform values ​​based on the Hough transform values ​​and the Hough correction coefficients. The present invention introduces a quasi-Hough characteristic transform function, and uses the Hough correction coefficient to correct the Hough transform value to obtain two types of quasi-Hough characteristic transform values. Through the constraint of the Hough correction coefficient, the influence of regional anomalies with large numerical anomalies but small abnormal areas on fracture identification is reduced; the data distribution after the two types of quasi-Hough transforms has obvious Gaussian-like distribution characteristics, making it possible to directly perform occurrence analysis through statistical methods in the future.

Description

Technical Field

[0001] The invention relates to the technical field of imaging logging fracture identification, in particular to an electric imaging logging fracture automatic identification and occurrence tracking method, device and electronic equipment. Background Art

[0002] Microresistivity imaging (FMI) logging technology is of great significance for identifying the fracture distribution, geometric characteristics, lithology and formation flatness of key underground formations. In order to reduce the result deviation caused by artificial methods and improve the accuracy and intelligence of imaging logging fracture identification, many scholars have done a lot of research on image information enhancement, fracture feature extraction, and automatic fracture identification in recent years, and have achieved many results.

[0003] The first existing public patent document, with the publication number CN108181665B, discloses a method and device for determining fractures, the method comprising: obtaining electrical imaging logging data of a target area; determining an electrical imaging logging image of the target area according to the electrical imaging logging data; performing morphological processing on the electrical imaging logging image in the vertical direction to obtain a first angle fracture image; performing morphological processing on the electrical imaging logging image in the horizontal direction to obtain a second angle fracture image; and determining fractures in the target area according to the first angle fracture image and the second angle fracture image. Since the scheme takes into account the specific characteristics of the fracture direction and morphology in the electrical imaging logging image, morphological processing is first performed in the vertical direction and the horizontal direction respectively, and then the fracture is identified, thereby reducing the interference signal in the image and improving the fracture accuracy in the image, thereby solving the technical problems of low accuracy and large error in identifying fractures in the existing methods, and achieving the technical effect of accurately identifying micro fractures.

[0004] The existing public patent document 2, with the publication number of CN115829990B, discloses a natural fracture identification method based on imaging logging image processing, which relates to the technical field of fracture identification in metamorphic rock reservoirs, and the steps are as follows: Step 1, identify the static image and dynamic image of imaging logging, and take their intersection to draw a dynamic and static combined image; Step 2, first perform binary image closing operation processing on the dynamic and static combined image, and then quantitatively mark the connected domain of the binary image according to the connected domain, and calculate the connected domain attribute parameters; Step 3, based on the connected domain attributes of the entire well, the machine learning fitting model is used to calculate the surface rate and identify the fractures. Through the above method, the present invention takes into account the imaging logging characteristics of pores, natural fractures and induced fractures in metamorphic rock reservoirs, and uses parameters such as area, aspect ratio, eccentricity, side length and convex hull ratio to perform multi-parameter logistic regression operations, effectively eliminate pores and induced fractures, and finally realize accurate identification and prediction of natural fractures in metamorphic rock reservoirs.

[0005] The existing imaging logging fracture identification technology has the following main problems: 1) Some imaging logging data preprocessing methods lead to large information loss. Among them, although methods such as binarization and edge extraction reduce the complexity of the data, a large amount of amplitude information and image lateral change information are discarded, which increases the possibility of misidentification of fractures. Wavelet transform technology not only consumes a lot of calculations, but also has poor high-frequency information preservation capabilities and signal enhancement effects under high noise. 2) Hough transform and improved Hough transform feature extraction methods ignore the numerical distribution characteristics along the characteristic curve, which may cause additional information loss and misleading in non-binarization preprocessing. 3) The existing fracture identification effect is not ideal, and the workload of manual identification is large. For example, the K-Means cluster analysis method relies on the introduction of prior fracture distribution information, while the unsupervised learning method requires a large number of training samples and is not universal. 4) Many of the existing occurrence statistics are performed after fracture identification, which will introduce information bias in the identification process into the statistics, resulting in information interference and sample dilution. Summary of the invention

[0006] The present invention provides a method, device and electronic equipment for automatic identification and occurrence tracking of electrical imaging logging fractures, which overcome the shortcomings of the above-mentioned prior art and can effectively solve the problem that the Hough transform and improved Hough transform feature extraction methods in the existing imaging logging fracture identification technology ignore the numerical distribution characteristics along the characteristic curve, which may cause additional information loss in non-binarization preprocessing.

[0007] One of the technical solutions of the present invention is achieved by the following measures: a method for automatic identification and occurrence tracking of electrical imaging logging fractures, comprising:

[0008] Preprocess the FMI dynamic logging image data to be analyzed;

[0009] Perform Hough transform on the preprocessed FMI dynamic logging image data to obtain the Hough transform value, and obtain the Hough correction coefficient from the error function constructed by the deviation distribution relative to the Hough mean;

[0010] Based on the Hough transform value and the Hough correction coefficient, two types of quasi-Hough feature transform values ​​are obtained, wherein the two types of quasi-Hough feature transform values ​​are obtained using two types of quasi-Hough feature transform functions as shown below;

[0011]

[0012]

[0013] in, It is a type of quasi-Hough transform value, which is used for the identification and analysis of high-conductivity cracks; for The second type of quasi-Hough transform value is used for the identification and analysis of high-resistance cracks; is the Hough transform value; is the Hough correction coefficient.

[0014] The following are further optimizations and / or improvements to the above technical solutions:

[0015] The above-mentioned Hough transform is performed on the pre-processed FMI dynamic logging image data to obtain the Hough transform value, and the Hough correction coefficient is obtained by the error function constructed by the deviation distribution relative to the Hough mean, including:

[0016] Determine the effective pixel number of the pre-processed FMI dynamic logging image data;

[0017]

[0018] in, is the number of effective pixels; It is a counting function, which indicates the total number of points that meet the condition; It is a non-number, that is, an invalid pixel; It is the pre-processed FMI dynamic logging image data; x is the horizontal pixel coordinate; is the vertical pixel coordinate; is the maximum amplitude of the trigonometric function of the inclination angle; is the phase shift value reflecting the crack azimuth; is the total number of horizontal pixels;

[0019] Obtain the Hough transform value and the Hough correction coefficient based on the number of effective pixels;

[0020]

[0021]

[0022] in, is the Hough transform value; is the Hough correction coefficient; is the number of effective pixels; It is a counting function, which indicates the total number of points that meet the condition; It is the pre-processed FMI dynamic logging image data; x is the horizontal pixel coordinate; is the vertical pixel coordinate; is the maximum amplitude of the trigonometric function that characterizes the inclination angle; is the phase shift value reflecting the crack azimuth; is the total number of horizontal pixels; is the mean value of the preprocessed FMI dynamic logging image data; is the error function; is a controllable parameter that controls the shape of the error function.

[0023] The above preprocessing of the FMI dynamic logging image data to be analyzed includes:

[0024] Convert the FMI dynamic logging image data to be analyzed into grayscale data;

[0025] Perform transverse periodic extension on the grayscale data, and use the improved iterative double spherical harmonic function interpolation algorithm to complete the grayscale data after periodic extension, complete the blank part of the data, and obtain the completed FMI dynamic logging image data;

[0026] The completed FMI dynamic logging image data is subjected to de-averaging, edge truncation and high-pass filtering in sequence;

[0027] The complementary data part in the filtered FMI dynamic logging image data is deleted, and the lateral periodic extension part is cropped.

[0028] The above also includes the occurrence analysis of the pre-processed FMI dynamic logging image data based on two types of quasi-Hough eigentransform values, including:

[0029] Determine the deviation coefficient and the standard deviation of the two types of quasi-Hough feature transform values, and bring them into the segmentation conditions to perform FMI dynamic logging image data segmentation, where the segmentation conditions include:

[0030] If satisfied It is a suspected high-conductivity seam;

[0031] If satisfied It is a suspected high-resistance seam;

[0032] in, is a type of quasi-Hough transform value; is the second type of quasi-Hough transform value; is the coefficient of deviation; is the standard deviation of the two types of quasi-Hough eigentransform values;

[0033] For the segmented suspected crack array, according to its Coordinates and corresponding well diameter curves , calculate the depth of the fracture ,inclination and normal azimuth , the calculation formula is as follows:

[0034]

[0035] in, for The corresponding well depth; is the length corresponding to each pixel interval in the vertical direction; The depth of the well is Well diameter at the time; is the vertical pixel coordinate; is the maximum amplitude of the trigonometric function of the inclination angle; is the phase shift value reflecting the crack azimuth.

[0036] The above also includes fracture identification of pre-processed FMI dynamic logging image data, including:

[0037] Divide the suspected fracture groups obtained by the occurrence analysis into several small depth segments;

[0038] Expand all arrays of suspected fracture groups obtained from the occurrence analysis outside the azimuth domain;

[0039] Sort the array for each depth segment, including:

[0040] (1) For a certain depth segment The array is processed sequentially for Z-Score data standardization and decorrelation dimensionality reduction;

[0041] (2) Determine the minimum number of cluster points and use the density-based clustering algorithm to perform cluster analysis on the array after decorrelation and dimensionality reduction processing;

[0042] Eliminate all unclassified arrays, and obtain the average value of all arrays in the same class based on the classification results, and use it as the microcrack occurrence corresponding to this class;

[0043] Eliminate the extended part outside the azimuth domain.

[0044] The above also includes the production statistics of the pre-processed FMI dynamic logging image data, including:

[0045] The suspected fracture groups obtained by the occurrence analysis are classified according to the dip threshold interval. For each category, multiple small units are divided according to the azimuth, and the weighted fracture number in each small unit is counted to form a weighted superposition azimuth distribution, where the weighted fracture number is the sum of the quasi-Hough function values ​​of the suspected fracture groups falling in the small unit.

[0046] The azimuth distribution is deflected 90° clockwise / counterclockwise and added together to form the strike-occurrence distribution;

[0047] The weighted superposition azimuth distribution and strike-occurrence distribution are normalized to draw rose diagrams of azimuth and strike-occurrence distribution;

[0048] The suspected crack group is divided into multiple equally spaced units according to the maximum amplitude, and the number of corrected cracks in each unit is counted. The calculation formula for the number of corrected cracks is as follows:

[0049] ;

[0050] in, is the inclination effect factor;

[0051] The dip angle is further divided into multiple small units, and the sum of the number of corrected cracks falling in each small unit is counted;

[0052] The sum of the corrected fracture numbers in each small unit is normalized and a dip angle distribution histogram is drawn.

[0053] The above-mentioned obtaining of the tilt effect factor includes:

[0054] Preset a crack thickness range, numerically simulate the number of cracks that can be identified at different maximum amplitudes for each crack thickness within the preset crack thickness range, and obtain a numerical simulation result curve of the inclination effect factor;

[0055] The average value of the numerical simulation result curve of the inclination effect factor is processed and the average value is used as the corrected inclination effect factor. .

[0056] The second technical solution of the present invention is achieved by the following measures: an automatic identification and occurrence tracking device for electrical imaging logging fractures, comprising:

[0057] A preprocessing unit, for preprocessing the FMI dynamic logging image data to be analyzed;

[0058] Feature transformation unit, including:

[0059] Perform Hough transform on the preprocessed FMI dynamic logging image data to obtain the Hough transform value, and obtain the Hough correction coefficient from the error function constructed by the deviation distribution relative to the Hough mean;

[0060] Based on the Hough transform value and the Hough correction coefficient, two types of quasi-Hough feature transform values ​​are obtained, wherein the two types of quasi-Hough feature transform values ​​are obtained using two types of quasi-Hough feature transform functions as shown below;

[0061]

[0062]

[0063] in, It is a type of quasi-Hough transform value, which is used for the identification and analysis of high-conductivity cracks; It is the second type of quasi-Hough transform value, which is used for the identification and analysis of high-resistance cracks; is the Hough transform value; is the Hough correction coefficient.

[0064] The following are further optimizations and / or improvements to the above technical solutions:

[0065] The above also includes:

[0066] The occurrence analysis unit determines the deviation coefficient and the standard deviation of the two types of quasi-Hough characteristic transform values, and brings them into the segmentation conditions to perform FMI dynamic logging image data segmentation, where the segmentation conditions include:

[0067] If satisfied It is a suspected high-conductivity seam;

[0068] If satisfied It is a suspected high-resistance seam;

[0069] in, is a type of quasi-Hough transform value; is the second type of quasi-Hough transform value; is the coefficient of deviation; is the standard deviation of the two types of quasi-Hough eigentransform values;

[0070] For the segmented suspected crack array, according to its Coordinates and corresponding well diameter curves , calculate the depth of the fracture ,inclination and normal azimuth , the calculation formula is as follows:

[0071]

[0072] in, for The corresponding well depth; is the length corresponding to each pixel interval in the vertical direction; The depth of the well is Well diameter at the time; is the vertical pixel coordinate; is the maximum amplitude of the trigonometric function of the inclination angle; is the phase shift value reflecting the crack azimuth.

[0073] The fracture identification unit performs fracture identification on the pre-processed FMI dynamic logging image data, including:

[0074] Divide the suspected fracture groups obtained by the occurrence analysis into several small depth segments;

[0075] Expand all arrays of suspected fracture groups obtained from the occurrence analysis outside the azimuth domain;

[0076] Sort the array for each depth segment, including:

[0077] (1) For a certain depth segment The array is processed sequentially for Z-Score data standardization and decorrelation dimensionality reduction;

[0078] (2) Determine the minimum number of cluster points and use the density-based clustering algorithm to perform cluster analysis on the array after decorrelation and dimensionality reduction processing;

[0079] Eliminate all unclassified arrays, and obtain the average value of all arrays in the same class based on the classification results, and use it as the microcrack occurrence corresponding to this class;

[0080] Eliminate the extended part outside the azimuth domain;

[0081] The production statistics unit performs production statistics on the pre-processed FMI dynamic logging image data, including:

[0082] The suspected fracture groups obtained by the occurrence analysis are classified according to the dip threshold interval. For each category, multiple small units are divided according to the azimuth, and the weighted fracture number in each small unit is counted to form a weighted superposition azimuth distribution, where the weighted fracture number is the sum of the quasi-Hough function values ​​of the suspected fracture groups falling in the small unit.

[0083] The azimuth distribution is deflected 90° clockwise / counterclockwise and added together to form the strike-occurrence distribution;

[0084] The weighted superposition azimuth distribution and strike-occurrence distribution are normalized to draw rose diagrams of azimuth and strike-occurrence distribution;

[0085] The suspected crack group is divided into multiple equally spaced units according to the maximum amplitude, and the number of corrected cracks in each unit is counted. The calculation formula for the number of corrected cracks is as follows:

[0086] ;

[0087] in, is the inclination effect factor;

[0088] The dip angle is further divided into multiple small units, and the sum of the number of corrected cracks falling in each small unit is counted;

[0089] The sum of the corrected fracture numbers in each small unit is normalized and a dip angle distribution histogram is drawn.

[0090] The third technical solution of the present invention is achieved through the following measures: an electronic device, including a processor and a memory, wherein the memory stores a computer program, and the computer program is loaded and executed by the processor to implement the steps in the method of automatic identification of electrical imaging logging fractures and production tracking.

[0091] The beneficial effects of the present invention are as follows:

[0092] The present invention is based on improved iterative double spherical harmonic function interpolation to fill in the blank parts of the data, and performs de-averaging, edge truncation and high-pass filtering. Compared with mainstream preprocessing methods (such as binary image segmentation, wavelet transform feature extraction, image edge extraction, etc.), it retains a large amount of amplitude information and avoids the deviation caused by artificial threshold segmentation. The data values ​​after preprocessing are Gaussian distributed, ensuring the feasibility of subsequent statistical analysis.

[0093] The present invention introduces a quasi-Hough characteristic transformation function, uses the Hough correction coefficient to correct the Hough transformation value, and obtains two types of quasi-Hough characteristic transformation values. Through the constraint of the Hough correction coefficient, the influence of regional anomalies with large numerical anomalies but small abnormal areas on crack identification is reduced; and the data distribution after the two types of quasi-Hough transformation has obvious Gaussian distribution characteristics, making it possible to directly perform occurrence analysis through statistical methods in the subsequent process.

[0094] The present invention realizes data cleaning through two steps of Z-Score standardization and singular value decomposition dimensionality reduction, eliminates the dimensional effect and intrinsic coherence of the occurrence data, improves the reliability of the subsequent clustering algorithm, and uses a density-based clustering algorithm (DBSCAN clustering). When the minimum number of cluster points is artificially limited, the neighborhood radius is automatically determined, and efficient and accurate crack identification is achieved. There is no need to predetermine the number of clusters, which simplifies the workload and reduces the false identification rate of cracks.

[0095] In the present invention, the occurrence analysis is no longer a subsequent product of crack identification, but can be carried out simultaneously with crack identification. Therefore, its reliability is no longer affected by the deviation of the crack identification process, and the two types of quasi-Hough transform values ​​obtained in the two steps of data preprocessing and feature transformation in the previous step both show good Gaussian-like distribution characteristics. Therefore, in the present invention, the occurrence analysis can be directly completed using image threshold segmentation and statistical analysis based on standard deviation. Compared with the existing technical process of first identifying and then doing statistics, the recognition deviation caused by the crack identification algorithm is avoided, and the time spent on the occurrence statistics process is effectively reduced. BRIEF DESCRIPTION OF THE DRAWINGS

[0096] Attached Figure 1 A flow chart of an analysis method provided in one embodiment of the present invention.

[0097] Attached Figure 2 A flow chart of a preprocessing method provided in one embodiment of the present invention.

[0098] Attached Figure 3 A flow chart of a crack identification method provided in one embodiment of the present invention.

[0099] Attached Figure 4 A flow chart of an occurrence statistics method provided in one embodiment of the present invention.

[0100] Attached Figure 5 This is a preprocessing example diagram provided for an embodiment of the present invention.

[0101] Attached Figure 6 For attachment Figure 5 A comparison of local details.

[0102] Attached Figure 7 This is an example diagram of feature transformation provided by one embodiment of the present invention.

[0103] Attached Figure 8 This is a comparative example diagram of quasi-Hough transform crack identification (b) and traditional Hough transform crack identification (a) provided in an embodiment of the present invention.

[0104] Attached Fig. 9 A schematic diagram of occurrence statistics provided for an embodiment of the present invention.

[0105] Attached Fig.10 A schematic diagram of numerical simulation results of the inclination effect factor provided in one embodiment of the present invention.

[0106] Attached Fig.11 An embodiment of the present invention provides a distribution map of fracture occurrence between small cells obtained based on the "first identification and then statistics" process.

[0107] Attached Fig.12 This is an example diagram of high-conductivity crack identification provided by an embodiment of the present invention.

[0108] Attached Fig.13 A schematic diagram of a device structure provided for an embodiment of the present invention.

[0109] Attached Fig.14 A schematic diagram of another device structure provided for an embodiment of the present invention. DETAILED DESCRIPTION

[0110] The present invention is not limited by the following embodiments, and specific implementation methods can be determined based on the technical solution of the present invention and actual conditions.

[0111] Those skilled in the art will understand that, unless otherwise stated, a "module" or "unit" in the embodiments of the present invention refers to a computer program or a part of a computer program with a predetermined function, and works together with other related parts to achieve a predetermined goal, and can be implemented in whole or in part by using software, hardware (such as a processing circuit or memory) or a combination thereof. Similarly, a processor (or multiple processors or memories) can be used to implement one or more modules or units. In addition, each module or unit can be part of an overall module or unit that includes the function of the module or unit.

[0112] In addition, the term “plurality” in the embodiments of the present invention refers to two or more than two, and “first” and “second” are used to distinguish descriptions and should not be understood as implying relative importance.

[0113] The embodiment of the present invention provides a method, device and electronic device for automatic identification and occurrence tracking of electrical imaging logging fractures. The electrical imaging logging fracture automatic identification and occurrence tracking device can be integrated in a computer device, which can be a server or a terminal or other device; it can also be executed by a terminal and a server together, and the above examples should not be understood as limiting the present invention.

[0114] The above-mentioned terminal may include a mobile phone, a wearable smart device, a tablet computer, a notebook computer, a personal computer (PC), and a vehicle-mounted computer, etc., and the present invention does not limit this. The present invention does not limit the number of terminal devices.

[0115] The above-mentioned server can be an independent physical server, or a server cluster or distributed system composed of multiple physical servers. It can also be a cloud server that provides basic cloud computing services such as cloud services, cloud databases, cloud computing, cloud functions, cloud storage, network services, cloud communications, middleware services, domain name services, security services, CDN, as well as big data and artificial intelligence platforms. The present invention does not impose any restrictions on this.

[0116] For example, the computer equipment uses the FMI dynamic logging image data to be analyzed to perform preprocessing, feature transformation, occurrence statistics, and fracture identification. Specifically, in the preprocessing stage, an improved iterative double spherical harmonic function interpolation algorithm is used to fill in the blank part of the image, and the boundary information of high / low resistivity anomalies is enhanced by high-pass filtering, while the lateral amplitude transformation information is retained while enhancing the image information of smaller wavelengths; in the feature transformation stage, two types of quasi-Hough transforms are introduced to convert each (depth , azimuth , the maximum amplitude of trigonometric function ) Two types of quasi-Hough values ​​corresponding to coordinates and , are all derived from the original Hough transform values and a Hough correction factor defined by an error function constructed from the distribution of deviations from the Hough mean Joint control greatly reduces the impact of regional anomalies with large numerical anomalies but small abnormal areas on crack identification; in the occurrence statistics stage, the Gaussian distribution characteristics of the quasi-Hough transform data are used to perform image segmentation of the cracks based on the standard deviation, calculate the weighted occurrence distribution, and eliminate the inclination effect, so as to directly obtain the occurrence distribution; in the crack identification stage, the feature transformation results are standardized, decorrelated, and reduced in dimension, and the DBSCAN clustering algorithm with adaptive neighborhood radius is used for initial screening of crack clustering, and then similar occurrence items are merged through secondary classification to finally obtain the micro-crack occurrence and its error sequence.

[0117] Based on this, the technical solution of the present invention will be introduced and explained in combination with several examples below.

[0118] Embodiment 1: As attached Figure 1 As shown, the embodiment of the present invention discloses a method for automatic identification and occurrence tracking of electrical imaging logging fractures, comprising:

[0119] Step S110, preprocessing the FMI dynamic logging image data to be analyzed;

[0120] The relationship between images and two-dimensional arrays is reflected in the representation and storage of images in computers. Due to the need for image display and processing, images are stored in the form of two-dimensional arrays in computers, where each element represents a pixel in the image, and its value usually represents the color or grayscale value of the pixel. Therefore, the FMI dynamic logging image data to be analyzed obtained in this embodiment can also be the FMI dynamic logging two-dimensional array to be analyzed.

[0121] The FMI dynamic logging image data to be analyzed is preprocessed to eliminate noise and irrelevant information in the image and enhance the image's resistivity anomaly boundary information, thereby improving the feasibility of subsequent statistical-based fracture occurrence analysis and identification.

[0122] Step S120, performing Hough transform on the pre-processed FMI dynamic logging image data to obtain a Hough transform value, and obtaining a Hough correction coefficient from an error function constructed by a deviation distribution relative to the Hough mean;

[0123] Step S130, obtaining two types of quasi-Hough feature transformation values ​​based on the Hough transformation value and the Hough correction coefficient, wherein the two types of quasi-Hough feature transformation values ​​are obtained using two types of quasi-Hough feature transformation functions as shown below;

[0124]

[0125]

[0126] in, It is a type of quasi-Hough transform value, which is used for the identification and analysis of high-conductivity cracks; It is the second type of quasi-Hough transform value, which is used for the identification and analysis of high-resistance cracks; is the Hough transform value; is the Hough correction coefficient.

[0127] Steps S120 to S130 implement feature transformation, which can reduce interference in the image and the influence of local feature abnormal points, making subsequent crack identification and occurrence statistics more accurate.

[0128] Since the original Hough transform (Hough) and the improved Hough transform have obvious defects in identifying characteristic curves, that is, the Hough transform is prone to misidentification of cracks near areas with large numerical anomalies but small abnormal areas, especially when the feature transformation is performed on data that has not been binarized or high-quality filtered, this defect is particularly obvious. Therefore, this embodiment introduces a quasi-Hough feature transformation function (Quasi-Hough Transform, abbreviated as QHT), that is, the Hough transform value is corrected using the Hough correction coefficient to obtain two types of quasi-Hough feature transformation values. Through the constraint of the Hough correction coefficient, the influence of regional anomalies with large numerical anomalies but small abnormal areas on crack identification is reduced. At the same time, the data distribution after the two types of quasi-Hough transformations has obvious Gaussian distribution characteristics, making it possible to directly perform occurrence analysis through statistical methods in the future.

[0129] Embodiment 2: The embodiment of the present invention is a further optimization of the above embodiment, wherein the feature transformation is performed on the pre-processed FMI dynamic logging image data, specifically including:

[0130] Step S210, performing Hough transform on the pre-processed FMI dynamic logging image data to obtain a Hough transform value, and obtaining a Hough correction coefficient from an error function constructed by a deviation distribution relative to the Hough mean, including:

[0131] (1) Determine the number of effective pixels of the preprocessed FMI dynamic logging image data;

[0132]

[0133] in, is the number of effective pixels; It is a counting function, which indicates the total number of points that meet the condition; It is a non-number, that is, an invalid pixel; It is the pre-processed FMI dynamic logging image data; x is the horizontal (azimuth) pixel coordinate; is the vertical (depth) pixel coordinate; The maximum amplitude of the trigonometric function of the inclination angle, in pixels; It is the phase shift value reflecting the crack azimuth, in degrees; is the total number of pixels in the horizontal direction (azimuth);

[0134] (2) Obtaining the Hough transform value and Hough correction coefficient based on the number of effective pixels;

[0135]

[0136] in, is the Hough transform value; is the Hough correction coefficient; is the number of effective pixels; It is a counting function, which indicates the total number of points that meet the condition; It is the pre-processed FMI dynamic logging image data; x is the horizontal pixel coordinate; is the vertical pixel coordinate; is the maximum amplitude of the trigonometric function that characterizes the inclination angle; is the phase shift value reflecting the crack azimuth; is the total number of horizontal pixels; is the mean value of the preprocessed FMI dynamic logging image data; is the error function; is a controllable parameter that controls the shape of the error function ( , hour Approaching a step function, hour Approaching a linear function);

[0137] in, is the error function, and its functional formula is:

[0138]

[0139] Step S220, obtaining two types of quasi-Hough feature transformation values ​​based on the Hough transformation value and the Hough correction coefficient, wherein the two types of quasi-Hough feature transformation values ​​are obtained using two types of quasi-Hough feature transformation functions as shown below;

[0140]

[0141]

[0142] in, It is a type of quasi-Hough transform value, which is used for the identification and analysis of high-conductivity cracks; It is the second type of quasi-Hough transform value, which is used for the identification and analysis of high-resistance cracks; is the Hough transform value; is the Hough correction coefficient.

[0143] Embodiment 3: As attached Figure 2 As shown, the embodiment of the present invention is a further optimization of the above embodiment, wherein the FMI dynamic logging image data to be analyzed is preprocessed, specifically including:

[0144] Step S310, converting the FMI dynamic logging image data to be analyzed into grayscale data;

[0145] In order to facilitate subsequent processing, the FMI dynamic logging image data to be analyzed is converted into grayscale data, that is, converted into a two-dimensional array and recorded as ,in and They are the total number of pixels in the horizontal direction (azimuth) and the total number of pixels in the vertical direction (depth direction); x and are the horizontal pixel coordinates and the vertical pixel coordinates respectively. It should be noted that if the FMI dynamic logging two-dimensional array to be analyzed is obtained, there is no need to perform grayscale data conversion.

[0146] Step S320, performing transverse periodic extension on the grayscale data, and using an improved iterative double spherical harmonic interpolation algorithm to complete the grayscale data after periodic extension, to complete the blank part of the data, and obtain the completed FMI dynamic logging image data;

[0147] The above specifically include:

[0148] (1) The grayscale data is periodically expanded horizontally, that is, the azimuth range is expanded. The expanded grayscale data is recorded as .

[0149] For example, The azimuth range extends to Scope (which extends The data is data duplication).

[0150] (2) Use the improved iterative double spherical harmonic interpolation algorithm to complete the grayscale data after periodic external expansion, fill in the blank part of the data, and obtain the completed FMI dynamic logging image data, which includes:

[0151] (2.1) Constructing the initial mask function , Data blank part , the rest are 0 (if the image is read as an image containing manual interpretation lines, the manual interpretation part is also 1);

[0152] (2.2) Searching for the boundary part of the mask function , use double spherical harmonic interpolation to complete the boundary data;

[0153] (2.3) Calculate the new mask function ,like If not all are 0, return to step 2.2;

[0154] (2.4) When the final mask function is all 0, the loop is exited and the completed image is obtained, which is recorded as .

[0155] The reason for the lateral extension in this step is to ensure the lateral periodicity of the image when the data is completed, and after edge truncation, high-pass filtering, and data cropping. In addition, the improved iterative double spherical harmonic function interpolation algorithm can well complete the information when there are fewer data gaps, and when there are more data gaps, it can ensure the continuity of the image at the edge of the data gap, so that subsequent processing will not produce artifact interference at the edge of the data gap.

[0156] Step S330, performing de-averaging, edge truncation and high-pass filtering in sequence on the completed FMI dynamic logging image data;

[0157] The above specifically include:

[0158] (1) De-meaning the FMI dynamic logging image data after interpolation;

[0159] (2) After removing the mean, a tapered cosine window function with a cosine fraction is added in the horizontal dimension to achieve truncation of the horizontal signal edge, that is:

[0160]

[0161] in, is the image mean; is the total number of vertical pixels;

[0162] Among them, the tapered cosine window function As shown below:

[0163]

[0164] in, is the cosine fraction; is the total number of horizontal pixels after periodic expansion;

[0165] (4) Perform high-pass filtering on the data after signal truncation to suppress low-wavelength information, retain the high-frequency part, and enhance the relative intensity of the crack boundary information, denoted as ;

[0166] The high-pass filter is selected according to the requirements, and may be, but is not limited to, a two-dimensional fourth-order Butterworth high-pass filter.

[0167] In this step, the completed FMI dynamic logging image data is subjected to de-averaging, edge truncation and high-pass filtering in sequence, which effectively suppresses the background amplitude information of longer wavelengths while retaining the short-wavelength information and enhancing the information amplitude of the fracture (numerical anomaly) area.

[0168] Step S340, deleting the supplementary data part in the filtered FMI dynamic logging image data, and cutting out the lateral periodic extension part.

[0169] Compared with the existing mainstream preprocessing methods (such as binary image segmentation, wavelet transform feature extraction, image edge extraction, etc.), this embodiment has the following beneficial effects: (1) the preprocessed data retains a large amount of amplitude information while avoiding the deviation caused by artificial threshold segmentation; (2) the values ​​of the preprocessed data are Gaussian distributed, ensuring the feasibility of subsequent statistical analysis.

[0170] Embodiment 4: As attached Figure 3 As shown, the embodiment of the present invention is a further optimization of the above embodiment, and also includes performing fracture identification on the pre-processed FMI dynamic logging image data, specifically including:

[0171] Step S410, based on two types of quasi-Hough eigentransform values, the pre-processed FMI dynamic logging image data is subjected to occurrence analysis, and the suspected fracture groups obtained by the occurrence analysis are divided into several small depth segments, including:

[0172] (1) Determine the deviation coefficient and the standard deviation of the two types of quasi-Hough eigentransform values , and bring it into the segmentation conditions for preprocessing FMI dynamic logging image data segmentation, where the segmentation conditions include:

[0173] If satisfied It is a suspected high-conductivity seam;

[0174] If satisfied It is a suspected high-resistance seam;

[0175] in, is a type of quasi-Hough transform value; is the second type of quasi-Hough transform value; is the coefficient of deviation; is the standard deviation of the two types of quasi-Hough eigentransform values;

[0176] (2) For the suspected crack arrays (suspected high-conductivity crack group and suspected high-resistance crack group) after segmentation, Coordinates and corresponding well diameter curves , calculate the depth of the fracture ,inclination and normal azimuth , the calculation formula is as follows:

[0177]

[0178] in, for The corresponding well depth; is the length corresponding to each pixel interval in the vertical direction (in meters); The depth of the well is Well diameter at the time, (in m); is the vertical pixel coordinate; is the maximum amplitude of the trigonometric function of the inclination angle; is the phase shift value reflecting the crack azimuth.

[0179] Step S420, performing azimuth extension on all arrays of suspected fracture groups obtained from the occurrence analysis.

[0180] For example, the azimuth range is , then expand it to ,in Copied from data, Copied from data.

[0181] Step S430, classifying the array of each depth segment, including:

[0182] (1) For a certain depth segment The array is processed sequentially for Z-Score data standardization and decorrelation dimensionality reduction;

[0183] This step includes:

[0184] (1.1) For each depth segment Array is Z-Score data standardized to obtain the standardized array , removing the dimension effect, the standardized formula is:

[0185]

[0186] in, , , , For each depth segment The mean of , , , The variation of standard error;

[0187] (1.2) Decorrelate and reduce the dimension of the standardized array, that is, perform singular value decomposition, dimension reduction and calculate the right matrix of the standardized array, including:

[0188] (1.2.1) Calculate the 4x4 covariance matrix of the normalized array ;

[0189] (1.2.2) Calculation The eigenvector of (sorted by eigenvalue from large to small);

[0190] (1.2.3) Reduce to three dimensions and calculate the output array ;

[0191] (2) Determine the minimum number of cluster points and use the density-based clustering algorithm to perform cluster analysis on the array after decorrelation and dimensionality reduction processing;

[0192] This step includes:

[0193] (2.1) Determine the minimum number of cluster points minPts, and based on minPts, use the nearest neighbor algorithm (k-NearestNeighbors, k-NN method) to grid search the optimal neighborhood radius of the output array (corresponding to the inflection point of the k-NN image);

[0194] (2.2) According to minPts and the corresponding optimal neighborhood radius, the DBSCAN method is used to classify the array after decorrelation dimensionality reduction processing;

[0195] Step S440, remove all unclassified arrays, and obtain the average value of all arrays in the same class according to the classification results, and use it as the microcrack occurrence corresponding to the class.

[0196] Step S450, removing the extended portion outside the azimuth domain.

[0197] It should be noted that the classified fracture occurrence data can also be subjected to secondary screening, including but not limited to secondary clustering to merge similar items, secondary data segmentation, and manual screening in combination with other data.

[0198] This embodiment implements data cleaning through two steps of Z-Score standardization and singular value decomposition dimensionality reduction, eliminates the dimensional effect and intrinsic coherence of the occurrence data, improves the reliability of the subsequent clustering algorithm, and uses a density-based clustering algorithm (DBSCAN clustering). When the minimum number of cluster points is artificially limited, the neighborhood radius is automatically determined, and efficient and accurate crack identification is achieved. Compared with the cluster analysis method of the K-Means algorithm, there is no need to predetermine the number of clusters, which simplifies the workload; abnormal point groups that do not meet the minimum number of cluster points can also be effectively identified, reducing the false recognition rate of cracks.

[0199] Embodiment 5: As attached Figure 4 As shown, the embodiment of the present invention is a further optimization of the above embodiment, and also includes performing production statistics on the pre-processed FMI dynamic logging image data, specifically including:

[0200] Step S510, performing production analysis on the pre-processed FMI dynamic logging image data based on two types of quasi-Hough feature transform values, including:

[0201] (1) Determine the deviation coefficient and the standard deviation of the two types of quasi-Hough eigentransform values , and bring it into the segmentation conditions for preprocessing FMI dynamic logging image data segmentation, where the segmentation conditions include:

[0202] If satisfied It is a suspected high-conductivity seam;

[0203] If satisfied It is a suspected high-resistance seam;

[0204] Among them, is a type of quasi-Hough transform value; is the second type of quasi-Hough transform value; is the coefficient of deviation; is the standard deviation of the two types of quasi-Hough eigentransform values;

[0205] (2) For the suspected crack arrays (suspected high-conductivity crack group and suspected high-resistance crack group) after segmentation, Coordinates and corresponding well diameter curves , calculate the depth of the fracture ,inclination and normal azimuth , the calculation formula is as follows:

[0206]

[0207] in, for The corresponding well depth; is the length corresponding to each pixel interval in the vertical direction; The depth of the well is Well diameter at the time; is the vertical pixel coordinate; is the maximum amplitude of the trigonometric function of the inclination angle; is the phase shift value reflecting the crack azimuth.

[0208] Step S520, classifying the suspected fracture groups (suspected high-conductivity fracture groups and suspected high-resistance fracture groups) obtained by the occurrence analysis according to the dip threshold interval, dividing multiple small units of each category according to the azimuth, and counting the weighted fracture number in each small unit to form a weighted superposition azimuth distribution, wherein the weighted fracture number is the sum of the quasi-Hough function values ​​of the suspected fracture groups falling in the small unit;

[0209] Among them, the inclination threshold interval is set according to the actual situation;

[0210] Among them, if the suspected crack group is a suspected high-conductivity crack group, the weighted crack number in each small unit is the sum of the first-class quasi-Hough transform values ​​of the suspected crack group falling in the small unit; if the suspected crack group is a suspected high-resistance crack group, the weighted crack number in each small unit is the sum of the second-class quasi-Hough transform values ​​of the suspected crack group falling in the small unit.

[0211] Step S530, deflect the azimuth distribution by 90° clockwise / counterclockwise and add them together to form a strike-occurrence distribution.

[0212] Step S540, normalizing the weighted superposition azimuth distribution and strike-occurrence distribution, and drawing rose diagrams of azimuth and strike-occurrence distribution.

[0213] Step S550, divide the suspected crack group into a plurality of equally spaced units according to the maximum amplitude, and count the number of corrected cracks in each unit, wherein the calculation formula for the number of corrected cracks is as follows:

[0214] ;

[0215] in, is the inclination effect factor, and the corresponding obtaining process specifically includes:

[0216] (1) Preset the crack thickness range, numerically simulate the number of cracks that can be identified at different maximum amplitudes for each crack thickness within the preset crack thickness range, and obtain a numerical simulation result curve of the inclination effect factor;

[0217] (2) Perform mean processing on the numerical simulation result curve of the inclination effect factor and use the mean as the corrected inclination effect factor .

[0218] Step S560, dividing the inclination angle into a plurality of small units, and counting the sum of the number of corrected cracks falling within each small unit.

[0219] Step S570, normalize the sum of the corrected fracture numbers in each small unit and draw a dip angle distribution histogram.

[0220] In this embodiment, the occurrence analysis is no longer a subsequent product of fracture identification, but can be performed simultaneously with fracture identification, so its reliability is no longer affected by the deviation of the fracture identification process. Since the two types of quasi-Hough transform values ​​obtained from the two steps of data preprocessing and feature transformation in the previous embodiment both show good Gaussian-like distribution characteristics, the occurrence analysis can be directly completed using image threshold segmentation and statistical analysis based on standard deviation in this embodiment. Further, by introducing the inclination effect factor, the inherent influence of different inclinations on the number of fracture identification is calculated and corrected, thereby improving the accuracy of the longitudinal (inclination) occurrence analysis; through weighted superposition calculation, the quasi-Hough function value that physically characterizes the possibility of fractures is included in the occurrence analysis, thereby improving the accuracy of the lateral (azimuth / strike) analysis. Therefore, compared with the existing technical process of first identifying and then doing statistics, this embodiment avoids the recognition deviation caused by the fracture identification algorithm, and the occurrence statistics process is effectively shortened, thereby achieving "fast and accurate" occurrence statistics.

[0221] Embodiment 6: The present invention introduces an example to verify the present invention, which specifically includes:

[0222] 1. Pretreatment

[0223] The image of the X well with a depth range of 3413.2-3423.4 m is intercepted as an example, and contains 6000 and 746 pixels in the vertical and horizontal directions respectively (see Appendix Figure 5 (a), (b)). For the two-dimensional array or image that has been read in (denoted as ,in and They are the total number of horizontal pixels and the total number of vertical pixels respectively; x and are horizontal pixel coordinates and vertical pixel coordinates respectively), and preprocessing is performed according to the steps of Embodiment 3 of the present invention. The preprocessing example is shown in the attached Figure 5 As shown (attached Figure 5 (c) to (h) correspond to each step of preprocessing respectively), and some of the parameter settings are as follows:

[0224] In step S310, the FMI dynamic logging image data to be analyzed is converted into grayscale data, with a range of [0,255], see the attached Figure 5 (c);

[0225] In step S320, the grayscale data is horizontally periodically expanded. The azimuth range extends to Range (i.e. The data is duplication of data);

[0226] In step S320, an initial mask function is constructed. , Data blank part , the rest are 0 (if the image is read as an image containing manual interpretation lines, the manual interpretation part is also 1, see Appendix Figure 5 (d));

[0227] The tapered cosine window function in step S330 In the cosine fraction Pick .

[0228] It should also be noted that Figure 5 The three blue rectangles in the figure are used for local detail comparison, as shown in the attached figure. Figure 6 As shown in the figure, the iterative double spherical harmonics method can be used to interpolate the data when there are few gaps (such as the attached Figure 6 (a) and Figure 6 (b), red dotted box); when there are many data blanks, the continuity of the image at the edge of the data blank is guaranteed, so that subsequent processing will not produce artifact interference at the edge of the data blank (compared with the attached Figure 6 (a) and Figure 6 (c), blue dashed box).

[0229] The completed FMI dynamic logging image data is subjected to de-averaging, edge truncation and high-pass filtering in sequence, which effectively suppresses the background amplitude information of longer wavelengths while retaining the short-wavelength information, effectively enhancing the information amplitude of the fracture (numerical anomaly) area (compared with the attached Figure 6 (a) and Figure 6 (c)).

[0230] 2. Feature Transformation

[0231] After preprocessing, the input data The following steps are performed to calculate two types of quasi-Hough feature transformation functions. In this embodiment, the well depth of 3421m is selected as the center point as an example, and feature transformation is performed according to the steps of embodiment 2 of the present invention. The feature transformation example is shown in the attached figure. Figure 7 As shown, Figure 7 (a) to (e) correspond to each step of feature transformation respectively.

[0232] As attached Figure 7 As can be seen from the black solid box in the figure, a type of quasi-Hough transform requires not only a small Hough transform value, but also a sufficient number of small values ​​along the trigonometric function curve, as shown in the following figure. Figure 7 As can be seen from the black dotted box in the middle, the amplitude of the part with abnormal local Hough transform value and poor coherence with the correction coefficient is compressed, thereby reducing the influence of local singular points.

[0233] As attached Figure 8As shown in the figure, the crack identification results based on Hough feature transform are compared with those based on quasi-Hough feature transform. Among them, (a) is an example of crack identification based on Hough transform. The three sub-images on the left are the results of abnormal point image segmentation based on different Hough domain standard deviations ("3sigma" means 3 Hough domain standard deviations away from the center point). The brown-red scattered points are abnormal points representing the suspected crack occurrence. The green dotted line in the right sub-image is the suspected crack occurrence after image segmentation corresponding to the left middle image; (b) is an example of crack identification based on quasi-Hough feature transform, which is similar to (a).

[0234] Attached Figure 8 In the paper, the results based on the quasi-Hough feature transform are subsequently applied to image segmentation (see Appendix Figure 8 (a), (b) left sub-image) and crack identification (Attachment Figure 8 (a), (b) right sub-figure), the interference is obviously less, and the accuracy of crack judgment is significantly improved (compared with the attached Figure 8 (a) and Figure 8 (b)).

[0235] (III) Occurrence statistics

[0236] In this embodiment, a 5 m well depth section (3435-3440 m) and a 30 m well depth section (3415-3445 m) are selected to perform high conductivity fracture occurrence statistics using the steps of Embodiment 5 of the present invention. The high conductivity fracture occurrence statistics are shown in the attached figure. Fig. 9 shown.

[0237] Attached Fig. 9 (a), (b), and (c) are the production statistics of the 3435-3440 m well depth section; (d), (e), and (f) are the production statistics of the 3415-3445 m well depth section. Fig. 9 In the middle left four (a)(b)(d)(e), the red one is the low angle seam ( ) distribution, green is oblique seams ( ) distribution, black is high angle seam ( ) Occurrence distribution.

[0238] The inclination effect factors used are as follows: Fig.10 As shown, the colored lines represent the number of identifiable cracks simulated according to different crack thicknesses; the black line is the mean (i.e. ). The horizontal axis has been converted to inclination.

[0239] Attached Fig.11 The fracture occurrence distribution map between small areas obtained based on the “first identification and then statistics” process is used for the same well depth section and compared with the attached Fig. 9By comparison, even in rhyolite formations with poor formation integration, many types of faults and many interfering phases, the present invention can still well restore the occurrence characteristics of the dominant high-angle induced fractures and oblique faults in the segment within a small depth interval (compare with the attached Fig. 9 (a)(b)(c) and attached Fig.11 (a)(b)(c)); and in a slightly larger depth range, it can more accurately reflect the overall fracture characteristics and stress state of the formation, such as the overall trend of the induced fracture (compared with the attached Fig. 9 (e) and Fig.11 (e) black rose diagram) and the overall near east-west geostress (see Appendix Fig. 9 (e) Blue double arrow).

[0240] (IV) Crack identification

[0241] In this embodiment, the 3435-3440 m well depth section of Well X is selected, and the steps of Embodiment 4 of the present invention are used to identify high-conductivity fractures. An example of the high-conductivity fracture identification process is shown in the attached figure. Fig.12 The crack identification results are shown in Table 1.

[0242] In this embodiment, when the azimuth angle is extended outside the domain, the azimuth angle is extended outside the domain to ,in Copied from data, Copied from data.

[0243] Attached Fig.12 In the figure, (a) preprocessed image; (b) cross-plot of the suspected fracture; (c) the normalized three-dimensional fracture Intersection diagram; (d) SVD right matrix The intersection diagram, the black "X" is the center point of each cluster; different scattered point colors in (b), (c), and (d) correspond to different clusters; (e) Schematic diagram of the neighborhood radius search of the k-NN method; (f) Crack identification result image, the green dotted line is the crack; (g) Crack identification result after secondary data screening; (h) Identification result provided by the data provider.

[0244] By attaching Fig.12 (f) (g) and attached Fig.12 (h) Comparison shows that the crack identification of the present invention is significantly better than the existing method regardless of whether it has undergone secondary screening or not.

[0245] Appendix 1 Crack identification results (occurrence and uncertainty)

[0246] .

[0247] Embodiment 7: As attached Fig.13As shown, the embodiment of the present invention discloses an automatic identification and occurrence tracking device of electrical imaging logging fractures, comprising:

[0248] A preprocessing unit, for preprocessing the FMI dynamic logging image data to be analyzed;

[0249] Feature transformation unit, including:

[0250] Perform Hough transform on the preprocessed FMI dynamic logging image data to obtain the Hough transform value, and obtain the Hough correction coefficient from the error function constructed by the deviation distribution relative to the Hough mean;

[0251] Based on the Hough transform value and the Hough correction coefficient, two types of quasi-Hough feature transform values ​​are obtained, wherein the two types of quasi-Hough feature transform values ​​are obtained using two types of quasi-Hough feature transform functions as shown below;

[0252]

[0253]

[0254] in, It is a type of quasi-Hough transform value, which is used for the identification and analysis of high-conductivity cracks; It is the second type of quasi-Hough transform value, which is used for the identification and analysis of high-resistance cracks; is the Hough transform value; is the Hough correction coefficient.

[0255] The specific process of the above preprocessing and feature transformation is the same as that of Embodiments 2 and 3 and will not be repeated here.

[0256] Embodiment 8: As attached Fig.14 As shown, the embodiment of the present invention discloses an automatic identification and occurrence tracking device of electrical imaging logging fractures, comprising:

[0257] The preprocessing unit preprocesses the FMI dynamic logging image data to be analyzed.

[0258] Feature transformation unit, including:

[0259] Perform Hough transform on the preprocessed FMI dynamic logging image data to obtain the Hough transform value, and obtain the Hough correction coefficient from the error function constructed by the deviation distribution relative to the Hough mean;

[0260] Based on the Hough transform value and the Hough correction coefficient, two types of quasi-Hough feature transform values ​​are obtained, wherein the two types of quasi-Hough feature transform values ​​are obtained using two types of quasi-Hough feature transform functions as shown below;

[0261]

[0262]

[0263] in, It is a type of quasi-Hough transform value, which is used for the identification and analysis of high-conductivity cracks; It is the second type of quasi-Hough transform value, which is used for the identification and analysis of high-resistance cracks; is the Hough transform value; is the Hough correction coefficient.

[0264] The occurrence analysis unit determines the deviation coefficient and the standard deviation of the two types of quasi-Hough characteristic transform values, and brings them into the segmentation conditions to perform FMI dynamic logging image data segmentation, where the segmentation conditions include:

[0265] If satisfied It is a suspected high-conductivity seam;

[0266] If satisfied It is a suspected high-resistance seam;

[0267] in, is a type of quasi-Hough transform value; is the second type of quasi-Hough transform value; is the coefficient of deviation; is the standard deviation of the two types of quasi-Hough eigentransform values;

[0268] For the segmented suspected crack array, according to its Coordinates and corresponding well diameter curves , calculate the depth of the fracture ,inclination and normal azimuth , the calculation formula is as follows:

[0269]

[0270] in, for The corresponding well depth; is the length corresponding to each pixel interval in the vertical direction; The depth of the well is Well diameter at the time; is the vertical pixel coordinate; is the maximum amplitude of the trigonometric function of the inclination angle; is the phase shift value reflecting the crack azimuth.

[0271] The fracture identification unit performs fracture identification on the pre-processed FMI dynamic logging image data, including:

[0272] Divide the suspected fracture groups obtained by the occurrence analysis into several small depth segments;

[0273] Expand all arrays of suspected fracture groups obtained from the occurrence analysis outside the azimuth domain;

[0274] Sort the array for each depth segment, including:

[0275] (1) For a certain depth segment The array is processed sequentially for Z-Score data standardization and decorrelation dimensionality reduction;

[0276] (2) Determine the minimum number of cluster points and use the density-based clustering algorithm to perform cluster analysis on the array after decorrelation and dimensionality reduction processing;

[0277] Eliminate all unclassified arrays, and obtain the average value of all arrays in the same class based on the classification results, and use it as the microcrack occurrence corresponding to this class;

[0278] Eliminate the extended part outside the azimuth domain;

[0279] The production statistics unit performs production statistics on the pre-processed FMI dynamic logging image data, including:

[0280] The suspected fracture groups obtained by the occurrence analysis are classified according to the dip threshold interval. For each category, multiple small units are divided according to the azimuth, and the weighted fracture number in each small unit is counted to form a weighted superposition azimuth distribution, where the weighted fracture number is the sum of the quasi-Hough function values ​​of the suspected fracture groups falling in the small unit.

[0281] The azimuth distribution is deflected 90° clockwise / counterclockwise and added together to form the strike-occurrence distribution;

[0282] The weighted superposition azimuth distribution and strike-occurrence distribution are normalized to draw rose diagrams of azimuth and strike-occurrence distribution;

[0283] The suspected crack group is divided into multiple equally spaced units according to the maximum amplitude, and the number of corrected cracks in each unit is counted. The calculation formula for the number of corrected cracks is as follows:

[0284] ;

[0285] in, is the inclination effect factor;

[0286] The dip angle is further divided into multiple small units, and the sum of the number of corrected cracks falling in each small unit is counted;

[0287] The sum of the corrected fracture numbers in each small unit is normalized and a dip angle distribution histogram is drawn.

[0288] The specific steps of the above-mentioned occurrence analysis, fracture identification, and occurrence statistics are the same as those of Examples 4 and 5, and will not be repeated here.

[0289] Embodiment 9: The embodiment of the present invention discloses an electronic device, including a processor and a memory, wherein the memory stores a computer program, and the computer program is loaded and executed by the processor to implement a method for automatic identification of electrical imaging logging fractures and production tracking.

[0290] The processor may be a central processing unit (CPU), a general purpose processor, a digital signal processor (DSP), an ASIC, an FPGA or other programmable logic device, a transistor logic device, a hardware component or any combination thereof. It may implement or execute various exemplary logic blocks, modules and circuits described in conjunction with the disclosure of the present invention. It may also be a combination that implements computing functions, such as a combination of one or more microprocessors, a combination of a DSP and a microprocessor, and the like. The memory may include, but is not limited to, various media that can store computer programs, such as a USB flash drive, a read-only memory, a mobile hard disk, a magnetic disk or an optical disk.

[0291] It should be understood by those skilled in the art that the embodiments of the present invention may be provided as methods, systems, or computer program products. Therefore, the present invention may take the form of a complete hardware embodiment, a complete software embodiment, or an embodiment combining software and hardware. Moreover, the present invention may take the form of a computer program product implemented on one or more computer-usable storage media (including but not limited to disk storage, CD-ROM, optical storage, etc.) containing computer-usable program codes. The solutions in the embodiments of the present invention may be implemented in various computer languages, for example, object-oriented programming language Java and interpreted scripting language JavaScript, etc.

[0292] The present invention is described with reference to flowcharts and / or block diagrams of methods, devices (systems), and computer program products according to embodiments of the present invention. It should be understood that each process and / or block in the flowchart and / or block diagram, as well as the combination of processes and / or blocks in the flowchart and / or block diagram, can be implemented by computer program instructions. These computer program instructions can be provided to a processor of a general-purpose computer, a special-purpose computer, an embedded processor, or other programmable data processing device to produce a machine, so that the instructions executed by the processor of the computer or other programmable data processing device generate instructions for implementing the processes in the flowchart and / or block diagram. Figure 1 A process or multiple processes and / or boxes Figure 1 A device that provides the functions specified in a block or multiple blocks.

[0293] These computer program instructions may also be stored in a computer-readable memory capable of directing a computer or other programmable data processing device to operate in a specific manner, so that the instructions stored in the computer-readable memory produce an article of manufacture comprising an instruction device, which implements the process Figure 1 A process or multiple processes and / or boxes Figure 1A function specified in one or more boxes.

[0294] The above content is only a specific implementation mode of the present invention, which has strong adaptability and implementation effect, but the protection scope of the present invention is not limited to this. Any technician familiar with the technical field can easily think of changes or substitutions within the technical scope disclosed by the present invention, which should be covered within the protection scope of the present invention. Therefore, equivalent changes made according to the claims of the present invention are still within the scope covered by the present invention.

Claims

1. A method for automatic identification and occurrence tracking of electrical imaging logging fractures, characterized in that: include: Preprocess the FMI dynamic logging image data to be analyzed; The pre-processed FMI dynamic logging image data is subjected to Hough transform to obtain the Hough transform value, and the Hough correction coefficient is obtained by the error function constructed by the deviation distribution relative to the Hough mean, which specifically includes: Determine the effective pixel number of the pre-processed FMI dynamic logging image data; in, is the number of effective pixels; It is a counting function, which indicates the total number of points that meet the condition; It is a non-number, that is, an invalid pixel; It is the pre-processed FMI dynamic logging image data; x is the horizontal pixel coordinate; is the vertical pixel coordinate; is the maximum amplitude of the trigonometric function of the inclination angle; is the phase shift value reflecting the crack azimuth; is the total number of horizontal pixels; Obtain the Hough transform value and the Hough correction coefficient based on the number of effective pixels; in, is the Hough transform value; is the Hough correction coefficient; is the number of effective pixels; It is a counting function, which indicates the total number of points that meet the condition; It is the pre-processed FMI dynamic logging image data; x is the horizontal pixel coordinate; is the vertical pixel coordinate; is the maximum amplitude of the trigonometric function that characterizes the inclination angle; is the phase shift value reflecting the crack azimuth; is the total number of horizontal pixels; is the mean value of the preprocessed FMI dynamic logging image data; is the error function; is a controllable parameter that controls the shape of the error function; Based on the Hough transform value and the Hough correction coefficient, two types of quasi-Hough feature transform values ​​are obtained, wherein the two types of quasi-Hough feature transform values ​​are obtained using two types of quasi-Hough feature transform functions as shown below; in, It is a type of quasi-Hough transform value, which is used for the identification and analysis of high-conductivity cracks; It is the second type of quasi-Hough transform value, which is used for the identification and analysis of high-resistance cracks; is the Hough transform value; is the Hough correction coefficient.

2. The method for automatic identification and occurrence tracking of electrical imaging logging fractures according to claim 1, characterized in that: The preprocessing of the FMI dynamic logging image data to be analyzed includes: Convert the FMI dynamic logging image data to be analyzed into grayscale data; Perform transverse periodic extension on the grayscale data, and use the improved iterative double spherical harmonic function interpolation algorithm to complete the grayscale data after periodic extension, complete the blank part of the data, and obtain the completed FMI dynamic logging image data; The completed FMI dynamic logging image data is subjected to de-averaging, edge truncation and high-pass filtering in sequence; The complementary data part in the filtered FMI dynamic logging image data is deleted, and the lateral periodic extension part is cropped.

3. The method for automatic identification and occurrence tracking of electrical imaging logging fractures according to claim 1 or 2, characterized in that: It also includes the production analysis of the pre-processed FMI dynamic logging image data based on two types of quasi-Hough characteristic transform values, including: Determine the deviation coefficient and the standard deviation of the two types of quasi-Hough feature transform values, and bring them into the segmentation conditions to perform FMI dynamic logging image data segmentation, where the segmentation conditions include: If satisfied It is a suspected high-conductivity seam; If satisfied It is a suspected high-resistance seam; in, is a type of quasi-Hough transform value; is the second type of quasi-Hough transform value; is the coefficient of deviation; is the standard deviation of the two types of quasi-Hough eigentransform values; For the segmented suspected crack array, according to its Coordinates and corresponding well diameter curves , calculate the depth of the fracture ,inclination and normal azimuth , the calculation formula is as follows: in, for The corresponding well depth; is the length corresponding to each pixel interval in the vertical direction; The depth of the well is Well diameter at the time; is the vertical pixel coordinate; is the maximum amplitude of the trigonometric function of the inclination angle; is the phase shift value reflecting the crack azimuth.

4. The method for automatic identification and occurrence tracking of electrical imaging logging fractures according to claim 3 is characterized in that: It also includes fracture identification of pre-processed FMI dynamic logging image data, including: Divide the suspected fracture groups obtained by the occurrence analysis into several small depth segments; Expand all arrays of suspected fracture groups obtained from the occurrence analysis outside the azimuth domain; Sort the array for each depth segment, including: (1) For a certain depth segment The array is processed sequentially for Z-Score data standardization and decorrelation dimensionality reduction; (2) Determine the minimum number of cluster points and use the density-based clustering algorithm to perform cluster analysis on the array after decorrelation and dimensionality reduction processing; Eliminate all unclassified arrays, and obtain the average value of all arrays in the same class based on the classification results, and use it as the microcrack occurrence corresponding to this class; Eliminate the extended part outside the azimuth domain.

5. The method for automatic identification and occurrence tracking of electrical imaging logging fractures according to claim 3, characterized in that: It also includes the production statistics of the pre-processed FMI dynamic logging image data, including: The suspected fracture groups obtained by the occurrence analysis are classified according to the dip threshold interval. For each classified category, it is divided into multiple small units according to the azimuth. The weighted fracture number in each small unit is counted to form a weighted superimposed azimuth distribution, where the weighted fracture number is the sum of the quasi-Hough function values ​​of the suspected fracture groups falling in the small unit. The azimuth distribution is deflected 90° clockwise / counterclockwise and added together to form the strike-occurrence distribution; The weighted superposition azimuth distribution and strike-occurrence distribution are normalized to draw rose diagrams of azimuth and strike-occurrence distribution; The suspected crack group is divided into multiple equally spaced units according to the maximum amplitude, and the number of corrected cracks in each equally spaced unit is counted. The calculation formula for the number of corrected cracks is as follows: ； in, is the inclination effect factor; Divide the suspected fracture group into multiple small units according to the inclination angle, and count the sum of the corrected fractures in each small unit; The sum of the corrected fracture numbers in each small unit is normalized and a dip angle distribution histogram is drawn.

6. The method for automatic identification and occurrence tracking of electrical imaging logging fractures according to claim 5, characterized in that: Obtaining the tilt effect factor includes: Preset a crack thickness range, numerically simulate the number of cracks that can be identified at different maximum amplitudes for each crack thickness within the preset crack thickness range, and obtain a numerical simulation result curve of the inclination effect factor; The average value of the numerical simulation result curve of the inclination effect factor is processed and the average value is used as the corrected inclination effect factor. .

7. An automatic identification and occurrence tracking device for electrical imaging logging fractures using the method as claimed in any one of claims 1 to 6, characterized in that: include: A preprocessing unit, for preprocessing the FMI dynamic logging image data to be analyzed; Feature transformation unit, including: The pre-processed FMI dynamic logging image data is subjected to Hough transform to obtain the Hough transform value, and the Hough correction coefficient is obtained by the error function constructed by the deviation distribution relative to the Hough mean, which specifically includes: Determine the effective pixel number of the pre-processed FMI dynamic logging image data; in, is the number of effective pixels; It is a counting function, which indicates the total number of points that meet the condition; It is a non-number, that is, an invalid pixel; It is the pre-processed FMI dynamic logging image data; x is the horizontal pixel coordinate; is the vertical pixel coordinate; is the maximum amplitude of the trigonometric function of the inclination angle; is the phase shift value reflecting the crack azimuth; is the total number of horizontal pixels; Obtain the Hough transform value and the Hough correction coefficient based on the number of effective pixels; in, is the Hough transform value; is the Hough correction coefficient; is the number of effective pixels; It is a counting function, which indicates the total number of points that meet the condition; It is the pre-processed FMI dynamic logging image data; x is the horizontal pixel coordinate; is the vertical pixel coordinate; is the maximum amplitude of the trigonometric function that characterizes the inclination angle; is the phase shift value reflecting the crack azimuth; is the total number of horizontal pixels; is the mean value of the preprocessed FMI dynamic logging image data; is the error function; is a controllable parameter that controls the shape of the error function; Based on the Hough transform value and the Hough correction coefficient, two types of quasi-Hough feature transform values ​​are obtained, wherein the two types of quasi-Hough feature transform values ​​are obtained using two types of quasi-Hough feature transform functions as shown below; in, It is a type of quasi-Hough transform value, which is used for the identification and analysis of high-conductivity cracks; It is the second type of quasi-Hough transform value, which is used for the identification and analysis of high-resistance cracks; is the Hough transform value; is the Hough correction coefficient.

8. The electrical imaging logging fracture automatic identification and occurrence tracking device according to claim 7, characterized in that: Also includes: The occurrence analysis unit determines the deviation coefficient and the standard deviation of the two types of quasi-Hough characteristic transform values, and brings them into the segmentation conditions to perform FMI dynamic logging image data segmentation, where the segmentation conditions include: If satisfied It is a suspected high-conductivity seam; If satisfied It is a suspected high-resistance seam; in, is a type of quasi-Hough transform value; is the second type of quasi-Hough transform value; is the coefficient of deviation; is the standard deviation of the two types of quasi-Hough eigentransform values; For the segmented suspected crack array, according to its Coordinates and corresponding well diameter curves , calculate the depth of the fracture ,inclination and normal azimuth , the calculation formula is as follows: in, for The corresponding well depth; is the length corresponding to each pixel interval in the vertical direction; The depth of the well is Well diameter at the time; is the vertical pixel coordinate; is the maximum amplitude of the trigonometric function of the inclination angle; is the phase shift value reflecting the crack azimuth; The fracture identification unit performs fracture identification on the pre-processed FMI dynamic logging image data, including: Divide the suspected fracture groups obtained by the occurrence analysis into several small depth segments; Expand all arrays of suspected fracture groups obtained from the occurrence analysis outside the azimuth domain; Sort the array for each depth segment, including: (1) For a certain depth segment The array is processed sequentially for Z-Score data standardization and decorrelation dimensionality reduction; (2) Determine the minimum number of cluster points and use the density-based clustering algorithm to perform cluster analysis on the array after decorrelation and dimensionality reduction processing; Eliminate all unclassified arrays, and obtain the average value of all arrays in the same class based on the classification results, and use it as the microcrack occurrence corresponding to this class; Eliminate the extended part outside the azimuth domain; The production statistics unit performs production statistics on the pre-processed FMI dynamic logging image data, including: The suspected fracture groups obtained by the occurrence analysis are classified according to the dip threshold interval. For each classified category, it is divided into multiple small units according to the azimuth. The weighted fracture number in each small unit is counted to form a weighted superimposed azimuth distribution, where the weighted fracture number is the sum of the quasi-Hough function values ​​of the suspected fracture groups falling in the small unit. The azimuth distribution is deflected 90° clockwise / counterclockwise and added together to form the strike-occurrence distribution; The weighted superposition azimuth distribution and strike-occurrence distribution are normalized to draw rose diagrams of azimuth and strike-occurrence distribution; The suspected crack group is divided into multiple equally spaced units according to the maximum amplitude, and the number of corrected cracks in each equally spaced unit is counted. The calculation formula for the number of corrected cracks is as follows: ； in, is the inclination effect factor; Divide the suspected fracture group into multiple small units according to the inclination angle, and count the sum of the corrected fractures in each small unit; The sum of the corrected fracture numbers in each small unit is normalized and a dip angle distribution histogram is drawn.

9. An electronic device, characterized in that: It comprises a processor and a memory, wherein a computer program is stored in the memory, and the computer program is loaded and executed by the processor to implement the steps in the method for automatic identification of electrical imaging well logging fractures and occurrence tracking as described in any one of claims 1 to 6.

Images