An em and porosity based digital core pore segmentation method

Abstract

The application discloses a digital core pore segmentation method based on EM and porosity, and comprises the following steps: selecting a preset number of samples meeting preset requirements from a digital core data volume to be segmented; obtaining initial probability parameters of an EM algorithm in a current iteration; obtaining prior probability and probability distribution parameters of pores and backgrounds; then obtaining updated probability parameters in the current iteration; judging whether an iteration stopping condition is met; if not, calculating porosity by using the updated probability distribution parameters in the current iteration; determining whether iteration is accelerated or decelerated based on the porosity in the current iteration, the porosity in the last iteration and measured porosity; then correcting the updated probability distribution parameters in the current iteration to obtain initial probability parameters in the next iteration for continuous iteration; and if yes, obtaining pore segmentation results of each pixel in the digital core data volume by using the updated probability distribution parameters in the current iteration. The application has strong universality, high segmentation precision and high efficiency.

Description

TECHNICAL FIELD

[0001] The application belongs to the field of digital core research, and particularly relates to a digital core pore segmentation method based on EM and porosity. BACKGROUND

[0002] With the continuous development of the oil industry, the demand for oil and gas resources is also increasing, and at present, the exploitation of oil and gas layers has gradually transformed from conventional reservoirs to unconventional reservoirs. Therefore, the role of rock physics research in oil and gas reservoir evaluation is becoming more and more important. Rock physics experiments have limitations in quantitative characterization and microscopic parameter analysis when analyzing physical properties, and at the same time, these experiments take a long time to carry out and are difficult to carry out, and the structure of unconventional reservoirs is more complex, which leads to low core success rate.

[0003] Digital core technology is developed based on CT technology and information processing technology, and in recent years, it has been gradually applied to rock physics. The more popular research direction is to use numerical simulation to calculate and analyze porosity, particle size distribution, permeability and saturation. Digital core technology can extract the internal pore structure of the core without damage and provide visual three-dimensional results, which is of great help to the simulation and analysis of reservoir properties and has important significance for oil reservoir exploration and development. Compared with traditional rock physics experiments, digital core technology has the advantages of fast speed, high efficiency and strong repeatability. Compared with thin section analysis, digital core technology can reconstruct the three-dimensional structure of the core and can more deeply analyze the microscopic characteristics of the rock, which is conducive to establishing the relationship between reservoir parameters and rock physics properties. Therefore, CT technology is becoming an advantageous means for reservoir property analysis, and the research on digital core has important significance.

[0004] Digital core is composed of image sequences collected by CT scanning technology, which are generally grayscale images. The pores correspond to the darker part of the image, i.e. the lower grayscale value. The content of a single digital core is simple, and the pore characteristics are obvious. However, the pore structure is complex and has no fixed shape and characteristics, the pore size range spans microns to nanometers, and the contrast and saturation in different CT data are also different. The extraction and reconstruction of digital core pores depend on the pore segmentation technology, which also affects the subsequent simulation based on pores. Therefore, accurate pore segmentation is the primary problem of digital core research.

[0005] Traditional digital core segmentation methods are mainly based on pixel brightness, morphological features and texture features, etc. The main methods include gray threshold method, multi-threshold method, region growing method, etc. Moreover, scholars have also developed more advanced processing methods, for example, Zhipeng Xu et al. proposed a new determination method for three-dimensional CT scanning images of pores, calculated the geometric and topological parameters of the image pore structure, and used the commonly used mathematical morphological operation for pore detection. Deng et al. generated a custom mask through the crack propagation process, amplified the signal on the intensity histogram to generate cracks for local threshold extraction, and Qiu Maoxin proposed to use image gray and covariance to establish two thresholds to segment the gray value in the CT image to extract pores.

[0006] With the development of information processing technology and computer technology, scholars have also begun to try to use machine learning-based methods to realize the segmentation of rock pores. For example, Cheng Guojian et al. used the mean clustering combined with BP neural network method to extract pores in the rock image. Liu Chun et al. improved the seed algorithm, which automatically segmented pores from the cast thin section image and calculated the geometric features of the pores. Jiang Feng improved the entire pore recognition process based on machine learning of the cast thin section, used transfer learning and multi-angle image augmentation to learn network parameters, and used Rock Net network to automatically extract features.

[0007] The above methods have effectively helped the pore segmentation, recognition, detection and classification of digital cores, but also have different degrees of drawbacks. Mainly reflected in, the traditional method is simple to implement, and is also the commonly used method in current engineering applications, but this method needs to be completed with human interaction, the segmentation accuracy is also low, and is affected by the CT scanning effect, the porosity calculated based on the pore segmentation is often greatly deviated from the experiment. The method based on machine learning can automatically complete the segmentation process, but the calculation amount of a single CT image is large, and at present it is mainly realized for cast thin sections, and there are few pore extraction methods for CT images. At the same time, the gray effect of digital cores obtained by different cores and different scanning processes is different, and there is no unified training sample library at present, and the generalization ability of the method is not high.

[0008] Therefore, how to analyze the characteristics of the digital core itself and design a pore segmentation method with strong universality and high segmentation accuracy and segmentation efficiency is a problem to be solved in the field. SUMMARY

[0009] The purpose of the embodiment of the present application is to provide a digital core pore segmentation method based on EM and porosity, to realize the purpose of improving the universality of the pore segmentation method, and improving the segmentation accuracy and segmentation efficiency. The specific technical solutions are as follows:

[0010] selecting a preset number of samples meeting preset requirements from the digital core data volume to be pore segmented; wherein the preset requirements are that a mean value and a variance of each sample correspond to a mean value and a variance of the digital core data volume respectively;

[0011] obtaining initial probability parameters of the EM algorithm in the current iteration; wherein the probability parameters include prior probabilities of pores and backgrounds, and probability distribution parameters of the pores and the backgrounds;

[0012] based on the preset number of samples and the initial probability parameters in the current iteration, obtaining updated probability parameters in the current iteration; wherein the updated probability parameters in the current iteration include updated prior probabilities and updated probability distribution parameters in the current iteration;

[0013] judging whether the current iteration meets an iteration stopping condition;

[0014] if not, calculating a porosity of the current iteration by using the updated probability distribution parameters in the current iteration; determining a to-be-changed category of iteration speed based on the porosity of the current iteration, a porosity of a previous iteration, and a measured porosity obtained by actually measuring a core sample corresponding to the digital core data volume; and correcting the updated probability distribution parameters in the current iteration according to the to-be-changed category to obtain initial probability parameters in a next iteration to continue iteration; wherein the to-be-changed category includes acceleration and deceleration;

[0015] if yes, obtaining a pore segmentation result of each pixel in the digital core data volume by using the updated probability distribution parameters in the current iteration.

[0016] In an embodiment of the present application, the selecting of the preset number of samples meeting the preset requirements from the digital core data volume to be pore segmented comprises:

[0017] randomly selecting a plurality of digital core slice images from the digital core data volume;

[0018] respectively calculating a mean value and a variance for each of the plurality of digital core slice images;

[0019] taking an average result of all obtained mean values as the mean value of the digital core data volume, and taking an average result of all obtained variances as the variance of the digital core data volume;

[0020] Each time, an unused digital core slice image is selected from the digital core data volume, a part of sub-images is selected from the digital core slice image as a candidate sample image by using a preset size, the mean and variance of the candidate sample image are calculated, the candidate sample image is determined as a sample when the mean and variance of the candidate sample image meet the preset requirements, and the sample selection process is repeated until the number of samples reaches the preset number.

[0021] In an embodiment of the present application, the probability distribution parameters of the pore and the background include:

[0022] The mean and variance of the pore class under the Gaussian probability distribution, and the mean and variance of the background class under the Gaussian probability distribution.

[0023] In an embodiment of the present application, the updated probability parameters in the current iteration are obtained based on the preset number of samples and the initial probability parameters in the current iteration, and the method includes:

[0024] According to the gray value of each pixel in the preset number of samples, the initial probability parameters in the current iteration, the classification soft labels of the pore and the background corresponding to each gray value representing the posterior probability are calculated;

[0025] The prior probabilities of the pore and the background are recalculated by using the classification soft labels of the pore and the background corresponding to each gray value, and the updated prior probabilities in the current iteration are obtained;

[0026] The probability distribution parameters of the pore and the background are recalculated by using the updated prior probabilities in the current iteration, and the updated probability distribution parameters in the current iteration are obtained.

[0027] In an embodiment of the present application, the porosity of the current iteration is calculated by using the updated probability distribution parameters in the current iteration, and the method includes:

[0028] According to the updated probability distribution parameters in the current iteration, the probability that each pixel in the preset number of samples belongs to the pore and the background is calculated;

[0029] The digital core data volume is subjected to pore segmentation by using the probability that each pixel in the preset number of samples belongs to the pore and the background, and a binary result of the pore segmentation is obtained;

[0030] According to the binary result of the pore segmentation, the porosity of the digital core data volume in the current iteration is calculated.

[0031] In an embodiment of the present application, the category to be changed of the iteration speed is determined based on the porosity of the current iteration, the porosity of the last iteration, and the measured porosity obtained by actually measuring the core sample corresponding to the digital core data volume.

[0032] The initial discrimination value of the to-be-changed category in the current iteration is calculated by using the porosity of the current iteration, the porosity of the last iteration, and the measured porosity, and the initial discrimination value of the to-be-changed category is normalized to obtain the discrimination value of the to-be-changed category in the current iteration.

[0033] If the discrimination value of the to-be-changed category is greater than 0.5, it is determined that the to-be-changed category of the iteration speed is acceleration; otherwise, it is determined that the to-be-changed category of the iteration speed is deceleration.

[0034] In an embodiment of the present application, the calculation formula of the initial discrimination value of the to-be-changed category in the current iteration comprises:

[0035]

[0036] wherein, ΔPS i represents the initial discrimination value of the to-be-changed category in the current iteration; PS i represents the porosity of the current iteration; PS i-1 represents the porosity of the last iteration; PS test represents the measured porosity.

[0037] In the process of normalizing the initial discrimination value of the to-be-changed category to obtain the discrimination value of the to-be-changed category in the current iteration, the calculation formula used comprises:

[0038]

[0039] wherein, ΔPS' i represents the discrimination value of the to-be-changed category in the current iteration.

[0040] In an embodiment of the present application, the process of correcting the updated probability distribution parameter in the current iteration according to the to-be-changed category comprises:

[0041] calculating a corresponding correction weight according to the discrimination value of the to-be-changed category in the current iteration;

[0042] based on the correction weight, using the adjustment formula corresponding to the to-be-changed category to correct the mean value of the pore class and the background class in the Gaussian probability distribution in the updated probability distribution parameter in the current iteration.

[0043] In an embodiment of the present application, in the process of calculating a corresponding correction weight according to the discrimination value of the to-be-changed category in the current iteration, the calculation formula used comprises:

[0044]

[0045] wherein, ω represents the correction weight.

[0046] In an embodiment of the present application, the adjustment formula corresponding to the to-be-changed category comprises:

[0047] m cj '=(1+ω)×m cj

[0048] wherein, m cj represents the mean value of the pore class and the background class in the Gaussian probability distribution in the updated probability distribution parameter in the current iteration; j=1, 2 respectively represent the pore class and the background class; m cj ' represents the value after the correction of m cj .

[0049] In the method for digital core pore segmentation based on EM and porosity provided by the embodiment of the present application, first, in order to avoid the influence of high gray value minerals, a preset number of samples meeting preset requirements are selected from the digital core data volume to be segmented based on the gray distribution characteristics, so that the mean value and the variance of the samples correspond to the mean value and the variance of the digital core data volume respectively. At the same time, this kind of sample screening method makes the samples required in the EM iteration process collected from the digital core data volume to be segmented itself, avoiding the problem of sample library data mismatch, and can improve the universality of the algorithm. Then, the EM algorithm is used for unsupervised learning on the preset number of samples. For each iteration, the initial probability parameter is used to obtain the updated probability parameter in the current iteration. When the iteration stopping condition is met, the updated probability distribution parameter in the current iteration is used to obtain the pore segmentation result of each pixel in the digital core data volume. When the iteration stopping condition is not met, the updated probability distribution parameter in the current iteration is used to calculate the porosity in the current iteration, and the measured porosity obtained through experiments is combined to determine whether the iteration is accelerated or decelerated, and the updated probability distribution parameter in the current iteration is corrected accordingly to obtain the initial probability parameter in the next iteration for continuing the iteration. In the present application, the variation of the probability distribution parameter is finely adjusted by combining the measured porosity in the parameter learning process to determine the final probability distribution parameter, which can avoid excessive parameter adjustment and accelerate the iteration process and reduce the operation time. Therefore, the accuracy of the automatic pore segmentation and the convergence speed can be improved when the pores of the digital core are segmented. BRIEF DESCRIPTION OF DRAWINGS

[0050] Figure 1 The flowchart of the method for digital core pore segmentation based on EM and porosity provided by the embodiment of the present application is shown in the figure;

[0051] Figure 2 The mathematical waveform diagram of the correction weight of the embodiment of the present application is shown in the figure;

[0052] Fig. 3 is a result graph of the influence comparison of different sample learning segmentation on the same digital core slice image to be segmented in the experiment of the embodiment of the present application; wherein, Fig. 3(a) is a digital core slice image to be segmented; Fig. 3(b) is the image corresponding to sample 1, sample 2, sample 3 and sample 4 from left to right; Figures 3(c) to 3(g) are the segmentation results of sample 1, sample 2, sample 3 and sample 4 after 1 time, 6 times, 11 times, 16 times and 21 times of iteration of the classic EM algorithm respectively;

[0053] Fig. 4 is a sample screening example given by the experiment of the embodiment of the present application; wherein, Fig. 4(a) is three typical learning samples screened by the method of the embodiment of the present application and meeting the preset requirements; Fig. 4(b) is three unqualified samples screened by the method of the embodiment of the present application;

[0054] Fig. 5 is an effect graph of pore extraction by using the Otus method, the classic EM algorithm and the method of the embodiment of the present application for the same sample in the experiment of the embodiment of the present application; wherein, Fig. 5(a) is the image corresponding to sample a-1, a-2 and a-3. Fig. 5(b) is the effect graph of pore extraction by using the Otus method for sample a-1, a-2 and a-3; Fig. 5(c) is the effect graph of pore extraction by using the classic EM algorithm for sample a-1, a-2 and a-3; Fig. 5(d) is the effect graph of pore extraction by using the method of the embodiment of the present application for sample a-1, a-2 and a-3. DETAILED DESCRIPTION

[0055] The technical solutions in the embodiments of the present application will be clearly and completely described below with reference to the drawings in the embodiments of the present application. Obviously, the described embodiments are only a part of the embodiments of the present application, rather than all the embodiments of the present application. Based on the embodiments in the present application, all other embodiments obtained by those skilled in the art without creative work fall within the scope of protection of the present application.

[0056] How to analyze the characteristics of the digital core itself to design a pore segmentation method with strong universality and capable of comprehensively segmenting the precision and the time-consuming of calculation is the research purpose of the present application. Therefore, the embodiment of the present application provides a digital core pore segmentation method based on EM and porosity.

[0057] It should be noted that the execution subject of the digital core pore segmentation method based on EM and porosity provided by the embodiment of the present application can be a digital core pore segmentation device based on EM and porosity, which can run in an electronic device. Wherein, the electronic device can be a server or a terminal device, of course, but is not limited to this.

[0058] As Figure 1As shown, the method for digital core pore segmentation based on EM and porosity provided by the embodiment of the present application can include the following steps:

[0059] S1, selecting a preset number of samples meeting preset requirements from a digital core data volume to be segmented.

[0060] This step is to collect unlabeled learning samples from the digital core data volume to be segmented as the basis for learning and training. However, it can be understood that the digital core data volume is composed of multiple digital core slice images, and the digital core slice image generally contains two components, i.e., pores with low gray value and minerals with other gray values. Since the types of minerals are complex, the gray values shown in the digital core slice image are also uneven. If the learning samples contain a large number of minerals with high brightness values, the overall gray mean value of the samples will be increased, and therefore the probability distribution model will be biased towards high brightness values, which is likely to cause minerals with relatively low gray values to be divided into pores, and thus it is difficult to effectively segment the pores. Therefore, when selecting samples, minerals with high gray values should be avoided as much as possible.

[0061] In addition, the inventors found in the research process that when the digital core slice image contains a large number of minerals with high gray values, the mean value and variance of the digital core data volume as a whole are large, and the quality of the samples will affect the selection of the probability distribution parameters. Therefore, the embodiment of the present application performs sample screening based on the statistical characteristics of the samples on the basis of randomly selecting samples to ensure the quality of the samples and the learning effect.

[0062] Based on the above reasons, the embodiment of the present application selects a preset number of samples meeting preset requirements from a digital core data volume to be segmented. The preset requirements are that the mean value and variance of each sample correspond to the mean value and variance of the digital core data volume, respectively, which are less than the mean value and variance of the digital core data volume.

[0063] In an optional embodiment, S1 can include the following steps:

[0064] S11, randomly selecting a plurality of digital core slice images from the digital core data volume.

[0065] Specifically, M digital core slice images are randomly selected from the plurality of digital core slice images of the digital core data volume, where M is a natural number greater than 0.

[0066] S12, calculating the mean value and variance of each of the plurality of digital core slice images, respectively.

[0067] Specifically, the mean value and variance of the M randomly selected digital core slice images are calculated, respectively. The mean value and variance of the i-th digital core slice image are represented by and respectively.

[0068] For the calculation process of the mean and the variance, please refer to the related technology, which will not be specifically described here.

[0069] S13, the mean value of all the obtained mean values is taken as the mean value of the digital core data volume, and the mean value of all the obtained variances is taken as the variance of the digital core data volume.

[0070] In this step, the mean value of all the mean values obtained from the M digital core slice images is summed and then averaged to obtain the mean value of the digital core data volume, and the variance of all the variances obtained from the M digital core slice images is summed and then averaged to obtain the variance of the digital core data volume, which is to approximate the gray mean value and the variance of the whole sample by the statistical results of the M digital core slice images. Specifically, the mean value and the variance of the digital core data volume are represented by and respectively, and the calculation method is as follows:

[0071]

[0072] S14, each time, an unused digital core slice image is selected from the digital core data volume, and a part of the sub-image is selected as a candidate sample image by using a preset size, the mean value and the variance of the candidate sample image are calculated, and when the mean value and the variance of the candidate sample image meet the preset requirement, the candidate sample image is determined as a sample, and the sample selection process is repeated until the number of samples reaches the preset number.

[0073] Wherein, the preset size can be represented as n x m, and n x m is less than the size of a single digital core slice image. Therefore, it can be understood that multiple candidate sample images can be obtained from the same digital core slice image of the digital core data volume.

[0074] For this step, each time, an unused, that is, non-repeated digital core slice image can be selected in a certain order from each layer of digital core slice images of the digital core data volume, and then a candidate sample image is selected from the digital core slice image, or each time, a non-repeated digital core slice image can be randomly selected from each layer of digital core slice images of the digital core data volume, and then a candidate sample image is selected therefrom.

[0075] For each selected candidate sample image, the mean value and the variance thereof are calculated, and if the mean value of the candidate sample image is less than the mean value of the digital core data volume and the variance of the candidate sample image is less than the variance of the digital core data volume If yes, it is determined that the candidate sample image can be used as a sample; otherwise, the next candidate sample image is selected for judgment, until k samples are selected, wherein k represents the preset number, and a specific value greater than 0 can be set according to needs. The selected k samples are used as samples in subsequent steps in unsupervised learning.

[0076] In an optional embodiment, Q unused candidate sample images can be randomly selected from each layer of digital core slice images of the digital core data body at a time, as a batch of candidate sample images, and then the mean and variance of each candidate sample image in the batch of candidate sample images are calculated respectively to determine whether it can be reserved as a sample. If the number of samples obtained from the batch of candidate sample images does not meet k, the next batch of Q candidate sample images is randomly selected from the unselected positions of the digital core data body for screening again, until the cumulative number of samples obtained meets k. Wherein Q≥k.

[0077] S2, obtaining initial probability parameters of the EM algorithm in the current iteration.

[0078] When the CT image is subjected to pore segmentation, it is also a problem of whether the probability that a pixel belongs to a pore is greater than the probability that it belongs to a background according to the sample distribution characteristics, so if the probability distribution of the pore and the background can be determined, the segmentation can be realized. Based on the Bayesian theory, scholars introduce the maximum expectation and establish the EM algorithm (Expectation-maximization algorithm, maximum expectation algorithm or expectation maximization algorithm), and calculate the parameters of the probability distribution to determine the specific model through iteration, but the process has certain randomness, and the operation amount in the iteration process is large. Therefore, the embodiment of the application starts from the characteristics of the CT scanned digital core, and proposes a new digital core pore segmentation method from the aspects of effect and convergence.

[0079] In the embodiment of the application, the probability parameters include the prior probabilities of the pore and the background, and the probability distribution parameters of the pore and the background.

[0080] In the initial iteration of the pore segmentation based on the EM algorithm, whether a pixel at a certain gray value belongs to a pore or a background is unknown, so it is necessary to initialize a pore-gray probability distribution carrying parameters and a background-gray probability distribution. In an optional embodiment, both can be realized by using Gaussian distribution.

[0081] For the convenience of description, c1 and c2 represent the pore class and the background class respectively. The prior probabilities of the pore and the background in the probability parameters can be represented by P c1 and P c2The probability distribution parameters of the pores and background include: the mean m of the pore class under a Gaussian probability distribution. c1 variance s c1 And the mean m of the background class under the Gaussian probability distribution. c2 variance s c2 In other words, for k samples, in each iteration of the EM algorithm, the initial probability parameters include six values: m c1 s c1 m c2 s c2 P c1 P c2 .

[0082] For the first iteration, the initial probability parameters can be preset to a value as needed, for example, m. c1 =0.05; s c1 =0.1; m c2 =0.5; s c2 =0.3; P c1 =0.5; P c2 =0.5.

[0083] S3. Based on the preset number of samples and the initial probability parameters in the current iteration, obtain the updated probability parameters in the current iteration.

[0084] The probability parameters updated in the current iteration include the prior probability and the updated probability distribution parameters updated in the current iteration.

[0085] In one optional implementation, S3 may include the following steps:

[0086] S31, based on the gray values ​​of each pixel in the preset number of samples and the initial probability parameters in the current iteration, calculate the classification soft labels for the apertures and background corresponding to each gray value representing the posterior probability.

[0087] Based on Bayesian theory, this invention provides soft labels for classifying pores and background corresponding to gray values ​​in samples according to formulas (2) and (3).

[0088]

[0089] P c2_soft (i)=1-P c1_soft (i) (3)

[0090] Where i is a pixel, x(i) represents the pixel grayscale value, which takes values ​​between [0-255]; PDF() represents the Gaussian probability distribution function; P c1_soft (i) represents the soft label for classifying pores corresponding to grayscale values; Pc2_soft (i) a classification soft label representing the background corresponding to the gray value; P c1_soft (i) and P c2_soft (i) represents the probability that the gray value belongs to the pore or the sample, that is, the posterior probability. j = 1, 2 respectively represent the pore and the background.

[0091] S32, re-computing the prior probability of the pore and the background by using the classification soft label of the pore and the background corresponding to each gray value, to obtain the updated prior probability in the current iteration.

[0092] The calculation formula used in this step is as follows:

[0093]

[0094] Wherein, N represents the number of pixel points contained in the sample, and j = 1, 2 respectively represent the pore and the background.

[0095] S33, re-computing the probability distribution parameter of the pore and the background by using the updated prior probability in the current iteration, to obtain the updated probability distribution parameter in the current iteration.

[0096] The calculation formula used in this step is as follows:

[0097]

[0098] Wherein, m cj , s cj in formula (5) represent the updated probability distribution parameter; j = 1, 2 respectively represent the pore and the background; N cj represents the number of pixel points belonging to the j class in the sample.

[0099] S4, judging whether the current iteration meets the iteration stopping condition.

[0100] Wherein, the iteration stopping condition includes that the current iteration reaches a preset iteration number, or the updated probability distribution parameter in the current iteration achieves convergence.

[0101] The preset iteration number can be set according to needs, such as 30 times, etc. The updated probability distribution parameter in the current iteration achieves convergence means that compared with the updated probability distribution parameter in the last iteration, the difference of the corresponding parameter is less than a certain preset value, such as 0.0001, etc.

[0102] Specifically, if any of the above conditions is met, it means that the current iteration meets the iteration stopping condition, and then S6 is executed; if neither of the two conditions is met, it means that the current iteration does not meet the iteration stopping condition, and then S5 is executed.

[0103] If no, execute S5, using the updated probability distribution parameters in the current iteration, calculate the porosity of the current iteration; based on the porosity of the current iteration, the porosity of the last iteration, and the measured porosity obtained by measuring the core sample corresponding to the digital core data body, determine the category to be changed of the iteration speed; according to the category to be changed, modify the updated probability distribution parameters in the current iteration to obtain the initial probability parameters in the next iteration to continue iteration.

[0104] The category to be changed includes acceleration and deceleration.

[0105] The pore space of the core is the channel of oil and gas accumulation and transportation, and the porosity is an important parameter reflecting the property of the core, which is the ratio of the pore volume to the core volume. The conventional method is to measure the core sample by using a double-chamber porosity tester, and the measurement accuracy is also high. However, this method can only obtain the total porosity of the core, that is, the ratio of the total pore volume to the core volume, and cannot obtain the microscopic pore morphology of the core. The microscopic pore morphology will affect the migration of oil and gas. After the digital core extraction pore is performed on the digital core slice image obtained by CT scanning, the visual pore structure can be established, all the pore morphologies in the core are displayed, and the pore morphology parameters such as porosity are quantitatively calculated. However, the prerequisite is that the accuracy of pore extraction is high, that is, the pore segmentation is accurate. Therefore, the embodiment of the present application considers combining the measured porosity with the pore extraction process of the digital core, and comparing with the experimentally measured porosity to improve the accuracy of segmentation.

[0106] In an optional embodiment, the calculation of the porosity of the current iteration using the updated probability distribution parameters in the current iteration comprises:

[0107] A1, according to the updated probability distribution parameters in the current iteration, the probability that each pixel in the preset number of samples belongs to the pore and the background is calculated.

[0108] Specifically, the updated probability distribution parameters m cj and s cj , j = 1, 2, can be calculated. The probability that each pixel in the k samples belongs to the pore and the background, respectively, such as the probability that a pixel belongs to the pore and the background is 0.2, 0.8, etc. The specific process of this step is understood in combination with related technologies, and will not be described in detail here.

[0109] A2, using the probability that each pixel in the preset number of samples belongs to the pore and the background, the pore of the digital core data body is segmented to obtain the binary result of the pore segmentation.

[0110] For each pixel in the k samples, the larger one of the two probability values corresponding to the class of the pixel belonging to the pore and the background respectively is determined as the final class of the pixel, i.e. the final class of the pixel is pore or background, so as to perform pore segmentation on the digital core data volume to obtain a binary result of the pore segmentation. This process belongs to the prior art and will not be described in detail here. It can be understood that the pore in the obtained binary result of the pore segmentation can be represented by 0 and displayed as black, and the background can be represented by 1 and displayed as white.

[0111] A3, according to the binary result of the pore segmentation, the porosity of the digital core data volume in the current iteration is calculated.

[0112] Specifically, according to the binary result of the pore segmentation, the porosity of the digital core data volume in the current iteration can be calculated by using the number of pixels divided into pores in the digital core data volume / the total number of pixels in the digital core data volume.

[0113] It can be understood that the probability distribution parameters obtained after each iteration are different, and therefore the calculated porosities are different.

[0114] In an optional embodiment, the determination of the to-be-changed class of the iteration speed based on the porosity in the current iteration, the porosity in the last iteration, and the measured porosity of the core sample corresponding to the digital core data volume comprises:

[0115] B1, the to-be-changed class initial discriminant value in the current iteration is calculated by using the porosity in the current iteration, the porosity in the last iteration, and the measured porosity, and the to-be-changed class initial discriminant value is normalized to obtain the to-be-changed class discriminant value in the current iteration.

[0116] When the current iteration is the first iteration, the porosity in the last iteration can be calculated by using the initial probability parameters in the first iteration by the EM algorithm, or a preset value can be given.

[0117] In this step, the calculation formula of the to-be-changed class initial discriminant value in the current iteration comprises:

[0118]

[0119] wherein, ΔPS i represents the to-be-changed class initial discriminant value in the current iteration; PS i represents the porosity in the current iteration; PS i-1 represents the porosity in the last iteration; wherein, i in formula (6) represents the iteration number of the current iteration; PS testrepresents the measured porosity. Wherein, the measured porosity is obtained by using a dual-chamber porosity tester in advance on a core sample corresponding to the digital core data volume to be pore segmented.

[0120] The following briefly describes the process of conceiving the formula (6).

[0121] The inventor of the embodiment of the present application finds that the iteration process is time-consuming when using the EM algorithm, which leads to low efficiency of the overall processing of the digital core and is prone to over-segmentation, affecting the accuracy of the pore segmentation. Therefore, the embodiment of the present application proposes a method based on porosity parameter adjustment in combination with the porosity, to improve the accuracy of pore segmentation and accelerate the convergence speed. Since the probability distribution parameters obtained in each iteration process can reflect different segmentation effects, and further reflect the changes of the segmented porosity. Specifically, in the case where the iteration has not ended, compared with the last iteration, if the change of the probability distribution parameters leads to a small change in the porosity of the slice to be segmented and a large difference from the measured porosity, it means that the iteration convergence process should be accelerated, and the change range of the probability distribution parameters should be increased, otherwise, it means that the iteration convergence process should be slowed down, and the change range of the probability distribution parameters should be reduced. In short, the probability distribution parameters are adjusted based on the porosity in the iteration, so as to adjust the speed of the iteration convergence.

[0122] Specifically, in the process of considering the change of the porosity of the slice to be segmented and the difference from the measured porosity, the relationship among the porosity PS i of the current iteration, the porosity PS i-1 of the last iteration and the measured porosity PS test should be considered. Therefore, the embodiment of the present application first establishes the change relationship of the porosity among the three as formula (7) and (8).

[0123]

[0124]

[0125] Formula (7) represents the porosity change rate of the i-th iteration and the i-1-th iteration; formula (8) represents the approximation degree of the porosity after the i-th iteration to the measured porosity. When the result of formula (7) is greater than the result of formula (8), it means that the porosity change rate is higher than the rate of the porosity tending to the measured value, and the convergence process should be accelerated; when the result of formula (7) is less than the result of formula (8), it means that the porosity change rate is lower than the rate of the porosity tending to the measured value, and the convergence process should be slowed down; of course, if the result of formula (7) is equal to the result of formula (8), it means that the speed of the convergence process can not be changed.

[0126] Therefore, by comparing formula (7) and formula (8), the following formula (9) is obtained In order to prevent the denominator from being 0, 1 is added to the denominator, and finally the initial discriminant value of the category to be changed in the current iteration is obtained as shown in formula (6) i .

[0127] When |ΔPS i |>1, the adjustment speed of the probability distribution parameter is increased; when |ΔPS i |<1, the adjustment speed of the probability distribution parameter is decreased. However, the range of ΔPS i is not the same in the scanning data of different cores, so the initial discriminant value of the category to be changed ΔPS i is normalized to obtain the discriminant value of the category to be changed ΔPS i ' in the current iteration.

[0128] Specifically, in the process of normalizing the initial discriminant value of the category to be changed to obtain the discriminant value of the category to be changed in the current iteration, the calculation formula used includes:

[0129]

[0130] Wherein, ΔPS i ' represents the discriminant value of the category to be changed in the current iteration.

[0131] It can be understood that through formula (9), the numerical range of the discriminant value of the category to be changed ΔPS i ' in the current iteration is normalized to be between [0.1, 0.9].

[0132] B2, if the discriminant value of the category to be changed is greater than 0.5, it is determined that the category to be changed of the iteration speed is acceleration; otherwise, it is determined that the category to be changed of the iteration speed is deceleration.

[0133] Specifically, if the discriminant value of the category to be changed ΔPS i ' in the current iteration is greater than 0.5, it is determined that the category to be changed of the iteration speed is acceleration; if the discriminant value of the category to be changed ΔPS i ' in the current iteration is less than 0.5, it is determined that the category to be changed of the iteration speed is deceleration.

[0134] In an optional embodiment, correcting the updated probability distribution parameter in the current iteration according to the category to be changed includes:

[0135] C1, calculating the corresponding correction weight according to the discriminant value of the category to be changed in the current iteration.

[0136] In an optional embodiment, the calculation formula used in this step includes:

[0137]

[0138] wherein ω represents the correction weight.

[0139] The correction weight of formula (10) is derived from the coordinate transformation of the sigmoid function, and the mathematical waveform is as shown in formula (11). Figure 2 Figure 2 The horizontal axis of the middle part is ΔPS i , and the vertical axis is ω.

[0140] C2, based on the correction weight, uses an adjustment formula corresponding to the to-be-changed category to correct the mean value of the pore class and the background class under the Gaussian probability distribution in the updated probability distribution parameter in the current iteration.

[0141] wherein the adjustment formula corresponding to the to-be-changed category comprises:

[0142] m cj '=(1+ω)×m cj (11)

[0143] wherein m cj represents the mean value of the pore class and the background class under the Gaussian probability distribution in the updated probability distribution parameter in the current iteration; j=1,2 respectively represent the pore class and the background class; m cj ' represents the value after m cj is corrected.

[0144] It can be seen that, after the updated probability distribution parameter in the current iteration is corrected according to the to-be-changed category, the corrected result is taken as the initial probability parameter in the next iteration, and the iteration is continued in S2. The corrected result comprises the variance s cj of the pore class and the background class under the Gaussian probability distribution in the updated probability distribution parameter in the current iteration, the corrected m cj ', and the prior probability P c1 ,P c2 of the pore and the background in the current iteration.

[0145] It can be seen that, based on the porosity change rate and the approximation degree of the segmented porosity and the measured porosity, the adjustment formula corresponding to the to-be-changed category is established, the convergence speed of the EM algorithm can be adjusted, and the segmentation precision is improved.

[0146] If yes, S6 is executed, and the pore segmentation result of each pixel in the digital core data volume is obtained by using the updated probability distribution parameter in the current iteration.

[0147] Specifically, the updated probability distribution parameter m c1 ,s​c1 , m c2 , s c2 After that, the probability of each pixel in each slice of the digital core data volume belonging to the pore and the background can be calculated by using the above parameters. For each pixel, its pixel type is determined by the higher probability value determined, and the pore segmentation is realized. The specific process can be understood in combination with the related prior art, and will not be described in detail here.

[0148] The digital core is generally a gray-scale image, wherein the pore corresponds to a darker part, i.e., a lower gray-scale value, in the image. In a single digital core, the content is simple, and the pore characteristics are relatively obvious. However, the pore structure is complex and has no fixed shape and characteristics, the pore diameter range spans microns to nanometers, the contrast and saturation in different CT data are also different, and at the same time, considering the actual engineering problem, the time consumption of the entire digital core data processing is also a necessary consideration factor.

[0149] In the digital core pore segmentation method based on EM and porosity provided by the embodiment of the application, first, in order to avoid the influence of high gray-scale value minerals, a preset number of samples satisfying a preset requirement are selected from the digital core data volume to be segmented based on the gray-scale distribution characteristics, so that the mean and variance of the samples correspond to the mean and variance of the digital core data volume respectively. At the same time, this kind of sample screening method makes the samples required in the EM iteration process collected from the digital core data volume to be segmented itself, avoiding the problem of sample library data mismatch, and can improve the universality of the algorithm. Then, the EM algorithm is used for unsupervised learning on the preset number of samples. For each iteration, the updated probability parameter in the current iteration is obtained by using the initial probability parameter. When the iteration stopping condition is met, the pore segmentation result of each pixel in the digital core data volume is obtained by using the updated probability distribution parameter in the current iteration. When the iteration stopping condition is not met, the porosity of the current iteration is calculated by using the updated probability distribution parameter in the current iteration, and whether the iteration is accelerated or decelerated is determined in combination with the measured porosity obtained in the experiment, and the updated probability distribution parameter in the current iteration is corrected accordingly to obtain the initial probability parameter in the next iteration for continuing the iteration. In the embodiment of the application, the change of the probability distribution parameter is finely adjusted by combining the measured porosity in the parameter learning process to determine the final probability distribution parameter, which can avoid excessive parameter adjustment and accelerate the iteration process and reduce the operation time consumption. Therefore, the accuracy of the pore automatic segmentation and the convergence speed can be improved when the pores of the digital core are segmented.

[0150] In order to verify the effectiveness of the method of the embodiment of the application, the experimental results are described below.

[0151] The method is used for extracting pores of collected digital cores, and is compared with EM and Otus algorithms. Otus (also known as Otsu method / maximum between-cluster variance method) is a classical algorithm in the field of image segmentation, and is one of the most commonly used pore extraction methods in the processing of digital core engineering at present. In the experiment, the method is compared with two typical algorithms from the perspective of engineering realization.

[0152] 1) Learning sample collection effect on pore extraction

[0153] The method collects learning samples from a digital core data volume to be segmented, and selects samples without minerals with high gray values as learning samples. When the samples are not selected properly, the segmentation result is affected. FIG. 3 is a comparison of the effects of different sample learning segmentation on the same digital core slice image to be segmented. FIG. 3(a) is a digital core slice image to be segmented; FIG. 3(b) is an image corresponding to sample 1, sample 2, sample 3 and sample 4 from left to right. Sample 1 and sample 2 are samples selected by the method and meeting the preset requirements; sample 3 contains minerals with high gray values and has fewer pores; and sample 4 has obvious minerals with medium gray values and has fewer pores.

[0154] Figures 3(c) to 3(g) are segmentation results of the four samples after 1, 6, 11, 16 and 21 iterations of the classical EM algorithm. Figures 3(c) to 3(g) Each row in FIG. 3 corresponds to the segmentation results of sample 1, sample 2, sample 3 and sample 4 under the corresponding iteration number.

[0155] As can be seen from FIG. 3, when the training samples are sample 1 and sample 2, the segmentation results after iteration of the classical EM algorithm are similar. Sample 3 contains minerals with high gray values and has fewer pores, so the pore segmentation is seriously distorted when the classical EM algorithm is iterated for 6 times. According to the classical EM algorithm theory, even if the result is wrong, the iteration operation will continue until the iteration stopping condition is met. Sample 4 has obvious minerals with medium gray values and has fewer pores. From the segmentation effect, after 11 iterations, the segmentation result also starts to appear obvious distortion, and the following iteration operation does not correct the error, but rather aggravates the distortion.

[0156] It can be seen that the above experiment verifies the necessity of selecting samples meeting the preset requirements from the digital core data volume to be segmented by the method. The sample selection of the method can improve the segmentation effect.

[0157] The method of the embodiment of the present application automatically extracts samples from the digital core data volume to be segmented, and filters out truly effective samples, and Figure 4 is a sample filtering example given by the embodiment of the present application; wherein Figure 4(a) is three typical learning samples satisfying preset requirements filtered out by the method of the embodiment of the present application; and Figure 4(b) is three unqualified samples filtered out by the method of the embodiment of the present application.

[0158] 2) Pore extraction result

[0159] For the three samples selected by the embodiment of the present application, Otus, classical EM algorithm and the method of the embodiment of the present application are respectively used for processing, and the pore extraction effects obtained are respectively shown in Figures 5(b), 5(c) and 5(d).

[0160] Among them, Figure 5(a) is three samples, that is, three to-be-segmented slices, which are respectively represented as a-1, a-2 and a-3. Figure 5(b) is an effect diagram of pore extraction of samples a-1, a-2 and a-3 using the Otus method; Figure 5(c) is an effect diagram of pore extraction of samples a-1, a-2 and a-3 using the classical EM algorithm; and Figure 5(d) is an effect diagram of pore extraction of samples a-1, a-2 and a-3 using the method of the embodiment of the present application.

[0161] As can be seen from Figure 5, a-1 of Figure 5(a) contains more large-size pores with obvious features, and the high-brightness mineral distribution is relatively uniform, and the three methods can all extract pores well; a-2 of Figure 5(a) contains obvious pores and unevenly distributed high-brightness minerals, and the Otus method extracts more small pores and noises in the extraction result, compared with the classical EM method and the method of the embodiment of the present application, which can better overcome the influence of noise; a-3 of Figure 5(a) contains fewer pores and complex morphological features, and the brightness and distribution of the mineral are also extremely uneven, and from the extraction effect, the classical EM and the method of the embodiment of the present application can effectively realize the segmentation of the pores, while the Otus method identifies the high-brightness mineral and fails to extract the pores.

[0162] That is, when high gray value minerals appear in the slice and the pores are relatively few, the Otus method will lead to obvious distortion in segmentation, while the classical EM algorithm and the method of the embodiment of the present application can effectively avoid this problem for such samples.

[0163] The three core samples are processed by the above three methods, and the porosity of a single slice is shown in Table 1.

[0164] Table 1. Comparison of porosity

[0165]

[0166] The measured porosity is used as a reference standard in Table 1. The measured porosity is obtained by experimental testing of the whole core sample, and the digital core used for processing is a part obtained by removing the scanning influence on the basis of the scanned data body, so the porosity after porosity extraction of the digital core is inherently deviated from the experimental porosity of the core, but the overall trend is consistent when multiple samples are compared. As can be seen from Table 1, the porosity change rate calculated by the Otus method is too high, and compared with the EM algorithm and the method of the embodiment of the application, the trend after segmentation of the latter two methods is consistent, and can match the test results of the actual porosity.

[0167] The comparison of the iteration process of the EM algorithm and the method of the embodiment of the application can be seen from Tables 2 to 4. Tables 2 to 4 are respectively for three digital core slice images, the iteration times and the corresponding porosity when the EM algorithm and the method of the embodiment of the application are iterated, wherein the decimal value under the iteration times in the table is the porosity obtained by the iteration.

[0168] Table 2 slice1 iteration process comparison

[0169]

[0170] Table 3 slice2 iteration process comparison

[0171]

[0172] Table 4 slice3 iteration process comparison

[0173]

[0174] As can be seen from Tables 2 to 4, the method of the embodiment of the application can converge at a faster speed while ensuring the segmentation effect, and from the processing data of the three slices, compared with the EM algorithm, the method of the embodiment of the application can reduce an average of 2 iterations, and the processing rate is increased by an average of 28%.

[0175] By comparing the processing effect of the method of the embodiment of the application with the classical EM algorithm and the Otus algorithm, it is shown that the sample collection method based on statistical characteristics of the embodiment of the application does not need to establish a sample library in advance, can effectively improve the accuracy of segmentation, avoid the influence of high-light minerals on the sample, and can obtain good segmentation effect for different types of samples. At the same time, compared with the traditional EM algorithm, the iteration efficiency of the embodiment of the application is improved by approximately 28%, which is very beneficial to the engineering application of digital core analysis, and promotes the segmentation of digital core and the remaining research based on digital core.

[0176] The above merely provides the preferred embodiments of the application, and not intended to limit the protection scope of the application. Any modification, equivalent replacement, improvement, etc. made within the spirit and principle of the application shall fall within the protection scope of the application.

Claims

1. A method for EM and porosity based digital core pore segmentation, characterized in that, The method comprises the following steps: selecting a preset number of samples meeting preset requirements from a digital core data volume to be pore segmented, wherein the preset requirements are that the mean and variance of each sample correspond to the mean and variance of the digital core data volume respectively; obtaining initial probability parameters in the current iteration of the EM algorithm, wherein the probability parameters include the prior probabilities of pores and backgrounds and the probability distribution parameters of pores and backgrounds; based on the preset number of samples and the initial probability parameters in the current iteration, obtaining updated probability parameters in the current iteration, wherein the updated probability parameters in the current iteration include updated prior probabilities and updated probability distribution parameters in the current iteration; determining whether the current iteration meets an iteration stopping condition; if not, calculating the porosity of the current iteration by using the updated probability distribution parameters in the current iteration, determining the category to be changed of the iteration speed based on the porosity of the current iteration, the porosity of the last iteration and the measured porosity obtained by actually measuring the core sample corresponding to the digital core data volume, and correcting the updated probability distribution parameters in the current iteration according to the category to be changed to obtain the initial probability parameters in the next iteration to continue the iteration, wherein the category to be changed includes acceleration and deceleration; if yes, obtaining the pore segmentation result of each pixel in the digital core data volume by using the updated probability distribution parameters in the current iteration; wherein, based on the preset number of samples and the initial probability parameters in the current iteration, obtaining the updated probability parameters in the current iteration comprises: calculating the classification soft labels of pores and backgrounds corresponding to each gray value representing the posterior probability according to the gray values of each pixel in the preset number of samples and the initial probability parameters in the current iteration; re-calculating the prior probabilities of pores and backgrounds by using the classification soft labels of pores and backgrounds corresponding to each gray value to obtain the updated prior probabilities in the current iteration; re-calculating the probability distribution parameters of pores and backgrounds by using the updated prior probabilities in the current iteration to obtain the updated probability distribution parameters in the current iteration; wherein, based on the porosity of the current iteration, the porosity of the last iteration and the measured porosity obtained by actually measuring the core sample corresponding to the digital core data volume, determining the category to be changed of the iteration speed comprises: calculating the initial discrimination value of the category to be changed in the current iteration by using the porosity of the current iteration, the porosity of the last iteration and the measured porosity, and performing normalization processing on the initial discrimination value of the category to be changed to obtain the discrimination value of the category to be changed in the current iteration; if the discrimination value of the category to be changed is greater than 0.5, it is determined that the category to be changed of the iteration speed is acceleration; otherwise, it is determined that the category to be changed of the iteration speed is deceleration; wherein, according to the category to be changed, the updated probability distribution parameters in the current iteration are corrected, which comprises: calculating the corresponding correction weight according to the discrimination value of the category to be changed in the current iteration; Based on the correction weight, using the adjustment formula corresponding to the to-be-changed category, the mean value of the pore class and the background class in the updated probability distribution parameter in the current iteration under the Gaussian probability distribution is corrected.

2. The EM and porosity based digital core pore segmentation method of claim 1, wherein, The selecting the preset number of samples meeting the preset requirement from the digital core data volume to be segmented includes: Randomly selecting a plurality of digital core slice images from the digital core data volume; For each of the plurality of digital core slice images, the mean value and the variance are calculated respectively; The average of all obtained mean values is taken as the mean value of the digital core data volume; and the average of all obtained variances is taken as the variance of the digital core data volume; Each time, an unused digital core slice image is selected from the digital core data volume, and a part of sub-images of a preset size is selected therefrom as a candidate sample image, the mean value and the variance of the candidate sample image are calculated, and when the mean value and the variance of the candidate sample image meet the preset requirement, the candidate sample image is determined as a sample, and the sample selection process is repeated until the number of samples reaches the preset number.

3. The EM and porosity based digital core pore segmentation method of claim 1, wherein, The probability distribution parameters of the pore and the background include: The mean value and the variance of the pore class under the Gaussian probability distribution, and the mean value and the variance of the background class under the Gaussian probability distribution.

4. The EM and porosity based digital core pore segmentation method of claim 2 or 3, wherein, The calculation of the porosity of the current iteration using the updated probability distribution parameter in the current iteration includes: According to the updated probability distribution parameter in the current iteration, the probability that each pixel in the preset number of samples belongs to the pore and the background is calculated; Using the probability that each pixel in the preset number of samples belongs to the pore and the background, the digital core data volume is segmented to obtain a binary result of pore segmentation; According to the binary result of pore segmentation, the porosity of the digital core data volume in the current iteration is calculated.

5. The EM and porosity based digital core pore segmentation method of claim 4, wherein, The calculation formula of the initial discrimination value of the to-be-changed category in the current iteration includes: wherein, ΔPS i represents the initial discriminant value of the category to be changed in the current iteration; PS i represents the porosity in the current iteration; PS i-1 represents the porosity in the last iteration; PS test represents the measured porosity; In the process of normalizing the initial discrimination value of the to-be-changed category to obtain the discrimination value of the to-be-changed category in the current iteration, the calculation formula used includes: where ΔPS ' i represents the to-be-changed category discrimination value of the current iteration.

6. The EM and porosity based digital core pore segmentation method of claim 5, wherein, In the process of calculating the correction weight corresponding to the discrimination value of the to-be-changed category in the current iteration, the calculation formula used includes: Wherein, ω represents the correction weight.

7. The EM and porosity based digital core pore segmentation method of claim 6, wherein, The adjustment formula corresponding to the to-be-changed category includes: m cj ' = (1 + ω) x m cj wherein m cj denotes the mean of the pore class and the background class under the Gaussian probability distribution in the updated probability distribution parameter in the current iteration; j = 1, 2 respectively denotes the pore class and the background class; m cj ' denotes the value after the correction of m cj .

Images