A spectral CT image reconstruction method and storage medium
By decomposing the spectral CT image reconstruction problem into low-dimensional sub-problems and introducing prior information and constraints, and utilizing the properties of the Jacobian matrix for transformation and optimization, the reconstruction accuracy and efficiency problems caused by noise and sparse data are solved, and efficient and accurate reconstruction of spectral CT images at low doses is achieved.
Patent Information
- Application Number
- CN202510919445.8
- Authority / Receiving Office
- CN · China
- Patent Type
- Patents(China)
- Current Assignee / Owner
- Filing Date
- 2025-07-04
- Publication Date
- 2025-09-23
- Estimated Expiration
- 2045-07-04
AI Technical Summary
During spectral CT image reconstruction, reducing the tube voltage or tube current will result in more noise or reducing the sampling angle will result in sparse projection data, which will affect the image reconstruction accuracy and efficiency.
The high-dimensional spectral CT image reconstruction problem is decomposed into low-dimensional sub-problems, and prior information and constraints are introduced for each sub-problem. The reconstruction efficiency and accuracy are improved through parallel computing. The special properties of the Jacobian matrix are used for transformation and optimization, and image reconstruction is performed by combining the solution algorithms of the optimization problem and the sub-problems.
Efficient and accurate reconstruction of spectral CT images is achieved under low-dose conditions, improving reconstruction accuracy and efficiency.
Smart Images

Figure CN120411295B_ABST
Abstract
Description
Technical Field
[0001] The embodiments of the present invention relate to the technical field of medical image processing, and in particular to a spectral CT image reconstruction method and a storage medium. Background Art
[0002] Compared with traditional computed tomography (CT), spectral CT, as an advanced imaging technology, has obvious advantages in suppressing artifacts and distinguishing different substances. Therefore, it has been widely used in medical imaging and industrial imaging.
[0003] In spectral CT, the target object is scanned using X-rays of varying energy spectra. The detected projection data is then used to reconstruct a spectral CT image, or a base material decomposition image (hereinafter referred to as the base image). Furthermore, in clinical applications, X-ray dose is often reduced by lowering the tube voltage or current, or by reducing the sampling angle, to ensure the safety of the scanned object.
[0004] However, reducing the tube voltage or current can lead to more noise in the projection data, while reducing the sampling angle can lead to sparser projection data. Practical experience has shown that noisy or sparse projection data can affect the reconstruction accuracy and efficiency of spectral CT images, and this issue needs to be addressed urgently. Summary of the Invention
[0005] The embodiments of the present invention provide a spectral CT image reconstruction method and storage medium, which solve the problem of low reconstruction accuracy and reconstruction efficiency of spectral CT images, especially solve the problem of low reconstruction accuracy and reconstruction efficiency of low-dose spectral CT images.
[0006] According to one aspect of the present invention, a method for reconstructing a spectral CT image is provided, comprising:
[0007] During the current iterative reconstruction process of the energy spectrum CT image, for each energy spectrum of at least two energy spectra, the transformed projection data is determined based on the projection matrix and original projection data corresponding to the energy spectrum, as well as the previous basis image reconstructed in the previous iterative reconstruction process; the transformed projection data and the projection matrix are substituted into a sub-problem pre-constructed for the energy spectrum, and the substituted sub-problem is solved based on the first constraint preset for the sub-problem to obtain a transformed basis image that satisfies the first constraint and has the smallest error with the transformed projection data; the transformed basis images corresponding to each energy spectrum are combined to obtain a combined basis image, and the combined basis image is substituted into the pre-constructed optimization problem to reconstruct a current basis image with the smallest error with the combined basis image, thereby realizing the reconstruction of the energy spectrum CT image.
[0008] According to another aspect of the present invention, a computer-readable storage medium is provided, on which computer instructions are stored. The computer instructions are used to enable a processor to implement the energy spectrum CT image reconstruction method provided by any embodiment of the present invention when executed.
[0009] The technical solution of the embodiment of the present invention can efficiently and accurately reconstruct spectral CT images, especially in low-dose situations.
[0010] It should be understood that the content described in this section is not intended to identify the key or important features of the embodiments of the present invention, nor is it intended to limit the scope of the present invention. Other features of the present invention will become readily understood through the following description. BRIEF DESCRIPTION OF THE DRAWINGS
[0011] In order to more clearly illustrate the technical solutions in the embodiments of the present invention, the following briefly introduces the drawings required for use in the description of the embodiments. Obviously, the drawings described below are only some embodiments of the present invention. For ordinary technicians in this field, other drawings can be obtained based on these drawings without creative work.
[0012] Figure 1 This is a flow chart of a spectral CT image reconstruction method provided according to an embodiment of the present invention;
[0013] Figure 2 is a flowchart of another spectral CT image reconstruction method provided according to an embodiment of the present invention;
[0014] Figure 3 This is a flowchart of another spectral CT image reconstruction method provided according to an embodiment of the present invention;
[0015] Figure 4a This is a graph showing how the relative error of a base image changes with the number of iterations under sparse noise-free data in another spectral CT image reconstruction method provided by an embodiment of the present invention, corresponding to Experiment 1;
[0016] Figure 4b This is a graph showing how the relative error of projection data changes with the number of iterations under sparse noise-free data in another spectral CT image reconstruction method provided by an embodiment of the present invention, corresponding to Experiment 1;
[0017] Figure 5 1 is a schematic diagram of a head phantom image reconstructed using different algorithms under sparse noise-free data in another spectral CT image reconstruction method provided by an embodiment of the present invention, corresponding to Experiment 1;
[0018] Figure 6aThis is a graph showing how the relative error of a base image changes with the number of iterations under sparse noise-free data in another spectral CT image reconstruction method provided by an embodiment of the present invention, corresponding to Experiment 2;
[0019] Figure 6b This is a graph showing how the relative error of projection data changes with the number of iterations under sparse noise-free data in another spectral CT image reconstruction method provided by an embodiment of the present invention, corresponding to Experiment 2;
[0020] Figure 7 Schematic diagram of XCAT phantom images reconstructed by different algorithms under sparse noise-free data in another spectral CT image reconstruction method provided by an embodiment of the present invention, corresponding to Experiment 2;
[0021] Figure 8a This is a graph showing how the relative error of a base image changes with the number of iterations under sparse noisy data in another spectral CT image reconstruction method provided by an embodiment of the present invention, corresponding to Experiment 3;
[0022] Figure 8b This is a graph showing how the relative error of projection data changes with the number of iterations under sparse noisy data in another spectral CT image reconstruction method provided by an embodiment of the present invention, corresponding to Experiment 3;
[0023] Figure 9a This is a graph showing how the relative error of a base image changes over time under sparse noisy data in another spectral CT image reconstruction method provided by an embodiment of the present invention, corresponding to Experiment 3;
[0024] Figure 9b This is a graph showing how the relative error of projection data changes over time under sparse and noisy data in another spectral CT image reconstruction method provided by an embodiment of the present invention, corresponding to Experiment 3;
[0025] Figure 10a This is a graph showing how the relative error between two adjacent basis images changes with the number of iterations in another spectral CT image reconstruction method provided by an embodiment of the present invention under sparse noisy data, corresponding to Experiment 3;
[0026] Figure 10b This is a graph showing how the relative error between two adjacent projection data changes with the number of iterations under sparse noisy data in another spectral CT image reconstruction method provided by an embodiment of the present invention, corresponding to Experiment 3;
[0027] Figure 11 Schematic diagram of lung CT images reconstructed using different algorithms under sparse noisy data in another spectral CT image reconstruction method provided by an embodiment of the present invention, corresponding to Experiment 3;
[0028] Figure 12 This is a structural block diagram of a spectral CT image reconstruction device provided by an embodiment of the present invention;
[0029] Figure 13 It is a structural block diagram of an electronic device for implementing the energy spectrum CT image reconstruction method according to an embodiment of the present invention. DETAILED DESCRIPTION
[0030] In order to enable those skilled in the art to better understand the solutions of the present invention, the technical solutions in the embodiments of the present invention will be clearly and completely described below in conjunction with the drawings in the embodiments of the present invention. Obviously, the embodiments described are only part of the embodiments of the present invention, not all of the embodiments. Based on the embodiments of the present invention, all other embodiments obtained by ordinary technicians in this field without making creative efforts should fall within the scope of protection of the present invention.
[0031] It should be noted that the terms "first", "second", etc. in the description and claims of the present invention and the above-mentioned drawings are used to distinguish similar objects, and are not necessarily used to describe a specific order or sequence. It should be understood that the data used in this way can be interchangeable where appropriate, so that the embodiments of the present invention described herein can be implemented in an order other than those illustrated or described herein. The situations of "target", "original", etc. are similar and will not be repeated here. In addition, the terms "including" and "having" and any variations thereof are intended to cover non-exclusive inclusions, for example, a process, method, system, product or device that includes a series of steps or units is not necessarily limited to those steps or units clearly listed, but may include other steps or units that are not clearly listed or inherent to these processes, methods, products or devices.
[0032] Before introducing the embodiments of the present invention, we first provide an illustrative explanation of the specific reasons why the application scenarios and related solutions have low reconstruction accuracy and efficiency, so as to better understand the energy spectrum CT image reconstruction solution proposed in the embodiments of the present invention and how to improve reconstruction accuracy and efficiency.
[0033] For example, when X-rays enter a target object, they absorb some of the X-ray energy. By detecting the change in X-ray intensity attenuation before and after penetration, the linear attenuation coefficient distribution of the material within the target object can be calculated, thereby reconstructing a CT image. From a physical perspective, the linear attenuation coefficient of a material is correlated with its composition and the energy of the incident photons.
[0034] To facilitate processing, researchers often need to decompose the linear attenuation coefficient based on physical effects (such as the photoelectric effect and Compton scattering) or material composition. Unlike traditional CT, which treats X-rays as a single-energy photon beam, spectral CT fully considers the energy spectrum of X-rays. The rays emitted by an X-ray tube actually contain photons of multiple energies, and this energy distribution is called the X-ray energy spectrum. Specifically, dual-energy CT uses two X-rays with different energy spectra (referred to as the low-energy spectrum and the high-energy spectrum) to scan the target object.
[0035] Available Representative energy spectrum, where It is The ratio of the number of photons in an energy interval (or energy node) to the total number of photons. The data model of energy spectrum CT can be expressed as:
[0036] ,
[0037] in, is the data noise generation operator, represents all projection data or observation data, represents the full spectrum CT imaging operator, which is defined as:
[0038] ,
[0039] in, represents the number of energy spectra, Indicates the The projection data of all detector units and all projection angles under the energy spectrum, For the The imaging operator of the energy spectrum is:
[0040] .
[0041] According to Beer's law, The energy spectrum Imaging operator corresponding to X-ray It can be modeled as:
[0042] ,
[0043] in, Corresponding to The index of all X-rays under the energy spectrum, The energy spectrum exists X-rays. is the index of the energy node, is the number of energy nodes. is an indicator of the decomposed base material, is the amount of base material. For the The base material is The mass attenuation coefficient at each energy node should be noted that this is a known physical quantity that reflects the absorption ability of a substance to X-rays of a specific energy. For the The energy spectrum X-rays pass through the The path length of a pixel block (or voxel) (e.g., length in 2D imaging, area in 3D imaging). For the Zhang Ji's image pixel values, there are pixel values, then For the The vector corresponding to the Zhang basis image is ,in, represents the transpose of a vector. In physics, Zhang Ji's image represents the The discrete spatial distribution of the base material. Represents the vector formed by concatenating all the basis images. is the discrete X-ray transformation vector.
[0044] In the field of spectral CT image reconstruction, the core issue is based on the physical model , by designing a reconstruction algorithm from the projection data The basis image is accurately solved in Furthermore, in clinical applications, as described above, X-ray dose is often reduced by lowering the tube voltage or current, or by reducing the sampling angle, to ensure the safety of the scanned object. However, the projection data obtained in this way is noisy or sparse, which not only results in low reconstruction accuracy of spectral CT images, but also significantly increases the modeling complexity due to the sparse projection data, resulting in low reconstruction efficiency of spectral CT images.
[0045] On this basis, to address the issues of low reconstruction accuracy and efficiency, particularly low reconstruction accuracy, related solutions have introduced prior information about the base image to be reconstructed. This allows for the reconstruction of a more accurate base image guided by this prior information. However, the original problem to be solved in these reconstruction processes is already high-dimensional. Introducing prior information further increases the dimensionality of the original problem, making it more difficult to solve. This means that while these solutions can improve reconstruction accuracy to a certain extent, the improvement is limited and further reduces reconstruction efficiency.
[0046] In response to this, after careful study, it was found that based on certain characteristics of the original problem, the original problem can be decomposed into two or more independent sub-problems, especially decomposed into sub-problems corresponding to each energy spectrum. Therefore, a reconstruction idea is proposed to decompose the high-dimensional original problem into low-dimensional sub-problems, introduce corresponding prior information for each sub-problem, and reconstruct on this basis. Compared with the high-dimensional original problem, the low-dimensional sub-problems are easier to solve, thereby improving the efficiency and accuracy of solving the overall problem. Furthermore, the sub-problems are structurally independent of each other and are naturally suitable for parallel computing, which can greatly shorten the reconstruction time and meet the clinical demand for fast imaging. Next, this reconstruction idea will be elaborated in detail.
[0047] Before reconstruction, known data may include Projection data of energy spectrum 、 Discrete X-ray transformation matrix of energy spectra 、 Energy spectrum distribution And the used The mass attenuation coefficient of the base material .in, Also called the projection matrix or system matrix.
[0048] On this basis, the X-ray energy spectrum distribution matrix is defined as:
[0049] ,
[0050] And define the mass attenuation coefficient matrix:
[0051] .
[0052] In the reconstruction process, what needs to be solved is the unknown Zhang Ji Image On this basis, it should be noted that represents the number of energy spectra used, Indicates the amount of base material, regardless of and Regardless of the situation, the embodiments of the present invention can achieve the decomposition of the base material.
[0053] Specifically, the embodiment of the present invention first uses the energy spectrum CT imaging operator An important property of: its Jacobian matrix at certain special points The special form of is:
[0054] ,
[0055] in, , indicating that The diagonal matrix formed by is a constant matrix, and It can be calculated in advance, for example, Specifically,
[0056] ,
[0057] in, , It can be defined as a point that satisfies the following conditions:
[0058] , for any ,
[0059] in, is a set of real constants. In particular, when When all are set to 0, .
[0060] Utilizing this property, the embodiment of the present invention considers the following optimization problem (i.e., the original problem described above):
[0061] ,
[0062] in, is the set of constraints corresponding to the prior information. The physical meaning of this optimization problem is: in each iterative reconstruction process, by measuring exist The approximate first-order Taylor expansion at ) and projection data The error between them is gradually updated to approximate the true base image. By minimizing the above formula, Under the premise of , the reconstruction error is effectively reduced and the image quality is improved.
[0063] However, since the base image is located in the space The dimension (i.e. ) is usually high, and it is difficult to directly solve the above optimization problem. Using the Jacobian matrix The special form of , after derivation, the objective function of the above optimization problem can be restated as:
[0064] ,
[0065] here, Since they are transformed based on the corresponding base images, they can be called transformed base images. These transformed base images are concatenated into vectors, which are consistent with the base images There are the following relationships:
[0066] ,
[0067] in, is block matrix multiplication.
[0068] also, in Because it is projecting data It is obtained by transformation based on , so it can be called transformed projection data, which is defined as:
[0069] .
[0070] Further, in the definition and Afterwards, we can consider the following A low-dimensional sub-problem:
[0071] ,
[0072] in, for The corresponding constraint set. The above sub-problems are used to measure and The error between .
[0073] In update Afterwards, according to its , we can consider the following optimization problem to update the base image :
[0074] ,
[0075] in, Represents the base image The constraint set corresponding to the prior information of It can be expressed as a positive definite weight matrix. The physical meaning of this optimization problem is: by minimizing and transform base image The error between them ensures that the two are consistent in physics and mathematics.
[0076] Based on the above derivation process, a spectral CT image reconstruction method described in the following embodiment is proposed.
[0077] Figure 1The following is a flow chart of a spectral CT image reconstruction method provided by an embodiment of the present invention. This embodiment is applicable to spectral CT image reconstruction, particularly in low-dose scenarios. This method can be performed by a spectral CT image reconstruction device provided by an embodiment of the present invention. This device can be implemented using software and / or hardware and can be integrated into an electronic device, such as a user terminal or server.
[0078] See also Figure 1 The method of the embodiment of the present invention specifically includes the following steps:
[0079] S110. During the current iterative reconstruction process of the energy spectrum CT image, for each energy spectrum of the at least two energy spectra, the transformed projection data is determined based on the projection matrix and original projection data corresponding to the energy spectrum, and the previous basis image reconstructed in the previous iterative reconstruction process.
[0080] The energy spectrum CT image is obtained through multiple iterative reconstructions. Here, the iterative reconstruction of the current time is taken as an example. In the iterative reconstruction process of the current time, the transformed projection data and the transformed basis image corresponding to each energy spectrum in the entire energy spectrum are determined in turn. Take the energy spectrum as an example, the specific implementation is as follows:
[0081] Get the The projection matrix corresponding to the energy spectrum and projection data It should be noted that in order to obtain the projection data obtained by subsequent transformation To distinguish, the projection data can be The projection data is called the original projection data and the projection data It is called transformed projection data. And obtain the base image reconstructed in the previous iterative reconstruction process of the current iterative reconstruction process , similarly, in order to reconstruct the base image in the iterative reconstruction process To distinguish, here the base image is called the previous base image and the base image is called the current base image. Based on this, according to the above explanation, the previous base image It can be considered as the base image corresponding to each base material reconstructed in the previous iterative reconstruction process The base image is obtained after stitching, so according to the previous base image The base images corresponding to each base material can be deduced Of course, if the iterative reconstruction process is the first iterative reconstruction process, that is, there is no previous iterative reconstruction process, the preset initial base image As the previous base image application.
[0082] Further, according to the obtained projection matrix and the original projection data , and the previous base image , it can be determined that The transformed projection data corresponding to the energy spectrum , which can be achieved through the projection matrix and the previous base image , for the original projection data Perform the transformation to obtain the transformed projection data , which transforms the projection data It can be understood as the first Transformation basis image under energy spectrum The corresponding projection data.
[0083] For example, here is a method for transforming projection data Determined plan:
[0084] Obtain the projection matrix, original projection data and imaging operator corresponding to the energy spectrum; obtain the previous base image reconstructed in the previous iterative reconstruction process; obtain the second projection data based on the imaging operator and the previous base image, and obtain the third projection data based on the projection matrix and the previous base image; obtain the transformed projection data based on the second projection data, the original projection data and the third projection data. , original projection data , imaging operator and the previous base image Afterwards, based on the imaging operator Process the previous base image , and obtain the second projection data. For example, the two can be multiplied to obtain the second projection data. ; and according to the projection matrix and the previous base image , obtain the third projection data, for example, based on the previous base image Get the base images corresponding to each base material , then according to the projection matrix And each base image , and obtain the third projection data. Further, according to the second projection data, the original projection data and the third projection data, the transformed projection data is obtained. For example, the addition and subtraction results of the three projection data can be used as the transformed projection data. .
[0085] On this basis, in order to better understand the above technical solution, the transformed projection data can be obtained by the following formula: :
[0086] ,
[0087] in, represents the second projection data; represents the third projection data; is a constant in the constant matrix, and its specific meaning has been explained above;
[0088] On this basis, the transformed projection data corresponding to each energy spectrum are , which can be expressed by the following formula:
[0089] .
[0090] S120. Substitute the transformed projection data and the projection matrix into the sub-problem pre-constructed for the energy spectrum, and solve the sub-problem based on the first constraint preset for the sub-problem to obtain a transformed basis image that satisfies the first constraint and has the smallest error with the transformed projection data.
[0091] Among them, here we continue with As an example, a corresponding sub-problem is constructed in advance for the energy spectrum, which is used to measure the transformation basis image After the projection matrix Post-projection and transformation of projection data The error between them can be solved by solving this subproblem and obtaining the projection matrix Post-projection and transformation of projection data The transformation basis image with the smallest error between , that is, the transformed basis image After the projection matrix Post-projection and transformation of projection data The error between them is the smallest. And for the sake of easy distinction, the transformed base image with the smallest error can be Transformation basis image , that is, it is the transformation basis image The optimal solution of .
[0092] In an embodiment of the present invention, optionally, the above sub-problem can be expressed by minimizing the fourth equation, which can be obtained according to the fifth equation, which can express the sixth equation and the transformed projection data. The error between the two, the sixth formula can be obtained by the projection matrix and the transformation basis image to be solved Alternatively, the sub-problems pre-constructed for each energy spectrum can be expressed as follows:
[0093] ,
[0094] in, Represents the sixth formula, Represents the fifth formula, It can be seen from this that the subproblems are independent of each other and have a natural parallel solution structure, which improves the reconstruction efficiency.
[0095] Furthermore, in order to improve the reconstruction accuracy, the first constraint can be pre-set for the sub-problem , the first constraint Can be used to transform the basis image in the subproblem Constraints are imposed to ensure that the final transformed basis image is In accordance with the actual situation. In the embodiment of the present invention, optionally, the first constraint The number can be one, two or more; optionally, the first constraint It can be at least one of the non-negative constraint and the total variation (TV) constraint, among which the non-negative constraint can ensure that the transformation basis image The total variation constraint helps to suppress noise and enhance image edges, further improving reconstruction accuracy. The above contents can be set according to actual needs and are not specifically limited here. On this basis, the first constraint is used here. Including non-negative constraints and total variation constraints as an example, the first constraint is It can be expressed by the following formula:
[0096] ,
[0097] in, represents the gradient operator in discrete form; For the The energy spectrum is preset; Represents the set of all non-negative vectors.
[0098] After obtaining the subproblem and the corresponding first constraint, the transformed projection data can be and the projection matrix Substitute into the subproblem and based on the first constraint Solve the substituted subproblem to obtain the solution that satisfies the first constraint. And after the projection matrix Post-projection and transformation of projection data The transformation basis image with the smallest error between , that is, the transformation basis image is obtained In an embodiment of the present invention, optionally, this solution process can be implemented based on any one of the primal-dual hybrid gradient algorithm, the alternating direction multiplier method, and the semi-smooth Newton method, etc., which can be set according to actual needs and is not specifically limited here.
[0099] On this basis, as an example, combined with the above examples, here is an example of solving a sub-problem:
[0100] For the candidate transformation base images that satisfy the first constraint, the first projection data is calculated based on the candidate transformation base images and the projection matrix, and the error between the first projection data and the transformed projection data is determined to solve the candidate transformation base image corresponding to the minimum error among the errors corresponding to the candidate transformation base images as the transformation base image. Among them, the candidate transformation base image can be understood as satisfying the first constraint. Transformation basis image , the number can be one, two or more, which depends on the actual situation and is not specifically limited here. For each candidate transformation base image, the candidate transformation base image and the projection matrix Calculate the first projection data. For example, according to the above example, the candidate transformation base image and the projection matrix The product result is taken as the first projection data; then, the first projection data and the transformed projection data are determined For example, the error between the first projection data and the transformed projection data The difference between them is taken as the error, and the difference can also be processed to obtain the error, and so on. After obtaining the errors corresponding to each candidate transformation base image, the candidate transformation base image corresponding to the minimum error can be used as the transformation base image The optimal solution of , realizing the transformation of the base image The accurate solution of .
[0101] In addition, it should be noted that in spectral CT, the scanning geometric parameters of X-rays under different energy spectra are often inconsistent or mismatched, which will also affect the reconstruction accuracy of spectral CT images. In the embodiment of the present invention, a corresponding sub-problem is constructed for each energy spectrum, that is, only one projection matrix is involved in each sub-problem. , to transform each projection matrix The separate processing solves the problem of low reconstruction accuracy caused by inconsistent or mismatched scanning geometric parameters.
[0102] S130. Combine the transformed basis images corresponding to each energy spectrum to obtain a combined basis image, and substitute the combined basis image into a pre-constructed optimization problem to reconstruct a current basis image with the smallest error between the current basis image and the combined basis image, thereby realizing the reconstruction of the energy spectrum CT image.
[0103] Among them, after solving the transformation basis images corresponding to each energy spectrum Afterwards, these transformed base images can be Combine and get the combined base image , that is, through The combined basis image is obtained by the iterative process. It can be expressed as .
[0104] The optimization problem can be understood as a pre-built method for measuring the combination of base images and the basis image to be solved The error problem between the two is solved by combining the basis images The base image to be solved with the smallest error between , that is, solving the current basis image , the current base image is the basis image to be solved In the embodiment of the present invention, the optimization problem can be expressed in a variety of ways. Here is an example of expression: the optimization problem can be expressed by minimizing the first formula; the first formula is at least expressed by the weight matrix The second formula is obtained by weighted norm, such as directly obtained by weighted norm, or obtained by secondary processing based on the weighted norm. The weight matrix It can be obtained by default, for example, the unit matrix can be taken or when Q=D, it can be taken This can be set according to actual needs and is not specifically limited here. is a constant matrix; the second formula can express the third formula and the combined basis image The error between the three, the third formula is obtained by the preset constant matrix and the basis image to be solved Optionally, the optimization problem can be expressed as follows:
[0105] ,
[0106] in, Represents the first formula, Represents the second formula, Represents the third formula, Can be defined as , Indicates the The transpose of the vector of the basis image of the basis material, represents block matrix multiplication, Represents the basis image to be solved The second constraint that should be satisfied will be described in detail later.
[0107] Further optionally, for the above constant matrix , assuming that the above energy spectrum CT image reconstruction method involves Energy spectrum, Energy nodes and Base materials, of which 、 and are all integers greater than 1, then the constant matrix According to the X-ray energy spectrum distribution matrix And the mass attenuation coefficient matrix Get, for example Among them, the X-ray energy spectrum distribution matrix May include The ratio represents the ratio of the number of photons at the corresponding energy spectrum and energy node to the total number of photons. For example, it means that Energy spectrum The ratio of the number of photons at each energy node to the total number of photons. include The mass attenuation coefficient represents the mass attenuation coefficient of the corresponding base material at the corresponding energy node. Here, the mass attenuation coefficient is For example, it means The base material is The mass attenuation coefficient at each energy node.
[0108] After obtaining the optimization problem, the combined base image Substitute into the optimization problem to reconstruct the combined base image The current base image with the smallest error between , especially the reconstruction of the constant matrix Transformed and combined base images The current base image with the smallest error between , that is, the current base image obtained After the constant matrix Transformed and combined base images The error between them is the smallest. To put it another way, we can update the previous base image by solving the optimization problem. , get the current base image , so as to realize the reconstruction of the energy spectrum CT image. For example, after obtaining the current base image Afterwards, if the current number of iterations has reached the maximum number of iterations or the current base image With the previous base image The error between When the error is less than the preset threshold, the current base image can be output. , the energy spectrum CT image reconstruction is completed; otherwise, the above steps can be repeated for iterative reconstruction to achieve energy spectrum CT image reconstruction. On this basis, optionally, the above error The current base image With the previous base image The difference between the two can also be expressed by the formula This can be set according to actual needs and is not specifically limited here.
[0109] The technical solution of the embodiment of the present invention is that in the iterative reconstruction process of the energy spectrum CT image, for each energy spectrum of at least two energy spectra, the projection matrix corresponding to the energy spectrum is used. and the original projection data , and the previous base image reconstructed in the previous iterative reconstruction process , determine the transformed projection data ; Transform the projection data and the projection matrix Substitute into the subproblem pre-constructed for the energy spectrum, and based on the first constraint preset for the subproblem Solve the substituted subproblem to obtain the solution that satisfies the first constraint And after the projection matrix Post-projection and transformation of projection data The transformation basis image with the smallest error between Then, the transformation basis images corresponding to each energy spectrum are Combine and get the combined base image , and the combined base image Substitute into the pre-built optimization problem to reconstruct the constant matrix Transformed and combined base images The current base image with the smallest error between , to achieve the reconstruction of spectral CT images. The above technical solution, by defining the transformation basis image , decompose the high-dimensional original problem into low-dimensional sub-problems, and set the corresponding first constraint for each sub-problem (i.e., prior information or regularization information), thus by utilizing the first constraint On the basis of improving the reconstruction accuracy, the original problem can be decomposed to reduce the problem dimension, and a low-dimensional, especially a low-dimensional and parallel-solvable sub-problem can be obtained. By solving the sub-problem, energy spectrum CT image reconstructed, thereby improving the reconstruction efficiency and further improving the reconstruction accuracy, achieving the effect of efficient and accurate reconstruction, especially in low-dose situations.
[0110] An optional technical solution, solving the substituted subproblem based on the first constraint preset for the subproblem, includes:
[0111] Obtain a preset primal-dual hybrid gradient algorithm, and use the primal-dual hybrid gradient algorithm to solve the subproblem after substitution based on the first constraint preset for the subproblem.
[0112] In this technical solution, the subproblem can be solved based on the primal-dual hybrid gradient algorithm. sub-problem (i.e. subproblem) as an example, the solution steps are as follows:
[0113] ,
[0114] Among them, given the step size parameter and extrapolation parameters ; Given a regularization parameter ; Let the iteration count index in the subproblem be ,Right now represents the number of iterations of the subproblem, Indicates assignment; initialization ,like Based on its The relationship is initialized, and You can use The dual variables recorded after the last iteration of the subproblem As the initial value, You can use The dual variables recorded after the last iteration of the subproblem As the initial value, Variables used to solve subproblems are continuously updated during iterations.
[0115] It can be performed using either a cold start or a hot start. The core idea of the hot start method is to use the results or related information of the previous iteration as the initial value for the current subproblem, thereby accelerating convergence and improving computational efficiency. The cold start method, on the other hand, does not rely on any existing information at all and reinitializes the variables each time a subproblem is solved.
[0116] Specifically, given the initial dual variables , hot start initialization corresponds to:
[0117] .
[0118] Cold start initialization corresponds to:
[0119] .
[0120] Sub-step 1: Update
[0121] ,
[0122] in, Indicates that the vector Orthogonal projection to the set of all non-negative vectors . This is a gradient descent type update, is the gradient information from the data fidelity term, Is the gradient information from the TV constraint term, passed through the dual variable. Finally, the updated vector Projecting onto the non-negative set ensures the non-negativity of pixels in the image.
[0123] Sub-step 2: Let
[0124] ,
[0125] Among them, this is an acceleration technique that uses the difference between the previous two iterations to perform "advance" updates, which helps to accelerate convergence.
[0126] Sub-step 3: Update
[0127] ,
[0128] where this is the dual variable Updates to the data fidelity term "Related. When and When the difference is large, There are usually larger adjustments. is the step size parameter.
[0129] Sub-step 4: Update
[0130] ,
[0131] in, Indicates that the vector Projection to a collection This is the dual variable The update of Related.
[0132] Sub-step 5: If Reach the maximum number of inner iterations, then let:
[0133] ,
[0134] Otherwise, the iteration count is increased by 1, i.e. , continue iterating and jump to sub-step 1.
[0135] The above technical solution achieves efficient and accurate solution of subproblems through the primal-dual hybrid gradient algorithm.
[0136] Figure 2 This is a flowchart of another energy spectrum CT image reconstruction method provided in an embodiment of the present invention. This embodiment is optimized based on the above-mentioned technical solutions. In this embodiment, optionally, the combined base image is substituted into a pre-constructed optimization problem to reconstruct a current base image with the minimum error between the combined base image and the image, including: obtaining a pre-constructed optimization problem and a second constraint preset for the optimization problem; substituting the combined base image into the optimization problem, and solving the substituted optimization problem based on the second constraint to reconstruct a current base image that satisfies the second constraint and has the minimum error between the combined base image and the image. The explanations of the terms that are the same as or corresponding to the above-mentioned embodiments are not repeated here.
[0137] See also Figure 2 The method of this embodiment may specifically include the following steps:
[0138] S210. During the current iterative reconstruction process of the energy spectrum CT image, for each energy spectrum of the at least two energy spectra, the transformed projection data is determined based on the projection matrix and original projection data corresponding to the energy spectrum, and the previous basis image reconstructed in the previous iterative reconstruction process.
[0139] S220. Substitute the transformed projection data and the projection matrix into the sub-problem pre-constructed for the energy spectrum, and solve the sub-problem based on the first constraint preset for the sub-problem to obtain a transformed basis image that satisfies the first constraint and has the smallest error with the transformed projection data.
[0140] Among them, this step can satisfy the first constraint And after the projection matrix Post-projection and transformation of projection data The transformation basis image with the smallest error between .
[0141] S230. Combine the transformed basis images corresponding to the respective energy spectra to obtain a combined basis image.
[0142] S240. Substitute the combined base image into the pre-constructed optimization problem, and solve the optimization problem based on the second constraint preset for the optimization problem to reconstruct a current base image that satisfies the second constraint and has the smallest error with the combined base image, thereby realizing reconstruction of the energy spectrum CT image.
[0143] Among them, in order to further improve the reconstruction accuracy, a second constraint can be preset for the optimization problem , the second constraint Used to solve the basis image in the optimization problem Constraints are made to take advantage of this second constraint Solve the optimization problem to obtain the second constraint And through the constant matrix Transformed and combined base images The current base image with the smallest error .
[0144] In an embodiment of the present invention, optionally, the second constraint It can be at least one of a real vector space, a non-negative constraint, and a total variation constraint, etc. This can be set according to actual needs and is not specifically limited here.
[0145] The technical solution of the embodiment of the present invention is to preset a second constraint for the optimization problem , and then using this second constraint Solving the optimization problem will help to further improve the reconstruction accuracy of spectral CT images.
[0146] An optional technical solution, the optimization problem is expressed by minimizing the first formula;
[0147] The first formula is obtained by performing a weighted norm on the second formula at least through the weight matrix;
[0148] The second equation represents the error between the third equation and the combined basis image;
[0149] The third equation is expressed by a preset constant matrix and a basis image to be solved, wherein the current basis image is a basis image obtained after solving the basis image to be solved.
[0150] The representation scheme of the above optimization problem has been explained above and will not be repeated here.
[0151] On this basis, the second constraint is given and the weight matrix Several optional selection methods and solutions to optimization problems, especially in constant matrices reversible case is given, whereby one can adjust and , flexibly determine the specific update method of the base image, and through Adjust the importance of different base materials or spatial regions to better adapt to different application scenarios and data characteristics.
[0152] Optionally, in the second constraint is in the real vector space , , and the weight matrix When it is a positive definite matrix, the current basis image According to the inverse matrix of the constant matrix and the combined base image get.
[0153] For example (i.e. Example 1), when Take , , and the constant matrix When the column is full rank, no matter For what positive definite matrix, the current basis image Both are:
[0154] ,
[0155] At this time, to solve the basis image Without imposing any prior constraints, the current base image The combined base image can be directly The transformation basis images corresponding to different energy spectra This represents an idealized, direct material decomposition process. The optimization problem can be simplified to the unconstrained linear least squares problem shown in the following formula. The above linear combination can be considered as the standard solution to the unconstrained linear least squares problem. The weight matrix W does not change the solution:
[0156] .
[0157] Optionally, in the second constraint Is a non-negative constraint , , constant matrix Reversible, and the weight matrix is the inverse matrix of the above Transposed matrix and inverse matrix In the case of the product of It can be obtained by performing non-negative projection on the estimated basis image, which is obtained according to the inverse matrix and the combined base image get.
[0158] For example (i.e. Example 2), when Take , , and the constant matrix Reversible, that is
[0159] ,
[0160] Take When , the above optimization problem has the following closed-form solution:
[0161] ,
[0162] in, Represents the estimated base image. The above example 2 can be understood as first performing direct material decomposition with example 1 to obtain a preliminary estimated base image ; Then, estimate the base image for this estimate Apply a non-negative projection and force all calculated negative pixel values to zero. This ensures that the reconstructed current base image In line with physical reality (ie non-negative). Example 2 above uses a specially selected weight matrix This makes the optimization problem equivalent to:
[0163] ,
[0164] The optimal solution to this optimization problem is to convert the unconstrained solution Orthogonal projection to feasible set In general, this can be considered a computationally efficient way to introduce non-negativity constraints.
[0165] On this basis, the non-negative constraints explained above are For further explanation:
[0166] ,
[0167] This constraint reflects the non-negative property of the image pixel value and has a clear physical meaning, that is, in the clinical scenario, the base image to be solved The pixel value of corresponds to a non-negative quantity physically. Therefore, by introducing the non-negative constraint , which not only better conforms to the physical characteristics of the actual problem, but also effectively narrows the scope of the solution space, thereby improving the solution efficiency and reconstruction accuracy.
[0168] In addition, when In order to ensure the existence and uniqueness of the solution to the optimization problem, it is necessary to solve the basis image Additional prior information is introduced. For example, when When for
[0169] ,
[0170] This is a physical assumption that is consistent with reality.
[0171] Optionally, the second constraint includes a non-negativity constraint and the total variation constraint, and the weight matrix When the current basis image is the identity matrix It can be solved based on the preset primal-dual hybrid gradient algorithm.
[0172] For example (i.e. Example 3), when and There can be any size relationship between them and Take
[0173] ,
[0174] When the matrix is the identity matrix, the above optimization problem has no closed-form solution. Here, it can be solved based on the primal-dual hybrid gradient algorithm. Of course, it can also be solved based on the alternating direction multiplier method or the semi-smooth Newton method. There is no specific limitation here. Here, we seek a set of non-negative basis images. These basis images should not only be able to correspond well to the transformed basis images Data fidelity items:
[0175] ,
[0176] Moreover, each basis image itself should have piecewise constant properties (i.e., satisfy the total variation constraint) to suppress noise and protect edges.
[0177] On this basis, the total variation constraint described above is further explained:
[0178] ,
[0179] in, represents the discrete form of the gradient operator, represents the 1-norm of the vector, represents the regularization parameter, which is preset for each energy spectrum and controls the strength of the total variation constraint. The above total variation constraint corresponds to the basis image to be solved The piecewise constant characteristic, that is, in the basis image to be solved Most areas in the image are characterized by piecewise constant changes, and large gradient changes are allowed at the boundaries or edges. In practical applications, the total variation constraint can effectively characterize the basis image to be solved. The spatial distribution characteristics of are particularly suitable for describing images with sparse gradients (such as tissue structures in medical images). By limiting the 1-norm of the gradient, the total variation constraint can suppress the impact of noise on the reconstruction results while preserving important edge information in the image.
[0180] Furthermore, the implementation process of solving the optimization problem based on the primal-dual hybrid gradient algorithm described above is explained:
[0181] Given a step size parameter , and the extrapolated parameters . Initialize the dual variable When solving the optimization problem here, you can also refer to the two initialization methods of hot start or cold start described above, and will not repeat them again. Sub-step 1: Update :
[0182] ,
[0183] Among them, in the current dual variable Under given conditions, update the base image by gradient descent , and immediately projected onto the non-negative set.
[0184] Sub-step 2: Let
[0185] .
[0186] Sub-step 3: Update
[0187] ,
[0188] Among them, according to the extrapolated value of the current base image With the transformed base image The residual between , to adjust the dual variable .
[0189] Sub-step 4: Update
[0190] ,
[0191] Among them, according to the extrapolated value of the current base image gradient To adjust the dual variables associated with the total variation constraints .
[0192] Sub-step 5: If the maximum number of inner iterations is reached, then output ; Otherwise, the iteration count is increased by 1, i.e. , continue iterating and jump to sub-step 1.
[0193] The above examples achieve efficient and accurate solutions to optimization problems.
[0194] Figure 3 This is a flow chart of another spectral CT image reconstruction method provided in an embodiment of the present invention. This embodiment is optimized based on the above-mentioned technical solutions. Explanations of terms that are identical or corresponding to those in the above-mentioned embodiments are not repeated here.
[0195] See also Figure 3 The method of this embodiment may specifically include the following steps:
[0196] S310. Obtaining an initial base image ,and The original projection data corresponding to the energy spectrum , X-ray energy spectrum distribution matrix and the mass attenuation coefficient matrix , and let the constant matrix and the current iteration number .
[0197] S320. Calculate the transformed projection data corresponding to each energy spectrum , that is, we get .
[0198] in, It is to construct the corresponding sub-problem in S330 below The projection data.
[0199] S330. Parallel solution subproblems to solve the predefined transformation basis image , and after all sub-problems are solved, let , get the transformed base image .
[0200] S340. Update the previous base image by solving the optimization problem , get the current base image .
[0201] S350. If the current number of iterations The maximum number of iterations or the current base image has been reached With the previous base image The error between If the error is less than the preset threshold, the current base image is output. ; Otherwise, the current iteration number Add 1, that is , and jump to S320 to continue iteration.
[0202] The technical solution of the embodiment of the present invention achieves the effect of accurate and efficient reconstruction of low-dose spectral CT images.
[0203] The above technical solution has been verified by experiments. Under the condition of sparse and noise-free data, the method described in the embodiment of the present invention can quickly converge to a result that is highly close to the real image (see Figure 4a 、 Figure 4b 、 Figure 6a and Figure 6b ), there is no visual difference (see Figure 5 and Figure 7 ); significantly outperformed the results of related schemes (i.e., spectral CT image reconstruction using a non-convex primal-dual algorithm) in numerical metrics such as relative error (see Tables 1 and 2 for details). In particular, the relative error between the reconstructed image and the true image achieved by the method described in this embodiment of the present invention approaches machine precision.
[0204] In addition, under the condition of sparse and high noise data, the method described in the embodiment of the present invention can quickly converge to a result that is closer to the real image (see Figure 8a 、 Figure 8b 、 Figure 9a 、 Figure 9b 、 Figure 10a and Figure 10b The obtained image is highly consistent with the real image in terms of details (see Figure 11 ), and has strong robustness to noise and is not easily affected by noise; it is significantly better than the results of related schemes in terms of numerical indicators such as relative error (see Table 3 for details).
[0205] In addition, the method described in the embodiment of the present invention can be considered as an algorithm framework, and its specific implementation module can flexibly adopt a variety of design schemes, has a wide range of applications, and has good scalability. and the weight matrix Different forms can be selected according to needs; the number of inner iterations of the subproblems can also be flexibly adjusted. When the number of inner iterations is large, different subproblems can be solved simultaneously in parallel, which can improve the efficiency of the overall algorithm.
[0206] Here, a dual-energy CT simulation experiment under sparse projection angles is used to verify the advantages of the method described in the embodiment of the present invention. In order to make a fair comparison with related solutions, considering that the method described in the embodiment of the present invention needs to solve the sub-problem through internal iteration, while the related solutions do not set internal iteration, the number of internal iterations for the sub-problem in the method described in the embodiment of the present invention can be set to 1. In fact, the method described in the embodiment of the present invention allows for flexible setting of different numbers of internal iterations and can effectively reconstruct the base image under different settings. No additional experiments will be designed to illustrate this point. It should be pointed out that in the method described in the embodiment of the present invention, multiple independent sub-problems have a natural parallel solution structure. When the number of internal iterations is large, parallel computing can be used to accelerate the solution to further improve the computational efficiency. In addition, the running time is related to the number of iterations, while the memory usage is not affected by the number of iterations.
[0207] Furthermore, the methods described in the embodiments of the present invention contain some hyperparameters. For example, the total variation regularization parameter (TVR) represents the total variation of the base image being solved. In practice, appropriate parameter ranges can be empirically set based on statistical analysis of similar image datasets, based on the target imaging task (e.g., organ location such as the head or chest) and the imaging device type.
[0208] Here we define the relative error about the basis image and the relative error about the projection data:
[0209] ,
[0210] in, is a real image; is the basis image determined after k iterations, such as each basis image and the virtual monoenergetic image at the corresponding energy; Represents the 2-norm of the vector. When reconstructing noisy data, in order to verify the numerical convergence of the method described in the embodiment of the present invention, the relative error of the basis images between two adjacent steps (i.e., the residual described above) and the relative error of the projection data between two adjacent steps can be considered, which are defined as:
[0211] .
[0212] In addition, the following quantitative indicators can be used to measure the quality of the reconstructed image. Given two image vectors and , their normalized root mean-squared error (NRMSE) is defined as:
[0213] ,
[0214] And structural similarity (SSIM) is defined as:
[0215] ,
[0216] in, correspond No. indivual window, correspond No. indivual window, represents the number of windows, then the local similarity index is defined as:
[0217] ,
[0218] in, express The mean of express The variance of express The mean of express The variance of express and In addition, The peak signal-to-noise ratio (PSNR) is defined as:
[0219] .
[0220] On this basis, the following experiments 1, 2 and 3 were carried out (the following experiments are based on As an example):
[0221] Experiment 1: Reconstructing a head phantom using sparse noise-free data
[0222] Experiment 1 uses a head phantom as the real image. The head phantom is made of water-based materials and bone-based materials, distributed in an area of [-15cm, 15cm] × [-15cm, 15cm], and the image size is 256 × 256. Figure 5As shown in the first row. The experiment used a low-energy spectrum of 80 kVp and a high-energy spectrum of 140 kVp. The X-rays from the high-energy spectrum were processed through a 1 mm thick copper filter. Specifically, the low-energy spectrum used 30 projection angles evenly distributed between 0 and 180 degrees, with each projection angle setting 301 X-ray beams evenly spaced within the detection range of [-15 cm, 15 cm]. The high-energy spectrum maintained the same number of rays, but its 30 projection angles were evenly distributed between 3 and 183 degrees. The energy spectrum distribution, mass attenuation coefficient data, and real images were substituted into the imaging model to generate noise-free dual-energy projection data.
[0223] In this experiment, the initial basis image is set to zero and all dual variables are also initialized to zero, i.e.
[0224] ,
[0225] Pick Identity matrix, regularization parameter are all taken as the total variation value of the real transformation basis image, that is,
[0226] ,
[0227] Step size parameter Both are taken as 0.08, and the extrapolated parameters is taken as 1.0. In addition, in the non-convex primal-dual algorithm in the related scheme, the initial basis image is set to zero, all dual variables are also initialized to zero, and the parameter of the total variation regularization term is taken as the total variation value of the real image, that is, , the step size parameter is also taken as 0.08, and the extrapolation parameter is also taken as 1.0.
[0228] In the experiment, the method described in the embodiment of the present invention and the non-convex primal-dual algorithm in the related scheme were iterated 200,000 times, and the relative error of the base image and the relative error of the projection data were plotted as a function of the number of iterations. Figure 4a and Figure 4b As shown in the graph of relative error versus iteration number, it can be seen that during the iterative reconstruction process of the method described in the embodiment of the present invention, both the relative error of the base image and the relative error of the projection data decrease significantly with increasing iteration number. In contrast, although the relative error of the base image and the relative error of the projection data of the related scheme also decrease with increasing iteration number, the rate of error decrease is much lower than that of the method described in the embodiment of the present invention. Ultimately, the relative error level achieved by the method described in the embodiment of the present invention is significantly lower than that of the related scheme. This fully demonstrates the advantages of the method described in the embodiment of the present invention in speed and accuracy.
[0229] Furthermore, we present the reconstruction results of the method and related solutions described in the embodiment of the present invention in Figure 5Among them, the first line is the real image, the second line is the reconstruction result of the method described in the embodiment of the present invention, and the third line is the reconstruction result of the related scheme; from left to right, they are water-based image, bone-based image, 60keV virtual monoenergetic image and 100keV virtual monoenergetic image. Figure 5 It can be seen that the image reconstructed by the method described in this embodiment of the present invention is almost completely consistent with the real image. However, in the reconstruction results of related solutions, the base image has obvious stripe artifacts caused by undersampling, and the virtual monoenergetic image is blurred, resulting in poor overall image quality. This comparison further verifies the advantages of the method described in this embodiment of the present invention in terms of reconstruction accuracy and image quality.
[0230] Table 1 shows the evaluation metrics provided by the method described in the embodiment of the present invention. The cells in Table 1, from top to bottom, represent the normalized root mean square error, structural similarity, and peak signal-to-noise ratio. As can be seen from Table 1, under the same experimental conditions, for the reconstruction of sparse, noise-free dual-energy CT images of a head phantom of the same scale, the method described in the embodiment of the present invention significantly outperforms the non-convex primal-dual algorithm in related solutions in all numerical metrics.
[0231] Table 1 Quantitative evaluation of image reconstruction quality (head phantom)
[0232]
[0233] Experiment 2: Reconstructing the XCAT Phantom Using Sparse Noise-Free Data
[0234] Experiment 2 uses the XCAT phantom to construct a digital phantom as the real image. The phantom is made of water-based materials and bone-based materials, distributed in an area of [-15cm, 15cm] × [-15cm, 15cm], and the image size is 256 × 256. Figure 7 As shown in the first row. The imaging system uses an 80kVp low-energy spectrum and a 140kVp high-energy spectrum. The X-rays in the high-energy spectrum are processed through a 1mm thick copper filter. Specifically, the low-energy spectrum uses 60 projection angles evenly distributed between 0 and 180 degrees, and each projection angle is equipped with 301 X-ray beams evenly spaced within the detection range of [-15cm, 15cm]. The high-energy spectrum maintains the same number of rays, but the 60 projection angles are evenly distributed between 1.5 and 181.5 degrees. Substituting the energy spectrum distribution, mass attenuation coefficient data, and real images into the imaging model generates noise-free dual-energy projection data.
[0235] In this experiment, the initial basis image is set to zero and all dual variables are also initialized to zero, i.e.
[0236] ,
[0237] Pick Identity matrix, regularization parameter are all taken as the total variation value of the real transformation basis image, that is,
[0238] ,
[0239] Step size parameter Both are taken as 0.08, and the extrapolated parameters are all set to 1.0. In addition, in the non-convex primal-dual algorithm in the related scheme, the initial basis image is set to zero, all dual variables are also initialized to zero, and the parameter of the total variation regularization term is taken as the total variation value of the real image, that is, In all relevant schemes, the step size parameter is taken as 0.1 and the extrapolation parameter is taken as 1.0.
[0240] In the experiment, the method described in the embodiment of the present invention and the non-convex primal-dual algorithm (including the total variation regularization term) in the related scheme were iterated 150,000 times, and the relative error of the base image and the relative error of the projection data were plotted as a function of the number of iterations. Figure 6a and Figure 6b shown.
[0241] As can be seen from the graph of relative error versus iteration number, during the iterative reconstruction process, the relative error of the base image and the relative error of the projection data decrease as the number of iterations increases in the method according to the embodiment of the present invention. In contrast, in the iterative reconstruction process of the related scheme, although the relative error of the base image and the relative error of the projection data also decrease with increasing iterations, the rate of error reduction is much slower than that of the method according to the embodiment of the present invention. The relative error level of the method according to the embodiment of the present invention is significantly lower than that of the related scheme, demonstrating its advantages in speed and accuracy.
[0242] We present the results of the reconstruction of the method and related solutions described in the embodiment of the present invention in Figure 7 , where the first row represents the real image, the second row represents the reconstruction result of the method described in the embodiment of the present invention, and the third row represents the reconstruction result of the related scheme; from left to right they are water-based image, bone-based image, 60keV virtual monoenergetic image and 100keV virtual monoenergetic image. Figure 7 It can be seen that the image reconstructed by the method described in the embodiment of the present invention is almost the same as the real image, while the base image reconstructed by the related scheme has obvious artifacts.
[0243] Table 2 shows the evaluation metrics provided by the method described in the embodiment of the present invention. The cells in Table 2, from top to bottom, represent the normalized root mean square error, structural similarity, and peak signal-to-noise ratio. As can be seen from Table 2, under the same experimental conditions, for the sparse, noise-free dual-energy CT image reconstruction problem of the XCAT phantom of the same scale, the method described in the embodiment of the present invention significantly outperforms the non-convex primal-dual algorithm in related solutions in all numerical metrics.
[0244] Table 2 Quantitative evaluation of image reconstruction quality (XCAT phantom)
[0245]
[0246] Experiment 3: Reconstructing real patient lung CT images using sparse noisy data
[0247] Experiment 3 uses real patient lung CT images from the LIDC-IDRI (The Lung Image Database Consortium and Image Database Resource Initiative) dataset as real images. This image is composed of water-based materials and bone-based materials, distributed in an area of [-15cm, 15cm] × [-15cm, 15cm], and has an image size of 256 × 256. Figure 11 As shown in the first row. The imaging system uses an 80kVp low-energy spectrum and a 140kVp high-energy spectrum. The high-energy spectrum's X-rays are filtered through a 1mm-thick copper filter. Specifically, the low-energy spectrum uses 120 projection angles evenly distributed between 0 and 180 degrees, with each projection angle configured with 301 X-ray beams equally spaced within the detection range of [-15cm, 15cm]. The high-energy spectrum maintains the same number of rays, but the 120 projection angles are evenly distributed between 0.75 and 180.75 degrees. Substituting the energy spectrum distribution, mass attenuation coefficient data, and real images into the imaging model generates noise-free dual-energy projection data. Subsequently, noisy dual-energy projection data is generated by adding Gaussian white noise. In this experiment, the signal-to-noise ratio of the dual-energy projection data is 23.0dB.
[0248] In this experiment, the initial basis image is set to zero and all dual variables are also initialized to zero, i.e.
[0249] .
[0250] In step 4, the weight matrix Take it as the identity matrix, the second constraint Take
[0251] ,
[0252] in, . Step size parameter Both are taken as 0.04, and the extrapolated parameters is taken as 1.0, and 50 iterations are used to update the basis image. In addition, in the non-convex primal-dual algorithm in the related scheme, the initial basis image is set to zero, all dual variables are also initialized to zero, and the parameter of the total variation regularization term is taken as the total variation value of the real image, that is, The step size parameter in the related schemes is also set to 0.04, and the extrapolation parameter is also set to 1.0.
[0253] In this experiment, the method described in the embodiment of the present invention and the non-convex primal-dual algorithm in the related scheme were iterated 100,000 times, and a curve chart showing the relative error of the base image and the relative error of the projection data as a function of the number of iterations was drawn. Figure 8a and Figure 8b As shown in the graph of relative error versus iteration number, it can be seen that during the iterative reconstruction process, the base image relative error and the projection data relative error gradually decrease with increasing iteration number in the method described in the embodiment of the present invention, and converge quickly. In contrast, although the base image relative error and the projection data relative error show a trend of gradually decreasing over time in the early stages of the iteration, after a certain number of iterations, the base image relative error actually begins to increase. In addition, the relative error level of the method described in the embodiment of the present invention is significantly lower than that of the related solutions.
[0254] In addition, a graph showing the relative error of the base image and the relative error of the projection data over time is also drawn, such as Figure 9a and Figure 9b As shown in the figure, it can be seen that within the time range shown, the relative error of the base image and the relative error of the projection data of the method described in the embodiment of the present invention gradually decrease with the running time, and quickly converge to a lower level. Although the relative error of the related scheme also gradually decreases with the running time, the rate of decrease is slower. In addition, the relative error level achieved by the method described in the embodiment of the present invention is significantly lower than that of the related scheme. These results fully demonstrate that the method described in the embodiment of the present invention has obvious advantages in convergence speed and reconstruction accuracy, can achieve high-quality image reconstruction in a shorter time, and effectively suppress noise, showing higher computational efficiency and robustness.
[0255] Figure 10a It shows the change of the relative error of the two adjacent step basis images. Figure 10b The figure shows the change in the relative error between two adjacent projection steps. As can be seen from the figure, the error between two adjacent steps of the method described in this embodiment of the present invention is close to machine precision. In contrast, the non-convex primal-dual algorithm still requires more iterations to reach convergence.
[0256] Figure 11Comparisons of the reconstruction results obtained using the method described in the embodiment of the present invention and related solutions are shown. The first row shows the real image, the second row shows the reconstruction results using the method described in the embodiment of the present invention, and the third row shows the reconstruction results using the related solution. From left to right, they are the water-based image, the bone-based image, the 60keV virtual monoenergetic image, and the 100keV virtual monoenergetic image. The figures clearly show that the image reconstructed using the method described in the embodiment of the present invention is virtually identical to the real image, while the base image reconstructed using the related solution exhibits significant artifacts caused by noise.
[0257] also, Figure 11 A magnified comparison of local details between the real image and the reconstructed result is also provided. These details further demonstrate that the water-based image, bone-based image, 60keV virtual monoenergetic image, and 100keV virtual monoenergetic image reconstructed by the method described in this embodiment of the present invention are highly consistent with the real image in detail, with minimal differences. In contrast, the water-based image, 60keV virtual monoenergetic image, and 100keV virtual monoenergetic image reconstructed by related schemes exhibit significant noise interference, resulting in blurred details of some tissues and excessive noise in the images, making it difficult to distinguish fine structures. These results further demonstrate the advantages of the method described in this embodiment of the present invention in terms of image quality and noise immunity.
[0258] Table 3 shows the evaluation metrics provided by the method described in the embodiment of the present invention. The cells in Table 3, from top to bottom, represent the normalized root mean square error, structural similarity, and peak signal-to-noise ratio. As can be seen from Table 3, under the same experimental conditions, for the reconstruction of sparse and noisy dual-energy CT images of the same scale, the method described in the embodiment of the present invention significantly outperforms the non-convex primal-dual algorithm in related solutions in all numerical metrics.
[0259] Table 3 Quantitative evaluation of image reconstruction quality (lung CT images of real patients)
[0260]
[0261] The above three experiments all prove that the method described in the embodiment of the present invention can efficiently and accurately reconstruct energy spectrum CT images.
[0262] Figure 12 This is a block diagram of the structure of the energy spectrum CT image reconstruction device provided in an embodiment of the present invention. The device is used to execute the energy spectrum CT image reconstruction method provided in any of the above embodiments. The device and the energy spectrum CT image reconstruction method of each of the above embodiments belong to the same inventive concept. For details not fully described in the embodiment of the energy spectrum CT image reconstruction device, please refer to the embodiment of the energy spectrum CT image reconstruction method. Figure 12 The device may specifically include: a transformed projection data determination module 410, a transformed basis image acquisition module 420 and an energy spectrum CT image reconstruction module 430.
[0263] The transformed projection data determining module 410 is configured to determine, during the current iterative reconstruction process of the spectral CT image, for each of the at least two energy spectra, the transformed projection data based on the projection matrix and original projection data corresponding to the energy spectrum, and the previous base image reconstructed in the previous iterative reconstruction process;
[0264] A transformation basis image obtaining module 420 is configured to substitute the transformed projection data and the projection matrix into a pre-constructed sub-problem for the energy spectrum, and solve the substituted sub-problem based on a first constraint preset for the sub-problem to obtain a transformation basis image that satisfies the first constraint and minimizes the error between the transformed projection data and the transformed basis image;
[0265] The energy spectrum CT image reconstruction module 430 is used to combine the transformed basis images corresponding to each energy spectrum to obtain a combined basis image, and substitute the combined basis image into a pre-constructed optimization problem to reconstruct a current basis image with the minimum error between the combined basis image and the current basis image, thereby realizing the reconstruction of the energy spectrum CT image.
[0266] Optionally, the first constraint is used to constrain the transformation basis image to be solved by the subproblem, and the subproblem is solved by the following units:
[0267] The transformation basis image obtaining unit is used to calculate the first projection data according to the candidate transformation basis images and the projection matrix for the candidate transformation basis images that satisfy the first constraint, and determine the error between the first projection data and the transformed projection data, so as to solve the candidate transformation basis image corresponding to the minimum error among the errors corresponding to the candidate transformation basis images as the transformation basis image.
[0268] Optionally, the base image transformation module 420 may include:
[0269] The subproblem solving unit is used to obtain a preset primal-dual hybrid gradient algorithm, and use the primal-dual hybrid gradient algorithm to solve the substituted subproblem based on the first constraint preset for the subproblem.
[0270] Optionally, the first constraint includes a non-negativity constraint and / or a total variation constraint.
[0271] Optionally, the spectral CT image reconstruction module 430 may include:
[0272] A second constraint acquisition unit, configured to acquire a pre-constructed optimization problem and a second constraint preset for the optimization problem;
[0273] The current base image reconstruction unit is used to substitute the combined base image into the optimization problem and solve the optimization problem based on the second constraint to reconstruct the current base image that satisfies the second constraint and has the smallest error with the combined base image.
[0274] Based on any of the above devices, optionally, the optimization problem is expressed by minimizing the first formula;
[0275] The first formula is obtained by performing a weighted norm on the second formula at least through the weight matrix;
[0276] The second equation represents the error between the third equation and the combined basis image;
[0277] The third equation is expressed by a preset constant matrix and a basis image to be solved, wherein the current basis image is a basis image obtained after solving the basis image to be solved.
[0278] On this basis, an optional embodiment of the above device involves Energy spectrum and Base materials, among which and are all integers greater than 1. If the current base image is reconstructed based on the optimization problem and the second constraint preset for the optimization problem, then:
[0279] The second constraint is located in real vector space, , when the constant matrix is invertible and the weight matrix is a positive definite matrix, the current basis image is obtained according to the inverse matrix of the constant matrix and the combined basis image; and / or,
[0280] The second constraint is a non-negative constraint, , when the constant matrix is invertible and the weight matrix is the product of the inverse matrix and the transposed matrix of the inverse matrix, the current basis image is obtained by performing non-negative projection on the estimated basis image, and the estimated basis image is obtained according to the inverse matrix and the combined basis image; and / or,
[0281] When the second constraint includes a non-negative constraint and a total variation constraint, and the weight matrix is an identity matrix, the current basis image is solved based on a preset primal-dual hybrid gradient algorithm.
[0282] Alternatively, the above device involves Energy spectrum, Energy nodes and Base materials, among which 、 and are all integers greater than 1; the constant matrix is obtained according to the X-ray energy spectrum distribution matrix and the mass attenuation coefficient matrix;
[0283] Among them, the X-ray energy spectrum distribution matrix includes A ratio, which represents the ratio of the number of photons at the corresponding energy spectrum and energy node to the total number of photons;
[0284] The mass attenuation coefficient matrix includes A mass attenuation coefficient, which represents the mass attenuation coefficient of the corresponding base material at the corresponding energy node.
[0285] Optionally, the transformed projection data determination module 410 may include:
[0286] An imaging operator acquisition unit, used to acquire a projection matrix, original projection data, and an imaging operator corresponding to the energy spectrum;
[0287] A previous base image acquisition unit, configured to acquire a previous base image reconstructed in a previous iterative reconstruction process;
[0288] a third projection data obtaining unit, configured to obtain second projection data according to the imaging operator and the previous base image, and to obtain third projection data according to the projection matrix and the previous base image;
[0289] The transformed projection data obtaining unit is used to obtain the transformed projection data according to the second projection data, the original projection data and the third projection data.
[0290] The spectral CT image reconstruction device provided in the embodiment of the present invention can achieve efficient and accurate spectral CT image reconstruction through the mutual cooperation of various modules, especially in low-dose situations.
[0291] The spectral CT image reconstruction device provided by the embodiment of the present invention can execute the spectral CT image reconstruction method provided by any embodiment of the present invention, and has the corresponding functional modules and beneficial effects of the execution method.
[0292] It is worth noting that in the embodiment of the above-mentioned energy spectrum CT image reconstruction device, the various units and modules included are only divided according to functional logic, but are not limited to the above-mentioned division, as long as the corresponding functions can be achieved; in addition, the specific names of the various functional units are only for the convenience of distinguishing each other and are not used to limit the scope of protection of the present invention.
[0293] Figure 13 A schematic diagram of an electronic device 10 that can be used to implement an embodiment of the present invention is shown. The electronic device is intended to represent various forms of digital computers, such as laptop computers, desktop computers, workstations, personal digital assistants, servers, blade servers, mainframe computers, and other suitable computers. The electronic device can also represent various forms of mobile devices, such as personal digital assistants, cellular phones, smartphones, wearable devices (such as helmets, glasses, watches, etc.), and other similar computing devices. The components shown herein, their connections and relationships, and their functions are merely examples and are not intended to limit the implementation of the present invention described and / or claimed herein.
[0294] like Figure 13 As shown, electronic device 10 includes at least one processor 11 and memory, such as read-only memory (ROM) 12 and random access memory (RAM) 13, communicatively connected to at least one processor 11. The memory stores computer programs executable by the at least one processor. Processor 11 can perform various appropriate actions and processes based on the computer programs stored in ROM 12 or loaded from storage unit 18 into RAM 13. RAM 13 can also store various programs and data required for the operation of electronic device 10. Processor 11, ROM 12, and RAM 13 are interconnected via bus 14. An input / output (I / O) interface 15 is also connected to bus 14.
[0295] Multiple components in the electronic device 10 are connected to the I / O interface 15, including an input unit 16, such as a keyboard, a mouse, etc.; an output unit 17, such as various types of displays, speakers, etc.; a storage unit 18, such as a magnetic disk, an optical disk, etc.; and a communication unit 19, such as a network card, a modem, a wireless communication transceiver, etc. The communication unit 19 allows the electronic device 10 to exchange information / data with other devices via a computer network such as the Internet and / or various telecommunication networks.
[0296] The processor 11 can be any general-purpose and / or specialized processing component with processing and computing capabilities. Some examples of the processor 11 include, but are not limited to, a central processing unit (CPU), a graphics processing unit (GPU), various specialized artificial intelligence (AI) computing chips, various processors running machine learning model algorithms, a digital signal processor (DSP), and any other suitable processor, controller, microcontroller, etc. The processor 11 executes the various methods and processes described above, such as the spectral CT image reconstruction method.
[0297] In some embodiments, the spectral CT image reconstruction method can be implemented as a computer program tangibly embodied in a computer-readable storage medium, such as storage unit 18. In some embodiments, part or all of the computer program can be loaded and / or installed on electronic device 10 via ROM 12 and / or communication unit 19. When the computer program is loaded into RAM 13 and executed by processor 11, one or more steps of the spectral CT image reconstruction method described above can be performed. Alternatively, in other embodiments, processor 11 can be configured to perform the spectral CT image reconstruction method via any other suitable means (e.g., via firmware).
[0298] Various embodiments of the systems and techniques described herein can be implemented in digital electronic circuit systems, integrated circuit systems, field programmable gate arrays (FPGAs), application specific integrated circuits (ASICs), application specific standard products (ASSPs), systems on chips or systems on chips (SOCs), complex programmable logic devices (CPLDs), computer hardware, firmware, software, and / or combinations thereof. These various embodiments can include being implemented in one or more computer programs that are executable and / or interpreted on a programmable system that includes at least one programmable processor, which can be a special purpose or general purpose programmable processor that can receive data and instructions from a storage system, at least one input device, and at least one output device, and transmit data and instructions to the storage system, the at least one input device, and the at least one output device.
[0299] Computer programs for implementing the methods of the present invention may be written in any combination of one or more programming languages. These computer programs may be provided to a processor of a general-purpose computer, a special-purpose computer, or other programmable data processing device, such that when the computer program is executed by the processor, the functions / operations specified in the flowcharts and / or block diagrams are implemented. The computer program may be executed entirely on the machine, partially on the machine, as a stand-alone software package, partially on the machine and partially on a remote machine, or entirely on a remote machine or server.
[0300] In the context of the present invention, a computer-readable storage medium may be a tangible medium that may contain or store a computer program for use by or in conjunction with an instruction execution system, device, or apparatus. A computer-readable storage medium may include, but is not limited to, an electronic, magnetic, optical, electromagnetic, infrared, or semiconductor system, device, or apparatus, or any suitable combination of the foregoing. Alternatively, a computer-readable storage medium may be a machine-readable signal medium. More specific examples of machine-readable storage media may include an electrical connection based on one or more wires, a portable computer disk, a hard disk, a random access memory (RAM), a read-only memory (ROM), an erasable programmable read-only memory (EPROM or flash memory), an optical fiber, a portable compact disk read-only memory (CD-ROM), an optical storage device, a magnetic storage device, or any suitable combination of the foregoing.
[0301] To provide interaction with a user, the systems and techniques described herein can be implemented on an electronic device that has: a display device (e.g., a CRT (cathode ray tube) or LCD (liquid crystal display) monitor) for displaying information to the user; and a keyboard and pointing device (e.g., a mouse or trackball) through which the user can provide input to the electronic device. Other types of devices can also be used to provide interaction with the user; for example, the feedback provided to the user can be any form of sensory feedback (e.g., visual feedback, auditory feedback, or tactile feedback); and input from the user can be received in any form (including acoustic input, voice input, or tactile input).
[0302] The systems and techniques described herein can be implemented in a computing system that includes back-end components (e.g., as a data server), or a computing system that includes middleware components (e.g., an application server), or a computing system that includes front-end components (e.g., a user computer with a graphical user interface or web browser through which a user can interact with implementations of the systems and techniques described herein), or a computing system that includes any combination of such back-end components, middleware components, or front-end components. The components of the system can be interconnected by any form or medium of digital data communication (e.g., a communication network). Examples of communication networks include: a local area network (LAN), a wide area network (WAN), a blockchain network, and the Internet.
[0303] A computing system may include clients and servers. The clients and servers are typically remote from each other and typically interact via a communication network. This client-server relationship arises through computer programs running on the respective computers, creating a client-server relationship. The server may be a cloud server, also known as a cloud computing server or cloud host. This server is a hosting product within the cloud computing service ecosystem that addresses the management difficulties and limited scalability of traditional physical hosting and VPS services.
[0304] In particular, according to an embodiment of the present invention, the process described above with reference to the flowchart can be implemented as a computer software program. For example, an embodiment of the present invention includes a computer program product that includes a computer program carried on a non-transitory computer-readable medium, the computer program containing program code for executing the method shown in the flowchart. In such an embodiment, the computer program can be downloaded and installed from a network via the communication unit 19, or installed from the storage unit 18, or installed from the ROM 12. When the computer program is executed by the processor 11, the above-mentioned functions defined in the method of the embodiment of the present invention are performed.
[0305] It should be understood that the various forms of the processes shown above can be used to reorder, add, or delete steps. For example, the steps described in the present invention can be performed in parallel, sequentially, or in a different order, as long as the desired results of the technical solution of the present invention can be achieved. This is not limited herein.
[0306] The above specific embodiments do not limit the scope of protection of the present invention. Those skilled in the art will appreciate that various modifications, combinations, sub-combinations, and substitutions may be made based on design requirements and other factors. Any modifications, equivalent substitutions, and improvements made within the spirit and principles of the present invention are intended to be included within the scope of protection of the present invention.
Claims
1. A spectral CT image reconstruction method, characterized in that: include: During the current iterative reconstruction process of the energy spectrum CT image, for each of the at least two energy spectra, determining transformed projection data based on the projection matrix and original projection data corresponding to the energy spectrum, and a previous basis image reconstructed in a previous iterative reconstruction process; Substituting the transformed projection data and the projection matrix into a pre-constructed sub-problem for the energy spectrum, and solving the sub-problem based on a first constraint preset for the sub-problem to obtain a transformed basis image that satisfies the first constraint and has a minimum error with the transformed projection data; Combining the transformed basis images corresponding to the energy spectra to obtain a combined basis image, and substituting the combined basis image into a pre-constructed optimization problem to reconstruct a current basis image with the minimum error between the combined basis image and the current basis image, thereby reconstructing the energy spectrum CT image; The step of determining the transformed projection data according to the projection matrix and the original projection data corresponding to the energy spectrum, and the last base image reconstructed in the last iterative reconstruction process, includes: Acquiring a projection matrix, raw projection data, and an imaging operator corresponding to the energy spectrum; Obtain the last base image reconstructed in the last iterative reconstruction process; Obtaining second projection data according to the imaging operator and the previous basis image, and obtaining third projection data according to the projection matrix and the previous basis image; Transformed projection data is obtained according to the second projection data, the original projection data and the third projection data.
2. The spectral CT image reconstruction method according to claim 1, characterized in that: The first constraint is used to constrain the transformation basis image to be solved by the subproblem, and the subproblem is solved in the following manner: For the candidate transformation basis images that satisfy the first constraint, first projection data is calculated based on the candidate transformation basis images and the projection matrix, and the error between the first projection data and the transformation projection data is determined to solve the candidate transformation basis image corresponding to the smallest error among the errors corresponding to each of the candidate transformation basis images as the transformation basis image.
3. The spectral CT image reconstruction method according to claim 1, characterized in that: Solving the substituted subproblem based on the first constraint preset for the subproblem includes: A preset primal-dual hybrid gradient algorithm is obtained, and the primal-dual hybrid gradient algorithm is used to solve the subproblem after substitution based on the first constraint preset for the subproblem.
4. The spectral CT image reconstruction method according to claim 1, characterized in that: The first constraint includes a non-negativity constraint and / or a total variation constraint.
5. The spectral CT image reconstruction method according to claim 1, characterized in that: Substituting the combined base image into a pre-constructed optimization problem to reconstruct a current base image having a minimum error with the combined base image includes: Obtaining a pre-constructed optimization problem and a second constraint preset for the optimization problem; Substitute the combined base image into the optimization problem, and solve the optimization problem based on the second constraint to reconstruct a current base image that satisfies the second constraint and has the smallest error with the combined base image.
6. The spectral CT image reconstruction method according to any one of claims 1 to 5, characterized in that: The optimization problem is expressed by minimizing the first equation: The first formula is obtained by performing a weighted norm on the second formula at least through a weight matrix; The second equation represents the error between the third equation and the combined basis image; The third equation is represented by a preset constant matrix and a basis image to be solved, wherein the current basis image is a basis image obtained by solving the basis image to be solved.
7. The spectral CT image reconstruction method according to claim 6, characterized in that: The method involves Energy spectrum and Base materials, among which and are all integers greater than 1. If the current base image is reconstructed based on the optimization problem and the second constraint preset for the optimization problem, then: The second constraint is located in the real vector space, , when the constant matrix is invertible and the weight matrix is a positive definite matrix, the current basis image is obtained according to the inverse matrix of the constant matrix and the combined basis image; and / or, The second constraint is a non-negativity constraint, , when the constant matrix is invertible and the weight matrix is the product of the inverse matrix and the transposed matrix of the inverse matrix, the current basis image is obtained by performing non-negative projection on an estimated basis image, and the estimated basis image is obtained according to the inverse matrix and the combined basis image; and / or, When the second constraint includes a non-negative constraint and a total variation constraint, and the weight matrix is an identity matrix, the current basis image is solved based on a preset primal-dual hybrid gradient algorithm.
8. The spectral CT image reconstruction method according to claim 6, characterized in that: The method involves Energy spectrum, Energy nodes and Base materials, among which 、 and are all integers greater than 1; The constant matrix is obtained according to the X-ray energy spectrum distribution matrix and the mass attenuation coefficient matrix; Wherein, the X-ray energy spectrum distribution matrix includes A ratio, wherein the ratio represents the ratio of the number of photons at the corresponding energy spectrum and energy node to the total number of photons; The mass attenuation coefficient matrix includes A mass attenuation coefficient is a mass attenuation coefficient, where the mass attenuation coefficient represents the mass attenuation coefficient of the corresponding base material at the corresponding energy node.
9. A computer-readable storage medium, characterized in that The computer-readable storage medium stores computer instructions, and the computer instructions are used to enable a processor to implement the energy spectrum CT image reconstruction method according to any one of claims 1 to 8 when executed.
Citation Information
Patent Citations
X ray multi-energy-spectrum CT (Computed Tomography) finite angle scanning and image iterative reconstruction method
CN108010099A
Image representation method based on practical robust PCA
CN111340120A