Three-dimensional X-ray low-dose imaging method and device

Abstract

The present invention discloses a three-dimensional X-ray low-dose imaging method and device. The method includes collecting air exposure data and estimating parameters within a three-dimensional X-ray system; using a three-dimensional X-ray system constructed based on the three-dimensional X-ray system parameters to collect original projection data of a patient, and calibrating the original projection data of the patient to obtain an integral image; performing data rearrangement on the integral image to obtain a noise-contaminated low-dose projection chord diagram; estimating an ideal projection chord diagram without noise interference based on the low-dose projection chord diagram using a projection recovery model; performing un-rearrangement on the ideal projection chord diagram to obtain a restored integral image; and reconstructing the restored integral image to obtain a target CBCT tomographic image. The present invention has low noise and high image resolution, can effectively solve the technical pain points of high noise and poor image resolution in current low-dose three-dimensional images, and can be widely used in the field of computed tomography.

Description

Technical Field

[0001] The present invention relates to the field of computer tomography technology, and in particular to a three-dimensional X-ray low-dose imaging method and device. Background Art

[0002] Three-dimensional X-ray imaging, exemplified by cone-beam computed tomography (CBCT), allows for extremely rapid data acquisition and delivers 3D images in a single scan. It has wide applications in clinical fields such as orthodontics, breast tomosynthesis, precise image guidance, and 4D-CBCT for cardiovascular and cerebrovascular examinations. In line with the low-dose principle (As Low As Reasonably Achievable, ALARA), low-dose X-ray imaging methods have become a research hotspot in recent years.

[0003] Low-dose imaging primarily reduces radiation dose by lowering X-ray tube voltage or current, improving detector energy response, reducing radiation duration, and sparse projection sampling. According to relevant research, the most critical impact of low-dose X-rays on imaging is a significant decrease in the signal-to-noise ratio (SNR). This causes some or much tissue information (particularly low-contrast soft tissue information) to be drowned out by noise, rendering diagnostic requirements unattainable. To address this issue, numerous experts and researchers have conducted in-depth research in the field of CT. Strategies for low-dose CT imaging primarily encompass projection data filtering, iterative reconstruction, and post-processing denoising. Projection filtering includes nonlinear filtering, such as the common mean and median filters, and statistical model-based projection filtering, such as projection denoising based on the Poisson distribution model. However, these methods often ignore the effects of beam energy polychromaticity, detector response, and neighborhood correlation, and therefore fail to accurately reflect the noise characteristics of projections. Iterative reconstruction methods, due to their inherent mathematical optimization, exhibit excellent noise immunity and can still reconstruct high-quality images even with relatively few projection data. However, iterative reconstruction suffers from high time complexity, and most current iterative reconstruction algorithms fail to meet the real-time requirements of clinical imaging. Post-processing denoising methods can also significantly reduce the noise in reconstructed images and improve the signal-to-noise ratio of reconstructed images. However, most post-processing methods will reduce the spatial resolution of the image to varying degrees, resulting in blurred edges of tissue structures, and the design of post-processing filters will greatly affect the denoising effect.

[0004] While there has been extensive research on low-dose CT imaging, there has been relatively little research on low-dose CBCT three-dimensional imaging. A key reason for this is that CT detectors are linear arrays, while CBCT detectors are flat-panel detectors. CBCT detectors lack the ability to use grids to filter out scattered radiation, further compromising the already poor low-dose CBCT images. The challenge of three-dimensional low-dose imaging remains a technical pain point that urgently needs to be addressed. Summary of the Invention

[0005] In view of this, an embodiment of the present invention provides a three-dimensional X-ray low-dose imaging method and apparatus with low noise and high image resolution.

[0006] One aspect of the present invention provides a three-dimensional X-ray low-dose imaging method, comprising:

[0007] Collect air exposure data and estimate the internal parameters of the 3D X-ray system;

[0008] The three-dimensional X-ray system constructed according to the three-dimensional X-ray system parameters collects original projection data of the patient, and calibrates the original projection data of the patient to obtain an integral image;

[0009] performing data rearrangement on the integral image to obtain a noise-contaminated low-dose projection chord diagram;

[0010] estimating an ideal projection chord diagram without noise interference by using a projection recovery model according to the low-dose projection chord diagram;

[0011] Rearranging the ideal projection chord diagram to obtain a restored integral image;

[0012] The restored integral image is reconstructed to obtain a target CBCT tomographic image.

[0013] Optionally, collecting air exposure data and estimating internal parameters of the three-dimensional X-ray system includes:

[0014] Configuring geometric parameters of a three-dimensional low-dose imaging device, including the distance from the radiation source to the rotation center, the distance from the rotation center to the detector, and the distance from the laser center to the imaging center;

[0015] Repeatedly collecting multiple air exposure data using the three-dimensional low-dose imaging device at a preset angle;

[0016] Determining a two-dimensional adjustment function matrix and a scale coefficient matrix of the same size as the projection according to the projection sequence in the air exposure data;

[0017] The expressions of the adjustment function and the scaling coefficient are:

[0018]

[0019] in, is the mean of the array p(u,v,:); σ 2 (u,v) is the variance of the array p(u,v,:); g(u,v) is the adjustment function; a n is the scale coefficient; N is the order of the scale coefficient.

[0020] Optionally, the three-dimensional X-ray system constructed according to the three-dimensional X-ray system parameters collects original projection data of the patient, and calibrates the original projection data of the patient to obtain an integral image, including:

[0021] performing afterglow calibration on the original projection data;

[0022] Performing scattering calibration on the afterglow-calibrated data using a Monte Carlo method;

[0023] performing hardening calibration on the scattering-calibrated data;

[0024] An integral image is determined according to the results of the persistence calibration, the scatter calibration, and the hardening calibration.

[0025] Optionally, in the step of rearranging the integral image to obtain a noise-contaminated low-dose projection chord diagram,

[0026] The result of the data rearrangement includes three-dimensional pseudo fan-beam data or three-dimensional pseudo parallel beam data;

[0027] The expression for data rearrangement is:

[0028]

[0029] Where u represents the row of the projection image in the integral image; v represents the column of the projection image in the integral image; s represents the sth projection image in the integral image; d is the detector pixel size; D is the distance from the ray source to the detector; p(s,v,u) represents the integral image data; q(u′,v′,s′) represents the rearranged three-dimensional pseudo fan-beam data.

[0030] Optionally, estimating an ideal projection chord diagram without noise interference by using a projection recovery model according to the low-dose projection chord diagram includes:

[0031] The projection recovery model is described by the maximum a posteriori probability;

[0032] Derivative the projection restoration model to obtain an optimized iterative solution formula;

[0033] Determine the projection restoration model by calculating the product of the target filter window and the image gradient under the target filter window according to the iterative solution formula;

[0034] Calculating an ideal projection chord diagram without noise interference according to the weighted mean of the filter window and the projection restoration model;

[0035] The expression of the projection restoration model is:

[0036]

[0037] in, is the lossless projection image to be restored; q i is the rearranged projected chord diagram; is the one-dimensional representation of the noise model; m is the total number of pixels; β is the relaxation factor; Ω is the neighborhood of pixel i; γ ij is the texture weight function.

[0038] Optionally, in the step of performing de-rearrangement on the ideal projection chord diagram to obtain a restored integral image, an expression for the de-rearrangement is:

[0039]

[0040] Where p(u,v,s) represents the integral image recovered after de-rearrangement; represents the ideal projection chord diagram without noise interference; u represents the row of the projection image in the integral image; v represents the column of the projection image in the integral image; s represents the s-th projection image in the integral image; d is the detector pixel size.

[0041] Optionally, reconstructing the restored integral image to obtain a target CBCT tomographic image includes:

[0042] reconstructing the restored integral image using an iterative reconstruction algorithm that has undergone sharpening modification;

[0043] determining an objective function based on the reconstructed image;

[0044] Wherein, the expression of the objective function is:

[0045]

[0046] Where f represents the reconstructed image obtained by iterative optimization; Pf represents the forward projection of the reconstructed image f; p is the de-permutated integral image; μ is the relaxation adjustment factor; R(f) is the regularization term;

[0047] The reconstructed image is accelerated by GPU parallel computing; the regularization term is constructed through gradient variance; and the iterative reconstruction algorithm is solved by the CGLS conjugate gradient descent algorithm.

[0048] Another aspect of the present invention provides a three-dimensional X-ray low-dose imaging device, comprising:

[0049] The first module is used to collect air exposure data and estimate the internal parameters of the 3D X-ray system;

[0050] A second module is configured to acquire original projection data of a patient using a three-dimensional X-ray system constructed according to the three-dimensional X-ray system parameters, and to calibrate the original projection data of the patient to obtain an integral image;

[0051] A third module is used to rearrange the data of the integral image to obtain a noise-contaminated low-dose projection chord diagram;

[0052] A fourth module is configured to estimate an ideal projection chord diagram without noise interference through a projection recovery model based on the low-dose projection chord diagram;

[0053] A fifth module is used to perform de-arrangement on the ideal projection chord diagram to obtain a restored integral image;

[0054] The sixth module is used to reconstruct the restored integral image to obtain a target CBCT tomographic image.

[0055] Another aspect of an embodiment of the present invention further provides an electronic device, including a processor and a memory;

[0056] The memory is used to store programs;

[0057] The processor executes the program to implement the method described above.

[0058] Another aspect of the embodiments of the present invention further provides a computer-readable storage medium, wherein the storage medium stores a program, and the program is executed by a processor to implement the method described above.

[0059] The present invention also discloses a computer program product or computer program, which includes computer instructions stored in a computer-readable storage medium. A processor of a computer device can read the computer instructions from the computer-readable storage medium and execute the computer instructions, causing the computer device to perform the above method.

[0060] The embodiment of the present invention first collects air exposure data and estimates the internal parameters of the three-dimensional X-ray system; the three-dimensional X-ray system constructed according to the three-dimensional X-ray system parameters collects the patient's original projection data, and calibrates the patient's original projection data to obtain an integral image; the integral image is data rearranged to obtain a noise-contaminated low-dose projection chord diagram; based on the low-dose projection chord diagram, an ideal projection chord diagram without noise interference is estimated through a projection recovery model; the ideal projection chord diagram is de-arranged to obtain a restored integral image; the restored integral image is reconstructed to obtain a target CBCT tomographic image. The present invention has low noise and high image resolution, and can effectively solve the technical pain points of high noise and poor image resolution of current low-dose three-dimensional images. BRIEF DESCRIPTION OF THE DRAWINGS

[0061] In order to more clearly illustrate the technical solutions in the embodiments of the present application, 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 application. For ordinary technicians in this field, other drawings can be obtained based on these drawings without any creative work.

[0062] Figure 1 An overall step flow chart provided for an embodiment of the present invention;

[0063] Figure 2 A schematic diagram of the histogram distribution of a certain projected pixel point sequence provided by an embodiment of the present invention;

[0064] Figure 3 Schematic diagram of projections before and after rearrangement in an embodiment of the present invention;

[0065] Figure 4 A schematic diagram of a method for calculating pixel mean in a projection restoration model provided by an embodiment of the present invention;

[0066] Figure 5 A comparison chart of image effects before and after processing using the low-dose reconstruction method provided by an embodiment of the present invention;

[0067] Figure 6 A comparison of images reconstructed using the low-dose reconstruction method provided by an embodiment of the present invention and other methods. DETAILED DESCRIPTION

[0068] In order to make the purpose, technical solutions and advantages of this application more clearly understood, the present application is further described in detail below with reference to the accompanying drawings and embodiments. It should be understood that the specific embodiments described herein are only used to explain this application and are not intended to limit this application.

[0069] In view of the problems existing in the prior art, one aspect of the present invention provides a three-dimensional X-ray low-dose imaging method, comprising:

[0070] Collect air exposure data and estimate the internal parameters of the 3D X-ray system;

[0071] The three-dimensional X-ray system constructed according to the three-dimensional X-ray system parameters collects original projection data of the patient, and calibrates the original projection data of the patient to obtain an integral image;

[0072] performing data rearrangement on the integral image to obtain a noise-contaminated low-dose projection chord diagram;

[0073] estimating an ideal projection chord diagram without noise interference by using a projection recovery model according to the low-dose projection chord diagram;

[0074] Rearranging the ideal projection chord diagram to obtain a restored integral image;

[0075] The restored integral image is reconstructed to obtain a target CBCT tomographic image.

[0076] Optionally, collecting air exposure data and estimating internal parameters of the three-dimensional X-ray system includes:

[0077] Configuring geometric parameters of a three-dimensional low-dose imaging device, including the distance from the radiation source to the rotation center, the distance from the rotation center to the detector, and the distance from the laser center to the imaging center;

[0078] Repeatedly collecting multiple air exposure data using the three-dimensional low-dose imaging device at a preset angle;

[0079] Determining a two-dimensional adjustment function matrix and a scale coefficient matrix of the same size as the projection according to the projection sequence in the air exposure data;

[0080] The expressions of the adjustment function and the scaling coefficient are:

[0081]

[0082] in, is the mean of the array p(u,v,:); σ 2 (u,v) is the variance of the array p(u,v,:); g(u,v) is the adjustment function; a n is the scale coefficient; N is the order of the scale coefficient.

[0083] Optionally, the three-dimensional X-ray system constructed according to the three-dimensional X-ray system parameters collects original projection data of the patient, and calibrates the original projection data of the patient to obtain an integral image, including:

[0084] performing afterglow calibration on the original projection data;

[0085] Performing scattering calibration on the afterglow-calibrated data using a Monte Carlo method;

[0086] performing hardening calibration on the scattering-calibrated data;

[0087] An integral image is determined according to the results of the persistence calibration, the scatter calibration, and the hardening calibration.

[0088] Optionally, in the step of rearranging the integral image to obtain a noise-contaminated low-dose projection chord diagram,

[0089] The result of the data rearrangement includes three-dimensional pseudo fan-beam data or three-dimensional pseudo parallel beam data;

[0090] The expression for data rearrangement is:

[0091]

[0092] Where u represents the row of the projection image in the integral image; v represents the column of the projection image in the integral image; s represents the sth projection image in the integral image; d is the detector pixel size; D is the distance from the ray source to the detector; p(s,v,u) represents the integral image data; q(u′,v′,s′) represents the rearranged three-dimensional pseudo fan-beam data.

[0093] Optionally, estimating an ideal projection chord diagram without noise interference by using a projection recovery model according to the low-dose projection chord diagram includes:

[0094] The projection recovery model is described by the maximum a posteriori probability;

[0095] Derivative the projection restoration model to obtain an optimized iterative solution formula;

[0096] Determine the projection restoration model by calculating the product of the target filter window and the image gradient under the target filter window according to the iterative solution formula;

[0097] Calculating an ideal projection chord diagram without noise interference according to the weighted mean of the filter window and the projection restoration model;

[0098] The expression of the projection restoration model is:

[0099]

[0100] in, is the lossless projection image to be restored; q i is the rearranged projected chord diagram; is the one-dimensional representation of the noise model; m is the total number of pixels; β is the relaxation factor; Ω is the neighborhood of pixel i; γ ij is the texture weight function.

[0101] Optionally, in the step of performing de-rearrangement on the ideal projection chord diagram to obtain a restored integral image, an expression for the de-rearrangement is:

[0102]

[0103] Where p(u,v,s) represents the integral image recovered after de-rearrangement; represents the ideal projection chord diagram without noise interference; u represents the row of the projection image in the integral image; v represents the column of the projection image in the integral image; s represents the s-th projection image in the integral image; d is the detector pixel size.

[0104] Optionally, reconstructing the restored integral image to obtain a target CBCT tomographic image includes:

[0105] reconstructing the restored integral image using an iterative reconstruction algorithm that has undergone sharpening modification;

[0106] determining an objective function based on the reconstructed image;

[0107] Wherein, the expression of the objective function is:

[0108]

[0109] Where f represents the reconstructed image obtained by iterative optimization; Pf represents the forward projection of the reconstructed image f; p is the de-permutated integral image; μ is the relaxation adjustment factor; R(f) is the regularization term;

[0110] The reconstructed image is accelerated by GPU parallel computing; the regularization term is constructed through gradient variance; and the iterative reconstruction algorithm is solved by the CGLS conjugate gradient descent algorithm.

[0111] Another aspect of the present invention provides a three-dimensional X-ray low-dose imaging device, comprising:

[0112] The first module is used to collect air exposure data and estimate the internal parameters of the 3D X-ray system;

[0113] A second module is configured to acquire original projection data of a patient using a three-dimensional X-ray system constructed according to the three-dimensional X-ray system parameters, and to calibrate the original projection data of the patient to obtain an integral image;

[0114] A third module is used to rearrange the data of the integral image to obtain a noise-contaminated low-dose projection chord diagram;

[0115] A fourth module is configured to estimate an ideal projection chord diagram without noise interference through a projection recovery model based on the low-dose projection chord diagram;

[0116] A fifth module is used to perform de-arrangement on the ideal projection chord diagram to obtain a restored integral image;

[0117] The sixth module is used to reconstruct the restored integral image to obtain a target CBCT tomographic image.

[0118] Another aspect of an embodiment of the present invention further provides an electronic device, including a processor and a memory;

[0119] The memory is used to store programs;

[0120] The processor executes the program to implement the method described above.

[0121] Another aspect of the embodiments of the present invention further provides a computer-readable storage medium, wherein the storage medium stores a program, and the program is executed by a processor to implement the method described above.

[0122] The present invention also discloses a computer program product or computer program, which includes computer instructions stored in a computer-readable storage medium. A processor of a computer device can read the computer instructions from the computer-readable storage medium and execute the computer instructions, causing the computer device to perform the above method.

[0123] The specific implementation principle of the present invention is described in detail below with reference to the accompanying drawings:

[0124] The specific operation steps of the present invention are as shown in the attached Figure 1 The invention is divided into 6 steps, which are explained in detail below:

[0125] 101: Collect air exposure data and estimate the system internal parameters. The geometric parameters of the three-dimensional low-dose imaging device used in the embodiment of the present invention are as follows: the distance from the ray source to the rotation center is 100 cm, the distance from the rotation center to the detector is 50 cm, and the distance from the laser center to the imaging center is 55 cm. In order to estimate the internal parameters of the system, it is necessary to repeatedly collect several air exposure data at a fixed angle. Figure 2 As shown, Figure 2 (a) is a schematic diagram of the projection sequence collected in this embodiment. Each projection pixel in this sequence can be represented by the symbol I(u,v,s), where u represents the row of the projection image, v represents the column of the projection image, and s represents the sth projection image, with s≥500. To accurately estimate the internal parameters of the system, the projection sequence for repeated collection in this embodiment is set to 1000. Since the scan is of air, it can be assumed that the rays will not experience the physical phenomena of afterglow, scattering, and hardening, but will only be affected by noise. Therefore, the integral image of air can be calculated using the following formula:

[0126]

[0127] All pixels at the (u, v) position can be extracted to form an array p(u,v,:) of 1000 elements, where: represents all the pixels. According to the knowledge of ray physics, particles have random properties when penetrating objects. The numerical distribution of the array p(u,v,:) conforms to the non-stationary Gaussian distribution, such as Figure 2 (b). Based on this physical characteristic, this embodiment uses the following formula to model the noise:

[0128]

[0129] in, is the mean of the array p(u,v,:); σ 2 (u,v) is the variance of the array p(u,v,:); g(u,v) is the adjustment function; a n is the scale coefficient; N is the order of the scale coefficient.

[0130] In this embodiment, N is selected as a 2nd-order scale. Calculate g(u,v) and a n The process can be described as follows: Taking the logarithm of both sides of the noise model formula, we can get:

[0131]

[0132] By fitting the above equation with a polynomial, we can get logg(u,v) and a n The value of g(u,v) and a n By traversing the pixels at all positions, we can obtain a two-dimensional adjustment function matrix and a scale coefficient matrix of the same size as the projection.

[0133] 102: Perform afterglow, scattering, and beam hardening calibrations on the patient's raw projection data to obtain an integral image. 3D X-ray systems typically use flat-panel detectors, which have a low photon capture rate, are prone to afterglow, and cannot use grids to filter out scattered radiation. Consequently, the quality of the resulting projection images is significantly inferior to that of CT-type linear array detectors. Furthermore, beam hardening in 3D X-ray systems is much more severe than in CT equipment. The afterglow, scattering, and beam hardening described above can be considered noise, causing offset and interference in the noise model established in 101, reducing accuracy. Therefore, it is necessary to calibrate these offsets before applying the projection recovery model. The calibration sequence is to perform afterglow calibration first, then scattering calibration, and finally hardening calibration on the data corrected by the first two calibrations. Because afterglow and hardening calibration techniques are well-established in the industry, they will not be further described in detail in this embodiment. This embodiment primarily describes scattering correction techniques based on Monte Carlo simulation. The Monte Carlo method can accurately simulate the entire physical process of X-ray photons, after being emitted from the ray source, penetrating objects and interacting (Compton, Rayleigh and photoelectric effects), and finally reaching the detector and depositing energy. It is the gold standard for scattering distribution estimation. The Monte Carlo simulation algorithm used in this embodiment is an independently developed algorithm, and the experimental verification results are consistent with the open source Monte Carlo code package (such as PENELOPE, EGSnrc, MCNP and GEANT4, etc.). The algorithm requires the input of the patient's prior three-dimensional volume image as the object to interact. The source of the prior data in this embodiment is the patient's spiral CT volume data. If the patient has no spiral CT volume data, the three-dimensional volume data reconstructed by FDK analysis is selected. Assume that S(u,v) is the scattering distribution obtained by Monte Carlo simulation, Ip (u,v) is the projection image after afterglow calibration, and the projection after scattering calibration can be described as:

[0134] I′ p (u,v)=I p (u,v)-λS(u,v)

[0135] Among them, I′ p (u, v) is the calibrated projection, and λ is the magnification factor used to match the amplitude of the Monte Carlo simulation data with the amplitude of the real data.

[0136] 103: Rearrange the integral image data to obtain a noise-contaminated low-dose projection chord diagram. In order to more efficiently recover the low-dose projection data, it is necessary to rearrange the integral image obtained in 102. The data rearrangement formula can be described as:

[0137]

[0138] Where d is the pixel size of the detector and D is the distance from the ray source to the detector. The images before and after data rearrangement are as follows: Figure 3 (a) and (b) show the projection data. The projection data is rearranged into 3D pseudo-fan-beam data, using linear interpolation as the interpolation algorithm. The projection data can also be rearranged into 3D pseudo-parallel-beam data, but this results in data redundancy loss.

[0139] 104: Apply the projection restoration model to estimate the ideal projection chord diagram without noise interference. According to Bayesian theory, the process of recovering a lossless image from a noisy image can be described by the following optimization model using the maximum a posteriori probability:

[0140]

[0141] in, is the lossless projection image to be restored. i is the rearranged projected chord diagram, represented here by a one-dimensional index. is a one-dimensional representation of the noise model defined in step 101. m is the total number of pixels, β is the relaxation factor. Ω is the neighborhood of pixel i, which is 4 in this embodiment. ij It is a texture weight function that can adaptively adjust the parameter value according to the texture of the image to remove noise while retaining the detailed information of the projection as much as possible. is called the data fidelity term of the model, It is called the penalty term of the model. By deriving the above projection recovery model, we can obtain the iterative solution formula for the above optimization:

[0142]

[0143] Among them, γ ijIt can be decomposed into the product of a filter window and the image gradient under the window. The filter window used in this embodiment is: The image gradient is calculated as: because The original calculation of requires multiple data acquisition to obtain the data mean, but it is impossible to acquire data multiple times in actual diagnosis, so the weighted mean of the filter window is used. To replace the definition in the 101 noise model The calculation is as follows Figure 4 As shown, it can be described that this mean not only considers the influence of the current filter window, but also the influence between the chord diagrams of different layers. It can be expressed as follows using mathematical formula:

[0144]

[0145] Among them, ω(s′)=[0.25 0.5 1 0.5 0.25] is the weight of the chord graph at different levels. As close to the ideal mean as possible

[0146] 105: De-rearrange the ideal projection chord diagram to obtain a restored integral image. After step 104, the noisy chord diagram acquired at low doses has been freed of most random noise interference, and the current chord diagram can be considered to be close to the ideal projection chord diagram. De-rearrange this chord diagram to obtain a restored integral image. The mathematical formula for de-rearrangement can be described as follows, which is the inverse operation of the rearrangement operation in step 103:

[0147]

[0148] The definitions of symbols d and D are consistent with those in step 103 .

[0149] 106: An improved and accelerated iterative reconstruction algorithm is used to reconstruct the restored integral image, resulting in a high-quality CBCT tomographic image with low noise and high resolution. Due to the excellent noise immunity of iterative reconstruction, the image reconstruction process does not introduce noise caused by the reconstruction algorithm, as is the case with analytical reconstruction. However, due to multiple iterations, the reconstructed image will tend to be blurred, and the edges of tissue structures will be weakened to varying degrees. To maximize the image resolution of the reconstructed image, this embodiment uses an iterative reconstruction algorithm that has undergone sharpening modification. The algorithm optimization process can be described by the following formula:

[0150]

[0151] Wherein, f represents the reconstructed image obtained by iterative optimization, Pf represents the forward projection of the reconstructed image f, and p is the integral image of the solution rearrangement. μ is the relaxation adjustment factor, and R(f) is the regularization term. Usually, R(f) is calculated by the total variation (TV) of the reconstructed image, but the TV term is still not enough to improve the weakened edge. In this embodiment, gradient variance is used to construct R(f), and its mathematical expression is described as follows:

[0152]

[0153] in, is the mean value of G. The above optimization function is solved by using the CGLS conjugate gradient descent algorithm, which can converge faster than the conventional GDLS gradient descent algorithm. Since the CGLS algorithm has been disclosed in many documents, the present invention will not repeat it and will not be included in the scope of protection. In order to further accelerate the iterative reconstruction operation, this embodiment uses GPU parallel computing to perform iterative reconstruction, and can complete the calculation of the 512x512x300 reconstruction matrix within 1 minute. After projection recovery and the improved iterative reconstruction algorithm, a high-quality CBCT tomographic image with low noise and high resolution can be obtained. As shown in the attached figure Figure 5 As shown in FIG. 1 , FIG. (a) is the reconstructed image obtained without using the solution of the present invention, and FIG. (b) is the reconstructed image obtained by using the solution of the present invention. It can be clearly seen that the improvement effect of the method of the present invention on the image effect is obvious. Figure 6 Figure (a) is an iteratively reconstructed image based on TV, and Figure (b) is an iteratively reconstructed image of the method of the present invention. It can also be seen that the reconstructed image of the method of the present invention has better uniformity and sharpness in the flat area, which shows that the method of the present invention can not only effectively solve the noise problem caused by low-dose imaging, but also maintain the texture characteristics of the tissue structure itself as much as possible.

[0154] In summary, compared with the prior art, the low-dose imaging method of the three-dimensional X-ray system of the present invention has the following advantages:

[0155] (1) There is currently no good CBCT low-dose imaging solution. The method provided by the present invention can effectively solve the technical pain points of current low-dose three-dimensional images, such as high noise and poor image resolution.

[0156] (2) The technical solution provided by the present invention can fully consider the interference of the afterglow effect, ray scattering and hardening effect of the flat-panel detector on the projection distribution model, making the projection recovery model more accurate.

[0157] (3) The iterative reconstruction algorithm provided by the present invention is a parallel iterative reconstruction algorithm that has undergone sharpening transformation and acceleration optimization. The spatial resolution of the image is higher, and the reconstruction speed is faster than other iterative algorithms.

[0158] In some optional embodiments, the function / operation mentioned in the block diagram may not occur in the order mentioned in the operation diagram. For example, depending on the function / operation involved, the two boxes shown in succession can actually be executed substantially simultaneously or the boxes can sometimes be executed in reverse order. In addition, the embodiment presented and described in the flow chart of the present invention is provided in an exemplary manner for the purpose of providing a more comprehensive understanding of the technology. The disclosed method is not limited to the operation and logic flow presented herein. Optional embodiments are contemplated in which the order of the various operations is changed and the sub-operations described as a part of a larger operation are performed independently.

[0159] Furthermore, although the present invention is described in the context of functional modules, it should be understood that, unless otherwise indicated, one or more of the functions and / or features described may be integrated into a single physical device and / or software module, or one or more functions and / or features may be implemented in separate physical devices or software modules. It will also be understood that a detailed discussion of the actual implementation of each module is not necessary for understanding the present invention. More specifically, given the properties, functions, and internal relationships of the various functional modules in the devices disclosed herein, the actual implementation of the module will be understood within the ordinary skill of an engineer. Therefore, a person skilled in the art using ordinary skill will be able to implement the present invention set forth in the claims without undue experimentation. It will also be understood that the specific concepts disclosed are merely illustrative and are not intended to limit the scope of the present invention, which is determined by the full scope of the appended claims and their equivalents.

[0160] If the functions are implemented in the form of software functional units and sold or used as independent products, they can be stored in a computer-readable storage medium. Based on this understanding, the technical solution of the present invention, or the part that contributes to the prior art, or the part of the technical solution, can be embodied in the form of a software product. The computer software product is stored in a storage medium and includes several instructions for enabling a computer device (which can be a personal computer, a server, or a network device, etc.) to execute all or part of the steps of the method described in each embodiment of the present invention. The aforementioned storage medium includes various media that can store program codes, such as a USB flash drive, a mobile hard disk, a read-only memory (ROM), a random access memory (RAM), a magnetic disk, or an optical disk.

[0161] The logic and / or steps represented in the flowcharts or otherwise described herein, for example, can be considered as an ordered list of executable instructions for implementing the logical functions, and can be embodied in any computer-readable medium for use by, or in conjunction with, an instruction execution system, apparatus, or device (e.g., a computer-based system, a system including a processor, or other system that can fetch and execute instructions from an instruction execution system, apparatus, or device). For purposes of this specification, a "computer-readable medium" can be any device that can contain, store, communicate, propagate, or transport a program for use by, or in conjunction with, an instruction execution system, apparatus, or device.

[0162] More specific examples (a non-exhaustive list) of computer-readable media include the following: an electrical connection with one or more wires (electronic devices), a portable computer disk cartridge (magnetic devices), a random access memory (RAM), a read-only memory (ROM), an erasable and programmable read-only memory (EPROM or flash memory), a fiber optic device, and a portable compact disc read-only memory (CDROM). In addition, the computer-readable medium may even be paper or other suitable medium on which the program is printed, since the program may be obtained electronically, for example, by optically scanning the paper or other medium, followed by editing, deciphering, or processing in another suitable manner as necessary, and then stored in a computer memory.

[0163] It should be understood that various parts of the present invention can be implemented using hardware, software, firmware, or a combination thereof. In the above-described embodiments, multiple steps or methods can be implemented using software or firmware stored in a memory and executed by a suitable instruction execution system. For example, if implemented using hardware, as in another embodiment, any one of the following technologies known in the art or a combination thereof can be used: a discrete logic circuit having a logic gate circuit for implementing a logic function on a data signal, an application-specific integrated circuit having a suitable combination of logic gate circuits, a programmable gate array (PGA), a field programmable gate array (FPGA), etc.

[0164] Throughout this specification, reference to terms such as "one embodiment," "some embodiments," "examples," "specific examples," or "some examples" means that a specific feature, structure, material, or characteristic described in conjunction with that embodiment or example is included in at least one embodiment or example of the present invention. In this specification, schematic representations of the above terms do not necessarily refer to the same embodiment or example. Furthermore, the specific features, structures, materials, or characteristics described may be combined in any suitable manner in any one or more embodiments or examples.

[0165] While embodiments of the present invention have been shown and described, it will be appreciated by those skilled in the art that various changes, modifications, substitutions, and variations may be made to the embodiments without departing from the principles and spirit of the invention, and that the scope of the invention is defined by the claims and their equivalents.

[0166] The above is a specific description of the preferred implementation of the present invention, but the present invention is not limited to the embodiments. Those skilled in the art can make various equivalent modifications or substitutions without violating the spirit of the present invention. These equivalent modifications or substitutions are all included in the scope defined by the claims of this application.

Claims

1. A three-dimensional X-ray low-dose imaging method, characterized in that: include: Collect air exposure data and estimate the internal parameters of the 3D X-ray system; The three-dimensional X-ray system constructed according to the three-dimensional X-ray system parameters collects original projection data of the patient, and calibrates the original projection data of the patient to obtain an integral image; performing data rearrangement on the integral image to obtain a noise-contaminated low-dose projection chord diagram; estimating an ideal projection chord diagram without noise interference by using a projection recovery model according to the low-dose projection chord diagram; Rearranging the ideal projection chord diagram to obtain a restored integral image; Reconstructing the restored integral image to obtain a target CBCT tomographic image, The three-dimensional X-ray system constructed according to the three-dimensional X-ray system parameters collects original projection data of the patient and calibrates the original projection data of the patient to obtain an integral image, including: performing afterglow calibration on the original projection data; Performing scattering calibration on the afterglow-calibrated data using a Monte Carlo method; performing hardening calibration on the scattering-calibrated data; determining an integral image according to the results of the persistence calibration, the scattering calibration and the hardening calibration, In the step of rearranging the integral image to obtain a noise-contaminated low-dose projection chord diagram, The result of the data rearrangement includes three-dimensional pseudo fan-beam data or three-dimensional pseudo parallel beam data; The expression for data rearrangement is: Where u represents the row of the projection image in the integral image; v represents the column of the projection image in the integral image; s represents the sth projection image in the integral image; d is the detector pixel size; D is the distance from the ray source to the detector; p(s,v,u) represents the integral image data; q(u′,v′,s′) represents the rearranged three-dimensional pseudo fan beam data. The method of estimating an ideal projection chord diagram without noise interference by using a projection recovery model according to the low-dose projection chord diagram comprises: The projection recovery model is described by the maximum a posteriori probability; Derivative the projection restoration model to obtain an optimized iterative solution formula; Determine the projection restoration model by calculating the product of the target filter window and the image gradient under the target filter window according to the iterative solution formula; Calculating an ideal projection chord diagram without noise interference according to the weighted mean of the filter window and the projection restoration model; The expression of the projection restoration model is: in, is the lossless projection image to be restored; qi is the rearranged projection chord diagram; is the one-dimensional representation of the noise model; m is the total number of pixels; β is the relaxation factor; Ω is the neighborhood of pixel i; γij is the texture weight function.

2. The three-dimensional X-ray low-dose imaging method according to claim 1, characterized in that: The reconstructing the restored integral image to obtain a target CBCT tomographic image includes: reconstructing the restored integral image using an iterative reconstruction algorithm that has undergone sharpening modification; determining an objective function based on the reconstructed image; Wherein, the expression of the objective function is: Where f represents the reconstructed image obtained by iterative optimization; Pf represents the forward projection of the reconstructed image f; p is the de-permutated integral image; μ is the relaxation adjustment factor; R(f) is the regularization term; The reconstructed image is accelerated by GPU parallel computing; the regularization term is constructed by gradient variance; The iterative reconstruction algorithm is solved by using the CGLS conjugate gradient descent algorithm.

Images