A method and apparatus for X-ray finite-angle CT image reconstruction based on generalized total variation.

Abstract

This invention discloses a method and apparatus for X-ray finite-angle CT image reconstruction based on generalized total variation, comprising: Step 1, using a finite-angle CT scan dataset p, and employing an image reconstruction operator R related to the scan geometry parameter set G. G Update estimated image u (k) Obtain the intermediate image u (k+1 / 2) k is the iteration number; Step 2, use the sparsity of the gradient to constrain the intermediate image u (k+1 / 2) The estimated image u is obtained. (k+1) Step 3: Determine if the iteration requirement has been met. If so, terminate the iteration and output the estimated image u. (k+1) Conversely, if the boundary is not clear, proceed to step 1. This invention overcomes the problem of blurred boundaries in existing finite-angle CT image reconstruction algorithms.

Description

Technical Field

[0001] This invention relates to the field of X-ray CT image reconstruction technology, and in particular to a method and apparatus for X-ray finite-angle CT image reconstruction based on generalized total variation. Background Technology

[0002] In industrial nondestructive testing, the inspection of plate-shaped objects such as circuit boards and wings is a common challenge. Because the object under test has a long boundary in a certain direction, the following situations may arise: 1. To reconstruct a high-resolution image, the distance between the X-ray source and the object must be small. To avoid collisions between the X-ray source and the object, the turntable cannot rotate a full circle during the scanning process; 2. Even if the turntable can rotate a full circle, some X-rays are severely attenuated due to passing through the long boundary of the object, resulting in poor signal-to-noise ratio of the detection data at the corresponding angle, rendering it unusable. These situations all lead to the acquisition of effective projection data only within a certain angular range.

[0003] The problem of reconstructing an image using data collected within a certain angular range is called finite-angle CT image reconstruction, which belongs to the category of incomplete data image reconstruction problems. Images obtained using traditional image reconstruction algorithms (such as FDK, SART, etc.) often suffer from severely blurred boundaries, affecting the usability of CT images.

[0004] Therefore, it is desirable to have a technical solution to overcome or at least mitigate at least one of the aforementioned defects of the prior art. Summary of the Invention

[0005] The purpose of this invention is to provide a method for X-ray finite-angle CT image reconstruction based on generalized total variation to overcome or at least mitigate at least one of the above-mentioned defects of the prior art, and improve the quality and usability of finite-angle CT images.

[0006] To achieve the above objectives, this invention provides a method for X-ray finite-angle CT image reconstruction based on generalized total variation, comprising:

[0007] Step 1: Using the finite-angle CT scan dataset p, the image reconstruction operator R is used in conjunction with the scan geometry parameter set G. G Update estimated image u (k) Obtain the intermediate image u (k+1 / 2) k is the number of iterations;

[0008] Step 2, according to equation (1), use the sparsity of the gradient to constrain the intermediate image u. (k+1 / 2) The estimated image u is obtained. (k+1) ;

[0009] u (k+1) =P(u (k+1 / 2) (1)

[0010] In equation (1), P represents the operator constraining the gradient sparsity of the image, obtained by solving the optimization problem provided in equation (2):

[0011]

[0012] In formula (2):

[0013] For regularization terms;

[0014] For data fidelity items;

[0015] These represent the gradients in the horizontal x-direction and the vertical y-direction of the input image u, respectively.

[0016] u represents the input image;

[0017] diag(·) represents the weight matrix, which assigns penalty weights to the horizontal and vertical gradients of each pixel.

[0018] ||·||2 represents the l2 norm;

[0019] Represents the square of the l2 norm;

[0020] ||·||1 represents the l1 norm;

[0021] ||·|| 2,1 The l1 norm represents the gradient magnitude;

[0022] |·| indicates performing absolute value operations on each pixel;

[0023] κ represents the parameter used to adjust the penalty weight;

[0024] G σ Represents the Gaussian convolution kernel;

[0025] * indicates a convolution operation;

[0026] T represents transpose;

[0027] λ represents the parameter used to control the weights of the regularization term and the data fidelity term;

[0028] Step 3: Determine if the iteration requirement has been met. If so, terminate the iteration and output the estimated image u. (k+1) Conversely, proceed to step 1.

[0029] Furthermore, the methods for calculating the penalty weights in the horizontal direction of pixels include:

[0030] Based on the gradient of the pixel in the horizontal direction, the first part of the weight matrix diag(·) in equation (2) is used. The penalty weight in the horizontal direction of the pixel is calculated.

[0031] Furthermore, the methods for calculating the penalty weights in the vertical direction of pixels include:

[0032] Based on the vertical gradient of the pixel, the second part of the weight matrix diag(·) in equation (2) is used. Calculate the penalty weight in the vertical direction of the pixel.

[0033] Furthermore, estimate image u (k) initial value u (0) All pixel values ​​are set to 0.

[0034] Furthermore, the scanning geometry parameter set G includes the distance from the X-ray source to the detector center, the distance from the X-ray source to the turntable center, the number of detector units, the detector unit size, the number of scanning angles, and the angle sampling interval.

[0035] The present invention also provides an X-ray finite-angle CT image reconstruction device based on generalized total variation, comprising:

[0036] The iterative processing estimation module is used to utilize a finite-angle CT scan dataset p, through an image reconstruction operator R associated with the scan geometry parameter set G. G Update estimated image u (k) Obtain the intermediate image u (k+1 / 2) k is the number of iterations;

[0037] The gradient sparsity constraint module is used to constrain the intermediate image u using the gradient sparsity according to equation (1). (k +1 / 2) The estimated image u is obtained. (k+1) ;

[0038] u (k+1) =P(u (k+1 / 2) (1)

[0039] In equation (1), P represents the operator constraining the gradient sparsity of the image, obtained by solving the optimization problem provided in equation (2):

[0040]

[0041] In formula (2):

[0042] For regularization terms;

[0043] For data fidelity items;

[0044] These represent the gradients in the horizontal x-direction and the vertical y-direction of the input image u, respectively.

[0045] u represents the input image;

[0046] diag(·) represents the weight matrix, which assigns penalty weights to the horizontal and vertical gradients of each pixel.

[0047] ||·||2 represents the l2 norm;

[0048] Represents the square of the l2 norm;

[0049] ||·||1 represents the l1 norm;

[0050] ||·|| 2,1 The l1 norm represents the gradient magnitude;

[0051] |·| indicates performing absolute value operations on each pixel;

[0052] κ represents the parameter used to adjust the penalty weight;

[0053] G σ Represents the Gaussian convolution kernel;

[0054] * indicates a convolution operation;

[0055] T represents transpose;

[0056] λ represents the parameter used to control the weights of the regularization term and the data fidelity term;

[0057] The judgment module determines whether the iteration requirement has been met. If so, the iteration terminates and the estimated image u is output. (k+1) Conversely, the iterative estimation module continues to update the estimated image u. (k+1) .

[0058] Furthermore, the methods for calculating the penalty weights in the horizontal direction of pixels include:

[0059] Based on the gradient of the pixel in the horizontal direction, the first part of the weight matrix diag(·) in equation (2) is used. The penalty weight in the horizontal direction of the pixel is calculated.

[0060] Furthermore, the methods for calculating the penalty weights in the vertical direction of pixels include:

[0061] Based on the vertical gradient of the pixel, the second part of the weight matrix diag(·) in equation (2) is used. Calculate the penalty weight in the vertical direction of the pixel.

[0062] Furthermore, estimate image u (k)initial value u (0) All pixel values ​​are set to 0.

[0063] Furthermore, the scanning geometry parameter set G includes the distance from the X-ray source to the detector center, the distance from the X-ray source to the turntable center, the number of detector units, the detector unit size, the number of scanning angles, and the angle sampling interval.

[0064] The present invention has the following advantages due to the adoption of the above technical solutions:

[0065] The method provided by this invention adaptively assigns different penalty weights to the gradients of the image u in the ideal placement state (horizontal direction) and the non-ideal placement state (non-horizontal direction), and utilizes the sparsity of the image gradient to gradually recover the boundary information of the image. This can effectively suppress image blurring and artifacts, improve the reconstruction quality of finite-angle CT images, and thus overcome the problem of boundary blurring in existing finite-angle CT image reconstruction algorithms. Attached Figure Description

[0066] Figure 1 This is a flowchart of an embodiment of the present invention;

[0067] Figure 2 The image shown is a scanned phantom image from an embodiment of the present invention.

[0068] Figure 3 The image is a reconstructed image of SART in an embodiment of the present invention;

[0069] Figure 4 The image is a reconstructed image obtained using an embodiment of the present invention. Detailed Implementation

[0070] The present invention will now be described in detail with reference to the accompanying drawings and embodiments.

[0071] like Figure 1 As shown, the X-ray finite-angle CT image reconstruction method based on generalized total variation provided in this embodiment of the invention includes:

[0072] Step 1: Using the finite-angle CT scan dataset p, the image reconstruction operator R is used in conjunction with the scan geometry parameter set G. G Update estimated image u (k) Obtain the intermediate image u (k+1 / 2) , where k is the iteration number, k = 0, 1, ..., N. The update method can be an existing update method such as ART, SART, or FBP. The CT scan geometric parameter set G includes the distance from the X-ray source to the detector center (SDD), the distance from the X-ray source to the turntable center (SOD), the number of detector elements, the detector element size, the number of scan angles, and the angle sampling interval. R GThe reconstruction is obtained by selecting either an iterative reconstruction algorithm or an analytical reconstruction algorithm. The iterative reconstruction algorithms are ART, SART, and EM; the analytical reconstruction algorithms are FBP and BPF.

[0073] Step 2, according to equation (1), use the sparsity of the gradient to constrain the intermediate image u. (k+1 / 2) The estimated image u is obtained. (k+1) ;

[0074] u (k+1) =P(u (k+1 / 2) (1)

[0075] In equation (1), P represents the operator constraining the gradient sparsity of the image, obtained by solving the optimization problem provided in equation (2), which can be solved using the ADMM method:

[0076]

[0077] In formula (2):

[0078] For regularization terms;

[0079] For data fidelity items;

[0080] Let x and y represent the gradients in the horizontal and vertical directions of the input image u, respectively. Taking the input image u as a reference, the horizontal direction (left-right direction) is considered the horizontal direction in this invention, while the vertical direction (up-down direction) is considered the vertical direction in this invention. Therefore… The value is the difference between the current pixel value and its adjacent right pixel value. The value is the difference between the current pixel value and its adjacent pixel value below it. When encountering a boundary, The value is the difference between the current pixel value and the pixel value in the first row. The value is the difference between the current pixel value and the pixel value in the first column;

[0081] u represents the input image, which must be compared with u (k+1 / 2) Approaching closely;

[0082] diag(·) represents the weight matrix, which depends on the local information of the input image u during the iteration process. Different penalty weights are assigned to the horizontal and vertical gradients of each pixel. The specific method for determining the penalty weights may include: based on the horizontal gradient of the pixel, using the first part of the weight matrix diag(·) in equation (2) Calculate the horizontal penalty weight of the pixel; similarly, based on the vertical gradient of the pixel, use the second part of the weight matrix diag(·) in equation (2). The vertical penalty weight for that pixel is calculated. Therefore, when the horizontal and vertical gradients of the same pixel are different, the two numbers in diag(·) are different, and thus the horizontal and vertical penalty weights for that pixel are different. For different pixels, the values ​​of the horizontal and vertical gradients are different, so the penalty weights for different pixels are also different.

[0083] ||·||2 represents the l2 norm, which is the square root of the sum of the squares of the parameters, i.e., the Euclidean distance. A smaller l2 norm makes each element of the parameters very small, close to 0, and the regularization term is differentiable everywhere. For example, u represents an M×M input image, u... i,j Represents the pixel in the i-th row and j-th column of u, ||u||2=

[0084] Represents the square of the l2 norm;

[0085] ||·||1 represents the L1 norm, which is the sum of the absolute values ​​of the parameters. It can prevent overfitting. The L1 norm does not make the parameters equal to 0 but close to 0, which has the effect of smoothing the parameters. The optimized parameter vector is often relatively sparse. For example: u = (u1, u2), where u1 and u2 represent matrices of size M × M.

[0086] ||·|| 2,1 The l1 norm represents the magnitude of the gradient, and the magnitude of the gradient is defined as follows: The square root of the sum of squares is obtained;

[0087] |·| indicates performing absolute value operations on each pixel;

[0088] κ represents the parameter used to adjust the penalty weight. This parameter can better distinguish between smooth regions and boundary regions of the image. The specific value was obtained through experiments.

[0089] G σ Represents the Gaussian convolution kernel;

[0090] * indicates a convolution operation;

[0091] T represents transpose;

[0092] λ represents the parameter used to control the weights of the regularization term and the data fidelity term; the specific value was obtained through experiments.

[0093] Equation (2) can solve the problem that existing finite-angle CT image reconstruction algorithms reconstruct images that are blurred along a specific direction and have finite-angle artifacts.

[0094] Of course, the operator P that constrains the gradient sparsity of an image can also be obtained by applying existing methods such as applying boundary-preserving diffusion and boundary-preserving smoothness regularization terms, and applying curvature constraint regularization terms.

[0095] Step 3: Determine if the iteration requirement has been met. If so, terminate the iteration and output the estimated image u. (k+1) Conversely, proceed to step 1.

[0096] The method provided in this embodiment, through the... Figure 2 The input image u is adaptively assigned different penalty weights in the horizontal and vertical directions. By utilizing the sparsity of the image gradient, the boundary information of the image is gradually recovered, thereby effectively suppressing image blurring and artifacts, improving the reconstruction quality of finite-angle CT images, and overcoming the boundary blurring problem in existing finite-angle CT image reconstruction algorithms.

[0097] In one embodiment, the estimated image u (k) initial value u (0) All pixel values ​​are set to 0, but other values ​​are also possible.

[0098] In one embodiment, an iteration termination threshold ε and an upper limit for the number of iterations N can also be set. Then, when determining whether the iteration requirement has been met in step 3, it can be determined whether the difference between two adjacent iteration images is less than a given threshold, i.e., ||u (k+1) -u (k) The value ||≤ε can also be used to determine whether the number of iterations has reached the upper limit N of the number of iterations.

[0099] This invention also provides an X-ray finite-angle CT image reconstruction device based on generalized total variation, which includes an iterative processing estimation module, a gradient sparsity constraint module, and a judgment module, wherein:

[0100] The iterative processing estimation module is used to utilize a finite-angle CT scan dataset p, through an image reconstruction operator R associated with the scan geometry parameter set G. G Update estimated image u (k) Obtain the intermediate image u (k+1 / 2) k is the number of iterations.

[0101] The gradient sparsity constraint module is used to constrain the intermediate image u using the gradient sparsity according to equation (1). (k+1 / 2) The estimated image u is obtained. (k+1) The operator P used to constrain the sparsity of the image gradient can be solved by the optimization problem provided by equation (2).

[0102] The decision module determines whether the iteration requirement has been met. If so, the iteration terminates and the estimated image u is output. (k+1) Conversely, the iterative estimation module continues to update the estimated image u. (k+1) .

[0103] To better demonstrate the advantages of the finite-angle CT image reconstruction method provided by this invention in terms of reconstruction effect, the method described in this invention is compared with the existing typical algorithm SART in conjunction with a specific embodiment below.

[0104] In this embodiment, data was acquired from an industrial CT system. The scanned sample was a PCB phantom. The PCB board was not ideally placed on the scanning stage, and the image is as follows. Figure 2 As shown. The experimental voltage and current were 140kV and 100uA, respectively, and projection data were collected every 0.5 degrees within a 120-degree range.

[0105] The PCB board scan data were reconstructed using the SART method and the method of this invention, respectively. The reconstruction results are as follows: Figure 3 , Figure 4 As shown. Wherein: Figure 3 This is the reconstruction result of the SART method; Figure 4 This is the reconstruction result of the method of this invention. It can be seen that the image reconstructed by the SART method suffers from severe boundary blurring and stripe artifacts due to incomplete scanning data; the boundary information in the PCB mold is almost completely lost. Figure 4 As can be seen, the reconstructed image obtained by the embodiment of the present invention has clear boundaries and accurately reconstructs long strip-shaped boundaries that are difficult to reconstruct.

[0106] Finally, it should be noted that the above embodiments are only used to illustrate the technical solutions of the present invention, and not to limit them. Those skilled in the art should understand that modifications can be made to the technical solutions described in the foregoing embodiments, or equivalent substitutions can be made to some of the technical features; these modifications or substitutions do not cause the essence of the corresponding technical solutions to deviate from the spirit and scope of the technical solutions of the embodiments of the present invention.

Claims

1. A method for reconstructing X-ray finite-angle CT images based on generalized total variation, characterized in that, include: Step 1: Using the finite-angle CT scan dataset p, the image reconstruction operator R is used in conjunction with the scan geometry parameter set G. G Update estimated image u (k) Obtain the intermediate image u (k+1 / 2) k is the iteration number, k = 0, 1, ..., N; where the update method is the existing update method of ART, SART, or FBP; the CT scan geometric parameter set G includes the distance from the X-ray source to the detector center (SDD), the distance from the X-ray source to the turntable center (SOD), the number of detector units, the detector unit size, the number of scan angles, and the angle sampling interval; R G The reconstruction can be obtained by selecting either an iterative reconstruction algorithm or an analytical reconstruction algorithm. The iterative reconstruction algorithms are ART, SART, and EM; the analytical reconstruction algorithms are FBP and BPF. Step 2, according to equation (1), use the sparsity of the gradient to constrain the intermediate image u. (k+1 / 2) The estimated image u is obtained. (k+1) ; in (k+1) =P(in (k+1 / 2) ) (1) In equation (1), P represents the operator that constrains the gradient sparsity of the image, obtained by solving the optimization problem provided in equation (2): In formula (2): For regularization terms; For data fidelity items; These represent the gradients in the horizontal x-direction and the vertical y-direction of the input image u, respectively; u represents the input image; diag(·) represents the weight matrix, which assigns penalty weights to the horizontal and vertical gradients of each pixel. ||·||2 represents the l2 norm; Denotes the square of the l2 norm; ||·||1 represents the l1 norm; ||·|| 2,1 The 1-norm represents the gradient magnitude. |·| indicates performing absolute value operations on each pixel; κ represents the parameter used to adjust the penalty weight; G σ Represents the Gaussian convolution kernel; * indicates a convolution operation; T represents transpose; λ represents the parameter used to control the weights of the regularization term and the data fidelity term; Step 3: Determine if the iteration requirement has been met. If so, terminate the iteration and output the estimated image u. (k+1) Conversely, proceed to step 1.

2. The X-ray finite-angle CT image reconstruction method based on generalized total variation as described in claim 1, characterized in that, Methods for calculating the horizontal penalty weights of pixels include: Based on the gradient of the pixel in the horizontal direction, the first part of the weight matrix diag(·) in equation (2) is used. The penalty weight in the horizontal direction of the pixel is calculated.

3. The X-ray finite-angle CT image reconstruction method based on generalized total variation as described in claim 1 or 2, characterized in that, Methods for penalizing pixel weights in the vertical direction include: Based on the vertical gradient of the pixel, the second part of the weight matrix diag(·) in equation (2) is used. Calculate the penalty weight in the vertical direction of the pixel.

4. The X-ray finite-angle CT image reconstruction method based on generalized total variation as described in claim 3, characterized in that, Estimated image u (k) initial value u (0) All pixel values ​​are set to 0.

5. A device for X-ray finite-angle CT image reconstruction based on generalized total variation, characterized in that, include: The iterative processing estimation module is used to utilize a finite-angle CT scan dataset p, through an image reconstruction operator R associated with the scan geometry parameter set G. G Update estimated image u (k) Obtain the intermediate image u (k+1 / 2) k is the iteration number, k = 0, 1, ..., N; where the update method is the existing update method of ART, SART, or FBP; the CT scan geometric parameter set G includes the distance from the X-ray source to the detector center (SDD), the distance from the X-ray source to the turntable center (SOD), the number of detector units, the detector unit size, the number of scan angles, and the angle sampling interval; R G The reconstruction can be obtained by selecting either an iterative reconstruction algorithm or an analytical reconstruction algorithm. The iterative reconstruction algorithms are ART, SART, and EM; the analytical reconstruction algorithms are FBP and BPF. The gradient sparsity constraint module is used to constrain the intermediate image u using the gradient sparsity according to equation (1). (k+1 / 2) The estimated image u is obtained. (k+1) ; in (k+1) =P(in (k+1 / 2) ) (1) In equation (1), P represents the operator that constrains the gradient sparsity of the image, obtained by solving the optimization problem provided in equation (2): In formula (2): For regularization terms; For data fidelity items; These represent the gradients in the horizontal x-direction and the vertical y-direction of the input image u, respectively; u represents the input image; diag(·) represents the weight matrix, which assigns penalty weights to the horizontal and vertical gradients of each pixel. ||·||2 represents the l2 norm; Denotes the square of the l2 norm; ||·||1 represents the l1 norm; ||·|| 2,1 The l1 norm represents the gradient magnitude; |·| indicates performing absolute value operations on each pixel; κ represents the parameter used to adjust the penalty weight; G σ Represents the Gaussian convolution kernel; * indicates a convolution operation; T represents transpose; λ represents the parameter used to control the weights of the regularization term and the data fidelity term; The judgment module determines whether the iteration requirement has been met. If so, the iteration terminates and the estimated image u is output. (k+1) Conversely, the iterative estimation module continues to update the estimated image u. (k+1) .

6. The X-ray finite-angle CT image reconstruction device based on generalized total variation as described in claim 5, characterized in that, Methods for calculating the horizontal penalty weights of pixels include: Based on the gradient of the pixel in the horizontal direction, the first part of the weight matrix diag(·) in equation (2) is used. The penalty weight in the horizontal direction of the pixel is calculated.

7. The X-ray finite-angle CT image reconstruction device based on generalized total variation as described in claim 5 or 6, characterized in that, Methods for penalizing pixel weights in the vertical direction include: Based on the vertical gradient of the pixel, the second part of the weight matrix diag(·) in equation (2) is used. Calculate the penalty weight in the vertical direction of the pixel.

8. The X-ray finite-angle CT image reconstruction device based on generalized total variation as described in claim 7, characterized in that, Estimated image u (k) initial value u (0) All pixel values ​​are set to 0.

Images