Hardening artifact correction method, system, device and media for images

Abstract

The application discloses a kind of hardening artifact correction method, system, equipment and medium of image, which comprises: obtaining the energy spectrum projection data to be corrected, according to the energy spectrum projection original data obtained by the standard water model of scanning eccentricity, the effective energy spectrum value of the detector channel of energy spectrum CT system in different energy interval is determined;According to the energy spectrum projection data to be corrected and effective energy spectrum value, the correction compensation amount is calculated;The energy spectrum projection data to be corrected is corrected based on the correction compensation amount.The application effectively removes the hardening artifact phenomenon caused by the existence of high atomic number material in the scanned object when the CT device scans the object;Avoid the problem of large amount of calculation and low accuracy caused by image segmentation and forward and backward projection involved in high-order hardening artifact correction method;Reduce the calculation amount of hardening artifact correction, improve the accuracy of hardening artifact correction.

Description

TECHNICAL FIELD

[0001] The present application relates to the technical field of medical devices, and in particular to a method, system, device and medium for correcting a hardened artifact of an image. BACKGROUND

[0002] The X-ray generated by a ball tube in a clinical CT (Computed Tomography) imaging system has a wide energy spectrum. When interacting with a substance, the attenuation amount of low-energy photons is greater than that of high-energy photons. When polychromatic X-rays penetrate an object, low-energy photons are more easily absorbed than high-energy photons. The share of low-energy photons in the rays becomes less and less, the X-ray energy spectrum gradually shifts to the high-energy end, the average energy of the rays becomes high, the energy spectrum becomes "hard", and the beam hardening phenomenon occurs. For a CT imaging system, if hardening correction is not performed, hardening artifacts in the form of uneven bands will appear in the soft tissue region. In order to obtain a CT image that can be used for clinical diagnosis, beam hardening correction needs to be performed during image reconstruction, especially when the scanned object contains high-atomic-number substances such as bone tissue or contrast agent iodine. If hardening correction is not performed, it will cause obvious hardening artifacts in the soft tissue region.

[0003] Traditional hardening correction includes a high-order hardening correction method, which uses water-bone or water-iodine as reference substances for beam hardening correction. The image is segmented to obtain a base substance image, and then forward projection is performed to calculate the deviation caused by hardening by using the projection of the base substance, thereby performing hardening correction. However, this method has a large amount of calculation and low accuracy. SUMMARY

[0004] The technical problem to be solved by the present application is to overcome the defects of large amount of calculation and low accuracy of the high-order hardening correction method in the prior art, and to provide a method, system, device and medium for correcting a hardened artifact of an image.

[0005] The present application solves the above technical problems by the following technical solutions:

[0006] In a first aspect, the present application provides a method for correcting a hardened artifact of an image, which includes:

[0007] Obtaining energy spectrum projection data to be corrected;

[0008] Determining effective energy spectrum values of a detector channel of an energy spectrum CT system in different energy intervals according to energy spectrum projection original data obtained from a standard water phantom placed in a scanning eccentric position;

[0009] Calculating a correction compensation amount according to the energy spectrum projection data to be corrected and the effective energy spectrum values;

[0010] The energy spectrum projection data to be corrected is corrected based on the correction compensation amount.

[0011] Preferably, the step of calculating the correction compensation amount based on the energy spectrum projection data to be corrected and the effective energy spectrum value includes:

[0012] The line integral data of the density of the base material is determined based on the energy spectrum projection data to be corrected;

[0013] The correction compensation amount is calculated based on the effective energy spectrum value, the line integral data, and the linear decay coefficient of each of the base materials.

[0014] Preferably, the step of determining the line integral data of the density of the base material based on the energy spectrum projection data to be corrected includes:

[0015] Determine the quantity of the base material and the quantity of energy ranges; the quantity of energy ranges is greater than the quantity of the base material.

[0016] The line integral data is calculated based on the amount of the base material, the number of energy ranges, and the projection data to be corrected.

[0017] Preferably, the step of calculating the line integral data based on the energy spectrum projection data to be corrected and the effective attenuation coefficient includes:

[0018] The line integral data is calculated using the least squares method or matrix inversion based on the energy spectrum projection data to be corrected and the effective attenuation coefficient.

[0019] Preferably, the step of determining the effective energy spectrum values ​​of the detector channel of the energy spectrum CT system in different energy ranges based on the raw energy spectrum projection data obtained from scanning an eccentrically placed standard water phantom includes:

[0020] Obtain the correspondence between the energy spectrum projection data and the effective energy spectrum data;

[0021] The effective energy spectrum value is obtained by iterative calculation based on the original energy spectrum projection data and the corresponding relationship.

[0022] Preferably, the step of correcting the energy spectrum projection data to be corrected based on the correction compensation amount includes:

[0023] The energy spectrum projection data to be corrected is iteratively calculated based on the correction compensation amount.

[0024] Secondly, the present invention provides an image hardening artifact correction system, the hardening artifact correction system comprising:

[0025] The acquisition module is used to acquire the energy spectrum projection data to be corrected.

[0026] a determination module, configured to determine effective energy spectrum values of detector channels of the spectral CT system in different energy intervals according to energy spectrum projection original data obtained by the standard water phantom placed in the scanning eccentricity;

[0027] a calculation module, configured to calculate a correction compensation according to the to-be-corrected energy spectrum projection data and the effective energy spectrum values;

[0028] a correction module, configured to correct the to-be-corrected energy spectrum projection data based on the correction compensation.

[0029] In a third aspect, the present application provides an electronic device, comprising a processor, a memory, and a computer program stored in the memory and executable on the processor, and when the computer program is executed by the processor, the image hardening artifact correction method is realized.

[0030] In a fourth aspect, the present application provides a computer readable storage medium, and a computer program stored in the computer readable storage medium, and when the computer program is executed by a processor, the image hardening artifact correction method is realized.

[0031] The positive progress effect of the present application is that: the to-be-corrected energy spectrum projection data is obtained; the effective energy spectrum values of the detector channels of the spectral CT system in different energy intervals are determined according to the energy spectrum projection original data obtained by the standard water phantom placed in the scanning eccentricity; the correction compensation is calculated according to the to-be-corrected energy spectrum projection data and the effective energy spectrum values; and the to-be-corrected energy spectrum projection data is corrected based on the correction compensation. The present application effectively removes the hardening artifact phenomenon caused by the existence of high atomic number materials in the scanned object when the CT device scans the object; avoids the problem of large calculation amount and low accuracy caused by the involvement of image segmentation and forward and backward projection in the high-order hardening artifact correction method; reduces the calculation amount of the hardening artifact correction, and improves the accuracy of the hardening artifact correction. BRIEF DESCRIPTION OF DRAWINGS

[0032] Figure 1 The flowchart of the image hardening artifact correction method of the present application embodiment 1.

[0033] Figure 2 The flowchart of step S3 of the image hardening artifact correction method of the present application embodiment 1.

[0034] Figure 3 The flowchart of step S31 of the image hardening artifact correction method of the present application embodiment 1.

[0035] Figure 4 The flowchart of step S2 of the image hardening artifact correction method of the present application embodiment 1.

[0036] Figure 5 Fig. 1 is a first application example diagram of the hardening artifact correction method of the image of embodiment 1 of the present application.

[0037] Figure 6 Fig. 2 is a second application example diagram of the hardening artifact correction method of the image of embodiment 1 of the present application.

[0038] Figure 7 Fig. 3 is a first module schematic diagram of the hardening artifact correction system of the image of embodiment 2 of the present application.

[0039] Figure 8 Fig. 4 is a second module schematic diagram of the hardening artifact correction system of the image of embodiment 2 of the present application.

[0040] Figure 9 Fig. 5 is a structure schematic diagram of the electronic device for implementing the hardening artifact correction method of the image of embodiment 3 of the present application. DETAILED DESCRIPTION

[0041] The present application will be further described below by way of examples, but the present application is not limited in the scope of the described examples.

[0042] Embodiment 1

[0043] The hardening artifact correction method of the image of the present embodiment, as shown in Fig. 1, comprises: Figure 1

[0044] S1, acquiring the to-be-corrected spectral projection data.

[0045] S2, determining the effective spectral values of the detector channels of the spectral CT system in different energy intervals according to the spectral projection original data obtained by the standard water phantom placed in the scanning eccentric position.

[0046] S3, calculating the correction compensation amount according to the to-be-corrected spectral projection data and the effective spectral values.

[0047] S4, correcting the to-be-corrected spectral projection data based on the correction compensation amount.

[0048] ​The hardening artifact correction method provided by the application is applicable to a cone beam CT device, i.e., a CBCT (Cone beam Computed Tomography) device. The CT device mainly comprises an X-ray source, a detector and a rotating device. The hardening artifact correction method is applied to a scenario of imaging by scanning an object by using the CT device. When scanning a high-atomic-number substance such as bone tissue or contrast agent iodine, the X-rays emitted by the X-ray source are wide energy bands, and the X-ray attenuation coefficients of substances with the same density change with the intensity of X-rays, resulting in the phenomenon of hardening artifacts in soft tissue regions or the phenomenon of uneven strip substances between dense substances. In the CT imaging process, the X-ray source and the detector rotate relative to the object, thereby obtaining CT projection data, also referred to as CT projection values, at different rotation angles.

[0049] For the above step S1, the target object is placed on a scanning bed of the cone beam CT device, and the target object is scanned by using multi-energy rays. In order to reduce the total scanning time, a "spiral" scan can be performed on the target object, so that the angle at which the radiation beam intersects the target object changes constantly to obtain initial energy spectrum projection data. After water hardening correction processing is performed on the initial energy spectrum projection data, the energy spectrum projection data to be corrected is obtained, which is marked as p b Water hardening correction equivalent to the scanned object is equivalent to water, and the projection value of X-rays with a specific energy spectrum is corrected into the projection value of ideal single-energy X-rays. This correction is prior art known to those skilled in the art, and will not be described here.

[0050] For the above step S2, it is considered that the cone beam CT device is mainly used for scanning human body parts, and the density of human soft tissue is close to the density of water. Therefore, the material of the phantom can be selected to be close to the density of water, and the pre-set density liquid can be selected to be water or a liquid with a density close to water for water phantom experiments. It should be noted that in the present embodiment, when designing a standard water phantom, the size of the flat panel detector of the cone beam CT device and the actual required scanning thickness and other parameters are designed to ensure that the flat panel detector can receive all the X-rays penetrating the phantom.

[0051] The standard uniform cylindrical phantom of a preset size is placed eccentrically on a scanning bed of a cone beam CT device, and a multi-energy ray is used to perform a rotating scan on the standard phantom filled with a liquid of a preset density to obtain energy spectrum projection original data, which is multi-energy projection data. After the standard phantom placed eccentrically is scanned, a water phantom image is reconstructed by using full-energy data. The water phantom image is divided into a water phantom wall and water by threshold segmentation. The water phantom wall is made of polymethyl methacrylate (PMMA), a high-molecular polymer. The segmented image is assigned with density values of PMMA and water, respectively. After 3D cone beam forward projection is performed on the eccentrically placed cylindrical water phantom by using the system geometry of the CT device during the scan, a penetration thickness of the water phantom wall and a penetration thickness of the water through the base material are obtained, and linear attenuation coefficients of the water phantom wall and the water are determined. The linear attenuation coefficient of an object can represent the relationship between the X-ray intensity incident on the object and the X-ray intensity out of the object, and can be related to the thickness, area, shape and the like of the object. The correspondence between the energy spectrum projection original data and the effective energy spectrum data is determined according to the energy spectrum projection original data, the penetration thicknesses of the water phantom wall and the water through the base material, and the linear attenuation coefficients. This way reduces the calculation amount of the hardening artifact correction and improves the accuracy of the hardening artifact correction.

[0052] The energy spectrum projection original data can be estimated by using an iterative algorithm, and the effective energy spectrum values of each energy interval of each detector unit of the cone beam CT device are calculated, but the effective energy spectrum values of the edge channels at different energy intervals need to be obtained by extension based on the effective energy spectrum of the middle region. Since the effective energy spectrum values of each energy interval of the detector unit need to be determined by using the eccentric water phantom, but the eccentric water phantom has certain limitations, the edge channels cannot be calibrated, and therefore the effective energy spectrum of the edge channels can be transformed and extended based on the middle region.

[0053] For the above step S3, when there are multiple base substances, the effective attenuation coefficients of each base substance are determined, and the energy spectrum projection data p i is corrected according to the effective attenuation coefficients of each base substance and the energy spectrum projection data p b , the subscript b is an energy interval index, and the line integral data of each base substance and the effective energy spectrum values of the detector channels at different energy intervals are iteratively calculated to obtain a correction compensation amount, which is the calculated projection value after water hardening correction minus the projection value in the previous step.

[0054] For the above step S4, the energy spectrum projection data p b is iteratively calculated according to the correction compensation amount, and the value of the bth energy interval is denoted as pb (0) , p b (0) As the initial value of the energy spectrum projection data to be corrected, R w (·) represents the process of water hardening correction. The water hardening correction is a prior art known to those skilled in the art, which corrects the cup-shaped artifact of the water phantom by correcting the nonlinear relationship of the projection value p with the thickness L of the water passing through into a linear relationship, and the specific water hardening correction process is not described here. The iterative correction process is as follows:

[0055]

[0056] Specifically, when the initial value of the energy spectrum projection data to be corrected p is calculated by the first iteration, when the second iteration is calculated, each time is based on the initial value of the energy spectrum projection data to be corrected p b (0) The iterative calculation is performed, and so on to realize multiple corrections of the energy spectrum projection data to be corrected. This method effectively removes the hardening artifact phenomenon caused by the existence of high atomic number materials in the scanned object when the CT device scans the object.

[0057] In an optional embodiment, as Figure 2 shown, step S3 specifically includes:

[0058] S31, determining the line integral data of the density of the base material according to the energy spectrum projection data to be corrected.

[0059] S32, calculating the correction compensation according to the effective energy spectrum value, the line integral data and the linear attenuation coefficient of each base material.

[0060] For the above steps S31-S32, the following formula is used to calculate the correction compensation projection value of the selected X-ray energy interval lower limit value E min and the selected X-ray energy interval upper limit value E max according to the effective energy spectrum value, the line integral data of the density of each base material, and the linear attenuation coefficient of each base material. This method avoids the problem of large calculation amount and low accuracy caused by image segmentation and forward and backward projection involved in high-order hardening artifact correction method.

[0061]

[0062] Wherein, j represents the energy sampling point in the energy interval from Emin to Emax, S jrepresents the effective energy spectrum value at the energy sampling point j, L1 (n) represents the line integral data of the first kind of base material in the n th iteration calculation, L2 (n) represents the line integral data of the second kind of base material in the n th iteration calculation, E j represents the set energy interval, μ1(E j represents the linear attenuation coefficient of the first kind of base material at the energy sampling point j, which is a function of the energy value of the energy sampling point j, μ2(E j represents the linear attenuation coefficient of the second kind of base material at the energy sampling point j, p cal (n) represents the calculated projection value.

[0063] In an optional embodiment, as shown in Figure 3 step S31 specifically includes:

[0064] S311, determine the number of base materials and the number of energy intervals; the number of energy intervals is greater than the number of base materials.

[0065] S312, calculate the line integral data according to the number of base materials, the number of energy intervals, and the projection data to be corrected.

[0066] For the above steps S311-S312, M represents the number of selected base materials, which can be water / iodine decomposition, water / bone decomposition, or water / calcium decomposition, etc. After the upper and lower limit values of the energy range are determined, the energy range is divided into a preset B number of energy intervals. In this embodiment, the data value of B is greater than the data value of M. The projection value of each energy bin and the model of the base material line integral are as follows:

[0067]

[0068] wherein p represents a vector matrix form of [Bx1], and the component p b represents the projection value of the b th energy interval, represents a matrix form of [BxM], L represents a vector form of [Mx1], and the component L m represents the line integral of the m th base material; represents the effective attenuation coefficient, and the component represents the effective attenuation coefficient of the m th base material in the b th energy interval, which can be obtained by phantom calibration. For example, when the number of base materials is two and there are two energy intervals, the above formula can be represented as follows:

[0069]

[0070] In an alternative embodiment, step S312 specifically comprises:

[0071] S3121, determining the effective attenuation coefficient according to the number of base substances and the range of energy intervals;

[0072] S3122, calculating the line integral data according to the energy spectrum projection data to be corrected and the effective attenuation coefficient.

[0073] For steps S3121-S3122 above, the line integral data of the density of each base substance can be calculated by the following formula, n represents the number of iterations, and the matrix inversion is used to calculate the line integral data according to the energy spectrum projection data to be corrected and the effective attenuation coefficient of each base substance in each energy interval.

[0074]

[0075] In an alternative embodiment, the line integral data can be calculated by the least square method or matrix inversion according to the energy spectrum projection data to be corrected and the effective attenuation coefficient, and other algorithms can also be used to calculate the line integral data, which is not specifically limited in the present embodiment.

[0076] In an alternative embodiment, as shown in Figure 4 step S2 specifically comprises:

[0077] S21, obtaining the corresponding relationship between the energy spectrum projection data and the effective energy spectrum data.

[0078] S22, performing iterative calculation based on the energy spectrum projection original data and the corresponding relationship to obtain the effective energy spectrum value.

[0079] For step S21 above, when the effective energy spectrum data is calibrated by scanning the standard water phantom placed eccentrically to establish the corresponding relationship between the energy spectrum projection data and the effective energy spectrum data, for a certain detector unit of the CT device, the incident rays will pass through the water phantom wall and water of different thicknesses in the scanning process. After projection of the base material water phantom wall and water from multiple viewing angles, the projection data of each viewing angle is sampled, p i represents the energy spectrum projection data of a certain energy interval of the i-th sampling point, S j represents the effective energy spectrum data of the CT device, and j represents the energy spectrum index, which can be reasonably sampled. For example, the energy interval can be 5 keV or other values. According to the effective energy spectrum data S j and the energy spectrum projection data p i the corresponding relationship is established as shown in the following formula:

[0080]

[0081] wherein, E​j represents the energy interval, μ P represents the linear attenuation coefficient of the water phantom wall, μ w represents the linear attenuation coefficient of water, L P represents the path integral value (penetration thickness value) of the water phantom wall, L w represents the path integral value (penetration thickness value) of water, A ji represents the attenuation effect value of the selected energy interval E j under the water phantom wall of the base material and the water combination of the base material. E min represents the lower limit value of the energy interval of the selected X-ray, E max represents the upper limit value of the energy interval of the selected X-ray, which is composed of E min and E max , the range of the selected X-ray energy interval, p i , L P and L w represent known quantities.

[0082] For the above step S22, the energy spectrum projection original data and the attenuation effect value A ji are mathematically operated by the above formula to obtain the effective energy spectrum data S j After that, the EM algorithm is used to iteratively calculate the obtained energy spectrum projection original data P i and the effective energy spectrum data S j , k represents the iteration number, and the initial effective energy spectrum value is set to is updated by the following iterative formula until the iteration stopping condition is reached, and the effective energy spectrum value is obtained according to the final iteration result. For example, after integrating calculation is performed at an accuracy of every 1 keV in the energy interval of 30 keV-60 keV, the current effective energy spectrum value is obtained as a function of the result of the previous step . After the effective energy spectrum data S j is iterated for multiple times by the following formula, more accurate energy spectrum data, i.e., the effective energy spectrum value, is obtained. It can be understood that the iteration number can be set according to actual conditions, and preferably, the effective energy spectrum data S j is iterated for an even number of times for correction.

[0083] It should be noted that the EM algorithm is a prior art known to those skilled in the art, and the energy spectrum value of the incident X-ray is estimated by measuring the X-ray attenuation data, which will not be described here. Other algorithms can also be used to iteratively calculate the energy spectrum projection original data p i , which is not specifically limited in the present embodiment.

[0084]

[0085] Where n represents the selected energy range S j The energy subscript of the effective energy spectrum value represents the energy spectrum precision, and k represents the number of iterations.

[0086] In one alternative implementation, such as Figure 5 As shown, during the effective energy spectrum calibration process using an eccentric water phantom, a CT energy spectrum scan is performed on the eccentric water phantom to obtain the raw energy spectrum projection data. After reconstructing the water phantom image obtained after scanning using full-bin data, an appropriate threshold is selected to segment the water phantom image. For example, when the CT value of the water phantom image is greater than the first threshold, the water phantom wall image is segmented, while when the CT value of the water phantom image is less than the second threshold, the water phantom image is segmented. The line integral data of X-rays passing through the water phantom wall and water density of each detector channel of the CT equipment are calculated. Using the EM algorithm, iterative calculations are performed based on the line integral data and the raw energy spectrum projection data. When the iteration conditions are met, the effective energy spectrum value of each energy range of each detector channel is output.

[0087] In one alternative implementation, such as Figure 6 As shown, in the process of multi-material artifact hardening correction, pre-processed energy spectrum scanning data is acquired. The line integral data of each base material is calculated using the least squares method or matrix inversion. Then, the line integral data of the X-ray density passing through the water phantom wall and water density at each detector channel of the CT scanner are determined. The projection value is calculated using the calibrated effective energy spectrum value and the line integral data of each base material. The hardening artifact error is subtracted from the original energy spectrum projection data, and this process is iterated multiple times until the iteration conditions are met. Finally, the corrected energy spectrum projection data is output. This method effectively removes the hardening artifact phenomenon caused by the presence of high atomic number materials in the scanned object when the CT scanner scans it.

[0088] This embodiment provides a method for correcting hardening artifacts in images. The method involves acquiring the energy spectrum projection data to be corrected; determining the effective energy spectrum values ​​of the detector channels in different energy ranges of the energy spectrum CT system based on the original energy spectrum projection data obtained from scanning a standard water phantom placed off-center; calculating a correction compensation amount based on the energy spectrum projection data to be corrected and the effective energy spectrum values; and correcting the energy spectrum projection data to be corrected based on the correction compensation amount. This invention effectively removes hardening artifacts caused by the presence of high atomic number materials in the scanned object during CT scanning; avoids the problems of high computational load and low accuracy caused by image segmentation and orthographic / reverse projection required by higher-order hardening artifact correction methods; reduces the computational load of hardening artifact correction; and improves the accuracy of hardening artifact correction.

[0089] Example 2

[0090] The hardened artifact correction system of the image of the embodiment, as shown in Figure 7 includes an acquisition module 110, a determination module 120, a calculation module 130 and a correction module 140.

[0091] The acquisition module 110 is configured to acquire the energy spectrum projection data to be corrected.

[0092] The determination module 120 is configured to determine the effective energy spectrum values of the detector channels of the energy spectrum CT system in different energy intervals according to the energy spectrum projection original data obtained by the standard water phantom placed in the scanning eccentricity.

[0093] The calculation module 130 is configured to calculate the correction compensation according to the energy spectrum projection data to be corrected and the effective energy spectrum values.

[0094] The correction module 140 is configured to correct the energy spectrum projection data to be corrected based on the correction compensation.

[0095] The hardened artifact correction system provided by the present application is applicable to a cone beam CT device, i.e. a CBCT (Cone beam Computed Tomography) device. The CT device mainly includes an X-ray source, a detector and a rotating device. The hardened artifact correction method is applied to a scene of imaging an object by using the CT device to scan. When scanning a high atomic number substance such as bone tissue or contrast agent iodine, the X-rays emitted by the X-ray source are wide energy bands, and the X-ray attenuation coefficients of substances with the same density change with the X-ray intensity, resulting in the phenomenon of hardened artifacts in soft tissue regions or the phenomenon of uneven strip substances between dense substances. In the CT imaging process, the X-ray source and the detector rotate relative to the object, thereby obtaining CT projection data, also known as CT projection values, at different rotation angles.

[0096] The target object is placed on the scanning bed of the cone beam CT device, and the target object is scanned by using multi-energy rays. In order to reduce the total scanning time, a "spiral" scan can be performed on the target object, so that the angle at which the radiation beam intersects the target object changes constantly to obtain energy spectrum projection initial data. The acquisition module 110 performs water hardening correction processing on the energy spectrum projection initial data to obtain the energy spectrum projection data to be corrected, denoted as p b Water hardening correction equivalent to the scanned object is equivalent to water, and the projection value of X-rays with a specific energy spectrum is corrected into the projection value of ideal single-energy X-rays. This correction is a prior art known to those skilled in the art, and will not be described here.

[0097] It is considered that the cone beam CT device is mainly used for scanning human body parts, and the density of human soft tissue is close to the density of water, so the material of the phantom can be selected to be close to the density of water, and the preset density liquid can be selected to be water or a liquid with a density close to water for water phantom experiment. It should be noted that in this embodiment, when designing the standard water phantom, the size of the flat panel detector of the cone beam CT device and the actual required scanning thickness and other parameters are designed to ensure that the flat panel detector can receive all the X-rays penetrating the phantom.

[0098] A standard uniform cylindrical phantom of a preset size is placed eccentrically on the scanning bed of the cone beam CT device, and a multi-energy ray is used to rotate scan the standard phantom filled with a liquid of a preset density to obtain energy spectrum projection original data, which is multi-energy projection data. After scanning the eccentrically placed standard phantom, a water phantom image is reconstructed using full-energy data. The water phantom image is divided into a water phantom wall and water through threshold segmentation. The material of the water phantom wall is polymethyl methacrylate, a high molecular polymer, which is referred to as PMMA. The segmented image is assigned a density value of PMMA and water, respectively. After 3D cone beam forward projection of the eccentrically placed cylindrical water phantom using the system geometry of the CT device during scanning, the penetration thickness of the base material water phantom wall and the penetration thickness of the base material water are obtained, and the linear attenuation coefficients of the water phantom wall and the water are determined. The linear attenuation coefficient of an object can represent the relationship between the incident X-ray intensity of the object and the intensity of the X-ray exiting the object. The linear attenuation coefficient of an object can be related to the thickness, area, shape, etc. of the object. The correspondence between the energy spectrum projection original data and the effective energy spectrum data is determined according to the energy spectrum projection original data, the penetration thickness of the base material water phantom wall and the base material water, and the linear attenuation coefficients. This method reduces the calculation amount of the hardening artifact correction and improves the accuracy of the hardening artifact correction.

[0099] The iterative algorithm can be used to estimate the energy spectrum projection original data, and the determination module 120 calculates the effective energy spectrum value of each energy interval of each detector unit of the cone beam CT device, but the effective energy spectrum value of the edge channel under different energy intervals needs to be obtained by extension based on the effective energy spectrum of the middle region. Since the effective energy spectrum value of each energy interval of the detector unit needs to be calibrated by the eccentric water phantom, but the eccentric water phantom has certain limitations, the edge channel cannot be calibrated, so the effective energy spectrum of the edge channel can be transformed and extended based on the middle region.

[0100] When there are multiple base substances, the effective attenuation coefficient of each base substance is determined, and the effective attenuation coefficient of each base substance and the energy spectrum projection data p iThe line integral data of each base substance is calculated by matrix inversion. The calculation module 130 calculates the correction compensation quantity according to the energy spectrum projection data p i to be corrected and the line integral data of each base substance and the effective energy spectrum value of the detector channel in different energy intervals.

[0101] The correction module 140 corrects the energy spectrum projection data p b to be corrected according to the correction compensation quantity by using the following formula: b (0) p b (0) The initial value of the energy spectrum projection data to be corrected is p w (·) represents the process of water hardening correction. The calculated projection value measured in different energy intervals approaches the calculated projection value of the ideal monochromatic spectrum and iodine water. The water hardening correction is a prior art known to those skilled in the art. By correcting the nonlinear relationship between the projection value p and the thickness L of the water through which the projection value p passes into a linear relationship, the cup-shaped artifact of the water phantom is corrected. The specific water hardening correction process is not described here. The iterative correction process is as follows:

[0102]

[0103] Specifically, when the initial value p of the energy spectrum projection data to be corrected is subjected to the first iteration calculation by using the above formula, when the initial value p of the energy spectrum projection data to be corrected is subjected to the second iteration calculation, each time is based on the initial value p b (0) of the energy spectrum projection data to be corrected. Iterative calculation is performed, and so on to realize multiple corrections of the energy spectrum projection data to be corrected. This method effectively removes the hardening artifact phenomenon caused by the existence of high atomic number materials in the scanned object when the CT device scans the object.

[0104] In an optional embodiment, as shown in Figure 8 , the calculation module 130 specifically includes:

[0105] The determination unit 131 is configured to determine the line integral data of the density of the base substance according to the energy spectrum projection data to be corrected.

[0106] The calculation unit 132 is configured to calculate the correction compensation quantity according to the effective energy spectrum value, the line integral data, and the linear attenuation coefficient of each base substance.

[0107] The calculation unit 132 calculates the projection value of the correction compensation amount according to the effective energy spectrum value, the linear integral data of the density of each base material, and the linear attenuation coefficient of each base material by using the following formula, to obtain the projection value of the correction compensation amount in the lower limit value E min and the upper limit value E max of the energy range of the selected X-ray. In this way, the problem of large amount of calculation and low accuracy caused by image segmentation and forward and back projection involved in the high-order hardening artifact correction method is avoided.

[0108]

[0109] wherein j represents the energy sampling point in the energy range of Emin to Emax, S j represents the effective energy spectrum value at the energy sampling point j, L1 (n) represents the linear integral data of the first base material in the n-th iteration calculation, L2 (n) represents the linear integral data of the second base material in the n-th iteration calculation, E j represents the set energy range, μ1(E j ) represents the linear attenuation coefficient of the first base material at the energy sampling point j, which is a function of the energy value of the energy sampling point j, μ2(E j ) represents the linear attenuation coefficient of the second base material at the energy sampling point j, which is a function of the energy value of the energy sampling point j, p cal (n) represents the calculated projection value.

[0110] In an optional embodiment, the determination unit 131 is specifically configured to:

[0111] determine the number of base materials and the number of energy ranges; the number of energy ranges is greater than the number of base materials.

[0112] calculate the linear integral data according to the number of base materials, the number of energy ranges, and the projection data to be corrected.

[0113] M represents the number of selected base materials, which can be water / iodine decomposition, water / bone decomposition, or water / calcium decomposition, etc. After the upper limit value and the lower limit value of the energy range are determined, the energy range is divided into a preset B number of energy ranges. In this embodiment, the data value of B is greater than the data value of M, and the determination unit 131 calculates the linear integral data of the density of each base material according to the following formula:

[0114]

[0115] wherein p represents a vector matrix form of [B x 1], and the component p b is the projection value of the b-th energy range, The matrix form of [B x M] represents, L represents the vector form of [M x 1], and the components of L are m denotes the line integral of the m-th base material; denotes the effective attenuation coefficient, and the components of denotes the effective attenuation coefficient of the m-th base material in the b-th energy interval, which can be obtained by phantom calibration. Exemplarily, when the number of base materials is two and there are two energy intervals, the above formula can be expressed as follows:

[0116]

[0117] In an optional embodiment, step S312 specifically comprises:

[0118] S3121, determining the effective attenuation coefficient according to the number of base materials and the range of energy intervals;

[0119] S3122, calculating the line integral data according to the energy spectrum projection data to be corrected and the effective attenuation coefficient.

[0120] For steps S3121-S3122, the line integral data of the density of each base material can be calculated by the following formula, n represents the number of iterations, and the inverse matrix is used to calculate the inverse of the energy spectrum projection data to be corrected and the effective attenuation coefficient of each base material in each energy interval

[0121]

[0122] In an optional embodiment, the line integral data can be calculated by using the least square method or matrix inversion according to the energy spectrum projection data to be corrected and the effective attenuation coefficient, and other algorithms can also be used to calculate the line integral data, which is not specifically limited in the present embodiment.

[0123] In an optional embodiment, as shown in Figure 8 the determining module 120 specifically comprises:

[0124] The acquisition unit 121 is configured to acquire the corresponding relationship between the energy spectrum projection data and the effective energy spectrum data.

[0125] The iterative calculation unit 122 is configured to perform iterative calculation based on the energy spectrum projection original data and the corresponding relationship to obtain the effective energy spectrum value.

[0126] ​In the process of calibrating the effective energy spectrum data by scanning the standard water phantom placed eccentrically to establish the correspondence between the energy spectrum projection data and the effective energy spectrum data, for a certain detector unit of the CT device, the incident rays will pass through the water phantom wall and water with different thicknesses during the scanning process. After projection of the base material water phantom wall and water from multiple viewing angles, the projection data of each viewing angle is sampled correspondingly, p i represents the energy spectrum projection data of the i-th sampling point, S j represents the effective energy spectrum data of the CT device, j represents the energy spectrum index, and j can be reasonably sampled. For example, the energy interval can be 5 keV or other values, and the acquisition unit 121 obtains the effective energy spectrum data S j and the energy spectrum projection data p i The correspondence is established as shown in the following formula:

[0127]

[0128] wherein, E j represents the energy interval, μ P represents the linear attenuation coefficient of the water phantom wall, μ w represents the linear attenuation coefficient of the water, L P represents the path integral value (penetration thickness value) of the water phantom wall, L w represents the path integral value (penetration thickness value) of the water. A ji represents the attenuation effect value of the base material water phantom wall and the base material water combination under the selected energy interval E j E min represents the lower limit value of the energy interval of the selected X-ray, E max represents the upper limit value of the energy interval of the selected X-ray, which constitutes the range of the selected X-ray energy interval with E min and E max , wherein p i , L P and L w represent known quantities.

[0129] The iterative calculation unit 122 performs mathematical operations on the energy spectrum projection raw data and the attenuation effect value A ji through the foregoing formula to obtain the effective energy spectrum data S j After that, the EM algorithm is used to perform iterative calculation on the obtained energy spectrum projection raw data P i and the effective energy spectrum data S j , wherein k represents the number of iterations, and the initial effective energy spectrum value is set as The updating is performed by using the following iterative formula until the iteration stopping condition is reached, and the effective energy spectrum value is obtained according to the final iteration result. For example, after integral calculation is performed in the energy interval of 30 keV-60 keV with the accuracy of 1 keV, the current effective energy spectrum value is obtained The result is a function of the result of the previous step After the effective energy spectrum data S j is iterated for multiple times by using the following formula, more accurate energy spectrum data, i.e., the effective energy spectrum value, is obtained. It can be understood that the iteration times can be set according to actual conditions, and preferably, the effective energy spectrum data S j is iterated for even times for correction.

[0130] It should be noted that the EM algorithm is a prior art known to those skilled in the art, the energy spectrum value of incident X-rays is estimated by measuring the X-ray attenuation data, and details are not described herein. Other iterative algorithms can also be used to iteratively calculate the energy spectrum projection original data P i , and the present embodiment is not limited in this regard.

[0131]

[0132] wherein n represents the energy index of the effective energy spectrum value in the selected energy range S j , represents the energy spectrum accuracy, and k represents the iteration times.

[0133] In the present embodiment, a system for correcting hardening artifacts of an image is provided, an acquisition module acquires energy spectrum projection data to be corrected; a determination module determines effective energy spectrum values of a detector channel of an energy spectrum CT system in different energy intervals according to energy spectrum projection original data obtained by a standard water phantom placed in a scanning eccentric position; a calculation module calculates a correction compensation amount according to the energy spectrum projection data to be corrected and the effective energy spectrum values; and a correction module corrects the energy spectrum projection data to be corrected based on the correction compensation amount. The present application effectively removes the hardening artifact phenomenon caused by the existence of high atomic number materials in the scanned object when the CT device scans the object; avoids the problem of large calculation amount and low accuracy caused by image segmentation and forward and backward projection involved in the high-order hardening artifact correction method; reduces the calculation amount of hardening artifact correction, and improves the accuracy of hardening artifact correction.

[0134] Embodiment 3

[0135] Figure 9 A structural schematic diagram of an electronic device is provided for the present embodiment. The electronic device includes a memory, a processor, and a computer program stored on the memory and executable on the processor, and the processor implements the image hardening artifact correction method of embodiment 1 when executing the program. Figure 9The electronic device 90 shown is merely one example and should not be taken as limiting the scope of functionality of various embodiments of the present application.

[0136] As Figure 9 shown, the electronic device 90 can be a form of a general computing device, such as a server device. The components of the electronic device 90 can include, but are not limited to, the at least one processor 91 described above, the at least one memory 92 described above, and a bus 93 that connects the various system components, including the memory 92 and the processor 91.

[0137] The bus 93 includes a data bus, an address bus, and a control bus.

[0138] The memory 92 can include volatile memory, such as random access memory (RAM) 921 and / or cache memory 922, and can further include non-volatile memory, such as read only memory (ROM) 923.

[0139] The memory 92 can also include a program / utility 925 having a set of program modules 924, including but not limited to an operating system, one or more application programs, other program modules, and program data, each of which can include implementation of a network environment, individually or in some combination.

[0140] The processor 91 can execute the various functions and data processing of the embodiments of the present application by running computer programs stored in the memory 92.

[0141] The electronic device 90 can also communicate with one or more external devices 94 such as a keyboard or a pointing device, through an input / output (I / O) interface 95. The I / O interface 95 can also include a visual display device, such as a monitor or screen, for displaying information to a user. Figure 9 As shown, the network adapter 96 communicates with the other modules of the model generation device 90 through the bus 93. It should be understood that although not shown, other hardware and / or software components could be used in conjunction with the model generation device 90. These include, but are not limited to, microcode, device drivers, redundant processing units, external disk drive arrays, RAID systems, tape drives, and data archival storage systems, etc.

[0142] It should be noted that although several units / modules or sub-units / modules of the electronic device are mentioned in the above detailed description, such division is merely exemplary and not mandatory. Indeed, according to embodiments of the application, the features and functionalities of two or more units / modules described above can be embodied in one unit / module. Conversely, the features and functionalities of one unit / module described above can be further divided into units / modules embodied by several units / modules.

[0143] Embodiment 4

[0144] The embodiment provides a computer readable storage medium, and a computer program is stored on the computer readable storage medium. The program is executed by a processor to implement steps in the image hardening artifact correction method of embodiment 1.

[0145] More specifically, the readable storage medium can include, but is not limited to, a portable disc, a hard disk, a random access memory, a read-only memory, an erasable programmable read-only memory, an optical storage device, a magnetic storage device, or any suitable combination of the above.

[0146] In possible embodiments, the application can also be implemented in the form of a program product, which includes program codes for causing a terminal device to execute steps in the image hardening artifact correction method of embodiment 1 when the program product is run on the terminal device.

[0147] The program codes for executing the application can be written in any combination of one or more programming languages, and can be executed completely on a user device, partially on a user device, as a separate software package, partially on a user device and partially on a remote device, or completely on a remote device.

[0148] Although the specific embodiments of the application are described above, those skilled in the art should understand that this is only an example, and the protection scope of the application is defined by the appended claims. Those skilled in the art can make various changes or modifications to these embodiments without departing from the principles and essence of the application, and such changes and modifications fall within the protection scope of the application.

Claims

1. A method of hardening artifact correction of an image, characterized by, The hardening artifact correction method includes: Obtain the energy spectrum projection data to be corrected; Based on the raw data of the energy spectrum projection obtained from scanning the eccentrically placed standard water model, the effective energy spectrum values ​​of the detector channel of the energy spectrum CT system in different energy ranges are determined. The correction compensation amount is calculated based on the energy spectrum projection data to be corrected and the effective energy spectrum value. The energy spectrum projection data to be corrected is corrected based on the correction compensation amount; The step of calculating the correction compensation amount based on the energy spectrum projection data to be corrected and the effective energy spectrum value includes: The line integral data of the density of the base material is determined based on the energy spectrum projection data to be corrected; The correction compensation amount is calculated based on the effective energy spectrum value, the line integral data, and the linear attenuation coefficient of each of the base materials; The step of determining the line integral data of the density of the base material based on the energy spectrum projection data to be corrected includes: The quantity of the base material and the quantity of energy ranges are determined; the quantity of energy ranges is greater than the quantity of the base material; once the upper and lower limits of the energy range are determined, the energy range is divided into a preset number of energy ranges. The line integral data is calculated based on the amount of the base material, the number of energy ranges, and the projection data to be corrected.

2. The method of hardening artifact correction of an image according to claim 1, wherein, The step of calculating the line integral data based on the quantity of the base material, the quantity of the energy range, and the projection data to be corrected includes: The effective attenuation coefficient is determined based on the quantity of the base material and the range of the energy range; The line integral data is calculated based on the energy spectrum projection data to be corrected and the effective attenuation coefficient.

3. The method of hardening artifact correction of an image according to claim 2, wherein, The step of calculating the line integral data based on the energy spectrum projection data to be corrected and the effective attenuation coefficient includes: The line integral data is calculated using the least squares method or matrix inversion based on the energy spectrum projection data to be corrected and the effective attenuation coefficient.

4. The method of hardening artifact correction of an image according to claim 1, wherein, The step of determining the effective energy spectrum values ​​of the detector channel of the energy spectrum CT system in different energy ranges based on the raw energy spectrum projection data obtained from scanning an eccentrically placed standard water phantom includes: Obtain the correspondence between the energy spectrum projection data and the effective energy spectrum data; The effective energy spectrum value is obtained by iterative calculation based on the original energy spectrum projection data and the corresponding relationship.

5. The image hardening artifact correction method as described in claim 1, characterized in that, The step of correcting the energy spectrum projection data to be corrected based on the correction compensation amount includes: The energy spectrum projection data to be corrected is iteratively calculated based on the correction compensation amount.

6. A system for correcting hardening artifacts in images, characterized in that, The hardening artifact correction system includes: The acquisition module is used to acquire the energy spectrum projection data to be corrected. The determination module is used to determine the effective energy spectrum values ​​of the detector channels of the energy spectrum CT system in different energy ranges based on the raw energy spectrum projection data obtained from scanning the eccentrically placed standard water model. The calculation module is used to calculate the correction compensation amount based on the energy spectrum projection data to be corrected and the effective energy spectrum value; A correction module is used to correct the energy spectrum projection data to be corrected based on the correction compensation amount; The calculation module specifically includes: The determination unit is used to determine the line integral data of the density of the base material based on the energy spectrum projection data to be corrected. The calculation unit is used to calculate the correction compensation amount based on the effective energy spectrum value, the line integral data, and the linear attenuation coefficient of each of the base materials; The determining unit is specifically used for: The quantity of the base material and the quantity of energy ranges are determined; the quantity of energy ranges is greater than the quantity of the base material; once the upper and lower limits of the energy range are determined, the energy range is divided into a preset number of energy ranges. Linear integral data is calculated based on the amount of the base material, the number of energy ranges, and the projection data to be corrected.

7. An electronic device, characterized in that, The system includes a processor, a memory, and a computer program stored in the memory and executable on the processor, wherein the computer program, when executed by the processor, implements the image hardening artifact correction method as described in any one of claims 1-5.

8. A computer-readable storage medium, characterized in that, A computer program is stored on the computer-readable storage medium, which, when executed by a processor, implements the image hardening artifact correction method as described in any one of claims 1-5.

Images