CN105957118A - Shale pore imaging method and device - Google Patents

Abstract

The invention provides a shale pore imaging method and a device, wherein the method comprises the steps: X-ray scanning is performed on a shale sample to obtain projection data of X-ray scanning; noise elimination correction processing is performed on the projection data; light intensity of X-ray and transport of intensity equation TIE are adopted to construct a phase interference model of the above projection data; discretization process is performed on the interference model based on a spatial domain, and a spatial domain discretization operator equation is obtained; the processed projection data serves as input data of the spatial domain discretization operator equation, and phase correction projection data is obtained; a filter back projection algorithm is utilized to process the above phase correction projection data, and an image of the shale sample is obtained. According to the invention, a method for solving the spatial domain is adopted, the interference of the phase information is reduced, resistance effect for noise is also presented, and meticulous expression can be performed on nanoscale pore throat and millimeter-scale and micron-scale pores of the shale sample.

Description

Shale pore imaging method and device

Technical Field

The invention relates to the technical field of geophysical exploration, in particular to a shale pore imaging method and device.

Background

In recent years, shale has attracted wide attention due to the fact that rich oil and gas resources are stored, shale oil and gas storage pore structures are complex and various, and a nanoscale pore throat system is taken as a main part, so that the shale oil and gas storage pore structures can be used for realizing micro-nano scale imaging of shale, and the shale oil and gas storage pore structure has great significance for researching geological problems such as oil and gas occurrence states, migration modes and the like. The traditional method is difficult to realize the research of the micro-nano hole cracks of the shale without damage, even if X-ray scanning is used, the problems of insufficient resolution, low signal-to-noise ratio, difficult imaging of a weakly-absorbed object and the like exist, and when the object is away from a detector by a certain distance, the problem of serious phase interference exists. The introduction of third generation synchrotron radiation light will provide hardware possibilities for improving resolution at the light source level, and the following problem is how to process information contained in the phase after the sample is scaled down, and commonly used methods are coaxial phase contrast imaging (Lee,2015), phase contrast imaging based on grating interference (Nesterets,2008), phase contrast imaging based on analyzer (Bravin,2003), and the like.

Through patent retrieval and domestic and foreign literature research, the current method for processing shale pore X-ray projection data by utilizing phase information mostly adopts a frequency domain filtering method, the method eliminates the phase influence through a filtering function, the problems of limited imaging resolution, unstable result and the like caused by easy interference of local values exist, and the research on the shale micro-nano pore cracks is difficult to accurately realize.

Aiming at the problems that the imaging resolution is limited and the result is unstable due to the fact that the shale pore imaging is carried out by adopting a frequency domain filtering method, an effective solution is not provided at present.

Disclosure of Invention

In view of this, an object of the embodiments of the present invention is to provide a method and an apparatus for shale pore imaging, so as to improve imaging resolution, reduce interference of phase information, and reduce relative errors.

In a first aspect, an embodiment of the present invention provides a shale pore imaging method, where the method includes: carrying out X-ray scanning on the shale sample to obtain projection data of the X-ray scanning; carrying out denoising correction processing on the projection data; constructing an interference model of the phase to the projection data by using the light intensity of the X-ray and a light intensity propagation equation TIE; carrying out discretization processing based on a space domain on the interference model to obtain a discretization operator equation of the space domain; taking the projection data after the denoising correction processing as input data of a space domain discretization operator equation to obtain phase correction projection data, wherein the phase correction projection data are the projection data without phase information; and processing the phase correction projection data by using a filtering back projection algorithm to obtain an image of the shale sample.

With reference to the first aspect, an embodiment of the present invention provides a first possible implementation manner of the first aspect, wherein constructing a model of interference of a phase on projection data by using an intensity of X-rays and an intensity propagation equation TIE includes: selecting a phase shift absorption ratio according to the light intensity of the X-ray and the prior information of the shale sample, wherein the prior information comprises the linear absorption coefficient and the absorption edge information of the main components of the shale sample; and making a unity assumption on the phase shift absorption ratio to obtain an interference model, wherein the interference model is a TIE continuous equation taking the projection thickness as an independent variable and is expressed as follows:

$(-\frac{d\mathrm{\δ}}{\mathrm{\μ}}{\▿}^2+1)e^{-\mathrm{\μ}T\left(r\right)}=\frac{I_{\mathrm{\θ}}^d}{I^{in}}$

is the residual light intensity after X-ray transmission through the shale sample, wherein the residual light intensity is recorded by a detector, IⁱⁿIs the incident light intensity of the X-ray, d is the distance between the shale sample and the detector, is the phase factor of the shale sample,for laplacian, μ is the linear absorption coefficient of the shale sample, and t (r) represents the projected thickness of the shale sample.

With reference to the first possible implementation manner of the first aspect, an embodiment of the present invention provides a second possible implementation manner of the first aspect, wherein constructing a model of interference of a phase on projection data by using an intensity of the X-ray and an intensity propagation equation TIE further includes:

order tof＝e^-μT(r)，The interference model is represented by the observation data obtained by recording the residual light intensity of the detector under the noisy condition of a laboratory as follows:

where error is laboratory noise.

With reference to the second possible implementation manner of the first aspect, an embodiment of the present invention provides a third possible implementation manner of the first aspect, where performing discretization processing based on a spatial domain on the interference model to obtain a spatial domain discretization operator equation includes:

to pairSecond order difference operator inCarrying out weighted expression on the space domain by adopting five surrounding points to obtain a space domain discretization expression form of the second-order difference operator:

$\left\{\begin{array}{c}\frac{{\∂}^2{f}_{i,j}}{\∂x^2}\≈\frac{{\mathrm{\δ}}^{2}f}{{\mathrm{\δx}}^2}=\frac{1}{{\left(\mathrm{\Δ}x\right)}^2}(a_1{f}_{i,j}+a_2{f}_{i+1,j}+a_3{f}_{i+2,j}+a_4{f}_{i-1,j}+a_5{f}_{i-2,j})\\ \frac{{\∂}^2{f}_{i,j}}{\∂y^2}\≈\frac{{\mathrm{\δ}}^{2}f}{{\mathrm{\δy}}^2}=\frac{1}{{\left(\mathrm{\Δ}y\right)}^2}(a_1{f}_{i,j}+a_2{f}_{i,j+1}+a_3{f}_{i,j+2}+a_4{f}_{i,j-1}+a_5{f}_{i,j-2})\end{array}\right.$

wherein i, j represent grid points in different directions;

second order difference operator in interference modelAnd setting the two-order difference operator into a space domain discretization expression form of the second-order difference operator to obtain a space domain discretization operator A of the interference model.

With reference to the second possible implementation manner of the first aspect, an embodiment of the present invention provides a fourth possible implementation manner of the first aspect, where taking the projection data after the denoising and correction processing as input data of a spatial domain discretization operator equation, and obtaining phase-corrected projection data includes: establishing a Gihonov regularization model by using a regularization method, and solving the space domain discretization operator equation by using an iteration method, wherein the Gihonov regularization model is expressed as follows:

wherein, min represents the minimum value,f is phase corrected projection data, u is projection data containing phase information, the mathematical notation: -represents the definition, a represents a spatial domain discretization operator derived from the interference model, α is a regularization factor (α)>0)，Is represented by₂And (4) norm.

In a second aspect, an embodiment of the present invention further provides a shale void imaging apparatus, where the apparatus includes: the projection data acquisition module is used for carrying out X-ray scanning on the shale sample to obtain projection data of the X-ray scanning; the de-noising processing module is used for carrying out de-noising correction processing on the projection data; the interference model building module is used for building an interference model of the phase to the projection data by using the light intensity of the X-ray and a light intensity propagation equation TIE; the discretization processing module is used for carrying out discretization processing on the interference model based on a space domain to obtain a space domain discretization operator equation; the equation solving module is used for taking the projection data after the denoising correction processing as input data of a space domain discretization operator equation to obtain phase correction projection data, wherein the phase correction projection data are projection data without phase information; and the result display module is used for processing the phase correction projection data by utilizing a filtering back projection algorithm to obtain an image of the shale sample.

With reference to the second aspect, an embodiment of the present invention provides a first possible implementation manner of the second aspect, where the interference model building module includes: the phase shifting absorption ratio selecting unit is used for selecting a phase shifting absorption ratio according to the light intensity of the X-ray and the prior information of the shale sample, wherein the prior information comprises the linear absorption coefficient and the absorption edge information of the main components of the shale sample; and the interference model characterization unit is used for making a unity assumption on the phase shift absorption ratio to obtain an interference model, wherein the interference model is a TIE continuous equation taking the projection thickness as an independent variable and is expressed as:

$(-\frac{d\mathrm{\δ}}{\mathrm{\μ}}{\▿}^2+1)e^{-\mathrm{\μ}T\left(r\right)}=\frac{I_{\mathrm{\θ}}^d}{I^{in}}$

is the residual light intensity of X-rays after transmission through the shale sample, wherein the residual light intensity is recorded by a detector, IⁱⁿIs the incident light intensity of the X-ray, d is the distance between the shale sample and the detector, is the phase factor of the shale sample,for laplacian, μ is the linear absorption coefficient of the shale sample, and t (r) represents the projected thickness of the shale sample.

With reference to the first possible implementation manner of the second aspect, an embodiment of the present invention provides a second possible implementation manner of the second aspect, where the interference model building module further includes: a noise interference model characterization unit for characterizing the interference model in the noisy condition of the laboratory, and orderf＝e^-μT(r)，The interference model is represented by the observation data obtained by recording the residual light intensity by a detector under the noisy condition in a laboratory as follows:

where error is laboratory noise.

With reference to the second possible implementation manner of the second aspect, an embodiment of the present invention provides a third possible implementation manner of the second aspect, where the discretization processing module includes: a second order difference operator discretization unit for discretizing the second order difference operatorSecond order difference operator inCarrying out weighted expression on the space domain by adopting five surrounding points to obtain a space domain discretization expression form of a second-order difference operator:

$\left\{\begin{array}{c}\frac{{\∂}^2{f}_{i,j}}{\∂x^2}\≈\frac{{\mathrm{\δ}}^{2}f}{{\mathrm{\δx}}^2}=\frac{1}{{\left(\mathrm{\Δ}x\right)}^2}(a_1{f}_{i,j}+a_2{f}_{i+1,j}+a_3{f}_{i+2,j}+a_4{f}_{i-1,j}+a_5{f}_{i-2,j})\\ \frac{{\∂}^2{f}_{i,j}}{\∂y^2}\≈\frac{{\mathrm{\δ}}^{2}f}{{\mathrm{\δy}}^2}=\frac{1}{{\left(\mathrm{\Δ}y\right)}^2}(a_1{f}_{i,j}+a_2{f}_{i,j+1}+a_3{f}_{i,j+2}+a_4{f}_{i,j-1}+a_5{f}_{i,j-2})\end{array}\right.$

wherein i, j represent grid points in different directions;

with reference to the second possible implementation manner of the second aspect, an embodiment of the present invention provides a fourth possible implementation manner of the second aspect, where the equation solving module includes: the regularization unit is used for establishing a Gihonov regularization model by using a regularization method, and solving a spatial domain discretization operator equation by using an iteration method, wherein the Gihonov regularization model is expressed as follows:

wherein, min represents the minimum value,f is phase corrected projection data, u is projection data containing phase information, the mathematical notation: -represents the definition, a represents the spatial domain discretization operator derived from the interference model, α is the regularization factor (α)>0)，Is represented by₂And (4) norm.

According to the shale pore imaging method and device provided by the embodiment of the invention, a shale sample is scanned by utilizing X-rays to obtain projection data scanned by the X-rays, an interference model of phases to the projection data is constructed, the interference model is subjected to discretization processing based on a space domain to obtain a space domain discretization operator equation, the space domain discretization operator equation is solved by a regularization method to obtain phase correction projection data, and an image of a shale pore based on the space domain is obtained according to the phase correction projection data by utilizing a filtering back projection algorithm. The method for imaging the shale pores by adopting the space domain solving method reduces the interference of phase information while improving the imaging resolution, reduces relative errors and also shows the noise resistance effect, thereby being easier to identify the small-scale pores related to the space connectivity of the shale gas reservoir and being capable of finely depicting the nano-scale pore throats and the millimeter-micron-scale pores of the shale samples.

Further, the shale pore imaging method and device provided by the embodiment of the invention establish a solving model based on a spatial domain and solve the shale pore imaging by using an iterative Gihonov regularization method to obtain a relatively stable imaging result under a relatively high resolution. Compared with a general imaging method based on a frequency domain, the method has the advantages that the relative error is smaller in numerical experiments, the imaging result is not influenced by the surrounding pixel values in details, and the like. Meanwhile, compared with a direct solution method of a space domain, the Gihonov iterative regularization method has obvious advantages for resisting noise, and the multi-solution property can be effectively reduced.

In order to make the aforementioned and other objects, features and advantages of the present invention comprehensible, preferred embodiments accompanied with figures are described in detail below.

Drawings

In order to more clearly illustrate the technical solutions of the embodiments of the present invention, the drawings needed to be used in the embodiments will be briefly described below, it should be understood that the following drawings only illustrate some embodiments of the present invention and therefore should not be considered as limiting the scope, and for those skilled in the art, other related drawings can be obtained according to the drawings without inventive efforts.

FIG. 1 is a flow chart illustrating a shale pore imaging method according to an embodiment of the present invention;

fig. 2 shows a structural block diagram of a shale void imaging apparatus provided in an embodiment of the present invention.

Detailed Description

In order to make the objects, technical solutions and advantages of the embodiments of the present invention clearer, the technical solutions in the embodiments of the present invention will be clearly and completely described below with reference to the drawings in the embodiments of the present invention, and it is obvious that the described embodiments are only a part of the embodiments of the present invention, and not all of the embodiments. The components of embodiments of the present invention generally described and illustrated in the figures herein may be arranged and designed in a wide variety of different configurations. Thus, the following detailed description of the embodiments of the present invention, presented in the figures, is not intended to limit the scope of the invention, as claimed, but is merely representative of selected embodiments of the invention. All other embodiments, which can be derived by a person skilled in the art from the embodiments of the present invention without making any creative effort, shall fall within the protection scope of the present invention.

In view of the problems that in the related art, the method of directly utilizing phase imaging through a filter function is mostly focused on a frequency domain, the imaging resolution is limited, and the method is easily interfered by local values, so that the result is unstable, embodiments of the present invention provide a shale pore imaging method and apparatus, which are described below by way of embodiments.

Example 1

Referring to the flow chart of the shale pore imaging method shown in fig. 1, the method comprises the following steps:

step S102, carrying out X-ray scanning on the shale sample to obtain projection data of the X-ray scanning.

The above process of performing X-ray scanning on the shale sample is performed under a laboratory condition, the X-ray light source may be a parallel light source, and the shale sample is scanned to obtain projection data, where the projection data includes phase information of the shale sample.

And step S104, carrying out denoising correction processing on the projection data.

The projection data of the shale sample is the projection data which is obtained by carrying out X-ray scanning on the shale sample under the laboratory condition and contains phase information, noise and errors which are generated under the laboratory condition are considered, the noise comprises bright and dark field noise caused by an X-ray light source, device position offset and some unpredictable noise, and the noise of the bright and dark field and the errors caused by the device position offset can be corrected. Based on this, the step S104 may include: the projection data is subjected to bright dark field noise removal and device position offset correction processing.

And S106, constructing an interference model of the phase to the projection data by using the light intensity of the X-ray and a light intensity propagation equation TIE.

Considering that the pore structure of the shale sample is complex and various, in the specific operation, the phase shift absorption ratio is selected according to the light intensity of the X-ray and the prior information of the shale sample, wherein the prior information comprises the linear absorption coefficient and the absorption edge information of the main components of the shale sample, the phase shift absorption ratio is subjected to the unity assumption, and then an interference model of the phase to the projection data is constructed according to the light intensity of the X-ray and the light intensity propagation equation TIE, and the interference model is expressed as:

$(-\frac{d\mathrm{\δ}}{\mathrm{\μ}}{\▿}^2+1)e^{-\mathrm{\μ}T\left(r\right)}=\frac{I_{\mathrm{\θ}}^d}{I^{in}}$

as described aboveThe residual light intensity after the X-ray penetrates through the shale sample and the light intensity is attenuated can be recorded by a detector, IⁱⁿIs the incident light intensity of the X-ray, d is the distance between the shale sample and the detector, is the phase factor of the shale sample,for laplacian, μ is the linear absorption coefficient of the shale sample, and t (r) represents the projected thickness of the shale sample.

Order tof＝e^-μT(r)，And taking into account the unpredictable noise effects under laboratory conditions, the interference model is expressed as:

where error is laboratory noise, and the noise includes unpredictable noise u after the denoising correction processing in step S104_eCan be obtained by recording the residual intensity of the X-ray by the detector. Based on this, the step S106 may include: selecting a phase shift absorption ratio according to the light intensity of the X-ray and the prior information of the shale sample, wherein the prior information comprises the main components of the shale sampleAnd (3) making a unity assumption on the phase shift absorption ratio by the linear absorption coefficient and the absorption side information, and considering the noise influence under the laboratory condition to obtain an interference model.

And S108, carrying out discretization processing on the interference model based on the space domain to obtain a discretization operator equation of the space domain.

The interference model obtained in step S106 is a continuous equation including a second order difference operatorThis operator represents the second derivative of the function in a continuous equation. Considering that the model needs a computer to solve, discretizing the second-order difference operator, and performing weighted expression on the spatial domain by adopting five surrounding points to obtain a spatial domain discrete expression form of the second-order difference operator:

$\left\{\begin{array}{c}\frac{{\∂}^2{f}_{i,j}}{\∂x^2}\≈\frac{{\mathrm{\δ}}^{2}f}{{\mathrm{\δx}}^2}=\frac{1}{{\left(\mathrm{\Δ}x\right)}^2}(a_1{f}_{i,j}+a_2{f}_{i+1,j}+a_3{f}_{i+2,j}+a_4{f}_{i-1,j}+a_5{f}_{i-2,j})\\ \frac{{\∂}^2{f}_{i,j}}{\∂y^2}\≈\frac{{\mathrm{\δ}}^{2}f}{{\mathrm{\δy}}^2}=\frac{1}{{\left(\mathrm{\Δ}y\right)}^2}(a_1{f}_{i,j}+a_2{f}_{i,j+1}+a_3{f}_{i,j+2}+a_4{f}_{i,j-1}+a_5{f}_{i,j-2})\end{array}\right.$

where i, j denote grid points in different directions, willSecond order difference operator inIs arranged as aboveAnd obtaining a space domain discretization operator A in a space domain discretization expression form, and substituting the space domain discretization operator A into a two-dimensional or three-dimensional space domain discretization operator equation for calculation in specific implementation. The discretization treatment based on the space domain, which is carried out on the interference model, enables the actual scale corresponding to each pixel point to reach 50nm in the imaging process, and further realizes the imaging of the shale pore structure with the nanoscale.

And step S110, taking the projection data after the denoising correction processing as input data of a space domain discretization operator equation to obtain phase correction projection data.

Inputting the projection data subjected to denoising correction into a space domain discretization operator equation, solving the space domain discretization operator equation by using a computer, solving by using a regularization method in consideration of the uncertainty of the solving process, establishing a Gihonov regularization model, and completing the calculating process by using an iterative Gihonov regularization method, wherein the Gihonov regularization model is expressed as:

wherein, min represents the minimum value,f is phase corrected projection data, u is projection data containing phase information, the mathematical notation represents definition, A represents the space domain discretization operator obtained according to the interference model, α is a regularization factor (α)>0)，Is represented by₂And (4) norm. In a specific implementation process, the solving process includes:

(1) inputting regularization parameters α (α)>0) Recording the total number of X-rays k_maxLet k: ═ 1, where the regularization parameterAnd selecting according to the discretization degree of the shale sample and the prior information of the main materials.

(2) If k is>k_maxJumping to the step (5); otherwise, setting i to 0, representing the result of the 0 th iteration in the presence of noise.

(3) Is performed by Gauss elimination i is i + 1; wherein I is an identity matrix, A^*Is a transpose of the matrix a,

and (4) returning to the step (3) until,wherein,₀for terminating the parameters, given by the user, the general settings₀＝10^-3。

(4) And execute itk＝k+1；WhereinRepresenting the k ray projection after i iterations, and the phase recovery result after the iteration with α as a regular parameter, and then returning to the step (2).

(5) And outputting the final data f_α,e。

Wherein the above final numberAccording to f_α,eThe projection data is phase corrected. Based on this, the step S110 may include: and taking the projection data after denoising as input data of a space domain discretization operator equation, establishing a Gihonov regularization model by adopting a regularization method, and finishing a calculation process by an iteration Gihonov regularization method to obtain phase correction projection data.

And step S112, processing the phase correction projection data by using a filtering back projection algorithm to obtain an image of the shale sample.

The method of example 1 above, which scans the shale sample with X-rays to obtain the projection data of the X-ray scan, constructs the interference model of the phase with respect to the projection data, discretizing the interference model based on a space domain to obtain a discretization operator equation of the space domain, solving by a regularization method to obtain phase correction projection data, obtaining an image of the shale pore based on the space domain according to the phase correction projection data by using a filtering back projection algorithm, and carrying out imaging research on the shale pore by adopting a method of solving the space domain, the imaging resolution is improved, the interference of phase information is reduced, the relative error is reduced, the noise resistance effect is also shown, therefore, small-scale pores related to the space connectivity of the shale gas reservoir can be more easily identified, and the nanometer pore throat and the millimeter-micron pores of the shale sample can be more finely depicted.

Example 2

Corresponding to the method provided by the above embodiment, the embodiment of the present invention further provides an apparatus for imaging shale voids, referring to fig. 2, the apparatus includes the following modules:

the projection data obtaining module 202 is configured to perform X-ray scanning on the shale sample to obtain projection data of the X-ray scanning, where the projection data includes phase information of the shale sample.

And the denoising processing module 204 is used for denoising and correcting the projection data, and can remove bright and dark field noise and correct the position offset of the device in the experiment.

And an interference model constructing module 206, configured to construct an interference model of the phase to the projection data by using the light intensity of the X-ray and the light intensity propagation equation TIE, where the interference model is an interference model obtained by considering noise influence under laboratory conditions.

The discretization processing module 208 is configured to perform discretization processing based on a spatial domain on the interference model to obtain a spatial domain discretization operator equation, and when the discretization processing is specifically implemented, the discretization processing is performed on the second-order difference operator in the interference model by considering that the model needs to be solved by a computer, and the spatial domain discretization expression form of the second-order difference operator is obtained by performing weighted expression on five surrounding points in the spatial domain.

An equation solving module 210, configured to use the projection data after the denoising correction processing as input data of a space domain discretization operator equation to obtain phase correction projection data, where the phase correction projection data is projection data without phase information; the spatial domain discretization operator equation of the embodiment is solved by a regularization method, in the specific implementation, a Gihonov regularization model is established, and a calculation process is completed by an iteration Gihonov regularization method to obtain phase correction projection data.

And the result display module 212 is configured to process the phase correction projection data by using a filtered back projection algorithm to obtain an image of the shale sample.

The device scans the shale sample by using X-rays to obtain projection data scanned by the X-rays, constructs an interference model of the phase to the projection data, discretizing the interference model based on a space domain to obtain a discretization operator equation of the space domain, solving by a regularization method to obtain phase correction projection data, obtaining an image of the shale pore based on the space domain according to the phase correction projection data by using a filtering back projection algorithm, and carrying out imaging research on the shale pore by adopting a method of solving the space domain, the imaging resolution is improved, the interference of phase information is reduced, the relative error is reduced, the noise resistance effect is also shown, therefore, small-scale pores related to the space connectivity of the shale gas reservoir can be more easily identified, and the nanometer pore throat and the millimeter-micron pore fracture of the shale sample can be more finely depicted.

In a specific implementation, the interference model building module 206 further includes the following units:

the phase shifting absorption ratio selecting unit is used for selecting a phase shifting absorption ratio according to the light intensity of the X-ray and the prior information of the shale sample, wherein the prior information comprises the linear absorption coefficient and the absorption edge information of the main components of the shale sample;

and the interference model characterization unit is used for making a unity assumption on the phase shift absorption ratio to obtain an interference model, wherein the interference model is a TIE continuous equation taking the projection thickness as an independent variable and is expressed as:

$(-\frac{d\mathrm{\δ}}{\mathrm{\μ}}{\▿}^2+1)e^{-\mathrm{\μ}T\left(r\right)}=\frac{I_{\mathrm{\θ}}^d}{I^{in}}$

is the residual light intensity after X-ray transmission through the shale sample, wherein the residual light intensity is recorded by a detector, IⁱⁿIs the incident light intensity of the X-ray, d is the distance between the shale sample and the detector, is the phase factor of the shale sample,for laplacian, μ is the linear absorption coefficient of the shale sample, and t (r) represents the projected thickness of the shale sample.

A noise interference model characterization unit for characterizing the interference model in the noisy condition of the laboratory, and orderf＝e^-μT(r)，The interference model is represented by the observation data obtained from the residual light intensity recorded by the detector under the noisy condition in the laboratory as:

where error is laboratory noise.

In a specific implementation, the discretization processing module 208 further includes the following units:

a second order difference operator discretization unit for discretizing the second order difference operatorSecond order difference operator inCarrying out weighted expression on the space domain by adopting five surrounding points to obtain a space domain discretization expression form of a second-order difference operator:

$\left\{\begin{array}{c}\frac{{\∂}^2{f}_{i,j}}{\∂x^2}\≈\frac{{\mathrm{\δ}}^{2}f}{{\mathrm{\δx}}^2}=\frac{1}{{\left(\mathrm{\Δ}x\right)}^2}(a_1{f}_{i,j}+a_2{f}_{i+1,j}+a_3{f}_{i+2,j}+a_4{f}_{i-1,j}+a_5{f}_{i-2,j})\\ \frac{{\∂}^2{f}_{i,j}}{\∂y^2}\≈\frac{{\mathrm{\δ}}^{2}f}{{\mathrm{\δy}}^2}=\frac{1}{{\left(\mathrm{\Δ}y\right)}^2}(a_1{f}_{i,j}+a_2{f}_{i,j+1}+a_3{f}_{i,j+2}+a_4{f}_{i,j-1}+a_5{f}_{i,j-2})\end{array}\right.$

where i, j denote grid points in different directions.

In a specific implementation, the equation solving module 210 further includes the following units:

the regularization unit is used for establishing a Gihonov regularization model by using a regularization method, and solving a spatial domain discretization operator equation by using an iteration method, wherein the Gihonov regularization model is expressed as follows:

wherein, min represents the minimum value,f is phase corrected projection data, u is projection data containing phase information, the mathematical notation: -represents the definition, a represents a spatial domain discretization operator derived from an interference model, α is a regularization factor (α)>0)，Is represented by₂And (4) norm.

The device provided by the embodiment of the present invention has the same implementation principle and technical effect as the method embodiments, and for the sake of brief description, reference may be made to the corresponding contents in the method embodiments without reference to the device embodiments.

The shale pore imaging device provided in embodiment 2 obtains a relatively stable and relatively high-resolution imaging result by establishing a spatial domain-based solution model and solving by using an iterative givenov regularization method, and can perform relatively fine depiction on the nano-scale pore throats and the millimeter-micron-scale pores of the shale samples.

The flowchart and block diagrams in the figures illustrate the architecture, functionality, and operation of possible implementations of methods, apparatus, and computer program products according to various embodiments of the present invention. In this regard, each block in the flowchart or block diagrams may represent a module, segment, or portion of code, which comprises one or more executable instructions for implementing the specified logical function(s). It should also be noted that, in some alternative implementations, the functions noted in the block may occur out of the order noted in the figures. For example, two blocks shown in succession may, in fact, be executed substantially concurrently, or the blocks may sometimes be executed in the reverse order, depending upon the functionality involved. It will also be noted that each block of the block diagrams and/or flowchart illustration, and combinations of blocks in the block diagrams and/or flowchart illustration, can be implemented by special purpose hardware-based systems which perform the specified functions or acts, or combinations of special purpose hardware and computer instructions.

In the several embodiments provided in the present application, it should be understood that the disclosed method and apparatus may be implemented in other ways. The above-described embodiments of the apparatus are merely illustrative, and for example, the division of the units is only one logical division, and there may be other divisions when actually implemented, and for example, a plurality of units or components may be combined or integrated into another system, or some features may be omitted, or not executed. In addition, the shown or discussed mutual coupling or direct coupling or communication connection may be an indirect coupling or communication connection of devices or units through some communication interfaces, and may be in an electrical, mechanical or other form.

In addition, functional units in the embodiments of the present invention may be integrated into one processing unit, or each unit may exist alone physically, or two or more units are integrated into one unit.

The functions, if implemented in the form of software functional units and sold or used as a stand-alone product, may be stored in a computer readable storage medium. Based on such understanding, the technical solution of the present invention may be embodied in the form of a software product, which is stored in a storage medium and includes instructions for causing a computer device (which may be a personal computer, a server, or a network device) to execute all or part of the steps of the method according to the embodiments of the present invention. And the aforementioned storage medium includes: a U-disk, a removable hard disk, a Read-Only Memory (ROM), a Random Access Memory (RAM), a magnetic disk or an optical disk, and other various media capable of storing program codes.

The above description is only for the specific embodiments of the present invention, but the scope of the present invention is not limited thereto, and any person skilled in the art can easily conceive of the changes or substitutions within the technical scope of the present invention, and all the changes or substitutions should be covered within the scope of the present invention. Therefore, the protection scope of the present invention shall be subject to the protection scope of the claims.

Claims (10)

1. A method of shale pore imaging, comprising:

carrying out X-ray scanning on the shale sample to obtain projection data of the X-ray scanning;

carrying out denoising correction processing on the projection data;

constructing an interference model of the phase on the projection data by using the light intensity of the X-ray and a light intensity propagation equation TIE;

carrying out discretization processing based on a space domain on the interference model to obtain a discretization operator equation of the space domain;

taking the projection data after the denoising correction processing as input data of the space domain discretization operator equation to obtain phase correction projection data, wherein the phase correction projection data are the projection data without phase information;

and processing the phase correction projection data by using a filtering back projection algorithm to obtain an image of the shale sample.

2. The method of claim 1, wherein constructing a model of the disturbance of phase to the projection data using the intensity of the X-rays and an intensity propagation equation TIE comprises:

selecting a phase shift absorption ratio according to the light intensity of the X-ray and the prior information of the shale sample, wherein the prior information comprises the linear absorption coefficient and the absorption edge information of the main components of the shale sample;

and making a unity assumption on the phase shift absorption ratio to obtain the interference model, wherein the interference model is a TIE continuous equation with the projection thickness as an independent variable and is expressed as follows:

$(-\frac{d\mathrm{\δ}}{\mathrm{\μ}}{\▿}^2+1)e^{-\mathrm{\μ}T\left(r\right)}=\frac{I_{\mathrm{\θ}}^d}{I^{in}}$

the residual light intensity of the X-ray after penetrating the shale sample is recorded by a detector, IⁱⁿIs the incident light intensity of the X-ray, d is the distance between the shale sample and the detector, is the phase factor of the shale sample,for laplacian, μ is the linear absorption coefficient of the shale sample, and t (r) represents the projected thickness of the shale sample.

3. The method of claim 2, wherein constructing a model of the disturbance of phase to the projection data using the intensity of the X-rays and an intensity propagation equation TIE further comprises:

order toThe interference model gives the observation data obtained by the residual light intensity recorded by the detector under the noisy laboratory condition as:

where error is laboratory noise.

4. The method of claim 3, wherein performing a spatial domain-based discretization process on the interference model to obtain a spatial domain discretization operator equation comprises:

to pairSecond order difference operator inIn the skyThe inter-domain adopts the surrounding five points to carry out weighted expression, and the space domain discretization expression form of the second-order difference operator is obtained:

$\left\{\begin{array}{c}\frac{{\∂}^2{f}_{i,j}}{\∂x^2}\≈\frac{{\mathrm{\δ}}^{2}f}{{\mathrm{\δx}}^2}=\frac{1}{{\left(\mathrm{\Δ}x\right)}^2}(a_1{f}_{i,j}+a_2{f}_{i+1,j}+a_3{f}_{i+2,j}+a_4{f}_{i-1,j}+a_5{f}_{i-2,j})\\ \frac{{\∂}^2{f}_{i,j}}{\∂y^2}\≈\frac{{\mathrm{\δ}}^{2}f}{{\mathrm{\δy}}^2}=\frac{1}{{\left(\mathrm{\Δ}y\right)}^2}(a_1{f}_{i,j}+a_2{f}_{i,j+1}+a_3{f}_{i,j+2}+a_4{f}_{i,j-1}+a_5{f}_{i,j-2})\end{array}\right.$

wherein i, j represent grid points in different directions;

second order difference operator in the interference modelAnd setting the discretization expression form of the space domain of the second-order difference operator to obtain a discretization operator A of the space domain of the interference model.

5. The method of claim 3, wherein using the de-noised and corrected projection data as input data to the spatial domain discretization operator equation to derive phase corrected projection data comprises:

establishing a Gihono regularization model by using a regularization method, and solving the space domain discretization operator equation by using an iteration method, wherein the Gihono regularization model is expressed as follows:

wherein, min represents the minimum value,f is phase corrected projection data, u is projection data containing phase information, the mathematical notation represents definition, A represents the space domain discretization operator obtained according to the interference model, α is a regularization factor (α)>0)，Is represented by₂And (4) norm.

6. A shale pore imaging apparatus, comprising:

the projection data acquisition module is used for carrying out X-ray scanning on the shale sample to obtain projection data of the X-ray scanning;

the de-noising processing module is used for carrying out de-noising correction processing on the projection data;

the interference model building module is used for building an interference model of the phase to the projection data by utilizing the light intensity of the X-ray and a light intensity propagation equation TIE;

the discretization processing module is used for carrying out discretization processing on the interference model based on a space domain to obtain a space domain discretization operator equation;

the equation solving module is used for taking the projection data after the denoising correction processing as input data of the space domain discretization operator equation to obtain phase correction projection data, wherein the phase correction projection data are the projection data without phase information;

and the result display module is used for processing the phase correction projection data by utilizing a filtering back projection algorithm to obtain an image of the shale sample.

7. The apparatus of claim 6, wherein the interference model building module comprises:

the phase shifting absorption ratio selecting unit is used for selecting a phase shifting absorption ratio according to the light intensity of the X-ray and the prior information of the shale sample, wherein the prior information comprises a linear absorption coefficient and absorption side information of main components of the shale sample;

an interference model characterization unit, configured to make a unity assumption on the phase shift absorption ratio to obtain the interference model, where the interference model is a TIE continuous equation using projection thickness as an argument, and is expressed as:

$(-\frac{d\mathrm{\δ}}{\mathrm{\μ}}{\▿}^2+1)e^{-\mathrm{\μ}T\left(r\right)}=\frac{I_{\mathrm{\θ}}^d}{I^{in}}$

the residual light intensity of the X-ray after penetrating the shale sample is recorded by a detector, IⁱⁿIs the incident light intensity of the X-ray, d is the distance between the shale sample and the detector, is the phase factor of the shale sample,for laplacian, μ is the linear absorption coefficient of the shale sample, and t (r) represents the projected thickness of the shale sample.

8. The apparatus of claim 7, wherein the interference model building module further comprises:

a noise interference model characterization unit for characterizing the interference model in the noisy condition of the laboratory, and orderThe interference model gives the observation data obtained by the residual light intensity recorded by the detector under the noisy laboratory condition as:

where error is laboratory noise.

9. The apparatus of claim 8, wherein the discretization processing module comprises:

a second order difference operator discretization unit for discretizing the second order difference operatorSecond order difference operator inPerforming weighted expression on the spatial domain by adopting five surrounding points to obtain a spatial domain discretization expression form of the second-order difference operator:

$\left\{\begin{array}{c}\frac{{\∂}^2{f}_{i,j}}{\∂x^2}\≈\frac{{\mathrm{\δ}}^{2}f}{{\mathrm{\δx}}^2}=\frac{1}{{\left(\mathrm{\Δ}x\right)}^2}(a_1{f}_{i,j}+a_2{f}_{i+1,j}+a_3{f}_{i+2,j}+a_4{f}_{i-1,j}+a_5{f}_{i-2,j})\\ \frac{{\∂}^2{f}_{i,j}}{\∂y^2}\≈\frac{{\mathrm{\δ}}^{2}f}{{\mathrm{\δy}}^2}=\frac{1}{{\left(\mathrm{\Δ}y\right)}^2}(a_1{f}_{i,j}+a_2{f}_{i,j+1}+a_3{f}_{i,j+2}+a_4{f}_{i,j-1}+a_5{f}_{i,j-2})\end{array}\right.$

where i, j denote grid points in different directions.

10. The apparatus of claim 8, wherein the equation solving module comprises:

the regularization unit is used for establishing a Gihonov regularization model by using a regularization method, and solving the spatial domain discretization operator equation by using an iteration method, wherein the Gihonov regularization model is expressed as:

wherein, min represents the minimum value,f is phase corrected projection data, u is projection data containing phase information, the mathematical notation represents definition, A represents the space domain discretization operator obtained according to the interference model, α is a regularization factor (α)>0)，Is represented by₂And (4) norm.