Dual-energy-spectrum material decomposition calibration phantom and method

Abstract

The present invention relates to the technical field of CT, and discloses a dual-energy-spectrum material decomposition calibration phantom and method. The dual-energy-spectrum material decomposition calibration phantom comprises an organic glass cylinder having a length of 16 cm, two Teflon cylinders having different diameters are embedded in the organic glass cylinder, and the diameters of the two Teflon cylinders are respectively one third and one fourth of the diameter of the organic glass cylinder. According to the dual-energy-spectrum material decomposition calibration phantom and method, calibration is performed on a projection, without relying on an image reconstruction process, so that energy spectrums and original projection data information can be more fully utilized, so as to obtain a more accurate image and a more accurate material decomposition result.

Description

A dual-energy spectroscopy material decomposition calibration phantom and method Technical Field

[0001] This invention relates to the field of CT technology, specifically to a dual-energy spectroscopy material decomposition calibration phantom and method. Background Technology

[0002] A typical CT scanner uses only one energy spectrum. The resulting image's grayscale values ​​correspond to a linear attenuation coefficient at a specific energy level, or are derived from this CT value. Such CT equipment lacks material resolution capability. Spectral CT, on the other hand, acquires scan data from multiple energy spectra simultaneously and performs corresponding multi-spectral reconstruction to obtain images with material resolution.

[0003] Spectral calibration in spectral CT is crucial for obtaining accurate quantitative scan images and physical resolution. A common practice is to calibrate and reconstruct spectral CT as multiple single-energy scans, obtaining images at two different energies, and then performing material decomposition on both sets of images.

[0004] A common dual-energy decomposition method involves first obtaining two sets of single-energy images using traditional methods, and then performing decomposition based on these images. This approach is relatively simple and intuitive, but it is limited by image quality. If there are artifacts in the image or the CT values ​​are inaccurate, the material decomposition results will be significantly compromised. In addition, errors can easily occur during the calculation process, affecting the calculation results. Therefore, this invention provides a dual-energy spectroscopy material decomposition calibration phantom and method. Summary of the Invention

[0005] To address the aforementioned shortcomings of existing technologies, this invention provides a dual-energy spectroscopy material decomposition calibration phantom and method, which can reduce errors, improve calculation accuracy, and obtain more accurate images and material decomposition results.

[0006] The present invention provides the following technical solution: a dual-energy spectroscopy material decomposition calibration phantom, comprising a 16cm long plexiglass cylinder, and two Teflon cylinders of different diameters are inlaid on the plexiglass cylinder, the diameters of the two Teflon cylinders being one-third and one-quarter of the diameter of the plexiglass cylinder, respectively.

[0007] A dual-energy spectroscopy material decomposition calibration method involves placing a dual-energy spectroscopy material decomposition calibration phantom at different positions and performing several scans. Then, the two sets of X-ray source spectra are directly calculated from the obtained projection data. The specific operation is as follows:

[0008] S1. The emission spectrum of an X-ray tube at one kVp can be expressed as the superposition of characteristic peaks and Gaussian peaks, as given by the following formula:

[0009]

[0010] Where, tri() represents a triangular peak, E0 is the characteristic peak position of the X-ray tube target material, W0 is the peak width, A0 is the peak amplitude, and the second term represents the superposition of n Gaussian peaks, {A n E n and W n}n>=1 represents the amplitude, position, and width of the Gaussian peak;

[0011] S2. Based on the data P0 obtained from the scanned model, extract the two trajectories of the model edge from the sine curve, and determine the coordinates of the model center accordingly.

[0012] S3. Using the actual projection data P0 and the theoretical projection data P obtained through scanning, we can obtain {A} using the following formula. n E n and W n The value of n>=1 is used to obtain the spectrum of X-rays, and the formula is:

[0013] {A n E n W n} = argmin(Projection(A n E n W n ;d1,d2)-P0) 2 ;

[0014] S4. Repeat this process for the data scanned at different energies to obtain two sets of spectra at different energies. Calculate the projection data P. Use the formula P = Projection(An, En and Wn; d1, d2) to obtain the projection data at two different energies, which are represented by P1 and P2 respectively.

[0015] S5. Sort the theoretical projection data at different angles according to the length of the plexiglass and Teflon along the ray path to obtain the data pairs {path length d1 of the ray through the plexiglass, path length d2 of the ray through the Teflon} => {projection value P1, projection value P2}.

[0016] S6. Using two-dimensional interpolation, the above data pairs are upsampled to obtain the projection data corresponding to the equally spaced combinations of the lengths of plexiglass and Teflon within a certain range. This set of data is denoted as R: {d1,d2}=>{P1,P2}.

[0017] S7. Perform corresponding air correction and detector geometric correction on the actual collected phantom data. The obtained data is recorded as {P1', P2'}. Based on {P1', P2'}, perform a table lookup operation in R to obtain the ideal material length {d1', d2'} corresponding to {P1', P2'}.

[0018] S8. Using a nonlinear model F to fit the mapping relationship from {P1,P2} to {d1,d2}, the data scanned in S4 can be processed using the nonlinear least squares method to obtain the ten coefficients required for model F. The specific form of this nonlinear model is as follows:

[0019]

[0020] S9. During the process, the range of {d1,d2} should be ≥0. The system scanning range is defined as an equally spaced distribution array between 0 and 200 mm. d1 is selected from 21 points with 10 mm equal intervals between 0 and 200 mm, and d2 is selected from 21 points with 3 mm equal intervals between 0 and 60 mm for sampling for subsequent calculations.

[0021] S10. Supplement the projection data between d1 and d2 (0 to 60 mm) corresponding to d1 (0 to -50 mm), and the projection data between d1 and 200 mm corresponding to d2 (0 to -15 mm). When d1 or d2 is negative, treat it as 0 for projection data to complete the data, so that the corresponding P1 and P2 results are positive.

[0022] S11. After obtaining the coefficients, save them as correction vectors. When the actual scan data arrives, model F can be used to directly obtain a set of length values ​​from a set of projection values.

[0023] S12. Based on the length values ​​at all angles, images of both plexiglass and Teflon can be reconstructed. After reconstructing the images of plexiglass (A(x,y)) and Teflon (B(x,y)) for d1 and d2 respectively, virtual monoenergetic images can be used to display images at different energies according to actual needs. The formula for any given energy E is:

[0024] μ(E, x, y) = μ PMMA (E)·A(x,y)+μ Teflon (E)·B(x,y);

[0025] According to the definition of CT value, the attenuation coefficient image is further converted into CT value through the attenuation of water under this energy for diagnostic use. When E is small, the obtained image has higher soft tissue contrast; when E is large, the metal artifacts on the image will be reduced.

[0026] Preferably, given the positions of the scanning light source, detector, and rotation center in S2, the path lengths {d1, d2} of each ray passing through the phantom and the plexiglass and Teflon, respectively, can be calculated using the obtained phantom position. This calculation is performed using a set of parameters {A}. n E n and W nWith n>=1, the theoretical projection data P on each detector pixel can be obtained, and the calculation formula is P=Projection(A n E n and W n ;d1,d2), where Projection represents the orthographic projection based on spectrum and geometry, and the projection data P is calculated as each ray passes through plexiglass and Teflon.

[0027] Compared with the prior art, the present invention has the following beneficial effects:

[0028] This dual-energy spectroscopy material decomposition calibration phantom and method performs calibration on the projection, without relying on the image reconstruction process. This allows for more full utilization of the energy spectrum and the original projection data, resulting in more accurate images and material decomposition results. Furthermore, since the d value is prone to negative values ​​during the calculation process, and d represents the physical length without any negative values, the dataset is expanded to ensure that d is positive, thereby reducing experimental errors and improving experimental accuracy. Attached Figure Description

[0029] Figure 1 is a schematic diagram of the calibration phantom structure of the present invention;

[0030] Figure 2 is a shape diagram of the emission spectrum from the X-ray tube of the present invention:

[0031] Figure 3 shows the theoretical projection data of the two sets of spectra in this invention;

[0032] Figure 4 is a data graph of d1 and d2 of the present invention;

[0033] Figure 5 shows the truncation effect of the dataset of this invention;

[0034] Figure 6 is a data projection diagram after the expansion of the present invention;

[0035] Figure 7 is a graph showing the calculation of positive values ​​of d1 and d2 in this invention. Detailed Implementation

[0036] To make the objectives, technical solutions, and advantages of the embodiments of this disclosure clearer, the technical solutions of the embodiments of this disclosure will be clearly and completely described below with reference to the accompanying drawings. In order to keep the following description of the embodiments of this disclosure clear and concise, detailed descriptions of known functions and known components are omitted to avoid unnecessarily obscuring the concept of the present invention.

[0037] Please refer to Figure 1. A dual-energy spectroscopy material decomposition calibration phantom includes a 16cm long plexiglass cylinder, on which two Teflon cylinders of different diameters are embedded. The diameters of the two Teflon cylinders are one-third and one-quarter of the diameter of the plexiglass cylinder, respectively.

[0038] A dual-energy spectroscopy material decomposition calibration method involves placing a dual-energy spectroscopy material decomposition calibration phantom at different positions and performing several scans. Then, the two sets of X-ray source spectra are directly calculated from the obtained projection data. The specific operation is as follows:

[0039] S1. The emission spectrum of an X-ray tube at one kVp can be represented as the superposition of characteristic peaks and Gaussian peaks, as shown in Figure 2. The formula is:

[0040]

[0041] Where tri() represents a triangular peak, E0 is the characteristic peak position of the X-ray tube target material, W0 is the peak width, A0 is the peak amplitude, and the second term represents the superposition of n Gaussian peaks, {An, En and Wn}n>=1 is the amplitude, peak position and peak width of the Gaussian peak. For common tungsten targets, E0 is taken as a fixed value of 69.525keV. According to the experimental test results, n can be an integer between 3 and 6.

[0042] S2. Based on the data P0 obtained from the scanned model, extract the two trajectories of the model edge from the sine curve, and determine the coordinates of the model center accordingly.

[0043] Given the positions of the scanning light source, detector, and rotation center, the path lengths {d1, d2} of each ray passing through the phantom and the plexiglass and Teflon, respectively, can be calculated using the obtained phantom position. This calculation is performed using a set of parameters {A}. n E n and W n With n>=1, the theoretical projection data P on each detector pixel can be obtained, and the calculation formula is P=Projection(A n E n and W n ;d1,d2), where Projection represents the orthographic projection based on spectrum and geometry, and the projection data P obtained by formula when each ray passes through plexiglass and Teflon.

[0044] S3. Using the actual projection data P0 and the theoretical projection data P obtained through scanning, we can obtain {A} using the following formula. n E n and W n The value of n>=1 is used to obtain the spectrum of X-rays, and the formula is:

[0045] {A n E n W n} = argmin(Projection(A n E n W n;d1,d2)-P0) 2 .

[0046] S4. Repeat this process for the scanned data at different energies to obtain two sets of spectra at different energies. Calculate the projection data P and use the formula P = Projection(An, En and Wn; d1, d2) to obtain the projection data at the two different energies, which are represented by P1 and P2 respectively.

[0047] S5. Sort the theoretical projection data at different angles according to the length of the plexiglass and Teflon along the ray path to obtain the data pairs {path length d1 of the ray through the plexiglass, path length d2 of the ray through the Teflon} => {projection value P1, projection value P2}.

[0048] S6. Using two-dimensional interpolation, the above data pairs are upsampled to obtain the projection data corresponding to the equally spaced combinations of the lengths of plexiglass and Teflon within a certain range. This set of data is denoted as R: {d1,d2}=>{P1,P2}.

[0049] S7. Perform corresponding air correction and detector geometric correction on the actual collected phantom data. The obtained data is recorded as {P1', P2'}. Based on {P1', P2'}, perform a table lookup operation in R to obtain the ideal material length {d1', d2'} corresponding to {P1', P2'}.

[0050] S8. Using a nonlinear model F to fit the mapping relationship from {P1,P2} to {d1,d2}, the data scanned in S4 can be processed using the nonlinear least squares method to obtain the ten coefficients required for model F. The specific form of this nonlinear model is as follows:

[0051]

[0052] S9. During the process, the range of {d1,d2} should be ≥0. The system scanning range is defined as an equally spaced distribution array between 0 and 200 mm. {P1,P2} are the theoretical projection data of different lengths of {d1,d2} under two sets of spectra, with a numerical range of 0 to 1. The correspondence between {d1,d2} and P under one set of spectra is shown in Figure 3. d1 is sampled from 21 points with 10 mm equal intervals between 0 and 200 mm, and d2 is sampled from 21 points with 3 mm equal intervals between 0 and 60 mm for subsequent calculations. The calculation data is shown in Figure 4.

[0053] If such a dataset is used directly for curve fitting, the fitting curve will fluctuate greatly around 0 due to the cutoff effect at the positions where d1 and d2 are 0, and may even result in negative values ​​for d1 and d2, as shown in Figure 5.

[0054] S10. Supplement the projection data between d1 (0 to -50 mm) and d2 (0 to 60 mm), and the projection data between d1 (0 to 200 mm) corresponding to d2 (0 to -15 mm). When d1 or d2 is negative, treat it as 0 for projection data to complete the data, as shown in Figure 6, so that the corresponding P1 and P2 results are positive, as shown in Figure 7. At the same time, even when there is system noise in p1 and p2, the fitting result can still be guaranteed to be positive.

[0055] S11. After obtaining the coefficients, save them as correction vectors. When the actual scan data arrives, model F can be used to directly obtain a set of length values ​​from a set of projection values.

[0056] S10. Based on the length values ​​at all angles, images of both plexiglass and Teflon can be reconstructed. After reconstructing the images of plexiglass (A(x,y)) and Teflon (B(x,y)) for d1 and d2 respectively, virtual monoenergetic images can be used to display images at different energies according to actual needs. The formula for any given energy E is:

[0057] μ(E, x, y) = μ PMMA (E)·A(x,y)+μ Teflon (E)·B(x,y);

[0058] According to the definition of CT value, the attenuation coefficient image is further converted into CT value through the attenuation of water under this energy for diagnostic use. When E is small, the obtained image has higher soft tissue contrast; when E is large, the metal artifacts on the image will be reduced.

[0059] The above embodiments are merely exemplary embodiments of the present invention and are not intended to limit the present invention. The scope of protection of the present invention is defined by the claims. Those skilled in the art can make various modifications or equivalent substitutions to the present invention within its spirit and scope of protection, and such modifications or equivalent substitutions should also be considered to fall within the scope of protection of the present invention.

Claims

1. A dual-energy spectroscopy material decomposition and calibration method, characterized in that: The calibration phantom is a 16cm long acrylic cylinder, and two Teflon cylinders of different diameters are inlaid on the acrylic cylinder. The diameters of the two Teflon cylinders are one-third and one-quarter of the diameter of the acrylic cylinder, respectively. The calibration phantom was placed in different positions and scanned several times. Then, the two sets of X-ray source spectra were directly calculated from the obtained projection data. The specific operation is as follows: S1. The emission spectrum of an X-ray tube at one kVp can be expressed as the superposition of characteristic peaks and Gaussian peaks, as given by the following formula: Where tri() represents a triangular peak, E0 is the characteristic peak position of the X-ray tube target material, W0 is the peak width, A0 is the peak amplitude, and the second term represents the superposition of n Gaussian peaks, n>=1, A n E n and W n This represents the amplitude, position, and width of the Gaussian peak. S2. Based on the data P0 obtained from scanning the model, extract the two trajectories of the model edge from the sine curve, and calculate the coordinates of the model center accordingly. S3. Using the actual projection data P0 and the theoretical projection data P obtained through scanning, we can obtain {A} using the following formula. n E n W n The value of} is used to obtain the spectrum of X-rays, and the formula is: {AT n ,IT n ,W n }=argmin(Projection(A n ,IT n ,W n ;d1,d2)-P0) 2 ; S4. Repeat this process for the scanned data at different energies to obtain two sets of spectra at different energies. Calculate the projection data P. Calculate the projection data at the two different energies using the formula P = Projection(An, En, Wn; d1, d2), and represent them as P1 and P2 respectively. S5. Sort the theoretical projection data at different angles according to the length of the plexiglass and Teflon along the ray path to obtain the data pairs {path length d1 of the ray through the plexiglass, path length d2 of the ray through the Teflon} => {projection value P1, projection value P2}. S6. Using two-dimensional interpolation, the above data pairs are upsampled to obtain the projection data corresponding to the equally spaced combinations of the lengths of plexiglass and Teflon within a certain range. This set of data is denoted as R: {d1,d2}=>{P1,P2}. S7. Perform corresponding air correction and detector geometric correction on the actual collected phantom data. The obtained data is recorded as {P1', P2'}. Based on {P1', P2'}, perform a table lookup operation in R to obtain the ideal material length {d1', d2'} corresponding to {P1', P2'}. S8. Use a nonlinear model F to fit the mapping relationship from {P1,P2} to {d1,d2}. Using the data scanned in S4, apply the nonlinear least squares method to obtain the ten coefficients required for model F. The specific form of this nonlinear model is as follows: S9. During the process, the range of {d1,d2} is greater than or equal to 0. The system scanning range is defined as an equally spaced distribution array between 0 and 200 mm. d1 is selected from 21 points with equal intervals of 10 mm between 0 and 200 mm, and d2 is selected from 21 points with equal intervals of 3 mm between 0 and 60 mm for sampling for subsequent calculations. S10. Supplement the projection data between d1 and d2, which are between 0 and 60 mm, corresponding to d1 being between 0 and -50 mm, and the projection data between d1 and 200 mm, corresponding to d2 being between 0 and -15 mm. When d1 or d2 is negative, treat it as 0 and use it as projection data to fill in the data, so that the corresponding P1 and P2 results are positive. S11. After obtaining the coefficients, save them as correction vectors. When the actual scan data arrives, use model F to directly obtain a set of length values ​​from a set of projection values. S12. Based on the length values ​​at all angles, reconstruct images of both plexiglass and Teflon. After reconstructing images of plexiglass (A(x,y)) and Teflon (B(x,y)) for d1 and d2 respectively, use a virtual monoenergetic image to display images at different energies as needed. The formula for any given energy E is: μ(E,x,y)=μ PMMA (E)·A(x,y)+μ Teflon (E)·B(x,y); According to the definition of CT value, the attenuation coefficient image is further converted into CT value through the attenuation of water under this energy for diagnostic use. When E is small, the obtained image has higher soft tissue contrast; when E is large, the metal artifacts on the image will be reduced.

2. The dual-energy spectroscopy material decomposition and calibration method according to claim 1, characterized in that: Given the positions of the scanning light source, detector, and rotation center in S2, the path lengths {d1, d2} of each ray passing through the phantom and plexiglass and Teflon, respectively, are calculated using the obtained phantom position. This calculation is performed using a set of parameters {A}. n E n W n Under the following conditions, the theoretical projection data P on each detector pixel is calculated, and the calculation formula is P = Projection(A) / P(A) ... n E n W n ;d1,d2), where Projection represents the orthographic projection based on spectrum and geometry, and the projection data P obtained when each ray passes through plexiglass and Teflon.

Images