CN113391341B - X-ray energy spectrum estimation method considering influence of scattered photons - Google Patents
X-ray energy spectrum estimation method considering influence of scattered photons Download PDFInfo
- Publication number
- CN113391341B CN113391341B CN202110568625.8A CN202110568625A CN113391341B CN 113391341 B CN113391341 B CN 113391341B CN 202110568625 A CN202110568625 A CN 202110568625A CN 113391341 B CN113391341 B CN 113391341B
- Authority
- CN
- China
- Prior art keywords
- energy spectrum
- ray energy
- ray
- value
- estimated
- Prior art date
- Legal status (The legal status is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the status listed.)
- Active
Links
- 238000001228 spectrum Methods 0.000 title claims abstract description 122
- 238000000034 method Methods 0.000 title claims abstract description 42
- 238000004422 calculation algorithm Methods 0.000 claims abstract description 14
- 238000002083 X-ray spectrum Methods 0.000 claims description 8
- 238000005070 sampling Methods 0.000 claims description 4
- XAGFODPZIPBFFR-UHFFFAOYSA-N aluminium Chemical compound [Al] XAGFODPZIPBFFR-UHFFFAOYSA-N 0.000 description 15
- 229910052782 aluminium Inorganic materials 0.000 description 15
- 238000010586 diagram Methods 0.000 description 12
- 238000012545 processing Methods 0.000 description 6
- 238000012360 testing method Methods 0.000 description 4
- 238000002591 computed tomography Methods 0.000 description 3
- 238000012937 correction Methods 0.000 description 3
- 238000005516 engineering process Methods 0.000 description 3
- 238000010606 normalization Methods 0.000 description 3
- NAWXUBYGYWOOIX-SFHVURJKSA-N (2s)-2-[[4-[2-(2,4-diaminoquinazolin-6-yl)ethyl]benzoyl]amino]-4-methylidenepentanedioic acid Chemical compound C1=CC2=NC(N)=NC(N)=C2C=C1CCC1=CC=C(C(=O)N[C@@H](CC(=C)C(O)=O)C(O)=O)C=C1 NAWXUBYGYWOOIX-SFHVURJKSA-N 0.000 description 2
- 238000013170 computed tomography imaging Methods 0.000 description 2
- 238000002474 experimental method Methods 0.000 description 2
- 239000000126 substance Substances 0.000 description 2
- 238000004458 analytical method Methods 0.000 description 1
- 238000004364 calculation method Methods 0.000 description 1
- 230000007547 defect Effects 0.000 description 1
- 238000001514 detection method Methods 0.000 description 1
- 238000003745 diagnosis Methods 0.000 description 1
- 230000000694 effects Effects 0.000 description 1
- 238000003384 imaging method Methods 0.000 description 1
- 230000036039 immunity Effects 0.000 description 1
- 238000004519 manufacturing process Methods 0.000 description 1
- 239000000463 material Substances 0.000 description 1
- 238000005259 measurement Methods 0.000 description 1
- 238000012986 modification Methods 0.000 description 1
- 230000004048 modification Effects 0.000 description 1
- 230000011218 segmentation Effects 0.000 description 1
- 238000004088 simulation Methods 0.000 description 1
- 238000006467 substitution reaction Methods 0.000 description 1
- 230000000007 visual effect Effects 0.000 description 1
Classifications
-
- G—PHYSICS
- G01—MEASURING; TESTING
- G01T—MEASUREMENT OF NUCLEAR OR X-RADIATION
- G01T1/00—Measuring X-radiation, gamma radiation, corpuscular radiation, or cosmic radiation
- G01T1/36—Measuring spectral distribution of X-rays or of nuclear radiation spectrometry
-
- G—PHYSICS
- G01—MEASURING; TESTING
- G01N—INVESTIGATING OR ANALYSING MATERIALS BY DETERMINING THEIR CHEMICAL OR PHYSICAL PROPERTIES
- G01N23/00—Investigating or analysing materials by the use of wave or particle radiation, e.g. X-rays or neutrons, not covered by groups G01N3/00 – G01N17/00, G01N21/00 or G01N22/00
- G01N23/02—Investigating or analysing materials by the use of wave or particle radiation, e.g. X-rays or neutrons, not covered by groups G01N3/00 – G01N17/00, G01N21/00 or G01N22/00 by transmitting the radiation through the material
- G01N23/04—Investigating or analysing materials by the use of wave or particle radiation, e.g. X-rays or neutrons, not covered by groups G01N3/00 – G01N17/00, G01N21/00 or G01N22/00 by transmitting the radiation through the material and forming images of the material
- G01N23/046—Investigating or analysing materials by the use of wave or particle radiation, e.g. X-rays or neutrons, not covered by groups G01N3/00 – G01N17/00, G01N21/00 or G01N22/00 by transmitting the radiation through the material and forming images of the material using tomography, e.g. computed tomography [CT]
Landscapes
- Physics & Mathematics (AREA)
- Health & Medical Sciences (AREA)
- Spectroscopy & Molecular Physics (AREA)
- Life Sciences & Earth Sciences (AREA)
- General Physics & Mathematics (AREA)
- Nuclear Medicine, Radiotherapy & Molecular Imaging (AREA)
- Chemical & Material Sciences (AREA)
- Engineering & Computer Science (AREA)
- High Energy & Nuclear Physics (AREA)
- Pulmonology (AREA)
- Radiology & Medical Imaging (AREA)
- Theoretical Computer Science (AREA)
- Molecular Biology (AREA)
- Analytical Chemistry (AREA)
- Biochemistry (AREA)
- General Health & Medical Sciences (AREA)
- Immunology (AREA)
- Pathology (AREA)
- Analysing Materials By The Use Of Radiation (AREA)
Abstract
The embodiment of the invention discloses an X-ray energy spectrum estimation method considering scattered photon influence, which comprises the following steps: constructing a nonlinear equation set about the X-ray energy spectrum and the scattering constant; setting an X-ray energy spectrum S (E) in the nonlinear equation set and an initial value S (E) of the scattering constant Sc (0) And Sc (Sc) (0) The method comprises the steps of carrying out a first treatment on the surface of the At each iteration point (S (E) (n) ，Sc (n) ) Performing first-order taylor expansion on the nonlinear equation set to obtain a linear equation set; converting a linear equation set into a form ax=b, and solving the linear equation set ax=b by using an EM algorithm to obtain estimated values of the X-ray energy spectrum and the scattering constant; judging whether the estimated value converges or not, if not, executing the next iteration; and when the estimated values of the X-ray energy spectrum and the scattering constant are converged, obtaining a final X-ray energy spectrum according to the estimated values of the X-ray energy spectrum and the scattering constant.
Description
Technical Field
The invention relates to the technical field of X-ray CT imaging, in particular to an X-ray energy spectrum estimation method considering the influence of scattered photons.
Background
The X-ray computed tomography (X ray Computed Tomography, abbreviated as X-ray CT) imaging technology can display the internal structural information of the object to be measured without damage, and is increasingly widely used in medical diagnosis and treatment, industrial defect detection, manufacturing error analysis, reverse engineering and other fields. The accuracy of CT imaging depends on the accuracy of X-ray polychromatic orthographic modeling, and the X-ray energy spectrum is an important factor affecting the accuracy of X-ray polychromatic orthographic modeling. The X-ray energy spectrum plays an important role in applications such as energy spectrum CT image reconstruction, CT image hardening correction, and the like. The X-ray flow generated by an X-ray source in the CT system is very strong, and the X-ray energy spectrum is difficult to be measured by directly utilizing equipment such as an X-ray spectrometer and the like; the scattering size is related to the attribute of the measured object and the environment of the CT system, so that accurate estimation is difficult, and scattered photons have a certain influence on the accuracy of energy spectrum estimation. Therefore, to accurately model the X-ray polychromatic orthographic projection process, accurate estimation of X-ray energy spectrum and scattering is required.
Disclosure of Invention
It is an object of embodiments of the present invention to provide a method of estimating an X-ray energy spectrum taking into account the influence of scattered photons, which overcomes or at least alleviates at least one of the above-mentioned drawbacks of the prior art.
To achieve the above object, an embodiment of the present invention provides an X-ray energy spectrum estimation method considering influence of scattered photons, including:
step 1, setting a scattered photon signal as a constant, and constructing a nonlinear equation set about an X-ray energy spectrum and a scattering constant:
wherein ,pi Is the X-ray polychromatic projection value obtained under the ith ray path; i represents the number of the ray path, the range of values starts from 1, and the maximum value is the product of the number of units of the detector and the number of scanning angles of the object;an index set for an X-ray path; s (E) = (S) 1 ,S 2 ,…，S M ) Is a discrete form of an unknown X-ray energy spectrum, where M represents the maximum energy of the X-ray energy spectrum; delta is the energy spectrum discrete interval; phi (E) = (phi) 1 ,φ 2 ,…,φ M ) Is a discrete form of the mass attenuation coefficient of the measured object; r is R i ＝(r i1 ,r i2 ,…,r iJ ) Is the projection vector of the ith ray path, r ij Representing the projection contribution of the jth pixel to the ith ray path; f= (f 1 ,f 2 ,…,f J ) T Representing a discretized image, f j A sampling value at the j-th pixel for the image f;
step 2, setting an X-ray energy spectrum S (E) in the nonlinear equation set and an initial value S (E) of the scattering constant Sc (0) and Sc(0) The method comprises the steps of carrying out a first treatment on the surface of the Wherein, the S (E) (0) The X-ray energy spectrum under the preset voltage is simulated by the X-ray energy spectrum software to obtain,the Sc is (0) A preset value within a range of 0.001-0.1;
step 3, at each iteration point (S (E) (n) ,Sc (n) ) Performing first-order taylor expansion on the nonlinear equation set (1) to obtain a linear equation set as follows:
wherein ：
wherein S (E) (n) and Sc(n) Values of S (E) and Sc after the end of the nth iteration are respectively represented;
step 4, converting the equation of the linear equation set (2) into a form ax=b, and solving the linear equation set ax=b by using an EM (Expectation-Maximization) algorithm to obtain estimated values of the X-ray energy spectrum and the scattering constant;
wherein ,
x＝(x 1 ,x 2 ,…,x N ) T ＝(S(E),Sc) T (7)
solving the linear system of equations ax=b using the EM algorithm using the formula:
wherein x represents the unknown number to be solved, x j Represents the j-th unknown number, x j (k) X representing the kth iteration j M represents the number of equations in the linear system of equations, N represents the number of unknowns in the linear system of equations, a ij Representing the coefficient corresponding to the jth unknown of the ith equation, b i A numerical value representing the right-hand term of the i-th equation;
step 5, judging whether the estimated values of the X-ray energy spectrum and the scattering constant obtained in the step 4 are converged, if not, executing the next iteration, and repeating the steps 3 and 4; when the estimated values of the X-ray energy spectrum and the scattering constant are converged, executing a step 6;
and 6, obtaining a final X-ray energy spectrum according to the X-ray energy spectrum and the estimated value of the scattering constant.
Optionally, the initial value of the scattering constant Sc is 0.01.
Optionally, S (E) (0) ＞0，Sc (0) > 0, and:
optionally, determining whether the estimated values of the X-ray energy spectrum and the scattering constant obtained in the step 4 converge includes:
calculating an estimated polychromatic projection value p according to the formula (3) using the X-ray energy spectrum and the estimated value of the scattering constant i (n) The method comprises the steps of carrying out a first treatment on the surface of the If p i (n) And the measured polychromatic projection value p i The square of the two norms of the difference is smaller than a given threshold epsilon 1 Judging convergence, namely:
wherein ε1 ＞0。
Optionally, determining whether the estimated values of the X-ray energy spectrum and the scattering constant obtained in the step 4 converge includes:
the square of the two norms of the difference between the X-ray spectra estimated in two adjacent iterations is less than a given threshold ε 2 Judging convergence, namely:
wherein ε2 ＞0。
Compared with the prior art, the invention has at least the following advantages:
because the method considers the low-frequency scattering information when estimating the X-ray energy spectrum, the accuracy of the X-ray energy spectrum in a high-energy section is improved; when estimating the X-ray energy spectrum and the scattering value, linearizing the nonlinear equation set, and then carrying out gradual iterative updating by using an EM algorithm, so that the dependence of the EM algorithm on an initial value is reduced; the method has low requirement on the precision of the calibration die body, can meet the requirement by adopting a cylinder made of single material, and is simple and easy to use.
Drawings
FIG. 1 shows a flow chart of an X-ray energy spectrum estimation method taking into account scattered photon effects provided by an embodiment of the invention;
FIG. 2a is a schematic diagram of an energy spectrum to be estimated and an initial value set according to an embodiment of the present invention;
FIG. 2b is a schematic diagram of a circular aluminum block mold for testing according to an embodiment of the present invention;
FIG. 3 is a schematic diagram of an X-ray energy spectrum obtained by an X-ray energy spectrum estimation method considering the influence of scattered photons according to an embodiment of the present invention;
FIG. 4a is a photograph of a cylindrical aluminum block phantom of 4cm diameter for estimating an X-ray energy spectrum according to another embodiment of the present invention;
FIG. 4b is a reconstructed image obtained by directly reconstructing the aluminum block of FIG. 4a using the SART algorithm after scanning in accordance with another embodiment of the present invention;
FIG. 4c is a view of the reconstructed image of FIG. 4b after thresholding and binarizing operations according to another embodiment of the present invention;
FIG. 4d is a scatter diagram showing a correspondence between the intersection line length of the aluminum block mold body and the polychromatic projection obtained after performing the CT orthographic projection operation on FIG. 4c according to another embodiment of the present invention;
FIG. 4e is a graph showing a relationship between the intersection line length of the aluminum block mold body and the polychromatic projection obtained by performing a weighted average operation on the polychromatic projection values corresponding to the same intersection line length in FIG. 4d according to another embodiment of the present invention;
FIG. 5a is a schematic diagram showing an initial value of an X-ray energy spectrum to be estimated according to another embodiment of the present invention and an estimated result of the X-ray energy spectrum by the method of the present invention;
FIG. 5b is a schematic diagram showing how the polychromatic projection curve calculated by using the estimated X-ray energy spectrum in FIG. 5a fits to the polychromatic projection curve obtained in practice according to another embodiment of the present invention;
FIG. 6a is a photograph of a cylindrical aluminum block mold body of 8cm diameter provided in an example of the present invention for a hardening artifact removal experiment;
FIG. 6b is a reconstructed image obtained by directly reconstructing the aluminum block of FIG. 6a using the SART algorithm after scanning, which is also an image before correction of hardening artifacts, according to an example of the present invention;
FIG. 6c is a graph showing the result of the present invention after the X-ray energy spectrum estimated in FIG. 5a has been used to cure;
FIG. 6d is a line density curve image obtained by centering the image of FIG. 6c in accordance with an example of the present invention.
Detailed Description
The present invention will be described in detail with reference to the accompanying drawings and examples.
Fig. 1 shows a flowchart of an X-ray energy spectrum estimation method taking into account the influence of scattered photons provided by an embodiment of the present invention. As shown in fig. 1, the method for estimating an X-ray energy spectrum taking into account influence of scattered photons provided by the embodiment of the invention includes the following steps:
step 110, setting scattered photon signals as constants, and constructing a nonlinear equation set about the X-ray energy spectrum and the scattering constants by combining a nonlinear model of the X-ray polychromatic projection forward process.
The method comprises the following steps of setting a low-frequency signal as a scattered photon signal, setting a constant Sc, and constructing a nonlinear equation set about an X-ray energy spectrum and a scattering constant by combining a nonlinear model of an X-ray polychromatic projection positive process:
wherein ,pi Is the polychromatic projection value of X-ray obtained under the ith ray path, and the value is obtained through measurement. i represents the number of the ray path, the range of values starts from 1, and the maximum value is the number of units of the detector multiplied by the number of object scanning angles.For the index set of the X-ray path, the value of which is determined by the CT system, in the case of a scan parameter determination,/-level>Is fixed. S (E) = (S) 1 ,S 2 ,…,S M ) Is a discrete form of an unknown X-ray energy spectrum, where M represents the maximum energy of the X-ray energy spectrum, its value is determined by the voltage of the CT system, and assuming that the current CT system voltage is 140kV, then the value of M is 140. Delta is the energy spectrum discrete interval, and the value is fixed, and is 1 in the embodiment. Phi (E) = (phi) 1 ,φ 2 ,…,φ M ) Is a discrete form of the mass attenuation coefficient of the object under test, and the form of phi (E) is determined when the object under test is determined. R is R i ＝(r i1 ,r i2 ,…,r iJ ) Is the projection vector of the ith ray path, r ij Representing the projection contribution of the jth pixel to the ith ray path, R i and rij Is determined by the CT system, R in the case of a scan parameter determination i and rij The values are fixed. f= (f 1 ,f 2 ,…,f J ) T Representing a discretized image, f j For the sample value of the image f at the j-th pixel, f= (f in the case of scan parameter determination 1 ,f 2 ,…,f J ) T The value of (2) is fixed.
By the above formula (1), a nonlinear relationship between the X-ray energy spectrum and the scattering constant is established.
Step 120, initializing the X-ray energy spectrum and scatter in the nonlinear system of equations.
In one embodiment, let the initial value of the X-ray spectrum be S (E) (0) The initial value of scattering is Sc (0) The method comprises the steps of carrying out a first treatment on the surface of the Wherein S (E) (0) ＞0，Sc (0) > 0, and
the simulation Spectrum at the corresponding voltage can be generated as S (E) by using X-ray Spectrum open source software Spectrum GUI (0) Is set to be a constant value. S (E) (0) The initial value of (2) may be selected from empirical values, with approximate values ranging from 0.001 to 0.1.
And 130, performing first-order Taylor expansion on the nonlinear equation set to obtain a corresponding linear equation set.
In one embodiment, at each iteration point (S (E) (n) ,Sc (n) ) The nonlinear equation set (1) is subjected to first-order Taylor expansion, wherein S (E) (n) and Sc(n) The values of S (E) and Sc after the end of the nth iteration are shown, respectively. The linear equation set obtained after taylor expansion is:
wherein ：
and 140, solving the linear equation set by using an EM (estimation-Maximization) algorithm, and updating the estimated values of the X-ray energy spectrum and the scattering.
In one embodiment, the equation of the system of linear equations (2) is converted to the form ax=b,
x＝(x 1 ,x 2 ,…,x N ) T ＝(S(E),Sc) T (7)
the linear equation set ax=b is solved by using the EM algorithm, and the formula is:
wherein x represents the unknown number to be solved, x j Represents the j-th unknown number, x j (k) X representing the kth iteration j M represents the number of equations in the linear system of equations, N represents the number of unknowns in the linear system of equations, a ij Representing the coefficient corresponding to the jth unknown of the ith equation, b i The value of the right-hand term representing the i-th equation.
In one embodiment, in this step, the system of linear equations ax=b is solved according to (9) to obtain S (E) (n) and Sc(n) And updating the estimated value of the X-ray energy spectrum and the scattering by using the value obtained by solving, and then carrying out normalization processing on the estimated X-ray energy spectrum and the scattering value, namely:
to constrain the magnitude of the scatter value Sc.
Step 150, judging whether the estimated X-ray energy spectrum and the scattering value are converged, if not, repeating steps 130 to 140, and if so, executing step 160.
In one embodiment, it is determined whether the estimated X-ray energy spectrum and the scatter value converge on two conditions, either of which is satisfied, and the iteration may be stopped.
Condition 1, calculating an estimated polychromatic projection value p according to formula (3) using the estimated X-ray energy spectrum and scatter values i (n) The method comprises the steps of carrying out a first treatment on the surface of the If p i (n) And the measured polychromatic projection value p i The square of the two norms of the difference is smaller than a given threshold epsilon 1 Then it is considered to have converged, i.e.:
wherein ε1 The value of > 0 in this example is 0.001.
Condition 2, the square of the two norms of the difference between the X-ray energy spectra estimated in two adjacent iterations is smaller than a given threshold epsilon 2 Then it is considered to have converged, i.e.:
wherein ε2 The value of > 0 in this example is 0.001.
It will be readily appreciated that other convergence conditions may be selected, including, but not limited to, a maximum number of iterations, and a degree of visual fit of the estimated polychromatic projection curve generated using the estimated X-ray energy spectrum to the measured polychromatic projection curve, etc., to determine whether the convergence condition is met.
Step 160, obtaining a final X-ray energy spectrum according to the converged X-ray energy spectrum and the scattering value.
The following describes a specific implementation procedure of the X-ray energy spectrum estimation method provided by the invention, which considers the influence of scattered photons, through a specific embodiment.
Fig. 2a is a schematic diagram of the energy spectrum to be estimated and the set initial value provided in the present embodiment. The energy Spectrum to be estimated is generated by utilizing an X-ray energy Spectrum open source software Spectrum GUI, the energy Spectrum of the HiRay 7 ray source under 140kV voltage is simulated, and normalization processing is carried out on the energy Spectrum, as shown by a dotted line in fig. 2 a; the initial values of the spectra were also generated using a spectra GUI, simulating the spectra of an Oxford Series6000 source at 140kV voltage, and normalized as shown by the solid line in fig. 2 a. Fig. 2b is a schematic diagram of a circular aluminum block phantom for testing provided in this embodiment, the mass attenuation coefficient of aluminum being available from, for example, the National Institute of Standards and Technology (NIST) website. The sampling interval of the mass attenuation coefficient of the X-ray energy spectrum and the mass attenuation coefficient of the substance is 1kV.
CT scanning is performed on the circular aluminum block of fig. 2b by using the energy spectrum to be estimated, and projection data p of 720 angles under 1536 detector units is calculated according to formula (1), i.e. i=720×1536, where the scattering constant in formula (1) is set to sc=0.01. The scan parameters were set as follows: the distance (SOD) from the source to the center of the turntable was 500mm and the distance (SDD) from the source to the detector was 1000mm. The number of line detector units is 1024, and the size of each detector unit is 0.2mm. Adding an initial incident intensity I to the obtained projection data p 0 ＝10 5 The noise addition formula is as follows:
wherein ,pnoisy To add noisy projection data, poissrnd is a function of the poisson distribution random number.
The specific implementation process of the X-ray energy spectrum estimation method taking the influence of scattered photons into consideration provided in the example comprises the following steps:
1) Applying the initial value to the step 110 in the formula (1)Is given an initial value of S (E) (0) The initial value of the scattering is set to Sc (0) ＝0.001。
2) Assume that n iterations have been performed, at which time the estimated values of the X-ray spectrum and scatter are (S (E) (n) ,Sc (n) ). Regarding the formula (1) as (S (E), sc) is (S (E) (n) ,Sc (n) ) Performing first-order taylor expansion to obtain a linear equation set as the formula (2) in the step 130; wherein p is i (n) 、q i (n) And Φ, refer to the above formulas (3) - (5).
3) And (3) according to formulas (6) - (8), the linear equation set (2) is correspondingly in the form of ax=b, the iteration format, namely formula (9), is utilized to solve and update the X-ray energy spectrum and the scattering value, and normalization processing is carried out on the updated X-ray energy spectrum and the updated scattering value according to formula (10).
4) If the X-ray energy spectrum and the scattering value obtained in the step 3) are converged, stopping iteration, and obtaining a final X-ray energy spectrum according to the converged X-ray energy spectrum and the scattering value; otherwise, the next iteration (i.e. n+1 iterations) is performed, and the operations in steps 2) and 3) are repeated until the X-ray energy spectrum and the scatter value converge.
Fig. 3 is a graph of the result of the X-ray spectrum obtained by the method for estimating the X-ray spectrum considering the influence of scattered photons in this embodiment. As shown in fig. 3, the algorithm herein can effectively recover the X-ray energy spectrum to be estimated from the projection data containing noise and scattering, and has a certain noise immunity.
The use of the invention in an actual scanning process and its use for removing stiffening artifacts is described below by way of a specific embodiment.
Fig. 4a is a photograph of a cylindrical aluminum block phantom of 4cm diameter for estimating X-ray energy spectrum provided in this example, the mass attenuation coefficient of aluminum obtained from National Institute of Standards and Technology (NIST) website. The Spectrum of the Oxford Series6000 ray source under 140kV voltage is simulated by using the Spectrum open source software Spectrum GUI generation, and normalized as an initial value S (E) (0) . The sampling interval of the mass attenuation coefficient of the X-ray energy spectrum and the mass attenuation coefficient of the substance is 1kV. The industrial CT apparatus is used for the image in FIG. 4aThe aluminum block is scanned, and the scanning parameters are configured as follows: the distance (SOD) from the ray source to the center of the turntable is 315.6mm, and the distance (SDD) from the ray source to the detector is 754.0mm; the voltage of the ray source is set to 140kV, the current is set to 0.14mA, and an aluminum filter plate with the thickness of 1.5mm is added in front of the ray source; the number of the line detector units is 1920, and the size of each detector unit is 0.2mm; collecting 1440 angles of data, wherein the exposure time of each angle is 0.2s, and finally obtaining projection data p with the size of 1920 multiplied by 1440 [M] 。
Before applying the method provided by the invention, the projection data p obtained by scanning is needed [M] And (5) processing. First, the SART algorithm is used to directly apply the projection data p [M] Reconstructing to obtain a reconstructed image, as shown in fig. 4 b; then, carrying out threshold segmentation processing on the reconstructed image, and carrying out binarization processing on the reconstructed image, wherein the obtained binarized image is shown in fig. 4 c; performing CT orthographic projection operation on the binarized image according to the same scanning configuration to obtain intersection line length data of the die body of FIG. 4 a; pairing the intersection line length data with the multicolor projection data to obtain a scatter diagram of the corresponding relationship between the intersection line length of the die body and the multicolor projection shown in fig. 4 d; the polychromatic projection data values under the same intersection line length are weighted and averaged, and each intersection line length can obtain a unique polychromatic projection value p, and finally a corresponding relation curve of the intersection line length of the die body and the polychromatic projection is obtained, as shown in fig. 4 e.
Then, the method for estimating the X-ray energy spectrum taking the influence of scattered photons into consideration is adopted for energy spectrum estimation, and the specific implementation steps are as follows:
1) The initial value is adopted to give an initial value of S (E) to the X-ray energy spectrum in the formula (1) (0) The initial value of the scattering is set to Sc (0) ＝0.001。
2) Assume that n iterations have been performed, at which time the estimated values of the X-ray spectrum and scatter are (S (E) (n) ,Sc (n) ). Regarding the formula (1) as (S (E), sc) is (S (E) (n) ,Sc (n) ) Performing first-order taylor expansion, and obtaining a linear equation set according to the processed polychromatic projection value p, wherein the equation set is shown as a formula (2) in the step 130; wherein p is i (n) 、q i (n) The calculation formula of sum phi can refer to formulas (3) - (5)。
3) And (3) according to formulas (6) - (8), the linear equation set (2) is correspondingly in the form of ax=b, the iterative format, namely formula (9), is utilized to solve and update the X-ray energy spectrum and the scattering value, and the updated X-ray energy spectrum and the updated scattering value are normalized according to formula (10).
4) If the X-ray energy spectrum and the scattering value obtained in the step 3) are converged, stopping iteration, and obtaining a final X-ray energy spectrum according to the converged X-ray energy spectrum and the scattering value; otherwise, the next iteration (i.e. n+1 iterations) is performed, and the operations in steps 2) and 3) are repeated until the X-ray energy spectrum and the scatter value converge.
Fig. 5a is a schematic diagram of an initial value of an X-ray energy spectrum to be estimated and an X-ray energy spectrum result obtained by the method according to the present invention in this embodiment. The true X-ray energy spectrum of the actual acquisition data is unknown, and the estimated X-ray energy spectrum result can be used for generating corresponding polychromatic projection values according to the formula (1) and then the polychromatic projection values are matched with the polychromatic projection p obtained by actual acquisition [M] Comparison was performed. Fig. 5b is a schematic diagram of fitting degree between a polychromatic projection curve calculated by using the X-ray energy spectrum estimated in fig. 5a and an actually obtained polychromatic projection curve, and it can be seen that the polychromatic projection curve generated by using the X-ray energy spectrum estimated by the X-ray energy spectrum estimation method considering the influence of scattered photons provided by the invention can be well fitted with actual acquired data, thereby proving feasibility of the method of the invention.
The application of the X-ray energy spectrum estimated by the X-ray energy spectrum estimation method taking the influence of scattered photons into consideration in the removal of hardening artifacts is further described below by way of an example.
Fig. 6a is a photograph of a cylindrical aluminum block phantom of 8cm diameter provided in this example for a hardening artifact removal experiment. The same scanning arrangement is used to scan the phantom of fig. 5a, and the reconstructed image obtained by direct reconstruction using the SART algorithm is shown in fig. 6b, and fig. 6b is also the image before correction of the hardening artifact. The result obtained after the X-ray energy spectrum estimated in fig. 5a is used for curing is shown in fig. 6c, and the line density curve obtained by taking the center line of the image of fig. 6c is shown in fig. 6 d. It can be seen that the hardening artifact in the center part of fig. 6c is significantly improved, and the curve recess originally caused by the ray hardening in fig. 6d is corrected. Fig. 6c and 6d illustrate that the X-ray energy spectrum estimated using the method provided by the present invention may remove hardening artifacts to some extent.
Those of ordinary skill in the art will appreciate that: the drawing is a schematic diagram of one embodiment and the modules or flows in the drawing are not necessarily required to practice the invention.
Finally, it should be pointed out that: the above embodiments are only for illustrating the technical solution of the present invention, and are not limiting. Those of ordinary skill in the art will appreciate that: the technical schemes described in the foregoing embodiments may be modified or some of the technical features may be replaced equivalently; such modifications and substitutions do not depart from the spirit and scope of the technical solutions of the embodiments of the present invention.
Claims (5)
1. An X-ray energy spectrum estimation method considering the influence of scattered photons, comprising:
step 1, setting a scattered photon signal as a constant, and constructing a nonlinear equation set about an X-ray energy spectrum and a scattering constant:
wherein ,pi Is the X-ray polychromatic projection value obtained under the ith ray path; i represents the number of the ray path, the range of values starts from 1, and the maximum value is the product of the number of units of the detector and the number of scanning angles of the object;an index set for an X-ray path; s (E) = (S) 1 ,S 2 ,…,S M ) Is a discrete form of an unknown X-ray energy spectrum, where M represents the maximum energy of the X-ray energy spectrum; delta is the energy spectrum discrete interval; phi (E) = (phi) 1 ,φ 2 ,…,φ M ) Is a discrete form of the mass attenuation coefficient of the measured object; r is R i ＝(r i1 ,r i2 ,…,r iJ ) Is the projection vector of the ith ray path, where r ij Representing the projection contribution of the jth pixel to the ith ray path; f= (f 1 ,f 2 ,…,f J ) T Representing a discretized image, f j A sampling value at the j-th pixel for the image f;
step 2, setting an X-ray energy spectrum S (E) in the nonlinear equation set and an initial value S (E) of the scattering constant Sc (0) and Sc(0) The method comprises the steps of carrying out a first treatment on the surface of the Wherein, the S (E) (0) The Sc is obtained by simulating an X-ray energy spectrum under a preset voltage through X-ray energy spectrum software (0) A preset value within a range of 0.001-0.1;
step 3, at each iteration point (S (E) (n) ,Sc (n) ) Performing first-order taylor expansion on the nonlinear equation set (1) to obtain a linear equation set as follows:
wherein ：
wherein S (E) (n) and Sc(n) Values of S (E) and Sc after the end of the nth iteration are respectively represented;
step 4, converting the equation of the linear equation set (2) into a form ax=b, and solving the linear equation set ax=b by using an EM (Expectation-Maximization) algorithm to obtain estimated values of the X-ray energy spectrum and the scattering constant;
wherein ,
x＝(x 1 ,x 2 ,…,x N ) T ＝(S(E),Sc) T (7)
solving the linear system of equations ax=b using the EM algorithm using the formula:
wherein x represents the unknown number to be solved, x j Represents the j-th unknown number, x j (k) X representing the kth iteration j M represents the number of equations in the linear system of equations, N represents the number of unknowns in the linear system of equations, a ij Representing the coefficient corresponding to the jth unknown of the ith equation, b i A numerical value representing the right-hand term of the i-th equation;
step 5, judging whether the estimated values of the X-ray energy spectrum and the scattering constant obtained in the step 4 are converged, if not, executing the next iteration, and repeating the steps 3 and 4; when the estimated values of the X-ray energy spectrum and the scattering constant are converged, executing a step 6;
and 6, obtaining a final X-ray energy spectrum according to the X-ray energy spectrum and the estimated value of the scattering constant.
2. A method according to claim 1, wherein the initial value of the scattering constant Sc is a preset value in the range 0.001-0.1.
3. The method according to claim 1 or 2, wherein S (E) (0) ＞0，Sc (0) > 0, and:
4. the method according to claim 1 or 2, wherein determining whether the estimated values of the X-ray energy spectrum and the scattering constant obtained in step 4 converge comprises:
calculating an estimated polychromatic projection value p according to the formula (3) using the X-ray energy spectrum and the estimated value of the scattering constant i (n) The method comprises the steps of carrying out a first treatment on the surface of the If p i (n) And the measured polychromatic projection value p i The square of the two norms of the difference is smaller than a given threshold epsilon 1 Judging convergence, namely:
wherein ε1 ＞0。
5. The method according to claim 1 or 2, wherein determining whether the estimated values of the X-ray energy spectrum and the scattering constant obtained in step 4 converge comprises:
the square of the two norms of the difference between the X-ray spectra estimated in two adjacent iterations is less than a given threshold ε 2 Judging convergence, namely:
wherein ε2 ＞0。
Priority Applications (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CN202110568625.8A CN113391341B (en) | 2021-05-25 | 2021-05-25 | X-ray energy spectrum estimation method considering influence of scattered photons |
Applications Claiming Priority (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CN202110568625.8A CN113391341B (en) | 2021-05-25 | 2021-05-25 | X-ray energy spectrum estimation method considering influence of scattered photons |
Publications (2)
| Publication Number | Publication Date |
|---|---|
| CN113391341A CN113391341A (en) | 2021-09-14 |
| CN113391341B true CN113391341B (en) | 2023-08-11 |
Family
ID=77618883
Family Applications (1)
| Application Number | Title | Priority Date | Filing Date |
|---|---|---|---|
| CN202110568625.8A Active CN113391341B (en) | 2021-05-25 | 2021-05-25 | X-ray energy spectrum estimation method considering influence of scattered photons |
Country Status (1)
| Country | Link |
|---|---|
| CN (1) | CN113391341B (en) |
Citations (2)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN101606082A (en) * | 2006-10-27 | 2009-12-16 | 洛斯阿拉莫斯国家安全股份有限公司 | Statistical tomographic reconstruction based on the charged particle measurement |
| EP3223002A1 (en) * | 2016-03-23 | 2017-09-27 | Carl Zeiss X-Ray Microscopy, Inc. | Method, tool and computer program product for measurement and estimation of an x-ray source sprectrum in computed tomography, and x-ray ct system |
-
2021
- 2021-05-25 CN CN202110568625.8A patent/CN113391341B/en active Active
Patent Citations (2)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN101606082A (en) * | 2006-10-27 | 2009-12-16 | 洛斯阿拉莫斯国家安全股份有限公司 | Statistical tomographic reconstruction based on the charged particle measurement |
| EP3223002A1 (en) * | 2016-03-23 | 2017-09-27 | Carl Zeiss X-Ray Microscopy, Inc. | Method, tool and computer program product for measurement and estimation of an x-ray source sprectrum in computed tomography, and x-ray ct system |
Non-Patent Citations (1)
| Title |
|---|
| 基于计算层析成像扫描数据的X射线能谱估计方法;张慧滔;张朋;;光学学报(11);全文 * |
Also Published As
| Publication number | Publication date |
|---|---|
| CN113391341A (en) | 2021-09-14 |
Similar Documents
| Publication | Publication Date | Title |
|---|---|---|
| CN109712212B (en) | Industrial CT artifact correction method | |
| CN107356615B (en) | Method and system for dual-energy X-ray CT | |
| CN111899188B (en) | Neural network learning cone beam CT noise estimation and suppression method | |
| US20200187893A1 (en) | Beam hardening correction in x-ray dark-field imaging | |
| CN113284210B (en) | TGS transmission image reconstruction method, device, medium and electronic equipment | |
| CN104166962A (en) | Cone beam CT scattering correction method by use of scattering nucleus method | |
| TW201801674A (en) | Post processing system and post processing method of electrical impedance tomography images | |
| Yang et al. | Cupping artifacts correction for polychromatic X-ray cone-beam computed tomography based on projection compensation and hardening behavior | |
| CN104783819B (en) | Scatter correction method and device | |
| JP5675978B2 (en) | Computer tomography method, computer program, computer device, and computer tomography system | |
| CN109916933B (en) | X-ray computed tomography energy spectrum estimation method based on convolutional neural network | |
| CN113391341B (en) | X-ray energy spectrum estimation method considering influence of scattered photons | |
| Lifton12 et al. | The application of beam hardening correction for industrial x-ray computed tomography | |
| CN110246199B (en) | Projection domain data noise removing method for energy spectrum CT | |
| CN116609819A (en) | X-ray energy spectrum estimation method with scattering correction | |
| CN107563972B (en) | CT data correction method and step test piece | |
| JP2017221339A (en) | X-ray CT image reconstruction method and computer program | |
| CN113284209B (en) | TGS transmission image reconstruction method, device, medium and electronic equipment | |
| Yu et al. | Heel effect adaptive flat field correction of digital x‐ray detectors | |
| CN113284208A (en) | TGS transmission image reconstruction method, device, medium and electronic equipment | |
| Rahman et al. | Deep learning based workflow for accelerated industrial X-ray Computed Tomography | |
| Duncan | MARS imaging for geology | |
| Roberts et al. | An Overview of 3D X-Ray Reconstruction Algorithms for PCB Inspection | |
| Krumm et al. | Referenceless beam hardening correction in 3d computed tomography images of multi-material objects | |
| Yahaghi | Using total variation denoising for detecting defects in industrial radiography |
Legal Events
| Date | Code | Title | Description |
|---|---|---|---|
| PB01 | Publication | ||
| PB01 | Publication | ||
| SE01 | Entry into force of request for substantive examination | ||
| SE01 | Entry into force of request for substantive examination | ||
| GR01 | Patent grant | ||
| GR01 | Patent grant |