CT image reconstruction method and system based on maximum residual extension oblique projection

Abstract

A CT image reconstruction method and system based on maximum residual extended oblique projection belongs to the field of medical image processing technology. The method involves splitting the CT image matrix to be reconstructed column-wise to form several linear equation systems; and setting auxiliary vectors. z k Sum of solutions x k The initial value is z 0 and x 0; In the first iteration, the working column and working row are selected based on the maximum residual index, and the iteration is performed using the orthogonal projection method. z 0 and x 0; then continue to select working columns / rows according to the maximum residual index, and calculate the new column / row projection direction vector using the working columns / rows of the two adjacent steps, and iterate the auxiliary vector using the new column projection direction vector and the appropriate step size. z k If the new row projection direction vector is not zero, then the solution vector is obtained by iterating through the new row projection direction vector and a suitable step size. x k Otherwise, solve the vector by iteratively projecting the orthogonal vector. x k Output the final solution vector when the termination condition is met. x k+1 ;Arrange all final solution vectors by column x k+1 By combining the data, a reconstructed CT image is obtained. This disclosure improves the accuracy and efficiency of CT image reconstruction.

Description

Technical Field

[0001] This disclosure belongs to the field of medical image processing technology, and in particular relates to a method and system for CT image reconstruction based on maximum residual extended oblique projection. Background Technology

[0002] Computed tomography (CT) technology has become an indispensable imaging method in clinical medical diagnosis. CT technology utilizes the attenuation properties of X-rays after they penetrate the human body to generate high-resolution cross-sectional images of the human body. Although CT plays a vital role in clinical medical diagnosis, the ionizing radiation emitted by X-rays during the scanning process poses a potential carcinogenic risk. Therefore, how to generate high-quality images while minimizing CT radiation is a key research direction in the field of CT.

[0003] As CT radiation dose decreases, image noise increases significantly, leading to a reduction in image quality. Therefore, advanced image reconstruction methods are essential to reduce image noise. Algebraic reconstruction, also known as the Kaczmarz method, is a core iterative method applied in CT image reconstruction. Its basic principle is to model the physical process of imaging as a large system of linear equations and solve iteratively using an orthogonal iterative format of line-by-line projection. This method is widely used due to its simple structure, ease of programming, and low memory requirements.

[0004] The classic Kaczmarz method for CT image reconstruction has the following problems:

[0005] 1. The classic Kaczmarz method iteratively selects only a single row of the system matrix for calculation, sequentially choosing working rows and performing orthogonal projection iterations on the solution vectors. This method has a slow convergence speed, resulting in a long image reconstruction time. To improve convergence efficiency, researchers have proposed several improved methods: the stochastic Kaczmarz method, which significantly improves convergence speed by probabilistically selecting working rows; and the maximum residual Kaczmarz method, which further accelerates iteration efficiency by selecting the row containing the maximum residual as the working row.

[0006] 2. When the CT radiation dose is reduced, the noise in the reconstructed image increases significantly, and linear systems often exhibit incompatible characteristics, making it difficult for the above methods to achieve ideal denoising results. Therefore, the Randomized Extended Kaczmarz (REK) method introduces auxiliary vectors and combines row and column double projection to solve for the minimum norm least squares solution, further improving numerical stability; the Maximum Residual Extended Kaczmarz (MREK) method uses the maximum residual row and column as the working unit to improve the convergence speed. Both of these methods use the traditional orthogonal projection method for iteration to achieve denoising effects.

[0007] 3. For the image reconstruction problem of low-dose CT, the REK method has a large number of iterations and a long computation time. Even if the working cell selection strategy is optimized, the MREK method still uses the traditional orthogonal projection direction, and there is still room for improvement in convergence efficiency. It is difficult to simultaneously achieve both reconstruction accuracy and imaging efficiency. Summary of the Invention

[0008] To address the aforementioned technical problems and overcome the shortcomings of existing methods in simultaneously achieving both reconstruction accuracy and imaging efficiency in low-dose CT images, this disclosure provides a CT image reconstruction method and system based on maximum-residual extended oblique projection, also known as the Maximum-Residual Extended Kaczmarz method with Oblique projection (MREKO). This disclosure uses oblique projection instead of traditional orthogonal projection, improving reconstruction accuracy while ensuring imaging efficiency.

[0009] The purpose of this disclosure is achieved through the following technical solution:

[0010] This disclosure presents a CT image reconstruction method based on maximum residual extended oblique projection, comprising the following steps:

[0011] (1) Load the CT image matrix and system matrix to be reconstructed, and set the maximum number of iterations and error range;

[0012] (2) Split the CT image matrix to be reconstructed into columns to form several linear equations;

[0013] (3) Set auxiliary vectors Sum of solutions The initial values ​​are respectively and The working column and working row are selected according to the maximum residual criterion, and the initial values ​​of the auxiliary vector and solution vector are iterated through the orthogonal projection method.

[0014] (4) Select the working column according to the maximum residual index, calculate the new column projection direction through the working columns of two adjacent steps, and calculate the appropriate step size through the working columns of two adjacent steps. Iterate the auxiliary vector according to the new column projection direction and the appropriate step size. ;

[0015] (5) Select the working row according to the maximum residual index. Calculate the new row projection direction using the working rows of the two adjacent steps. If the new row projection direction is not a zero vector, then use the auxiliary vector obtained from the two adjacent working rows and step (4). Calculate the appropriate step size, and iterate the solution vector based on the new row projection direction and the appropriate step size. If the new row projection direction is the zero vector, then the vector is solved iteratively through orthogonal projection. ;

[0016] (6) Determine the current solution vector Check if the termination condition is met. If so, change the current solution vector. Output the final solution vector and, starting from step (3), solve the next system of linear equations until the solution is complete; otherwise, let the iteration number be... Repeat steps (4)-(6);

[0017] (7) Combine all the final solution vectors obtained by solving the problem. By combining columns, the reconstructed CT images are obtained.

[0018] Furthermore, in step (2), the process of reconstructing the CT image matrix is ​​regarded as solving a linear system. The process, in which , , ,in For the system matrix, For the ideal target CT image matrix, It is a known matrix of CT images to be reconstructed; the known matrix of CT images to be reconstructed Split by column Then solve Transform into solving , ,in yes The The column, i.e., the first column of the ideal target CT image matrix. Column pixel values; yes The The column, i.e., the first column of the CT image matrix to be reconstructed. Column pixel values.

[0019] Furthermore, in steps (3) and (4), the working columns and working rows are selected according to the maximum residual criterion. The selection index for the working columns is in the system matrix. In all columns, a total Columns, selected to make transpose and auxiliary vector The column index with the largest absolute value of the inner product. And use it as a target for the next step of work. The selection of working row indicators is based on the system matrix. In all rows, a total Line, select the residual term Largest row index And use it as a target for the next step of work. ;in and They represent The Line and number List, for The item.

[0020] Furthermore, in step (3),

[0021] Iterate the initial value of the auxiliary vector using the orthogonal projection method. Initial value of the solution vector The auxiliary vectors after the first iteration are obtained respectively. Sum of solutions Specifically:

[0022] Auxiliary vector It is determined by the initial auxiliary vector To the hyperplane The result obtained by orthogonal projection onto the top;

[0023] Solution vector It is based on the initial solution vector To the hyperplane The result obtained by orthogonal projection;

[0024] in and They represent and The item, and Representing the system matrix respectively The Column and number OK, It is a system matrix The transpose of a column and The column and row indicators selected for the first time.

[0025] Furthermore, in step (4), the working column is selected according to the maximum residual index, based on the first... The work column selected in the next iteration With the The work column selected in the next iteration Calculate a new column projection direction vector that matches the projection noise suppression requirements of CT images:

[0026] ,

[0027] Among the symbols Represents the dot product of two vectors. The Euclidean norm of a vector;

[0028] And through the obtained new column projection direction vector The first Auxiliary vector for the next iteration Project onto the following hyperplane to eliminate noise:

[0029] ,

[0030] Let the auxiliary vector be at this time. The step size is ,Right now ,

[0031] The step size is obtained as follows: ;

[0032] Therefore, the auxiliary vector The iterative formula is:

[0033] .

[0034] Further, in step (5), the working row is selected according to the maximum residual index, and the working row is selected according to the maximum residual index, based on the first... The working row selected in the next iteration With the The working row selected in the next iteration Calculate a new row projection direction to match the pixel reconstruction requirements of the CT image:

[0035] ,

[0036] When the new line projection direction When it is a non-zero vector, the first... Solution vector of the next iteration Project onto the following hyperplane to eliminate artifacts:

[0037] ,

[0038] in, For vectors The item, The auxiliary vector obtained in step (4) is The item;

[0039] Let the step size of the solution vector be at this point. ,Right now ,

[0040] Obtain the solution vector The step size is: ;

[0041] Therefore, the solution vector The iterative formula is:

[0042] ,

[0043] Otherwise, when the new row projection direction When it is a zero vector, for the solution vector Perform orthographic projection, at this time The iterative formula is:

[0044] .

[0045] Further, in step (6), the termination condition is: the Euclidean norm of the current error is less than the set value, or the number of iterations reaches the maximum number of iterations.

[0046] This disclosure discloses a CT image reconstruction system based on maximum residual extended oblique projection, used to perform the method, including:

[0047] The data loading module is used to load the CT image matrix and system matrix to be reconstructed, and to set the maximum number of iterations and the error range;

[0048] The matrix splitting module is used to split the CT image matrix to be reconstructed into columns, forming several systems of linear equations;

[0049] The initialization module is used to set the initial values ​​of the auxiliary vector and solution vector of the linear equation system, select the working column and working row according to the maximum residual criterion, and perform the first iteration through the orthogonal projection method to obtain the auxiliary vector and solution vector after the first iteration.

[0050] The column iteration module is used to select the working column based on the maximum residual index, calculate the new column projection direction and column direction step size through two adjacent working columns, and iterate the auxiliary vector.

[0051] The row iteration module is used to select working rows based on the maximum residual index, calculate new row projection directions through two adjacent working rows, and select different methods to iterate the solution vector based on whether the row projection direction is a zero vector.

[0052] The termination judgment module is used to determine whether the current solution vector meets the termination condition, and output the final solution vector or continue iteration;

[0053] The image combination module is used to combine all the final solution vectors obtained by solving the problem column by column to obtain the reconstructed CT image.

[0054] A computer-readable storage medium having a computer program stored thereon, which, when executed by a processor, implements the CT image reconstruction method based on maximum residual extended oblique projection.

[0055] An electronic device, comprising:

[0056] processor;

[0057] Memory, used to store computer programs;

[0058] When the processor executes the computer program, it implements the CT image reconstruction method based on maximum residual extended oblique projection as described in this disclosure.

[0059] The beneficial effects of this disclosure are as follows:

[0060] The method disclosed herein calculates a new projection vector and an appropriate step size in each iteration, and then performs oblique projection on the current solution vector. Compared with traditional orthogonal projection, this method effectively improves imaging efficiency and enhances image reconstruction accuracy. Attached Figure Description

[0061] Figure 1 The images shown are a comparison of CT images reconstructed using the method of this disclosure and existing methods. Image (a) is an ideal target CT image, and images (b), (c), and (d) are images reconstructed using the REK, MREK, and MREKO methods of this disclosure, respectively. The PSNR values ​​corresponding to the three methods are: REK (27.9654dB), MREK (44.5375dB), and MREKO (48.9740dB).

[0062] Figure 2The images shown are a comparison of the CT images reconstructed using the method of this disclosure and existing methods. Image (a) is an ideal target CT image, and images (b), (c), and (d) are images reconstructed using the REK, MREK, and MREKO methods of this disclosure, respectively. The PSNR values ​​corresponding to the three methods are: REK (26.402 dB), MREK (42.3788 dB), and MREKO (45.6958 dB).

[0063] Figure 3 The convergence curves of different methods are shown when the system matrix in this disclosure is rand(10000,2000) (elements ∈ (0,1)).

[0064] Figure 4 The convergence curves of different methods are given when the system matrix is ​​rand(10000,2000) (elements ∈ (0.5,1)).

[0065] Figure 5 The convergence curves of different methods are given when the system matrix is ​​rand(10000,2000) (elements ∈ (0.75,1)).

[0066] Figure 6 This is a system block diagram of the present disclosure. Detailed Implementation

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

[0068] Example: This disclosure describes a CT image reconstruction method based on maximum residual extended oblique projection, comprising the following steps:

[0069] (1) Load the CT image matrix and system matrix to be reconstructed, and set the maximum number of iterations and error range;

[0070] (2) Split the CT image matrix to be reconstructed into columns to form several linear equations;

[0071] (3) Set auxiliary vectors Sum of solutions The initial values ​​are respectively and The working column and working row are selected according to the maximum residual criterion, and the initial value of the auxiliary vector is iterated through the orthogonal projection method. Initial value of the solution vector , This represents the number of iterations.

[0072] (4) Select the working column according to the maximum residual index, calculate the new column projection direction through the working columns of two adjacent steps, and calculate the appropriate step size through the working columns of two adjacent steps. Finally, iterate the auxiliary vector according to the new column projection direction and the obtained appropriate step size. This aims to maximize the elimination of noise in the corresponding columns of the CT image to be reconstructed;

[0073] (5) Select a working row based on the maximum residual index. Calculate the new row projection direction using the working rows from the two adjacent steps. If the calculated new row projection direction is not a zero vector, then use the auxiliary vector obtained from the two adjacent working rows and step (4). Calculate the appropriate step size, and then iterate the solution vector according to the new row projection direction and step size. If the resulting new row projection direction is the zero vector, then the vector is solved iteratively through orthogonal projection. This aims to maximize the optimization of pixels in the corresponding columns of the reconstructed CT image;

[0074] (6) Determine the current solution vector Check if the termination condition is met. If so, change the current solution vector. Output the final solution vector and, starting from step (3), solve the next system of linear equations until the solution is complete; otherwise, let the iteration number be... , and repeat steps (4)-(6);

[0075] (7) Combine all the final solution vectors obtained by solving the problem. By combining columns, the reconstructed CT images are obtained.

[0076] In step (2), the process of reconstructing the CT image matrix is ​​regarded as solving a linear system. The process, in which , , ,in For the system matrix; For the ideal target CT image matrix, It is a known matrix of CT images to be reconstructed; the known matrix of CT images to be reconstructed Split by column Then solve Transform into solving , ,in yes The The column, i.e., the first column of the ideal target CT image matrix. Column pixel values; yes The The column, i.e., the first column of the CT image matrix to be reconstructed. Column pixel values.

[0077] In step (3), an auxiliary vector is set. Sum of solutions The initial values ​​are respectively and However, the two are different in this example.

[0078] Selecting the system matrix based on the maximum residual criterion Work column and work row indicators:

[0079] Selecting working column indicators is done in the system matrix. In all columns, a total Columns, selected to make transpose and auxiliary vector The column index with the largest absolute value of the inner product. And use it as a target for the next step of work. The specific method is as follows:

[0080] ,

[0081] Selecting working row indicators is done in the system matrix. In all rows, a total Line, select the residual term Largest row index And use it as a target for the next step of work. The specific method is as follows:

[0082] ,

[0083] in and They represent The Line and number List, for The Item; using the maximum residual criterion, the auxiliary vector can be made Sum of solutions Fast convergence.

[0084] During the first iteration, orthogonal projection is performed as follows: the initial values ​​of the auxiliary vector are iterated using the orthogonal projection method. Initial value of the solution vector The auxiliary vectors after the first iteration are obtained respectively. and :

[0085] , ,

[0086] in The Euclidean norm of a vector. and They represent and The item, and Representing the system matrix respectively The Column and number OK, It is a system matrix The transpose of a column and For the first selection of column and row indicators;

[0087] That is: during the first iteration, the auxiliary vector It is determined by the initial auxiliary vector To the hyperplane The result obtained by orthogonal projection onto the top;

[0088] Solution vector It is based on the initial solution vector To the hyperplane The result obtained by orthogonal projection;

[0089] Oblique projection cannot be used in the first iteration because oblique projection requires not only the current auxiliary vector Sum of solutions Located on the hyperplane of the previous iteration column (row), a new projection direction needs to be calculated using two adjacent working columns and rows. The auxiliary vector selected during the first iteration... Sum of solutions initial value and Generally, it is not located on the corresponding hyperplane, and only one column and one row are selected, so it is impossible to calculate the new projection direction.

[0090] Specifically, in step (4), the working column is selected according to the maximum residual index, based on the first... The work column selected in the next iteration With the The work column selected in the next iteration A new column projection direction was calculated to adapt to the system matrix characteristics and match the CT projection noise suppression requirements:

[0091] ,

[0092] Among the symbols Represents the dot product of two vectors. The Euclidean norm of a vector.

[0093] And through the obtained new column projection direction vector The first Auxiliary vector for the next iteration Project onto the following hyperplane to eliminate noise:

[0094] ,

[0095] Let the auxiliary vector be at this time. The step size is ,Right now ,

[0096] The step size is obtained as follows: ;

[0097] Therefore, the auxiliary vector The iterative formula is:

[0098] .

[0099] In step (5), the working row is selected according to the maximum residual index. The working row selected in the next iteration With the The working row selected in the next iteration Calculate a new row projection direction to match the pixel reconstruction requirements of the CT image:

[0100] ,

[0101] When the new line projection direction When it is a non-zero vector, the first... Solution vector of the next iteration Project onto the following hyperplane to eliminate artifacts:

[0102] ,

[0103] in, For vectors The item, The auxiliary vector obtained in step (4) is The item;

[0104] At this point, the step size of the solution vector is ,Right now ,

[0105] Obtain the solution vector The step size is: ;

[0106] Therefore, the solution vector The iterative formula is:

[0107] ,

[0108] Otherwise, when the new row projection direction When it is a zero vector, for the solution vector Perform orthographic projection, at this time The iterative formula is:

[0109] .

[0110] In step (5), to ensure the stability of CT image reconstruction, oblique projection can only be used for iteration when the calculated new row projection direction vector is non-zero. This is because the calculated step size needs to be divided by the L2 norm of the corresponding projection vector. If the projection vector is a zero vector, it will result in a division by zero operation. After analysis, the working column selection index ensures that the new column projection direction vector selected in step (4) is non-zero. Therefore, it is only necessary to determine whether the row projection direction obtained in step (5) is non-zero. If it is zero, the solution vector in this step is changed to use orthogonal projection for iteration.

[0111] Determine the current solution vector Check if the termination condition is met. If so, change the current solution vector. Output the final solution vector and, starting from step (3), solve the next system of linear equations until the solution is complete; otherwise, let the iteration number be... Repeat steps (4)-(6).

[0112] The termination condition is: the Euclidean norm of the current error is less than the set value, or the number of iterations reaches the maximum number of iterations.

[0113] The final solution vector of all outputs The images are combined in columns to form the reconstructed CT images.

[0114] This disclosure presents the maximum residual extended oblique projection Kaczmarz method, which is implemented as follows:

[0115] With system matrix CT image matrix to be reconstructed Maximum number of iterations For input; with For output;

[0116] Will Split into column pixel values ​​to be reconstructed. ,make ,

[0117] set up and ,

[0118] Determine work column ,

[0119] calculate: ,

[0120] Determine work line ,

[0121] calculate: ,

[0122] ,

[0123] Determine work column ,

[0124] calculate: ,

[0125] renew: ,

[0126] Determine work line ,

[0127] calculate: ,

[0128] if ,

[0129] renew: ,

[0130] Otherwise, update ,

[0131] Determine if the termination condition is met; if yes, output the result. ,Will Set as the pixel value for the next column, and reset. and No, let the number of iterations be... And repeat the loop.

[0132] In its specific implementation, this disclosure does not entirely employ oblique projection to iteratively solve the vector, but instead introduces a method for... The check for whether it is 0 is because when... If it continues to follow the next oblique projection iteration... This would lead to a division by zero operation, and oblique projection would not be able to update the solution vector. Therefore, orthogonal projection is used for iteration.

[0133] right The judgment is made based on whether it is a zero vector, but not on whether it is a zero vector. Because It is possible to obtain the zero vector, and No. In fact, A zero vector is equivalent to the working row vector selected in the current step. Compared with the working row vector in the previous step Linear correlation. Because... In the hyperplane Above, therefore always Valid. Assume a certain row... and Linear dependence, i.e., the existence of real numbers satisfy We cannot guarantee According to the row indicator selection rules, it is still possible to select those that match... Linearly correlated rows are used as working rows, thus It is possible to obtain the zero vector. And a vector... If and only if and Collinear. Due to ,and Collinear vectors Must be satisfied Therefore, according to the work column selection rules, it is known that... Collinear column vectors will not be included, therefore .

[0134] The maximum residual extended oblique projection (MREKO) method proposed in this disclosure is compared with the existing REK and MREK methods in CT image reconstruction. The CT images used in the following experiments are all from the professional medical image database Radiopaedia. Figure 1 and Figure 2 As shown, (a) is the ideal target CT image, and (b), (c), and (d) are the images reconstructed using the REK, MREK, and MREKO methods of this disclosure, respectively. Peak signal-to-noise ratio (PSNR) is used to measure the similarity between the reconstructed CT image and the accurate image. Figure 1 In the data, the PSNR values ​​corresponding to the three methods are: REK (27.9654dB), MREK (44.5375dB), and MREKO (48.9740dB) of this disclosure; Figure 2 The corresponding PSNR values ​​are: REK (26.4021dB), MREK (42.3788dB), and MREKO (45.6958dB) of this disclosure.

[0135] The comparison results show that the method proposed in this disclosure is superior to the other two comparative methods in terms of the accuracy of CT image reconstruction, demonstrating better CT image reconstruction performance.

[0136] Three sets of examples are given below to illustrate the effectiveness of the MREKO method disclosed herein. All experiments were implemented on MATLAB (R2024b) with an Intel(R) Core(TM) i9-14900HX (2.20 GHz) processor, 32 GB RAM, and Windows 11 operating system. The three sets of system matrices in the experiments were obtained in the following two ways: (1) uniform random matrices generated by the MATLAB function rand(m, n). (2) sparse matrices selected from the SuiteSparse matrix library.

[0137] ,

[0138] Experimental results include the number of iterations (IT) and time consumed (CPU / second) for solving linear systems using relevant methods. The IT and CPU values ​​are the median results after repeating the corresponding calculations 10 times. The performance of the REK method, MREK method, and the MREKO method of this disclosure are compared in the examples. The experimental stopping criteria are as follows:

[0139] ,

[0140] Or the maximum number of iterations reaches ,in This represents the least-norm least-squares solution of the corresponding linear system of equations. If the number of iterations exceeds the maximum number of iterations, it is denoted as "—". Additionally, the speedup ratio is defined. and This is used to record the ratio of MREKO's CPU time to that of REK and MREK, respectively.

[0141] , .

[0142] Tables 1 and 2 show the number of iterations (IT) and computation time (CPU) for the REK, MREK, and MREKO methods in overdetermined and underdetermined linear systems, respectively. As can be seen from the tables, the MREKO method of this disclosure requires fewer iterations and has less computation time than the traditional REK and MREK methods in overdetermined and underdetermined linear systems. In particular, the speedup of the MREKO method of this disclosure increases with the increase in the system matrix size.

[0143] Table 1. Experimental results of different methods for the system matrix rand(m,n) in overdetermined linear systems.

[0144]

[0145] Table 2. Experimental results of different methods for the system matrix rand(m,n) in underdetermined linear systems.

[0146]

[0147] For a sparse matrix selected from the SuiteSparse matrix library, its density is defined as:

[0148] ,

[0149] Table 3. Experimental results of different methods when the system matrix in a linear system is a SuiteSparse sparse matrix.

[0150]

[0151] Table 4 shows the properties of the five sparse matrices selected from the SuiteSparse matrix library; Table 3 shows the number of iterations (IT) and computation time (CPU) for the REK, MREK, and MREKO methods of this disclosure when the system matrix is ​​the corresponding sparse matrix. It can be seen from the tables that the MREKO method of this disclosure has fewer iterations and less computation time than the traditional REK and MREK methods for sparse system matrices with different properties. In particular, the speedup of the MREKO method of this disclosure remains high for positive feedback across system matrices with different properties.

[0152] Table 4 Properties of different sparse matrices

[0153]

[0154] Figure 3-5 The graphs show the convergence curves of three methods when the system matrix is ​​a uniform random matrix of 10000×2000, with all its elements located in the ranges (0,1), (0.5,1), and (0.75,1), respectively. The more concentrated the elements in the system matrix, the stronger the correlation between its rows (or columns), and the more ill-conditioned the matrix. The graphs show that the more concentrated the distribution of the system matrix elements, the greater the impact on the convergence performance of traditional methods (REK, MREK). The MREKO method disclosed in this paper effectively overcomes the interference of matrix element distribution, maintaining high efficiency and strong stability in matrices with varying element concentrations.

[0155] This disclosure applies the proposed maximum residual extended oblique projection Kaczmarz method to CT image reconstruction. Since the core problem of CT reconstruction can be reduced to solving a large-scale incompatible linear system, this method, in the iterative solution process, selects the row (column) with the maximum residual as the working row (column) and uses oblique projection for iteration, thus effectively converging to the solution of the incompatible system. In specific reconstruction experiments, the method is validated by downloading image data from the professional medical image database Radiopaedia, simulating sampling using the system matrix, and then applying the MREKO method of this disclosure to reconstruct CT images. Then, the MREKO method of this disclosure is used iteratively to finally restore a high-quality reconstructed CT image. Compared with traditional REK and MREK methods, the MREKO method of this disclosure exhibits superior computational efficiency, robustness, and effectiveness in solving such incompatible reconstruction problems, giving it significant advantages in both accuracy and speed for CT image reconstruction.

[0156] The parts not described in detail in this application are all existing conventional technologies and will not be elaborated here.

[0157] It is understood that the above specific description of this disclosure is only for illustrating this disclosure and is not limited to the technical solutions described in the embodiments of this disclosure. Those skilled in the art should understand that modifications or equivalent substitutions can still be made to this disclosure to achieve the same technical effect; as long as the usage needs are met, they are all within the protection scope of this disclosure.

Claims

1. A CT image reconstruction method based on maximum residual extended oblique projection, characterized in that: Includes the following steps: (1) Load the CT image matrix and system matrix to be reconstructed, and set the maximum number of iterations and error range; (2) Split the CT image matrix to be reconstructed into columns to form several linear equations; (3) Set the initial values of the auxiliary vector and the solution vector respectively as and , select the working column and the working row according to the maximum residual criterion, and iterate the initial values of the auxiliary vector and the solution vector through the orthogonal projection method; (4) Select the working column according to the maximum residual index, calculate the new column projection direction through the adjacent two working columns, and calculate the appropriate step size through the adjacent two working columns, and iterate the auxiliary vector according to the new column projection direction and the appropriate step size ; Specifically: Select the working column based on the maximum residual index, according to the first... The work column selected in the next iteration With the The work column selected in the next iteration Calculate a new column projection direction vector that matches the projection noise suppression requirements of CT images: ， where the notation represents the inner product of two vectors, represents the Euclidean norm of a vector; and through the new column projection direction vector the first iteration's helper vector is projected to the hyperplane that eliminates the noise: ， Let the step of the auxiliary vector at this time be i.e. , The step size is obtained as: ; Thus the auxiliary vector is given by the iteration formula ； (5) Select the working row according to the maximum residual index. Calculate the new row projection direction using the working rows of the two adjacent steps. If the new row projection direction is not a zero vector, then use the auxiliary vector obtained from the two adjacent working rows and step (4). Calculate the appropriate step size, and iterate the solution vector based on the new row projection direction and the appropriate step size. If the new row projection direction is the zero vector, then the vector is solved iteratively through orthogonal projection. ; (6) Determine the current solution vector Check if the termination condition is met. If so, change the current solution vector. Output the final solution vector and, starting from step (3), solve the next system of linear equations until the solution is complete; otherwise, let the iteration number be... Repeat steps (4)-(6); (7) Combine all the final solution vectors obtained by solving the problem. By combining columns, the reconstructed CT images are obtained.

2. The method according to claim 1, characterized in that: In step (2), the process of reconstructing the CT image matrix is ​​regarded as solving a linear system. The process, in which , , ,in For the system matrix, For the ideal target CT image matrix, It is a known matrix of CT images to be reconstructed; the known matrix of CT images to be reconstructed Split by column Then solve Transform into solving , ,in yes The The column, i.e., the first column of the ideal target CT image matrix. Column pixel values; yes The The column, i.e., the first column of the CT image matrix to be reconstructed. Column pixel values.

3. The method according to claim 1, characterized in that: In steps (3) and (4), the working columns and working rows are selected according to the maximum residual criterion. The selection index for the working columns is in the system matrix. In all columns, a total Columns, selected to make transpose and auxiliary vector The column index with the largest absolute value of the inner product. And use it as a target for the next step of work. ; Selecting working row indicators is done in the system matrix. In all rows, a total Line, select the residual term Largest row index And use it as a target for the next step of work. ;in and They represent The Line and number List, for The item.

4. The method according to claim 1, characterized in that: In step (3), Iterate the initial value of the auxiliary vector using the orthogonal projection method. Initial value of the solution vector The auxiliary vectors after the first iteration are obtained respectively. Sum of solutions Specifically: Auxiliary vector It is determined by the initial auxiliary vector To the hyperplane The result obtained by orthogonal projection onto the top; Solution vector It is based on the initial solution vector To the hyperplane The result obtained by orthogonal projection; in and They represent and The item, and Representing the system matrix respectively The Column and number OK, It is a system matrix The transpose of a column and The column and row indicators selected for the first time.

5. The method according to claim 1, characterized in that: In step (5), the working row is selected according to the maximum residual index, based on the first... The working row selected in the next iteration With the The working row selected in the next iteration Calculate a new row projection direction to match the pixel reconstruction requirements of the CT image: ， When the new line projection direction When it is a non-zero vector, the first... Solution vector of the next iteration Project onto the following hyperplane to eliminate artifacts: ， in, For vectors The item, Auxiliary vector obtained through iteration The item; Let the step size of the solution vector be at this point. ,Right now , Obtain the solution vector The step size is: ; Therefore, the solution vector The iterative formula is: ， Otherwise, when the new row projection direction When it is a zero vector, for the solution vector Perform orthographic projection, at this time The iterative formula is: 。 6. The method according to claim 1, characterized in that: In step (6), the termination condition is: the Euclidean norm of the current error is less than the set value, or the number of iterations reaches the maximum number of iterations.

7. A CT image reconstruction system based on maximum residual extended oblique projection, used to perform the method according to any one of claims 1-6, characterized in that: include: The data loading module is used to load the CT image matrix and system matrix to be reconstructed, and to set the maximum number of iterations and the error range; The matrix splitting module is used to split the CT image matrix to be reconstructed into columns, forming several systems of linear equations; The initialization module is used to set the initial values ​​of the auxiliary vector and solution vector of the linear equation system, select the working column and working row according to the maximum residual criterion, and perform the first iteration through the orthogonal projection method to obtain the auxiliary vector and solution vector after the first iteration. The column iteration module is used to select the working column based on the maximum residual index, calculate the new column projection direction and column direction step size through two adjacent working columns, and iterate the auxiliary vector. The row iteration module is used to select working rows based on the maximum residual index, calculate new row projection directions through two adjacent working rows, and select different methods to iterate the solution vector based on whether the row projection direction is a zero vector. The termination judgment module is used to determine whether the current solution vector meets the termination condition, and output the final solution vector or continue iteration; The image combination module is used to combine all the final solution vectors obtained by solving the problem column by column to obtain the reconstructed CT image.

8. A computer-readable storage medium having a computer program stored thereon, characterized in that, When the computer program is executed by the processor, it implements the CT image reconstruction method based on maximum residual extended oblique projection as described in any one of claims 1-6.

9. An electronic device, characterized in that, include: processor; Memory, used to store computer programs; When the processor executes the computer program, it implements the CT image reconstruction method based on maximum residual extended oblique projection as described in any one of claims 1-6.

Images