Image phase recovery method and device based on X-ray coaxial phase contrast imaging

Abstract

The invention discloses an image phase recovery method and device based on X-ray coaxial phase contrast imaging, and the method comprises the steps: carrying out the first processing of an initial intensity projection image through a convolution block, and obtaining feature data X1; performing iteration processing on the image data X1 for four times through maximum pooling operation and a convolution block to sequentially obtain feature data X2, X3, X4 and X5; processing the feature data through an up-sampling convolution block and a convolution block in sequence according to a set sequence to obtain a phase-shifted image after phase recovery, and performing Fresnel propagation processing on the phase-shifted image to obtain a calculated intensity projection image; performing calculation through a loss function calculation strategy to obtain an L1 loss function value; and iteratively optimizing and updating parameters of the recovery network through an L1 loss function value and a back propagation algorithm, and when the loss function value converges, outputting a final phase shift image and obtaining a reconstructed image through filtering back projection reconstruction. According to the invention, a reconstructed image which clearly shows internal material information can be obtained.

Description

Technical Field

[0001] The present invention relates to the technical field of medical image processing, and in particular to an image phase recovery method and device based on X-ray coaxial phase-contrast imaging. Background Art

[0002] The principle of X-ray CT imaging technology is mainly that the density differences of different parts of an object result in different X-ray absorption characteristics, so that an intensity distribution map of the object's absorption information can be obtained on a detector. However, for samples of low-Z elements such as biological soft tissues, the intensity contrast images obtained on the detector are not sufficient to show the absorption information differences of the parts. When X-rays pass through a substance, in addition to the absorption characteristics of the object to X-rays, the object will produce a phase shift effect in imaging due to the wave nature of X-rays. The amplitude of the X-ray phase change caused by low-Z element samples is at least 3 to 5 orders of magnitude higher than its absorption value of X-rays. Phase-contrast imaging technology can better present the internal structure of low-Z element samples compared with traditional absorption imaging technology. X-ray coaxial phase-contrast CT uses the phase change caused by X-rays passing through a tissue body for imaging. Since X-ray detection can only record intensity information and cannot directly measure phase information, it is necessary to obtain phase change information from intensity data through phase recovery technology to achieve accurate imaging of the internal structure of the sample.

[0003] Currently, traditional X-ray coaxial phase-contrast CT phase recovery methods include analytical methods, iterative methods, and deep learning-based methods. Among them, the Transport of Intensity Equation (TIE) method is a traditional phase recovery analytical method. The TIE method is usually based on small-angle approximation. However, for samples with large optical thickness or strong phase changes, this may make the approximation invalid, resulting in inaccurate recovered phase information. The Contrast Transfer Function (CTF) algorithm is a traditional phase recovery iterative method. The CTF algorithm needs to use phase-contrast images at multiple distances to accurately recover the phase of an object. In practical applications, it is limited by data acquisition and may have the problem of smoothing the phase. When facing different imaging conditions in practical applications, the above methods cannot have good generalization ability. Phase-GAN is a relatively new deep learning-based phase recovery method. Phase-GAN is based on a generative adversarial network and can achieve a relatively high-quality phase recovery effect through adversarial training. However, it highly depends on the quality and quantity of training data, and its generalization ability is limited.

[0004] Therefore, how to invent a new deep learning-based method that does not require the quantity of data, and at the same time has generalization and high-quality phase recovery capabilities has become an urgent problem to be solved. Summary of the Invention

[0005] To this end, the present invention provides an image phase recovery method and device based on X-ray coaxial phase-contrast imaging. For a projected image of a sample with different materials as input, the intensity projection image of Fresnel diffraction is simulated through the Fresnel propagation function to guide the recovery network, and finally a reconstructed image that can clearly display the internal material information is achieved. The present invention can achieve an ideal multi-material image phase recovery effect.

[0006] To achieve the above object, the present invention provides the following technical solutions: An image phase recovery method based on X-ray coaxial phase-contrast imaging, including:

[0007] Input the initial intensity projection image into the recovery network; perform the first processing on the initial intensity projection image through a convolutional block to obtain feature data X 1 ; perform four iterative processes on the image data X 1 through max pooling operation and the convolutional block, and sequentially obtain feature data X 2 、feature data X 3 、feature data X 4 and feature data X 5 ;

[0008] Input the feature data X 1 、the feature data X 2 、the feature data X 3 、the feature data X 4 and the feature data X 5 in a set order, and sequentially process them through an upsampling convolutional block and the convolutional block to obtain a phase-shifted image after phase recovery;

[0009] Perform Fresnel propagation processing on the phase-shifted image after phase recovery to obtain a calculated intensity projection image;

[0010] Calculate the L1 loss function value between the calculated intensity projection image and the initial intensity projection image through a loss function calculation strategy to obtain the L1 loss function value;

[0011] Iteratively optimize and update the parameters of the recovery network through the L1 loss function value and the backpropagation algorithm; when the L1 loss function value reaches a set convergence value, output the final phase-shifted image;

[0012] Reconstruct the final phase-shifted image through a filtered back-projection reconstruction algorithm to obtain a reconstructed image.

[0013] As a preferred solution of the image phase recovery method based on X-ray coaxial phase-contrast imaging, in the process of obtaining the phase-shifted image after phase recovery, the steps of sequentially processing through the upsampling convolutional block and the convolutional block are:

[0014] Input the feature data X 5 into the upsampling convolutional block for processing to obtain a first output result; connect the first output result with the feature data X 4 to obtain a first connection result;

[0015] Input the first connection result into the convolutional block and the upsampling convolutional block for processing to obtain a second output result; connect the second output result with the feature data X 3 to obtain a second connection result;

[0016] Input the second connection result into the convolutional block and the upsampling convolutional block for processing to obtain a third output result; connect the third output result with the feature data X 2 to obtain a third connection result;

[0017] Input the third connection result into the convolutional block and the upsampling convolutional block for processing to obtain a fourth output result; connect the fourth output result with the feature data X 1 to obtain a fourth connection result;

[0018] Input the fourth connection result into the convolutional block for 1×1 convolution processing to obtain the phase-shifted image after phase recovery.

[0019] As a preferred solution of the image phase recovery method based on X-ray in-line phase-contrast imaging, during the process of processing the feature data by the convolutional block and the upsampling convolutional block, perform 3×3 two-dimensional convolution processing on the feature data by the convolutional block and the upsampling convolutional block to obtain two-dimensional convolution data; perform normalization processing on the two-dimensional convolution data by setting a normalization strategy to obtain normalized data; perform processing on the normalized data by using a Leaky ReLu calculation strategy to obtain processed data; perform the next-round processing on the processed data by the convolutional block and the upsampling convolutional block;

[0020] The expression of the set normalization strategy is:

[0021]

[0022] wherein, is the jth feature of the kth sample; μ is the mean of each feature; σ 2 is the variance of each feature; ε` is a very small constant used to prevent the denominator from being zero;

[0023] The expression of the Leaky ReLu calculation strategy is:

[0024]

[0025] Wherein, x is any gray value of the image; α is a small positive number, defaulting to 0.01.

[0026] As a preferred solution of the image phase recovery method based on X-ray in-line phase contrast imaging, in the process of performing Fresnel propagation processing on the phase-shifted image after phase recovery, the calculation formula of Fresnel propagation is:

[0027]

[0028] Wherein, T 0 is the light intensity immediately after the sample, that is, the light intensity when the propagation distance is 0; is the phase information; η is the absorption information; ε is the phase recovery parameter obtained based on absorption-phase duality; e is the natural constant; P D (x, y) is the Fresnel propagation factor; I D (x, y) is the planar light intensity distribution; λ is the light wavelength; represents convolution; D is the propagation distance.

[0029] As a preferred solution of the image phase recovery method based on X-ray in-line phase contrast imaging, in the process of calculating the L1 loss function value between the calculated intensity projection image and the initial intensity projection image through the loss function calculation strategy, the expression of the loss function calculation strategy is:

[0030]

[0031] Wherein, n is the number of projection image pixels; y i is the gray value of the i-th pixel of the input projection image; is the gray value of the i-th pixel of the calculated projection image.

[0032] The present invention also provides an image phase recovery device based on X-ray in-line phase contrast imaging. Based on the above image phase recovery method based on X-ray in-line phase contrast imaging, it includes:

[0033] A feature data acquisition module, configured to input an initial intensity projection image into a recovery network; perform a first processing on the initial intensity projection image through a convolution block to obtain feature data X 1 ; perform four iterative processes on the image data X 1 through a max pooling operation and the convolution block, and sequentially obtain feature data X 2 , feature data X 3 , feature data X 4 and feature data X 5 ;

[0034] The restored phase-shifted image acquisition module is used to process the feature data X 1 、the feature data X 2 、the feature data X 3 、the feature data X 4 and the feature data X 5 in a set order, and sequentially process them through the upsampling convolutional block and the convolutional block to obtain a phase-restored phase-shifted image;

[0035] The calculated intensity projection image acquisition module is used to perform Fresnel propagation processing on the phase-restored phase-shifted image to obtain a calculated intensity projection image;

[0036] The L1 loss function value calculation module is used to calculate the L1 loss function value between the calculated intensity projection image and the initial intensity projection image through a loss function calculation strategy to obtain the L1 loss function value;

[0037] The final phase-shifted image acquisition module is used to iteratively optimize and update the parameters of the restoration network through the L1 loss function value and the backpropagation algorithm; when the L1 loss function value reaches a set convergence value, output the final phase-shifted image;

[0038] The reconstructed image acquisition module is used to reconstruct the final phase-shifted image through a filtered back-projection reconstruction algorithm to obtain a reconstructed image.

[0039] As a preferred solution of the image phase restoration device based on X-ray in-line phase-contrast imaging, in the restored phase-shifted image acquisition module, the sub-module for processing through the upsampling convolutional block and the convolutional block includes:

[0040] The first connection result acquisition sub-module is used to input the feature data X 5 into the upsampling convolutional block for processing to obtain a first output result; connect the first output result with the feature data X 4 to obtain a first connection result;

[0041] The second connection result acquisition sub-module is used to input the first connection result into the convolutional block and the upsampling convolutional block for processing to obtain a second output result; connect the second output result with the feature data X 3 to obtain a second connection result;

[0042] The third connection result acquisition sub-module is used to input the second connection result into the convolutional block and the upsampling convolutional block for processing to obtain a third output result; connect the third output result with the feature data X 2 to obtain a third connection result;

[0043] The fourth connection result acquisition sub-module is used to input the third connection result into the convolutional block and the upsampling convolutional block for processing to obtain a fourth output result; connect the fourth output result with the feature data X 1 for connection processing to obtain a fourth connection result;

[0044] The phase-shifted image acquisition sub-module after phase recovery is used to input the fourth connection result into the convolutional block for 1×1 convolutional processing to obtain the phase-shifted image after phase recovery.

[0045] As a preferred solution of the image phase recovery device based on X-ray in-line phase-contrast imaging, in the recovered phase-shifted image acquisition module, during the process of processing the feature data through the convolutional block and the upsampling convolutional block, 3×3 two-dimensional convolutional processing is performed on the feature data through the convolutional block and the upsampling convolutional block to obtain two-dimensional convolutional data; the two-dimensional convolutional data is normalized by setting a normalization strategy to obtain normalized data; the normalized data is processed through a LeakyReLu calculation strategy to obtain processed data; the processed data is subjected to the next round of processing through the convolutional block and the upsampling convolutional block;

[0046] The expression of the set normalization strategy is:

[0047]

[0048] where, is the j-th feature of the k-th sample; μ is the mean of each feature; σ 2 is the variance of each feature; ε` is a very small constant used to prevent the denominator from being zero;

[0049] The expression of the Leaky ReLu calculation strategy is:

[0050]

[0051] where, x is any gray value of the image; α is a small positive number, defaulting to 0.01.

[0052] As a preferred solution of the image phase recovery device based on X-ray in-line phase-contrast imaging, in the calculation intensity projection image acquisition module, during the process of performing Fresnel propagation processing on the phase-shifted image after phase recovery, the calculation formula of Fresnel propagation is:

[0053]

[0054] where, T 0 is the light intensity immediately after the sample, that is, the light intensity at a propagation distance of 0; is the phase information; η is the absorption information; ε is the phase recovery parameter obtained based on absorption-phase duality; e is the natural constant; P D (x, y) is the Fresnel propagation factor; I D (x, y) is the planar light intensity distribution; λ is the light wavelength; denotes convolution; D is the propagation distance.

[0055] As a preferred solution of the image phase recovery device based on X-ray in-line phase-contrast imaging, in the L1 loss function value calculation module, during the process of calculating the L1 loss function value between the calculated intensity projection image and the initial intensity projection image through the loss function calculation strategy, the expression of the loss function calculation strategy is:

[0056]

[0057] In the formula, n is the number of projection image pixels; y i is the gray value of the i-th pixel of the input projection image; is the gray value of the i-th pixel of the calculated projection image.

[0058] The present invention has the following advantages: The present invention inputs the initial intensity projection image into the recovery network; the initial intensity projection image is first processed through a convolution block to obtain the feature data X 1 ; through the max pooling operation and the convolution block, the image data X 1 is iteratively processed four times to successively obtain the feature data X 2 , the feature data X 3 , the feature data X 4 and the feature data X 5 ; the feature data X 1 , the feature data X 2 , the feature data X 3 , the feature data X 4 and the feature data X 5Process the input successively through the upsampling convolutional block and the convolutional block according to the set order to obtain the phase-shifted image after phase recovery; perform Fresnel propagation processing on the phase-shifted image after phase recovery to obtain the calculated intensity projection image; calculate the L1 loss function value between the calculated intensity projection image and the initial intensity projection image through the loss function calculation strategy to obtain the L1 loss function value; perform iterative optimization and update of the parameters of the recovery network through the L1 loss function value and the backpropagation algorithm; when the L1 loss function value reaches the set convergence value, output the final phase-shifted image; perform reconstruction on the final phase-shifted image through the filtered back-projection reconstruction algorithm to obtain the reconstructed image. The present invention uses the intensity projection images of samples with different materials to train the phase recovery network for X-ray coaxial phase-contrast imaging guided by physical information, and obtains relatively ideal results during the testing process. For the input projection images of samples with different materials, the present invention simulates the intensity projection map of Fresnel diffraction through the Fresnel propagation function to guide the network, and finally realizes a reconstructed image that can clearly display the internal material information, and can achieve relatively ideal phase recovery effects for multi-material images. In view of the characteristics of CT intensity images, the present invention uses the Leaky ReLu function as the activation function of the network model. The Leaky ReLu function ensures that there will be a small amount of gradient propagation even in the case of negative inputs, thereby alleviating gradient disappearance and improving the network training effect; the present invention solves the practical application limitations of the traditional deep learning process that requires a large amount of training data sets and their corresponding label information. The recovery network of the present invention does not require manual annotation, adopts a self-supervised learning method, and embeds the physical imaging mechanism into the loss function of the deep learning network to train and guide the model, so as to accurately simulate and predict the target information, which will greatly save manpower and material resources. The present invention calculates the L1 loss function between the real projection image before phase recovery and the projection image obtained after performing Fresnel propagation on the phase-contrast image after phase recovery, and backpropagates the L1 loss function value to optimize the network model. The gradient of the L1 loss function is relatively constant, and the penalty is linearly related, and will not change sharply due to the increase of errors or extreme values, which helps to enhance the stability of the network model and at the same time can retain the edges and details of the image. Description of the Drawings

[0059] In order to more clearly illustrate the embodiments of the present invention or the technical solutions in the prior art, the following will briefly introduce the drawings required for the description of the embodiments or the prior art. Obviously, the drawings in the following description are only exemplary, and for those of ordinary skill in the art, other implementation drawings can be obtained according to the provided drawings without creative efforts.

[0060] The structures, proportions, sizes, etc. shown in this specification are only used to match the content disclosed in the specification for those familiar with this technology to understand and read, and are not used to limit the implementation conditions of the present invention. Therefore, they do not have substantial technical significance. Any modification of the structure, change in the proportional relationship, or adjustment of the size, without affecting the effects that the present invention can produce and the purposes that can be achieved, should still fall within the scope covered by the technical content disclosed in the present invention.

[0061] Figure 1 It is a schematic flowchart of the image phase recovery method based on X-ray in-line phase-contrast imaging provided in Embodiment 1 of the present invention;

[0062] Figure 2 It is a schematic flowchart of the specific implementation process of the image phase recovery method based on X-ray in-line phase-contrast imaging provided in Embodiment 1 of the present invention;

[0063] Figure 3 It is a schematic flowchart of the convolution block processing process in the image phase recovery method based on X-ray in-line phase-contrast imaging provided in Embodiment 1 of the present invention;

[0064] Figure 4 It is a schematic flowchart of the upsampling convolution block processing process in the image phase recovery method based on X-ray in-line phase-contrast imaging provided in Embodiment 1 of the present invention;

[0065] Figure 5 It is a schematic diagram of the Ground Truth CT reconstruction of the simulation sample in a possible embodiment provided in Embodiment 1 of the present invention;

[0066] Figure 6 It is a schematic diagram of the CT reconstruction of the real sample in a possible embodiment provided in Embodiment 1 of the present invention;

[0067] Figure 7 It is a schematic diagram of the simulation projection image before phase recovery and its reconstruction effect in a possible embodiment provided in Embodiment 1 of the present invention;

[0068] Figure 8 It is a schematic diagram of the simulation projection image after phase recovery and its reconstruction effect in a possible embodiment provided in Embodiment 1 of the present invention;

[0069] Figure 9 It is a schematic diagram of the real projection image before phase recovery and its reconstruction effect in a possible embodiment provided in Embodiment 1 of the present invention;

[0070] Figure 10 It is a schematic diagram of the real projection image after phase recovery and its reconstruction effect in a possible embodiment provided in Embodiment 1 of the present invention;

[0071] Figure 11It is a schematic architecture diagram of the image phase recovery device based on X-ray coaxial phase-contrast imaging provided in Embodiment 2 of the present invention. Detailed implementation manners

[0072] The following specific embodiments illustrate the implementation manners of the present invention. Those skilled in the art can easily understand the other advantages and effects of the present invention from the content disclosed in this specification. Obviously, the described embodiments are part of the embodiments of the present invention, rather than all of the embodiments. All other embodiments obtained by those of ordinary skill in the art based on the embodiments of the present invention without creative efforts shall fall within the protection scope of the present invention.

[0073] Embodiment 1

[0074] Refer to Figure 1 and Figure 2 Embodiment 1 of the present invention provides an image phase recovery method based on X-ray coaxial phase-contrast imaging, including the following steps:

[0075] S1. Input the initial intensity projection image into the recovery network; perform the first processing on the initial intensity projection image through a convolutional block to obtain feature data X 1 ; perform four iterative processes on the image data X 1 through a max pooling operation and the convolutional block to sequentially obtain feature data X 2 , feature data X 3 , feature data X 4 and feature data X 5 ;

[0076] S2. Process the feature data X 1 , the feature data X 2 , the feature data X 3 , the feature data X 4 and the feature data X 5 in a set order through an upsampling convolutional block and the convolutional block to obtain a phase-shifted image after phase recovery;

[0077] S3. Perform Fresnel propagation processing on the phase-shifted image after phase recovery to obtain a calculated intensity projection image;

[0078] S4. Calculate the L1 loss function value between the calculated intensity projection image and the initial intensity projection image through a loss function calculation strategy to obtain the L1 loss function value;

[0079] S5. Iteratively optimize and update the parameters of the recovery network through the L1 loss function value and the backpropagation algorithm; when the L1 loss function value reaches a set convergence value, output the final phase-shifted image;

[0080] S6. Reconstruct the final phase-shifted image through a filtered back-projection reconstruction algorithm to obtain a reconstructed image.

[0081] In this embodiment, in step S1, the initial intensity projection image is input into the restoration network; the initial intensity projection image is first processed through a convolutional block to obtain feature data X 1 ; the image data X 1 is iteratively processed four times through a max-pooling operation and the convolutional block to sequentially obtain feature data X 2 、feature data X 3 、feature data X 4 and feature data X 5 ;

[0082] Specifically, the initial intensity projection images of samples with different materials are input into the restoration network, and after passing through a convolutional block once, feature data X 1 is obtained; the feature data X 1 is sequentially repeated through four max-pooling operations and convolutional block processing to respectively obtain feature data X 2 、feature data X 3 、feature data X 4 and feature data X 5 ;

[0083] Among them, the obtained feature data is downsampled to reduce the spatial dimension of the data.

[0084] In this embodiment, in step S2, the feature data X 1 、the feature data X 2 、the feature data X 3 、the feature data X 4 and the feature data X 5 are processed in a set order through an upsampling convolutional block and the convolutional block in sequence to obtain a phase-shifted image after phase recovery;

[0085] Specifically, the steps of processing through the upsampling convolutional block and the convolutional block in sequence are as follows:

[0086] S21. Input the feature data X 5 into the upsampling convolutional block for processing to obtain a first output result; the first output result is concatenated with the feature data X 4 to obtain a first concatenation result;

[0087] S22. Input the first concatenation result into the convolutional block and the upsampling convolutional block for processing to obtain a second output result; the second output result is concatenated with the feature data X 3Perform connection processing to obtain a second connection result;

[0088] S23. Input the second connection result into the convolutional block and the upsampling convolutional block for processing to obtain a third output result; Connect the third output result with the feature data X 2 Perform connection processing to obtain a third connection result;

[0089] S24. Input the third connection result into the convolutional block and the upsampling convolutional block for processing to obtain a fourth output result; Connect the fourth output result with the feature data X 1 Perform connection processing to obtain a fourth connection result;

[0090] S25. Input the fourth connection result into the convolutional block for 1×1 convolution processing to obtain the phase-shifted image after phase recovery.

[0091] In this embodiment, as Figure 3 and Figure 4 shown, in the process of processing the feature data through the convolutional block and the upsampling convolutional block, perform 3×3 two-dimensional convolution processing on the feature data through the convolutional block and the upsampling convolutional block to obtain two-dimensional convolution data; perform normalization processing on the two-dimensional convolution data by setting a normalization strategy to obtain normalized data; perform processing on the normalized data through a Leaky ReLu calculation strategy to obtain processed data; perform the next-round processing on the processed data through the convolutional block and the upsampling convolutional block;

[0092] Among them, batch normalization is mainly to perform normalization processing on a batch. Given a batch with a size of m and a feature dimension of n, the batch contains multiple samples.

[0093] Calculate the mean μ and variance σ of each feature 2 :

[0094]

[0095] Then perform normalization processing on each feature:

[0096]

[0097] In the formula, is the j-th feature of the k-th sample; μ is the mean of each feature; σ 2 is the variance of each feature; ε` is a very small constant used to prevent the denominator from being zero;

[0098] After the data is normalized, a distribution with a mean of 0 and a variance of unit error is obtained.

[0099] In this embodiment, in order to ensure that the expression ability of the model does not decrease due to the normalization process, two learnable parameters γ j and β j are introduced to scale and shift the normalized values respectively:

[0100]

[0101] In the formula, γ j is a scaling coefficient used to enlarge or reduce the scale of the normalized features; β j is an offset coefficient used to increase or decrease the offset of the normalized features.

[0102] In this embodiment, the expression of the LeakyReLu calculation strategy is:

[0103]

[0104] In the formula, x is any gray value of the image; α is a small positive number, defaulting to 0.01.

[0105] In this embodiment, in step S3, the phase-shifted image after phase recovery is subjected to Fresnel propagation processing to obtain a calculated intensity projection image;

[0106] Specifically, the phase-shifted image after phase recovery is subjected to Fresnel propagation to simulate the diffraction and measurement processes, and a calculated intensity projection image is obtained.

[0107] Among them, when the default incident light intensity is 1, the calculation formula of Fresnel propagation is:

[0108]

[0109] In the formula, T 0 is the light intensity immediately after the sample, that is, the light intensity at a propagation distance of 0; is the phase information; η is the absorption information; ε is the phase recovery parameter obtained based on the absorption-phase duality; e is the natural constant; P D (x, y) is the Fresnel propagation factor; I D (x, y) is the planar light intensity distribution; λ is the light wavelength; denotes convolution; D is the propagation distance.

[0110] In this embodiment, ε is 0.0005 for the simulation sample and 0.0125 for the real sample.

[0111] In this embodiment, in step S4, the L1 loss function value between the calculated intensity projection image and the initial intensity projection image is calculated through the loss function calculation strategy to obtain the L1 loss function value;

[0112] Specifically, calculate the value of the L1 loss function between the initial intensity projection image and the intensity projection image calculated by the network.

[0113] Among them, the expression of the loss function calculation strategy is:

[0114]

[0115] In the formula, n is the number of projection image pixels; y i is the gray value of the i-th pixel of the input projection image; is the gray value of the i-th pixel of the calculated projection image.

[0116] In this embodiment, in step S5, the parameters of the recovery network are iteratively optimized and updated through the L1 loss function value and the backpropagation algorithm; when the L1 loss function value reaches the set convergence value, the final phase-shifted image is output;

[0117] Specifically, the parameters of the recovery network are optimized and updated multiple times through the L1 loss function value and the backpropagation algorithm. When the parameters are updated 1000 times, that is, when the L1 loss value reaches the set convergence value, the final phase-shifted image is output.

[0118] In this embodiment, in step S6, the final phase-shifted image is reconstructed through the filtered back-projection reconstruction algorithm to obtain a reconstructed image.

[0119] Specifically, the final phase-shifted image is reconstructed using the filtered back-projection reconstruction algorithm (FBP algorithm) to obtain a reconstructed image.

[0120] In a possible embodiment, examples of reconstructed images of simulation samples and real samples are provided as follows:

[0121] Perform image reconstruction on the simulation sample and the real sample respectively. The Ground Truth CT reconstruction image of the simulation sample is as Figure 5 shown; the CT reconstruction image of the real sample is as Figure 6 shown.

[0122] Visual evaluation:

[0123] By comparing Figure 7 and Figure 8 , Figure 9 and Figure 10 respectively, it can be seen that the quality of the reconstructed image after phase recovery is better, and it can show the multi-material structure information inside the sample more clearly, and can well distinguish different material information inside the sample.

[0124] Quantitative evaluation:

[0125] Calculate Figure 5 andFigure 8 The SSIM index of the final reconstructed image, with a value of 0.9261, can reflect that the image reconstructed after phase recovery has a high structural similarity with the original image.

[0126] Select Figure 5 the average gray value of the size areas of each detailed part in the simulation samples marked by different colored boxes in

[0127]

[0128] Table 1 Comparison of CNR index values of different simulated materials and air in the images before and after phase recovery

[0129] The higher the CNR value, the better the contrast of the image after phase recovery and the higher the image quality. It can be seen from the quantitative CNR index of the simulation experiment in Table 1 that the CNR value of the image reconstructed after phase recovery is significantly higher than that of the image reconstructed before phase recovery, indicating that the image after phase recovery can more intuitively and realistically display and distinguish different material information.

[0130] Select Figure 6 the average gray value of the size areas of each detailed part in the teeth marked by different colored boxes in

[0131]

[0132] Table 1 Comparison of CNR index values of different real materials and air in the pulp chamber in the images before and after phase recovery

[0133] It can be seen from the quantitative CNR index of the real experiment in Table 2 that the CNR value of the image reconstructed after phase recovery is also significantly higher than that of the image reconstructed before phase recovery.

[0134] In summary, the present invention inputs the initial intensity projection image into the recovery network; the initial intensity projection image is first processed through a convolutional block to obtain feature data X 1 ; through a max-pooling operation and the convolutional block, the image data X 1 is iteratively processed four times to successively obtain feature data X 2 、feature data X 3 、feature data X 4 and feature data X5 ; Process the feature data X 1 , the feature data X 2 , the feature data X 3 , the feature data X 4 and the feature data X 5 in a set order, successively through an upsampling convolutional block and the convolutional block for processing to obtain a phase-shifted image after phase recovery; perform Fresnel propagation processing on the phase-shifted image after phase recovery to obtain a calculated intensity projection image; calculate the L1 loss function value between the calculated intensity projection image and the initial intensity projection image through a loss function calculation strategy to obtain the L1 loss function value; perform iterative optimization and update on the parameters of the recovery network through the L1 loss function value and the backpropagation algorithm; when the L1 loss function value reaches a set convergence value, output the final phase-shifted image; perform reconstruction on the final phase-shifted image through a filtered back-projection reconstruction algorithm to obtain a reconstructed image. The present invention uses intensity projection images of samples with different materials to train a phase recovery network for X-ray coaxial phase-contrast imaging guided by physical information, and obtains relatively ideal results during the testing process. For the input projection images of samples with different materials, the present invention simulates the intensity projection diagram of Fresnel diffraction through a Fresnel propagation function to guide the network, and finally realizes a reconstructed image that can clearly display the internal material information, and can achieve relatively ideal phase recovery effects for multi-material images. In view of the characteristics of CT intensity images, the present invention uses the Leaky ReLu function as the activation function of the network model. The Leaky ReLu function ensures that there will be a small amount of gradient propagation even in the case of negative input, thereby alleviating gradient disappearance and improving the network training effect; the present invention solves the practical application limitations of the traditional deep learning process that requires a large amount of training data sets and their corresponding label information. The recovery network of the present invention does not require manual annotation, adopts a self-supervised learning method, embeds the physical imaging mechanism into the loss function of the deep learning network, and uses this to guide the training of the model, so as to accurately simulate and predict the target information, which will greatly save manpower and material resources. The present invention calculates the L1 loss function between the true projection image before phase recovery and the projection image obtained after performing Fresnel propagation on the phase-contrast image after phase recovery, and the L1 loss function value is backpropagated to optimize the network model. The gradient of the L1 loss function is relatively constant, and the penalty is linearly related, and it will not change sharply due to the increase of errors or extreme values, which helps to enhance the stability of the network model and can retain the edges and details of the image.

[0135] Note that the method of the embodiments of the present disclosure can be executed by a single device, such as a computer or a server. The method of this embodiment can also be applied to a distributed scenario, and multiple devices cooperate with each other to complete it. In such a distributed scenario, one of the multiple devices can only execute one or more steps of the method of the embodiments of the present disclosure, and these multiple devices will interact with each other to complete the described method.

[0136] Note that some embodiments of the present disclosure have been described above. Other embodiments are within the scope of the appended claims. In some cases, the actions or steps recited in the claims can be executed in a different order than in the above embodiments and still achieve the desired results. Additionally, the processes depicted in the drawings do not necessarily require the specific order or sequential order shown to achieve the desired results. In certain embodiments, multitasking and parallel processing are also possible or may be advantageous.

[0137] Embodiment 2

[0138] See Figure 11 , Embodiment 2 of the present invention also provides an image phase recovery device based on X-ray coaxial phase-contrast imaging, including:

[0139] A feature data acquisition module 001, configured to input an initial intensity projection image into a recovery network; perform a first processing on the initial intensity projection image through a convolutional block to obtain feature data X 1 ; perform four iterative processes on the image data X 1 through a max-pooling operation and the convolutional block, and sequentially obtain feature data X 2 , feature data X 3 , feature data X 4 and feature data X 5 ;

[0140] A recovered phase shift image acquisition module 002, configured to sequentially process the feature data X 1 , the feature data X 2 , the feature data X 3 , the feature data X 4 and the feature data X 5 in a set order through an upsampling convolutional block and the convolutional block to obtain a phase shift image after phase recovery;

[0141] A calculated intensity projection image acquisition module 003, configured to perform Fresnel propagation processing on the phase shift image after phase recovery to obtain a calculated intensity projection image;

[0142] The L1 loss function value calculation module 004 is configured to calculate the L1 loss function value between the calculated intensity projection image and the initial intensity projection image through a loss function, so as to obtain the L1 loss function value.

[0143] The final phase shift image acquisition module 005 is configured to iteratively optimize and update the parameters of the recovery network through the L1 loss function value and the backpropagation algorithm; when the L1 loss function value reaches a set convergence value, output the final phase shift image.

[0144] The reconstructed image acquisition module 006 is configured to reconstruct the final phase shift image through a filtered back-projection reconstruction algorithm to obtain a reconstructed image.

[0145] In this embodiment, in the recovered phase shift image acquisition module 002, the sub-module for processing through the upsampling convolutional block and the convolutional block includes:

[0146] The first connection result acquisition sub-module 021 is configured to input the feature data X 5 into the upsampling convolutional block for processing to obtain a first output result; connect the first output result with the feature data X 4 to obtain a first connection result.

[0147] The second connection result acquisition sub-module 022 is configured to input the first connection result into the convolutional block and the upsampling convolutional block for processing to obtain a second output result; connect the second output result with the feature data X 3 to obtain a second connection result.

[0148] The third connection result acquisition sub-module 023 is configured to input the second connection result into the convolutional block and the upsampling convolutional block for processing to obtain a third output result; connect the third output result with the feature data X 2 to obtain a third connection result.

[0149] The fourth connection result acquisition sub-module 024 is configured to input the third connection result into the convolutional block and the upsampling convolutional block for processing to obtain a fourth output result; connect the fourth output result with the feature data X 1 to obtain a fourth connection result.

[0150] The phase shift image acquisition sub-module 025 after phase recovery is configured to input the fourth connection result into the convolutional block for 1×1 convolutional processing to obtain the phase shift image after phase recovery.

[0151] In this embodiment, in the restored phase-shifted image acquisition module 002, during the process of processing the feature data through the convolutional block and the upsampling convolutional block, the feature data is subjected to 3×3 two-dimensional convolutional processing through the convolutional block and the upsampling convolutional block to obtain two-dimensional convolutional data; the two-dimensional convolutional data is normalized through a set normalization strategy to obtain normalized data; the normalized data is processed through a Leaky ReLu calculation strategy to obtain processed data; the processed data is subjected to the next round of processing through the convolutional block and the upsampling convolutional block.

[0152] The expression of the set normalization strategy is:

[0153]

[0154] In the formula, is the jth feature of the kth sample; μ is the mean of each feature; σ 2 is the variance of each feature; ε` is a very small constant used to prevent the denominator from being zero;

[0155] The expression of the Leaky ReLu calculation strategy is:

[0156]

[0157] In the formula, x is any gray value of the image; α is a small positive number, defaulting to 0.01.

[0158] In this embodiment, in the calculation intensity projection image acquisition module 003, during the process of performing Fresnel propagation processing on the phase-restored phase-shifted image, the calculation formula for Fresnel propagation is:

[0159]

[0160] In the formula, T 0 is the light intensity immediately after the sample, that is, the light intensity when the propagation distance is 0; is the phase information; η is the absorption information; ε is the phase recovery parameter obtained based on the absorption-phase duality; e is the natural constant; P D (x,y) is the Fresnel propagation factor; I D (x,y) is the plane light intensity distribution; λ is the light wavelength; denotes convolution; D is the propagation distance.

[0161] In this embodiment, in the L1 loss function value calculation module 004, during the process of calculating the L1 loss function value between the calculated intensity projection image and the initial intensity projection image through the loss function calculation strategy, the expression of the loss function calculation strategy is:

[0162]

[0163] where n is the number of pixels in the projection image; y i is the gray value of the i-th pixel of the input projection image; is the gray value of the i-th pixel of the calculated projection image.

[0164] It should be noted that for the information interaction, execution process, etc. among the above-mentioned system modules, since they are based on the same concept as the method embodiment in Embodiment 1 of the present application, the technical effects brought by them are the same as those of the method embodiment of the present application. For specific content, reference can be made to the description in the method embodiment shown above in the present application, and details will not be repeated here.

[0165] Embodiment 3

[0166] Embodiment 3 of the present invention provides a non-transitory computer-readable storage medium, in which program codes of an image phase recovery method based on X-ray in-line phase contrast imaging are stored, and the program codes include instructions for executing the image phase recovery method based on X-ray in-line phase contrast imaging in Embodiment 1 or any possible implementation manner thereof.

[0167] The computer-readable storage medium can be any available medium that can be accessed by a computer or a data storage device such as a server or a data center that integrates one or more available media. The available medium can be a magnetic medium (for example, a floppy disk, a hard disk, a magnetic tape), an optical medium (for example, a DVD), or a semiconductor medium (for example, a solid-state disk (SSD)), etc.

[0168] Embodiment 4

[0169] Embodiment 4 of the present invention provides an electronic device, including: a memory and a processor;

[0170] The processor and the memory complete communication with each other through a bus; the memory stores program instructions executable by the processor, and the processor can execute the image phase recovery method based on X-ray in-line phase contrast imaging in Embodiment 1 or any possible implementation manner thereof by calling the program instructions.

[0171] Specifically, the processor can be implemented by hardware or by software. When implemented by hardware, the processor can be a logic circuit, an integrated circuit, etc.; when implemented by software, the processor can be a general-purpose processor, which is implemented by reading software codes stored in the memory. The memory can be integrated in the processor or can exist independently outside the processor.

[0172] In the above embodiments, it can be implemented in whole or in part by software, hardware, firmware, or any combination thereof. When implemented using software, it can be implemented in whole or in part in the form of a computer program product. The computer program product includes one or more computer instructions. When the computer program instructions are loaded and executed on a computer, the processes or functions described in the embodiments of the present invention are generated in whole or in part. The computer can be a general-purpose computer, a special-purpose computer, a computer network, or other programmable systems. The computer instructions can be stored in a computer-readable storage medium or transmitted from one computer-readable storage medium to another. For example, the computer instructions can be transmitted from a website, computer, server, or data center to another website, computer, server, or data center by wire (such as coaxial cable, optical fiber, digital subscriber line (DSL)) or wireless (such as infrared, wireless, microwave, etc.).

[0173] Obviously, those skilled in the art should understand that the above-mentioned modules or steps of the present invention can be implemented by a general-purpose computing system. They can be concentrated on a single computing system or distributed on a network composed of multiple computing systems. Optionally, they can be implemented by program code executable by the computing system. Thus, they can be stored in a storage system and executed by the computing system. And in some cases, the steps shown or described can be executed in a different order than here, or they can be separately fabricated into individual integrated circuit modules, or multiple modules or steps of them can be fabricated into a single integrated circuit module to be implemented. In this way, the present invention is not limited to any specific combination of hardware and software.

[0174] Although the present invention has been described in detail above with general descriptions and specific embodiments, based on the present invention, some modifications or improvements can be made, which are obvious to those skilled in the art. Therefore, these modifications or improvements made without departing from the spirit of the present invention all fall within the scope of protection required by the present invention.

Claims

1. An image phase recovery method based on X-ray coaxial phase contrast imaging, characterized in that: include: Feed the initial intensity projection image into the restoration network; The initial intensity projection image is processed for the first time by a convolution block to obtain feature data X1; the image data X1 is iteratively processed four times by a maximum pooling operation and the convolution block to obtain feature data X2, feature data X3, feature data X4 and feature data X5 in sequence; Processing the feature data X1, the feature data X2, the feature data X3, the feature data X4 and the feature data X5 in a set order through an upsampling convolution block and the convolution block in turn to obtain a phase-shifted image after phase recovery; Performing Fresnel propagation processing on the phase-recovered phase-shifted image to obtain a calculated intensity projection image; Calculating the L1 loss function value between the calculated intensity projection image and the initial intensity projection image using a loss function calculation strategy to obtain the L1 loss function value; Iteratively optimizing and updating the parameters of the restoration network through the L1 loss function value and the back propagation algorithm; when the L1 loss function value reaches the set convergence value, outputting the final phase-shifted image; The final phase-shifted image is reconstructed by a filtered back-projection reconstruction algorithm to obtain a reconstructed image.

2. The image phase recovery method based on X-ray coaxial phase contrast imaging according to claim 1, characterized in that: In the process of acquiring the phase-shifted image after phase recovery, the steps of sequentially processing through the upsampling convolution block and the convolution block are: Input the feature data X5 into the upsampling convolution block for processing to obtain a first output result; connect the first output result with the feature data X4 to obtain a first connection result; Inputting the first connection result into the convolution block and the upsampling convolution block for processing to obtain a second output result; Connecting the second output result with the feature data X3 to obtain a second connection result; Input the second connection result into the convolution block and the upsampling convolution block for processing to obtain a third output result; connect the third output result with the feature data X2 to obtain a third connection result; Inputting the third connection result into the convolution block and the upsampling convolution block for processing to obtain a fourth output result; Connecting the fourth output result with the feature data X1 to obtain a fourth connection result; The fourth connection result is input into the convolution block for 1×1 convolution processing to obtain the phase-recovered phase-shifted image.

3. The image phase recovery method based on X-ray coaxial phase contrast imaging according to claim 2, characterized in that: In the process of processing the feature data by the convolution block and the upsampling convolution block, the feature data is processed by 3×3 two-dimensional convolution by the convolution block and the upsampling convolution block to obtain two-dimensional convolution data; the two-dimensional convolution data is normalized by setting a normalization strategy to obtain normalized data; the normalized data is processed by the LeakyReLu calculation strategy to obtain processed data; the processed data is processed by the convolution block and the upsampling convolution block for the next round; The expression for setting the normalization strategy is: In the formula, is the jth feature of the kth sample; μ is the mean of each feature; σ 2 is the variance of each feature; ε` is a small constant used to prevent the denominator from being zero; The expression of the Leaky ReLu calculation strategy is: Where x is any grayscale value of the image; α is a small positive number, the default value is 0.

01.

4. The image phase recovery method based on X-ray coaxial phase contrast imaging according to claim 3 is characterized in that: In the process of performing Fresnel propagation processing on the phase-shifted image after phase recovery, the calculation formula of Fresnel propagation is: Where T0 is the light intensity immediately after the sample, that is, when the propagation distance is 0; is the phase information; η is the absorption information; ε is the phase recovery parameter based on the absorption-phase duality; e is the natural constant; P D (x,y) is the Fresnel propagation factor; I D (x, y) is the plane light intensity distribution; λ is the wavelength of light; represents convolution; D is the propagation distance.

5. The image phase recovery method based on X-ray coaxial phase contrast imaging according to claim 4, characterized in that: In the process of calculating the L1 loss function value between the calculated intensity projection image and the initial intensity projection image by the loss function calculation strategy, the expression of the loss function calculation strategy is: Where n is the number of pixels in the projection image; y i is the gray value of the i-th pixel of the input projection image; To calculate the gray value of the i-th pixel of the projection image.

6. An image phase recovery device based on X-ray coaxial phase contrast imaging, using the image phase recovery method based on X-ray coaxial phase contrast imaging according to any one of claims 1 to 5, characterized in that: include: A feature data acquisition module, used for inputting the initial intensity projection image into the restoration network; The initial intensity projection image is processed for the first time by a convolution block to obtain feature data X1; the image data X1 is iteratively processed four times by a maximum pooling operation and the convolution block to obtain feature data X2, feature data X3, feature data X4 and feature data X5 in sequence; A recovered phase-shift image acquisition module is used to process the feature data X1, the feature data X2, the feature data X3, the feature data X4 and the feature data X5 in a set order through an upsampling convolution block and the convolution block in turn to obtain a phase-recovered phase-shift image; A calculated intensity projection image acquisition module is used to perform Fresnel propagation processing on the phase-shifted image after phase recovery to obtain a calculated intensity projection image; An L1 loss function value calculation module, used to calculate the L1 loss function value between the calculated intensity projection image and the initial intensity projection image through a loss function calculation strategy to obtain the L1 loss function value; The final phase-shifted image acquisition module is used to iteratively optimize and update the parameters of the restoration network through the L1 loss function value and the back-propagation algorithm; when the L1 loss function value reaches the set convergence value, the final phase-shifted image is output; The reconstructed image acquisition module is used to reconstruct the final phase-shifted image by using a filtered back-projection reconstruction algorithm to obtain a reconstructed image.

7. The image phase recovery device based on X-ray coaxial phase contrast imaging according to claim 6, characterized in that: In the restored phase-shift image acquisition module, the submodules processed by the upsampling convolution block and the convolution block include: A first connection result acquisition submodule is used to input the feature data X5 into the upsampling convolution block for processing to obtain a first output result; and connect the first output result with the feature data X4 to obtain a first connection result; A second connection result acquisition submodule is used to input the first connection result into the convolution block and the upsampling convolution block for processing to obtain a second output result; and connect the second output result with the feature data X3 to obtain a second connection result; A third connection result acquisition submodule, used for inputting the second connection result into the convolution block and the upsampling convolution block for processing to obtain a third output result; and connecting the third output result with the feature data X2 to obtain a third connection result; A fourth connection result acquisition submodule is used to input the third connection result into the convolution block and the upsampling convolution block for processing to obtain a fourth output result; and connect the fourth output result with the feature data X1 to obtain a fourth connection result; The phase-shifted image acquisition submodule after phase recovery is used to input the fourth connection result into the convolution block for 1×1 convolution processing to obtain the phase-shifted image after phase recovery.

8. The image phase recovery device based on X-ray coaxial phase contrast imaging according to claim 7, characterized in that: In the restored phase-shift image acquisition module, in the process of processing the feature data by the convolution block and the upsampling convolution block, the feature data is subjected to 3×3 two-dimensional convolution processing by the convolution block and the upsampling convolution block to obtain two-dimensional convolution data; the two-dimensional convolution data is subjected to normalization processing by setting a normalization strategy to obtain normalized data; the normalized data is processed by the Leaky ReLu calculation strategy to obtain processed data; the processed data is subjected to the next round of processing by the convolution block and the upsampling convolution block; The expression for setting the normalization strategy is: In the formula, is the jth feature of the kth sample; μ is the mean of each feature; σ 2 is the variance of each feature; ε` is a small constant used to prevent the denominator from being zero; The expression of the Leaky ReLu calculation strategy is: Where x is any grayscale value of the image; α is a small positive number, the default value is 0.

01.

9. The image phase recovery device based on X-ray coaxial phase contrast imaging according to claim 8, characterized in that: In the calculation intensity projection image acquisition module, in the process of performing Fresnel propagation processing on the phase-recovered phase-shifted image, the calculation formula of Fresnel propagation is: Where T0 is the light intensity immediately after the sample, that is, when the propagation distance is 0; is the phase information; η is the absorption information; ε is the phase recovery parameter based on the absorption-phase duality; e is the natural constant; P D (x,y) is the Fresnel propagation factor; I D (x, y) is the plane light intensity distribution; λ is the wavelength of light; represents convolution; D is the propagation distance.

10. The image phase recovery device based on X-ray coaxial phase contrast imaging according to claim 9, characterized in that: In the L1 loss function value calculation module, in the process of calculating the L1 loss function value between the calculated intensity projection image and the initial intensity projection image by using the loss function calculation strategy, the expression of the loss function calculation strategy is: Where n is the number of pixels in the projection image; y i is the gray value of the i-th pixel of the input projection image; To calculate the gray value of the i-th pixel of the projection image.

Images